[Getfem-users] Reduction problems

[Torquil Macdonald Sørensen]

Hi!

I've started with the tests/laplacian.cc program, modified it a bit, and 
implemented an example reduction scheme. But whenever I use reduction, the 
solution is identical to zero apart from on the boundary, where I use Dirichlet 
conditions. In addition, I've also tried an approach using the Laplacian model, 
with the same incorrect results.

In more detail, I have a square 2d domain and the mesh is:

getfem::regular_unit_mesh(mesh, ref, bgeot::parallelepiped_geotrans(2, 1));

On this mesh, I'm using FEM_QK(2,2) for both parts of the equation, and 
IM_GAUSS_PARALLELEPIPED(2,k) for k=3, but I've also tried k=4 and 5.

The reduction scheme consists of using only the DOFs at the corners of the 
squares, and not the DOFs at the intermediate points along the edges or the 
central point of the square. I thought I'd try this as a simple example of 
reduction. The extention matrix E is the transpose of the reduction matrix R, 
and I do test that R*E is the identity, so that part should be OK.

When using reduction, I make sure to use the "general Dirichlet" method in 
laplacian.cc, since the simplified approach won't work in that case.

The program runs fine without error messages with and without reduction enabled, 
but as I mentioned, the solution is incorrect when using reduction (zero 
everywhere apart from the boundary).

Any ideas about why this is so? Perhaps my reduction scheme is mathematically 
invalid for some reason, there is some bug in connection with reduction, or I've 
done something wrong?

Best regards,
Torquil Sørensen