Implementation of Dirichlet boundary condition when tgv=-1

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.

1 Like