How to add zero-flux boundary conditions to a mesh generated in gmsh Fipy?

Question

I am trying to solve the Meinhart model in a mesh created by gmsh in Fipy. However, I am not sure how to add zero flux-boundary conditions. Below, you can find my code. I would like to know if the lines:

#u.constrain(0, where=mesh.exteriorFaces)
#u.faceGrad.constrain(mesh.faceCenters, where=mesh.exteriorFaces)

should be activated to assure zero-flux boundary conditions.

I read in: https://github.com/usnistgov/fipy/issues/674 that "FiPy's boundaries are no-flux by default". But, I don't know if it also applies for mesh created by gmsh options.

If there is not needed to activated or add additional lines to my code to set up zero flux boundary condition, it also works when I have more complex mesh, e.g. polygons, circles, or other irregular shapes?

Thanks,

"""
Emacs Editor

This is a temporary script file.
"""

# -*- coding: utf-8 -*-
"""
@author: Irbin B.
"""
'''Solving Meinhart model 2D'''


# 1. Libraries
import fipy as fi                   # Finite volume method's package

# 2. Building the domain (Gmsh square)
## 2.1. Domain lenght
nx = 100
ny = nx

## 2.2. Gmsh config
mesh = fi.Gmsh2D('''
Side = 100;
CellSize = 1;
Point(1) = {0, 0, 0, CellSize};   
Point(2) = {0, Side, 0, CellSize};
Point(3) = {Side, Side, 0, CellSize};
Point(4) = {Side, 0, 0, CellSize};
Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 1};
Line Loop(6) = {1, 2, 3, 4};
Plane Surface(7) = {6};
''' % locals())


## 2.3. Adding initial conditions values (Random)
noise_u = fi.GaussianNoiseVariable(mesh=mesh,
                                   mean=0.5,
                                   variance=0.05).value
noise_v = fi.GaussianNoiseVariable(mesh=mesh,
                                   mean=0.5,
                                   variance=0.05).value

# 3. Zero-Flux boundary conditions
BCs = (fi.FixedFlux(faces=mesh.facesRight, value=0.),
       fi.FixedFlux(faces=mesh.facesLeft, value=0.),
       fi.FixedFlux(faces=mesh.facesTop, value=0.),
       fi.FixedFlux(faces=mesh.facesBottom, value=0.))

# 4. Defining the variables
u = fi.CellVariable(name = "u",
                    mesh=mesh,
                    value=0.,
                    hasOld=True)

u[:] = noise_u
#u.constrain(0, where=mesh.exteriorFaces)
#u.faceGrad.constrain(mesh.faceCenters, where=mesh.exteriorFaces)

v = fi.CellVariable(name = "v",
                    mesh = mesh,
                    value = 0.,
                    hasOld = True)

v[:] = noise_v
#v.constrain(0, where=mesh.exteriorFaces)
#v.faceGrad.constrain(mesh.faceCenters, where=mesh.exteriorFaces)

# 5. Defining the parameters
Da = 1.
Db = 100
alpha = -0.005
beta = 10

# 6. Creating the system of PDEs
equ = fi.TransientTerm(var=u) == fi.DiffusionTerm(coeff=Da, var=u) + u - u**3 - v + alpha

eqv = fi.TransientTerm(var=v) == fi.DiffusionTerm(coeff=Db, var=v) + (u - v) * beta

eqn = (equ & eqv)

# 7. Solving the PDEs and showing the results in figures
timeStepDuration = .1
steps = 100


for step in range(steps):
    u.updateOld()
    v.updateOld()
    eqn.sweep(dt = timeStepDuration)
    print(step+1)

if __name__== '__main__':
    uviewer = fi.Viewer(vars=u, datamin=0., datamax=.7)
    vviewer = fi.Viewer(vars=v, datamin=0., datamax=.7)

Answer 1

FiPy's boundary conditions are no-flux by default for all meshes. It's a fundamental characteristic of the cell-centered finite volume method.