← Back to team overview

dolfin team mailing list archive

Re: pydolfin

 

It works if you replace

  Function(element, mesh, h_e)

by

  Function(element, mesh, Vector())

You should get an error message when trying to create a Function from
a numpy array (unless Martin has done some tricks for initialization
of vector constants that I have missed).

Martin?

-- 
Anders


On Thu, Aug 21, 2008 at 08:44:32PM +0200, Evan Lezar wrote:
> Hi there
> 
> I have implemented a simple test case and have attached it.  In writing it, I
> made a discovery.  The vector that is passed as a list of coefficients to the
> function constructor (h_e in the example) needs to be a multidimensional array
> with the first dimension 1.  If one swaps the dimensions (1, num_dofs) on line
> 26, then the test crashes when field.eval is called.
> 
> However, according to me, the eval function does not work as expected.  Since
> the list of coefficients is zero and the function should be evaluated as the
> weighted sum of the basis functions using the coefficients as weights, this
> means that the function value returned for any position within the mesh should
> be zero.  However, the first component (let us say x) of H is always -1 and the
> second component is not changed and retains the value that was assigned to it
> before the call to evaluate.
> 
> I will see about getting a C++ test up and running to test the same thing.
> 
> Please let me know if there are any other questions.
> Evan
> 
> On Thu, Aug 21, 2008 at 6:06 PM, Evan Lezar <evanlezar@xxxxxxxxx> wrote:
> 
> 
> 
>     On Thu, Aug 21, 2008 at 5:49 PM, Ola Skavhaug <skavhaug@xxxxxxxxx> wrote:
> 
>         Evan Lezar skrev den 21/08-2008 f lgende:
>         > Thanks for the speedy replies from all of you.
>         >
>         > Firstly platform information
>         >
>         > I am running ubuntu hardy, gcc version 4.2.3, swig 1.3.33, python 2.5
>         with
>         > release versions of numpy and scipy.  I am building dolfin directly
>         from
>         > source with mpi, but have not enabled petsc or any of the other
>         optional
>         > packages.
>         >
>         > What other information do you require?
>         >
>         > In terms of my specific problem, see the following code snippet:
>         >
>         > from dolfin import *
>         > import numpy as N
>         >
>         > mesh = UnitSquare(4,4)
>         > element = FiniteElement("Nedelec", "triangle", 0)
>         >
>         > ### Code to assemble system matrices and solve an eigenvalue problem
>         >
>         > # h_e is an array of coefficients associated with the degrees of
>         freedom in
>         > the system obtained from the eigenvectors V
>         > h_e = V[:,sorted_index[first_index]]
>         >
>         > # setup a function to evaluate - the solution of the magnetic field
>         > field = Function(element, self.mesh, h_e)
>         >
>         > # arbitrary point to evaluate the field
>         > points = N.array((0.5,0.25))
>         > # array to store the field value at a point
>         > H = N.array((0.0,0.0))
>         >
>         > print "pre_eval"
>         >
>         > #evaluate the function
>         > field.eval(H, points)
>         >
>         > # Ouput the field at the point
>         > print "H",
>         > print H
>         >
>         >
>         > When I run this code, I get the following output:
>         >
>         > pre_eval
>         > [labby:25160] *** Process received signal ***
>         > [labby:25160] Signal: Segmentation fault (11)
>         > [labby:25160] Signal code: Address not mapped (1)
>         > [labby:25160] Failing at address: 0xde8f58ec
>         > [labby:25160] [ 0] [0xb7f81440]
>         > [labby:25160] [ 1] /usr/bin/python(PyEval_EvalCodeEx+0x6e7)
>         [0x80cb0d7]
>         > [labby:25160] [ 2] /usr/bin/python(PyEval_EvalFrameEx+0x565e)
>         [0x80c92de]
>         > [labby:25160] [ 3] /usr/bin/python(PyEval_EvalCodeEx+0x6e7)
>         [0x80cb0d7]
>         > [labby:25160] [ 4] /usr/bin/python(PyEval_EvalFrameEx+0x565e)
>         [0x80c92de]
>         > [labby:25160] [ 5] /usr/bin/python(PyEval_EvalCodeEx+0x6e7)
>         [0x80cb0d7]
>         > [labby:25160] [ 6] /usr/bin/python(PyEval_EvalCode+0x57) [0x80cb227]
>         > [labby:25160] [ 7] /usr/bin/python(PyRun_FileExFlags+0xf8)
>         [0x80ea6d8]
>         > [labby:25160] [ 8] /usr/bin/python(PyRun_SimpleFileExFlags+0x199)
>         > [0x80ea979]
>         > [labby:25160] [ 9] /usr/bin/python(Py_Main+0xa35) [0x8059335]
>         > [labby:25160] [10] /usr/bin/python(main+0x22) [0x80587f2]
>         > [labby:25160] [11] /lib/tls/i686/cmov/libc.so.6
>         (__libc_start_main+0xe0)
>         > [0xb7de8450]
>         > [labby:25160] [12] /usr/bin/python [0x8058761]
>         > [labby:25160] *** End of error message ***
>         >
>         > With a previous build I have gotten different output and was able to
>         > actually display the contents of H, but it would shill crash once the
>         method
>         > where the code was being run was terminated.  I logged a bug against
>         it at
>         > http://www.fenics.org/cgi-bin/bugzilla/show_bug.cgi?id=55
>         >
>         > As I said in my mail, I see myself using fenics in the long run and
>         would
>         > like to get involved as much as possible, but some pointers as to
>         where to
>         > start would be much apprciated.
> 
>         Have you implemented the C++ version as well, such that you know that
>         what you
>         do is correct? 
> 
> 
>         Clearly, we are not catching what's going wrong here (missing check of
>         some
>         sort).
> 
>         My guess is that the error occurs in site-packages/dolfin/function.py,
>         when
>         the Function is initialized. Can you boil down the example to a shorter
>         script
>         and post it to this list? I'll have a look at it.
> 
>         Ola
> 
> 
>     Hi
> 
>     I have not implemented a C++ version, but I will look into it - I am a
>     little more comfortable with python at the moment, but I suppose I need to
>     jump into the deep end :)
> 
>     I will write a simple example script and then attach that.
> 
>     Evan
> 
> 
> 

