Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

code cleaning & cubit script updates #1711

Merged
merged 4 commits into from
Jul 8, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
96 changes: 46 additions & 50 deletions CUBIT_GEOCUBIT/geocubitlib/cubit2specfem3d.py
Original file line number Diff line number Diff line change
Expand Up @@ -652,14 +652,16 @@ def block_definition(self):
self.bc = bc
print('HEX Blocks:')
for m, f in zip(self.block_mat, self.block_flag):
print('block ', m, 'material flag ', f)
print(' block ', m, 'material flag ', f)
print('Absorbing Boundary Conditions:')
for m, f in zip(self.block_bc, self.block_bc_flag):
print('bc ', m, 'bc flag ', f)
print('Topography (free surface)')
print(self.topography)
print('Free surface')
print(self.free)
print(' bc ', m, 'bc flag ', f)
if self.topography:
print('Topography (free surface):')
print(' ',self.topography)
if self.free:
print('Free surface:')
print(' ',self.free)
except:
print('****************************************')
print('sorry, no blocks or blocks not properly defined')
Expand Down Expand Up @@ -765,7 +767,7 @@ def mat_parameter(self, properties):
else:
raise RuntimeError('Error: material id must be strictly positive or negative, but not equal to 0')
# info output
print("material: ",txt)
print(" material: ",txt)
return txt

def nummaterial_write(self, nummaterial_name, placeholder=True):
Expand Down Expand Up @@ -853,8 +855,7 @@ def create_hexnode_string(self, hexa, hexnode_string=True):
else:
map(int, txt.split())

