On Fri, Feb 22, 2013 at 12:10:51PM -0500, Daniel Wheeler wrote:
> On Wed, Feb 20, 2013 at 1:54 PM, Pierre de Buyl <[email protected]> wrote:
> > What I want is n(theta, r) where r is fixed. I can obtain this value by
> > calling
> > the variable "n" with a list of coordinates. I do this by adding the
> > following
> > code:
> >
> > def radial_value(th,r,order=1):
> > x = r*np.sin(th)
> > y = r*np.cos(th)
> > return n((x, y), order=order)
> >
> > and I plot
> > plot(theta, radial_value(theta, 1.5))
> > here, 1.5 is arbitrary. I need the value at different radii depending on
> > other parameters.
> >
> > Around theta=0 and theta=pi, which are close to the left of my domain
> > (cylindrical radius=0), I expect zero flux (-> flat function). However, the
> > slope is non-negligible. I attach a plot of n(theta,r=1.5), for sigma=9 and
> > v=-0.01).
> >
> > The result doesn't change much if I shift the mesh with
> > + ((1e-5,),(0.,))
> >
> >
> Does it get any better as you increase the perturbation?
For the range 1e-5 to 1e-3, there is not a lot of change, but the solution is
manifestly not appropriate. The contour lines close to r=0 are not perpendicular
to the boundary. At 1e-2, the solution does not display the expected asymmetry
anymore. (the top of the sphere shows an unnatural depletion of the
concentration with respect to what I expect from the solution).
> > Could it be the the multiplication of the cell volume and face area by the
> > radius cause of bad handling of the zero flux BC on the left of the domain
> > with
> > radius=0.
>
> Yes. I have experienced this issue before. Play with the shift and see if
> things approve. I'm curious as to why you aren't using the
> CylindricalGrid2D class? Try using that at least to test against.
I haven't used CylindricalGrid2D because the inner boundary in my problem is a
sphere (or a circle in polar coordinates) and that FiPy does not support
internal boundaries. I tried to use it with a constraint that "n=0 for r<1" but
it does produce very irregular results.
> > In that case, I couldn't solve satisfactorily my axisymmetric problem, is
> > this right?
> >
>
> I'm sure it is possible, but I haven't tried to reproduce your work yet and
> really figure out what is happening.
>
> Try the CylindricalGrid2D just to make sure that the various changes for
> the axisymmetric problem are being handled correctly in your code and then
> we can go from there.
> Also see if things improve as the shift is increased away from zero and
> that the solutions are adequately converged. That's the easy diagnostics.
> Confirm those first.
As written above, I cannot use CylindricalGrid2D. Except if I misunderstood
something and I can input an arbitrary mesh into a CylindricalGrid2D. The shift,
unfortunately, does not give nice results.
Thanks for your comments. If anyone has other suggestions, I am all ears.
Following my 2D experience, I went on to define a 3D mesh in gmsh that is made
of two concentric spheres, one at r=sigma and one at r=20*sigma ("far away").
The mesh is only defined on a slice of an eighth of a sphere. This gives very
good results, at the expense of more memory and more CPU time. As the precision
needs to be higher close to r=sigma, I defined the typical cell size in gmsh.
I use 0.1*sigma at r=sigma and sigma at r=3*sigma and r=20*sigma. I created an
extra shell at r=3*sigma just for the purpose of assigning different cell sizes.
I am just unsure of the precision of the solution. What I need, as a final
result, is
\int_0^\pi n(R,theta) sin(theta) cos(theta) dtheta
at some fixed value R. In the absence of flow, the solution should be perfeclty
symmetrical and the above integral should equal 0. I get, depending on the size
of the sphere, a significantly non-zero result (from 0.1 to 0.3).
Is there any typical way of checking wether the solution is precise enough for
this kind of computation?
Cheers,
Pierre de Buyl
PS: if needed, I'll send the new code and geometries.
>
> --
> Daniel Wheeler
_______________________________________________
fipy mailing list
[email protected]
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]