Merge pull request #28 from KVSlab/add-rigid-cylinder-mesh

Add cylinder mesh for (rigid) Cylinder problem

Author: hkjeldsberg

src/oasismove/problems/Cylinder.py

@@ -3,28 +3,7 @@
 __copyright__ = "Copyright (C) 2014 " + __author__
 __license__ = "GNU Lesser GPL version 3 or any later version"
 
-import os
-import platform
-
-from dolfin import Mesh, AutoSubDomain, near
-
-if not os.path.isfile("cylinder.xml"):
-    if platform.system() == "Linux":
-        os.system("wget -O cylinder.xml https://www.dropbox.com/s/d78g4cyjxl3ylay/cylinder.xml?dl=0")
-    elif platform.system() == "Darwin":
-        os.system("curl -L https://www.dropbox.com/s/d78g4cyjxl3ylay/cylinder.xml?dl=0 -o cylinder.xml")
-    else:
-        raise ImportError("Could not determine platform")
-
-    # try:
-    # os.system("gmsh mesh/cylinder.geo -2 -o mesh/cylinder.msh")
-    # os.system("dolfin-convert mesh/cylinder.msh mesh/cylinder.xml")
-    # os.system("rm mesh/cylinder.msh")
-    # except RuntimeError:
-    # os.system("wget -O cylinder.xml https://www.dropbox.com/s/d78g4cyjxl3ylay/cylinder.xml?dl=0")
-    # raise "Gmsh is required to run this demo"
-
-mesh = Mesh("cylinder.xml")
+from dolfin import AutoSubDomain, near
 
 H = 0.41
 L = 2.2
@@ -47,14 +26,18 @@
 
 
 # Overload post_import_problem to choose between the two cases
-def post_import_problem(NS_parameters, commandline_kwargs, **NS_namespace):
+def post_import_problem(NS_parameters, mesh, commandline_kwargs, **NS_namespace):
+    # If the mesh is a callable function, then create the mesh here.
+    if callable(mesh):
+        mesh = mesh(**NS_parameters)
+
     """ Choose case - case could be defined through command line."""
     NS_parameters.update(commandline_kwargs)
     case = NS_parameters['case'] if 'case' in NS_parameters else 1
     Um = cases[case]["Um"]
     Re = cases[case]["Re"]
     Umean = 2. / 3. * Um
     nu = Umean * D / Re
-    NS_parameters.update(nu=nu, Re=Re, Um=Um, Umean=Umean)
+    NS_parameters.update(nu=nu, Re=Re, Um=Um, Umean=Umean, mesh=mesh)
 
     return NS_parameters

src/oasismove/problems/NSCoupled/Cylinder.py

@@ -5,8 +5,8 @@
 __copyright__ = "Copyright (C) 2014 " + __author__
 __license__ = "GNU Lesser GPL version 3 or any later version"
 
-from ..Cylinder import *
-from ..NSCoupled import *
+from oasismove.problems.NSCoupled import *
+from oasismove.problems.Cylinder import *
 
 
 # Override some problem specific parameters
@@ -15,7 +15,9 @@ def problem_parameters(NS_parameters, scalar_components, **NS_namespace):
         omega=1.0,
         max_iter=100,
         plot_interval=10,
-        velocity_degree=2)
+        velocity_degree=2,
+        mesh_path="src/oasismove/mesh/cylinder.xml",
+    )
 
     scalar_components += ["c", "d"]
 
@@ -24,6 +26,15 @@ def scalar_source(c_, d_, **NS_namespace):
     return {"c": -Constant(0.1) * c_ * c_, "d": -Constant(0.25) * c_ * d_ * d_}
 
 
+def mesh(mesh_path, dt, **NS_namespace):
+    # Import mesh
+    print(mesh_path)
+    mesh = Mesh(mesh_path)
+
+    print_mesh_information(mesh, dt, dim=2)
+    return mesh
+
+
 def create_bcs(VQ, Um, CG, V, element, **NS_namespace):
     inlet = Expression(("4.*{0}*x[1]*({1}-x[1])/pow({1}, 2)".format(Um, H), "0"), element=V)
     ux = Expression(("0.00*x[1]", "-0.00*(x[0]-{})".format(center)), element=V)

src/oasismove/problems/NSfracStep/Cylinder.py

@@ -9,9 +9,8 @@
 from os import getcwd
 
 import matplotlib.pyplot as plt
-
-from ..Cylinder import *
-from ..NSfracStep import *
+from oasismove.problems.NSfracStep import *
+from oasismove.problems.Cylinder import *
 
 
 def problem_parameters(commandline_kwargs, NS_parameters, scalar_components,
@@ -35,7 +34,9 @@ def problem_parameters(commandline_kwargs, NS_parameters, scalar_components,
             plot_interval=10,
             velocity_degree=2,
             print_intermediate_info=100,
-            use_krylov_solvers=True)
+            use_krylov_solvers=True,
+            mesh_path="src/oasismove/mesh/cylinder.xml",
+        )
         NS_parameters['krylov_solvers'].update(dict(monitor_convergence=True))
         NS_parameters['velocity_krylov_solver'].update(dict(preconditioner_type='jacobi',
                                                             solver_type='bicgstab'))
@@ -44,6 +45,15 @@ def problem_parameters(commandline_kwargs, NS_parameters, scalar_components,
     Schmidt["alfa"] = 0.1
 
 
+def mesh(mesh_path, dt, **NS_namespace):
+    # Import mesh
+    print(mesh_path)
+    mesh = Mesh(mesh_path)
+
+    print_mesh_information(mesh, dt, dim=2)
+    return mesh
+
+
 def create_bcs(V, Q, Um, H, **NS_namespace):
     inlet = Expression(
         "4.*{0}*x[1]*({1}-x[1])/pow({1}, 2)".format(Um, H), degree=2)