def create_facenode_string(self, hexa, face, normal=None, cknormal=True,
facenode_string=True):
def create_facenode_string(self, hexa, face, normal=None, cknormal=True, facenode_string=True):
nodes = self.get_face_connectivity(face)
if cknormal:
nodes_ok = self.normal_check(nodes[0:4], normal)
Expand Down Expand Up @@ -886,7 +887,7 @@ def mesh_write(self, mesh_name):
meshfile.write(str(num_elems) + '\n')
for block, flag in zip(self.block_mat, self.block_flag):
hexes = cubit.get_block_hexes(block)
print('block ', block, ' hexes ', len(hexes))
print(' block ', block, ' hexes ', len(hexes))
for hexa in hexes:
txt = self.create_hexnode_string(hexa)
meshfile.write(txt)
Expand All @@ -897,7 +898,7 @@ def material_write(self, mat_name):
mat = open(mat_name, 'w')
print('Writing ' + mat_name + '.....')
for block, flag in zip(self.block_mat, self.block_flag):
print('block ', block, 'flag ', flag)
print(' block ', block, 'flag ', flag)
hexes = cubit.get_block_hexes(block)
for hexa in hexes:
mat.write(('%10i %10i\n') % (hexa, flag))
Expand Down Expand Up @@ -949,8 +950,7 @@ def free_write(self, freename=None):
for block, flag in zip(self.block_bc, self.block_bc_flag):
if block == self.topography:
name = cubit.get_exodus_entity_name('block', block)
print('free surface (topography) block name:', \
name, 'id:', block)
print(' free surface (topography): block name:', name, 'id:', block)
quads_all = cubit.get_block_faces(block)
print(' number of faces = ', len(quads_all))
dic_quads_all = dict(zip(quads_all, quads_all))
Expand All @@ -968,13 +968,13 @@ def free_write(self, freename=None):
faces = cubit.get_sub_elements('hex', h, 2)
for f in faces:
if f in s:
txt = self.create_facenode_string(
h, f, normal, cknormal=True)
# topography surface: checks face orientation against vertical normal (0,0,1)
txt = self.create_facenode_string(h, f, normal, cknormal=True)
freehex.write(txt)
print(' 100 %')
elif block == self.free:
name = cubit.get_exodus_entity_name('block', block)
print('free surface block name:', name, 'id:', block)
print(' free surface: block name:', name, 'id:', block)
quads_all = cubit.get_block_faces(block)
print(' number of faces = ', len(quads_all))
dic_quads_all = dict(zip(quads_all, quads_all))
Expand All @@ -992,8 +992,8 @@ def free_write(self, freename=None):
faces = cubit.get_sub_elements('hex', h, 2)
for f in faces:
if f in s:
txt = self.create_facenode_string(
h, f, normal, cknormal=False)
# free surface: no face orientation check against vertical normal
txt = self.create_facenode_string(h, f, normal, cknormal=False)
freehex.write(txt)
print(' 100 %')
freehex.close()
Expand Down Expand Up @@ -1122,60 +1122,60 @@ def abs_write(self, absname=None):
# loops through all block definitions
list_hex = cubit.parse_cubit_list('hex', 'all')
for block, flag in zip(self.block_bc, self.block_bc_flag):
if block != self.topography:
if block != self.topography and block != self.free:
name = cubit.get_exodus_entity_name('block', block)
print(' block name:', name, 'id:', block)
cknormal = True
abshex_local = False
# opens file
if re.search('xmin', name):
print("xmin")
print(" xmin")
abshex_local = open(absname + '_xmin', 'w')
normal = (-1, 0, 0)
elif re.search('xmax', name):
print("xmax")
print(" xmax")
abshex_local = open(absname + '_xmax', 'w')
normal = (1, 0, 0)
elif re.search('ymin', name):
print("ymin")
print(" ymin")
abshex_local = open(absname + '_ymin', 'w')
normal = (0, -1, 0)
elif re.search('ymax', name):
print("ymax")
print(" ymax")
abshex_local = open(absname + '_ymax', 'w')
normal = (0, 1, 0)
elif re.search('bottom', name):
print("bottom")
print(" bottom")
abshex_local = open(absname + '_bottom', 'w')
normal = (0, 0, -1)
elif re.search('abs', name):
print("abs all - experimental, check the output")
print(" abs all - experimental, check the output")
cknormal = False
abshex_local = open(absname, 'w')
normal = None
else:
if block == 1003:
print("xmin")
print(" xmin")
abshex_local = open(absname + '_xmin', 'w')
normal = (-1, 0, 0)
elif block == 1004:
print("ymin")
print(" ymin")
abshex_local = open(absname + '_ymin', 'w')
normal = (0, -1, 0)
elif block == 1005:
print("xmax")
print(" xmax")
abshex_local = open(absname + '_xmax', 'w')
normal = (1, 0, 0)
elif block == 1006:
print("ymax")
print(" ymax")
abshex_local = open(absname + '_ymax', 'w')
normal = (0, 1, 0)
elif block == 1002:
print("bottom")
print(" bottom")
abshex_local = open(absname + '_bottom', 'w')
normal = (0, 0, -1)
elif block == 1000:
print("custom")
print(" custom")
abshex_local = open(absname, 'w')
cknormal = False
normal = None
Expand All @@ -1199,8 +1199,7 @@ def abs_write(self, absname=None):
faces = cubit.get_sub_elements('hex', h, 2)
for f in faces:
if f in s:
txt = self.create_facenode_string(
h, f, normal=normal, cknormal=cknormal)
txt = self.create_facenode_string(h, f, normal=normal, cknormal=cknormal)
abshex_local.write(txt)
abshex_local.close()
print(' 100 %')
Expand All @@ -1215,19 +1214,23 @@ def surface_write(self, pathdir=None):
# > block 10 name 'moho_surface'
import re
for block in self.block_bc:
if block != self.topography:
if block != self.topography and block != self.free:
name = cubit.get_exodus_entity_name('block', block)
# skips block names like face_abs**, face_topo**
# skips block names like face_abs**, face_topo**, face_free**
if re.search('abs', name):
continue
elif re.search('topo', name):
continue
elif re.search('free', name):
continue
elif re.search('surface', name):
# e.g. moho_surface
filename = pathdir + name + '_file'
else:
continue
# gets face elements
print(' surface block name: ', name, 'id: ', block)
print('Optional surface:')
print(' surface: block name: ', name, 'id: ', block)
quads_all = cubit.get_block_faces(block)
print(' face = ', len(quads_all))
if len(quads_all) == 0:
Expand All @@ -1243,17 +1246,16 @@ def surface_write(self, pathdir=None):
Nhex = len(list_hex)
# Create a set to speed up "if f in dic_quads_all.keys()" big time
s = set(dic_quads_all.keys())
percentageDone = -1
for h in list_hex:
lastDisplayed = -1
for ih, h in enumerate(list_hex):
percentageDone = float(ih)/Nhex*100.0
if int(percentageDone)%10 == 0 and lastDisplayed != int(percentageDone):
print(' ', int(percentageDone), '%')
lastDisplayed = int(percentageDone)
faces = cubit.get_sub_elements('hex', h, 2)
for f in faces:
if f in s:
txt = self.create_facenode_string(
h, f, cknormal=False)
txt = self.create_facenode_string(h, f, cknormal=False)
surfhex_local.write(txt)
# closes file
surfhex_local.close()
Expand Down Expand Up @@ -1388,7 +1390,7 @@ def _get_bc_flag(self, block):
normal = (0, 0, -1)
elif re.search('abs', name):
label = "all"
print("abs all - experimental, check the output")
print(" abs all - experimental, check the output")
cknormal = False
else:
if block == 1003:
Expand Down Expand Up @@ -1424,9 +1426,7 @@ def _get_bc_flag(self, block):
# faces = cubit.get_sub_elements('hex', h, 2)
# for f in faces:
# if dic_quads_all.has_key(f):
# abs_array.append(self.create_facenode_string(
# h, f, normal=normal, cknormal=cknormal,
# facenode_string=False))
# abs_array.append(self.create_facenode_string(h, f, normal=normal, cknormal=cknormal, facenode_string=False))
# print('Ok')
# cubit.cd('set info on')
# cubit.cmd('set echo on')
Expand Down Expand Up @@ -1474,9 +1474,7 @@ def _get_bc_flag(self, block):
# for f in faces:
# if dic_quads_all.has_key(f):
# # print(f)
# free_array.append(self.create_facenode_string(
# h, f, normal, cknormal=True,
# facenode_string=False))
# free_array.append(self.create_facenode_string(h, f, normal, cknormal=True, facenode_string=False))
# elif block == self.free:
# name = cubit.get_exodus_entity_name('block', block)
# print('free surface block name:', name, 'id:', block)
Expand All @@ -1488,9 +1486,7 @@ def _get_bc_flag(self, block):
# faces = cubit.get_sub_elements('hex', h, 2)
# for f in faces:
# if dic_quads_all.has_key(f):
# free_array.append(self.create_facenode_string(
# h, f, normal, cknormal=False,
# facenode_string=False))
# free_array.append(self.create_facenode_string(h, f, normal, cknormal=False, facenode_string=False))
# print('Ok')
# cubit.cmd('set info on')
# cubit.cmd('set echo on')
Expand Down
2 changes: 1 addition & 1 deletion EXAMPLES/applications/Mount_StHelens/mesh_mount.py
Original file line number Diff line number Diff line change
Expand Up @@ -158,8 +158,8 @@
else:
print("#" + time.strftime("%a, %d %b %Y %H:%M:%S", time.localtime()))
print("# imprinting volume, this will take around 1 min, please be patience...")
cubit.cmd('merge all')
cubit.cmd('imprint all')
cubit.cmd('merge all')

