G-shape generation argument added

.gitignore

@@ -24,3 +24,4 @@ scadyn
 /Geometries/
 */T/**
 /T/
+*.h5

other/generate_gellip.py

@@ -82,7 +82,7 @@ def deform_mesh(mesh):
 
    return X
 
-def draw_mesh(meshname, mesh):
+def draw_mesh(meshname, mesh, refinement):
    fig = plt.figure(figsize=(6, 6),frameon=False)
    ax = mplot3d.Axes3D(fig)
 
@@ -104,7 +104,8 @@ def draw_mesh(meshname, mesh):
 
    # Generate the tetrahedral 3d mesh for the particle. Note that optimal 
    # refinement can be hard to automatize with tetgen, so quartet is used!
-   tetramesh = pymesh.tetrahedralize(mesh,0.15,engine='quartet')     
+   tetramesh = pymesh.tetrahedralize(mesh,refinement,engine='quartet')     
+   print('Number of tetras: ' + str(tetramesh.num_elements))
 
    param_r = eps_r*np.ones(tetramesh.voxels.shape[0])
    param_i = eps_i*np.ones(tetramesh.voxels.shape[0])
@@ -124,10 +125,11 @@ def args(argv):
    meshname = "mesh"
    gid   = 1      # G-ellipsoid id
    gotGid = False
+   refinement = 0.15
    try:
-      opts, args = getopt.getopt(argv,"i:o:f:")
+      opts, args = getopt.getopt(argv,"i:o:f:r:")
    except getopt.GetoptError:
-      print('generate_gellip.py -i <G-ID> -o <meshname> -f <saveOnlyPly>')
+      print('generate_gellip.py -i <G-ID> -o <meshname> -f <saveOnlyPly> -r <refinement_level>')
       sys.exit(2)
    for opt, arg in opts:
       if opt in ('-i'):
@@ -137,6 +139,12 @@ def args(argv):
          except ValueError:
             print('Argument of -i is an integer, so try again!')
             sys.exit(2)
+      if opt in ('-r'):
+         try:
+            refinement = float(arg)
+         except ValueError:
+            print('Argument of -r is a float, so try again!')
+            sys.exit(2)
       if opt in ('-o'):
          meshname = arg
       if opt in ('-f'):
@@ -149,10 +157,10 @@ def args(argv):
             saveOnlyPly = True
    if gotGid:
       meshname = meshname + str(gid)
-   return meshname, saveOnlyPly, gid
+   return meshname, saveOnlyPly, gid, refinement
 
 if __name__ == "__main__":
-   meshname, saveOnlyPly, gid = args(sys.argv[1:])
+   meshname, saveOnlyPly, gid, refinement = args(sys.argv[1:])
    seed(gid)
    ellipsoid = genellip()  
 
@@ -163,6 +171,6 @@ def args(argv):
    if(saveOnlyPly):
       pymesh.save_mesh(meshname+".ply", gellip)
    else:
-      draw_mesh(meshname, gellip)
+      draw_mesh(meshname, gellip, refinement)
 
 

other/generate_gsphere.py

@@ -95,7 +95,7 @@ def r_gsphere(a_lm, b_lm, z, phi, beta):
          s = s + legp[mm,ll]*(a_lm[ii]*CPHI[mm]+b_lm[ii]*SPHI[mm])
    return np.exp(s-0.5*beta**2)
 
-def draw_mesh(meshname,mesh):
+def draw_mesh(meshname, mesh, refinement):
    fig = plt.figure(figsize=(6, 6),frameon=False)
    ax = mplot3d.Axes3D(fig)
 
@@ -117,8 +117,9 @@ def draw_mesh(meshname,mesh):
 
    # Generate the tetrahedral 3d mesh for the particle. Note that optimal 
    # refinement can be hard to automatize with tetgen, so quartet is used!
-   tetramesh = pymesh.tetrahedralize(mesh,0.15,engine='quartet')   
-
+   tetramesh = pymesh.tetrahedralize(mesh,refinement,engine='quartet')     
+   print('Number of tetras: ' + str(tetramesh.num_elements))
+
    param_r = eps_r*np.ones(tetramesh.voxels.shape[0])
    param_i = eps_i*np.ones(tetramesh.voxels.shape[0])
 
@@ -137,10 +138,11 @@ def args(argv):
    meshname = "mesh"
    gid   = 1      
    gotGid = False
+   refinement = 0.15
    try:
-      opts, args = getopt.getopt(argv,"i:o:f:")
+      opts, args = getopt.getopt(argv,"i:o:f:r:")
    except getopt.GetoptError:
-      print('generate_gsphere.py -i <G-ID> -o <meshname> -f <saveOnlyPly>')
+      print('generate_gsphere.py -i <G-ID> -o <meshname> -f <saveOnlyPly> -r <refinement_level>')
       sys.exit(2)
    for opt, arg in opts:
       if opt in ('-i'):
@@ -150,6 +152,12 @@ def args(argv):
          except ValueError:
             print('Argument of -i is an integer, so try again!')
             sys.exit(2)
+      if opt in ('-r'):
+         try:
+            refinement = float(arg)
+         except ValueError:
+            print('Argument of -r is a float, so try again!')
+            sys.exit(2)
       if opt in ('-o'):
          meshname = arg
       if opt in ('-f'):
@@ -162,10 +170,10 @@ def args(argv):
             saveOnlyPly = True
    if gotGid:
       meshname = meshname + str(gid)
-   return meshname, saveOnlyPly, gid
+   return meshname, saveOnlyPly, gid, refinement
 
 if __name__ == "__main__":
-   meshname, saveOnlyPly, gid = args(sys.argv[1:])
+   meshname, saveOnlyPly, gid, refinement  = args(sys.argv[1:])
 
    gamma = gamma*np.pi/180
    ell = 2.0*np.sin(0.5*gamma)
@@ -189,5 +197,5 @@ def args(argv):
    if(saveOnlyPly):
       pymesh.save_mesh(meshname+".ply", gsphere)
    else:
-      draw_mesh(meshname,gsphere)
+      draw_mesh(meshname, gsphere, refinement)