Can’s sandbox

You can find here documentations, examples and test cases related to my implementation of a fictitious-domain method within basilisk.

Distributed Lagrange multiplier/fictitious domain method for rigid particle laden flows

We are looking for a numerical approximation for the solution of the following equations ρ(tu+u\nablau)=\nablap+\nabla(2μD)+ρa+\lambdaoverΩ (1ρρs)(M(tUg))=P(t)\lambdadvoverP(t) (1ρρs)(It\omega+\omega×I\omega)=P(t)r×\lambdadvoverP(t) u(U+\omega×r)=0overP(t) \nablau=0overΩ

with unknowns u,U,\omega and \lambda being respectively the fluid velocity field, the particle’s translational velocity, the particle’s rotational velocity and the Lagrange multipliers. The particle occupies the domain P(t) with density ρs, inertia tensor I and mass M. The vector g is the gravity acceleration here.

From a numerical point of view, this is a complicated problem to solve direclty. One has to deal with three main difficulties:

  • an advection-diffusion equation
  • the imcompressibility contraint with the pressure p as unknown
  • the particle’s solid body’s constraint with the Lagrange multipliers as unknow.

The classic strategy is to split (and decouple) the problem into subproblems and solve them successively. We chose here a two-steps time-spliting.

1st sub-problem: Navier-Stokes

We are using the centred(.h) solver for this problem.



2nd subproblem: fictitious-domain problem

We are solving a fictitious-domain problem given by:

ρtu=\lambdaoverΩ (1ρρs)(M(tUg))=P(t)\lambdadvoverP(t) (1ρρs)(It\omega+\omega×I\omega)=P(t)r×\lambdadvoverP(t) u(U+\omega×r)=0overP(t),

This leads to a saddle-point problem which is solved with an iterative solver (Uzawa, conjugate gradient algorithm).

Numerical scheme

Basilisk handles the resolution with its own temporal scheme. It reads:



Fictitious domain

The nunerical scheme for the fictitious domain problem reads: ρ(ufdn+1unsn+1Δt)\lambdan+1=0overΩ (1ρρs)(M([Un+1UnΔt]g))=P(t)\lambdan+1dvoverP(t)

(1ρρs)(I(\omegan+1\omeganΔt)+\omegan+1×I\omegan+1)=P(t)r×\lambdan+1dvoverP(t) ufdn+1(Un+1+\omegan+1×r)=0overP(t).

Test cases




Stencils and fictititous domains