Hello FreeFem community,

I am looking for how to impose conditions in vectorial PDE problems.

More preciesely, I am working on a NS problem. Assuming my variables are [u,v,p], in a domain, say [0,2]x[0,1], I look for how to impose conditions like:

- u=A in [1,2]x[0,1],
- v=B in [0,1]x[0,1],
- p=C on the segment [0,0]-[2,1],

where A, B, C are constants (or functions; these conditions are just for the purpose of presenting the question and we let aside any related mathematical issue).

For scalar problems. say wrt u, I can easily impose my condition - by identifying the index of nodes involved in the condition for u and modifying accordingly the stiffness matrix M and the Rhs vector.

For vectorial problems I can proceed similarly, but I need to know the precise structure of M. For example, to set

- u(X(i))=A,
- v(X(j))=B,
- p(X(k))=C.

Namely, is the matrix organised as N blocks of 3, or 3 blocs of N? So should I do:

- M(3
*i,3*i) = M(3*j+1,3*j+1) = M(3*k+2,3*k+2) = tgv, - M(i,i) = M(N+j,N+j) = M(2
*N+k,2*N+k) = tgv, - or ?

Here, N=number of nodes (so M is a (3*N)x(3*N) matrix).

Thank you for your help.

ps: For convenience, here is a piece of code I use

// domain, mesh, FE spaces

//

mesh Th = square(20,10, [2*x, y]); // rectangle [0,2]x[0,1]

fespace Vh(Th,P1b);

fespace Qh(Th,P1);

fespace Zh(Thf,[P1b,P1b,P1]);

//

varf a([u,v,p],[uh,vh,ph])

= int2d(Th)(…) + on(1, u=Cu) + on(2,v=Cb)+…;

//

matrix A(Zh,Zh);

//////

A(i,i)=tgv; // ?

// … ?