On Thu, Aug 30, 2012 at 3:42 PM, Scott Norris <[email protected]> wrote:
> > I'm particularly interested in using FiPy to study weakly-nonlinear >>> equations associated with pattern formation in multiple dimensions. A >>> classical problem there is the Kuramoto-Sivashinsky equation: >>> >>> h_t = - \nabla^2 h - \nabla^4 h - (1/2) | \nabla h |**2 >>> >>> where h is the dependent variable. >> >> >> One possible way to do this (with the coupled equation capability) would >> be >> >> h_t = - \nabla^2 h - \nabla^2 f - (1/2) \nabla \cdot \left(h \nabla >> h \right) + h f >> >> f = \nabla^2 h >> >> I can't say how well this will work, but it seems logical. The second to >> last term will be a convection term with the velocity "\nabla h" and the >> last term would be an implicit source term on f (h will be the >> coefficient), since it is based on the derivatives of h. >> > > > Thanks so much, and for the links to documentation on coupled equations. > So far, this single equation seems to work: > > mesh = PeriodicGrid1D(nx = 100, dx = 0.5) > u = CellVariable(name="u", mesh=mesh) > v = CellVariable(name="v=lapl u", mesh=mesh) > eqn0 = TransientTerm(var=u) == - DiffusionTerm(coeff=1., var=u) - > DiffusionTerm(coeff=(1., 1.), var=u) \ > + PowerLawConvectionTerm(coeff=.5*u.faceGrad, var=u) - > ImplicitSourceTerm(coeff=.5*u.faceGrad.divergence, var=u) > > But replacing the last line with this coupled equation fails: > > eqn1 = TransientTerm(var=u) == - DiffusionTerm(coeff=1., var=u) - > DiffusionTerm(coeff=1., var=v) \ > + PowerLawConvectionTerm(coeff=.5*u.faceGrad, var=u) - > ImplicitSourceTerm(coeff=.5*v, var=u) > eqn2 = ImplicitSourceTerm(var=v) == DiffusionTerm(1, var=u) > eqn = eqn1 & eqn2 > > Couldn't "PowerLawConvectionTerm(coeff=.5*u.faceGrad, var=u)" be "DiffusionTerm(coeff=.5 * u, var=u)"? That would make more sense. It removes the implicitness from second order to first order. It is mostly okay to have the explicitness in the coefficients of the diffusion term. That has less restriction on the time step. > It seems to run, but quickly fails with message alluding to: > > Warning: overflow encountered in square > Warning: overflow encountered in square > Warning: overflow encountered in square > Warning: overflow encountered in square > Warning: overflow encountered in square > (many more of this, then:) > ... > Traceback (most recent call last): > File "dKS-1D.py", line 76, in <module> > eqn.solve(dt=dt) > File > "/usr/local/lib/python2.7/dist-packages/FiPy-3.0_dev5303-py2.7.egg/fipy/terms/term.py", > line 211, in solve > solver._solve() > File > "/usr/local/lib/python2.7/dist-packages/FiPy-3.0_dev5303-py2.7.egg/fipy/solvers/pysparse/pysparseSolver.py", > line 103, in _solve > self._solve_(self.matrix, array, self.RHSvector) > File > "/usr/local/lib/python2.7/dist-packages/FiPy-3.0_dev5303-py2.7.egg/fipy/solvers/pysparse/linearLUSolver.py", > line 87, in _solve_ > b = b * (1 / maxdiag) > ValueError: cannot convert float NaN to integer > I don't think I can fiigure that out without debugging. It looks like a diagonal entry is zero, which you can't have. Obviously, the transient and implicit source should prevent zero diagonal, but then other terms are contributing. Make the change with the convection term and see if you still have the same error. I suspect it might deal with this issue otherwise I'll need the entire script to debug. > > I note either way I try to do it, I still seem to have to use faceGrad() > in the ConvectionTerm for $ \nabla (u \nabla u) $. Is this destroying > the implicit nature of the method? Is there a way around it? > Exactly, swap it for a diffusion term. Having it as a convection term is introducing second order explicitness in place of first order explicitness. -- Daniel Wheeler
_______________________________________________ fipy mailing list [email protected] http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
