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 
<vitesh.manches...@gmail.com<mailto:vitesh.manches...@gmail.com>> 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 
fipy+unsubscr...@list.nist.gov<mailto:fipy+unsubscr...@list.nist.gov>

View this message at https://list.nist.gov/fipy
---
To unsubscribe from this group and stop receiving emails from it, send an email 
to fipy+unsubscr...@list.nist.gov<mailto:fipy+unsubscr...@list.nist.gov>.
<KWC_with_stored_energy_example.py>

-- 
To unsubscribe from this group, send email to fipy+unsubscr...@list.nist.gov

View this message at https://list.nist.gov/fipy
--- 
To unsubscribe from this group and stop receiving emails from it, send an email 
to fipy+unsubscr...@list.nist.gov.

Reply via email to