Vitesh -
I get this warning from SciPy and PySparse, but not from PETSc or Trilinos.
Evidently, there’s something LU doesn’t like about it. I don’t think it
actually matters, though. The solution proceeds anyway. Often, the initial
matrix can raise objections that don’t cause the overall evolution of the
problem to fail.
A bigger issue is that you have several cases where you are casting away the
“Variable” nature of your variables. Everywhere where you write
`something.value`, you get the instantaneous value of `something`; it won’t
update as the problem involves. Just use `something`.
- return 1. +
beta_r*(1.-np.exp(-beta_s*phi.value*one_over(phi.grad.mag.value)))
+ return 1. + beta_r*(1.-np.exp(-beta_s*phi*one_over(phi.grad.mag)))
- diffusionCoeff = (g(phiFace,a).value * IGamma.value + 2.0*h(phiFace,d,e))
+ diffusionCoeff = (g(phiFace,a) * IGamma + 2.0*h(phiFace,d,e))
More subtly, appending `()` to a FiPy Variable also returns its instantaneous
value (`something()` is equivalent to `something.value`).
def build_E_st_equation(phi,theta,C_d,E_st):
- prefactor = -E_st*((phi - phi.old()) > 0.0)
+ prefactor = -E_st*((phi - phi.old) > 0.0)
prefactor = prefactor*localize_E(theta)
- return TransientTerm() == prefactor*(phi - phi.old())/dt
+ return TransientTerm() == prefactor*(phi - phi.old)/dt
Finally, you need to be careful with calling functions. You want to be sure
that your functions return a Variable expression, not just the instantaneous
value. Most of your helper functions are OK, but `one_over` won’t worked as
desired because `np.where` and `np.isclose` discard the Variable character and
just return numpy arrays. I recommend:
- return np.where(np.isclose(f,0.0,atol=gamma),1.0/gamma,np.tanh(f/gamma)/f)
+ return 1.0 / (f + gamma)
I suppose you could do
+ return 1.0 / (f + (abs(f) < gamma) * gamma)
but I doubt it makes much difference.
After that, the solution evolves, but it’s not stable.
Dropping the time step to 2e-4 or 2e-5 gets a smooth evolution toward something
like Fig. 2(b) in Abrivard et al. (I don’t know if it’s supposed to, as I
haven’t actually read that paper, but it does).
- Jon
On Feb 1, 2021, at 6:20 AM, vitesh shah
<[email protected]<mailto:[email protected]>> wrote:
Hello,
I am currently using Fipy, to solve a set of phase field equations. These
equations are based on the KWC model but with an additional extension made to
it, which enables it to account for stored energy of a material to be a driving
force.
Currently, I have three equations which I attempt to solve using Fipy:
1. Phase equation
2. Theta (misorientation) equation
3. Equation for the stored energy.
If I solve without the stored energy, the Fipy solver works fine and I get
expected behavior (in 1D at least). But when I also attempt to solve equation
for stored energy, I get a runtime warning:
“/nethome/v.shah/.local/lib/python3.8/site-packages/fipy/solvers/scipy/linearLUSolver.py:41:
RuntimeWarning: invalid value encountered in double_scalars if
(numerix.sqrt(numerix.sum(errorVector**2)) / error0) <= self.tolerance:”
So I suspect that solver is not handling these equations properly. And this
affects the solution that we get.
We expect that for 1D case, the phase boundary should move when stored energy
is introduced (in normal phase field it should remain stationary). But this
behavior is not observed.
Also when I compare the analytical rates and the rates that I observe from Fipy
and they have quite a different behavior.
I am attaching a jupyter notebook here, so that the whole code (with references
to the papers) can be looked at. (link: if the attachment doesnt work)
It would be great, if you can give me some suggestions or point out some
mistakes that I might have made.
With regards,
Vitesh Shah
--
To unsubscribe from this group, send email to
[email protected]<mailto:[email protected]>
View this message at https://list.nist.gov/fipy
---
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected]<mailto:[email protected]>.
<KWC_with_stored_energy_example.py>
--
To unsubscribe from this group, send email to [email protected]
View this message at https://list.nist.gov/fipy
---
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