# exports only surfaces which will create single volume
cubit.cmd('export acis "topo_2.acis" surface 3 10 12 14 15 9 ascii overwrite')
Expand Down
2 changes: 1 addition & 1 deletion EXAMPLES/applications/Mount_StHelens/mesh_mount_stl.py
Original file line number Diff line number Diff line change
Expand Up @@ -132,7 +132,7 @@
# imprints topography surfaces into brick volume, creates connected surfaces
# note: this is a develop feature
cubit.cmd('imprint all')

cubit.cmd('merge all')

cubit.cmd('version')

Expand Down
1 change: 1 addition & 0 deletions EXAMPLES/applications/fault_examples/tpv5/TPV5.py
Original file line number Diff line number Diff line change
Expand Up @@ -131,6 +131,7 @@

cubit.cmd("imprint all")
cubit.cmd("merge all")

cubit.cmd("surface 1 size "+str(elementsize))
cubit.cmd("volume 1 size "+str(elementsize))
cubit.cmd("surface 1 scheme pave")
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -94,8 +94,8 @@

cubit.cmd('delete volume 2 4')

cubit.cmd('merge all')
cubit.cmd('imprint all')
cubit.cmd('merge all')

cubit.cmd('volume 3 size '+str(elementsize))
cubit.cmd('mesh volume 3')
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -94,8 +94,8 @@

