examples.levelSet.distanceFunction package

Submodules

examples.levelSet.distanceFunction.circle module

Solve the level set equation in two dimensions for a circle.

The 2D level set equation can be written,

\abs{\nabla \phi} = 1

and the boundary condition for a circle is given by, \phi = 0 at (x - L / 2)^2 + (y - L / 2)^2 = (L / 4)^2.

The solution to this problem will be demonstrated in the following script. Firstly, setup the parameters.

>>> from fipy import CellVariable, Grid2D, DistanceVariable, TransientTerm, FirstOrderAdvectionTerm, AdvectionTerm, Viewer
>>> from fipy.tools import numerix, serialComm

>>> dx = 1.
>>> dy = 1.
>>> nx = 11
>>> ny = 11
>>> Lx = nx * dx
>>> Ly = ny * dy

Construct the mesh.

>>> mesh = Grid2D(dx=dx, dy=dy, nx=nx, ny=ny, communicator=serialComm)

Construct a distanceVariable object.

>>> var = DistanceVariable(name='level set variable',
...                        mesh=mesh,
...                        value=-1.,
...                        hasOld=1)

>>> x, y = mesh.cellCenters
>>> var.setValue(1, where=(x - Lx / 2.)**2 + (y - Ly / 2.)**2 < (Lx / 4.)**2)

>>> var.calcDistanceFunction(order=1) 

>>> if __name__ == '__main__':
...     viewer = Viewer(vars=var, datamin=-5., datamax=5.)
...     viewer.plot()

The result can be tested with the following commands.

>>> dY = dy / 2.
>>> dX = dx / 2.
>>> mm = min (dX, dY)
>>> m1 = dY * dX / numerix.sqrt(dY**2 + dX**2)
>>> def evalCell(phix, phiy, dx, dy):
...     aa = dy**2 + dx**2
...     bb = -2 * ( phix * dy**2 + phiy * dx**2)
...     cc = dy**2 * phix**2 + dx**2 * phiy**2 - dx**2 * dy**2
...     sqr = numerix.sqrt(bb**2 - 4. * aa * cc)
...     return ((-bb - sqr) / 2. / aa,  (-bb + sqr) / 2. / aa)
>>> v1 = evalCell(-dY, -m1, dx, dy)[0]
>>> v2 = evalCell(-m1, -dX, dx, dy)[0]
>>> v3 = evalCell(m1,  m1,  dx, dy)[1]
>>> v4 = evalCell(v3, dY, dx, dy)[1]
>>> v5 = evalCell(dX, v3, dx, dy)[1]
>>> MASK = -1000.
>>> trialValues = CellVariable(mesh=mesh, value= \
...     numerix.array((
...     MASK,  MASK, MASK, MASK, MASK, MASK, MASK, MASK, MASK, MASK, MASK,
...     MASK,  MASK, MASK, MASK, -3*dY, -3*dY, -3*dY, MASK, MASK, MASK, MASK,
...     MASK,  MASK, MASK,   v1,  -dY,  -dY,  -dY,   v1, MASK, MASK, MASK,
...     MASK,  MASK,   v2,  -m1,   m1,   dY,   m1,  -m1,   v2, MASK, MASK,
...     MASK, -dX*3,  -dX,   m1,   v3,   v4,   v3,   m1,  -dX, -dX*3, MASK,
...     MASK, -dX*3,  -dX,   dX,   v5, MASK,   v5,   dX,  -dX, -dX*3, MASK,
...     MASK, -dX*3,  -dX,   m1,   v3,   v4,   v3,   m1,  -dX, -dX*3, MASK,
...     MASK,  MASK,   v2,  -m1,   m1,   dY,   m1,  -m1,   v2, MASK, MASK,
...     MASK,  MASK, MASK,   v1,  -dY,  -dY,  -dY,   v1, MASK, MASK, MASK,
...     MASK,  MASK, MASK, MASK, -3*dY, -3*dY, -3*dY, MASK, MASK, MASK, MASK,
...     MASK,  MASK, MASK, MASK, MASK, MASK, MASK, MASK, MASK, MASK, MASK), 'd'))

>>> var[numerix.array(trialValues == MASK)] = MASK
>>> print(numerix.allclose(var, trialValues)) 
True

examples.levelSet.distanceFunction.interior module

Here we solve the level set equation in two dimension for an interior region. The equation is given by:

Do the tests:

>>> var.calcDistanceFunction(order=1) 

