> OK. My comments still stand. The discrepancies do not have the appearance of
> a "bad" solution, but rather that they are not solutions to the same equation
> and/or boundary conditions.
>
> Are you solving the same equation and boundary conditions? How do you know?
So the equation I am trying to solve is:
\frac{1}{D_r} * \frac{\partial \varphi(\alpha, t)}{\partial t} =
\frac{\partial^2 \varphi(\alpha, t)}{\partial \alpha^2} +
\frac{\partial}{\partial \alpha} (\frac{1}{K * T} * \frac{\partial H(\alpha,
t)}{\partial \alpha} * \varphi(\alpha, t))
Where D_r, K, T are constants. My FiPy equation is:
equation = TransientTerm(coeff = (1. / self.p.Dr)) == \
DiffusionTerm() + ExponentialConvectionTerm(coeff = (1 / (self.p.K *
self.p.temp)) * H.getFaceGrad())
\alpha has the domain [-pi / 2; pi / 2] so I have:
alphaStepsCount = int(math.pi / self.p.alphaStep)
alphaMesh = Grid1D(nx = alphaStepsCount, dx = self.p.alphaStep)
For calculating H I shift alpha by pi / 2 to get the correct domain.
The boundary conditions of my problem are:
\varphi(\alpha, 0) = p_0(\alpha)
\varphi(-\frac{\pi}{2}, 0) = \varphi(\frac{\pi}{2}, 0) = 0
Where p_0 is a known function.
In FiPy I represent them as:
prob.setValue(self._p_zero(alpha - math.pi / 2))
boundaries = (FixedValue(faces=alphaMesh.getFacesRight(), value=0.0), \
FixedValue(faces=alphaMesh.getFacesLeft(), value=0.0))
The reference I am comparing my result to is one computed my Mathematica which
has been "verified" by experiment. What is odd is that Mathematica uses only
455 gridpoints for alpha, while I use up to 30 000 and still get worse result.
I guess that really points to a non-numerical problem, but I have no idea
which...
Thanks for your help,
Matej
_______________________________________________
fipy mailing list
[email protected]
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