cubit.cmd('delete volume 2 4')

cubit.cmd('merge all')
cubit.cmd('imprint all')
cubit.cmd('merge all')

# Meshing the volumes
## middle volume
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -88,8 +88,8 @@

cubit.cmd('delete volume 2 4')

cubit.cmd('merge all')
cubit.cmd('imprint all')
cubit.cmd('merge all')

# Meshing the volumes
cubit.cmd('volume 3 size '+str(elementsize))
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -68,8 +68,8 @@

cubit.cmd('delete volume 2 4')

cubit.cmd('merge all')
cubit.cmd('imprint all')
cubit.cmd('merge all')

# Meshing the volumes
cubit.cmd('volume 3 size 3589.2')
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -64,8 +64,8 @@

cubit.cmd('delete volume 2 4')

cubit.cmd('merge all')
cubit.cmd('imprint all')
cubit.cmd('merge all')

# Meshing the volumes
#elementsize = 1196.4 #hi-resolution
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -95,8 +95,8 @@
cubit.cmd('delete volume 4')

# merges surfaces
cubit.cmd('merge all')
cubit.cmd('imprint all')
cubit.cmd('merge all')

# resets numbering
cubit.cmd('compress all')
Expand Down
4 changes: 2 additions & 2 deletions src/decompose_mesh/part_decompose_mesh.F90
Original file line number Diff line number Diff line change
Expand Up @@ -620,7 +620,7 @@ subroutine write_wavefield_discontinuity_database(iproc, outputpath_name, &
implicit none
integer, intent(in) :: iproc
integer, intent(in) :: nspec
integer, intent(in) :: nb_wd
integer, intent(in) :: nb_wd
integer, dimension(:), pointer :: glob2loc_elmnts
integer, dimension(1:nspec) :: part
integer, dimension(1:nb_wd), intent(in) :: boundary_to_ispec_wd, side_wd
Expand All @@ -632,7 +632,7 @@ subroutine write_wavefield_discontinuity_database(iproc, outputpath_name, &
character(len=MAX_STRING_LEN) :: prname
write(prname, "('proc',i6.6,'_')") iproc
open(unit=IFILE_WAVEFIELD_DISCONTINUITY, &
file=trim(outputpath_name)//'/'//trim(prname)//&
file = trim(outputpath_name)//'/'//trim(prname)//&
trim(FNAME_WAVEFIELD_DISCONTINUITY_MESH), &
form='unformatted', action='write')
local_nb_wd = 0
Expand Down
4 changes: 2 additions & 2 deletions src/decompose_mesh/read_mesh_files.F90
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ subroutine read_mesh_files()
use fault_scotch, only: ANY_FAULT,read_fault_files,save_nodes_coords,close_faults

use constants, only: FNAME_WAVEFIELD_DISCONTINUITY_INTERFACE

implicit none

! local parameters
Expand Down Expand Up @@ -592,7 +592,7 @@ subroutine read_mesh_files()
read(IIN_DB, *, iostat=ier) boundary_to_ispec_wd(ib), side_wd(ib)
enddo
close(IIN_DB)
endif
endif

! reads in absorbing boundary files
open(unit=IIN_DB, file=localpath_name(1:len_trim(localpath_name))//'/absorbing_surface_file_xmin', &
Expand Down
Loading
Loading