dolfin team mailing list archive

["Andy Ray Terrel" <aterrel@xxxxxxxxxxxx>]

Okay so I am trying to get the Scott Vogelius Stokes algorithm working
again, but I've hit a bit of a wall with the copying 4th ordered
functions. In the algorithm there is a term that gets updated every
iteration (w = w + r * u).  So I could do this with a form as Anders
has suggested before, but it seems silly to use the 10K lines to do
what essentially is a vector axpy.  Since some people at simula are
working on the function classes, I pose the question: Should we be
able to do this with functions? For example something like:

...
1   Function w;
2   w = Function(mesh, 0.0);
3   BilinearForm a(-1e3);
4   for (...) {
5     LinearForm L(f,w);
6     LinearPDE pde(a,L,mesh,bcs);
7     pde.solve(u);
8     w += 1e3 * u;
...
9   }
...

Okay so assignment doesn't work for Constant functions (simple fix I
put in my local repo).  Okay operaters +, +=, * are not defined, so I
expand line 8 to:

w.vector().axpy(1e3,u.vector());

Okay then line 8 won't work because w was previously a constant
function.  I tried adding something like:

if (iter == 0) w = u;

but that didn't work because it is 4th order and the ElementLibrary
doesn't define it.

I have tried monkeying around with creating a better assignment
operator and copy constructor but then I run into other problems.  For
example who owns the finite_element and DOF, in my current code they
are getting destroyed with the form and that kills later preprocessing
I do.

I have some hacks that basically work, but I wanted to poll the list
to see what is being planned and put my two cents in about it.