adding coeff

Author: Fimache

examples/Freefem/MMS/2D/incompressible.py

@@ -19,6 +19,7 @@
 from beam import (case_scene, run_reference_scene,
                   element_quad, element_tri,
                   quad_q1_rule, tri_p1_rule)
+coef=0.1
 
 
 class Incompressible(MMSCase2D):
@@ -31,15 +32,15 @@ class Incompressible(MMSCase2D):
 
     def u_ex(self, x, y, L):
         k = np.pi / L
-        return ( np.sin(k * x) * np.cos(k * y),
-                -np.cos(k * x) * np.sin(k * y))
+        return (coef* np.sin(k * x) * np.cos(k * y),
+                -coef* np.cos(k * x) * np.sin(k * y))
 
     def grad_u_ex(self, x, y, L):
         k = np.pi / L
-        dux_dx =  k * np.cos(k * x) * np.cos(k * y)
-        dux_dy = -k * np.sin(k * x) * np.sin(k * y)
-        duy_dx =  k * np.sin(k * x) * np.sin(k * y)
-        duy_dy = -k * np.cos(k * x) * np.cos(k * y)
+        dux_dx =  coef* k * np.cos(k * x) * np.cos(k * y)
+        dux_dy = - coef* k * np.sin(k * x) * np.sin(k * y)
+        duy_dx =  coef * k * np.sin(k * x) * np.sin(k * y)
+        duy_dy = - coef *k * np.cos(k * x) * np.cos(k * y)
         return np.array([[dux_dx, dux_dy],
                          [duy_dx, duy_dy]])
 
@@ -48,12 +49,12 @@ def source(self, x, y, E, nu, L, dim):
         k  = np.pi / L
         ux =  np.sin(k * x) * np.cos(k * y)
         uy = -np.cos(k * x) * np.sin(k * y)
-        d2ux_dxx = -k**2 * ux
-        d2ux_dyy = -k**2 * ux
-        d2ux_dxy = -k**2 * np.cos(k * x) * np.sin(k * y)
-        d2uy_dxx = -k**2 * uy
-        d2uy_dyy = -k**2 * uy
-        d2uy_dxy =  k**2 * np.sin(k * x) * np.cos(k * y)
+        d2ux_dxx = - coef* k**2 * ux
+        d2ux_dyy = - coef* k**2 * ux
+        d2ux_dxy = -coef* k**2 * np.cos(k * x) * np.sin(k * y)
+        d2uy_dxx = -coef* k**2 * uy
+        d2uy_dyy = -coef *k**2 * uy
+        d2uy_dxy =  coef *k**2 * np.sin(k * x) * np.cos(k * y)
         fx = -((lam + 2*mu) * d2ux_dxx + lam * d2uy_dxy
              + mu * (d2ux_dyy + d2uy_dxy))
         fy = -(mu * (d2ux_dxy + d2uy_dxx) + lam * d2ux_dxy
@@ -97,8 +98,8 @@ def apply_bcs(self, Beam, nodes_2d, L, dim):
 
 
 mms         = Incompressible()
-createScene = case_scene(mms, element_quad)
+createScene = case_scene(mms, element_tri)
 
 
 if __name__ == "__main__":
-    run_reference_scene(element_quad, mms)
+    run_reference_scene(element_tri, mms)