slim3D: add constant diff and visc + set z0

+ update Danilov test case

benchmarks/swe3d/danilov/square.geo

-//L = 9e5;
-L = 9e5/60;
+L = 9e5;
 l = 8e5;
 Point(1) = {-L, -l, 0};
 Point(2) = {L, -l, 0};
@@ -19,13 +18,9 @@ sf - > /2
 c - > /4
 */
 
-Transfinite Line{1} = 2;
-Transfinite Line{3} = 2;
-
-
-//Transfinite Line{1} = 6*20 /2;///8;
+Transfinite Line{1} = 6*20 /2;///8;
 Transfinite Line{2} = 7*15 /2;///8;
-//Transfinite Line{3} = 6*20 /2;///8;
+Transfinite Line{3} = 6*20 /2;///8;
 Transfinite Line{4} = 7*15 /2;///8;
 Transfinite Surface{1} = {1,2,3,4};
 

dgpy/scripts/slim3d.py

@@ -45,6 +45,7 @@ class Slim3d_equations :
         self._wind_stress = None
         self._boundary_coast = []
         self._boundary_open = []
+        self._z0 = [0.005, 0.02]
 
     def set_implicit_vertical_diffusion(self, flag):
         slimSolver = self._slimSolver
@@ -78,15 +79,17 @@ class Slim3d_equations :
         self._hor_diff_fact = factor
         self._hor_diff_max = maximum
 
-    def set_vertical_viscosity(self, viscosity):
-        if viscosity != "gotm":
-            dgpy.Msg.Fatal("Unknown viscosity for vertical_viscosity : "+viscosity)
-        self._vertical_viscosity = viscosity
+    def set_vertical_viscosity(self, mode, constant_value=1e-5):
+        if (mode != "gotm") and (mode != "constant"):
+            dgpy.Msg.Fatal("Unknown viscosity for vertical_viscosity : "+mode)
+        self._vertical_viscosity = mode
+        self._vertical_viscosity_value = constant_value
         
-    def set_vertical_diffusivity(self, diffusivity):
-        if diffusivity != "gotm":
-            dgpy.Msg.Fatal("Unknown diffusivity for vertical_diffusivity : "+diffusivity)
-        self._vertical_diffusivity = diffusivity
+    def set_vertical_diffusivity(self, mode, constant_value=1e-5):
+        if (mode != "gotm") and (mode != "constant"):
+            dgpy.Msg.Fatal("Unknown diffusivity for vertical_diffusivity : "+mode)
+        self._vertical_diffusivity = mode
+        self._vertical_diffusivity_value = constant_value
 
     def set_coriolis(self, coriolis):
         self._coriolis = coriolis
@@ -94,8 +97,12 @@ class Slim3d_equations :
     def set_atmospheric_pressure(self, atmPress):
         self._atmPress = atmPress
 
-    def set_bottom_friction(self, flag):
+    def set_bottom_friction(self, flag, z0B=None, z0S=None):
         self._slimSolver.setComputeBottomFriction(flag)
+        if z0B:
+            self._z0[0] = z0B
+        if z0S:
+            self._z0[1] = z0S
 
     def set_linear_density(self, salinity=False, temperature=False, constant_coefficient=0, linear_coefficient=0):
         if (temperature and salinity) or (not temperature and not salinity):

dgpy/scripts/slim3d_setup.py

@@ -99,9 +99,14 @@ def slim3d_setup(loop):
             eq._areaFunc = dgpy.functionPrecomputedExtrusion(slimSolver.extrusion(), 3, area2d.getFunction())
         f.kappahTotal = dgpy.okubo(eq._hor_diff_fact, eq._hor_diff_max, eq._areaFunc)
 
+    if eq._vertical_viscosity == 'constant':
+        f.nuvFunc = dgpy.functionConstant(eq._vertical_viscosity_value)
+    if eq._vertical_diffusivity == 'constant':
+        f.kappavFunc = dgpy.functionConstant(eq._vertical_diffusivity_value)
+    
     if slimSolver.getComputeBottomFriction():
-        f.z0BFunc = dgpy.functionConstant(0.005)
-        f.z0SFunc = dgpy.functionConstant(0.02)
+        f.z0BFunc = dgpy.functionConstant(eq._z0[0])
+        f.z0SFunc = dgpy.functionConstant(eq._z0[1])
     
     if eq._coriolis:
         coriolisFunc = slim_private._load_function(eq._coriolis, slimSolver.groups2d)