Skip to content

Commit f43e6fd

Browse files
authored
Merge pull request SPECFEM#1189 from danielpeter/devel
updates cubit2specfem2d python script (to be compatible w/ python3)
2 parents cea5867 + 71c33b6 commit f43e6fd

File tree

2 files changed

+93
-54
lines changed

2 files changed

+93
-54
lines changed

utils/cubit2specfem2d/cubit2specfem2d.py

Lines changed: 65 additions & 38 deletions
Original file line numberDiff line numberDiff line change
@@ -31,15 +31,26 @@
3131
#
3232
# The names of the block and the entities types must match the ones given during the definition of the class mesh on this file :
3333
# Below :
34-
# class mesh(object,mesh_tools):
34+
# class mesh(mesh_tools):
3535
# """ A class to store the mesh """
3636
# def __init__(self):
3737
#
3838
#!! Warning : a block in cubit != quad !! A block is a group of something (quads, edges, volumes, surfaces...)
3939
# On this case the blocks are used to gather faces corresponding to different materials and edges corresponding to free surfaces,
4040
# absorbing surfaces, topography or axis
4141
from __future__ import print_function
42-
import cubit
42+
import sys
43+
44+
try:
45+
import cubit
46+
except:
47+
print('error importing cubit')
48+
sys.exit()
49+
50+
try:
51+
set
52+
except NameError:
53+
from sets import Set as set
4354

4455
class mtools(object):
4556
"""docstring for mtools"""
@@ -70,7 +81,7 @@ def mesh_it(self):
7081
command = "mesh surf "+str(surf)
7182
cubit.cmd(command)
7283

73-
class block_tools:
84+
class block_tools(object):
7485
def __int__(self):
7586
pass
7687
def create_blocks(self,mesh_entity,list_entity = None,):
@@ -251,10 +262,12 @@ def sorter(x, y):
251262
print(self.ddt,self.dr)
252263
print('Deltat minimum => edge:'+str(self.ddt[0][0])+' dt: '+str(self.ddt[0][1]))
253264
print('Minimum frequency resolved => edge:'+str(self.dr[0][0])+' frequency: '+str(self.dr[0][1]))
265+
cubit.cmd('set info on')
266+
cubit.cmd('set echo on')
254267
return self.ddt[0],self.dr[0]
255268

256269

