Hi Michael,

Reply is  below.

On Fri, Feb 22, 2008 at 2:13 PM, Michael Waters <[EMAIL PROTECTED]> wrote:
>
>  I am trying to compute an interface angle in the same way that was done
>  in the paper "Extending phase field models of solidification to
>  polycrystalline materials"  by James Warren, Ryo Kobayashi, Alexander
>  Lobkovsky, Craig Carter. I used a snippet of code from the fipy
>  anisotropy example (section 9.4), but this returns an error that I'm not
>  sure how to correct. If anyone knows a better way to do this or how to
>  work around the error, I'd love to hear.
>
>  This is the code that I'm using,  I'd like it if I could keep theta as a
>  modular variable:
>
>  phase = CellVariable(name='phase field', mesh=mesh, value = 0.0)
>  theta = ModularVariable( name='theta', mesh=mesh, value=-pi + 0.0001,
>  hasOld=1 )
>
>  phaseY = phase.getFaceGrad().dot((0, 1))*1.0
>  phaseX = phase.getFaceGrad().dot((1, 0))*1.0
>  psi= arctan2(phaseY, phaseX)  + theta

I think this error is simply due to doing mixed operations with face
and cell values. It has
nothing to do with the fact that theta is a ModularVariable. The error
message is not particularly helpful though.

To fix it do the following. Make phaseY and phaseX cell values

    phaseY = phase.getGrad().dot((0, 1))
    phaseX = phase.getGrad().dot((1, 0))

Then

    psi= arctan2(phaseY, phaseX)  + theta

works. FYI is not a ModularVariable. Do you want it to be?

Looking into this I realize that ModularVariable does not have
addition defined. This doesn't work

    >>> print theta + theta
    [-6.28298531 -6.28298531 -6.28298531 -6.28298531] 1

This does work

    >>> print theta - (-theta)
   [ 0.0002  0.0002  0.0002  0.0002] 1

which I believe is correct. Also this works

    >>> print 4 * theta - theta
    [-3.14129265 -3.14129265 -3.14129265 -3.14129265] 1

I am not sure why the addition operator was never included. Probably
there is a reason, but I can't remember.
I am now worried that the other operators need to be fixed. What about this?

    >>> print theta * 4
    [-12.56597061 -12.56597061 -12.56597061 -12.56597061] 1

Anyhow, if you are only doing limited operations, then I think you can
work around it.

Is ModularVariable worth fixing up to work properly? Or is there a
better way such as just defining a calculation
and then just using a mod operator like this

    >>> a = CellVariable(..)
    >>> b = CellVariable(..)
    >>> c = (a + b) % numerix.pi

Seems more natural to me. The % operator should return a BinOp, which
is currently not the case (Is it?).
We would then have to do something about the weird face grads that
need to be done in some of
the examples in a modular fashion, but they probably could be defined
as special variables or even in the script.
I think this makes much more sense. Of course, the following,

    >>> d = ((a + b) % numerix.pi) + a

would not be a modular variable, which also makes sense to me.

Cheers

--
Daniel Wheeler

Reply via email to