← Back to team overview

dolfin team mailing list archive

Re: Value rank and geometrical dimension for Expressions

 

On Monday 05 October 2009 11:20:00 Anders Logg wrote:
> On Mon, Oct 05, 2009 at 07:05:40AM +0100, Garth N. Wells wrote:
> > Johan Hake wrote:
> > > On Monday 05 October 2009 07:17:48 Anders Logg wrote:
> > >> On Mon, Oct 05, 2009 at 12:13:30AM +0200, Johan Hake wrote:
> > >>> Hello!
> > >>>
> > >>> I am trying to work my way through the Python stuff that needs to be
> > >>> updated. I have now encountered a problem with the director typemaps
> > >>> for double* values and double* x.
> > >>>
> > >>> We need the value rank and geometrical dimension to be able to wrap
> > >>> the value to Python. Previously we just used the provided
> > >>> FunctionSpace. Now with this gone we need to come up with something
> > >>> else.
> > >>>
> > >>> I know Anders commented that we could add these as optional virtual
> > >>> functions. That could work, but that would cause another two
> > >>> callbacks to Python each time the expression is evaluated during
> > >>> assemble, and what even worse:
> > >>>
> > >>>    larger code when subclassing ;)
> > >>
> > >> Yes! We don't want that...
> > >>
> > >>> What if we defined an optional constructor taking the value rank and
> > >>> geometrical dimensions as arguments and added two access-methods?
> > >>>
> > >>>     /// Constructor
> > >>>     Expression():_value_rank(-1),_dim(-1);
> > >>>
> > >>>     /// Constructor
> > >>>     Expression(uint value_rank,
> > >>> uint::dim):_value_rank(value_rank),_dim(dim);
> > >>>
> > >>>     ...
> > >>>
> > >>>     uint value_rank()
> > >>>     {
> > >>>        if _value_rank == -1
> > >>>            error("value_rank is not initialized")
> > >>>        return _value_rank
> > >>>     }
> > >>>
> > >>>     uint dim()
> > >>>     {
> > >>>        if _dim == -1
> > >>>            error("dim is not initialized")
> > >>>        return _dim
> > >>>     }
> > >>>
> > >>> A user will not instantiate an Expression using the optional
> > >>> constructor, this will be done for him behind the scene. So user code
> > >>> would look like:
> > >>>
> > >>>   class MyExpression(Expression):
> > >>>       eval(value,x):
> > >>>          ...
> > >>>
> > >>>   f = MyExpression(v.ufl_element())
> > >>
> > >> Sounds good, but why do we need the two access functions? Shouldn't
> > >> the extra constructor be enough?
> > >
> > > In the typemap I can only access public variables and methods. I do not
> > > think that making these variables public is a good idea.
> > >
> > >> We could store the value rank and
> > >> dimension as std::vector<uint> (a shape thing, same as numpy) in the
> > >> Expression class:
> > >>
> > >> class Expression
> > >> {
> > >> public:
> > >>
> > >>   /// Create scalar expressioin
> > >>   Expression(); // scalar
> > >>
> > >>   /// Create vector-valued expression with given dimension
> > >>   Expression(uint dim);
> > >>
> > >>   /// Create tensor-valued expression with given shape
> > >>   Expression(const std::vector<uint>& shape);
> > >>
> > >>   std::uint rank() const;
> > >>   std::uint dim(uint i) const;
> > >>   const std::vector<uint>& shape() const;
> > >>
> > >> private:
> > >>
> > >>   std::vector<uint> _shape;
> > >>
> > >> };
> > >
> > > We can store the value_shape as std::vector, similare to ufl. But we
> > > need to store the (geometrical_)dimension too. Otherwise we cannot
> > > figure out what the size of 'x' should be. Note that the
> > > geometrical_dimension has nothing to do with the value rank.
> > >
> > > So I think we need one constructor for both dimension and value_rank.
> >
> > Agree.
> 
> Yes, I'll see what I can do.

Nice!

Johan

> --
> Anders
> 
> > >>> If this is too ugly for you I can probably just add such
> > >>> functionality into a PythonExpression class. The Expression class is
> > >>> not that large now, or would that interfere with some other logic in
> > >>> the new GenericFunction remake?
> > >>
> > >> This would be good to have and we could also check that things match
> > >> against the incoming element in the restrict() call (using assert).
> > >
> > > I assume you mean that this is good to have in the DOLFIN library?
> >
> > I assume this too, and agree that it would be good to have.
> >
> > Garth
> >
> > >> We could then also reimplement the checks that are currently commented
> > >> out in Assembler.cpp and SystemAssembler.cpp.
> > >
> > > Ok.
> > >
> > > Johan
> > >
> > >
> > > _______________________________________________
> > > DOLFIN-dev mailing list
> > > DOLFIN-dev@xxxxxxxxxx
> > > http://www.fenics.org/mailman/listinfo/dolfin-dev
> >
> > _______________________________________________
> > DOLFIN-dev mailing list
> > DOLFIN-dev@xxxxxxxxxx
> > http://www.fenics.org/mailman/listinfo/dolfin-dev
> 


References