adopting to Python standards

This commit is contained in:
Martin Diehl 2022-11-28 09:21:47 +01:00
parent 20f984fc87
commit 03a81b81e8
2 changed files with 195 additions and 195 deletions

View File

@ -11,195 +11,195 @@ import numpy as np
import damask import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0] script_name = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version]) script_id = ' '.join([script_name,damask.version])
# Convert .mfd file into a usable format # Convert .mfd file into a usable format
# Broken into labeled sections (eg. nodes, links, etc) # Broken into labeled sections (eg. nodes, links, etc)
# Each section has a list of labeled elements with formatted numerical data # Each section has a list of labeled elements with formatted numerical data
def parseMFD(dat): def parseMFD(dat):
formatted = [] formatted = []
section = 0 section = 0
formatted.append({'label': 'header', 'uid': -1, 'els': []}) formatted.append({'label': 'header', 'uid': -1, 'els': []})
# in between =beg= and =end= part of file # in between =beg= and =end= part of file
in_block = False in_block = False
for line in dat: for line in dat:
if in_block: # currently in a section if in_block: # currently in a section
# lines that start with a space are numerical data # lines that start with a space are numerical data
if line[0] == ' ': if line[0] == ' ':
formatted[section]['els'].append([]) formatted[section]['els'].append([])
# grab numbers # grab numbers
nums = re.split(r'\s+', line.strip()) nums = re.split(r'\s+', line.strip())
for num in nums: for num in nums:
# floating point has format ' -x.xxxxxxxxxxxxe+yy' # floating point has format ' -x.xxxxxxxxxxxxe+yy'
# scientific notation is used for float # scientific notation is used for float
if (len(num) >= 4) and (num[-4] == 'e'): if (len(num) >= 4) and (num[-4] == 'e'):
formatted[section]['els'][-1].append(float(num)) formatted[section]['els'][-1].append(float(num))
else: # integer else: # integer
formatted[section]['els'][-1].append(int(num)) formatted[section]['els'][-1].append(int(num))
else: # not numerical data, so it is a label for an element or section end else: # not numerical data, so it is a label for an element or section end
if line[0] == '=' and re.search(r'=end=$', line) is not None: # End of section, avoiding regex if possible if line[0] == '=' and re.search(r'=end=$', line) is not None: # End of section, avoiding regex if possible
in_block = False in_block = False
else: else:
formatted[section]['els'].append([]) formatted[section]['els'].append([])
formatted[section]['els'][-1] = line formatted[section]['els'][-1] = line
else: # Not in a section, we are looking for a =beg= now else: # Not in a section, we are looking for a =beg= now
search = re.search(r'=beg=\s+(\d+)\s\((.*?)\)', line) search = re.search(r'=beg=\s+(\d+)\s\((.*?)\)', line)
if search is not None: # found start of a new section if search is not None: # found start of a new section
section += 1 section += 1
in_block = True in_block = True
formatted.append({'label': search.group(2), 'uid': int(search.group(1)), 'els': []}) formatted.append({'label': search.group(2), 'uid': int(search.group(1)), 'els': []})
else: # No =beg= found, probably in the header else: # No =beg= found, probably in the header
# Either header or somthing we didn't plan for - just save the line so it isn't lost # Either header or somthing we didn't plan for - just save the line so it isn't lost
if formatted[section]['uid'] > 0: if formatted[section]['uid'] > 0:
section += 1 section += 1
formatted.append({'label': '', 'uid': -2, 'els': []}) # make dummy section to store unrecognized data formatted.append({'label': '', 'uid': -2, 'els': []}) # make dummy section to store unrecognized data
formatted[section]['els'].append(line) formatted[section]['els'].append(line)
return formatted return formatted
def asMFD(mfd_data): def asMFD(mfd_data):
result = '' result = ''
for section in mfd_data: for section in mfd_data:
if section['uid'] > 0: if section['uid'] > 0:
result += '=beg={0:5d} ({1})\n'.format(section['uid'], section['label']) result += '=beg={0:5d} ({1})\n'.format(section['uid'], section['label'])
for el in section['els']: for el in section['els']:
if type(el) == str: if type(el) == str:
result += el result += el
elif type(el) == list: elif type(el) == list:
for num in el: for num in el:
if type(num) == int: if type(num) == int:
result += '{:20d}'.format(num) result += '{:20d}'.format(num)
elif type(num) == float: elif type(num) == float:
result += '{:20.12e}'.format(num) result += '{:20.12e}'.format(num)
else: else:
print(f'WARNING: encountered unknown type: {type(el)}') print(f'WARNING: encountered unknown type: {type(el)}')
result += '\n' result += '\n'
else: else:
print(f'WARNING: encountered unknown type: {type(el)}') print(f'WARNING: encountered unknown type: {type(el)}')
if section['uid'] > 0: if section['uid'] > 0:
result += '=end=\n' result += '=end=\n'
return result.strip() return result.strip()
def add_servoLinks(mfd_data,active=[True,True,True]): # directions on which to add PBC def add_servoLinks(mfd_data,active=[True,True,True]): # directions on which to add PBC
base = ['x','y','z'] base = ['x','y','z']
box = {'min': np.zeros(3,dtype='d'), box = {'min': np.zeros(3,dtype='d'),
'max': np.zeros(3,dtype='d'), 'max': np.zeros(3,dtype='d'),
'delta': np.zeros(3,dtype='d'), 'delta': np.zeros(3,dtype='d'),
} }
mfd_dict = {} mfd_dict = {}
for i in range(len(mfd_data)): for i in range(len(mfd_data)):
mfd_dict[mfd_data[i]['label']] = i mfd_dict[mfd_data[i]['label']] = i
NodeCoords = np.array(mfd_data[mfd_dict['nodes']]['els'][1::4])[:,1:4] NodeCoords = np.array(mfd_data[mfd_dict['nodes']]['els'][1::4])[:,1:4]
Nnodes = NodeCoords.shape[0] Nnodes = NodeCoords.shape[0]
box['min'] = NodeCoords.min(axis=0) # find the bounding box box['min'] = NodeCoords.min(axis=0) # find the bounding box
box['max'] = NodeCoords.max(axis=0) box['max'] = NodeCoords.max(axis=0)
box['delta'] = box['max']-box['min'] box['delta'] = box['max']-box['min']
for coord in range(3): # calc the dimension of the bounding box for coord in range(3): # calc the dimension of the bounding box
if box['delta'][coord] != 0.0: if box['delta'][coord] != 0.0:
for extremum in ['min','max']: for extremum in ['min','max']:
rounded = round(box[extremum][coord]*1e+15/box['delta'][coord]) * \ rounded = round(box[extremum][coord]*1e+15/box['delta'][coord]) * \
1e-15*box['delta'][coord] # rounding to 1e-15 of dimension 1e-15*box['delta'][coord] # rounding to 1e-15 of dimension
box[extremum][coord] = 0.0 if rounded == 0.0 else rounded # get rid of -0.0 (negative zeros) box[extremum][coord] = 0.0 if rounded == 0.0 else rounded # get rid of -0.0 (negative zeros)
baseNode = {} baseNode = {}
linkNodes = [] linkNodes = []
#------------------------------------------------------------------------------------------------- #-------------------------------------------------------------------------------------------------
# loop over all nodes # loop over all nodes
for node in range(Nnodes): for node in range(Nnodes):
key = {} key = {}
maxFlag = [False, False, False] maxFlag = [False, False, False]
Nmax = 0 Nmax = 0
Nmin = 0 Nmin = 0
for coord in range(3): # for each direction for coord in range(3): # for each direction
if box['delta'][coord] != 0.0: if box['delta'][coord] != 0.0:
rounded = round(NodeCoords[node,coord]*1e+15/box['delta'][coord]) * \ rounded = round(NodeCoords[node,coord]*1e+15/box['delta'][coord]) * \
1e-15*box['delta'][coord] # rounding to 1e-15 of dimension 1e-15*box['delta'][coord] # rounding to 1e-15 of dimension
NodeCoords[node,coord] = 0.0 if rounded == 0.0 else rounded # get rid of -0.0 (negative zeros) NodeCoords[node,coord] = 0.0 if rounded == 0.0 else rounded # get rid of -0.0 (negative zeros)
key[base[coord]] = "%.8e"%NodeCoords[node,coord] # translate position to string key[base[coord]] = "%.8e"%NodeCoords[node,coord] # translate position to string
if (key[base[coord]] == "%.8e"%box['min'][coord]): # compare to min of bounding box (i.e. is on outer face?) if (key[base[coord]] == "%.8e"%box['min'][coord]): # compare to min of bounding box (i.e. is on outer face?)
Nmin += 1 # count outer (back) face membership Nmin += 1 # count outer (back) face membership
elif (key[base[coord]] == "%.8e"%box['max'][coord]): # compare to max of bounding box (i.e. is on outer face?) elif (key[base[coord]] == "%.8e"%box['max'][coord]): # compare to max of bounding box (i.e. is on outer face?)
Nmax += 1 # count outer (front) face membership Nmax += 1 # count outer (front) face membership
maxFlag[coord] = True # remember face membership (for linked nodes) maxFlag[coord] = True # remember face membership (for linked nodes)
if Nmin > 0: # node is on a back face if Nmin > 0: # node is on a back face
# prepare for any non-existing entries in the data structure # prepare for any non-existing entries in the data structure
if key['x'] not in baseNode.keys(): if key['x'] not in baseNode.keys():
baseNode[key['x']] = {} baseNode[key['x']] = {}
if key['y'] not in baseNode[key['x']].keys(): if key['y'] not in baseNode[key['x']].keys():
baseNode[key['x']][key['y']] = {} baseNode[key['x']][key['y']] = {}
if key['z'] not in baseNode[key['x']][key['y']].keys(): if key['z'] not in baseNode[key['x']][key['y']].keys():
baseNode[key['x']][key['y']][key['z']] = 0 baseNode[key['x']][key['y']][key['z']] = 0
baseNode[key['x']][key['y']][key['z']] = node+1 # remember the base node id baseNode[key['x']][key['y']][key['z']] = node+1 # remember the base node id
if Nmax > 0 and Nmax >= Nmin: # node is on at least as many front than back faces if Nmax > 0 and Nmax >= Nmin: # node is on at least as many front than back faces
if any([maxFlag[i] and active[i] for i in range(3)]): if any([maxFlag[i] and active[i] for i in range(3)]):
linkNodes.append({'id': node+1,'coord': NodeCoords[node], 'faceMember': [maxFlag[i] and active[i] for i in range(3)]}) linkNodes.append({'id':node+1,'coord':NodeCoords[node],'faceMember':[maxFlag[i] and active[i] for i in range(3)]})
mfd_data[mfd_dict['entities']]['els'][0][0] += len(linkNodes) * 3 mfd_data[mfd_dict['entities']]['els'][0][0] += len(linkNodes) * 3
baseCorner = baseNode["%.8e"%box['min'][0]]["%.8e"%box['min'][1]]["%.8e"%box['min'][2]] # detect ultimate base node baseCorner = baseNode["%.8e"%box['min'][0]]["%.8e"%box['min'][1]]["%.8e"%box['min'][2]] # detect ultimate base node
links = {'uid': 1705, 'label': 'links', 'els': [[7,0],[9,0]]} links = {'uid': 1705, 'label': 'links', 'els': [[7,0],[9,0]]}
linkID = 0 linkID = 0
for node in linkNodes: # loop over all linked nodes for node in linkNodes: # loop over all linked nodes
linkCoord = [node['coord']] # start list of control node coords with my coords linkCoord = [node['coord']] # start list of control node coords with my coords
for dir in range(3): # check for each direction for dir in range(3): # check for each direction
if node['faceMember'][dir]: # me on this front face if node['faceMember'][dir]: # me on this front face
linkCoord[0][dir] = box['min'][dir] # project me onto rear face along dir linkCoord[0][dir] = box['min'][dir] # project me onto rear face along dir
linkCoord.append(np.array(box['min'])) # append base corner linkCoord.append(np.array(box['min'])) # append base corner
linkCoord[-1][dir] = box['max'][dir] # stretch it to corresponding control leg of "dir" linkCoord[-1][dir] = box['max'][dir] # stretch it to corresponding control leg of "dir"
nLinks = len(linkCoord) nLinks = len(linkCoord)
for dof in [1,2,3]: for dof in [1,2,3]:
tied_node = node['id'] tied_node = node['id']
nterms = 1 + nLinks nterms = 1 + nLinks
linkID += 1 linkID += 1
# Link header # Link header
links['els'].append('link{0}\n'.format(linkID)) links['els'].append('link{0}\n'.format(linkID))
links['els'].append([linkID, 1]) links['els'].append([linkID, 1])
links['els'].append([0]) links['els'].append([0])
links['els'].append([0]) links['els'].append([0])
links['els'].append([0, 0, 0, tied_node]) links['els'].append([0, 0, 0, tied_node])
# these need to be put in groups of four # these need to be put in groups of four
link_payload = [dof, 0, nterms] link_payload = [dof, 0, nterms]
# Individual node contributions (node, dof, coef.) # Individual node contributions (node, dof, coef.)
for i in range(nterms): for i in range(nterms):
if i == nLinks: if i == nLinks:
link_payload.append(baseCorner) link_payload.append(baseCorner)
else: else:
link_payload.append(baseNode["%.8e"%linkCoord[i][0]]["%.8e"%linkCoord[i][1]]["%.8e"%linkCoord[i][2]]) link_payload.append(baseNode["%.8e"%linkCoord[i][0]]["%.8e"%linkCoord[i][1]]["%.8e"%linkCoord[i][2]])
for i in range(nterms): for i in range(nterms):
link_payload.append(dof) link_payload.append(dof)
for i in range(nterms): for i in range(nterms):
if i == nLinks: if i == nLinks:
link_payload.append(1.0 - nLinks) link_payload.append(1.0 - nLinks)
else: else:
link_payload.append(1.0) link_payload.append(1.0)
# Needs to be formatted 4 data points per row, character width of 20, so 80 total # Needs to be formatted 4 data points per row, character width of 20, so 80 total
for j in range(0, len(link_payload), 4): for j in range(0, len(link_payload), 4):
links['els'].append(link_payload[j:j+4]) links['els'].append(link_payload[j:j+4])
if j+4 < len(link_payload): if j+4 < len(link_payload):
links['els'].append(link_payload[j+4:]) links['els'].append(link_payload[j+4:])
i = 0 i = 0
while i < len(mfd_data) and mfd_data[i]['uid'] < 1705: i += 1 while i < len(mfd_data) and mfd_data[i]['uid'] < 1705: i += 1
if mfd_data[i]['uid'] == 1705: del mfd_data[i] if mfd_data[i]['uid'] == 1705: del mfd_data[i]
mfd_data.insert(i, links) mfd_data.insert(i, links)
#-------------------------------------------------------------------------------------------------- #--------------------------------------------------------------------------------------------------
@ -209,54 +209,54 @@ def add_servoLinks(mfd_data,active=[True,True,True]): # directions on which to
parser = OptionParser(usage='%prog options [file[s]]', description = """ parser = OptionParser(usage='%prog options [file[s]]', description = """
Set up servo linking to achieve periodic boundary conditions for a regular hexahedral mesh. Set up servo linking to achieve periodic boundary conditions for a regular hexahedral mesh.
Use *py_connection to operate on model presently opened in MSC.Mentat. Use *py_connection to operate on model presently opened in MSC.Mentat.
""", version = scriptID) """, version = script_id)
parser.add_option('-p', '--port', parser.add_option('-p', '--port',
type = int, metavar = 'int', default = None, type = int, metavar = 'int', default = None,
help = 'Mentat connection port') help = 'Mentat connection port')
parser.add_option('-x', parser.add_option('-x',
action = 'store_false', default = True, action = 'store_false', default = True,
help = 'no PBC along x direction') help = 'no PBC along x direction')
parser.add_option('-y', parser.add_option('-y',
action = 'store_false', default = True, action = 'store_false', default = True,
help = 'no PBC along y direction') help = 'no PBC along y direction')
parser.add_option('-z', parser.add_option('-z',
action = 'store_false', default = True, action = 'store_false', default = True,
help = 'no PBC along z direction') help = 'no PBC along z direction')
(options, filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
remote = options.port is not None remote = options.port is not None
if remote and filenames != []: if remote and filenames != []:
parser.error('file can not be specified when port is given.') parser.error('file can not be specified when port is given.')
if filenames == []: filenames = [None] if filenames == []: filenames = [None]
if remote: if remote:
sys.path.append(str(damask.solver.Marc().library_path)) sys.path.append(str(damask.solver.Marc().library_path))
import py_mentat import py_mentat
print(scriptName+': waiting to connect...') print(script_name+': waiting to connect...')
filenames = [os.path.join(tempfile._get_default_tempdir(), next(tempfile._get_candidate_names()) + '.mfd')] filenames = [os.path.join(tempfile._get_default_tempdir(), next(tempfile._get_candidate_names()) + '.mfd')]
try: try:
py_mentat.py_connect('',options.port) py_mentat.py_connect('',options.port)
py_mentat.py_send('*set_save_formatted on') py_mentat.py_send('*set_save_formatted on')
py_mentat.py_send('*save_as_model "{}" yes'.format(filenames[0])) py_mentat.py_send('*save_as_model "{}" yes'.format(filenames[0]))
py_mentat.py_get_int("nnodes()") py_mentat.py_get_int("nnodes()")
except py_mentat.InputError as err: except py_mentat.InputError as err:
print(f'{err}. Try Tools/Python/"Run as Separate Process" & "Initiate".') print(f'{err}. Try Tools/Python/"Run as Separate Process" & "Initiate".')
sys.exit(-1) sys.exit(-1)
print( 'connected...') print( 'connected...')
for name in filenames: for name in filenames:
while remote and not os.path.exists(name): time.sleep(0.5) while remote and not os.path.exists(name): time.sleep(0.5)
with open( name,'r') if name is not None else sys.stdin as fileIn: with open( name,'r') if name is not None else sys.stdin as file_in:
print(scriptName+': '+name) print(script_name+': '+name)
mfd = parseMFD(fileIn) mfd = parseMFD(file_in)
add_servoLinks(mfd,[options.x,options.y,options.z]) add_servoLinks(mfd,[options.x,options.y,options.z])
with open( name,'w') if name is not None else sys.stdout as fileOut: with open( name,'w') if name is not None else sys.stdout as file_out:
fileOut.write(asMFD(mfd)) file_out.write(asMFD(mfd))
if remote: if remote:
py_mentat.py_send('*open_model "{}"'.format(filenames[0])) py_mentat.py_send('*open_model "{}"'.format(filenames[0]))

View File

@ -7,8 +7,8 @@ from optparse import OptionParser
import damask import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0] script_name = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version]) script_id = ' '.join([script_name,damask.version])
#------------------------------------------------------------------------------------------------- #-------------------------------------------------------------------------------------------------
def outMentat(cmd,locals): def outMentat(cmd,locals):
@ -45,7 +45,7 @@ def output(cmds,locals,dest):
#------------------------------------------------------------------------------------------------- #-------------------------------------------------------------------------------------------------
def init(): def init():
return [ return [
"|"+' '.join([scriptID] + sys.argv[1:]), "|"+' '.join([script_id] + sys.argv[1:]),
"*draw_manual", # prevent redrawing in Mentat, should be much faster "*draw_manual", # prevent redrawing in Mentat, should be much faster
"*new_model yes", "*new_model yes",
"*reset", "*reset",
@ -170,7 +170,7 @@ def initial_conditions(material):
parser = OptionParser(usage='%prog options [file[s]]', description = """ parser = OptionParser(usage='%prog options [file[s]]', description = """
Generate MSC.Marc FE hexahedral mesh from geom file. Generate MSC.Marc FE hexahedral mesh from geom file.
""", version = scriptID) """, version = script_id)
parser.add_option('-p', '--port', parser.add_option('-p', '--port',
dest = 'port', dest = 'port',
@ -194,7 +194,7 @@ if options.port is not None:
if filenames == []: filenames = [None] if filenames == []: filenames = [None]
for name in filenames: for name in filenames:
print(scriptName+': '+name) print(script_name+': '+name)
geom = damask.Grid.load(StringIO(''.join(sys.stdin.read())) if name is None else name) geom = damask.Grid.load(StringIO(''.join(sys.stdin.read())) if name is None else name)
material = geom.material.flatten(order='F') material = geom.material.flatten(order='F')
@ -211,11 +211,11 @@ for name in filenames:
'*draw_automatic', '*draw_automatic',
] ]
outputLocals = {} output_locals = {}
if options.port: if options.port:
py_mentat.py_connect('',options.port) py_mentat.py_connect('',options.port)
output(cmds,outputLocals,'Mentat') output(cmds,output_locals,'Mentat')
py_mentat.py_disconnect() py_mentat.py_disconnect()
else: else:
with sys.stdout if name is None else open(os.path.splitext(name)[0]+'.proc','w') as f: with sys.stdout if name is None else open(os.path.splitext(name)[0]+'.proc','w') as f:
output(cmds,outputLocals,f) output(cmds,output_locals,f)