Anders Logg wrote: > Where does this show up? >
Not sure what's up. Demo seems to work now. 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(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 > _______________________________________________ Mailing list: https://launchpad.net/~dolfin Post to : [email protected] Unsubscribe : https://launchpad.net/~dolfin More help : https://help.launchpad.net/ListHelp
