If "psi = epsilon**2*laplace phi"(Split plan 2 in codes) is a better
choice than "psi = dfdphi - epsilon**2*laplace phi"(Split plan 1 in
codes)? Because "Split plan 2" runs much faster than "Split plan 1". Is
there any hidden trouble using "Split plan 2"?
import fipy as fp
nx = 64
ny = 64
mesh = fp.Grid2D(nx=nx, ny=ny, dx=1., dy=1.)
psi = fp.CellVariable( hasOld=1,mesh=mesh )
phi = fp.CellVariable( mesh = mesh, value =
fp.GaussianNoiseVariable(mesh=mesh, mean=0.5, variance=0.01).value,
hasOld=1 )
D = a = 1.
epsilon = 0.5
dfdphi_ = a**2 * (1 - phi) * (1 - 2 * phi)
d2fdphi2 = a**2 * (1 - 6 * phi * (1 - phi))
#Split plan 1
eq1 = (fp.TransientTerm(coeff=1.,var=phi) ==
fp.ImplicitDiffusionTerm(coeff=D, var=psi))
eq2 = (fp.ImplicitSourceTerm(coeff=1., var=psi) ==
fp.ImplicitSourceTerm(coeff=dfdphi_, var=phi) -
fp.DiffusionTerm(coeff=epsilon**2, var=phi))
#Split plan 2
#eq1 = (fp.TransientTerm(coeff=1.,var=phi) ==
fp.ImplicitDiffusionTerm(coeff=d2fdphi2, var=phi) +
fp.ImplicitDiffusionTerm(coeff=-D, var=psi))
#eq2 = (fp.ImplicitSourceTerm(coeff=1., var=psi) ==
fp.DiffusionTerm(coeff=epsilon**2, var=phi))
eq = eq1 & eq2
dt = 0.5
viewer = fp.Viewer(vars=(phi, psi))
for t in range(1000):
phi.updateOld()
psi.updateOld()
res = 1.e5
while res > 1.e-1:
res = eq.sweep(dt=dt)
print res
viewer.plot()
On 20/02/15 19:22, Guyer, Jonathan E. Dr. wrote:
On Feb 20, 2015, at 6:39 AM, Ronghai Wu <ronghai...@fau.de> wrote:
(1) The 4th-order is split by "psi = d2fdphi2(phi-phiold) + dfdphi -
epsilon**2*laplace phi". I do not understand why we need
d2fdphi2(phi-phiold)?
See the discussion on linearization of the source in
http://www.ctcms.nist.gov/fipy/examples/phase/generated/examples.phase.simple.html
You could write "psi = dfdphi - epsilon**2*laplace phi", but dfdphi would be a
fully explicit source. By expanding
source = source_old + d(source)/d(phi) (phi - phi_old)
you get a more implicit version of the same term, with better convergence
properties.
_______________________________________________
fipy mailing list
fipy@nist.gov
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
--
------------------------------------------
Ronghai Wu
Institute of Materials Simulation (WW8)
Department of Materials Science and Engineering
University of Erlangen-Nürnberg
Dr.-Mack-Str. 77, 90762 Fürth, Germany
Tel. +49 (0)911 65078-65064
_______________________________________________
fipy mailing list
fipy@nist.gov
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]