Where does this show up? -- Anders
On Mon, Nov 30, 2009 at 09:44:24AM +0000, Garth N. Wells wrote: > This change appears to have broken compiled expressions. > > Garth > > -------- Original Message -------- > Subject: [Branch ~dolfin-core/dolfin/main] Rev 4352: Work on > Expression in PyDOLFIN > Date: Mon, 30 Nov 2009 08:11:09 -0000 > From: [email protected] > Reply-To: [email protected] > To: Garth Wells <[email protected]> > > ------------------------------------------------------------ > revno: 4352 > committer: Johan Hake <[email protected]> > branch nick: dolfin > timestamp: Sun 2009-11-29 23:59:29 -0800 > message: > Work on Expression in PyDOLFIN > - Simplify expression.py > * I have split the logic in create_expression_class into on for > python Expressions and one for compiled ones > * Removed the possibilty of creating compiled classes that is > not instantiated. (never used?) > - Fixes in some typemaps > * director typemaps for std::vector<double> &values now works > * I have turned off director for eval(values,x). Could not get it > to work > It now works through eval(values,data) in python > * Got argout typemap for std::vector<double> & values to work. > > Bottom line we are up to par in features wrt Expression in PyDOLFIN > - Plotting works > - Subclassing in Python works for both eval and eval_data > - CompiledExpressions work: > * MultiLine expressions need to set its value shape if it is > different from 0 > > TODO: > - Cleanup in the SWIG interface files and docstrings in expression.py > - get_row, and set_row in GenericMatrix does not work as intended. > Needed to sacrify this one... > modified: > demo/pde/poisson/python/demo.py > dolfin/function/Expression.cpp > dolfin/function/Expression.h > dolfin/swig/function_post.i > dolfin/swig/function_pre.i > dolfin/swig/std_vector_typemaps.i > site-packages/dolfin/expression.py > > > === modified file 'demo/pde/poisson/python/demo.py' > --- demo/pde/poisson/python/demo.py 2009-11-29 21:17:18 +0000 > +++ demo/pde/poisson/python/demo.py 2009-11-30 07:59:29 +0000 > @@ -13,7 +13,7 @@ > """ > > __author__ = "Anders Logg ([email protected])" > -__date__ = "2007-08-16 -- 2009-11-24" > +__date__ = "2007-08-16 -- 2009-11-29" > __copyright__ = "Copyright (C) 2007-2009 Anders Logg" > __license__ = "GNU LGPL Version 2.1" > > @@ -25,7 +25,7 @@ > > class Source(Expression): > def eval(self, values, x): > - values[0] = > 4.0*DOLFIN_PI*DOLFIN_PI*DOLFIN_PI*DOLFIN_PI*sin(DOLFIN_PI*x[0])*sin(DOLFIN_PI*x[1]) > + values[0] = 10*exp(-((x[0] - 0.5)**2 + (x[1] - 0.5)**2) / 0.02) > > # Define Dirichlet boundary (x = 0 or x = 1) > def boundary(x): > @@ -38,8 +38,8 @@ > # Define variational problem > v = TestFunction(V) > u = TrialFunction(V) > -#f = Expression("10*exp(-(pow(x[0] - 0.5, 2) + pow(x[1] - 0.5, 2)) / 0.02)") > -f = Source(V) > +f = Expression("10*exp(-(pow(x[0] - 0.5, 2) + pow(x[1] - 0.5, 2)) / 0.02)") > +#f = Source() > g = Expression("sin(5*x[0])") > a = inner(grad(v), grad(u))*dx > L = v*f*dx - v*g*ds > @@ -53,4 +53,4 @@ > file << u > > # Plot solution > -#plot(u, interactive=True) > +plot(u, interactive=True) > > === modified file 'dolfin/function/Expression.cpp' > --- dolfin/function/Expression.cpp 2009-11-29 12:18:27 +0000 > +++ dolfin/function/Expression.cpp 2009-11-30 07:59:29 +0000 > @@ -2,7 +2,7 @@ > // Licensed under the GNU LGPL Version 2.1. > // > // First added: 2009-09-28 > -// Last changed: 2009-10-11 > +// Last changed: 2009-11-29 > // > // Modified by Johan Hake, 2009. > > @@ -114,8 +114,8 @@ > > //----------------------------------------------------------------------------- > void Expression::eval(std::vector<double>& values, const > std::vector<double>& x) const > { > - eval(&values[0], x); > - //error("Missing eval() for Expression (must be overloaded)."); > + //eval(&values[0], x); > + error("Missing eval() for Expression (must be overloaded)."); > } > > //----------------------------------------------------------------------------- > > > === modified file 'dolfin/function/Expression.h' > --- dolfin/function/Expression.h 2009-11-29 21:17:18 +0000 > +++ dolfin/function/Expression.h 2009-11-30 07:59:29 +0000 > @@ -2,7 +2,7 @@ > // Licensed under the GNU LGPL Version 2.1. > // > // First added: 2009-09-28 > -// Last changed: 2009-10-11 > +// Last changed: 2009-11-29 > > #ifndef __EXPRESSION_H > #define __EXPRESSION_H > @@ -75,13 +75,6 @@ > /// Evaluate expression, must be overloaded by user (simple version) > virtual void eval(std::vector<double>& values, const > std::vector<double>& x) const; > > - // Tempory fix while figuring out SWIG > - virtual void eval(double* values, const std::vector<double>& x) const > - { > - cout << "In eval " << endl; > - error("Missing eval() for Expression (must be overloaded)."); > - } > - > protected: > > // Value shape > > === modified file 'dolfin/swig/function_post.i' > --- dolfin/swig/function_post.i 2009-11-16 15:13:08 +0000 > +++ dolfin/swig/function_post.i 2009-11-30 07:59:29 +0000 > @@ -25,16 +25,16 @@ > > //----------------------------------------------------------------------------- > // Extend the Data class with an accessor function for the x coordinates > > //----------------------------------------------------------------------------- > -//%extend dolfin::Data { > -// PyObject* x_() { > -// npy_intp adims[1]; > -// adims[0] = self->cell().dim(); > -// PyArrayObject* array = > reinterpret_cast<PyArrayObject*>(PyArray_SimpleNewFromData(1, adims, > NPY_DOUBLE, (char *)(self->x))); > -// if ( array == NULL ) return NULL; > -// PyArray_INCREF(array); > -// return reinterpret_cast<PyObject*>(array); > -// } > -//} > +%extend dolfin::Data { > + PyObject* x_() { > + npy_intp adims[1]; > + adims[0] = self->cell().dim(); > + PyArrayObject* array = > reinterpret_cast<PyArrayObject*>(PyArray_SimpleNewFromData(1, adims, > NPY_DOUBLE, reinterpret_cast<char *>( &(const_cast<std::vector<double>& > >(self->x))[0] ))); > + if ( array == NULL ) return NULL; > + PyArray_INCREF(array); > + return reinterpret_cast<PyObject*>(array); > + } > +} > > > //----------------------------------------------------------------------------- > // Clear director typemaps > @@ -42,3 +42,4 @@ > %clear const double* x; > %clear double* y; > %clear const std::vector<double>& x; > +%clear std::vector<double>& values; > > === modified file 'dolfin/swig/function_pre.i' > --- dolfin/swig/function_pre.i 2009-11-29 21:17:18 +0000 > +++ dolfin/swig/function_pre.i 2009-11-30 07:59:29 +0000 > @@ -9,7 +9,7 @@ > // Modified by Kent-Andre Mardal, 2009 > // > // First added: 2007-08-16 > -// Last changed: 2009-10-07 > +// Last changed: 2009-11-29 > > // > =========================================================================== > // SWIG directives for the DOLFIN function kernel module (pre) > @@ -60,8 +60,7 @@ > > //----------------------------------------------------------------------------- > %ignore dolfin::Data::x; > %rename (x) dolfin::Data::x_(); > -%ignore dolfin::Expression::eval(std::vector<double>& values, const > std::vector<double>& x) const; > -//%ignore dolfin::Expression::eval(double* values, const > std::vector<double>& x) const; > +//%ignore dolfin::Expression::eval(std::vector<double>& values, const > std::vector<double>& x) const; > > > //----------------------------------------------------------------------------- > // Modifying the interface of Constant > @@ -118,19 +117,18 @@ > %feature("novaluewrapper") std::vector<double>; > > > //----------------------------------------------------------------------------- > -// Instantiate a dummy std::vector<dolfin::uint> so value wrapper is not used > +// Instantiate a dummy std::vector<dolfin::double> so value wrapper is not > used > > //----------------------------------------------------------------------------- > %template () std::vector<double>; > > > //----------------------------------------------------------------------------- > -// Typemap for std::vector<dolfin::uint> values > +// Typemap for std::vector<dolfin::double> values (used in Constant > constructor) > > //----------------------------------------------------------------------------- > -//%typecheck(SWIG_TYPECHECK_DOUBLE_ARRAY) std::vector<double> values > -//{ > -// $1 = PyList_Check($input) ? 1 : 0; > -//} > +%typecheck(SWIG_TYPECHECK_DOUBLE_ARRAY) std::vector<double> values > +{ > + $1 = PyList_Check($input) ? 1 : 0; > +} > > -/* > %typemap (in) std::vector<double> values > { > if (PyList_Check($input)) > @@ -153,7 +151,7 @@ > SWIG_exception(SWIG_TypeError, "expected list of floats"); > } > } > -*/ > + > > //----------------------------------------------------------------------------- > // Add director classes > > //----------------------------------------------------------------------------- > @@ -163,28 +161,14 @@ > %feature("nodirector") dolfin::Expression::gather; > %feature("nodirector") dolfin::Expression::value_dimension; > %feature("nodirector") dolfin::Expression::value_rank; > +%feature("nodirector") dolfin::Expression::eval(std::vector<double>& values, > const std::vector<double>& x) const; > > > //----------------------------------------------------------------------------- > // Director typemap for values in Expression > > //----------------------------------------------------------------------------- > -%typemap(directorin) double* values > -{ > - { > - // Compute size of value (number of entries in tensor value) > - dolfin::uint size = 1; > - for (dolfin::uint i = 0; i < this->value_rank(); i++) > - size *= this->value_dimension(i); > - > - npy_intp dims[1] = {size}; > - $input = PyArray_SimpleNewFromData(1, dims, NPY_DOUBLE, > - reinterpret_cast<char*>($1_name)); > - } > -} > - > %typemap(directorin) std::vector<double>& values > { > { > - std::cout << "In typemap " << std::endl; > // Compute size of x > npy_intp dims[1] = {$1_name.size()}; > $input = PyArray_SimpleNewFromData(1, dims, NPY_DOUBLE, > @@ -195,18 +179,11 @@ > > //----------------------------------------------------------------------------- > // Director typemap for coordinates in Expression > > //----------------------------------------------------------------------------- > -//%typemap(directorin) const double* x { > -// { > -// // Compute size of x > -// npy_intp dims[1] = {this->geometric_dimension()}; > -// $input = PyArray_SimpleNewFromData(1, dims, NPY_DOUBLE, > reinterpret_cast<char *>(const_cast<double*>($1_name))); > -// } > -//} > - > // FIXME: Is there a better way to map a std::vector to a numpy array? > %typemap(directorin) const std::vector<double>& x > { > { > + std::cout << "In typemap & x" << std::endl; > // Compute size of x > npy_intp dims[1] = {$1_name.size()}; > $input = PyArray_SimpleNewFromData(1, dims, NPY_DOUBLE, > @@ -214,3 +191,43 @@ > } > } > > +////----------------------------------------------------------------------------- > +//// In typemaps for std::vector<double> _array > +////----------------------------------------------------------------------------- > +//%typemap(in) std::vector<double> & _array (std::vector<double> vec_tmp) > +//{ > +// // Check arguments > +// if (!PyArray_Check($input)){ > +// PyErr_SetString(PyExc_TypeError, "numpy array of 'double' expected. > Make sure that the numpy array use dtype='d'."); > +// return NULL; > +// } > +// > +// PyArrayObject *xa = reinterpret_cast<PyArrayObject*>(input); > +// if (!PyArray_TYPE(xa) == NPY_DOUBLE ){ > +// PyErr_SetString(PyExc_TypeError, "numpy array of 'double' expected. > Make sure that the numpy array use dtype='d'."); > +// return NULL; > +// } > +// > +// // Get size and reserve the tmp vector > +// npy_int size = PyArray_Size(xa); > +// vec_tmp.reserve(size); > +// > +// // Get the data > +// double * data = static_cast<double*>(PyArray_DATA(xa)); > +// for (int i=0, i<size; i++) > +// vec_tmp[i] = data[i]; > +// > +// // Provide the out argument > +// $1 = &vec_tmp; > +//} > +// > +//%typecheck(SWIG_TYPECHECK_DOUBLE_ARRAY) std::vector<double> & values > +//{ > +// $1 = PyArray_Check($input) ? 1 : 0; > +//} > +// > +//%typecheck(SWIG_TYPECHECK_DOUBLE_ARRAY) const std::vector<double> & x > +//{ > +// $1 = PyArray_Check($input) ? 1 : 0; > +//} > +// > > === modified file 'dolfin/swig/std_vector_typemaps.i' > --- dolfin/swig/std_vector_typemaps.i 2009-11-29 21:17:18 +0000 > +++ dolfin/swig/std_vector_typemaps.i 2009-11-30 07:59:29 +0000 > @@ -3,7 +3,7 @@ > // Licensed under the GNU LGPL Version 2.1. > // > // First added: 2009-08-31 > -// Last changed: 2009-09-29 > +// Last changed: 2009-11-29 > > > //============================================================================= > // In this file we declare what types that should be able to be passed using > a > @@ -132,6 +132,7 @@ > > //----------------------------------------------------------------------------- > %typemap (in,numinputs=0) std::vector<TYPE>& ARG_NAME (std::vector<TYPE> > vec_temp) > { > + // ARGOUT ARG_NAME > $1 = &vec_temp; > } > > @@ -184,8 +185,11 @@ > IN_TYPEMAPS_STD_VECTOR_OF_POINTERS(FunctionSpace) > > ARGOUT_TYPEMAP_STD_VECTOR_OF_PRIMITIVES(dolfin::uint, INT32, cells, NPY_INT) > -ARGOUT_TYPEMAP_STD_VECTOR_OF_PRIMITIVES(dolfin::uint, INT32, columns, > NPY_INT) > -ARGOUT_TYPEMAP_STD_VECTOR_OF_PRIMITIVES(double, DOUBLE, values, NPY_DOUBLE) > +// FIXME: We cannot have argout typemaps for columns and values > +// FIXME: They work for get_row, but they interfere with eval interface. > +// FIXME: They should also _not_ kick in for const std::vector<TYPE>, but > they do > +//ARGOUT_TYPEMAP_STD_VECTOR_OF_PRIMITIVES(dolfin::uint, INT32, columns, > NPY_INT) > +//ARGOUT_TYPEMAP_STD_VECTOR_OF_PRIMITIVES(double, DOUBLE, values, NPY_DOUBLE) > > > //----------------------------------------------------------------------------- > // NumPy to const std::vector<double>& typemap. > @@ -194,28 +198,6 @@ > $1 = PyArray_Check($input) ? 1 : 0; > } > > -%typemap(in) std::vector<double>& x (std::vector<double> temp) > -{ > - { > - if (PyArray_Check($input)) > - { > - PyArrayObject *xa = reinterpret_cast<PyArrayObject*>($input); > - if ( PyArray_TYPE(xa) == NPY_DOUBLE ) > - { > - const unsigned int size = PyArray_DIM(xa, 0); > - temp.resize(size); > - double* array = static_cast<double*>(PyArray_DATA(xa)); > - std::copy(array, array + size, temp.begin()); > - $1 = &temp; > - } > - else > - SWIG_exception(SWIG_TypeError, "NumPy array expected"); > - } > - else > - SWIG_exception(SWIG_TypeError, "NumPy array expected"); > - } > -} > - > %typemap(in) const std::vector<double>& x (std::vector<double> temp) > { > { > @@ -237,19 +219,64 @@ > SWIG_exception(SWIG_TypeError, "NumPy array expected"); > } > } > -//----------------------------------------------------------------------------- > -// const std::vector<double>& to NumPy typemap. > -//----------------------------------------------------------------------------- > -%typemap(argout) const std::vector<double>& x > -{ > - // Do nothing > -} > -//----------------------------------------------------------------------------- > -// std::vector<double>& to NumPy typemap. > -//----------------------------------------------------------------------------- > -%typemap(argout) std::vector<double>& x > -{ > - SWIG_exception(SWIG_TypeError, "std::vector<double> to NumPy (argout) not > implemented"); > -} > -//----------------------------------------------------------------------------- > - > + > +//----------------------------------------------------------------------------- > +// std::vector<double>& to NumPy typemap. > +//----------------------------------------------------------------------------- > +%typemap(in) std::vector<double>& values (std::vector<double> temp, > PyArrayObject *out_array = 0) > +{ > + if (!PyArray_Check($input)) > + SWIG_exception(SWIG_TypeError, "numpy array of 'double' expected. Make > sure that the numpy array use dtype='d'."); > + out_array = reinterpret_cast<PyArrayObject*>($input); > + if ( !PyArray_TYPE(out_array) == NPY_DOUBLE ) > + SWIG_exception(SWIG_TypeError, "numpy array of 'double' expected. Make > sure that the numpy array use dtype='d'."); > + > + // Use the size of the incomming array to pass a correct sized vector to > the function > + const unsigned int size = PyArray_DIM(out_array, 0); > + temp.resize(size); > + $1 = &temp; > +} > + > +%typemap(argout) std::vector<double>& values > +{ > + // Get the pointer to the data in the passed NumPy array > + double* array = static_cast<double*>(PyArray_DATA(out_array$argnum)); > + > + // Copy the content from the temp array to the NumPy array > + std::copy(temp$argnum.begin(), temp$argnum.end(), array); > + > + // Tell the function to return 'None', which means not return value > + $result = Py_None; > + > +} > + > +//----------------------------------------------------------------------------- > +// std::vector<double>& to NumPy typemap. > +//----------------------------------------------------------------------------- > +%typemap(in) std::vector<double>& vertex_values (std::vector<double> temp, > PyArrayObject *out_array = 0, dolfin::uint init_size = 0) > +{ > + if (!PyArray_Check($input)) > + SWIG_exception(SWIG_TypeError, "numpy array of 'double' expected. Make > sure that the numpy array use dtype='d'."); > + out_array = reinterpret_cast<PyArrayObject*>($input); > + if ( !PyArray_TYPE(out_array) == NPY_DOUBLE ) > + SWIG_exception(SWIG_TypeError, "numpy array of 'double' expected. Make > sure that the numpy array use dtype='d'."); > + > + // Use the size of the incomming array to pass a correct sized vector to > the function > + const unsigned int size = PyArray_DIM(out_array, 0); > + init_size = size; > + temp.resize(size); > + $1 = &temp; > +} > + > +%typemap(argout) std::vector<double>& vertex_values > +{ > + // Get the pointer to the data in the passed NumPy array > + double* array = static_cast<double*>(PyArray_DATA(out_array$argnum)); > + > + // Copy the content from the temp array to the NumPy array > + std::copy(temp$argnum.begin(), temp$argnum.end(), array); > + > + // Tell the function to return 'None', which means not return value > + $result = Py_None; > + > +} > > === modified file 'site-packages/dolfin/expression.py' > --- site-packages/dolfin/expression.py 2009-11-28 21:18:15 +0000 > +++ site-packages/dolfin/expression.py 2009-11-30 07:59:29 +0000 > @@ -42,7 +42,6 @@ > __license__ = "GNU LGPL Version 2.1" > > # Modified by Anders Logg, 2008-2009. > -# Modified by Garth N. Wells, 2009. > > __all__ = ["Expression", "Expressions"] > > @@ -64,96 +63,93 @@ > from compile_expressions import compile_expressions > from functionspace import * > > -def create_expression_class(name, > - cpp_base, > - user_bases = None, > - user_dict = None, > - dim_needs_to_be_passed = False): > +def create_compiled_expression_class(cpp_base): > + # Check the cpp_base > + assert(isinstance(cpp_base, (types.ClassType, type))) > + > + def __init__(self, cpparg, defaults=None, element=None, degree=None): > + """ Initialize the Expression """ > + # Initialize the cpp base class first and extract value_shape > + cpp_base.__init__(self) > + value_shape = tuple(self.value_dimension(i) \ > + for i in range(self.value_rank())) > + > + # Store the dim > + if len(value_shape) == 0: > + self._dim = 0 > + elif len(value_shape) == 1: > + self._dim = value_shape[0] > + else: > + self._dim = value_shape > + > + # Select an appropriate element if not specified > + if element is None: > + element = _auto_select_element_from_shape(value_shape, degree) > + else: > + # Check that we have an element > + if not isinstance(element, ufl.FiniteElementBase): > + raise TypeError, "The 'element' argument must be a UFL > finite element." > + > + # Check same value shape of compiled expression and passed > element > + if not element.value_shape() == value_shape: > + raise ValueError, "The value shape of the passed 'element', > is not equal to the value shape of the compiled expression." > + > + # Initialize UFL base class > + self._ufl_element = element > + ufl.Function.__init__(self, self._ufl_element) > + > + # Create and return the class > + return type("CompiledExpression", (Expression, ufl.Function, cpp_base), > {"__init__":__init__}) > + > +def create_python_derived_expression_class(name, user_bases, user_dict): > """Return Expression class > > This function is used to create all the dynamically created Expression > classes. It takes a name, and a compiled cpp.Expression and returns > - a dolfin.Expression class. In addition to cpp.Expression and > dolfin.Expression > - it also inherits from ufl.Function. > + a class that inherits the compiled class together with dolfin.Expression > + and ufl.Function. > > @param name: > The name of the class > - @param cpp_base: > - The cpp.Expression base class which the created > - Expression class shall inherit. > @param user_bases: > - Optional user defined bases > + User defined bases > @param user_dict: > - Optional dict with user specified function or attributes > - @param dim_needs_to_be_passed: > - Optional if a simple expression using cpparg is created with no > - information about geometric dimensions > + Dict with user specified function or attributes > """ > > - # Check the name > + # Check args > assert(isinstance(name, str)) > - assert(name != "Expression"), "Cannot create a sub class of Expression > with the same name as Expression" > - > - assert(isinstance(cpp_base, (types.ClassType, type))) > + assert(isinstance(user_bases, list)) > + assert(isinstance(user_dict, dict)) > > # Define the bases > - user_bases = user_bases or [] > assert(all([isinstance(t, (types.ClassType, type)) for t in user_bases])) > - bases = tuple([Expression, ufl.Function, cpp_base] + user_bases) > - > - # Define the dictionary of the class > - dict_ = user_dict or {} > - assert isinstance(dict_, dict) > - > + bases = tuple([Expression, ufl.Function, cpp.Expression] + user_bases) > + > # If a user init is not provided create a dummy one > - if "__init__" not in dict_: > + if "__init__" not in user_dict: > def user_init(self, *arg, **kwargs): pass > else: > - user_init = dict_.pop("__init__") > - > - # If a user init is not provided create a dummy one > - if "dim" not in dict_: > - def user_dim(self): return 1 > - else: > - user_dim = dict_.pop("dim") > + user_init = user_dict.pop("__init__") > > def __init__(self, *args, **kwargs): > - # This is called if no user defined init function is provided > > # Get element and degree > - element = kwargs.pop("element", None) > - degree = kwargs.pop("degree", None) > + element = kwargs.get("element", None) > + degree = kwargs.get("degree", None) > > # Select an appropriate element if not specified > - if element is None and len(args) > 0: > - element = _auto_select_element_from_cpparg(args[0], degree) > - elif element is None: > - element = _auto_select_element_from_dim(user_dim(self), degree) > - > - # Check that we have have an element > - if not isinstance(element, ufl.FiniteElementBase): > + if element is None: > + element = _auto_select_element_from_dim(self.dim(), degree) > + elif not isinstance(element, ufl.FiniteElementBase): > raise TypeError, "The 'element' argument must be a UFL finite > element." > - > + > # Initialize UFL base class > self._ufl_element = element > ufl.Function.__init__(self, self._ufl_element) > > # Initialize cpp_base class > - > - # First check if we are instantiating a user-defined Python class > - if "eval" in dict_ or "eval_data" in dict_: > - assert cpp_base == cpp.Expression > - cpp_base.__init__(self, list(self._ufl_element.value_shape())) > - else: > - cpp_base.__init__(self) > - > - # Check that the value_shape of the ufl.FiniteElement corresponds > with the > - # created cpp.Expression > - shape = self._ufl_element.value_shape() > - if not (self.value_rank() == len(shape) and > - all(dim == self.value_dimension(i) for i, dim in > enumerate(shape))): > - exp_shape = tuple(self.value_dimension(i) for i in > xrange(self.value_rank())) > - raise ValueError, "value_shape of passed element does not match > value_shape of the Expression: %s != %s"%(str(shape), str(exp_shape)) > + cpp.Expression.__init__(self, list(self._ufl_element.value_shape())) > > # Calling the user defined_init > user_init(self, *args, **kwargs) > @@ -165,6 +161,7 @@ > __init__.__doc__ = """ Initialize the Expression""" > > # NOTE: Do not prevent the user to overload attributes "reserved" by > PyDOLFIN > + # NOTE: Why not? > > ## Collect reserved attributes from both cpp.Function and ufl.Function > #reserved_attr = dir(ufl.Function) > @@ -180,22 +177,22 @@ > # if attr in dict_: > # raise TypeError, "The Function attribute '%s' is reserved by > PyDOLFIN."%attr > > - # Fill the dict_ with constructed function > - dict_["__init__"] = __init__ > + # Add __init__ to the user_dict > + user_dict["__init__"] = __init__ > > # Create the class and return it > - return type(name, bases, dict_) > + return type(name, bases, user_dict) > > class ExpressionMetaClass(type): > - > + """ Meta Class for Expression""" > def __new__(cls, name, bases, dict_): > - """ Return a new Expression class """ > + """ Returns a new Expression class """ > > assert(isinstance(name, str)), "Expecting a 'str'" > assert(isinstance(bases, tuple)), "Expecting a 'tuple'" > assert(isinstance(dict_, dict)), "Expecting a 'dict'" > > - # First check if we are creating the Function class > + # First check if we are creating the Expression class > if name == "Expression": > # Assert that the class is _not_ a subclass of Expression, > # i.e., a user have tried to: > @@ -209,7 +206,9 @@ > # this module > return type.__new__(cls, name, bases, dict_) > > - # If subclassing Expression directly (used in specialfunctions.py) > + # If creating a fullfledged derived expression class, i.e, inheriting > + # dolfin.Expression, ufl.Function and cpp.Expression (or a subclass) > + # then just return the new class. > if len(bases) >= 3 and bases[0] == Expression and \ > bases[1] == ufl.Function and issubclass(bases[2], > cpp.Expression): > # Return the instantiated class > @@ -221,59 +220,28 @@ > # remove Expression, to be added later > user_bases.remove(Expression) > > - # Check the cppcode and eval attributes > - if 'cpparg' in dict_ and ('eval' in dict_ or 'eval_data' in dict_) : > - raise TypeError, "Cannot create class with both > 'eval'/'eval_data' and 'cpparg' attributes defined." > - > - # If the Expression class is a user defined python class, case 4. > from docstring > - if 'eval' in dict_ or 'eval_data' in dict_: > - # Get name of eval function > - eval_name = 'eval' if 'eval' in dict_ else 'eval_data' > - > - user_eval = dict_[eval_name] > - > - # Check type and number of arguments of user_eval function > - if not isinstance(user_eval, types.FunctionType): > - raise TypeError, "'%s' attribute must be a > 'function'"%eval_name > - if not user_eval.func_code.co_argcount == 3: > - raise TypeError, "The overloaded '%s' function must use > three arguments"%eval_name > - > - return create_expression_class(name, cpp.Expression, user_bases, > dict_) > - > - # If cpparg is provided, case 5. from docstring > - if 'cpparg' in dict_: > - > - # Check the handed attributes and return an args tuple > - cpparg = dict_.pop('cpparg') > - defaults = dict_.pop("defaults",None) > - > - # Check that the user has not provide any other attributes > - # than the allowed ones. > - if len(dict_) > 1: > - dict_.pop('__module__') > - raise TypeError, "Not allowed to provide user defined > attributes to a sub class of Expression when the compiled function interface > is used. Found: %s"%\ > - (", ".join(["'%s'"%key for key in dict_.iterkeys()])) > - > - # Check arguments > - _check_cpparg(cpparg) > - _check_defaults(defaults) > - > - # Compile the cppargs > - cpp_base = compile_expressions([cpparg], [defaults])[0] > - > - # Add back the cpparg as an attribute > - cpp_base.cpparg = cpparg > - > - # If defaults where handed add it back too > - if defaults is not None: > - cpp_base.defaults = defaults > - > - # Create the Expression class and return it > - return create_expression_class(name, cpp_base, user_bases, > - dim_needs_to_be_passed = not > _is_complex_expression(cpparg)) > - > - # If we have reached this stage raise error > - raise TypeError, "Error in subclassing Expression. For correct usage > see 4. and 5. in Expression docstring." > + # Check the user has provided either an eval or eval_data method > + if not ('eval' in dict_ or 'eval_data' in dict_): > + raise TypeError, "expected an overload 'eval' or 'eval_data' > method" > + > + # Get name of eval function > + eval_name = 'eval' if 'eval' in dict_ else 'eval_data' > + > + user_eval = dict_[eval_name] > + > + # Check type and number of arguments of user_eval function > + if not isinstance(user_eval, types.FunctionType): > + raise TypeError, "'%s' attribute must be a 'function'"%eval_name > + if not user_eval.func_code.co_argcount == 3: > + raise TypeError, "The overloaded '%s' function must use three > arguments"%eval_name > + # A nice hack to get around some SWIG director problems > + # In short: eval_data works but not eval... > + if eval_name == 'eval': > + def eval_data(self, values, data): > + user_eval(self, values, data.x()) > + dict_['eval_data'] = eval_data > + > + return create_python_derived_expression_class(name, user_bases, > dict_) > > #--- The user interface --- > > @@ -324,7 +292,6 @@ > return ufl.Function.__call__(self,*args) > > # Some help variables > - #dim = self.geometric_dimension() > value_size = ufl.common.product(self.ufl_element().value_shape()) > > # If values (return argument) is passed, check the type and length > @@ -349,19 +316,14 @@ > raise TypeError, "expected a scalar or an iterable as coordinate > argument" > # Check for scalar x > if isinstance(x[0], (int, float)): > - #if not dim == 1: > - # raise TypeError, "expected a coordinate argument of length > %d"%dim > x = numpy.fromiter(x, 'd') > else: > x = x[0] > - #if len(x) != dim: > - # raise TypeError, "expected an iterable of length %d as > coordinate argument"%dim > if isinstance(x, (list, tuple)): > x = numpy.fromiter(x, 'd') > > # If several x arguments have been provided > else: > - #if len(x) != dim or not all(isinstance(v,(int,float)) for v in x): > if not all(isinstance(v,(int,float)) for v in x): > raise TypeError, "expected different number of scalar arguments > for the coordinates" > x = numpy.fromiter(x,'d') > @@ -550,8 +512,8 @@ > Optional quadrature degree element. > """ > > - # If the __new__ function is called because we are instantiating a > sub > - # class of Expression, then use the object's __new__ function instead > + # If the __new__ function is called because we are instantiating a > python sub > + # class of Expression, then just return a new instant of the passed > class > if cls.__name__ != "Expression": > return object.__new__(cls) > > @@ -562,85 +524,111 @@ > # Compile module and get the cpp.Expression class > cpp_base = compile_expressions([cpparg], [defaults])[0] > > - # Store arguments for later use > + # Store compile arguments for later use > cpp_base.cpparg = cpparg > cpp_base.defaults = defaults > - cpp_base.element = element > - cpp_base.degree = degree > - return object.__new__(create_expression_class("CompiledExpression", > cpp_base)) > - > + > + return object.__new__(create_compiled_expression_class(cpp_base)) > + > + # This method is only included so a user can check what arguments > + # one should use in IPython using tab completion > + def __init__(self, cpparg=None, defaults=None, element=None, > degree=None):pass > + > + # Reuse the docstring from __new__ > + __init__.__doc__ = __new__.__doc__ > + > def ufl_element(self): > " Return the ufl FiniteElement." > return self._ufl_element > > def __str__(self): > + "x.__str__() <==> print(x)" > # FIXME: We might change this using rank and dimension instead > return "<Expression on a %s>" % str(self._ufl_element) > > def __repr__(self): > + "x.__repr__() <==> repr(x)" > return ufl.Function.__repr__(self) > > + # Default value of dim > + _dim = 0 > + > + def dim(self): > + """ Returns the dimension of the value""" > + return self._dim > + > + __call__ = expression__call__ > + > def Expressions(*args, **kwargs): > """ Batch-processed user-defined JIT-compiled expressions > ------------------------------------------------------- > > By specifying several cppargs one may compile more than one expression > - at a time. These may either be instantiated using a single FunctionSpace > - common for all expressions, using the optional kwarg 'V', or with > - a separate FunctionSpace for each cpparg: > + at a time: > > >>> f0, f1 = Expressions('sin(x[0]) + cos(x[1])', 'exp(x[1])', degree=3) > - >>> f0, f1 = Expressions('sin(x[0]) + cos(x[1])', 'exp(x[1])', > element=element) > - > - Here cppcode is a code snippet, which should be a string of C++ > + > + >>> f0, f1, f2 = Expressions((('A*sin(x[0])', 'B*cos(x[1])') > + ('0','1')), {'A':2.0,'B':3.0}, > + code, > + (('cos(x[0])','sin(x[1])'), > + ('sin(x[0])','cos(x[1])')), element=element) > + > + Here code is a C++ code snippet, which should be a string of C++ > code that implements a class that inherits from dolfin::Expression, > see user case 3. in Expression docstring > > Batch-processing of JIT-compiled expressions may significantly speed up > JIT-compilation at run-time. > - > -""" > + """ > > # Get the element/degree from kwarg > + if len(kwargs) > 1: > + raise TypeError, "Can only define one kwarg and that can only be > 'degree' or 'element'." > degree = kwargs.pop("degree", None) > element = kwargs.pop("element", None) > - if len(kwargs) > 1: > - raise TypeError, "Can only define one kwarg and that can only be > 'degree' or 'element'." > > # Iterate over the *args and collect input to compile_expressions > cppargs = []; defaults = []; i = 0; > - while i < len(args): > + while i < nargs: > + > + # Check type of cppargs > if not isinstance(args[i],(tuple, list, str)): > raise TypeError, "Expected either a 'list', 'tuple' or 'str' for > argument %d"%i > + > cppargs.append(args[i]) > - defaults.append(None) > i += 1 > > + # If we have more args and the next is a dict > + if i < nargs and isinstance(args[i], dict): > + # Append the dict to defaults > + _check_defaults(args[i]) > + defaults.append(args[i]) > + i += 1 > + else: > + # If not append None > + defaults.append(None) > + > # Compile the cpp.Expressions > cpp_bases = compile_expressions(cppargs, defaults) > > # Instantiate the return arguments > return_expressions = [] > - for i, cpp_base in enumerate(cpp_bases): > - > return_expressions.append(create_expression_class("CompiledExpression", > cpp_base) > - (cppargs[i], element=element, > degree=degree)) > + > + for cpp_base in cpp_bases: > + # If we only want the cpp.Expression > + > return_expressions.append(create_compiled_expression_class(cpp_base)(\ > + degree=degree, > + element=element)) > > # Return the instantiated Expressions > return tuple(return_expressions) > > -# Assign doc string > -expression__call__.__doc__ > - > -# Assign __call__ method > -Expression.__call__ = types.MethodType(expression__call__, None, Expression) > - > #--- Utility functions --- > > def _check_cpparg(cpparg): > "Check that cpparg makes sense" > > - if cpparg is None: return > - > # Check that we get a string expression or nested expression > if not isinstance(cpparg, (str, tuple, list)): > raise TypeError, "Please provide a 'str', 'tuple' or 'list' for the > 'cpparg' argument." > @@ -665,14 +653,6 @@ > "Check if cpparg is a complex expression" > return isinstance(cpparg, str) and "class" in cpparg and "Expression" in > cpparg > > -def _auto_select_element_from_cpparg(cpparg, degree=None): > - "Automatically select an appropriate element from cpparg." > - > - # Use numpy to get the shape > - shape = numpy.shape(cpparg) > - > - return _auto_select_element_from_shape(shape, degree) > - > def _auto_select_element_from_dim(dim, degree=None): > "Automatically select an appropriate element from dim." > > @@ -708,3 +688,12 @@ > cpp.info("Automatic selection of expression element: " + str(element)) > > return element > + > +def _check_name_and_base(name, cpp_base): > + # Check the name > + assert(isinstance(name, str)) > + assert(name != "Expression"), "Cannot create a sub class of Expression > with the same name as Expression" > + > + assert(isinstance(cpp_base, (types.ClassType, type))) > + > + > > > _______________________________________________ > Mailing list: https://launchpad.net/~dolfin > Post to : [email protected] > Unsubscribe : https://launchpad.net/~dolfin > More help : https://help.launchpad.net/ListHelp
signature.asc
Description: Digital signature
_______________________________________________ Mailing list: https://launchpad.net/~dolfin Post to : [email protected] Unsubscribe : https://launchpad.net/~dolfin More help : https://help.launchpad.net/ListHelp