> # A script to test the evaluation of a vector function from pydoflin
> 
> __author__ = "Evan Lezar (evanlezar@xxxxxxxxx)"
> __date__ = "2008-08-21"
> __copyright__ = "Copyright (C) 2008 Evan Lezar"
> __license__  = "GNU LGPL Version 2.1"
> 
> from dolfin import *
> import numpy as N
> 
> def test_zero_function(order = 1):
>     mesh = UnitSquare(2,2)
>     element = FiniteElement("Nedelec", "triangle", 0)
>     
>     v = TestFunction(element)
>     u = TrialFunction(element)
>             
>     L = dot(v, u)
>     # assemble T to get the size of the system
>     T = assemble(L*dx, mesh)
>     
>     num_dofs = T.size(1)
>     
>     
>     # set the coeficients of all the dofs to 0
>     h_e = N.zeros((1,num_dofs), dtype=N.float64)
>     
>     # initialize the function
>     field = Function(element, mesh, h_e)
>             
>     point = N.array((0.0,0.0))
>     H = N.array((0.0,0.0))
>             
>     print "About to evaluate the function"
>             
>     field.eval(H, point)
>             
>     print "H",
>     print H
>     
>     assert ( H[0] == 0 and H[1] == 0 )
> 
> if __name__ == "__main__":
>     order = 0
>     test_zero_function(order)
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/mailman/listinfo/dolfin-dev

Attachment: signature.asc
Description: Digital signature


Follow ups

References