Modelling Deforming Interfaces using Level Sets by Hans Mühlhaus, Laurent Bourgouin and

1. The Australian Computational Earth Systems Simulator (ACcESS) Modelling Deforming Interfaces using Level Sets by Hans Mühlhaus, Laurent Bourgouin and Lutz Gross

2. Overview Introducing Example What’s needed to model interfaces • Constitutive models • Surface tracking, level set, stress advection EScript & Finley Separation of Physics from Linear Algebra and Parallel Computing Applications • Advection, upwinding, implicit vs. explicit • Lava dome simulation,Subduction, Rayleigh-Taylor Instabilities

3. Moving Interface: a 1D Example v2 h=small j<0 x2 h=large j>0 h x1 We define: so that Update: Note that: Time integration:

4. Governing equations Stress Equilibrium Heat Equation Concentration advection: Temperature and concentration dependence of density:

5. Example for Rayleigh – Taylor Instabilities using level sets: Mantle Plumes

6. The General Case • Implicit representation of the interface by the zero level set of a smooth function φ • φis usually chosen as a • “signed” distance function ( ) • At each time step, φ is updated solving the (hyperbolic) advection equation:

7. Problems…… • Symmetric difference expressions like • (symm.) (non-symm.) • don’t work well in hyperbolic problems (upwinding etc!) • 2. Inhomogeneous velocity field causes loss • of distance function property ( ) of j

8. Problems……(cont.) • Symmetric difference expressions don’t • work well in hyperbolic problems (upwinding!) 1. Upwinding If v is constant:

9. Problems…… This can be transformed into a non-symmetric expression by adding…. We expect that the FE approx. of the PDE: is better conditioned than the original Hyperbolic problem

10. Generalisations….. Taylor-Galerkin: 2-step alternative to Taylor-Galerkin upwinding (very effective in the presence of diffusionterms….):

11. 2 Gaussians

12. 1 Gaussian

13. The Level Set Method: Solving the advection equation Test: A Gaussian is advected in a constant 1D velocity field. • Explicit • Implicit • Taylor Galerkin

14. Software can be downloaded fromwww.esscc.uq.edu.au, contact Ken Steube (esys@access.edu.au) If you need instructions re libraries etc Formulation Finley PDE: Example : Momentum and Heat equation

15. LinearPDE class General form (as relevant here): PDE: natural boundary condition Kronecker symbol: δij=0 for i=j and 0 otherwise

16. Helmholtz Class in mytools.py from esys.linearPDEs import LinearPDE import numarray class Helmholtz(LinearPDE): def setValue(self,kappa,omega,f,eta,g): ndim=self.getDim() # spatial dimension kronecker=numarray.identity(ndim) self._setValue(A=kappa*kronecker,\ D=omega,Y=f,d=eta,\ y=g)

17. Use the Helmholtz Class # Helmholtz class defined in mytools.py from mytools import Helmholtz mydomain=... mypde=Helmholtz(mydomain) mypde.setValue(kappa=10,omega=0.1,\ f=12,eta=0,g=0) u=mypde.getSolution()

18. 2. Problem: Inhomogeneous velocity field causes loss of distance function property of j Previous test: No topological change in the solution Need for a new test with: and New test: shear flow • Mesh: 100x100 • Courant Number: 0.25 • 1000 steps forward • 1000 steps with -v

19. The Level Set Method: Solving the advection equation The shape gets “noisy”… Problem: φ looses its distance function property Reinitialisation needed!

20. The Level Set Method: Reinitialisation • Idea: • Rebuild a “signed” distance function ψ from the distorted function φ • Requirements: • The interface must not be changed • ψ must represent a distance function • Solution: • Solve to steady state the equation: • Rewritten as: with Interpretation: The “distance information” is carried by w, a unit vector pointing away from the interface.

21. Remarks on re-initialisation….. • During iteration (pseudo time integration) the vector w is established once and then kept constant • In the explicit solution of the advection problem for y we found that only alumped mass matrix discretisation works

22. The Level Set Method: Reinitialisation (2/3) 1D 2D 3D

23. The Level Set Method: Reinitialisation Same test as before, with reinitialisation

24. Level set cont. : Merger of small and large bubbles Parameters: Surface tension: Calculation, includes inertia, Courant Number=0.5, msh:30 by 45 8 node quad’s

25. Level set cont. : Calculation of curvature for C_0 continuity Projection: and Representation of surface tension b.c. as volume force: =distance in the direction of the normal of l smoothing length, related to the element size

26. Level set cont. : Merger of small and large bubbles

27. Surface Tension: Benchmark

28. Level set: Surface membrane shell, surface tension where Inserting yields where at equilibrium.

29. Collapsing Cylinder

30. Lava Dome

31. Remarks • Escript & Finley: Rapid development of simulation software; parallelised assembly and solution phase; separation of physics from linear algebra • Level set modelling of interfaces: distance function property crucial • Modelling of surface tension; example of higher order b.c.’s • Upwinding strategy dependent on element type • Re-initialisation strategy has an (undesirable) element of mystique…..