Function assignment operator is very confusing and buggy

[Jan Blechta]

It is very confusing for users what Function::operator=, Function.assign do, especially in the context of subfunctions. The documentation is not very helpful.

Particular problems:

1. Function space is mutated. Arguably this should not be possible. Once Function object exists it should be defined on a space forever.

2. The implementation of the space reassignment feature is buggy though. See MWE2 below.

3. The behaviour of assignment with ufl.split()/Function.split() is very confusing. Basically subfunctions keep referencing the old vector while the function itself refers to a new vector after assignment.

4. New temporary vector is created even in the case when old one could be just modified efficiently in place.

5. Inefficent and confusing assignment of type u.assign(w.sub(0)). The space of u is replaced by space w.sub(0).function_space().collapse(). This is very inefficient because collapsed dofmap and a respective map are constructed on every such assignment, prone to bugs and confusing. The space might look the same (having same UFL element) but with different dofmap.

6. Python interface Function.assign() provides too many different features in a one method. Operand can be: function on same space, function on different space, subfunction, expression (just duplicates Function.interpolate()), function axpy (just duplicates Function::operator=(FunctionAXPY&)), UFL expressions with sum, product, division and component (emulates axpy).

The situation goes totally against Pythonic explicit is better than implicit. Function.assign() does something clever enough and "assigns".

Proposed first step to solution: remove Function.assign() completely and see how all needed features can be done more explicitly (and eventually efficiently) in demos/tests. As a second step, too magic Function::operator= (const Function&) and superfluous Function::operator= (const Expression&) can be removed.

The first step can be done within the pybind11 transition.

MWE2:

from dolfin import *

mesh = UnitSquareMesh(3, 3)
P1 = FunctionSpace(mesh, "P", 1)
P2 = FunctionSpace(mesh, "P", 2)

u1 = Function(P1)
u2 = Function(P2)

u1.assign(u2)

assert u1.function_space() == u2.function_space()
assert u1.ufl_element() == u2.ufl_element()  # Fails!

[Prof Garth Wells]

Agree. Another part of the fix is to introduce 'function views' for subfunctions.

A question is what should the syntax be for 'updating' a function, e.g. inside a time loop? I think settling this will guide changes.

[Jan Blechta]

A question is what should the syntax be for 'updating' a function, e.g. inside a time loop? I think settling this will guide changes.

I don't think that this is an issue of deciding syntax. This is rather matter of deciding which of the following features should be ditched/renamed (in pseudocode, = not meaning Python's =):

1. u = v (on same spaces)
2. u = v (on different spaces?!)
3. u[i] = v
4. u = v[i]

Option 1. is just u.vector()[:] = v.vector() which is nice because it's explicit. In the case that Function.assign() does this and only this (and otherwise raises), I wouldn't object having such a convenience method.

Options 2. should not be allowed.

Option 3. is not implemented. Code w.sub(0).assign(u1) does not change w at all although w.sub(0) is a view into component so one could it expect it does the thing. If options 2. is forbidden this one should raise automatically until/if it is implemented.

Option 4. should not be done unless really needed - it is assignment between two different spaces. Users should be discouraged to do such thing unless really in a need. The problem with this is that the space on the lhs is dropped and new one is every time computed as v.function_space().collapse(). It should rather assert that the space on the lhs is already matching and allow the assignment only in that case.

[Prof Garth Wells]

u.vector()[:] = v.vector() is Python syntax. We need C++ syntax too.

For subfunctions, I would like something along the lines of:

FunctionSpace V;
FunctionSpaceView V0_view = V.sub(0)  // 
FunctionSpace V0(V0_view); // Like 'collapse'

Function u(V); 
FunctionView u0_view = u.sub(0);  // View into u

Function u0(u0_view);  // 'collapses u0_view'

[Jan Blechta]

is Python syntax. We need C++ syntax too.

Syntax is not important at first step. Important is to delete flawed functionality (option 2. above) and fix inefficiency of option 4 (FunctionView looks like a good solution).

FunctionView subclassing GenericFunction is a good concept. Function.cpp is full of if is_view.

[Prof Garth Wells]

I think the syntax is important; we decide how it should work and then remove whatever is left over.

[Prof Garth Wells]

Related issue: #531

[Jan Blechta]

Ok, to be more constructive I want to see explicit

u.vector()[:] = v.vector()  # Python
*u.vector() = *v.vector();  // C++

in demos.

Eventually I would be ok with convenience method for that.

Function::operator=(Function& v)
{
  if (*_function_space != *v.function_space())
  {
    dolfin_error("Nonmatching space!");
  }

  *u.vector() = *v.vector();
}

Function::operator=(Expression&)
{
  error_deprecated("Use Function::interpolate")
}

Function::operator=(FunctionAXPY&)
{
  // Current impl is ok; delegates work to operator=(Function&);
  // will raise with non-matching spaces
}

i.e. Function::operator= raising with anything else than matching spaces assignment. Python Function.assign() can do only that, not any more magic. Summarizing Function.assign assigns only if spaces match.

Your subfunctions proposal looks good.

[Prof Garth Wells]

Sounds good. Will be a good change.

[Jan Blechta]

Have to figure out first whether this depends on first having FunctionSpaceView and FunctionView (#531).

[Tara Root]

Look who’s here! Our help essay company constantly upgrade this software to make certain it provides the most reliable results. Because of this, as a reputable academic assistance company, we guarantee each piece of writing is free of plagiarism against on-line based sources, all sample- or free-paper databases, and all work ever checked by our own software.