diff --git a/processing/legacy/seeds_fromDistribution.py b/processing/legacy/seeds_fromDistribution.py deleted file mode 100755 index 59d9f0213..000000000 --- a/processing/legacy/seeds_fromDistribution.py +++ /dev/null @@ -1,268 +0,0 @@ -#!/usr/bin/env python3 - -import threading -import time -import os -import sys -import random -from optparse import OptionParser -from io import StringIO - -import numpy as np - -import damask - -scriptName = os.path.splitext(os.path.basename(__file__))[0] -scriptID = ' '.join([scriptName,damask.version]) - -mismatch = None -currentSeedsName = None - -#--------------------------------------------------------------------------------------------------- -class myThread (threading.Thread): - """Perturb seed in seed file, performes Voronoi tessellation, evaluates, and updates best match.""" - - def __init__(self, threadID): - """Threading class with thread ID.""" - threading.Thread.__init__(self) - self.threadID = threadID - - def run(self): - global bestSeedsUpdate - global bestSeedsVFile - global nMaterials - global delta - global points - global target - global match - global baseFile - global maxSeeds - - s.acquire() - bestMatch = match - s.release() - - random.seed(options.randomSeed+self.threadID) # initializes to given seeds - knownSeedsUpdate = bestSeedsUpdate -1.0 # trigger update of local best seeds - randReset = True # aquire new direction - - myBestSeedsVFile = StringIO() # store local copy of best seeds file - perturbedSeedsVFile = StringIO() # perturbed best seeds file - -#--- still not matching desired bin class ---------------------------------------------------------- - while bestMatch < options.threshold: - s.acquire() # ensure only one thread acces global data - if bestSeedsUpdate > knownSeedsUpdate: # write best fit to virtual file - knownSeedsUpdate = bestSeedsUpdate - bestSeedsVFile.seek(0) - myBestSeedsVFile.close() - myBestSeedsVFile = StringIO() - i=0 - myBestSeedsVFile.writelines(bestSeedsVFile.readlines()) - s.release() - - if randReset: # new direction because current one led to worse fit - - randReset = False - - NmoveGrains = random.randrange(1,maxSeeds) - selectedMs = [] - direction = [] - for i in range(NmoveGrains): - selectedMs.append(random.randrange(1,nMaterials)) - - direction.append((np.random.random()-0.5)*delta) - - perturbedSeedsVFile.close() # reset virtual file - perturbedSeedsVFile = StringIO() - myBestSeedsVFile.seek(0) - - perturbedSeedsTable = damask.Table.load(myBestSeedsVFile) - coords = perturbedSeedsTable.get('pos') - i = 0 - for ms,coord in enumerate(coords): - if ms in selectedMs: - newCoords=coord+direction[i] - newCoords=np.where(newCoords>=1.0,newCoords-1.0,newCoords) # ensure that the seeds remain in the box - newCoords=np.where(newCoords <0.0,newCoords+1.0,newCoords) - coords[i]=newCoords - direction[i]*=2. - i+= 1 - perturbedSeedsTable.set('pos',coords).save(perturbedSeedsVFile,legacy=True) - -#--- do tesselation with perturbed seed file ------------------------------------------------------ - perturbedGeom = damask.Grid.from_Voronoi_tessellation(options.grid,np.ones(3),coords) - - -#--- evaluate current seeds file ------------------------------------------------------------------ - myNmaterials = len(np.unique(perturbedGeom.material)) - currentData = np.bincount(perturbedGeom.material.ravel())[1:]/points - currentError=[] - currentHist=[] - for i in range(nMaterials): # calculate the deviation in all bins per histogram - currentHist.append(np.histogram(currentData,bins=target[i]['bins'])[0]) - currentError.append(np.sqrt(np.square(np.array(target[i]['histogram']-currentHist[i])).sum())) - -# as long as not all grains are within the range of the target, use the deviation to left and right as error - if currentError[0]>0.0: - currentError[0] *=((target[0]['bins'][0]-np.min(currentData))**2.0+ - (target[0]['bins'][1]-np.max(currentData))**2.0)**0.5 # norm of deviations by number of usual bin deviation - s.acquire() # do the evaluation serially - bestMatch = match -#--- count bin classes with no mismatch ---------------------------------------------------------------------- - myMatch=0 - for i in range(nMaterials): - if currentError[i] > 0.0: break - myMatch = i+1 - - if myNmaterials == nMaterials: - for i in range(min(nMaterials,myMatch+options.bins)): - if currentError[i] > target[i]['error']: # worse fitting, next try - randReset = True - break - elif currentError[i] < target[i]['error']: # better fit - bestSeedsUpdate = time.time() # save time of better fit - damask.util.croak('Thread {:d}: Better match ({:d} bins, {:6.4f} --> {:6.4f})'\ - .format(self.threadID,i+1,target[i]['error'],currentError[i])) - damask.util.croak(' target: '+np.array_str(target[i]['histogram'])) - damask.util.croak(' best: '+np.array_str(currentHist[i])) - currentSeedsName = baseFile+'_'+str(bestSeedsUpdate).replace('.','-') # name of new seed file (use time as unique identifier) - perturbedSeedsVFile.seek(0) - bestSeedsVFile.close() - bestSeedsVFile = StringIO() - sys.stdout.flush() - with open(currentSeedsName+'.seeds','w') as currentSeedsFile: # write to new file - for line in perturbedSeedsVFile: - currentSeedsFile.write(line) - bestSeedsVFile.write(line) - for j in range(nMaterials): # save new errors for all bins - target[j]['error'] = currentError[j] - if myMatch > match: # one or more new bins have no deviation - damask.util.croak( 'Stage {:d} cleared'.format(myMatch)) - match=myMatch - sys.stdout.flush() - break - if i == min(nMaterials,myMatch+options.bins)-1: # same quality as before: take it to keep on moving - bestSeedsUpdate = time.time() - perturbedSeedsVFile.seek(0) - bestSeedsVFile.close() - bestSeedsVFile = StringIO() - bestSeedsVFile.writelines(perturbedSeedsVFile.readlines()) - for j in range(nMaterials): - target[j]['error'] = currentError[j] - randReset = True - else: #--- not all grains are tessellated - damask.util.croak('Thread {:d}: Material mismatch ({:d} material indices mapped)'\ - .format(self.threadID,myNmaterials)) - randReset = True - - - s.release() - - -# -------------------------------------------------------------------- -# MAIN -# -------------------------------------------------------------------- - -parser = OptionParser(usage='%prog options [file[s]]', description = """ -Monte Carlo simulation to produce seed file that gives same size distribution like given geometry file. - -""", version = scriptID) - -parser.add_option('-s','--seeds', dest='seedFile', metavar='string', - help='name of the intial seed file. If not found, a new one is created [%default]') -parser.add_option('-g','--grid', dest='grid', type='int', nargs=3, metavar='int int int', - help='a,b,c grid of hexahedral box [%default]') -parser.add_option('-t','--threads', dest='threads', type='int', metavar='int', - help='number of parallel executions [%default]') -parser.add_option('-r', '--rnd', dest='randomSeed', type='int', metavar='int', - help='seed of random number generator [%default]') -parser.add_option('--target', dest='target', metavar='string', - help='name of the geom file with target distribution [%default]') -parser.add_option('--tolerance', dest='threshold', type='int', metavar='int', - help='stopping criterion (bin number) [%default]') -parser.add_option('--scale', dest='scale',type='float', metavar='float', - help='maximum moving distance of perturbed seed in pixel [%default]') -parser.add_option('--bins', dest='bins', type='int', metavar='int', - help='bins to sort beyond current best fit [%default]') -parser.add_option('--maxseeds', dest='maxseeds', type='int', metavar='int', - help='maximum number of seeds to move simulateneously [number of seeds]') - -parser.set_defaults(seedFile = 'seeds', - grid = (64,64,64), - threads = 2, - randomSeed = None, - target = 'geom', - threshold = 20, - bins = 15, - scale = 1.0, - maxseeds = 0) - -options = parser.parse_args()[0] - -damask.util.report(scriptName,options.seedFile) - -if options.randomSeed is None: - options.randomSeed = int(os.urandom(4).hex(),16) -damask.util.croak(options.randomSeed) -delta = options.scale/np.array(options.grid) -baseFile = os.path.splitext(os.path.basename(options.seedFile))[0] -points = np.array(options.grid).prod().astype('float') - -# ----------- calculate target distribution and bin edges -targetGeom = damask.Grid.load_ASCII(os.path.splitext(os.path.basename(options.target))[0]+'.geom') -nMaterials = len(np.unique(targetGeom.material)) -targetVolFrac = np.bincount(targetGeom.material.flatten())/targetGeom.cells.prod().astype(np.float) -target = [] -for i in range(1,nMaterials+1): - targetHist,targetBins = np.histogram(targetVolFrac,bins=i) #bin boundaries - target.append({'histogram':targetHist,'bins':targetBins}) - -# ----------- create initial seed file or open existing one -bestSeedsVFile = StringIO() -if os.path.isfile(os.path.splitext(options.seedFile)[0]+'.seeds'): - initial_seeds = damask.Table.load(os.path.splitext(options.seedFile)[0]+'.seeds').get('pos') -else: - initial_seeds = damask.seeds.from_random(np.ones(3),nMaterials,options.grid,options.randomSeed) - -bestSeedsUpdate = time.time() - -# ----------- tessellate initial seed file to get and evaluate geom file -bestSeedsVFile.seek(0) -initialGeom = damask.Grid.from_Voronoi_tessellation(options.grid,np.ones(3),initial_seeds) - -if len(np.unique(targetGeom.material)) != nMaterials: - damask.util.croak('error. Material count mismatch') - -initialData = np.bincount(initialGeom.material.flatten())/points -for i in range(nMaterials): - initialHist = np.histogram(initialData,bins=target[i]['bins'])[0] - target[i]['error']=np.sqrt(np.square(np.array(target[i]['histogram']-initialHist)).sum()) - -# as long as not all grain sizes are within the range, the error is the deviation to left and right -if target[0]['error'] > 0.0: - target[0]['error'] *=((target[0]['bins'][0]-np.min(initialData))**2.0+ - (target[0]['bins'][1]-np.max(initialData))**2.0)**0.5 -match=0 -for i in range(nMaterials): - if target[i]['error'] > 0.0: break - match = i+1 - - -if options.maxseeds < 1: - maxSeeds = len(np.unique(initialGeom.material)) -else: - maxSeeds = options.maxseeds - -if match >0: damask.util.croak('Stage {:d} cleared'.format(match)) -sys.stdout.flush() - -# start mulithreaded monte carlo simulation -threads = [] -s = threading.Semaphore(1) - -for i in range(options.threads): - threads.append(myThread(i)) - threads[i].start() -for i in range(options.threads): - threads[i].join() diff --git a/python/damask/_rotation.py b/python/damask/_rotation.py index b68e25d49..4eb52eb87 100644 --- a/python/damask/_rotation.py +++ b/python/damask/_rotation.py @@ -569,6 +569,7 @@ class Rotation: array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]) + """ return Rotation._qu2om(self.quaternion)