Presentation of the computational code
This code, that I developed during my PhD in collaboration with Loïc Gouarin and Bertrand Maury, is written in C++ and uses the library PETSc. It works with the 4Q1/Q1 finite element to simulate the motion of immersed rigid bodies with the smooth extension method developed in my PhD. The library PETSc is used for the parallel structure of matrices and vectors, aswell as the parallel solvers for the linear systems. The code is also interfaced with the code SCoPI in order to handle contacts between the rigid bodies.
Algorithm
The idea of the method is to seek a smooth extension of the exact solution in the whole fictitious domain. For this purpose, the extension of the right-hand side within the rigid bodies is used. The problem amounts to look for an extension of the right-hand side such that the restriction of the solution of a classical Stokes problem, defined in the whole fictitious domain, to the fluid domain is the exact solution of the original problem (the fluid velocity satisfies the rigid-body motion constraint on the boundary of the rigid bodies). To find such an extension of the right-hand side, a cost function is minimized with a conjugate gradient like algorithm. This minimization process involves the resolution of two classical Stokes problems at each iteration of the conjugate gradient. The algorithm of the code interfaced with the contact code SCoPI is the following:
-
SCoPI sends the position of the rigid bodies to the code.
-
Resolution of a classical Stokes problem to initialize the conjugate gradient.
-
Conjugate gradient loop:
- Resolution of two classical Stokes problem to get the gradient of the cost function.
- Update of the extension of the right-hand side.
-
End of conjugate gradient.
-
Resolution of a classical Stokes problem with the extension of the right-hand side found by the minimization process.
-
Computation of the velocity of each rigid body from the fluid velocity.
-
The code sends these velocities to SCoPI that projects them on the space of admissible velocities.
A lot of Stokes problems have to be solved with this method but the only part that is changing from one Stokes problem to another is the right-hand side. The matrices are the one of a classical incompressible Stokes problem and are always the same. Moreover, because we are using cartesian meshes, fast iterative solvers like geometric multigrids are used to solve these Stokes problems.
Some simulations
| Description | Video |
|---|---|
| 7500 rigid bodies sedimenting in 3D. The rigid bodies are evenly distributed at initial time and have a radius of 0.0125. The cartesian meshes has the size 256x256x256 for the velocity. The simulation was ran on 24 cores and took 22 minutes by time iteration. | |
| 2000 rigid bodies sedimenting in 3D. The top of the cube is rotating. The rigid bodies have different radius. The cartesian mesh has the size 65x65x65 for the velocity. | |
| Sedimentation of 1008 rigid bodies with the same radius of 0.01. The cartesian mesh has a size of 1024x1024 for the velocity. |
Additional information
-
This code has been extended to solve other similar types of problem and rewritten in modern C++ by Loïc Gouarin. It is now named CAFES.
-
Related documents:
- B. Fabrèges, L. Gouarin and B. Maury, A smooth extension method, C. R. Acad. Sci. Paris, 2013.
- B. Fabrèges and B. Maury, Approximation of single layer distributions by Dirac masses in finite elements computations, J. Sci. Comput., 2014.
- The manuscript of the thesis and slides of the defense (in french) are available here.