One other point is that maybe phi.getGrad()[0] * TransientTerm(1.) works now. I'm not sure. Seems to work:
In [7]: m = Grid1D(nx=10) In [8]: v = CellVariable(mesh=m) In [9]: v.getGrad()[0] * TransientTerm() Out[9]: _gauss_grad[index] * TransientTerm(coeff=1.0) Maybe try that as well. On Fri, Mar 19, 2010 at 12:00 PM, Daniel Wheeler <[email protected]> wrote: > If alpha and beta are zero, which they seem to be in the equation that > you wrote down. Then you just have the three terms, which is easy to > set up. The source term will just have phi.getGrad()[0] - > phi.getOld().getGrad()[0] as the ImplictSourceTerm coefficient. It > doesn't have to be implicit either, maybe just start with an explicit > source at first, you'll have to sweep to convergence at each time step > and think carefully about the boundary conditions, which I haven't > discussed. > > On Fri, Mar 19, 2010 at 11:48 AM, Daniel Wheeler > <[email protected]> wrote: >> Hi Toni, >> >> I'm not sure about this problem. One thought is that you could try solving >> >> \left( \phi_x^* \phi \right)_t = \phi_{xx} + \left( \left[ \beta - >> \frac{\alpha}{2} \phi^*\right] \phi \right)_x + \phi_{xt}^* \phi >> >> where the * values are old sweeps or old time steps. The non-starred >> values are implicit. Basically multiply through by \phi_x. Get an >> awkward \phi_{xt} term, but maybe that will be okay. The value >> \phi_x^* will be part of the transient coefficient. I assume you >> understand that it is best to get coefficients inside the derivatives >> for FV schemes in general. >> >> Cheers >> >> On Thu, Mar 18, 2010 at 12:58 PM, Ivas Toni <[email protected]> wrote: >>> Hi, >>> >>> >>> >>> I would like to solve “infinite” diffusion PDE (3d): >>> >>> 1/(a\nabla\phi) a \otimes a \frac{\partial^2}{\partial x_i \partial x_j} >>> \phi -\alpha \phi + \beta = \partial \phi /\partial t >>> >>> >>> >>> where “a” is a vector and a \otimes a is outer product of this vector >>> (symmetric matrix a_ij). This equation reduces to: >>> >>> \frac{\partial^2 \phi / \partial x^2 }{\partial \phi / \partial x} -\alpha >>> \phi + \beta = \partial \phi /\partial t >>> >>> for one dimensional case and a scalar. >>> >>> >>> >>> Equilibrium solution for 1d case is reducing to Burgers equation which can >>> be solved analytically(Hopf-Cole transformation). The solution is >>> >>> -tanh(2\alpha/a(x-\delta)) where \delta is position of the interface >>> (boundary condition are +1 for x=-\inf and -1 for x=+\inf) . This represents >>> >>> a traveling wave solution with certain velocity \delta=velocity * time. >>> >>> >>> >>> I tried to solve 1d equation with Fipy something like: >>> >>> >>> >>> from fipy import * >>> >>> >>> >>> nx = 40 >>> >>> Lx = 2.5 * nx / 100. >>> >>> dx = Lx / nx >>> >>> mesh = Grid1D(dx=dx, nx=nx) >>> >>> >>> >>> timeStepDuration = 0.0002 >>> >>> phaseTransientCoeff = 1 >>> >>> thetaSmallValue = 1e-6 >>> >>> beta = 1e5 >>> >>> mu = 1e3 >>> >>> gamma = 1000 >>> >>> s = 1 >>> >>> >>> >>> phase = CellVariable(name='phase field',mesh=mesh, value=1.) >>> >>> x = mesh.getCellCenters()[0] >>> >>> phase.setValue(x*(1-x),where=mesh.getCellCenters()[0] < Lx ) >>> >>> #phase.setValue(1./3,where=mesh.getCellCenters()[0] > Lx /3.) >>> >>> #phase.setValue(1-x,where=mesh.getCellCenters()[0] > 2* Lx / 3.) >>> >>> print phase >>> >>> gradMag = phase.getFaceGrad().getMag() >>> >>> eps = 1./gamma/10. >>> >>> gradMag+=(gradMag < eps) * eps >>> >>> IGamma = (gradMag > 1./gamma) * (1/gradMag-gamma) + gamma >>> >>> #print IGamma >>> >>> diffusionCoeff = s*IGamma >>> >>> >>> >>> phaseEq = TransientTerm(phaseTransientCoeff) == >>> ImplicitDiffusionTerm(diffusionCoeff) >>> >>> #print diffusionCoeff >>> >>> if __name__ == '__main__': >>> >>> phaseViewer = Viewer(vars=phase, datamin=0., datamax=2) >>> >>> phaseViewer.plot() >>> >>> steps = 100 >>> >>> for i in range(steps): >>> >>> phase.updateOld(); >>> >>> phaseEq.solve(phase,dt=timeStepDuration) >>> >>> if __name__ =='__main__': >>> >>> phaseViewer.plot() >>> >>> >>> >>> If I fix the boundary condition in fipy to 1 an -1 at the domain boundaries, >>> the numerical solution converges to step like function. >>> >>> It seems very difficult to control the terms involving the 1/grad(Phi). Is >>> it risky to handle this term as a very high diffusion coefficient ??? >>> >>> The 1d solution of the previous equation has exact equilibrium solution as >>> tangent hyperbolic and not step function. >>> >>> >>> >>> I would be grateful if somebody can help me with this problem. >>> >>> >>> >>> Thank you, >>> >>> Toni Ivas >>> >>> >>> >>> >> >> >> >> -- >> Daniel Wheeler >> > > > > -- > Daniel Wheeler > -- Daniel Wheeler
