dolfin team mailing list archive

["Evan Lezar" <evanlezar@xxxxxxxxx>]

On Tue, Aug 19, 2008 at 6:18 PM, Anders Logg <logg@xxxxxxxxx> wrote:

> On Tue, Aug 19, 2008 at 05:33:21PM +0200, Evan Lezar wrote:
> >
> >
> > 2008/8/11 Anders Logg <logg@xxxxxxxxx>
> >
> >     On Sun, Aug 10, 2008 at 03:54:25PM +0200, Evan Lezar wrote:
> >     > Hi there
> >     >
> >     > After some time away I have once again been able to try and get do
> grips
> >     with
> >     > dolfin.  By a bit of a round about way I was able to solve a
> >     magnetic-field
> >     > based eigenvalue problem giving me the cutoff wavenumbers for a
> square
> >     > waveguide structure - it is still, however, a work in progress.
> >     >
> >     > One question I do have is regarding the visualization of the
> results I
> >     have
> >     > obtained.  To this end I have two questions - the first is simply
> the
> >     > visualization of the basis functions (vector Nedelec) themselves
> (see the
> >     > attached image taken from Finite Element Method for
> Electromagnetics by
> >     Volakis
> >     > et al) - how could this be achieved in dolfin?  Put another way -
> how can
> >     I set
> >     > the coefficient for a single basis function equal to 1 and zero the
> rest
> >     and
> >     > display the result?
> >
> >     Marie Rognes has code for doing just this (plotbases.py in the
> >     bdm-2007 repository). I don't want to post it myself, but I imagine
> >     she would be willing to post the code. Marie?
> >
> >     > Then the second question builds on this.  If I have the
> coefficients
> >     associated
> >     > with all the basis functions in a system, how would I set up a
> function
> >     to
> >     > calculate the weighted sum of the vector basis functions anywhere
> in the
> >     mesh?
> >     >
> >     > I am using pydolfin.
> >     >
> >     > If anything is unclear, please contact me so that I can clarify
> matters.
> >
> >     If you have a Function u, you can just call u.eval() to evaluate it
> at
> >     any point x. This will locate the cell to which x belongs and compute
> >     the linear combination. Take a look at the demo in
> demo/function/python/.
> >
> >
> > Anders
> >
> > I am communicating with Marie regarding the visualization of the basis
> > functions - I am still running into problems with the function
> evaluation.
> > Consider the following:
> >
> > mesh = UnitSquare(1,1)
> > element = FiniteElement("Nedelec", "triangle", 0)
> >
> > # obtain a list of coefficients for the basis fuctions in the finite
> element
> > system
> >
> > f = Function(element, mesh, coefficient_list)
> >
> > Now this is where the problems start.  As far as I am aware, f should be
> a
> > vector function which could be evaluated at any point (x,y), so for
> argument's
> > sake, let
> > p = array((0.5,0.5))
> >
> > then if I execute the following
> >
> > F_val = f.eval(p)
> >
> > Then F_val should contain the x and y components of the function as a
> weighted
> > sum of the basis functions?
> >
> > Please let me know if I am understanding this correctly, or if I am doing
> > something wrong.
> >
> > Evan
>
> When evaluating a Function in Python, you have to do something like
> this:
>
>  x = array((0.5, 0.5))
>   values = array((0.0, 0.0))
>  u.eval(values, x)
>
> Then values will be an array with the two values for the vector-valued
> function.
>

Ok, I think I have it working a little better now.  I am actually able to
get function values out.

I do however get an error for the following code

        points = N.array((0.5,0.5))
        H = N.array((0.0,0.0),)

        print "pre_eval"

        f.eval(H, points)

        print "H",
        print H

The output is:

pre_eval
H [ 0.02926315  0.31694778]
*** glibc detected *** /usr/bin/python: free(): invalid next size (fast):
0x08805f98 ***

And adding an additional print after print H causes the program to not print
the values of H.

I assume this is a bug?  I could have a look into fixing it ...

But I am already a lot further than what I was :)


>
> Ideally, one should be able to do
>
>  x = array((0.5, 0.5))
>   values = u.eval(x)
>
> or even
>
>  value = u(x)
>
> but this hasn't been fixed in the interface yet.
>
> --
> Anders
>
> -----BEGIN PGP SIGNATURE-----
> Version: GnuPG v1.4.6 (GNU/Linux)
>
> iD8DBQFIqvI7TuwUCDsYZdERAn8FAJ9nfLnoYiwiKTXnLGvWgQF35albHACgi9rM
> hLT5IHVEBfvVKbVrc4QNY0E=
> =3B7n
> -----END PGP SIGNATURE-----
>
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/mailman/listinfo/dolfin-dev
>
>