Hello everyone,
I wrote a FreeFEM++ code to solve the Poisson Problem -∆u=f in Ω u=g on ∂Ω using the Modified Weak Galerkin Finite Element Method based on the paper:
[1] “A modified weak Galerkin finite element method”. X. Wang, N.S. Malluwawadu, F. Gao, T.C. McMillan.
I modified a FreeFEM++ code for the Discontinuous Galerkin FEM to solve this problem to fit the Modified Weak Galerkin FEM, but the results are not the same as in the paper. I am wondering what’s wrong with my code. MWGFEM.edp (1.4 KB)
The weak gradient was not defined correctly in this code (it must not be the usual gradient). Here is a code with the weak gradient mwgfem.edp (1.4 KB)
However still it does not give the numerical values listed in the paper (Table 1). The results in Table 1 are doubtful because the error cannot be so small (even the L2 projection of uex gives an error to uex larger that stated!).
Respected sir, I have not so much idea about weak gradient but results coming if take a penaltity parameter 1000 in the intalledges term of the formulation using P1 element. I don’t know why it is happening. But i have not found any thing about penaltity parameter in the paper.
Here, is the Code> My Code_wgm.edp (2.2 KB)
The weak gradient is defined in the paper mentioned above.
If we put a large penalty parameter I agree that the method becomes second order, indeed it approaches the conforming P1 framework, which is known.
Here the paper claims that with P1 discontinuous it is better than classical, with the penalty parameter 1/h. We don’t get it.
Respected sir, I wrote a freefem code with the weak gradient, and the error can be small enough like the paper mentioned. MWG.edp (2.4 KB)
But When I try to change the exact solution’s boundary condition u=0
to u = u_{ex}, \partial_n u = (u_x, u_y)\cdot n ,
it doesn’t work. I don’t know why. TT MWG_2.edp (2.6 KB)
Hello Dear orange, your code is giving very good results. But which formulation you used in the code??. Why you are defining two spaces??.
Do you used the formulation in the following:
This two spaces Vh and Wh should be such that Wh has one less polynomial degree. Why you taken like this ??.
I need your clarification brother.
In the paper still they define with same degree>
Regarding boundary condition: +on(1,2,3,4,uh=0) . You have to take it as you used uex=sin(pi*x)sin(piy) which is zero on boundary of the domain D=[0,1]^2.
I have checked your code with P1, P2, and P3. All gives good results.
Indeed my code works also if I replace +on(1,2,3,4,uh=0) by +int1d(Th)(1.e10*uh*vh). It seems that there is a bug in the implementation of on() on P1dc.
The code is correct.
Apply Definition 1.1 to v\in V_h^0, taking q=\nabla_w u gives (\nabla_w v,\nabla_w u)_T=-(v,\nabla\cdot(\nabla_w u))_T+(\{v\},\nabla_w u\cdot n)_{\partial T}
Replacing this in eq. (3) gives the variational formulation -\sum_T\int_Tv\nabla\cdot(\nabla_w u) +\sum_T\int_{\partial T}\{v\}(\nabla_w u)\cdot n +\sum_e h^{-1}\int_e [[u]][[v]]-\int_\Omega fv=0
for all v\in V_h^0.
You write then Definition 1.1 for u, and call \nabla_wu=(\sigma_1,\sigma_2). Renaming the test function q as q=(\tau_1,\tau_2) you finally get the variational formulation of the code, with test functions v,\tau_1,\tau_2.
You may prefer this version mwgfem.edp (1.5 KB)
which is exactly the above formulas ((\sigma_1,\sigma_2) is replaced by (Gx,Gy)).
Thanks for quick response sir. So, can i apply this code ideas in any time dependent problem in modified weak Galerkin method and weak galerkin method sir??.
So, why the proof of stability, convergence and error estimate may be missing sir??. Numerically, we are getting optimal results. So, theoritically should come.
As i am understand from your ideas sir,
i am getting two set of equations , but i am not getting from where the term - int2d(Th)(Gxqx+Gyqy) comes. Can tell sir??. If i am understand well, we do not add after getting two set of equations otherwise the terms + intalledges(Th)( jump(uh)*jump(vh)/h/nTonEdge ) and
int2d(Th)(f*vh) will comes twice. Can you tell me please sir??.
My understanding: (q, (Gx,Gy)) giving the term - int2d(Th)(Gxqx+Gyqy). Am i right sir??. But this should be with “+” sign why “-” is in formulation sir??. (\tau_1=qx, \tau_2=qy in the code )
Also, in this way i think we are adding many extra terms in formulation sir.
I’m sorry for too late to reply.
In that artical, you can see (\nabla_w u, q)_T = -(u, \nabla\cdot q)_T + \langle\{u\}, q\rangle_{\partial T}.
So, for (\nabla_w u, \nabla_w v)_T, let \nabla_w u=\sigma, and we can get: (\nabla_w u, \nabla_w v)_T=(\sigma, \nabla_w v)_T= -(v, \nabla\cdot \sigma)_T + \langle\{v\}, \sigma\rangle_{\partial T}.
We can also get that \sigma satisfies following equation: for any \tau (\sigma, \tau)_T = (\nabla_w u, \tau)_T= -(u, \nabla\cdot \tau)_T + \langle\{u\}, \tau\rangle_{\partial T}.
Now we have following equations: (\nabla_w u, \nabla_w v)+s(u,v)=(f,v), (\sigma, \nabla_w v)= -(v, \nabla\cdot \sigma) + \langle\{v\}, \sigma\rangle_{\partial T_h}, (\sigma, \tau)= -(u, \nabla\cdot \tau)_T + \langle\{u\}, \tau\rangle_{\partial T_h}.
