CLSVOF

Introduction

CLSVOF is a coupled level set - volume of fluid method. The method corrects the level set to make it locally (for each element containing the interface) volume conserving.

Serial Algorithm

Loop over elements (find the elements containing interface pieces)

1.(I) Tag each element as incomplete

 Need a tag for elemental completeness
 Also need a tag for nodal completeness
 These tags get dimensioned (size numel and size numnp) in itrdrv?

2.(C) Truncate volume fractions so they’re 0 ≤ F ≤ 1

 Volume fraction is the second scalar
 

3.(C) If 0 < F < 1 (element contains interface)

3.a.(I) Use advected φ at element nodes to find interface slope for the element

 For tetrahedra, find the points at which the interface intersects the edges
 Use those three points (defining the plane) to find the interface normal
 <math> \bold  n = ( \bold p_2 - \bold p_1 ) \times ( \bold p_3 - \bold p_1 ) </math>
 Need a way to temporarily store those points (for volume computation)
 Need to temporarily store the normal

There are 2 possible cases for φ in a cell with 0 < VOF < 1:

1.(L) All nodal φ are the same sign - use smallest φ and assume interface plane is parallel to opposite face

2.(I) All φ are different signs

3.b.(I) Move interface along normal until computed F = advected F

Cases 3, 5, and 7 above have the interface plane passing through the element. (Cases 1, 2, 4, and 6 are degenerate cases.) When the element is divided by the interface plane, the resulting volumes are either 1 tetrahedron and 1 pentahedron (interface plane intersects 3 edges) or 2 pentahedra (interface plane intersects 4 edges).

(A) Figure out which side of the interface is a tetrahedron. F is either the volume of the interface tetrahedron or the volume of the element minus the interface tetrahedron.

 For a tetrahedron with vertices a, b, c, and d, the volume is:  <math>V = \frac { |(\mathbf{a}-\mathbf{d}) \cdot ((\mathbf{b}-\mathbf{d}) \times (\mathbf{c}-\mathbf{d}))| } {6}.</math>

When there are two pentahedra, the volume of either one can be found by dividing it into two tetrahedra and a 5-sided pyramid. The vertices of each tetrahedron are an element vertex, two points where the interface plane intersects element edges, and the midpoint of an element edge not intersected by the interface plane. The vertices of the pyramid are the four points where the interface plane intersects the element edges and the midpoint of an element edge not intersected by the interface plane.

 For a tetrahedron with vertices a, b, c, d, and e, the volume is:  <math>V = \frac { |(\mathbf{e}-\mathbf{a}) \cdot ((\mathbf{b}-\mathbf{a}) \times (\mathbf{d}-\mathbf{a}))| } {3}.</math>

Shashkov uses a line search algorithm to adjust the interface. -> Use a recursive bisection.

After each movement of the interface, the LS=0 points on the edges have to be found again, and the volume recomputed. When the computed volume fraction is within a tolerance of the actual volume fraction, the new interface has been found, and the nodes can be updates.

3.c.(I) Use the newly established interface to compute φ at each node and tag those nodes as reconstructed

 At each node, φ is the distance to the nearest point, unless two points are equidistant.
 Need a tag for nodal completeness/reconstruction

3.d.(A) If any node has already been reconstructed, φ is min (previous reconstruction, current reconstruction) (ensure that φ is min to interface everywhere)

3.e.(A) Mark the element as complete

3.f. Loop over elements adjacent by shared faces - If it’s not marked as complete and not an interface element, add it to the reconstruction list

 Need reconstruction list
 Need data structure of each element's faces
 Need data structure of which two elements are adjacent to each face

While reconstruction list is not empty, loop over the reconstruction list (establish and walk out φ)

1. N_fixed = 0

2. For each element

2.a. num_const = 0

2.b. loop over nodes

2.b.1. If a node has been reconstructed, num_const = num_const+1

2.c. If num_const <= 3

2.c.1 Compute gradient of φ for the element using the three smallest values of reconstructed φ

 if num_const = 3, then use the reconstructed nodes
 if num_const > 3, then loop over nodes to find the three smallest
 re-use PHASTA's gradient calculation method?

2.c.2. If computed gradient is very different from advected gradient (nearest interface is not the one walked out from)

 if computed gradient dotted with advected gradient is less than pi/4 (or some other suitable value)

2.c.2.a. The current element is done for now, and left in the list; move on to the next element

2.c.3. Else

2.c.3.a. Compute φ at all other nodes

 gradient gives normal to the planar field
 φ at a node is φ at another node plus the normal distance from the plane that other node is in to the current node

2.c.3.b.1. If a node is tagged as reconstructed, φ = min (current value, new computed value)

2.c.3.b.2. Else φ = computed value, node is tagged as reconstructed

2.c.3.c. Mark the element as complete

2.c.3.d. Loop over elements adjacent by shared faces

2.c.3.d.1. If the adjacent element is not marked as complete, and not in the reconstruction list, and reconstructed φ of shared nodes is < nε add it to the list (stop adding nodes when > nε from interface)

 check completeness
 check reconstruction list
 check φ of shared nodes, which requires some nodal adjacency information

2.c.3.e. Remove the element from the reconstruction list

2.c.3.f. N_fixed = N_fixed + 1

2.d. Else (element has only 1 or 2 reconstructed nodes), skip the element and move on to next in list

3. At the end of a cycle through the list, if N_fixed ≠ 0, go to start of loop

3.a. Else N_fixed = 0 and list is not empty, loop again through the list (less accurate method)

3.a.1. If reconstructed nodes ≤ 2

3.a.1.a. For each node not yet reconstructed, φ = min (distance to adjacent node + reconstructed φ of that node)

 for tetrahedra, all element nodes are adjacent

3.a.1.b. Mark the element as complete

3.a.1.c. Loop over elements adjacent by shared faces

3.a.1.c.1. If it’s not marked as complete, and not in the reconstruction list, and reconstructed φ of shared nodes is < nε add it to the list (stop adding nodes when > nε from interface)

3.a.1.d. Remove the element from the reconstruction list

3.a.1.e. N_fixed = N_fixed + 1

3.a.2. Else (at least 3 nodes are tagged as reconstructed)

3.a.2.a. Handle element as above (for at least three reconstructed nodes)

3.a.2.b. N_fixed = N_fixed + 1

When list is empty, for all elements not tagged as reconstructed, φ = nε for all nodes not tagged as reconstructed