Skip to content

Code: sandbox

William Jussiau edited this page Feb 24, 2025 · 73 revisions

Page with all code information, that will be split later into 3 pages

Under construction

Code: basic use

Convention

Use of the _prefix for methods that are not intended to be used outside of the body of a class.

Basic use

Define a new use-case: inherit FlowSolver abstract class

The simulation revolves around the abstract class FlowSolver that implements core features such as loading mesh, defining function spaces, trial/test functions, variational formulations, numerical schemes and solvers, handling the time-stepping and exporting information. The class is abstract as it does not implement a simulation case per se, but only provides utility for doing so. It features two abstract methods, that are redefined for each use-case:

  1. _make_boundaries provides a definition and naming of the boundaries of the mesh in a pandas DataFrame.
    @abstractmethod
    def _make_boundaries(self) -> pd.DataFrame:
        pass
  1. _make_bcs provides a description of the boundary conditions on the boundaries defined above, in a dictionary.
    @abstractmethod
    def _make_bcs(self) -> dict[str, Any]:
        pass

For the two aforementioned examples, these methods are reimplemented in the classes CylinderFlowSolver and CavityFlowSolver that inherit from FlowSolver.

You also need to be able to provide a mesh that is consistent with your definition of boundaries

Attach Sensors and Actuators to an instance of a FlowSolver subclass

In order to perform sensing and actuation (in order to close the loop), dedicated classes Sensor and Actuator are proposed. They are not aimed at being instantiated, but rather inherited.

The sensors and actuators are attached to a FlowSolver object as a list, embedded in the ParamControl object. The call to Sensors and Actuators is made automatically by FlowSolver. By attaching several sensors or actuators, it is possible to use Multiple-Input, Multiple-Output controllers in the loop.

In the example below (for the cylinder use-case), we are creating an actuator acting on boundary conditions and three point probes at different locations. They are gathered in a ParamControl object, which is passed as an argument when creating a CylinderFlowSolver.

# Actuator
actuator_bc = ActuatorBCParabolicV(angular_size_deg=10)
# Sensors
sensor_feedback = SensorPoint(sensor_type=SENSOR_TYPE.V, position=np.array([3, 0]))
sensor_perf_1 = SensorPoint(sensor_type=SENSOR_TYPE.V, position=np.array([3.1, 1]))
sensor_perf_2 = SensorPoint(sensor_type=SENSOR_TYPE.V, position=np.array([3.1, -1.3]))
# Gather actuators and sensors in ParamControl object
params_control = ParamControl(
    sensor_list=[sensor_feedback, sensor_perf_1, sensor_perf_2],
    actuator_list=[actuator_bc],
)

Run a closed-loop simulation

Once a use-case has been defined by implementing the corresponding class, the basic feedback syntax has the following philosophy:

# Instantiate and initialize FlowSolver object
fs = CylinderFlowSolver(...)
fs.compute_steady_state(...)
fs.initialize_time_stepping(...)

# Instantiate Controller (e.g. load from .mat file)
Kss = Controller.from_file(...)

# Time loop
y_meas = fs.y_meas
for _ in range(fs.params_time.num_steps):
    u_ctrl = Kss.step(y=-y_meas[0], dt=fs.params_time.dt)
    y_meas = fs.step(u_ctrl=u_ctrl)

See examples for a more detailed description.


More details


Initializations in general

Mesh

Mesh as xdmf. User is responsible for mesh coherence with respect to boundary conditions.

Initialization of a FlowSolver

Make all parameter structures from flowsolverparameters.py Timeline of function calls at __init__ for example.

Initialize time-stepping

initialize_time_stepping -> functions with two different ways, whether Tstart=0 (start from ic) or not (restart)

Perturbation of initial state:

_default_initial_perturbation(self, xloc: float = 0.0, yloc: float = 0.0, radius: float = 1.0) -> dolfin.Function

Base flow computation

Good practice is to do Picard iterations, then Newton

Initial guess for Picard only:

_default_steady_state_initial_guess(self) -> dolfin.UserExpression

Inlet flow profile

By default uniform u=1, can be changed in make_BCs. The perturbation bc is still 0 on the inlet.

FlowField, FlowFieldCollection, self.fields

Utility class to gather fields (velocity, pressure and mixed)


Actuators

Principle

Actuator is an abstract class that encapsulates a dolfin.Expression among other elements, and embeds it in the variational formulations or in boundary conditions.

Actuator can be FORCE or BC

  • FORCE: included in momentum equation (automatic in code)
  • BC: included in boundary conditions (by user) in _make_bcs() and should be handled by the user. An example can be found in examples/cylinder/cylinderflowsolver.py:
