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 ]

Reply via email to