Hello Sergio,
Instead of « for j in range(cells) » you can try « for j in
range(phi.mesh.numberOfCells) », it should match what you want to
accomplish.
Best,
Thibault Bridel-Bertomeu
—
MSc, PhD
Design and Development Engineer
CERFACS, CFD Team
42 Av. Gaspard Coriolis
31057 TOULOUSE CEDEX 1
FRANCE
2017-05-18 13:36 GMT+02:00 Sergio Manzetti <sergio.manze...@fjordforsk.no>:
> Hello, following the manual, at the first example with the 1D diffusion, I
> have tried to make the python script for this (python2.7) and I get a
> missing name, which is not defined in the example at all.
>
> This is the name "cell", which appears on p 58 in the manual.
>
> Please see this code which describes this example:
>
>
> Thanks!
>
>
>
> #Python script for solving a Partial Differential Equation in a diffusion
> case
>
>
>
> from fipy import *
> nx = 50
> dx = 1
> mesh = Grid1D(nx=nx, dx=dx)
> phi = CellVariable(name="Solution variable",
> mesh=mesh,
> value=0.)
>
>
> #We'll use the coefficient set to D=1
>
> D=1
>
> # We give then 2000 time-steps for the computation, which is defined by 90
> percent of the maximum stable timestep, which is given
>
> timeStepDuration = 0.9 * dx**2 / (2 * D)
> steps = 2000
>
>
> #Then we define the boundary conditions
>
> valueLeft = 1
>
> valueRight = 0
>
> #The boundary conditions are represented as faces around the exterior of
> the mesh, so we define the constraint by the given values on the left and
> right side of the boundary using the phi.constrain() command
>
> phi.constrain(valueLeft, mesh.facesRight)
> phi.constrain(valueRight, mesh.facesLeft)
>
> # We can also omit giving boundary conditions, then the default bc is
> equivalent to a zero gradient.
>
> #At this stage we define the partial differential equation as a numerical
> problem to be solved
>
> for step in range(steps):
> for j in range(cells):
> phi_new[j] = phi_old[j] \
> + (D * dt / dx**2) * (phi_old[j+1] - 2 * phi_old[j] +
> phi_old[j-1])
> time += dt
>
> #and additional code for the boundary conditions
>
> eqX = TransientTerm() == ExplicitDiffusionTerm(coeff=D)
>
>
> #We then want to view the results of the calculation and use the command
> Viewer() for this
>
> phiAnalytical = CellVariable(name="analytical value",
> mesh=mesh)
> if __name__ == '__main__':
> viewer = Viewer(vars=(phi, phiAnalytical),
> datamin=0., datamax=1.)
> viewer.plot()
>
> #If we have a semi-infinite domain, then the solution for this PDE is
> phi=1-erf(x/2(Dt)^(0.5)). This requires to import SciPy library, so we
> import that.
>
> x = mesh.cellCenters[0]
> t = timeStepDuration * steps
>
> try:
> from scipy.special import erf
> phiAnalytical.setValue(1-erf(x / (2 * numerix.sqrt(D * t))))
> except ImportError:
> print ("The Scipy Library is not avaliable to test the solution to the
> given trasient diffusion equation")
>
> # The equation is then solved by repeatedly looping
>
> for step in range(steps):
> eqX.solve(var=phi,
> dt=timeStepDuration)
> if __name__ == '__main__':
> viewer.plot()
> print (phi.allclose(phiAnalytical, atol = 7e-4))
> 1
>
> if __name__ == '__main__':
> raw_input("Explicit transient diffusion. Press <return> to proceed...")
>
>
> _______________________________________________
> fipy mailing list
> fipy@nist.gov
> http://www.ctcms.nist.gov/fipy
> [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
>
>
_______________________________________________
fipy mailing list
fipy@nist.gov
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
