Problem related to SUPG

General Discussion
1

I am trying to solve the Navier Stokes equation in a driven cavity but i am acctually having problem with adding the stabilization term related to SUPG.
Freefam can only solve weak formulation problem and in the integral related to the SUPG terms we are adding the strond formulation of the problem:

fespace Uh(Th,P2);
Uh ux,uu;
Uh uy,vv;

fespace Ph(Th,P1);
Ph p,pp;

//stabilization term related to SUPG
int2d(Th)(deltatau(uxdx(uu)+uydy(vv))(uxdx(ux)…)

The problem is related to the term: ux*dx(ux)
Can someone help me or does someone have some example to look at?

2

Your stabilization term is non linear so need a schema to solve this problem.

3

Thanks a lot for the answer, i’m am new to freefem and i didn’t understand that non linear terms are a problem for freefem. Can i ask you one more thing?

After running a first simulation (without stabilization) with lower reynolds i tried to run a second one with higher Reynolds, so i can use the first one as initial condition for the second. I’m having problem with the SUPG stabilization terms. This term is non linear, but i’am trying some escamotage to solve the problem without linearizing. (driven cavity problem).

Do you have any advice or raccomandation (or see any errors) because the code does not converge.

Immagine WhatsApp 2024-03-09 ore 19.38.27_9f20ac4d
Immagine WhatsApp 2024-03-09 ore 19.38.27_9f20ac4d869×162 14 KB

fespace Uh(Th,P2)
Uh ux, uu;
Uh uy,vv;

fespace Ph(Th,P1)
Ph p,pp;

unastro=1; //reference velocity
cout<<“Reynolds =”<<Re<<endl;

//SUPG PARAMETERS
real delta=1; //parametro libero
Uh U=sqrt(ux^2+uy^2); //modulo velocità
Uh Reu=U* hTriangle/(2*nu); //ReynoldsU
Uh tau; //parametro di stabilizzazione

tau=(Reu<3)* hTriangle* hTriangle/(12nu)+(Reu>=3) hTriangle/(2*U);

Uh pold;

macro div(u1,u2) (dx(u1)+dy(u2)) //
macro laplaciano(u1) (dxx(u1)+dyy(u1)) //
macro convettivo(u1,u2,u3) (u1dx(u3)+u2dy(u3)) //
macro prodint(u1,u2,v1,v2) (dx(u1)*dx(v1)+dy(u2)dy(v2)+0.5(dx(u2)*dx(v2)+
dx(u2)*dy(v1)+dy(u1)*dx(v2)+dy(u1)*dy(v1))) //

problem NavierStokesSUPG ([ux,uy,p],[uu,vv,pp],init=i) =
int2d(Th)(alpha* (uxuu + uyvv)
+ 2* invRe* prodint(ux,uy,uu,vv)
- p* pp* (0.000001) - p* div(uu,vv)
+ pp* div(ux,uy))
- int2d(Th)(alpha* (uoldx* uu + uoldy* vv))
+ int2d(Th)(convettivo(uoldx,uoldy,ux)* uu //termini convettivi trattati
+ convettivo(uoldx,uoldy,uy)* vv) //in modo semi-implicito
// Stabilizzazione SUPG
+ int2d(Th)(delta* tau* convettivo(uoldx,uoldy,uu)* alpha* ux)
+ int2d(Th)(delta* tau* convettivo(uoldx,uoldy,vv)* alpha* uy)
- int2d(Th)(delta* tau* convettivo(ux,uy,uu)(alpha uoldx+dx(pold)-invRe* laplaciano(uoldx)))
- int2d(Th)(delta* tau* convettivo(ux,uy,vv)(alpha uoldy+dy(pold)-invRe* laplaciano(uoldy)))
+ int2d(Th)(delta* tau* convettivo(uoldx,uoldy,uu)* convettivo(uoldx,uoldy,ux))
+ int2d(Th)(delta* tau* convettivo(uoldx,uoldy,vv)* convettivo(uoldx,uoldy,uy))
+ on(1,2,4,ux=0,uy=0) + on(3,ux=unastro,uy=0);

pold=p;
uoldx=ux; //simulazione precedente presa come condizione iniziale
uoldy=uy;
NavierStokesSUPG;
Th=adaptmesh(Th,[ux,uy]);
plot(Th, wait=1);
plot([ux,uy], wait=1);

for(i=0; i<100; i++){
uoldx=ux;
uoldy=uy;
NavierStokesSUPG;
if (!(i % 2)) // plot every 2 iteration
plot(p,[ux,uy]);}

plot(p, wait=1);
plot([ux,uy], wait=1, value=true);