# Problem in implementation of DG code for Stokes-Darcy interface problem **URL:** \ **Category:** General Discussion\ **Created:** [December 5, 2024, 6:46pm UTC](https://community.freefem.org/t/problem-in-implementation-of-dg-code-for-stokes-darcy-interface-problem/3616 "2024-12-05T18:46:55Z")\ **Posts on this page:** 8\ **Page:** 1
**Author:** ![Sweta](https://avatars.discourse-cdn.com/v4/letter/s/258eb7/32.png) [@Sweta](https://community.freefem.org/u/Sweta)\ **Post date:** [December 5, 2024, 6:46pm UTC](https://community.freefem.org/t/problem-in-implementation-of-dg-code-for-stokes-darcy-interface-problem/3616/1 "2024-12-05T18:46:55Z")
Dear all, I have implemented standard finite element code for coupled Stokes-Darcy problem (Here I am attaching the link of the pdf where the model problem has been described [[https://drive.google.com/drive/folders/1yjBxjZfF5rSFCa0V0SbETqHlfIeIje4Y?usp=drive\_link](https://drive.google.com/drive/folders/1yjBxjZfF5rSFCa0V0SbETqHlfIeIje4Y?usp=drive_link)]). It is running well. But Discontinuous Galerkin(DG) SIPG code for aforementioned problem doesn’t not providing desired results. The problem contains interface condition which occurs coupling between Stokes and Darcy equation. But I think the DG code is not approximating the interface properly and taking zero value at the interface. But in my model problem the specific boundary condition at the interface is not prescribed. Below, I have attached the DG code for the Stokes-Darcy problem. real nu = 1./1.; verbosity = 0; macro div(ux,uy) (dx(ux) + dy(uy))// int NN = 32; border D1(t=0, 1){x=t; y=0; label=1;} border D2(t=0, 1){x=1; y=t; label=2;} border D3(t=0, 1){x=1-t; y=1; label=3;}//interface of darcy region border D4(t=0, 1){x=0; y=1-t; label=4;} border N1(t=0, 1){x=t; y=1; label=5;} //interface of Stokes region border N2(t=1, 2){x=1; y=t; label=6;} border N3(t=0, 1){x=1-t; y=2; label=7;} border N4(t=0, 1){x=0; y=2-t; label=8;} mesh Th1 = buildmesh( N1(NN) + N2(NN)+ N3(NN) + N4(NN)); mesh Th2 = buildmesh( D1(NN) + D2(NN) + D3(NN) + D4(NN)); mesh Th = Th1+Th2; //mesh plot(Th); fespace Vh(Th1,P1dc); // Discontinous P1 finite element Vh ux, uy, uu1, uu2, dux,duy, vx, vy; fespace Ph1(Th1, P0); // Discontinous P0 finite element Ph1 p1, pp1,q1,dp1; fespace Ph2(Th2, P1dc);// Discontinous P1 finite element Ph2 p2, pp2,q2,dp2; real pena=50.; // a paramater to add penalisation Vh Uex1= x^2\*(y-1)^2+y; //exact solution Vh Uex2= -2_x_(y-1)^3/(3)+2-pi_sin(pi_x); //exact solution Ph1 p1ex= (2-pi_sin(pi_x))_sin(pi_y/(2)); //exact solution Ph2 p2ex= (2-pi_sin(pi_x))_(1-y-cos(pi_y));//exact solution func U1b= x^2\*(y-1)^2+y; //boundary condition func U2b= -2_x_(y-1)^3/(3)+2-pi_sin(pi_x); //boundary condition func P22b= (2-pi_sin(pi_x))_(1-y-cos(pi_y)); //boundary condition // Define right hand side for Stokes and Darcy func f1 = (-2_x^2-2_(y-1)^2-(pi^2)_cos(pi_x)_sin(pi_y/(2.))); func f2 = (4_x_(y-1)-(pi^3)_sin(pi_x)+(2-pi_sin(pi_x))_(pi/2)cos(piy/(2.))); func f3 = -(pi^3)sin(pix)_(1-y-cos(pi_y))-(2-pi_sin(pi_x))_(pi^2)_cos(pi_y); solve Stokes([dux,duy,dp1], [vx,vy,q1])= int2d(Th1)( nu\*(dx(dux)_dx(vx) + dy(dux)dy(vx) + dx(duy)dx(vy) + dy(duy)dy(vy)) - div(dux,duy) \* q1 - dp1div(vx,vy) - 1e-4dp1q1) + intalledges(Th1)( (nu_(mean(dx(dux))_N.x_jump(vx)+mean(dy(dux))_N.y_jump(vx)+mean(dx(duy))_N.x_jump(vy) +mean(dy(duy))_N.y_jump(vy) +mean(dx(vx))_N.x_jump(dux)+mean(dy(vx))_N.y_jump(dux)+mean(dx(vy))_N.x_jump(duy)+mean(dy(vy))_N.y_jump(duy) + (pena/lenEdge)_(jump(dux)jump(vx)+jump(duy)jump(vy))) -(jump(vx)N.x+jump(vy)N.y)mean(dp1) -(jump(dux)N.x+jump(duy)N.y)mean(q1) )) +int1d(Th1,6,7,8)(nu( (dx(vx))N.x(U1b)+(dy(vx))N.y(U1b)+(dx(vy))N.x(U2b)+(dy(vy))N.y(U2b) - (pena/lenEdge)((U1b)(vx)+(U2b)(vy))) -((U1b)N.x+(U2b)N.y)(q1) ) +int1d(Th1,5)((dp2(vx0+vy(-1)))) //interface +int1d(Th1,5)((dux(1)+duy0)_(vx\*(1)+vy_0)) //interface -int2d(Th1)([f1,f2]'_[vx,vy]); ``` solve DAR(dp2, q2)= int2d(Th2)( dx(dp2)*dx(q2) + dy(dp2)*dy(q2) ) + intalledges(Th2)(( ((mean(dx(dp2))*N.x*jump(q2)+mean(dy(dp2))*N.y*jump(q2) +mean(dx(q2))*N.x*jump(dp2)+mean(dy(q2))*N.y*jump(dp2) + (pena/lenEdge)*(jump(dp2)*jump(q2))) /nTonEdge) )) +int1d(Th2,1,2,4)(( (dx(q2))*N.x*(P22b)+(dy(q2))*N.y*(P22b) - (pena/lenEdge)*((P22b)*(q2))) ) -int1d(Th2,3)((dux*0+duy*(-1))*q2) //interface -int2d(Th2)(f3*q2) ; plot([Uex1,Uex2],wait=1, cmm="exact solution"); plot([dux,duy],wait=1, cmm="approximate solution"); plot(p2ex,wait=1, cmm="exact solution"); plot(dp2,wait=1, cmm="approximate solution"); ``` Please, help me out to resolve this problem. Thanks in advance. ---
**Author:** ![fb77](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/fb77/32/3796_2.png) [@fb77](https://community.freefem.org/u/fb77)\ **Post date:** [December 6, 2024, 7:45pm UTC](https://community.freefem.org/t/problem-in-implementation-of-dg-code-for-stokes-darcy-interface-problem/3616/2 "2024-12-06T19:45:05Z")
When you solve Stokes, you have an interface term involving dp2. At this time, Darcy is not yet solved, thus it takes dp2=0. You should try to solve both problems in a coupled way, using a composite space, and a variational formulation with Stokes and Darcy together. ---
**Author:** ![Sweta](https://avatars.discourse-cdn.com/v4/letter/s/258eb7/32.png) [@Sweta](https://community.freefem.org/u/Sweta)\ **Post date:** [December 8, 2024, 9:25am UTC](https://community.freefem.org/t/problem-in-implementation-of-dg-code-for-stokes-darcy-interface-problem/3616/3 "2024-12-08T09:25:42Z")
Thank you sir for your suggestion. I have tried accordingly but it is not working and again at interface they are taking zero value. Also, I would like to mention that the decoupled technique is nicely working for standard Finite element but not working for Discontinuous galerkin method. ---
**Author:** ![Sweta](https://avatars.discourse-cdn.com/v4/letter/s/258eb7/32.png) [@Sweta](https://community.freefem.org/u/Sweta)\ **Post date:** [December 19, 2024, 2:57pm UTC](https://community.freefem.org/t/problem-in-implementation-of-dg-code-for-stokes-darcy-interface-problem/3616/4 "2024-12-19T14:57:29Z")
I have tried to solve the coupled problem in this way but it is not working. Is there any other way to solve it? Any suggestions in this direction will be great help for me. ---
**Author:** ![fb77](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/fb77/32/3796_2.png) [@fb77](https://community.freefem.org/u/fb77)\ **Post date:** [December 20, 2024, 12:20pm UTC](https://community.freefem.org/t/problem-in-implementation-of-dg-code-for-stokes-darcy-interface-problem/3616/5 "2024-12-20T12:20:35Z")
I have defined a single mesh with two regions. The coupled problem is solved nicely with the spaces P2,P1,P2. [composite-space-DGStokesDarcy.edp](https://community.freefem.org/uploads/short-url/zArgpqI7yp5dEqRuzeIzSM5cM7B.edp) (4.1 KB) For DG it works with composite, but some corrections on the interface involve the exact solution (don’t know how to remove them) [composite-space-DGStokesDarcy.edp](https://community.freefem.org/uploads/short-url/bRJ7zoVupMTaqwdHntMHvWqyBBE.edp) (4.6 KB) ---
**Author:** ![fb77](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/fb77/32/3796_2.png) [@fb77](https://community.freefem.org/u/fb77)\ **Post date:** [December 21, 2024, 4:29pm UTC](https://community.freefem.org/t/problem-in-implementation-of-dg-code-for-stokes-darcy-interface-problem/3616/6 "2024-12-21T16:29:31Z")
Sorry for the DG implementation the code was wrong. A correct version is [composite-space-DGStokesDarcy.edp](https://community.freefem.org/uploads/short-url/nJQLZPcMuzN3uLBzurdFy3y66jO.edp) (4.4 KB) ---
**Author:** ![Sweta](https://avatars.discourse-cdn.com/v4/letter/s/258eb7/32.png) [@Sweta](https://community.freefem.org/u/Sweta)\ **Post date:** [December 24, 2024, 10:46am UTC](https://community.freefem.org/t/problem-in-implementation-of-dg-code-for-stokes-darcy-interface-problem/3616/7 "2024-12-24T10:46:02Z")
The plot is correct but the code is not giving us the correct convergence rate for error. Could you suggest a way to resolve the issue? ---
**Author:** ![fb77](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/fb77/32/3796_2.png) [@fb77](https://community.freefem.org/u/fb77)\ **Post date:** [December 27, 2024, 11:27am UTC](https://community.freefem.org/t/problem-in-implementation-of-dg-code-for-stokes-darcy-interface-problem/3616/8 "2024-12-27T11:27:42Z")
Dear Sweta, I don’t know how to improve it. From my point of view it converges nicely. I have showed you how to use the composite framework in the DG context, by imitation of the model code > Next, finding the best scheme is your job! Francois.