From 7dc8391c03d42425a806a51e22b6d8a18a104539 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 25 Sep 2020 07:40:12 +0200 Subject: [PATCH] not needed anymore --- .gitlab-ci.yml | 37 -- PRIVATE | 2 +- processing/{post => legacy}/vtk2ang.py | 0 processing/post/filterTable.py | 161 -------- .../pre/geom_fromVoronoiTessellation.py | 226 ------------ processing/pre/hybridIA_linODFsampling.py | 349 ------------------ processing/pre/seeds_fromGeom.py | 68 ---- processing/pre/seeds_fromRandom.py | 165 --------- 8 files changed, 1 insertion(+), 1007 deletions(-) rename processing/{post => legacy}/vtk2ang.py (100%) delete mode 100755 processing/post/filterTable.py delete mode 100755 processing/pre/geom_fromVoronoiTessellation.py delete mode 100755 processing/pre/hybridIA_linODFsampling.py delete mode 100755 processing/pre/seeds_fromGeom.py delete mode 100755 processing/pre/seeds_fromRandom.py diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 40662cb40..054e18ae2 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -253,27 +253,6 @@ grid_parsingArguments: - master - release -StateIntegration_compareVariants: - stage: grid - script: StateIntegration_compareVariants/test.py - except: - - master - - release - -nonlocal_densityConservation: - stage: grid - script: nonlocal_densityConservation/test.py - except: - - master - - release - -RGC_DetectChanges: - stage: grid - script: RGC_DetectChanges/test.py - except: - - master - - release - Nonlocal_Damage_DetectChanges: stage: grid script: Nonlocal_Damage_DetectChanges/test.py @@ -302,15 +281,6 @@ grid_all_loadCaseRotation: - master - release -grid_mech_MPI: - stage: grid - script: - - module load $IntelCompiler $MPICH_Intel $PETSc_MPICH_Intel - - grid_mech_MPI/test.py - except: - - master - - release - grid_all_restartMPI: stage: grid script: @@ -327,13 +297,6 @@ Plasticity_DetectChanges: - master - release -Homogenization: - stage: grid - script: Homogenization/test.py - except: - - master - - release - Phenopowerlaw_singleSlip: stage: grid script: Phenopowerlaw_singleSlip/test.py diff --git a/PRIVATE b/PRIVATE index b73dcfe74..b7c6e3a67 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit b73dcfe746eadce92ed0ab08a3bfbb19f5c8eb0a +Subproject commit b7c6e3a6742b6098530a6fa5121b23ec9487dd4f diff --git a/processing/post/vtk2ang.py b/processing/legacy/vtk2ang.py similarity index 100% rename from processing/post/vtk2ang.py rename to processing/legacy/vtk2ang.py diff --git a/processing/post/filterTable.py b/processing/post/filterTable.py deleted file mode 100755 index 4f4af088b..000000000 --- a/processing/post/filterTable.py +++ /dev/null @@ -1,161 +0,0 @@ -#!/usr/bin/env python3 - -import os -import sys -from optparse import OptionParser -import re -import fnmatch -import math # noqa - -import numpy as np - -import damask - - -scriptName = os.path.splitext(os.path.basename(__file__))[0] -scriptID = ' '.join([scriptName,damask.version]) - -def sortingList(labels,whitelistitems): - - indices = [] - names = [] - - for label in labels: - if re.match(r'^\d+_',label): - indices.append(int(label.split('_',1)[0])) - names.append(label.split('_',1)[1]) - else: - indices.append(0) - names.append(label) - - return [indices,names,whitelistitems] - - -# -------------------------------------------------------------------- -# MAIN -# -------------------------------------------------------------------- - -parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [ASCIItable(s)]', description = """ -Filter rows according to condition and columns by either white or black listing. - -Examples: -Every odd row if x coordinate is positive -- " #ip.x# >= 0.0 and #_row_#%2 == 1 ). -All rows where label 'foo' equals 'bar' -- " #s#foo# == 'bar' " - -""", version = scriptID) - -parser.add_option('-w','--white', - dest = 'whitelist', - action = 'extend', metavar = '', - help = 'whitelist of column labels (a,b,c,...)') -parser.add_option('-b','--black', - dest = 'blacklist', - action = 'extend', metavar='', - help = 'blacklist of column labels (a,b,c,...)') -parser.add_option('-c','--condition', - dest = 'condition', metavar='string', - help = 'condition to filter rows') - -parser.set_defaults(condition = None, - ) - -(options,filenames) = parser.parse_args() - -# --- loop over input files ------------------------------------------------------------------------- - -if filenames == []: filenames = [None] - -for name in filenames: - try: - table = damask.ASCIItable(name = name) - except IOError: - continue - damask.util.report(scriptName,name) - -# ------------------------------------------ assemble info --------------------------------------- - - table.head_read() - -# ------------------------------------------ process data --------------------------------------- - - specials = { \ - '_row_': 0, - } - labels = [] - positions = [] - - for position,label in enumerate(table.labels(raw = True)): - if (options.whitelist is None or any([ position in table.label_indexrange(needle) \ - or fnmatch.fnmatch(label,needle) for needle in options.whitelist])) \ - and (options.blacklist is None or not any([ position in table.label_indexrange(needle) \ - or fnmatch.fnmatch(label,needle) for needle in options.blacklist])): # a label to keep? - labels.append(label) # remember name... - positions.append(position) # ...and position - - if len(labels) > 0 and options.whitelist is not None and options.blacklist is None: # check whether reordering is possible - whitelistitem = np.zeros(len(labels),dtype=int) - for i,label in enumerate(labels): # check each selected label - match = [ positions[i] in table.label_indexrange(needle) \ - or fnmatch.fnmatch(label,needle) for needle in options.whitelist] # which whitelist items do match it - whitelistitem[i] = match.index(True) if np.sum(match) == 1 else -1 # unique match to a whitelist item --> store which - - order = range(len(labels)) if np.any(whitelistitem < 0) \ - else np.lexsort(sortingList(labels,whitelistitem)) # reorder if unique, i.e. no "-1" in whitelistitem - else: - order = range(len(labels)) # maintain original order of labels - -# --------------------------------------- evaluate condition --------------------------------------- - if options.condition is not None: - condition = options.condition # copy per file, since might be altered inline - breaker = False - - for position,(all,marker,column) in enumerate(set(re.findall(r'#(([s]#)?(.+?))#',condition))): # find three groups - idx = table.label_index(column) - dim = table.label_dimension(column) - if idx < 0 and column not in specials: - damask.util.croak('column "{}" not found.'.format(column)) - breaker = True - else: - if column in specials: - replacement = 'specials["{}"]'.format(column) - elif dim == 1: # scalar input - replacement = '{}(table.data[{}])'.format({ '':'float', - 's#':'str'}[marker],idx) # take float or string value of data column - elif dim > 1: # multidimensional input (vector, tensor, etc.) - replacement = 'np.array(table.data[{}:{}],dtype=float)'.format(idx,idx+dim) # use (flat) array representation - - condition = condition.replace('#'+all+'#',replacement) - - if breaker: continue # found mistake in condition evaluation --> next file - -# ------------------------------------------ assemble header --------------------------------------- - - table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:])) - table.labels_clear() - table.labels_append(np.array(labels)[order]) # update with new label set - table.head_write() - -# ------------------------------------------ process and output data ------------------------------------------ - - positions = np.array(positions)[order] - - atOnce = options.condition is None - if atOnce: # read full array and filter columns - try: - table.data_readArray(positions+1) # read desired columns (indexed 1,...) - table.data_writeArray() # directly write out - except Exception: - table.data_rewind() - atOnce = False # data contains items that prevent array chunking - - if not atOnce: # read data line by line - outputAlive = True - while outputAlive and table.data_read(): # read next data line of ASCII table - specials['_row_'] += 1 # count row - if options.condition is None or eval(condition): # valid row ? - table.data = [table.data[position] for position in positions] # retain filtered columns - outputAlive = table.data_write() # output processed line - -# ------------------------------------------ finalize output ----------------------------------------- - - table.close() # close input ASCII table (works for stdin) diff --git a/processing/pre/geom_fromVoronoiTessellation.py b/processing/pre/geom_fromVoronoiTessellation.py deleted file mode 100755 index aee79cc05..000000000 --- a/processing/pre/geom_fromVoronoiTessellation.py +++ /dev/null @@ -1,226 +0,0 @@ -#!/usr/bin/env python3 - -import os -import sys -import multiprocessing -from io import StringIO -from functools import partial -from optparse import OptionParser,OptionGroup - -import numpy as np -from scipy import spatial - -import damask - - -scriptName = os.path.splitext(os.path.basename(__file__))[0] -scriptID = ' '.join([scriptName,damask.version]) - - -def findClosestSeed(seeds, weights, point): - return np.argmin(np.sum((np.broadcast_to(point,(len(seeds),3))-seeds)**2,axis=1) - weights) - - -def Laguerre_tessellation(grid, size, seeds, weights, origin = np.zeros(3), periodic = True, cpus = 2): - - if periodic: - weights_p = np.tile(weights.squeeze(),27) # Laguerre weights (1,2,3,1,2,3,...,1,2,3) - seeds_p = np.vstack((seeds -np.array([size[0],0.,0.]),seeds, seeds +np.array([size[0],0.,0.]))) - seeds_p = np.vstack((seeds_p-np.array([0.,size[1],0.]),seeds_p,seeds_p+np.array([0.,size[1],0.]))) - seeds_p = np.vstack((seeds_p-np.array([0.,0.,size[2]]),seeds_p,seeds_p+np.array([0.,0.,size[2]]))) - coords = damask.grid_filters.cell_coord0(grid*3,size*3,-origin-size).reshape(-1,3) - else: - weights_p = weights.squeeze() - seeds_p = seeds - coords = damask.grid_filters.cell_coord0(grid,size,-origin).reshape(-1,3) - - if cpus > 1: - pool = multiprocessing.Pool(processes = cpus) - result = pool.map_async(partial(findClosestSeed,seeds_p,weights_p), [coord for coord in coords]) - pool.close() - pool.join() - closest_seed = np.array(result.get()).reshape(-1,3) - else: - closest_seed= np.array([findClosestSeed(seeds_p,weights_p,coord) for coord in coords]) - - if periodic: - closest_seed = closest_seed.reshape(grid*3) - return closest_seed[grid[0]:grid[0]*2,grid[1]:grid[1]*2,grid[2]:grid[2]*2]%seeds.shape[0] - else: - return closest_seed - - -def Voronoi_tessellation(grid, size, seeds, origin = np.zeros(3), periodic = True): - - coords = damask.grid_filters.cell_coord0(grid,size,-origin).reshape(-1,3) - KDTree = spatial.cKDTree(seeds,boxsize=size) if periodic else spatial.cKDTree(seeds) - devNull,closest_seed = KDTree.query(coords) - - return closest_seed - - -# -------------------------------------------------------------------- -# MAIN -# -------------------------------------------------------------------- - -parser = OptionParser(option_class=damask.extendableOption, usage='%prog option(s) [seedfile(s)]', description = """ -Generate geometry description and material configuration by tessellation of given seeds file. - -""", version = scriptID) - - -group = OptionGroup(parser, "Tessellation","") - -group.add_option('-l', - '--laguerre', - dest = 'laguerre', - action = 'store_true', - help = 'use Laguerre (weighted Voronoi) tessellation') -group.add_option('--cpus', - dest = 'cpus', - type = 'int', metavar = 'int', - help = 'number of parallel processes to use for Laguerre tessellation [%default]') -group.add_option('--nonperiodic', - dest = 'periodic', - action = 'store_false', - help = 'nonperiodic tessellation') - -parser.add_option_group(group) - -group = OptionGroup(parser, "Geometry","") - -group.add_option('-g', - '--grid', - dest = 'grid', - type = 'int', nargs = 3, metavar = ' '.join(['int']*3), - help = 'a,b,c grid of hexahedral box') -group.add_option('-s', - '--size', - dest = 'size', - type = 'float', nargs = 3, metavar=' '.join(['float']*3), - help = 'x,y,z size of hexahedral box [1.0 1.0 1.0]') -group.add_option('-o', - '--origin', - dest = 'origin', - type = 'float', nargs = 3, metavar=' '.join(['float']*3), - help = 'origin of grid [0.0 0.0 0.0]') - -parser.add_option_group(group) - -group = OptionGroup(parser, "Seeds","") - -group.add_option('-p', - '--pos', '--seedposition', - dest = 'pos', - type = 'string', metavar = 'string', - help = 'label of coordinates [%default]') -group.add_option('-w', - '--weight', - dest = 'weight', - type = 'string', metavar = 'string', - help = 'label of weights [%default]') -group.add_option('-m', - '--microstructure', - dest = 'microstructure', - type = 'string', metavar = 'string', - help = 'label of microstructures [%default]') -group.add_option('-e', - '--eulers', - dest = 'eulers', - type = 'string', metavar = 'string', - help = 'label of Euler angles [%default]') -group.add_option('--axes', - dest = 'axes', - type = 'string', nargs = 3, metavar = ' '.join(['string']*3), - help = 'orientation coordinate frame in terms of position coordinate frame') - -parser.add_option_group(group) - -group = OptionGroup(parser, "Configuration","") - -group.add_option('--without-config', - dest = 'config', - action = 'store_false', - help = 'omit material configuration header') -group.add_option('--phase', - dest = 'phase', - type = 'int', metavar = 'int', - help = 'phase index to be used [%default]') - -parser.add_option_group(group) - -parser.set_defaults(pos = 'pos', - weight = 'weight', - microstructure = 'microstructure', - eulers = 'euler', - phase = 1, - cpus = 2, - laguerre = False, - periodic = True, - config = True, - ) - -(options,filenames) = parser.parse_args() -if filenames == []: filenames = [None] - -for name in filenames: - damask.util.report(scriptName,name) - - table = damask.Table.load(StringIO(''.join(sys.stdin.read())) if name is None else name) - - size = np.ones(3) - origin = np.zeros(3) - for line in table.comments: - items = line.lower().strip().split() - key = items[0] if items else '' - if key == 'grid': - grid = np.array([ int(dict(zip(items[1::2],items[2::2]))[i]) for i in ['a','b','c']]) - elif key == 'size': - size = np.array([float(dict(zip(items[1::2],items[2::2]))[i]) for i in ['x','y','z']]) - elif key == 'origin': - origin = np.array([float(dict(zip(items[1::2],items[2::2]))[i]) for i in ['x','y','z']]) - if options.grid: grid = np.array(options.grid) - if options.size: size = np.array(options.size) - if options.origin: origin = np.array(options.origin) - - seeds = table.get(options.pos) - - grains = table.get(options.microstructure) if options.microstructure in table.labels else np.arange(len(seeds))+1 - grainIDs = np.unique(grains).astype('i') - - if options.eulers in table.labels: - eulers = table.get(options.eulers) - - if options.laguerre: - indices = grains[Laguerre_tessellation(grid,size,seeds,table.get(options.weight),origin, - options.periodic,options.cpus)] - else: - indices = grains[Voronoi_tessellation (grid,size,seeds,origin,options.periodic)] - - config_header = [] - if options.config: - - if options.eulers in table.labels: - config_header += [''] - for ID in grainIDs: - eulerID = np.nonzero(grains == ID)[0][0] # find first occurrence of this grain id - config_header += ['[Grain{}]'.format(ID), - '(gauss)\tphi1 {:.2f}\tPhi {:.2f}\tphi2 {:.2f}'.format(*eulers[eulerID]) - ] - if options.axes: config_header += ['axes\t{} {} {}'.format(*options.axes)] - - config_header += [''] - for ID in grainIDs: - config_header += ['[Grain{}]'.format(ID), - '(constituent)\tphase {}\ttexture {}\tfraction 1.0'.format(options.phase,ID) - ] - - config_header += [''] - - header = [scriptID + ' ' + ' '.join(sys.argv[1:])]\ - + config_header - geom = damask.Geom(indices.reshape(grid),size,origin, - comments=header) - damask.util.croak(geom) - - geom.save_ASCII(sys.stdout if name is None else os.path.splitext(name)[0]+'.geom',compress=False) diff --git a/processing/pre/hybridIA_linODFsampling.py b/processing/pre/hybridIA_linODFsampling.py deleted file mode 100755 index f99a2dd89..000000000 --- a/processing/pre/hybridIA_linODFsampling.py +++ /dev/null @@ -1,349 +0,0 @@ -#!/usr/bin/env python3 - -import os -import sys -import math -import random -from io import StringIO -from optparse import OptionParser - -import numpy as np - -import damask - - -scriptName = os.path.splitext(os.path.basename(__file__))[0] -scriptID = ' '.join([scriptName,damask.version]) - - -# --- helper functions --- -def integerFactorization(i): - - j = int(math.floor(math.sqrt(float(i)))) - while j>1 and int(i)%j != 0: - j -= 1 - return j - -def binAsBins(bin,intervals): - """Explode compound bin into 3D bins list.""" - bins = [0]*3 - bins[0] = (bin//(intervals[1] * intervals[2])) % intervals[0] - bins[1] = (bin//intervals[2]) % intervals[1] - bins[2] = bin % intervals[2] - return bins - -def binsAsBin(bins,intervals): - """Implode 3D bins into compound bin.""" - return (bins[0]*intervals[1] + bins[1])*intervals[2] + bins[2] - -def EulersAsBins(Eulers,intervals,deltas,center): - """Return list of Eulers translated into 3D bins list.""" - return [int((euler+(0.5-center)*delta)//delta)%interval \ - for euler,delta,interval in zip(Eulers,deltas,intervals) \ - ] - -def binAsEulers(bin,intervals,deltas,center): - """Compound bin number translated into list of Eulers.""" - Eulers = [0.0]*3 - Eulers[2] = (bin%intervals[2] + center)*deltas[2] - Eulers[1] = (bin//intervals[2]%intervals[1] + center)*deltas[1] - Eulers[0] = (bin//(intervals[2]*intervals[1]) + center)*deltas[0] - return Eulers - -def directInvRepetitions(probability,scale): - """Calculate number of samples drawn by direct inversion.""" - nDirectInv = 0 - for bin in range(len(probability)): # loop over bins - nDirectInv += int(round(probability[bin]*scale)) # calc repetition - return nDirectInv - - -# ---------------------- sampling methods ----------------------------------------------------------------------- - -# ----- efficient algorithm --------- -def directInversion (ODF,nSamples): - """ODF contains 'dV_V' (normalized to 1), 'center', 'intervals', 'limits' (in radians).""" - nOptSamples = max(ODF['nNonZero'],nSamples) # random subsampling if too little samples requested - - nInvSamples = 0 - repetition = [None]*ODF['nBins'] - - scaleLower = 0.0 - nInvSamplesLower = 0 - scaleUpper = float(nOptSamples) - incFactor = 1.0 - nIter = 0 - nInvSamplesUpper = directInvRepetitions(ODF['dV_V'],scaleUpper) - while (\ - (scaleUpper-scaleLower > scaleUpper*1e-15 or nInvSamplesUpper < nOptSamples) and \ - nInvSamplesUpper != nOptSamples \ - ): # closer match required? - if nInvSamplesUpper < nOptSamples: - scaleLower,scaleUpper = scaleUpper,scaleUpper+incFactor*(scaleUpper-scaleLower)/2.0 - incFactor *= 2.0 - nInvSamplesLower,nInvSamplesUpper = nInvSamplesUpper,directInvRepetitions(ODF['dV_V'],scaleUpper) - else: - scaleUpper = (scaleLower+scaleUpper)/2.0 - incFactor = 1.0 - nInvSamplesUpper = directInvRepetitions(ODF['dV_V'],scaleUpper) - nIter += 1 - damask.util.croak('%i:(%12.11f,%12.11f) %i <= %i <= %i'%(nIter,scaleLower,scaleUpper, - nInvSamplesLower,nOptSamples,nInvSamplesUpper)) - nInvSamples = nInvSamplesUpper - scale = scaleUpper - damask.util.croak('created set of %i samples (%12.11f) with scaling %12.11f delivering %i'%(nInvSamples, - float(nInvSamples)/nOptSamples-1.0, - scale,nSamples)) - repetition = [None]*ODF['nBins'] # preallocate and clear - - for bin in range(ODF['nBins']): # loop over bins - repetition[bin] = int(round(ODF['dV_V'][bin]*scale)) # calc repetition - -# build set - set = [None]*nInvSamples - i = 0 - for bin in range(ODF['nBins']): - set[i:i+repetition[bin]] = [bin]*repetition[bin] # fill set with bin, i.e. orientation - i += repetition[bin] # advance set counter - - orientations = np.zeros((nSamples,3),'f') - reconstructedODF = np.zeros(ODF['nBins'],'f') - unitInc = 1.0/nSamples - for j in range(nSamples): - if (j == nInvSamples-1): ex = j - else: ex = int(round(random.uniform(j+0.5,nInvSamples-0.5))) - bin = set[ex] - Eulers = binAsEulers(bin,ODF['interval'],ODF['delta'],ODF['center']) - orientations[j] = np.degrees(Eulers) - reconstructedODF[bin] += unitInc - set[ex] = set[j] # exchange orientations - - return orientations, reconstructedODF - - -# ----- trial and error algorithms --------- - -def MonteCarloEulers (ODF,nSamples): - """ODF contains 'dV_V' (normalized to 1), 'center', 'intervals', 'limits' (in radians).""" - countMC = 0 - maxdV_V = max(ODF['dV_V']) - orientations = np.zeros((nSamples,3),'f') - reconstructedODF = np.zeros(ODF['nBins'],'f') - unitInc = 1.0/nSamples - - for j in range(nSamples): - MC = maxdV_V*2.0 - bin = 0 - while MC > ODF['dV_V'][bin]: - countMC += 1 - MC = maxdV_V*random.random() - Eulers = [limit*random.random() for limit in ODF['limit']] - bins = EulersAsBins(Eulers,ODF['interval'],ODF['delta'],ODF['center']) - bin = binsAsBin(bins,ODF['interval']) - orientations[j] = np.degrees(Eulers) - reconstructedODF[bin] += unitInc - - return orientations, reconstructedODF, countMC - - -def MonteCarloBins (ODF,nSamples): - """ODF contains 'dV_V' (normalized to 1), 'center', 'intervals', 'limits' (in radians).""" - countMC = 0 - maxdV_V = max(ODF['dV_V']) - orientations = np.zeros((nSamples,3),'f') - reconstructedODF = np.zeros(ODF['nBins'],'f') - unitInc = 1.0/nSamples - - for j in range(nSamples): - MC = maxdV_V*2.0 - bin = 0 - while MC > ODF['dV_V'][bin]: - countMC += 1 - MC = maxdV_V*random.random() - bin = int(ODF['nBins'] * random.random()) - Eulers = binAsEulers(bin,ODF['interval'],ODF['delta'],ODF['center']) - orientations[j] = np.degrees(Eulers) - reconstructedODF[bin] += unitInc - - return orientations, reconstructedODF - - -def TothVanHoutteSTAT (ODF,nSamples): - """ODF contains 'dV_V' (normalized to 1), 'center', 'intervals', 'limits' (in radians).""" - orientations = np.zeros((nSamples,3),'f') - reconstructedODF = np.zeros(ODF['nBins'],'f') - unitInc = 1.0/nSamples - - selectors = [random.random() for i in range(nSamples)] - selectors.sort() - indexSelector = 0 - - cumdV_V = 0.0 - countSamples = 0 - - for bin in range(ODF['nBins']) : - cumdV_V += ODF['dV_V'][bin] - while indexSelector < nSamples and selectors[indexSelector] < cumdV_V: - Eulers = binAsEulers(bin,ODF['interval'],ODF['delta'],ODF['center']) - orientations[countSamples] = np.degrees(Eulers) - reconstructedODF[bin] += unitInc - countSamples += 1 - indexSelector += 1 - - damask.util.croak('created set of %i when asked to deliver %i'%(countSamples,nSamples)) - - return orientations, reconstructedODF - - -# -------------------------------------------------------------------- -# MAIN -# -------------------------------------------------------------------- -parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [file[s]]', description =""" -Transform linear binned ODF data into given number of orientations. -IA: integral approximation, STAT: Van Houtte, MC: Monte Carlo - -""", version = scriptID) - -algorithms = ['IA', 'STAT','MC'] -parser.add_option('-n', '--nsamples', - dest = 'number', - type = 'int', metavar = 'int', - help = 'number of orientations to be generated [%default]') -parser.add_option('-a','--algorithm', - dest = 'algorithm', - choices = algorithms, metavar = 'string', - help = 'sampling algorithm {%s} [IA]'%(', '.join(algorithms))) -parser.add_option('-p','--phase', - dest = 'phase', - type = 'int', metavar = 'int', - help = 'phase index to be used [%default]') -parser.add_option('-r', '--rnd', - dest = 'randomSeed', - type = 'int', metavar = 'int', \ - help = 'seed of random number generator [%default]') -parser.set_defaults(randomSeed = None, - number = 500, - algorithm = 'IA', - phase = 1, - ang = True, - ) - -(options,filenames) = parser.parse_args() -if filenames == []: filenames = [None] - -for name in filenames: - damask.util.report(scriptName,name) - - table = damask.Table.load(StringIO(''.join(sys.stdin.read())) if name is None else name) - - randomSeed = int(os.urandom(4).hex(),16) if options.randomSeed is None else options.randomSeed # random seed per file - random.seed(randomSeed) - -# --------------- figure out limits (left/right), delta, and interval ----------------------------- - ODF = {} - eulers = table.get('euler') - limits = np.array([np.min(eulers,axis=0),np.max(eulers,axis=0)]) # min/max euler angles in degrees - ODF['limit'] = np.radians(limits[1,:]) # right hand limits in radians - ODF['center'] = 0.0 if all(limits[0,:]<1e-8) else 0.5 # vertex or cell centered - - ODF['interval'] = np.array(list(map(len,[np.unique(eulers[:,i]) for i in range(3)])),'i') # steps are number of distict values - ODF['nBins'] = ODF['interval'].prod() - ODF['delta'] = np.radians(np.array(limits[1,0:3]-limits[0,0:3])/(ODF['interval']-1)) # step size - - if eulers.shape[0] != ODF['nBins']: - damask.util.croak('expecting %i values but got %i'%(ODF['nBins'],eulers.shape[0])) - continue - -# ----- build binnedODF array and normalize ------------------------------------------------------ - sumdV_V = 0.0 - ODF['dV_V'] = [None]*ODF['nBins'] - ODF['nNonZero'] = 0 - dg = ODF['delta'][0]*2.0*math.sin(ODF['delta'][1]/2.0)*ODF['delta'][2] - intensity = table.get('intensity') - for b in range(ODF['nBins']): - ODF['dV_V'][b] = \ - max(0.0,intensity[b,0]) * dg * \ - math.sin(((b//ODF['interval'][2])%ODF['interval'][1]+ODF['center'])*ODF['delta'][1]) - if ODF['dV_V'][b] > 0.0: - sumdV_V += ODF['dV_V'][b] - ODF['nNonZero'] += 1 - - for b in range(ODF['nBins']): - ODF['dV_V'][b] /= sumdV_V # normalize dV/V - - damask.util.croak(['non-zero fraction: %12.11f (%i/%i)'%(float(ODF['nNonZero'])/ODF['nBins'], - ODF['nNonZero'], - ODF['nBins']), - 'Volume integral of ODF: %12.11f\n'%sumdV_V, - 'Reference Integral: %12.11f\n'%(ODF['limit'][0]*ODF['limit'][2]*(1-math.cos(ODF['limit'][1]))), - ]) - - Functions = {'IA': 'directInversion', 'STAT': 'TothVanHoutteSTAT', 'MC': 'MonteCarloBins'} - method = Functions[options.algorithm] - - Orientations, ReconstructedODF = (globals()[method])(ODF,options.number) - -# calculate accuracy of sample - squaredDiff = {'orig':0.0,method:0.0} - squaredRelDiff = {'orig':0.0,method:0.0} - mutualProd = {'orig':0.0,method:0.0} - indivSum = {'orig':0.0,method:0.0} - indivSquaredSum = {'orig':0.0,method:0.0} - - for bin in range(ODF['nBins']): - squaredDiff[method] += (ODF['dV_V'][bin] - ReconstructedODF[bin])**2 - if ODF['dV_V'][bin] > 0.0: - squaredRelDiff[method] += (ODF['dV_V'][bin] - ReconstructedODF[bin])**2/ODF['dV_V'][bin]**2 - mutualProd[method] += ODF['dV_V'][bin]*ReconstructedODF[bin] - indivSum[method] += ReconstructedODF[bin] - indivSquaredSum[method] += ReconstructedODF[bin]**2 - indivSum['orig'] += ODF['dV_V'][bin] - indivSquaredSum['orig'] += ODF['dV_V'][bin]**2 - - damask.util.croak(['sqrt(N*)RMSD of ODFs:\t %12.11f'% math.sqrt(options.number*squaredDiff[method]), - 'RMSrD of ODFs:\t %12.11f'%math.sqrt(squaredRelDiff[method]), - 'rMSD of ODFs:\t %12.11f'%(squaredDiff[method]/indivSquaredSum['orig']), - 'nNonZero correlation slope:\t %12.11f'\ - %((ODF['nNonZero']*mutualProd[method]-indivSum['orig']*indivSum[method])/\ - (ODF['nNonZero']*indivSquaredSum['orig']-indivSum['orig']**2)), - 'nNonZero correlation confidence:\t %12.11f'\ - %((mutualProd[method]-indivSum['orig']*indivSum[method]/ODF['nNonZero'])/\ - (ODF['nNonZero']*math.sqrt((indivSquaredSum['orig']/ODF['nNonZero']-(indivSum['orig']/ODF['nNonZero'])**2)*\ - (indivSquaredSum[method]/ODF['nNonZero']-(indivSum[method]/ODF['nNonZero'])**2)))), - ]) - - if method == 'IA' and options.number < ODF['nNonZero']: - strOpt = '(%i)'%ODF['nNonZero'] - - formatwidth = 1+int(math.log10(options.number)) - - materialConfig = [ - '#' + scriptID + ' ' + ' '.join(sys.argv[1:]), - '# random seed %i'%randomSeed, - '#-------------------#', - '', - '#-------------------#', - ] - - for i,ID in enumerate(range(options.number)): - materialConfig += ['[Grain%s]'%(str(ID+1).zfill(formatwidth)), - '(constituent) phase %i texture %s fraction 1.0'%(options.phase,str(ID+1).rjust(formatwidth)), - ] - - materialConfig += [ - '#-------------------#', - '', - '#-------------------#', - ] - - for ID in range(options.number): - eulers = Orientations[ID] - - materialConfig += ['[Grain%s]'%(str(ID+1).zfill(formatwidth)), - '(gauss) phi1 {} Phi {} phi2 {} scatter 0.0 fraction 1.0'.format(*eulers), - ] - -#--- output finalization -------------------------------------------------------------------------- - - with (open(os.path.splitext(name)[0]+'_'+method+'_'+str(options.number)+'_material.config','w')) as outfile: - outfile.write('\n'.join(materialConfig)+'\n') diff --git a/processing/pre/seeds_fromGeom.py b/processing/pre/seeds_fromGeom.py deleted file mode 100755 index b8d74b651..000000000 --- a/processing/pre/seeds_fromGeom.py +++ /dev/null @@ -1,68 +0,0 @@ -#!/usr/bin/env python3 - -import os -import sys -from io import StringIO -from optparse import OptionParser - -import numpy as np - -import damask - -scriptName = os.path.splitext(os.path.basename(__file__))[0] -scriptID = ' '.join([scriptName,damask.version]) - -#-------------------------------------------------------------------------------------------------- -# MAIN -#-------------------------------------------------------------------------------------------------- - -parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [file[s]]', description = """ -Create seed file taking material indices from given geom file. -Indices can be black-listed or white-listed. - -""", version = scriptID) - -parser.add_option('-w', - '--white', - action = 'extend', metavar = '', - dest = 'whitelist', - help = 'whitelist of grain IDs') -parser.add_option('-b', - '--black', - action = 'extend', metavar = '', - dest = 'blacklist', - help = 'blacklist of grain IDs') - -parser.set_defaults(whitelist = [], - blacklist = [], - ) - -(options,filenames) = parser.parse_args() -if filenames == []: filenames = [None] - -options.whitelist = [int(i) for i in options.whitelist] -options.blacklist = [int(i) for i in options.blacklist] - -for name in filenames: - damask.util.report(scriptName,name) - - geom = damask.Geom.load_ASCII(StringIO(''.join(sys.stdin.read())) if name is None else name) - material = geom.material.reshape((-1,1),order='F') - - mask = np.logical_and(np.in1d(material,options.whitelist,invert=False) if options.whitelist else \ - np.full(geom.grid.prod(),True,dtype=bool), - np.in1d(material,options.blacklist,invert=True) if options.blacklist else \ - np.full(geom.grid.prod(),True,dtype=bool)) - - seeds = damask.grid_filters.cell_coord0(geom.grid,geom.size).reshape(-1,3,order='F') - - comments = geom.comments \ - + [scriptID + ' ' + ' '.join(sys.argv[1:]), - 'grid\ta {}\tb {}\tc {}'.format(*geom.grid), - 'size\tx {}\ty {}\tz {}'.format(*geom.size), - 'origin\tx {}\ty {}\tz {}'.format(*geom.origin), - ] - - damask.Table(seeds[mask],{'pos':(3,)},comments)\ - .add('material',material[mask].astype(int))\ - .save(sys.stdout if name is None else os.path.splitext(name)[0]+'.seeds',legacy=True) diff --git a/processing/pre/seeds_fromRandom.py b/processing/pre/seeds_fromRandom.py deleted file mode 100755 index 451e218aa..000000000 --- a/processing/pre/seeds_fromRandom.py +++ /dev/null @@ -1,165 +0,0 @@ -#!/usr/bin/env python3 - -import os -import sys -from optparse import OptionParser,OptionGroup - -import numpy as np -from scipy import spatial - -import damask - - -scriptName = os.path.splitext(os.path.basename(__file__))[0] -scriptID = ' '.join([scriptName,damask.version]) - - -# -------------------------------------------------------------------- -# MAIN -# -------------------------------------------------------------------- - -parser = OptionParser(option_class=damask.extendableOption, usage='%prog options', description = """ -Distribute given number of points randomly within rectangular cuboid. -Reports positions with random crystal orientations in seeds file format to STDOUT. - -""", version = scriptID) - -parser.add_option('-N', - dest = 'N', - type = 'int', metavar = 'int', - help = 'number of seed points [%default]') -parser.add_option('-s', - '--size', - dest = 'size', - type = 'float', nargs = 3, metavar = 'float float float', - help='size x,y,z of unit cube to fill %default') -parser.add_option('-g', - '--grid', - dest = 'grid', - type = 'int', nargs = 3, metavar = 'int int int', - help='min a,b,c grid of hexahedral box %default') -parser.add_option('-m', - '--microstructure', - dest = 'microstructure', - type = 'int', metavar = 'int', - help = 'first microstructure index [%default]') -parser.add_option('-r', - '--rnd', - dest = 'randomSeed', type = 'int', metavar = 'int', - help = 'seed of random number generator [%default]') - -group = OptionGroup(parser, "Laguerre Tessellation", - "Parameters determining shape of weight distribution of seed points" - ) -group.add_option( '-w', - '--weights', - action = 'store_true', - dest = 'weights', - help = 'assign random weights to seed points for Laguerre tessellation [%default]') -group.add_option( '--max', - dest = 'max', - type = 'float', metavar = 'float', - help = 'max of uniform distribution for weights [%default]') -group.add_option( '--mean', - dest = 'mean', - type = 'float', metavar = 'float', - help = 'mean of normal distribution for weights [%default]') -group.add_option( '--sigma', - dest = 'sigma', - type = 'float', metavar = 'float', - help='standard deviation of normal distribution for weights [%default]') -parser.add_option_group(group) - -group = OptionGroup(parser, "Selective Seeding", - "More uniform distribution of seed points using Mitchell's Best Candidate Algorithm" - ) -group.add_option( '--selective', - action = 'store_true', - dest = 'selective', - help = 'selective picking of seed points from random seed points') -group.add_option( '--distance', - dest = 'distance', - type = 'float', metavar = 'float', - help = 'minimum distance to next neighbor [%default]') -group.add_option( '--numCandidates', - dest = 'numCandidates', - type = 'int', metavar = 'int', - help = 'size of point group to select best distance from [%default]') -parser.add_option_group(group) - -parser.set_defaults(randomSeed = None, - grid = (16,16,16), - size = (1.0,1.0,1.0), - N = 20, - weights = False, - max = 0.0, - mean = 0.2, - sigma = 0.05, - microstructure = 1, - selective = False, - distance = 0.2, - numCandidates = 10, - ) - -(options,filenames) = parser.parse_args() -if filenames == []: filenames = [None] - -size = np.array(options.size) -grid = np.array(options.grid) -np.random.seed(int(os.urandom(4).hex(),16) if options.randomSeed is None else options.randomSeed) - -for name in filenames: - damask.util.report(scriptName,name) - - if options.N > np.prod(grid): - damask.util.croak('More seeds than grid positions.') - sys.exit() - if options.selective and options.distance < min(size/grid): - damask.util.croak('Distance must be larger than grid spacing.') - sys.exit() - if options.selective and options.distance**3*options.N > 0.5*np.prod(size): - damask.util.croak('Number of seeds for given size and distance should be < {}.'\ - .format(int(0.5*np.prod(size)/options.distance**3))) - - eulers = np.random.rand(options.N,3) # create random Euler triplets - eulers[:,0] *= 360.0 # phi_1 is uniformly distributed - eulers[:,1] = np.degrees(np.arccos(2*eulers[:,1]-1.0)) # cos(Phi) is uniformly distributed - eulers[:,2] *= 360.0 # phi_2 is uniformly distributed - - - if not options.selective: - coords = damask.grid_filters.cell_coord0(grid,size).reshape(-1,3,order='F') - seeds = coords[np.random.choice(np.prod(grid), options.N, replace=False)] \ - + np.broadcast_to(size/grid,(options.N,3))*(np.random.rand(options.N,3)*.5-.25) # wobble without leaving grid - else: - seeds = np.empty((options.N,3)) - seeds[0] = np.random.random(3) * size - - i = 1 - progress = damask.util._ProgressBar(options.N,'',50) - while i < options.N: - candidates = np.random.rand(options.numCandidates,3)*np.broadcast_to(size,(options.numCandidates,3)) - tree = spatial.cKDTree(seeds[:i]) - distances, dev_null = tree.query(candidates) - best = distances.argmax() - if distances[best] > options.distance: # require minimum separation - seeds[i] = candidates[best] # maximum separation to existing point cloud - i += 1 - progress.update(i) - - - comments = [scriptID + ' ' + ' '.join(sys.argv[1:]), - 'grid\ta {}\tb {}\tc {}'.format(*grid), - 'size\tx {}\ty {}\tz {}'.format(*size), - 'randomSeed\t{}'.format(options.randomSeed), - ] - - table = damask.Table(np.hstack((seeds,eulers)),{'pos':(3,),'euler':(3,)},comments)\ - .add('microstructure',np.arange(options.microstructure,options.microstructure + options.N,dtype=int)) - - if options.weights: - weights = np.random.uniform(low = 0, high = options.max, size = options.N) if options.max > 0.0 \ - else np.random.normal(loc = options.mean, scale = options.sigma, size = options.N) - table = table.add('weight',weights) - - table.save(sys.stdout if name is None else name,legacy=True)