Thanks !
I have changed the script into this now, and removed the non-classified parts
and it works. This is just for testing.
However, I was looking for a part in the script which refers to the original
PDE, the PDE of 1 D diffusion. However, I cant find any part in the manual,
where one can modify the PDE into other ones.
This is because I wanted to use fipy to solve a complex nonlinear PDE.
#Python script for solving a Partial Differential Equation in a diffusion case
#import multiprocessing
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...")
Sergio Manzetti
[ http://www.fjordforsk.no/logo_hr2.jpg ]
[ http://www.fjordforsk.no/ | Fjordforsk AS ] [ http://www.fjordforsk.no/ | ]
Midtun
6894 Vangsnes
Norge
Org.nr. 911 659 654
Tlf: +47 57695621
[ http://www.oekolab.com/ | Økolab ] | [ http://www.nanofact.no/ | Nanofactory
] | [ http://www.aq-lab.no/ | AQ-Lab ] | [ http://www.phap.no/ | FAP ]
From: "Thibault Bridel-Bertomeu" <[email protected]>
To: "fipy" <[email protected]>
Sent: Thursday, May 18, 2017 2:00:14 PM
Subject: Re: Problems running a simple PDE
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 < [
mailto:[email protected] | [email protected] ] > :
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
[ mailto:[email protected] | [email protected] ]
[ http://www.ctcms.nist.gov/fipy | http://www.ctcms.nist.gov/fipy ]
[ NIST internal ONLY: [ https://email.nist.gov/mailman/listinfo/fipy |
https://email.nist.gov/mailman/listinfo/fipy ] ]
_______________________________________________
fipy mailing list
[email protected]
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
_______________________________________________
fipy mailing list
[email protected]
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
