On Monday 30 November 2009 10:03:40 Garth N. Wells wrote: > Anders Logg wrote: > > Where does this show up? > > Not sure what's up. Demo seems to work now.
I thought you where aiming at Expressions with an 's'. That was broken until it was fixed simultaneously by me and Garth. Johan > Garth > > > -- > > 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(DOLF > >>IN_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 > > _______________________________________________ > Mailing list: https://launchpad.net/~dolfin > Post to : [email protected] > Unsubscribe : https://launchpad.net/~dolfin > More help : https://help.launchpad.net/ListHelp > _______________________________________________ Mailing list: https://launchpad.net/~dolfin Post to : [email protected] Unsubscribe : https://launchpad.net/~dolfin More help : https://help.launchpad.net/ListHelp
