[yt-users] generic computation of derivatives

[David Collins]

Hi, Wolfram--

I'm sure there are several ways to do this, but this is what I've done.

I defined a function that takes *data* and a list of field names, then
I call that function elsewhere.  So for instance for terms of the form
v1 dot grad v2
or
vi d (vj) /dxi

so I have wrote the following function (which I think might be poorly named)

def tensor_derivative(data,fields1,fields2):
    """performs (A \cdot \nabla) B,
    Ai di Bj
    where di = d/d[x,y,z]
    Ai in fields1, Bj in fields2
    Assumes centered difference"""
    iM1 = slice(None,-2)
    iP1 = slice(2,None)
    all = slice(1,-1)
    all_all=[all]*3
    out = [na.zeros(data[fields1[0]].shape)]*3
    dxi=1./(2*data.dds )
    for j,fj in enumerate(fields2):
        for i, fi in enumerate(fields1):
            Left,Right = [all]*3,[all]*3
            Right[i] = iP1
            Left[i] = iM1
            out[j][all_all] += data[fi][all_all]*(data[fj][ Right ]-
data[fj][ Left] )*dxi[i]
    return out

This returns an array of results, I typically take the norm.  For
instance, the norm of the magnetic curvature || Bi Di Bj || as

def _curvatureNorm(field,data):
    Bfields=['Bxhat','Byhat','Bzhat']
    curvature = tensor_derivative(data,Bfields,Bfields)
    return na.sqrt( curvature[0]**2 + curvature[1]**2 + curvature[2]**2 )

add_field('CurvatureNorm',function=_curvatureNorm,take_log=False,
        validators=[ValidateSpatial(1,['Bxhat','Byhat','Bzhat'])])

and Bxhat is the unit vector along Bx, though anything can be used here.

Hope that helps!
d.

On Mon, Sep 12, 2011 at 9:39 AM, Wolfram Schmidt
<schmidt at astro.physik.uni-goettingen.de> wrote:
> Hi everyone,
>
> how can I generically define the derivatives of fields in a yt script (or is
> there already something to do that in yt)?
>
> I want to define, say, second-order centered differences of a field without
> referring to data["particular field"] but with a place holder for just any
> field (like in a Fortran or C subroutine) so that I can use this single
> piece of code for several different fields.
>
> Thanks,
> Wolfram
>
> _______________________________________________
> yt-users mailing list
> yt-users at lists.spacepope.org
> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>



-- 
Sent from my computer.