Peter - I'm glad sweeping was able to get you a more reasonable solution (although I'm not seeing the non-linearity that would require it). `AdvectionTerm` only ever gets used in our level set examples, which I'm not intimately familiar with. I don't know of any fundamental reason it wouldn't work, but it's not something we've broadly tested.
As far as how you've rearranged terms, I think you've taken a reasonable approach. You could alternatively do something like (2 * L**2 * f.grad.dot([[1]])) (not tested), but that's probably not any different and possibly less implicit than the AdvectionTerm. Generally, though, when terms arise that are not in one of FiPy's canonical forms, it's good to ask why. `L^2 d/dL (L df/dL)` and, even more so, `L^2 d/dL (D/L^2 df/dL)` in Eq (1) of the paper you cite in the code is strongly reminiscent of a laplacian operator in spherical coordinates: `(1/r^2) d/dr (r^2 df/dr)`, but the r^2 and (1/r^2) are swapped. Far be it from me to question an eminence like Van Allen, but are you certain that equation is right? Is L some sort of scaling that is inverse to radius by any chance? We provide 1D and 2D CylindricalGrid meshes and it would not be hard to write a SphericalGrid1D class, if that's called for. The reason to do that is that our CylindricalGrid1D class accounts for the curvature of the faces and the expanding volume of the cells with increasing radius. The laplacian operator automatically takes the form `(1/r) d/dr (r df/dr)` when you do this. Even if the operator in question is not in spherical coordinates, I suspect the correct answer is to build `L^2 d/dL (1/L^2 ...)` into the mesh somehow. Once we understand why it looks like that, it shouldn't be too hard to create the appropriate mesh class. - Jon P.S. Thanks very much for your pull request! I've merged it to the develop branch. On Oct 19, 2015, at 7:05 PM, Peter Williams <[email protected]> wrote: > Hi, > > I'm trying to use fipy to solve some differential equations, and I'm not sure > how to work with the particular equation I've got at hand. I'm a total newbie > at this, so I'm probably doing something silly, but I can't figure out how to > proceed so I thought I'd ask for some help. > > The attached Python file demonstrates my problem and what I've done so far. > The key issue is that I have an equation with diffusion terms that have > spatially-dependent prefactors; in the simplest case, something like > > L^2 d/dL (L df/dL) > > where L is a spatial variable and I'm trying to determine f(L). As far as I > can tell, this kind of term is not directly allowed by fipy, and I get a > TermMultiplyError if I try to construct it. > > In the simplest case, demonstrated in the attachment, I can divide out the > prefactor and fipy works great. However, in my full problem, there will be > multiple diffusion terms (in multiple dimensions) with differing prefactors, > so I won't be able to make them all disappear. > > I'm pretty sure that I can generically re-express any prefactored diffusion > term as the sum of a prefactor-free diffusion term and an advection term. > However, I've tried to implement this with fipy, and it doesn't work. The > result from the solver is as if I didn't include the advection term at all. > Again this is demonstrated in the attached file. > > So, I have two questions: > > (1) I'm probably doing something wrong with the advection term in the > attached sample. Can anyone spot the problem? > > (2) It'd be nice not to have to do the error-prone elimination of my > diffusion prefactors. Are there any ways to avoid this? > > Thanks for any tips, > > Peter > <testcase.py>_______________________________________________ > fipy mailing list > [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 ]