def _make_bcs(self):
    ...
    bcu_actuation_up = dolfin.DirichletBC(
        self.W.sub(0),
        self.params_control.actuator_list[0].expression,
        self.get_subdomain["actuator_up"],
        )
    bcu_actuation_lo = dolfin.DirichletBC(
        self.W.sub(0),
        self.params_control.actuator_list[1].expression,
        self.get_subdomain["actuator_lo"],
        )
    ...
    return BoundaryConditions(bcu=bcu, bcp=[])

Actuators need to be loaded after the FlowSolver is instantiated (analytic Expression is projected onto FEM), which is handled by the code

Examples of actuators

  • ActuatorBCParabolicV: boundary condition actuator, 2nd component on velocity has parabolic profile

Mathematical expression:

$${v_{act}}({x}, t) = - \dfrac{(x_1-l)(x_1+l)}{l^2} u(t)$$, with $l = \frac{1}{2} D \sin \left( \frac{\delta}{2} \right)$ and $\delta$ is the tunable actuator opening in degrees.

FEniCS syntax:

    def load_expression(self, flowsolver):
        L = (
            1
            / 2
            * flowsolver.params_flow.user_data["D"]
            * np.sin(1 / 2 * self.angular_size_deg * dolfin.pi / 180)
        )
        expression = dolfin.Expression(
            [
                "0",
                "(x[0]>=L || x[0] <=-L) ? 0 : u_ctrl * -1*(x[0]+L)*(x[0]-L) / (L*L)",
            ],
            element=flowsolver.V.ufl_element(),
            L=L,
            u_ctrl=0.0,
        )

        self.expression = expression
  • ActuatorForceGaussianV: force actuator, gaussian-shaped on the 2nd component of velocity

Mathematical expression:

$$B({x})u(t)=\left[ 0, \eta \exp\left( \frac{\left(x_1 - x_1^0\right)^2 + \left(x_2 - x_2^0\right)^2}{2\sigma_0^2} \right)\right]^T u(t)$$ with $\eta$ such that $\int_\Omega B({x})^T B({x}) d\Omega = 1$.

FEniCS syntax:

    def load_expression(self, flowsolver):
        expression = dolfin.Expression(
            [
                "0",
                "u_ctrl * eta*exp(-0.5*((x[0]-x10)*(x[0]-x10)+(x[1]-x20)*(x[1]-x20))/(sig*sig))",
            ],
            element=flowsolver.V.ufl_element(),
            eta=1,
            sig=self.sigma,
            x10=self.position[0],
            x20=self.position[1],
            u_ctrl=1.0,
        )

        BtB = dolfin.norm(expression, mesh=flowsolver.mesh)
        expression.eta = 1 / BtB
        expression.u_ctrl = 0.0
        self.expression = expression

Defining new actuators

Define a subclass with a dedidacted expression.

⚠️ Do not forget ⚠️

  • Include a u_ctrl field in the Expression. Its value should be 0.0 when the method load_expression(self, flowsolver) exits. Its value may be set to something else in the body of the method (see ActuatorForceGaussianV),
  • Assign self.expression = expression at the end of def load_expression(self, flowsolver).

Sensors

Principle

Sensor is an abstract class that provides a method eval(self, up: dolfin.Function) -> float. Classes SensorPoint (point probe) and SensorIntegral (integration on a subdomain) are examples of subclasses that implement the eval method.

The evaluation of sensors is handled by FlowSolver in the step method.

Defining new sensors

User is responsible for parallel execution when implementing new sensors.

Saving files

save every

Controller

What is this class

control.StateSpace encapsulating a current state x and potentially a file

What is step

Wrapper around control.StateSpace.forced_response



Code: advanced use

  • Restart (show graph)
  • See other bullet points in Home
## Additional uses of the toolbox
The toolbox provides additional utility related to flow control:
* Compute dynamic operators A, B, C, D and mass matrix E,
* Restart a simulation from a saved file,
* Define an arbitrary number of actuators and sensors (e.g. feedback and performance sensors),
* Export time series (measurements from sensors, perturbation kinetic energy...) and fields for visualization,
* Modify the equations, the numerical schemes and the solvers used for the time simulation,
* Leverage parallel execution native to FEniCS,
* Use it as backend in an optimization tool (as in [Jussiau, W., Demourant, F., Leclercq, C., & Apkarian, P. (2025). Control of a Class of High-Dimensional Nonlinear Oscillators: Application to Flow Stabilization. IEEE Transactions on Control Systems Technology.](https://ieeexplore.ieee.org/abstract/document/10884641/)).


Code: utility

  • Operator computation
  • Frequency response computation
  • Export utils (spy, save)
  • Debug utils (export subdomains)
  • Mpi utils (= shortcuts)

Clone this wiki locally