Solve the steady-state convection-diffusion equation in one dimension.
This example solves the steady-state convection-diffusion equation given by
with coefficients and , or
>>> diffCoeff = 1. >>> convCoeff = (10.,)
We define a 1D mesh
>>> from fipy import CellVariable, Grid1D, DiffusionTerm, ExponentialConvectionTerm, Viewer >>> from fipy.tools import numerix
>>> L = 10. >>> nx = 10 >>> mesh = Grid1D(dx=L / nx, nx=nx)
>>> valueLeft = 0. >>> valueRight = 1.
The solution variable is initialized to
>>> var = CellVariable(mesh=mesh, name="variable")
and impose the boundary conditions
>>> var.constrain(valueLeft, mesh.facesLeft) >>> var.constrain(valueRight, mesh.facesRight)
The equation is created with the
ExponentialConvectionTerm. The scheme used by the convection term
needs to calculate a Péclet number and thus the diffusion term
instance must be passed to the convection term.
>>> eq = (DiffusionTerm(coeff=diffCoeff) ... + ExponentialConvectionTerm(coeff=convCoeff))
More details of the benefits and drawbacks of each type of convection
term can be found in Numerical Schemes.
both handle most types of convection-diffusion cases, with the
PowerLawConvectionTerm being more efficient.
We solve the equation
and test the solution against the analytical result
>>> axis = 0 >>> x = mesh.cellCenters[axis] >>> CC = 1. - numerix.exp(-convCoeff[axis] * x / diffCoeff) >>> DD = 1. - numerix.exp(-convCoeff[axis] * L / diffCoeff) >>> analyticalArray = CC / DD >>> print(var.allclose(analyticalArray)) 1
If the problem is run interactively, we can view the result:
>>> if __name__ == '__main__': ... viewer = Viewer(vars=var) ... viewer.plot()