You can do phi.grad.dot(phi.grad), but it's fully explicit, so best avoided.
On Jan 27, 2015, at 4:04 AM, "Ronghai Wu" <[email protected]<mailto:[email protected]>> wrote: Thank you so much. I just realized the solution this morning, which agrees with your answer. It is a stupid question. But I got stuck for days because I had one-track mind that the nabla terms must be fully expanded in following way: \nabla^2 \phi^2 = 2 (\nabla \phi) \cdot (\nabla \phi) + 2 \phi (\nabla^2 \phi) Then I stuck at implementing (\nabla \phi) \cdot (\nabla \phi) which is not applicable in Fipy. Thanks again. On 26/01/15 20:25, Guyer, Jonathan E. Dr. wrote: An _OperatorVariable is an internal FiPy object that represents a mathematical operation, such as phi**2. By saying 'var=phi**2' you are telling FiPy to solve for phi**2 rather than for phi. The error is telling you that FiPy doesn't know how to do that. DiffusionTerm(var=???, coeff=D) represents \nabla\cdot D \nabla ???, where ??? must be a solution variable, not a function of a solution variable. FiPy builds systems of linear equations (like any other finite-(difference|volume|element) code), so each of the terms must be a linear expression of a solution variable. The thing to do here is carry the derivatives through to put these expressions in a form FiPy understands: d\phi/d t = \nabla^2 \phi + \nabla^2 \phi^2 + \nabla^2 \phi^3 = \nabla \cdot \nabla \phi + \nabla\cdot(2\phi \nabla \phi) + \nabla\cdot(3\phi^2 \nabla \phi) There's no reason not to condense this to d\phi/d t = \nabla \cdot[(1 + 2\phi + 3\phi^2)\nabla \phi] or PHI = phi.faceValue TransientTerm(coeff=1., var=phi) == DiffusionTerm(coeff=1 + 2*PHI + 3*PHI**2, var=phi) On Jan 26, 2015, at 12:09 PM, Ronghai Wu <[email protected]><mailto:[email protected]> wrote: Dear Fipy developers and users, I have a strange PDE in 2D "dphi/dt = Laplacian(phi) + Laplacian(phi**2) + Laplacian(phi**3)", which I implement in the following way, but raise typeerror "The value of an `_OperatorVariable` cannot be assigned". Can anyone help me? Thanks. dx = 0.1 dy = dx nx = 100 ny = nx mesh = Grid2D(dx=dx, dy=dy, nx=nx, ny=ny) x = mesh.cellCenters[0] y = mesh.cellCenters[1] phi = CellVariable(mesh=mesh, hasOld=1, value=x*y) eq = TransientTerm(var=phi, coeff=1.) == ImplicitDiffusionTerm(var=phi, coeff=1.) + ImplicitDiffusionTerm(var=phi**2, coeff=1.) + ImplicitDiffusionTerm(var=phi**3, coeff=1.) viewer = Viewer(vars=phi) for step in range(10): phi.updateOld() res = 1.e5 while res > 1.e-1: res = eq.sweep(dt=1.e-3) viewer.plot() Regards Ronghai _______________________________________________ fipy mailing list [email protected]<mailto:[email protected]> http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ] _______________________________________________ fipy mailing list [email protected]<mailto:[email protected]> http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ] -- ------------------------------------------ Ronghai Wu Institute of Materials Simulation (WW8) Department of Materials Science and Engineering University of Erlangen-N?rnberg Dr.-Mack-Str. 77, 90762 F?rth, Germany Tel. +49 (0)911 65078-65064 _______________________________________________ fipy mailing list [email protected]<mailto:[email protected]> http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
_______________________________________________ fipy mailing list [email protected] http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