257-
class mesh(object,mesh_tools):
270+
class mesh(mesh_tools):
258271
""" A class to store the mesh """
259272
def __init__(self):
260273
super(mesh, self).__init__()
@@ -309,12 +322,13 @@ def block_definition(self):
309322
name = cubit.get_exodus_entity_name('block',block) # Contains the name of the blocks
310323
ty = cubit.get_block_element_type(block) # Contains the block element type (QUAD4...)
311324
if ty == self.face: # If we are dealing with a block containing faces
325+
print("block: ",name," contains faces")
312326
nAttributes = cubit.get_block_attribute_count(block)
313327
if (nAttributes != 1 and nAttributes != 6):
314328
print('Blocks not properly defined, 2d blocks must have one attribute (material id) or 6 attributes')
315329
return None,None,None,None,None,None,None,None
316330
flag=int(cubit.get_block_attribute_value(block,0)) # Fetch the first attribute value (containing material id)
317-
print("nAttributes : ",nAttributes)
331+
print(" nAttributes : ",nAttributes)
318332
if nAttributes == 6:
319333
self.write_nummaterial_velocity_file = True
320334
velP = cubit.get_block_attribute_value(block,1) # Fetch the first attribute value (containing P wave velocity)
@@ -343,21 +357,26 @@ def block_definition(self):
343357
# # (index 0 : pml_x_acoust, index 1 : pml_z_acoust, index 2 : pml_xz_acoust,
344358
# # index 3 : pml_x_elast, index 4 : pml_z_elast, index 5 : pml_xz_elast)
345359
elif ty == self.edge: # If we are dealing with a block containing edges
360+
print("block: ",name," contains edges")
346361
block_bc_flag.append(2) # Append "2" to block_bc_flag
347362
block_bc.append(block) # Append block id to block_bc
348363
bc[name] = 2 # Associate the name of the block with its connectivity : an edge has connectivity = 2
349364
if name == self.topo:
350365
self.topo_mesh = True
366+
print(" topo_mesh: ",self.topo_mesh)
351367
topography = block # If the block considered refered to topography store its id in "topography"
352368
if name in self.forcing_boun_name:
353369
self.forcing_mesh = True
370+
print(" forcing_mesh: ",self.forcing_mesh)
354371
forcing_boun[self.forcing_boun_name.index(name)] = block
355372
# -> Put it at the correct position in abs_boun (index 0 : bottom, index 1 : right, index 2 : top, index 3 : left)
356373
if name == self.axisname:
357374
self.axisymmetric_mesh = True
375+
print(" axisymmetric_mesh: ",self.axisymmetric_mesh)
358376
axisId = block # AXISYM If the block considered refered to the axis store its id in "axisId"
359377
if name in self.abs_boun_name : # If the block considered refered to one of the boundaries
360378
self.abs_mesh = True
379+
print(" abs_mesh: ",self.abs_mesh)
361380
abs_boun[self.abs_boun_name.index(name)] = block
362381
# -> Put it at the correct position in abs_boun (index 0 : bottom, index 1 : right, index 2 : top, index 3 : left)
363382
else:
@@ -415,6 +434,8 @@ def mesh_write(self,mesh_name):
415434
""" Write mesh (quads ids with their corresponding nodes ids) on file : mesh_name """
416435
meshfile = open(mesh_name,'w')
417436
print('Writing '+mesh_name+'.....')
437+
cubit.cmd('set info off') # Turn off return messages from Cubit commands
438+
cubit.cmd('set echo off') # Turn off echo of Cubit commands
418439
num_elems = cubit.get_quad_count() # Store the number of elements
419440
toWritetoFile = [""]*(num_elems+1)
420441
toWritetoFile[0] = str(num_elems)+'\n'
@@ -433,10 +454,14 @@ def mesh_write(self,mesh_name):
433454
meshfile.writelines(toWritetoFile)
434455
meshfile.close()
435456
print('Ok num elements/write =',str(num_elems), str(num_write))
457+
cubit.cmd('set info on')
458+
cubit.cmd('set echo on')
436459
def material_write(self,mat_name):
437460
""" Write quads material on file : mat_name """
438461
mat = open(mat_name,'w')
439462
print('Writing '+mat_name+'.....')
463+
cubit.cmd('set info off') # Turn off return messages from Cubit commands
464+
cubit.cmd('set echo off') # Turn off echo of Cubit commands
440465
num_elems = cubit.get_quad_count() # Store the number of elements
441466
toWritetoFile = [""]*num_elems
442467
print('block_mat:',self.block_mat)
@@ -450,12 +475,14 @@ def material_write(self,mat_name):
450475
mat.writelines(toWritetoFile)
451476
mat.close()
452477
print('Ok')
478+
cubit.cmd('set info on')
479+
cubit.cmd('set echo on')
453480
def pmls_write(self,pml_name):
454481
""" Write pml elements on file : mat_name """
455-
cubit.cmd('set info off') # Turn off return messages from Cubit commands
456-
cubit.cmd('set echo off') # Turn off echo of Cubit commands
457482
pml_file = open(pml_name,'w')
458483
print('Writing '+pml_name+'.....')
484+
cubit.cmd('set info off') # Turn off return messages from Cubit commands
485+
cubit.cmd('set echo off') # Turn off echo of Cubit commands
459486
npml_elements = 0
460487
#id_element = 0 # Global id
461488
indexFile = 1
@@ -493,6 +520,8 @@ def nodescoord_write(self,nodecoord_name):
493520
""" Write nodes coordinates on file : nodecoord_name """
494521
nodecoord = open(nodecoord_name,'w')
495522
print('Writing '+nodecoord_name+'.....')
523+
cubit.cmd('set info off') # Turn off return messages from Cubit commands
524+
cubit.cmd('set echo off') # Turn off echo of Cubit commands
496525
node_list = cubit.parse_cubit_list('node','all') # Import all the nodes of the model
497526
num_nodes = len(node_list) # Total number of nodes
498527
nodecoord.write('%10i\n' % num_nodes) # Write the number of nodes on the first line
@@ -502,26 +531,27 @@ def nodescoord_write(self,nodecoord_name):
502531
nodecoord.write(txt) # Write x and z coordinates on the file -> Model must be in x,z coordinates. TODO
503532
nodecoord.close()
504533
print('Ok')
534+
cubit.cmd('set info on') # Turn on return messages from Cubit commands
535+
cubit.cmd('set echo on') # Turn on echo of Cubit commands
505536
def free_write(self,freename): #freename = None):
506537
""" Write free surface on file : freename """
507-
cubit.cmd('set info off') # Turn off return messages from Cubit commands
508-
cubit.cmd('set echo off') # Turn off echo of Cubit commands
509-
cubit.cmd('set journal off') # Do not save journal file
510-
from sets import Set
511538
# if not freename: freename = self.freename
512539
freeedge = open(freename,'w')
513540
print('Writing '+freename+'.....')
541+
cubit.cmd('set info off') # Turn off return messages from Cubit commands
542+
cubit.cmd('set echo off') # Turn off echo of Cubit commands
543+
cubit.cmd('set journal off') # Do not save journal file
514544
if self.topo_mesh:
515545
for block,flag in zip(self.block_bc,self.block_bc_flag): # For each 1D block
516546
if block == self.topography: # If the block correspond to topography
517-
edges_all = Set(cubit.get_block_edges(block)) # Import all topo edges id as a Set
547+
edges_all = set(cubit.get_block_edges(block)) # Import all topo edges id as a Set
518548
toWritetoFile = [] #[""]*(len(edges_all)+1)
519549
toWritetoFile.append('%10i\n' % len(edges_all)) # Print the number of edges on the free surface
520550
for block,flag in zip(self.block_mat,self.block_flag): # For each 2D block
521551
print(block,flag)
522552
quads = cubit.get_block_faces(block) # Import quads id
523553
for quad in quads: # For each quad
524-
edges = Set(cubit.get_sub_elements("face", quad, 1)) # Get the lower dimension entities associated with a higher dimension entities.
554+
edges = set(cubit.get_sub_elements("face", quad, 1)) # Get the lower dimension entities associated with a higher dimension entities.
525555
# Here it gets the 1D edges associates with the face of id "quad". Store it as a Set
526556
intersection = edges & edges_all # Contains the edges of the considered quad that is on the free surface
527557
if len(intersection) != 0: # If this quad touch the free surface
@@ -550,19 +580,18 @@ def free_write(self,freename): #freename = None):
550580
cubit.cmd('set echo on') # Turn on echo of Cubit commands
551581
def forcing_write(self,forcname):
552582
""" Write forcing surfaces on file : forcname """
583+
forceedge = open(forcname,'w')
584+
print('Writing '+forcname+'.....')
553585
cubit.cmd('set info off') # Turn off return messages from Cubit commands
554586
cubit.cmd('set echo off') # Turn off echo of Cubit commands
555587
cubit.cmd('set journal off') # Do not save journal file
556-
from sets import Set
557-
forceedge = open(forcname,'w')
558-
print('Writing '+forcname+'.....')
559-
edges_forc = [Set()]*self.nforc # edges_forc[0] will be a Set containing the nodes describing the forcing boundary
588+
edges_forc = [set()]*self.nforc # edges_forc[0] will be a Set containing the nodes describing the forcing boundary
560589
# (index 0 : bottom, index 1 : right, index 2 : top, index 3 : left)
561590
nedges_all = 0 # To count the total number of forcing edges
562591
for block,flag in zip(self.block_bc,self.block_bc_flag): # For each 1D block
563592
for iforc in range(0, self.nforc): # iforc = 0,1,2,3 : for each forcing boundaries
564593
if block == self.forcing_boun[iforc]: # If the block considered correspond to the boundary
565-
edges_forc[iforc] = Set(cubit.get_block_edges(block)) # Store each edge on edges_forc
594+
edges_forc[iforc] = set(cubit.get_block_edges(block)) # Store each edge on edges_forc
566595
nedges_all = nedges_all+len(edges_forc[iforc]) # add the number of edges to nedges_all
567596
toWritetoFile = [""]*(nedges_all+1)
568597
toWritetoFile[0] = '%10i\n' % nedges_all # Write the total number of forcing edges to the first line of file
@@ -574,7 +603,7 @@ def forcing_write(self,forcname):
574603
quads = cubit.get_block_faces(block) # Import quads id
575604
for quad in quads: # For each quad
576605
#id_element = id_element+1 # id of this quad
577-
edges = Set(cubit.get_sub_elements("face", quad, 1)) # Get the lower dimension entities associated with a higher dimension entities.
606+
edges = set(cubit.get_sub_elements("face", quad, 1)) # Get the lower dimension entities associated with a higher dimension entities.
578607
# Here it gets the 1D edges associates with the face of id "quad". Store it as a Set
579608
for iforc in range(0,self.nforc): # iforc = 0,1,2,3 : for each forcing boundaries
580609
intersection = edges & edges_forc[iforc] # Contains the edges of the considered quad that is on the forcing boundary considered
@@ -605,20 +634,19 @@ def forcing_write(self,forcname):
605634
cubit.cmd('set echo on') # Turn on echo of Cubit commands
606635
def abs_write(self,absname):
607636
""" Write absorbing surfaces on file : absname """
608-
cubit.cmd('set info off') # Turn off return messages from Cubit commands
609-
cubit.cmd('set echo off') # Turn off echo of Cubit commands
610-
cubit.cmd('set journal off') # Do not save journal file.
611-
from sets import Set
612637
# if not absname: absname = self.absname
613638
absedge = open(absname,'w')
614639
print('Writing '+absname+'.....')
615-
edges_abs = [Set()]*self.nabs # edges_abs[0] will be a Set containing the nodes describing bottom adsorbing boundary
640+
cubit.cmd('set info off') # Turn off return messages from Cubit commands
641+
cubit.cmd('set echo off') # Turn off echo of Cubit commands
642+
cubit.cmd('set journal off') # Do not save journal file.
643+
edges_abs = [set()]*self.nabs # edges_abs[0] will be a Set containing the nodes describing bottom adsorbing boundary
616644
# (index 0 : bottom, index 1 : right, index 2 : top, index 3 : left)
617645
nedges_all = 0 # To count the total number of absorbing edges
618646
for block,flag in zip(self.block_bc,self.block_bc_flag): # For each 1D block
619647
for iabs in range(0, self.nabs): # iabs = 0,1,2,3 : for each absorbing boundaries
620648
if block == self.abs_boun[iabs]: # If the block considered correspond to the boundary
621-
edges_abs[iabs] = Set(cubit.get_block_edges(block)) # Store each edge on edges_abs
649+
edges_abs[iabs] = set(cubit.get_block_edges(block)) # Store each edge on edges_abs
622650
nedges_all = nedges_all+len(edges_abs[iabs]); # add the number of edges to nedges_all
623651
toWritetoFile = [""]*(nedges_all+1)
624652
toWritetoFile[0] = '%10i\n' % nedges_all # Write the total number of absorbing edges to the first line of file
@@ -630,7 +658,7 @@ def abs_write(self,absname):
630658
quads = cubit.get_block_faces(block) # Import quads id
631659
for quad in quads: # For each quad
632660
#id_element = id_element+1 # id of this quad
633-
edges = Set(cubit.get_sub_elements("face", quad, 1)) # Get the lower dimension entities associated with a higher dimension entities.
661+
edges = set(cubit.get_sub_elements("face", quad, 1)) # Get the lower dimension entities associated with a higher dimension entities.
634662
# Here it gets the 1D edges associates with the face of id "quad". Store it as a Set
635663
for iabs in range(0,self.nabs): # iabs = 0,1,2,3 : for each absorbing boundaries
636664
intersection = edges & edges_abs[iabs] # Contains the edges of the considered quad that is on the absorbing boundary considered
@@ -656,15 +684,14 @@ def abs_write(self,absname):
656684
cubit.cmd('set echo on') # Turn on echo of Cubit commands
657685
def axis_write(self,axis_name):
658686
""" Write axis on file """
687+
axisedge = open(axis_name,'w')
688+
print('Writing '+axis_name+'.....')
659689
cubit.cmd('set info off') # Turn off return messages from Cubit commands
660690
cubit.cmd('set echo off') # Turn off echo of Cubit commands
661691
cubit.cmd('set journal off') # Do not save journal file
662-
from sets import Set
663-
axisedge = open(axis_name,'w')
664-
print('Writing '+axis_name+'.....')
665692
for block,flag in zip(self.block_bc,self.block_bc_flag): # For each 1D block
666693
if block == self.axisId: # If the block correspond to the axis
667-
edges_all = Set(cubit.get_block_edges(block)) # Import all axis edges id as a Set
694+
edges_all = set(cubit.get_block_edges(block)) # Import all axis edges id as a Set
668695
toWritetoFile = [""]*(len(edges_all)+1)
669696
toWritetoFile[0] = '%10i\n' % len(edges_all) # Write the number of edges on the axis
670697
#axisedge.write('%10i\n' % len(edges_all)) # Write the number of edges on the axis
@@ -675,7 +702,7 @@ def axis_write(self,axis_name):
675702
quads = cubit.get_block_faces(block) # Import quads id
676703
for quad in quads: # For each quad
677704
#id_element = id_element+1 # id of this quad
678-
edges = Set(cubit.get_sub_elements("face", quad, 1)) # Get the lower dimension entities associated with a higher dimension entities.
705+
edges = set(cubit.get_sub_elements("face", quad, 1)) # Get the lower dimension entities associated with a higher dimension entities.
679706
# Here it gets the 1D edges associates with the face of id "quad". Store it as a Set
680707
intersection = edges & edges_all # Contains the edges of the considered quad that are on the axis
681708
if len(intersection) != 0: # If this quad touch the axis
@@ -710,11 +737,11 @@ def rec_write(self,recname):
710737
print('Ok')
711738
def write(self,path = ''):
712739
""" Write mesh in specfem2d format """
713-
print('Writing '+self.recname+'.....')
740+
print('Writing '+path+'.....')
714741
import os
715-
cubit.cmd('set info off') # Turn off return messages from Cubit commands
716-
cubit.cmd('set echo off') # Turn off echo of Cubit commands
717-
cubit.cmd('set journal off') # Do not save journal file
742+
#cubit.cmd('set info off') # Turn off return messages from Cubit commands
743+
#cubit.cmd('set echo off') # Turn off echo of Cubit commands
744+
#cubit.cmd('set journal off') # Do not save journal file
718745
if len(path) != 0: # If a path is supplied add a / at the end if needed
719746
if path[-1] != '/': path = path+'/'
720747
else:
@@ -736,8 +763,8 @@ def write(self,path = ''):
736763
if self.receivers:
737764
self.rec_write(path+self.recname) # If receivers has been set (as nodeset) write receiver file as well
738765
print('Mesh files has been writen in '+path)
739-
cubit.cmd('set info on') # Turn on return messages from Cubit commands
740-
cubit.cmd('set echo on') # Turn on echo of Cubit commands
766+
#cubit.cmd('set info on') # Turn on return messages from Cubit commands
767+
#cubit.cmd('set echo on') # Turn on echo of Cubit commands
741768

742769
def export2SPECFEM2D(path_exporting_mesh_SPECFEM2D='.'):
743770
# reads in mesh from cubit
@@ -747,7 +774,7 @@ def export2SPECFEM2D(path_exporting_mesh_SPECFEM2D='.'):
747774
print("END SPECFEM2D exporting process......")
748775

749776
# gets executed when run as script within cubit
750-
export2SPECFEM2D()
777+
#export2SPECFEM2D()
751778

752779

753780

0 commit comments

Comments
 (0)