>>> dX = dx / 2.
>>> dY = dy / 2.
>>> mm = dX * dY / numerix.sqrt(dX**2 + dY**2)
>>> def evalCell(phix, phiy, dx, dy):
...     aa = dy**2 + dx**2
...     bb = -2 * ( phix * dy**2 + phiy * dx**2)
...     cc = dy**2 * phix**2 + dx**2 * phiy**2 - dx**2 * dy**2
...     sqr = numerix.sqrt(bb**2 - 4. * aa * cc)
...     return ((-bb - sqr) / 2. / aa,  (-bb + sqr) / 2. / aa)
>>> v1 = evalCell(dY, dX, dx, dy)[1]
>>> v2 = max(-dY*3, -dX*3)
>>> values = numerix.array((  v1,   dY,   dY,  dY,  v1,
...                           dX,  -mm,   -dY,  -mm,   dX,
...                           dX,  -dX,   -v1,  -dX,   dX,
...                           dX,  -mm,   -dY,  -mm,   dX,
...                           v1,   dY,   dY,  dY,  v1  ))
>>> print(var.allclose(values, atol = 1e-10)) 
1

examples.levelSet.distanceFunction.mesh1D module

Create a level set variable in one dimension.

The level set variable calculates its value over the domain to be the distance from the zero level set. This can be represented succinctly in the following equation with a boundary condition at the zero level set such that,

\frac{\partial \phi}{\partial x} = 1

with the boundary condition, \phi = 0 at x = L / 2.

The solution to this problem will be demonstrated in the following script. Firstly, setup the parameters.

>>> from fipy import CellVariable, Grid1D, DistanceVariable, TransientTerm, FirstOrderAdvectionTerm, AdvectionTerm, Viewer
>>> from fipy.tools import numerix, serialComm

>>> dx = 0.5
>>> nx = 10

Construct the mesh.

>>> mesh = Grid1D(dx=dx, nx=nx, communicator=serialComm)

Construct a distanceVariable object.

>>> var = DistanceVariable(name='level set variable',
...                        mesh=mesh,
...                        value=-1.,
...                        hasOld=1)
>>> x = mesh.cellCenters[0]
>>> var.setValue(1, where=x > dx * nx / 2)

Once the initial positive and negative regions have been initialized the calcDistanceFunction() method can be used to recalculate var as a distance function from the zero level set.

>>> var.calcDistanceFunction()

The problem can then be solved by executing the solve() method of the equation.

>>> if __name__ == '__main__':
...     viewer = Viewer(vars=var, datamin=-5., datamax=5.)
...     viewer.plot()

The result can be tested with the following commands.

>>> print(numerix.allclose(var, x - dx * nx / 2)) 
1

examples.levelSet.distanceFunction.square module

Here we solve the level set equation in two dimensions for a square. The equation is given by:

\abs{\nabla \phi} &= 1 \\ \phi &= 0 \qquad \text{at} \qquad \begin{cases} x = \left( L / 3, 2 L / 3 \right) & \text{for $L / 3 \le y \le 2 L / 3$} \\ y = \left( L / 3, 2 L / 3 \right) & \text{for $L / 3 \le x \le 2 L / 3$} \end{cases}

Do the tests:

>>> var.calcDistanceFunction(order=1) 

>>> def evalCell(phix, phiy, dx, dy):
...     aa = dy**2 + dx**2
...     bb = -2 * ( phix * dy**2 + phiy * dx**2)
...     cc = dy**2 * phix**2 + dx**2 * phiy**2 - dx**2 * dy**2
...     sqr = numerix.sqrt(bb**2 - 4. * aa * cc)
...     return ((-bb - sqr) / 2. / aa,  (-bb + sqr) / 2. / aa)
>>> val = evalCell(-dy / 2., -dx / 2., dx, dy)[0]
>>> v1 = evalCell(val, -3. * dx / 2., dx, dy)[0]
>>> v2 = evalCell(-3. * dy / 2., val, dx, dy)[0]
>>> v3 = evalCell(v2, v1, dx, dy)[0]
>>> v4 = dx * dy / numerix.sqrt(dx**2 + dy**2) / 2
>>> arr = numerix.array((
...     v3, v2, -3. * dy / 2., v2, v3,
...     v1, val, -dy / 2., val, v1,
...     -3. * dx / 2., -dx / 2., v4, -dx / 2., -3. * dx / 2.,
...     v1, val, -dy / 2., val, v1,
...     v3, v2, -3. * dy / 2., v2, v3           ))
>>> print(var.allclose(arr)) 
1

examples.levelSet.distanceFunction.test module

Last updated on Jun 27, 2023. Created using Sphinx 6.2.1.