Implementation of Dirichlet boundary condition when tgv=-1

As far as I read, when using tgv=-1, freefem doesn’t apply dirichlet BC in a symmetric way.

The symmetric way would look like this:

image
image

The asymmetric way only modifies the third row without setting the column values to zero which freefem is actually doing now(?)

This could be good to be implemented, and let user select e.g tgv=-1 is asymmetric and tgv=0 symmetric strategy…

Dear @arousta,

According to About the functions LinearCG a symmetric elimination is possible with tgv=-2, although I never tried myself, and it is not available for the rhs, as far as I remember. On the other hand, the implementation is not that hard; at Exact Dirichlet BC I have presented a naive implementation at the level of FreeFEM itself.

Regards,
/Fotios Kasolis

1 Like

Dear Fotios,

Thanks for your info. So I tried tgv=-2 and it does what you said! Interesting, but this is no where mentioned in the documentation and I would assume it a hidden utility? (But don’t know why hidden :slight_smile:)

And tgv=-2 doesn’t do anything on rhs as you said but can do a fix. I don’t like to write too many nested loops in freefem so I do this to implement the rhs needed for the symmetric matrix:

border C(t=0,2*pi){ x=sin(t); y=cos(t); label=1; }
mesh Th = buildmesh(C(-6));
fespace Vh(Th,P1);
macro BC on(1,u=1) //
varf a(u,v) = int2d(Th)( dx(u)*dx(v)+dy(u)*dy(v) ) +int2d(Th)(v) +BC;
Vh uA, uS; // A: Asymmetric, S:Symmetric

matrix MA = a(Vh,Vh,tgv=-1);
real[int] rhsA = a(0,Vh,tgv=-1.);
uA[] = MA^-1*rhsA;
plot(uA,wait=1);

varf dirichletBC(u,v) = BC;
real[int] bc = dirichletBC(0,Vh,tgv=-1);
matrix MS = a(Vh,Vh,tgv=-2.);
matrix MD = MS-MA;
real[int] rhsS = MD*bc;
rhsS += rhsA;
uS[] = MS^-1*rhsS;
plot(uS);

If you look at the two figures in the opening post, MD is just the matrix (a13, a23) and we multiply it by the dirichlet bc values to get the final rhs.

But this is actually a horrible way of doing this, as I am using two extra matrices (MD & MS) to achieve such a simple thing. It would be much much better if the tgv=-2 becomes a built-in feature.

Thanks for your feedback!

I have opened an issue to add more infos about tgv values in the documentation (https://github.com/FreeFem/FreeFem-doc/issues/21).

2 Likes