On Feb 26, 2010, at 2:24 PM, Anjuli Appapillai wrote:
I'm working on a problem similar to the example 8.1 in the FiPy
manual. However
I don't quite understand the linearization of the source term using
a Taylor
expansion (p.108). The text seems to say that the Taylor expansion
is about the
previous value of the dependent variable, phase (phi), such that:
S = Sold + dS/dphi|old * (phi - phi_old)
However, the actual code for this source term does not involve
phase.getOld() as
I would expect:
I see the confusion. By "Taylor expansion of the source about the
previous value", we mean the value at the previous *sweep*, not the
previous *timestep*.
dmPhidPhi = 2 * W - 30 * (1 - 2 * phase) * enthalpy
S1 = dmPhidPhi * phase * (1 - phase) + mPhi * (1 - 2 * phase)
S0 = mPhi * phase * (1 - phase) - S1 * phase
eq = ImplicitDiffusionTerm(coeff=kappa) + S0 \
.. + ImplicitSourceTerm(coeff = S1)
So doesn't that mean that this source term still uses the current
value of phase
rather than the previous value?
It means that it uses the most recent value, but not the future
(implicit) value.
Or would it work better if I inserted
getOld() at each instance of phase?
No, this might work, but would probably have worse convergence
properties.
Even on p.105 where the source is added "explicitly", it seems to be
dependent
on the current value of phase rather than the previous value, making
it
implicit. Not sure if I am misunderstanding how this is supposed to
work.
Again, the current value is not the implicit value. With FiPy, you
never write the implicit value at all, it is taken care of
automatically by the form of the Term that you choose. For instance
\partial \phi / \partial t = \nabla\cdot\sin(\phi)\nabla \phi
is written in FiPy as
TransientTerm() == DiffusionTerm(coeff=sin(phase))
and is discretized at position j, timestep n, sweep k as something like
(phase_{j, n, k} - phase_{j, n-1}) / delta_t
= (sin(phase)_{j+1/2, n, k-1} (phase_{j+1, n, k} - phase_{j, n, k})
- sin(phase)_{j-1/2, n, k-1} (phase_{j, n, k} - phase_{j-1, n,
k})) / delta_x**2
so the only place that phase.getOld() appears is as phi_{j, n-1} in
the TransientTerm. When you write `phase`, you get phase_{j, n, k-1},
i.e., the previous sweep of the current timestep. You never write any
representation of phase_{j, n, k} at all; it is implicit in the
discretization and it is implicit in the code.
I hope this helps, but I concede that I may only have muddied the
waters. Please ask for further clarification if needed.
