extract the deformation anisotropic coefficients (r values) and strength anisotropic coefficients (normalized yield stress) automatically.
This commit is contained in:
Haiming Zhang
parent
47249bb52a
commit
61e3987bfa
|
|
@ -0,0 +1,251 @@
|
||||||
|
#!/usr/bin/python
|
||||||
|
# -*- coding: UTF-8 no BOM -*-
|
||||||
|
|
||||||
|
import threading,time,os,subprocess,shlex,string
|
||||||
|
import numpy as np
|
||||||
|
from optparse import OptionParser, OptionGroup
|
||||||
|
from shutil import copy2
|
||||||
|
from re import split
|
||||||
|
import damask
|
||||||
|
|
||||||
|
scriptID = string.replace('$Id: yieldSurface.py 3942 2015-02-25 17:24:33Z MPIE\hm.zhang $','\n','\\n')
|
||||||
|
scriptName = scriptID.split()[1][:-3]
|
||||||
|
|
||||||
|
def list_split(option, opt, value, parser):
|
||||||
|
setattr(parser.values, option.dest, value.split(','))
|
||||||
|
|
||||||
|
def execute(cmd,streamIn=None,wd='./'):
|
||||||
|
'''
|
||||||
|
executes a command in given directory and returns stdout and stderr for optional stdin
|
||||||
|
'''
|
||||||
|
initialPath=os.getcwd()
|
||||||
|
os.chdir(wd)
|
||||||
|
process = subprocess.Popen(shlex.split(cmd),stdout=subprocess.PIPE,stderr = subprocess.PIPE,stdin=subprocess.PIPE)
|
||||||
|
if streamIn != None:
|
||||||
|
out,error = process.communicate(streamIn.read())
|
||||||
|
else:
|
||||||
|
out,error = process.communicate()
|
||||||
|
os.chdir(initialPath)
|
||||||
|
return out,error
|
||||||
|
|
||||||
|
#---------------------------------------------------------------------------------------------------
|
||||||
|
class myThread (threading.Thread):
|
||||||
|
#---------------------------------------------------------------------------------------------------
|
||||||
|
'''
|
||||||
|
Runner class
|
||||||
|
'''
|
||||||
|
def __init__(self, threadID):
|
||||||
|
threading.Thread.__init__(self)
|
||||||
|
self.threadID = threadID
|
||||||
|
def run(self):
|
||||||
|
s.acquire()
|
||||||
|
conv=isFinished()
|
||||||
|
s.release()
|
||||||
|
while conv:
|
||||||
|
doSim(4.,self.name)
|
||||||
|
s.acquire()
|
||||||
|
conv=isFinished()
|
||||||
|
s.release()
|
||||||
|
|
||||||
|
def doSim(delay,thread):
|
||||||
|
# s.acquire() and s.release() are couple
|
||||||
|
#
|
||||||
|
global dirCurrent
|
||||||
|
s.acquire()
|
||||||
|
delta_angle = offsetPhi()
|
||||||
|
if len(str(delta_angle)) > 5:
|
||||||
|
file_angle = str(delta_angle)[:5]
|
||||||
|
else:
|
||||||
|
file_angle = str(delta_angle)
|
||||||
|
dire = dirCurrent+'/'+file_angle
|
||||||
|
|
||||||
|
if not os.path.isdir(dire):
|
||||||
|
os.mkdir(dire,0755)
|
||||||
|
for file in [options.geometry+'.geom',options.load+'.load','numerics.config']:
|
||||||
|
copy2(dirCurrent+'/'+file, dire)
|
||||||
|
newMatConfig = newMaterialConfig(dirCurrent,delta_angle)
|
||||||
|
|
||||||
|
os.chdir(dire)
|
||||||
|
if not os.path.isfile('%s_%s.spectralOut'%(options.geometry,options.load)):
|
||||||
|
print('starting uniaxial tension in direction of angle %s from %s'%(file_angle,thread))
|
||||||
|
s.release()
|
||||||
|
execute('DAMASK_spectral -g %s -l %s'%(options.geometry,options.load))
|
||||||
|
else: s.release()
|
||||||
|
|
||||||
|
s.acquire()
|
||||||
|
if not os.path.isfile('./%s/%s_%s.txt'%('Rvalues',options.geometry,options.load)):
|
||||||
|
print('starting post processing for angle %s from %s'%(file_angle,thread))
|
||||||
|
s.release()
|
||||||
|
execute('postResults --cr f,p -d %s %s_%s.spectralOut'%('Rvalues',options.geometry,options.load))
|
||||||
|
execute('addCauchy ./%s/%s_%s.txt'%('Rvalues',options.geometry,options.load))
|
||||||
|
execute('addStrainTensors -l -v ./%s/%s_%s.txt'%('Rvalues',options.geometry,options.load))
|
||||||
|
print('post processing for angle %s from %s is finished'%(file_angle,thread))
|
||||||
|
|
||||||
|
else:
|
||||||
|
s.release()
|
||||||
|
os.chdir(dirCurrent)
|
||||||
|
|
||||||
|
def isFinished():
|
||||||
|
global N_simulations
|
||||||
|
|
||||||
|
if N_simulations < options.number:
|
||||||
|
return True
|
||||||
|
else:
|
||||||
|
return False
|
||||||
|
|
||||||
|
def offsetPhi():
|
||||||
|
global N_simulations #, N_tensile
|
||||||
|
N_simulations+=1
|
||||||
|
return np.linspace(0,90,options.number)[N_simulations-1]
|
||||||
|
|
||||||
|
def newMaterialConfig(dire,angle):
|
||||||
|
filename = '/material.config'
|
||||||
|
if len(str(angle)) > 5:
|
||||||
|
file_angle = str(angle)[:5]
|
||||||
|
else:
|
||||||
|
file_angle = str(angle)
|
||||||
|
f = open(dire+'/'+file_angle+filename,'w')
|
||||||
|
data = open(dire+filename, 'r').readlines()
|
||||||
|
|
||||||
|
for line in data:
|
||||||
|
if '(gauss)' in line:
|
||||||
|
linesplit = split(r'[;,\s,\t]\s*',line)
|
||||||
|
index = linesplit.index('phi1')
|
||||||
|
phi = float(linesplit[index+1]) - angle
|
||||||
|
if phi > 360. :
|
||||||
|
phi = phi-360.0
|
||||||
|
elif phi < 0.0:
|
||||||
|
phi = phi+360.0
|
||||||
|
index2 = line.index(linesplit[index+1])
|
||||||
|
line2 = line[:index2]+str(phi)+line[index2+len(str(float(linesplit[index+1]))):]
|
||||||
|
else:
|
||||||
|
line2 = line
|
||||||
|
f.write(line2)
|
||||||
|
f.close()
|
||||||
|
return True
|
||||||
|
|
||||||
|
# --------------------------------------------------------------------
|
||||||
|
# MAIN
|
||||||
|
# --------------------------------------------------------------------
|
||||||
|
|
||||||
|
parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [file[s]]', description = """
|
||||||
|
Calculate the deformation anisotropic coefficients (r values) and
|
||||||
|
strength anisotropic coefficients (normalized yield stress)
|
||||||
|
|
||||||
|
""", version=string.replace(scriptID,'\n','\\n')
|
||||||
|
)
|
||||||
|
|
||||||
|
parser.add_option('-l','--load' , dest='load', type='string',
|
||||||
|
help='name of the load file [%default]', metavar='string')
|
||||||
|
parser.add_option('-g','--geometry', dest='geometry', type='string',
|
||||||
|
help='name of the geometry file [%default]', metavar='string')
|
||||||
|
parser.add_option('-s', '--strain', dest='strain', type='string', action='callback', callback=list_split,
|
||||||
|
help='the threshold strains, using comma to seperate multiple strains [%default]', metavar='string')
|
||||||
|
parser.add_option('-t','--threads', dest='threads', type='int',
|
||||||
|
help='number of parallel executions [%default]', metavar='int')
|
||||||
|
parser.add_option('-n','--number', dest='number', type='int',
|
||||||
|
help='Number of uni-axial tensile tests [%default]', metavar='int')
|
||||||
|
|
||||||
|
parser.set_defaults(geometry = '20grains16x16x16')
|
||||||
|
parser.set_defaults(load = 'tensionX')
|
||||||
|
parser.set_defaults(threads = 1)
|
||||||
|
parser.set_defaults(number = 7)
|
||||||
|
parser.set_defaults(strain = ('0.0025','0.025'))
|
||||||
|
|
||||||
|
options = parser.parse_args()[0]
|
||||||
|
|
||||||
|
if not os.path.isfile(options.geometry+'.geom'):
|
||||||
|
parser.error('geometry file %s.geom not found'%options.geometry)
|
||||||
|
if not os.path.isfile('material.config'):
|
||||||
|
parser.error('material.config file not found')
|
||||||
|
|
||||||
|
thresStrain = [float(i) for i in options.strain]
|
||||||
|
|
||||||
|
if not os.path.isfile(options.load+'.load'):
|
||||||
|
maxstrain = max(thresStrain)
|
||||||
|
f11 = str(1.2*maxstrain + 1.0)
|
||||||
|
if maxstrain < 0.05:
|
||||||
|
maxincs = '60'; maxtime = '30'
|
||||||
|
else:
|
||||||
|
maxincs = str(int(maxstrain*3000)); maxtime = str(maxstrain*600)
|
||||||
|
|
||||||
|
f=open('tensionX.load','w')
|
||||||
|
uniaxial_load = 'f ' +f11+ ' 0 0 0 * 0 0 0 * ' + \
|
||||||
|
'stress ' +'* * * * 0 * * * 0 ' + \
|
||||||
|
'time '+ maxtime + ' incs ' + maxincs + ' freq 3'
|
||||||
|
f.write(uniaxial_load)
|
||||||
|
f.close()
|
||||||
|
|
||||||
|
N_simulations=0
|
||||||
|
s=threading.Semaphore(1)
|
||||||
|
threads=[]
|
||||||
|
dirCurrent = os.getcwd()
|
||||||
|
for i in range(options.threads):
|
||||||
|
threads.append(myThread(i))
|
||||||
|
threads[i].start()
|
||||||
|
|
||||||
|
for i in range(options.threads):
|
||||||
|
threads[i].join()
|
||||||
|
|
||||||
|
anisotropy = {}
|
||||||
|
for i,delta_angle in enumerate(np.linspace(0,90,options.number)):
|
||||||
|
if len(str(delta_angle)) > 5:
|
||||||
|
file_angle = str(delta_angle)[:5]
|
||||||
|
else:
|
||||||
|
file_angle = str(delta_angle)
|
||||||
|
refFile = dirCurrent+'/'+file_angle+'/Rvalues/%s_%s.txt'%(options.geometry,options.load)
|
||||||
|
if not os.path.isfile(refFile):
|
||||||
|
print('post processing file is not found for the angle %s'%file_angle)
|
||||||
|
else:
|
||||||
|
File = open(refFile)
|
||||||
|
table = damask.ASCIItable(File)
|
||||||
|
table.head_read()
|
||||||
|
if not set(['1_ln(V)','5_ln(V)','9_ln(V)','1_Cauchy']).issubset(set(table.labels)):
|
||||||
|
print 'data missing in direction %s'%str(delta_angle)
|
||||||
|
table.data_readArray(['%i_Cauchy'%(i+1) for i in xrange(9)]+['%i_ln(V)'%(i+1) for i in xrange(9)])
|
||||||
|
line, lines = 0, np.shape(table.data)[0]
|
||||||
|
aniso_thres = {}
|
||||||
|
for threshold in thresStrain:
|
||||||
|
while line < lines:
|
||||||
|
if abs(table.data[line,9])>= threshold: # 1_ln(V), e_xx
|
||||||
|
upper,lower = abs(table.data[line,9]),abs(table.data[line-1,9]) # values for linear interpolation
|
||||||
|
interplate = lambda x_d, x_u : x_d * (upper-threshold)/(upper-lower) + x_u * (threshold-lower)/(upper-lower)
|
||||||
|
e_yy = interplate (table.data[line-1,13], table.data[line,13])
|
||||||
|
e_zz = interplate (table.data[line-1,17], table.data[line,17])
|
||||||
|
s_xx = interplate (table.data[line-1,0 ], table.data[line,0 ])
|
||||||
|
aniso_thres[str(threshold)] = [e_yy/e_zz, s_xx]
|
||||||
|
break
|
||||||
|
else:
|
||||||
|
line+=1
|
||||||
|
anisotropy[file_angle] = aniso_thres
|
||||||
|
|
||||||
|
f = open('./anisotropy.txt','w')
|
||||||
|
f.write('4 header \n')
|
||||||
|
f.write('#the first row mean the threshold strain e_xx \n#the first column means loading directions \n')
|
||||||
|
f.write('#none means the data is unavailable\n# \n')
|
||||||
|
title = ['*R values (Lankford coefficients) \n', '# \n*Normalized yield stress \n']
|
||||||
|
for i in xrange(2):
|
||||||
|
f.write(title[i])
|
||||||
|
f.write(' '*10+len(thresStrain)*'%-12.4f'%tuple(thresStrain)+'\n')
|
||||||
|
for j,delta_angle in enumerate(np.linspace(0,90,options.number)):
|
||||||
|
if len(str(delta_angle)) > 5:
|
||||||
|
file_angle = str(delta_angle)[:5]
|
||||||
|
else:
|
||||||
|
file_angle = str(delta_angle)
|
||||||
|
|
||||||
|
if file_angle in anisotropy.keys():
|
||||||
|
aniso_dic_ang = anisotropy[file_angle]
|
||||||
|
aniso_list_ang_strain = []; writeformat = ''
|
||||||
|
for threshold in thresStrain:
|
||||||
|
if str(threshold) in aniso_dic_ang.keys():
|
||||||
|
if i == 1: # extract the normalized stress
|
||||||
|
aniso = aniso_dic_ang[str(threshold)][i]/anisotropy['0.0'][str(threshold)][i]
|
||||||
|
else: # extract r value
|
||||||
|
aniso = aniso_dic_ang[str(threshold)][i]
|
||||||
|
aniso_list_ang_strain.append(aniso)
|
||||||
|
writeformat = writeformat+'%-12.6f'
|
||||||
|
else:
|
||||||
|
aniso_list_ang_strain.append('none')
|
||||||
|
writeformat = writeformat+'%-12s'
|
||||||
|
f.write('%-10s'%file_angle + writeformat%(tuple(aniso_list_ang_strain))+'\n')
|
||||||
|
f.close()
|
||||||
Loading…
Reference in New Issue