From c6863a612407b2e746ee54bb5638e27706b8b6fa Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 15 May 2019 23:32:23 +0200 Subject: [PATCH 01/29] also consider homogenization/materialpoint results --- python/damask/dadf5.py | 23 +++++++++++++++++++++-- 1 file changed, 21 insertions(+), 2 deletions(-) diff --git a/python/damask/dadf5.py b/python/damask/dadf5.py index 841029af8..eed5907e5 100644 --- a/python/damask/dadf5.py +++ b/python/damask/dadf5.py @@ -48,16 +48,24 @@ class DADF5(): for o in f['inc{:05}/constituent/{}'.format(self.increments[0]['inc'],c)].keys(): self.c_output_types.append(o) self.c_output_types = list(set(self.c_output_types)) # make unique + + self.m_output_types = [] + for m in self.materialpoints: + for o in f['inc{:05}/materialpoint/{}'.format(self.increments[0]['inc'],m)].keys(): + self.m_output_types.append(o) + self.m_output_types = list(set(self.m_output_types)) # make unique self.active= {'increments': self.increments, 'constituents': self.constituents, 'materialpoints': self.materialpoints, 'constituent': self.Nconstituents, - 'c_output_types': self.c_output_types} + 'c_output_types': self.c_output_types, + 'm_output_types': self.m_output_types} self.filename = filename self.mode = mode - + + def list_data(self): """Shows information on all datasets in the file""" with h5py.File(self.filename,'r') as f: @@ -73,6 +81,17 @@ class DADF5(): print(' {} ({})'.format(x,f[group_output_types+'/'+x].attrs['Description'].decode())) except: pass + for m in self.active['materialpoints']: + print('\n'+m) + group_materialpoint = group_inc+'/materialpoint/'+m + for t in self.active['m_output_types']: + print(' {}'.format(t)) + group_output_types = group_materialpoint+'/'+t + try: + for x in f[group_output_types].keys(): + print(' {} ({})'.format(x,f[group_output_types+'/'+x].attrs['Description'].decode())) + except: + pass def get_dataset_location(self,label): From 3c0c0a2cd10d38c472b92d84cb687f2268653b02 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 00:27:06 +0200 Subject: [PATCH 02/29] more flexible in selecting data --- processing/post/DADF5_postResults.py | 28 ++++++++++++++++++++++++---- processing/post/DADF5_vtk_cells.py | 9 +++++++-- python/damask/dadf5.py | 26 +++++++++++++++++++++----- 3 files changed, 52 insertions(+), 11 deletions(-) diff --git a/processing/post/DADF5_postResults.py b/processing/post/DADF5_postResults.py index caf09d536..5ad49604b 100755 --- a/processing/post/DADF5_postResults.py +++ b/processing/post/DADF5_postResults.py @@ -22,10 +22,15 @@ parser.add_argument('filenames', nargs='+', help='DADF5 files') parser.add_argument('-d','--dir', dest='dir',default='postProc',metavar='string', help='name of subdirectory to hold output') +parser.add_argument('--mat', nargs='+', + help='labels for materialpoint/homogenization',dest='mat') +parser.add_argument('--con', nargs='+', + help='labels for constituent/crystallite/constitutive',dest='con') options = parser.parse_args() -options.labels = ['Fe','Fp','xi_sl'] +if options.mat is None: options.mat=[] +if options.con is None: options.con=[] # --- loop over input files ------------------------------------------------------------------------ @@ -57,7 +62,7 @@ for filename in options.filenames: header+=' 1_pos 2_pos 3_pos' results.active['increments'] = [inc] - for label in options.labels: + for label in options.con: for o in results.c_output_types: results.active['c_output_types'] = [o] for c in results.constituents: @@ -67,12 +72,27 @@ for filename in options.filenames: continue label = x[0].split('/')[-1] array = results.read_dataset(x,0) - d = np.product(np.shape(array)[1:]) + d = int(np.product(np.shape(array)[1:])) + array = np.reshape(array,[np.product(results.grid),d]) + data = np.concatenate((data,array),1) + + header+= ''.join([' {}_{}'.format(j+1,label) for j in range(d)]) + + for label in options.mat: + for o in results.m_output_types: + results.active['m_output_types'] = [o] + for m in results.materialpoints: + results.active['materialpoints'] = [m] + x = results.get_dataset_location(label) + if len(x) == 0: + continue + label = x[0].split('/')[-1] + array = results.read_dataset(x,0) + d = int(np.product(np.shape(array)[1:])) array = np.reshape(array,[np.product(results.grid),d]) data = np.concatenate((data,array),1) header+= ''.join([' {}_{}'.format(j+1,label) for j in range(d)]) - dirname = os.path.abspath(os.path.join(os.path.dirname(filename),options.dir)) try: diff --git a/processing/post/DADF5_vtk_cells.py b/processing/post/DADF5_vtk_cells.py index dc1177488..0e79209e2 100755 --- a/processing/post/DADF5_vtk_cells.py +++ b/processing/post/DADF5_vtk_cells.py @@ -23,10 +23,15 @@ parser.add_argument('filenames', nargs='+', help='DADF5 files') parser.add_argument('-d','--dir', dest='dir',default='postProc',metavar='string', help='name of subdirectory to hold output') +parser.add_argument('--mat', nargs='+', + help='labels for materialpoint/homogenization',dest='mat') +parser.add_argument('--con', nargs='+', + help='labels for constituent/crystallite/constitutive',dest='con') options = parser.parse_args() -options.labels = ['Fe','Fp','xi_sl'] +if options.mat is None: options.mat=[] +if options.con is None: options.con=[] # --- loop over input files ------------------------------------------------------------------------ @@ -54,7 +59,7 @@ for filename in options.filenames: print('Output step {}/{}'.format(i+1,len(results.increments))) vtk_data = [] results.active['increments'] = [inc] - for label in options.labels: + for label in options.con: for o in results.c_output_types: results.active['c_output_types'] = [o] if o != 'generic': diff --git a/python/damask/dadf5.py b/python/damask/dadf5.py index eed5907e5..56082a1f4 100644 --- a/python/damask/dadf5.py +++ b/python/damask/dadf5.py @@ -82,7 +82,6 @@ class DADF5(): except: pass for m in self.active['materialpoints']: - print('\n'+m) group_materialpoint = group_inc+'/materialpoint/'+m for t in self.active['m_output_types']: print(' {}'.format(t)) @@ -108,6 +107,14 @@ class DADF5(): path.append(group_constituent+'/'+t+'/'+label) except: pass + for m in self.active['materialpoints']: + group_materialpoint = group_inc+'/materialpoint/'+m + for t in self.active['m_output_types']: + try: + f[group_materialpoint+'/'+t+'/'+label] + path.append(group_materialpoint+'/'+t+'/'+label) + except: + pass return path @@ -122,10 +129,19 @@ class DADF5(): dataset = np.full(shape,np.nan) for pa in path: label = pa.split('/')[2] - p = np.where(f['mapping/cellResults/constituent'][:,c]['Name'] == str.encode(label))[0] - u = (f['mapping/cellResults/constituent'][p,c]['Position']) - dataset[p,:] = f[pa][u,:] - + try: + p = np.where(f['mapping/cellResults/constituent'][:,c]['Name'] == str.encode(label))[0] + u = (f['mapping/cellResults/constituent'][p,c]['Position']) + dataset[p,:] = f[pa][u,:] # does not work for scalar datasets + except: + pass + try: + p = np.where(f['mapping/cellResults/materialpoint']['Name'] == str.encode(label))[0] + u = (f['mapping/cellResults/materialpoint'][p.tolist()]['Position']) + dataset[p,:] = f[pa][u,:] # does not work for scalar datasets + except: + pass + return dataset From 4599d1c34e219d9c1ab3d49971d4d37ab65049e3 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 00:42:22 +0200 Subject: [PATCH 03/29] does not match node does not make sense, is weirdly numbered --- processing/post/DADF5_postResults.py | 3 --- 1 file changed, 3 deletions(-) diff --git a/processing/post/DADF5_postResults.py b/processing/post/DADF5_postResults.py index 5ad49604b..b918388d9 100755 --- a/processing/post/DADF5_postResults.py +++ b/processing/post/DADF5_postResults.py @@ -53,9 +53,6 @@ for filename in options.filenames: data = np.array([inc['inc'] for j in range(np.product(results.grid))]).reshape([np.product(results.grid),1]) header+= 'inc' - - data = np.concatenate((data,np.array([j+1 for j in range(np.product(results.grid))]).reshape([np.product(results.grid),1])),1) - header+=' node' coords = coords.reshape([np.product(results.grid),3]) data = np.concatenate((data,coords),1) From 9f7fa5393a8af1d784433262dfbbe8bebded5ab0 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 09:31:13 +0200 Subject: [PATCH 04/29] fix for scalar datasets --- processing/post/DADF5_postResults.py | 12 +++++++++--- python/damask/dadf5.py | 11 +++++++++-- 2 files changed, 18 insertions(+), 5 deletions(-) diff --git a/processing/post/DADF5_postResults.py b/processing/post/DADF5_postResults.py index b918388d9..1af8c2aa9 100755 --- a/processing/post/DADF5_postResults.py +++ b/processing/post/DADF5_postResults.py @@ -73,8 +73,11 @@ for filename in options.filenames: array = np.reshape(array,[np.product(results.grid),d]) data = np.concatenate((data,array),1) - header+= ''.join([' {}_{}'.format(j+1,label) for j in range(d)]) - + if d>1: + header+= ''.join([' {}_{}'.format(j+1,label) for j in range(d)]) + else: + header+=' '+label + for label in options.mat: for o in results.m_output_types: results.active['m_output_types'] = [o] @@ -89,7 +92,10 @@ for filename in options.filenames: array = np.reshape(array,[np.product(results.grid),d]) data = np.concatenate((data,array),1) - header+= ''.join([' {}_{}'.format(j+1,label) for j in range(d)]) + if d>1: + header+= ''.join([' {}_{}'.format(j+1,label) for j in range(d)]) + else: + header+=' '+label dirname = os.path.abspath(os.path.join(os.path.dirname(filename),options.dir)) try: diff --git a/python/damask/dadf5.py b/python/damask/dadf5.py index 56082a1f4..0abde662b 100644 --- a/python/damask/dadf5.py +++ b/python/damask/dadf5.py @@ -126,19 +126,26 @@ class DADF5(): """ with h5py.File(self.filename,'r') as f: shape = (self.Nmaterialpoints,) + np.shape(f[path[0]])[1:] + if len(shape) == 1: shape = shape +(1,) dataset = np.full(shape,np.nan) for pa in path: label = pa.split('/')[2] try: p = np.where(f['mapping/cellResults/constituent'][:,c]['Name'] == str.encode(label))[0] u = (f['mapping/cellResults/constituent'][p,c]['Position']) - dataset[p,:] = f[pa][u,:] # does not work for scalar datasets + a = np.array(f[pa]) + if len(a.shape) == 1: + a=a.reshape([a.shape[0],1]) + dataset except: pass try: p = np.where(f['mapping/cellResults/materialpoint']['Name'] == str.encode(label))[0] u = (f['mapping/cellResults/materialpoint'][p.tolist()]['Position']) - dataset[p,:] = f[pa][u,:] # does not work for scalar datasets + a = np.array(f[pa]) + if len(a.shape) == 1: + a=a.reshape([a.shape[0],1]) + dataset[p,:] = a[u,:] except: pass From f2268d055ff9889e11aeafbc8d0a6c34e247a75f Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 09:42:37 +0200 Subject: [PATCH 05/29] first test relying on DADF5 --- PRIVATE | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PRIVATE b/PRIVATE index d96bfb329..df0b7cc97 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit d96bfb32920c96a8a43958f76a209d34c6bd841a +Subproject commit df0b7cc97963cd712f6e33397937b187635c99f4 From 84a82d0878f442067ac6c3fba30276c0c7523f21 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 09:46:50 +0200 Subject: [PATCH 06/29] rely on working HDF5 output for RGC homog --- PRIVATE | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PRIVATE b/PRIVATE index df0b7cc97..8cec21296 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit df0b7cc97963cd712f6e33397937b187635c99f4 +Subproject commit 8cec2129617f4a206193671e0e1070965c7c2e53 From 0958c4bb88bc441e598a460eb785344ca9edf9dd Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 09:48:30 +0200 Subject: [PATCH 07/29] phase out postResults --- PRIVATE | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PRIVATE b/PRIVATE index 8cec21296..639c6f4a5 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 8cec2129617f4a206193671e0e1070965c7c2e53 +Subproject commit 639c6f4a5eafc893c83c740c57f417eaaabc45ae From 39a75c201545327d62a37df8f3b8f89b21cf3a32 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 10:23:23 +0200 Subject: [PATCH 08/29] phasing out postResults starting with RGC because it is rarely used and removing it here allows to go ahead with the submodule structure for homogenization --- src/homogenization.f90 | 44 +---------------- src/homogenization_mech_RGC.f90 | 88 +++------------------------------ 2 files changed, 10 insertions(+), 122 deletions(-) diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 77a2e004f..1ff036c3a 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -120,7 +120,6 @@ subroutine homogenization_init character(len=32) :: outputName !< name of output, intermediate fix until HDF5 output is ready logical :: valid - if (any(homogenization_type == HOMOGENIZATION_NONE_ID)) call mech_none_init if (any(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID)) call mech_isostrain_init if (any(homogenization_type == HOMOGENIZATION_RGC_ID)) call homogenization_RGC_init @@ -139,35 +138,9 @@ subroutine homogenization_init call IO_write_jobFile(FILEUNIT,'outputHomogenization') do p = 1,size(config_homogenization) if (any(material_homogenizationAt == p)) then - i = homogenization_typeInstance(p) ! which instance of this homogenization type - valid = .true. ! assume valid - select case(homogenization_type(p)) ! split per homogenization type - case (HOMOGENIZATION_NONE_ID) - outputName = HOMOGENIZATION_NONE_label - thisOutput => null() - thisSize => null() - case (HOMOGENIZATION_ISOSTRAIN_ID) - outputName = HOMOGENIZATION_ISOSTRAIN_label - thisOutput => null() - thisSize => null() - case (HOMOGENIZATION_RGC_ID) - outputName = HOMOGENIZATION_RGC_label - thisOutput => homogenization_RGC_output - thisSize => homogenization_RGC_sizePostResult - case default - valid = .false. - end select write(FILEUNIT,'(/,a,/)') '['//trim(homogenization_name(p))//']' - if (valid) then - write(FILEUNIT,'(a)') '(type)'//char(9)//trim(outputName) - write(FILEUNIT,'(a,i4)') '(ngrains)'//char(9),homogenization_Ngrains(p) - if (homogenization_type(p) /= HOMOGENIZATION_NONE_ID .and. & - homogenization_type(p) /= HOMOGENIZATION_ISOSTRAIN_ID) then - do e = 1,size(thisOutput(:,i)) - write(FILEUNIT,'(a,i4)') trim(thisOutput(e,i))//char(9),thisSize(e,i) - enddo - endif - endif + write(FILEUNIT,'(a,i4)') '(ngrains)'//char(9),homogenization_Ngrains(p) + i = thermal_typeInstance(p) ! which instance of this thermal type valid = .true. ! assume valid select case(thermal_type(p)) ! split per thermal type @@ -837,8 +810,6 @@ end subroutine averageStressAndItsTangent function postResults(ip,el) use mesh, only: & mesh_element - use homogenization_mech_RGC, only: & - homogenization_RGC_postResults use thermal_adiabatic, only: & thermal_adiabatic_postResults use thermal_conduction, only: & @@ -861,17 +832,6 @@ function postResults(ip,el) postResults = 0.0_pReal - startPos = 1 - endPos = homogState(material_homogenizationAt(el))%sizePostResults - chosenHomogenization: select case (homogenization_type(mesh_element(3,el))) - - case (HOMOGENIZATION_RGC_ID) chosenHomogenization - instance = homogenization_typeInstance(material_homogenizationAt(el)) - of = mappingHomogenization(1,ip,el) - postResults(startPos:endPos) = homogenization_RGC_postResults(instance,of) - - end select chosenHomogenization - startPos = endPos + 1 endPos = endPos + thermalState(material_homogenizationAt(el))%sizePostResults chosenThermal: select case (thermal_type(mesh_element(3,el))) diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index d3feac1eb..d7b1b31bf 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -12,10 +12,6 @@ module homogenization_mech_RGC implicit none private - integer, dimension(:,:), allocatable,target, public :: & - homogenization_RGC_sizePostResult - character(len=64), dimension(:,:), allocatable,target, public :: & - homogenization_RGC_output ! name of each post result output enum, bind(c) enumerator :: & @@ -28,7 +24,7 @@ module homogenization_mech_RGC magnitudemismatch_ID end enum - type, private :: tParameters + type :: tParameters integer, dimension(:), allocatable :: & Nconstituents real(pReal) :: & @@ -43,7 +39,7 @@ module homogenization_mech_RGC outputID end type tParameters - type, private :: tRGCstate + type :: tRGCstate real(pReal), pointer, dimension(:) :: & work, & penaltyEnergy @@ -51,7 +47,7 @@ module homogenization_mech_RGC relaxationVector end type tRGCstate - type, private :: tRGCdependentState + type :: tRGCdependentState real(pReal), allocatable, dimension(:) :: & volumeDiscrepancy, & relaxationRate_avg, & @@ -62,12 +58,12 @@ module homogenization_mech_RGC orientation end type tRGCdependentState - type(tparameters), dimension(:), allocatable, private :: & + type(tparameters), dimension(:), allocatable :: & param - type(tRGCstate), dimension(:), allocatable, private :: & + type(tRGCstate), dimension(:), allocatable :: & state, & state0 - type(tRGCdependentState), dimension(:), allocatable, private :: & + type(tRGCdependentState), dimension(:), allocatable :: & dependentState public :: & @@ -75,16 +71,7 @@ module homogenization_mech_RGC homogenization_RGC_partitionDeformation, & homogenization_RGC_averageStressAndItsTangent, & homogenization_RGC_updateState, & - homogenization_RGC_postResults, & mech_RGC_results ! name suited for planned submodule situation - private :: & - relaxationVector, & - interfaceNormal, & - getInterface, & - grain1to3, & - grain3to1, & - interface4to1, & - interface1to4 contains @@ -111,7 +98,7 @@ subroutine homogenization_RGC_init() integer :: & Ninstance, & h, i, & - NofMyHomog, outputSize, & + NofMyHomog, & sizeState, nIntFaceTot character(len=65536), dimension(0), parameter :: emptyStringArray = [character(len=65536)::] @@ -139,9 +126,6 @@ subroutine homogenization_RGC_init() allocate(state0(Ninstance)) allocate(dependentState(Ninstance)) - allocate(homogenization_RGC_sizePostResult(maxval(homogenization_Noutput),Ninstance),source=0) - allocate(homogenization_RGC_output(maxval(homogenization_Noutput),Ninstance)) - homogenization_RGC_output='' do h = 1, size(homogenization_type) if (homogenization_type(h) /= HOMOGENIZATION_RGC_ID) cycle @@ -176,28 +160,20 @@ subroutine homogenization_RGC_init() case('constitutivework') outputID = constitutivework_ID - outputSize = 1 case('penaltyenergy') outputID = penaltyenergy_ID - outputSize = 1 case('volumediscrepancy') outputID = volumediscrepancy_ID - outputSize = 1 case('averagerelaxrate') outputID = averagerelaxrate_ID - outputSize = 1 case('maximumrelaxrate') outputID = maximumrelaxrate_ID - outputSize = 1 case('magnitudemismatch') outputID = magnitudemismatch_ID - outputSize = 3 end select if (outputID /= undefined_ID) then - homogenization_RGC_output(i,homogenization_typeInstance(h)) = outputs(i) - homogenization_RGC_sizePostResult(i,homogenization_typeInstance(h)) = outputSize prm%outputID = [prm%outputID , outputID] endif @@ -211,7 +187,7 @@ subroutine homogenization_RGC_init() + size(['avg constitutive work ','average penalty energy']) homogState(h)%sizeState = sizeState - homogState(h)%sizePostResults = sum(homogenization_RGC_sizePostResult(:,homogenization_typeInstance(h))) + homogState(h)%sizePostResults = 0 allocate(homogState(h)%state0 (sizeState,NofMyHomog), source=0.0_pReal) allocate(homogState(h)%subState0(sizeState,NofMyHomog), source=0.0_pReal) allocate(homogState(h)%state (sizeState,NofMyHomog), source=0.0_pReal) @@ -1033,54 +1009,6 @@ subroutine homogenization_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF, end subroutine homogenization_RGC_averageStressAndItsTangent -!-------------------------------------------------------------------------------------------------- -!> @brief return array of homogenization results for post file inclusion -!-------------------------------------------------------------------------------------------------- -pure function homogenization_RGC_postResults(instance,of) result(postResults) - - integer, intent(in) :: & - instance, & - of - - integer :: & - o,c - real(pReal), dimension(sum(homogenization_RGC_sizePostResult(:,instance))) :: & - postResults - - associate(stt => state(instance), dst => dependentState(instance), prm => param(instance)) - - c = 0 - - outputsLoop: do o = 1,size(prm%outputID) - select case(prm%outputID(o)) - - case (constitutivework_ID) - postResults(c+1) = stt%work(of) - c = c + 1 - case (magnitudemismatch_ID) - postResults(c+1:c+3) = dst%mismatch(1:3,of) - c = c + 3 - case (penaltyenergy_ID) - postResults(c+1) = stt%penaltyEnergy(of) - c = c + 1 - case (volumediscrepancy_ID) - postResults(c+1) = dst%volumeDiscrepancy(of) - c = c + 1 - case (averagerelaxrate_ID) - postResults(c+1) = dst%relaxationrate_avg(of) - c = c + 1 - case (maximumrelaxrate_ID) - postResults(c+1) = dst%relaxationrate_max(of) - c = c + 1 - end select - - enddo outputsLoop - - end associate - -end function homogenization_RGC_postResults - - !-------------------------------------------------------------------------------------------------- !> @brief writes results to HDF5 output file ! ToDo: check wheter units are correct From 339b86f784bdb31c1359b486203e12d1e1d0717a Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 11:44:03 +0200 Subject: [PATCH 09/29] bugfix + more verbose reporting --- processing/post/DADF5_postResults.py | 2 +- processing/post/DADF5_vtk_cells.py | 2 ++ python/damask/dadf5.py | 23 +++++++++++++---------- 3 files changed, 16 insertions(+), 11 deletions(-) diff --git a/processing/post/DADF5_postResults.py b/processing/post/DADF5_postResults.py index 1af8c2aa9..fa47805bb 100755 --- a/processing/post/DADF5_postResults.py +++ b/processing/post/DADF5_postResults.py @@ -78,7 +78,7 @@ for filename in options.filenames: else: header+=' '+label - for label in options.mat: + for label in options.mat: for o in results.m_output_types: results.active['m_output_types'] = [o] for m in results.materialpoints: diff --git a/processing/post/DADF5_vtk_cells.py b/processing/post/DADF5_vtk_cells.py index 0e79209e2..ef27a973c 100755 --- a/processing/post/DADF5_vtk_cells.py +++ b/processing/post/DADF5_vtk_cells.py @@ -59,7 +59,9 @@ for filename in options.filenames: print('Output step {}/{}'.format(i+1,len(results.increments))) vtk_data = [] results.active['increments'] = [inc] + for label in options.con: + for o in results.c_output_types: results.active['c_output_types'] = [o] if o != 'generic': diff --git a/python/damask/dadf5.py b/python/damask/dadf5.py index 0abde662b..4aa31ec56 100644 --- a/python/damask/dadf5.py +++ b/python/damask/dadf5.py @@ -99,22 +99,25 @@ class DADF5(): with h5py.File(self.filename,'r') as f: for i in self.active['increments']: group_inc = 'inc{:05}'.format(i['inc']) + for c in self.active['constituents']: group_constituent = group_inc+'/constituent/'+c for t in self.active['c_output_types']: try: f[group_constituent+'/'+t+'/'+label] path.append(group_constituent+'/'+t+'/'+label) - except: - pass - for m in self.active['materialpoints']: + except Exception as e: + print('unable to locate constituents dataset: '+ str(e)) + + for m in []: #self.active['materialpoints']: group_materialpoint = group_inc+'/materialpoint/'+m for t in self.active['m_output_types']: try: f[group_materialpoint+'/'+t+'/'+label] path.append(group_materialpoint+'/'+t+'/'+label) - except: - pass + except Exception as e: + print('unable to locate materialpoints dataset: '+ str(e)) + return path @@ -136,9 +139,9 @@ class DADF5(): a = np.array(f[pa]) if len(a.shape) == 1: a=a.reshape([a.shape[0],1]) - dataset - except: - pass + dataset[p,:] = a[u,:] + except Exception as e: + print('unable to read constituent: '+ str(e)) try: p = np.where(f['mapping/cellResults/materialpoint']['Name'] == str.encode(label))[0] u = (f['mapping/cellResults/materialpoint'][p.tolist()]['Position']) @@ -146,8 +149,8 @@ class DADF5(): if len(a.shape) == 1: a=a.reshape([a.shape[0],1]) dataset[p,:] = a[u,:] - except: - pass + except Exception as e: + print('unable to read materialpoint: '+ str(e)) return dataset From f8b335a3a485a1269e93670cd51e8f5db46d5c78 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 18:24:54 +0200 Subject: [PATCH 10/29] loop (forall) over integration points wrong this was done for each integration point, but this was not detected for the forall loop --- src/homogenization.f90 | 119 ++++++++++++++++++++--------------------- 1 file changed, 59 insertions(+), 60 deletions(-) diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 77a2e004f..605cf5094 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -375,43 +375,46 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) ! initialize restoration points of ... do e = FEsolving_execElem(1),FEsolving_execElem(2) myNgrains = homogenization_Ngrains(mesh_element(3,e)) - do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do g = 1,myNgrains + do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); + do g = 1,myNgrains + + plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = & + plasticState (phaseAt(g,i,e))%state0( :,phasememberAt(g,i,e)) + do mySource = 1, phase_Nsources(phaseAt(g,i,e)) + sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) = & + sourceState(phaseAt(g,i,e))%p(mySource)%state0( :,phasememberAt(g,i,e)) + enddo + + crystallite_partionedFp0(1:3,1:3,g,i,e) = crystallite_Fp0(1:3,1:3,g,i,e) + crystallite_partionedLp0(1:3,1:3,g,i,e) = crystallite_Lp0(1:3,1:3,g,i,e) + crystallite_partionedFi0(1:3,1:3,g,i,e) = crystallite_Fi0(1:3,1:3,g,i,e) + crystallite_partionedLi0(1:3,1:3,g,i,e) = crystallite_Li0(1:3,1:3,g,i,e) + crystallite_partionedF0(1:3,1:3,g,i,e) = crystallite_F0(1:3,1:3,g,i,e) + crystallite_partionedS0(1:3,1:3,g,i,e) = crystallite_S0(1:3,1:3,g,i,e) - plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = & - plasticState (phaseAt(g,i,e))%state0( :,phasememberAt(g,i,e)) - do mySource = 1, phase_Nsources(phaseAt(g,i,e)) - sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) = & - sourceState(phaseAt(g,i,e))%p(mySource)%state0( :,phasememberAt(g,i,e)) enddo - crystallite_partionedFp0(1:3,1:3,g,i,e) = crystallite_Fp0(1:3,1:3,g,i,e) ! ...plastic def grads - crystallite_partionedLp0(1:3,1:3,g,i,e) = crystallite_Lp0(1:3,1:3,g,i,e) ! ...plastic velocity grads - crystallite_partionedFi0(1:3,1:3,g,i,e) = crystallite_Fi0(1:3,1:3,g,i,e) ! ...intermediate def grads - crystallite_partionedLi0(1:3,1:3,g,i,e) = crystallite_Li0(1:3,1:3,g,i,e) ! ...intermediate velocity grads - crystallite_partionedF0(1:3,1:3,g,i,e) = crystallite_F0(1:3,1:3,g,i,e) ! ...def grads - crystallite_partionedS0(1:3,1:3,g,i,e) = crystallite_S0(1:3,1:3,g,i,e) ! ...2nd PK stress - enddo; enddo - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e)) - materialpoint_subF0(1:3,1:3,i,e) = materialpoint_F0(1:3,1:3,i,e) ! ...def grad + materialpoint_subF0(1:3,1:3,i,e) = materialpoint_F0(1:3,1:3,i,e) materialpoint_subFrac(i,e) = 0.0_pReal materialpoint_subStep(i,e) = 1.0_pReal/subStepSizeHomog ! <> materialpoint_converged(i,e) = .false. ! pretend failed step of twice the required size materialpoint_requested(i,e) = .true. ! everybody requires calculation - endforall - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - homogState(material_homogenizationAt(e))%sizeState > 0) & - homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - homogState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal homogenization state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - thermalState(material_homogenizationAt(e))%sizeState > 0) & - thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - thermalState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal thermal state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - damageState(material_homogenizationAt(e))%sizeState > 0) & - damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - damageState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal damage state + + if (homogState(material_homogenizationAt(e))%sizeState > 0) & + homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & + homogState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal homogenization state + + if (thermalState(material_homogenizationAt(e))%sizeState > 0) & + thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & + thermalState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal thermal state + + if (damageState(material_homogenizationAt(e))%sizeState > 0) & + damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & + damageState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal damage state + enddo enddo + NiterationHomog = 0 cutBackLooping: do while (.not. terminallyIll .and. & @@ -422,7 +425,7 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) myNgrains = homogenization_Ngrains(mesh_element(3,e)) IpLooping1: do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) - converged: if ( materialpoint_converged(i,e) ) then + converged: if (materialpoint_converged(i,e)) then #ifdef DEBUG if (iand(debug_level(debug_homogenization), debug_levelExtensive) /= 0 & .and. ((e == debug_e .and. i == debug_i) & @@ -443,22 +446,22 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) ! wind forward grain starting point of... crystallite_partionedF0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedF(1:3,1:3,1:myNgrains,i,e) ! ...def grads + crystallite_partionedF(1:3,1:3,1:myNgrains,i,e) crystallite_partionedFp0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_Fp (1:3,1:3,1:myNgrains,i,e) ! ...plastic def grads + crystallite_Fp (1:3,1:3,1:myNgrains,i,e) crystallite_partionedLp0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_Lp (1:3,1:3,1:myNgrains,i,e) ! ...plastic velocity grads + crystallite_Lp (1:3,1:3,1:myNgrains,i,e) crystallite_partionedFi0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_Fi (1:3,1:3,1:myNgrains,i,e) ! ...intermediate def grads + crystallite_Fi (1:3,1:3,1:myNgrains,i,e) crystallite_partionedLi0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_Li (1:3,1:3,1:myNgrains,i,e) ! ...intermediate velocity grads + crystallite_Li (1:3,1:3,1:myNgrains,i,e) crystallite_partionedS0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_S (1:3,1:3,1:myNgrains,i,e) ! ...2nd PK stress + crystallite_S (1:3,1:3,1:myNgrains,i,e) do g = 1,myNgrains plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = & @@ -469,23 +472,22 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) enddo enddo - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - homogState(material_homogenizationAt(e))%sizeState > 0) & + if(homogState(material_homogenizationAt(e))%sizeState > 0) & homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - homogState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) ! ...internal homogenization state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - thermalState(material_homogenizationAt(e))%sizeState > 0) & + homogState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) + if(thermalState(material_homogenizationAt(e))%sizeState > 0) & thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - thermalState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) ! ...internal thermal state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - damageState(material_homogenizationAt(e))%sizeState > 0) & + thermalState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) + if(damageState(material_homogenizationAt(e))%sizeState > 0) & damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - damageState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) ! ...internal damage state - materialpoint_subF0(1:3,1:3,i,e) = materialpoint_subF(1:3,1:3,i,e) ! ...def grad + damageState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) + + materialpoint_subF0(1:3,1:3,i,e) = materialpoint_subF(1:3,1:3,i,e) + endif steppingNeeded else converged - if ( (myNgrains == 1 .and. materialpoint_subStep(i,e) <= 1.0 ) .or. & ! single grain already tried internal subStepping in crystallite + if ( (myNgrains == 1 .and. materialpoint_subStep(i,e) <= 1.0 ) .or. & ! single grain already tried internal subStepping in crystallite subStepSizeHomog * materialpoint_subStep(i,e) <= subStepMinHomog ) then ! would require too small subStep ! cutback makes no sense !$OMP FLUSH(terminallyIll) @@ -514,16 +516,16 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) ! restore... if (materialpoint_subStep(i,e) < 1.0_pReal) then ! protect against fake cutback from \Delta t = 2 to 1. Maybe that "trick" is not necessary anymore at all? I.e. start with \Delta t = 1 crystallite_Lp(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedLp0(1:3,1:3,1:myNgrains,i,e) ! ...plastic velocity grads + crystallite_partionedLp0(1:3,1:3,1:myNgrains,i,e) crystallite_Li(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedLi0(1:3,1:3,1:myNgrains,i,e) ! ...intermediate velocity grads + crystallite_partionedLi0(1:3,1:3,1:myNgrains,i,e) endif ! maybe protecting everything from overwriting (not only L) makes even more sense crystallite_Fp(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedFp0(1:3,1:3,1:myNgrains,i,e) ! ...plastic def grads + crystallite_partionedFp0(1:3,1:3,1:myNgrains,i,e) crystallite_Fi(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedFi0(1:3,1:3,1:myNgrains,i,e) ! ...intermediate def grads + crystallite_partionedFi0(1:3,1:3,1:myNgrains,i,e) crystallite_S(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedS0(1:3,1:3,1:myNgrains,i,e) ! ...2nd PK stress + crystallite_partionedS0(1:3,1:3,1:myNgrains,i,e) do g = 1, myNgrains plasticState (phaseAt(g,i,e))%state( :,phasememberAt(g,i,e)) = & plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) @@ -532,18 +534,15 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) enddo enddo - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - homogState(material_homogenizationAt(e))%sizeState > 0) & + if(homogState(material_homogenizationAt(e))%sizeState > 0) & homogState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = & - homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) ! ...internal homogenization state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - thermalState(material_homogenizationAt(e))%sizeState > 0) & + homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) + if(thermalState(material_homogenizationAt(e))%sizeState > 0) & thermalState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = & - thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) ! ...internal thermal state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - damageState(material_homogenizationAt(e))%sizeState > 0) & + thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) + if(damageState(material_homogenizationAt(e))%sizeState > 0) & damageState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = & - damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) ! ...internal damage state + damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) endif endif converged From 6df563624dc858ed5204c7181006fcddeb0b956a Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 21:57:20 +0200 Subject: [PATCH 11/29] type needed for postResults --- src/homogenization.f90 | 18 ++++++++++++++++-- 1 file changed, 16 insertions(+), 2 deletions(-) diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 1ff036c3a..426b37cb8 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -138,9 +138,23 @@ subroutine homogenization_init call IO_write_jobFile(FILEUNIT,'outputHomogenization') do p = 1,size(config_homogenization) if (any(material_homogenizationAt == p)) then + i = homogenization_typeInstance(p) ! which instance of this homogenization type + valid = .true. ! assume valid + select case(homogenization_type(p)) ! split per homogenization type + case (HOMOGENIZATION_NONE_ID) + outputName = HOMOGENIZATION_NONE_label + case (HOMOGENIZATION_ISOSTRAIN_ID) + outputName = HOMOGENIZATION_ISOSTRAIN_label + case (HOMOGENIZATION_RGC_ID) + outputName = HOMOGENIZATION_RGC_label + case default + valid = .false. + end select write(FILEUNIT,'(/,a,/)') '['//trim(homogenization_name(p))//']' - write(FILEUNIT,'(a,i4)') '(ngrains)'//char(9),homogenization_Ngrains(p) - + if (valid) then + write(FILEUNIT,'(a)') '(type)'//char(9)//trim(outputName) + write(FILEUNIT,'(a,i4)') '(ngrains)'//char(9),homogenization_Ngrains(p) + endif i = thermal_typeInstance(p) ! which instance of this thermal type valid = .true. ! assume valid select case(thermal_type(p)) ! split per thermal type From 118c74a268ef9417903576fd7f886d2b4ab569c1 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 22:01:12 +0200 Subject: [PATCH 12/29] enable materialpoint output again accidently lost during a former commit --- python/damask/dadf5.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/python/damask/dadf5.py b/python/damask/dadf5.py index 4aa31ec56..f5fc81b16 100644 --- a/python/damask/dadf5.py +++ b/python/damask/dadf5.py @@ -109,7 +109,7 @@ class DADF5(): except Exception as e: print('unable to locate constituents dataset: '+ str(e)) - for m in []: #self.active['materialpoints']: + for m in self.active['materialpoints']: group_materialpoint = group_inc+'/materialpoint/'+m for t in self.active['m_output_types']: try: From dce4775c17abed2b2ca09dd8b1acd53a8e477ccb Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 17 May 2019 06:36:30 +0200 Subject: [PATCH 13/29] removal of RGC out led to undefined variable --- src/homogenization.f90 | 30 +++++++++--------------------- 1 file changed, 9 insertions(+), 21 deletions(-) diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 426b37cb8..c83322f26 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -134,27 +134,14 @@ subroutine homogenization_init !-------------------------------------------------------------------------------------------------- ! write description file for homogenization output - mainProcess2: if (worldrank == 0) then + mainProcess: if (worldrank == 0) then call IO_write_jobFile(FILEUNIT,'outputHomogenization') do p = 1,size(config_homogenization) if (any(material_homogenizationAt == p)) then - i = homogenization_typeInstance(p) ! which instance of this homogenization type - valid = .true. ! assume valid - select case(homogenization_type(p)) ! split per homogenization type - case (HOMOGENIZATION_NONE_ID) - outputName = HOMOGENIZATION_NONE_label - case (HOMOGENIZATION_ISOSTRAIN_ID) - outputName = HOMOGENIZATION_ISOSTRAIN_label - case (HOMOGENIZATION_RGC_ID) - outputName = HOMOGENIZATION_RGC_label - case default - valid = .false. - end select write(FILEUNIT,'(/,a,/)') '['//trim(homogenization_name(p))//']' - if (valid) then - write(FILEUNIT,'(a)') '(type)'//char(9)//trim(outputName) - write(FILEUNIT,'(a,i4)') '(ngrains)'//char(9),homogenization_Ngrains(p) - endif + write(FILEUNIT,'(a)') '(type) n/a' + write(FILEUNIT,'(a,i4)') '(ngrains)'//char(9),homogenization_Ngrains(p) + i = thermal_typeInstance(p) ! which instance of this thermal type valid = .true. ! assume valid select case(thermal_type(p)) ! split per thermal type @@ -184,6 +171,7 @@ subroutine homogenization_init enddo endif endif + i = damage_typeInstance(p) ! which instance of this damage type valid = .true. ! assume valid select case(damage_type(p)) ! split per damage type @@ -216,7 +204,7 @@ subroutine homogenization_init endif enddo close(FILEUNIT) - endif mainProcess2 + endif mainProcess call config_deallocate('material.config/homogenization') @@ -842,12 +830,12 @@ function postResults(ip,el) postResults integer :: & startPos, endPos ,& - of, instance, homog + homog postResults = 0.0_pReal - startPos = endPos + 1 - endPos = endPos + thermalState(material_homogenizationAt(el))%sizePostResults + startPos = 1 + endPos = thermalState(material_homogenizationAt(el))%sizePostResults chosenThermal: select case (thermal_type(mesh_element(3,el))) case (THERMAL_adiabatic_ID) chosenThermal From 0fe43c50abcee7a76982a07e7815b22192c23a51 Mon Sep 17 00:00:00 2001 From: Test User Date: Fri, 17 May 2019 09:06:28 +0200 Subject: [PATCH 14/29] [skip ci] updated version information after successful test of v2.0.3-304-g7b14263c --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index 69faed072..2d8c83361 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v2.0.3-301-g789420c9 +v2.0.3-304-g7b14263c From ed8af98d695a670f3a8df0e66975d804d0851723 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 18 May 2019 06:54:45 +0200 Subject: [PATCH 15/29] don't clutter with use statements --- src/homogenization.f90 | 111 +++----------------------- src/homogenization_mech_RGC.f90 | 75 +++-------------- src/homogenization_mech_isostrain.f90 | 13 +-- src/homogenization_mech_none.f90 | 8 +- 4 files changed, 29 insertions(+), 178 deletions(-) diff --git a/src/homogenization.f90 b/src/homogenization.f90 index c83322f26..d7b5b3fdc 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -6,7 +6,20 @@ !-------------------------------------------------------------------------------------------------- module homogenization use prec + use IO + use config + use debug + use math use material + use numerics + use constitutive + use crystallite + use mesh + use FEsolving +#if defined(PETSc) || defined(DAMASK_HDF5) + use results + use HDF5_utilities +#endif !-------------------------------------------------------------------------------------------------- ! General variables for the homogenization at a material point @@ -81,26 +94,6 @@ contains !> @brief module initialization !-------------------------------------------------------------------------------------------------- subroutine homogenization_init - use math, only: & - math_I3 - use debug, only: & - debug_level, & - debug_homogenization, & - debug_levelBasic, & - debug_e, & - debug_g - use mesh, only: & - theMesh, & - mesh_element - use constitutive, only: & - constitutive_plasticity_maxSizePostResults, & - constitutive_source_maxSizePostResults - use crystallite, only: & - crystallite_maxSizePostResults - use config, only: & - config_deallocate, & - config_homogenization, & - homogenization_name use homogenization_mech_RGC use thermal_isothermal use thermal_adiabatic @@ -108,9 +101,6 @@ subroutine homogenization_init use damage_none use damage_local use damage_nonlocal - use IO - use numerics, only: & - worldrank integer, parameter :: FILEUNIT = 200 integer :: e,i,p @@ -278,51 +268,6 @@ end subroutine homogenization_init !> @brief parallelized calculation of stress and corresponding tangent at material points !-------------------------------------------------------------------------------------------------- subroutine materialpoint_stressAndItsTangent(updateJaco,dt) - use numerics, only: & - subStepMinHomog, & - subStepSizeHomog, & - stepIncreaseHomog, & - nMPstate - use FEsolving, only: & - FEsolving_execElem, & - FEsolving_execIP, & - terminallyIll - use mesh, only: & - mesh_element - use crystallite, only: & - crystallite_F0, & - crystallite_Fp0, & - crystallite_Fp, & - crystallite_Fi0, & - crystallite_Fi, & - crystallite_Lp0, & - crystallite_Lp, & - crystallite_Li0, & - crystallite_Li, & - crystallite_S0, & - crystallite_S, & - crystallite_partionedF0, & - crystallite_partionedF, & - crystallite_partionedFp0, & - crystallite_partionedLp0, & - crystallite_partionedFi0, & - crystallite_partionedLi0, & - crystallite_partionedS0, & - crystallite_dt, & - crystallite_requested, & - crystallite_stress, & - crystallite_stressTangent, & - crystallite_orientations -#ifdef DEBUG - use debug, only: & - debug_level, & - debug_homogenization, & - debug_levelBasic, & - debug_levelExtensive, & - debug_levelSelective, & - debug_e, & - debug_i -#endif real(pReal), intent(in) :: dt !< time increment logical, intent(in) :: updateJaco !< initiating Jacobian update @@ -616,14 +561,6 @@ end subroutine materialpoint_stressAndItsTangent !> @brief parallelized calculation of result array at material points !-------------------------------------------------------------------------------------------------- subroutine materialpoint_postResults - use FEsolving, only: & - FEsolving_execElem, & - FEsolving_execIP - use mesh, only: & - mesh_element - use crystallite, only: & - crystallite_sizePostResults, & - crystallite_postResults integer :: & thePos, & @@ -673,10 +610,6 @@ end subroutine materialpoint_postResults !> @brief partition material point def grad onto constituents !-------------------------------------------------------------------------------------------------- subroutine partitionDeformation(ip,el) - use mesh, only: & - mesh_element - use crystallite, only: & - crystallite_partionedF use homogenization_mech_RGC, only: & homogenization_RGC_partitionDeformation @@ -710,13 +643,6 @@ end subroutine partitionDeformation !> "happy" with result !-------------------------------------------------------------------------------------------------- function updateState(ip,el) - use mesh, only: & - mesh_element - use crystallite, only: & - crystallite_P, & - crystallite_dPdF, & - crystallite_partionedF,& - crystallite_partionedF0 use homogenization_mech_RGC, only: & homogenization_RGC_updateState use thermal_adiabatic, only: & @@ -769,10 +695,6 @@ end function updateState !> @brief derive average stress and stiffness from constituent quantities !-------------------------------------------------------------------------------------------------- subroutine averageStressAndItsTangent(ip,el) - use mesh, only: & - mesh_element - use crystallite, only: & - crystallite_P,crystallite_dPdF use homogenization_mech_RGC, only: & homogenization_RGC_averageStressAndItsTangent @@ -810,8 +732,6 @@ end subroutine averageStressAndItsTangent !> if homogenization_sizePostResults(i,e) > 0 !! !-------------------------------------------------------------------------------------------------- function postResults(ip,el) - use mesh, only: & - mesh_element use thermal_adiabatic, only: & thermal_adiabatic_postResults use thermal_conduction, only: & @@ -868,14 +788,9 @@ end function postResults !-------------------------------------------------------------------------------------------------- subroutine homogenization_results #if defined(PETSc) || defined(DAMASK_HDF5) - use results - use homogenization_mech_RGC - use HDF5_utilities use config, only: & config_name_homogenization => homogenization_name ! anticipate logical name - use material, only: & - homogenization_typeInstance, & material_homogenization_type => homogenization_type integer :: p diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index d7b1b31bf..f27733eb5 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -8,10 +8,18 @@ !-------------------------------------------------------------------------------------------------- module homogenization_mech_RGC use prec + use IO + use config + use debug + use math use material + use numerics + use constitutive +#if defined(PETSc) || defined(DAMASK_HDF5) + use results +#endif implicit none - private enum, bind(c) enumerator :: & @@ -66,34 +74,12 @@ module homogenization_mech_RGC type(tRGCdependentState), dimension(:), allocatable :: & dependentState - public :: & - homogenization_RGC_init, & - homogenization_RGC_partitionDeformation, & - homogenization_RGC_averageStressAndItsTangent, & - homogenization_RGC_updateState, & - mech_RGC_results ! name suited for planned submodule situation - contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all necessary fields, reads information from material configuration file !-------------------------------------------------------------------------------------------------- -subroutine homogenization_RGC_init() - use debug, only: & -#ifdef DEBUG - debug_i, & - debug_e, & -#endif - debug_level, & - debug_homogenization, & - debug_levelBasic - use math, only: & - math_EulerToR, & - INRAD - use IO, only: & - IO_error - use config, only: & - config_homogenization +subroutine homogenization_RGC_init integer :: & Ninstance, & @@ -218,12 +204,6 @@ end subroutine homogenization_RGC_init !> @brief partitions the deformation gradient onto the constituents !-------------------------------------------------------------------------------------------------- subroutine homogenization_RGC_partitionDeformation(F,avgF,instance,of) -#ifdef DEBUG - use debug, only: & - debug_level, & - debug_homogenization, & - debug_levelExtensive -#endif real(pReal), dimension (:,:,:), intent(out) :: F !< partioned F per grain @@ -275,24 +255,6 @@ end subroutine homogenization_RGC_partitionDeformation ! "happy" with result !-------------------------------------------------------------------------------------------------- function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) -#ifdef DEBUG - use debug, only: & - debug_level, & - debug_homogenization,& - debug_levelExtensive -#endif - use math, only: & - math_invert2 - use numerics, only: & - absTol_RGC, & - relTol_RGC, & - absMax_RGC, & - relMax_RGC, & - pPert_RGC, & - maxdRelax_RGC, & - viscPower_RGC, & - viscModus_RGC, & - refRelaxRate_RGC real(pReal), dimension(:,:,:), intent(in) :: & P,& !< array of P @@ -712,10 +674,6 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) !> @brief calculate stress-like penalty due to deformation mismatch !-------------------------------------------------------------------------------------------------- subroutine stressPenalty(rPen,nMis,avgF,fDef,ip,el,instance,of) - use math, only: & - math_civita - use numerics, only: & - xSmoo_RGC real(pReal), dimension (:,:,:), intent(out) :: rPen !< stress-like penalty real(pReal), dimension (:,:), intent(out) :: nMis !< total amount of mismatch @@ -828,13 +786,6 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) !> @brief calculate stress-like penalty due to volume discrepancy !-------------------------------------------------------------------------------------------------- subroutine volumePenalty(vPen,vDiscrep,fAvg,fDef,nGrain,instance,of) - use math, only: & - math_det33, & - math_inv33 - use numerics, only: & - maxVolDiscr_RGC,& - volDiscrMod_RGC,& - volDiscrPow_RGC real(pReal), dimension (:,:,:), intent(out) :: vPen ! stress-like penalty due to volume real(pReal), intent(out) :: vDiscrep ! total volume discrepancy @@ -883,8 +834,6 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) ! deformation !-------------------------------------------------------------------------------------------------- function surfaceCorrection(avgF,instance,of) - use math, only: & - math_invert33 real(pReal), dimension(3) :: surfaceCorrection @@ -916,8 +865,6 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) !> @brief compute the equivalent shear and bulk moduli from the elasticity tensor !-------------------------------------------------------------------------------------------------- function equivalentModuli(grainID,ip,el) - use constitutive, only: & - constitutive_homogenizedC real(pReal), dimension(2) :: equivalentModuli @@ -1015,8 +962,6 @@ end subroutine homogenization_RGC_averageStressAndItsTangent !-------------------------------------------------------------------------------------------------- subroutine mech_RGC_results(instance,group) #if defined(PETSc) || defined(DAMASK_HDF5) - use results, only: & - results_writeDataset integer, intent(in) :: instance character(len=*) :: group diff --git a/src/homogenization_mech_isostrain.f90 b/src/homogenization_mech_isostrain.f90 index 7dd7bad7d..b17085b57 100644 --- a/src/homogenization_mech_isostrain.f90 +++ b/src/homogenization_mech_isostrain.f90 @@ -5,7 +5,10 @@ !> @brief Isostrain (full constraint Taylor assuption) homogenization scheme !-------------------------------------------------------------------------------------------------- submodule(homogenization) homogenization_mech_isostrain - + use config + use debug + use IO + implicit none enum, bind(c) @@ -30,14 +33,6 @@ contains !> @brief allocates all neccessary fields, reads information from material configuration file !-------------------------------------------------------------------------------------------------- module subroutine mech_isostrain_init - use debug, only: & - debug_HOMOGENIZATION, & - debug_level, & - debug_levelBasic - use IO, only: & - IO_error - use config, only: & - config_homogenization integer :: & Ninstance, & diff --git a/src/homogenization_mech_none.f90 b/src/homogenization_mech_none.f90 index e7a5a12e6..b8b74c267 100644 --- a/src/homogenization_mech_none.f90 +++ b/src/homogenization_mech_none.f90 @@ -5,6 +5,8 @@ !> @brief dummy homogenization homogenization scheme for 1 constituent per material point !-------------------------------------------------------------------------------------------------- submodule(homogenization) homogenization_mech_none + use config + use debug implicit none @@ -14,12 +16,6 @@ contains !> @brief allocates all neccessary fields, reads information from material configuration file !-------------------------------------------------------------------------------------------------- module subroutine mech_none_init - use debug, only: & - debug_HOMOGENIZATION, & - debug_level, & - debug_levelBasic - use config, only: & - config_homogenization integer :: & Ninstance, & From 2258bfb22161c60987043ef89bc2d4e4dfeac132 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 18 May 2019 07:23:46 +0200 Subject: [PATCH 16/29] RGC as submodule submodules inherit use-associated entities and implicit none/private statements --- src/homogenization.f90 | 71 ++++++++++++++++++++------- src/homogenization_mech_RGC.f90 | 70 +++++++++++--------------- src/homogenization_mech_isostrain.f90 | 5 -- src/homogenization_mech_none.f90 | 6 +-- 4 files changed, 83 insertions(+), 69 deletions(-) diff --git a/src/homogenization.f90 b/src/homogenization.f90 index d7b5b3fdc..cbc6471b4 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -60,11 +60,24 @@ module homogenization module subroutine mech_isostrain_init end subroutine mech_isostrain_init + module subroutine mech_RGC_init + end subroutine mech_RGC_init + + module subroutine mech_isostrain_partitionDeformation(F,avgF) real(pReal), dimension (:,:,:), intent(out) :: F !< partitioned deformation gradient real(pReal), dimension (3,3), intent(in) :: avgF !< average deformation gradient at material point end subroutine mech_isostrain_partitionDeformation + module subroutine mech_RGC_partitionDeformation(F,avgF,instance,of) + real(pReal), dimension (:,:,:), intent(out) :: F !< partitioned deformation gradient + real(pReal), dimension (3,3), intent(in) :: avgF !< average deformation gradient at material point + integer, intent(in) :: & + instance, & + of + end subroutine mech_RGC_partitionDeformation + + module subroutine mech_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) real(pReal), dimension (3,3), intent(out) :: avgP !< average stress at material point real(pReal), dimension (3,3,3,3), intent(out) :: dAvgPdAvgF !< average stiffness at material point @@ -73,7 +86,37 @@ module homogenization real(pReal), dimension (:,:,:,:,:), intent(in) :: dPdF !< partitioned stiffnesses integer, intent(in) :: instance end subroutine mech_isostrain_averageStressAndItsTangent + + module subroutine mech_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) + real(pReal), dimension (3,3), intent(out) :: avgP !< average stress at material point + real(pReal), dimension (3,3,3,3), intent(out) :: dAvgPdAvgF !< average stiffness at material point + + real(pReal), dimension (:,:,:), intent(in) :: P !< partitioned stresses + real(pReal), dimension (:,:,:,:,:), intent(in) :: dPdF !< partitioned stiffnesses + integer, intent(in) :: instance + end subroutine mech_RGC_averageStressAndItsTangent + + + module function mech_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) +logical, dimension(2) :: mech_RGC_updateState + real(pReal), dimension(:,:,:), intent(in) :: & + P,& !< array of P + F,& !< array of F + F0 !< array of initial F + real(pReal), dimension(:,:,:,:,:), intent(in) :: dPdF !< array of current grain stiffness + real(pReal), dimension(3,3), intent(in) :: avgF !< average F + real(pReal), intent(in) :: dt !< time increment + integer, intent(in) :: & + ip, & !< integration point number + el !< element number + end function + + module subroutine mech_RGC_results(instance,group) + + integer, intent(in) :: instance + character(len=*), intent(in) :: group + end subroutine end interface public :: & @@ -112,7 +155,7 @@ subroutine homogenization_init if (any(homogenization_type == HOMOGENIZATION_NONE_ID)) call mech_none_init if (any(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID)) call mech_isostrain_init - if (any(homogenization_type == HOMOGENIZATION_RGC_ID)) call homogenization_RGC_init + if (any(homogenization_type == HOMOGENIZATION_RGC_ID)) call mech_RGC_init if (any(thermal_type == THERMAL_isothermal_ID)) call thermal_isothermal_init if (any(thermal_type == THERMAL_adiabatic_ID)) call thermal_adiabatic_init @@ -610,8 +653,6 @@ end subroutine materialpoint_postResults !> @brief partition material point def grad onto constituents !-------------------------------------------------------------------------------------------------- subroutine partitionDeformation(ip,el) - use homogenization_mech_RGC, only: & - homogenization_RGC_partitionDeformation integer, intent(in) :: & ip, & !< integration point @@ -628,7 +669,7 @@ subroutine partitionDeformation(ip,el) materialpoint_subF(1:3,1:3,ip,el)) case (HOMOGENIZATION_RGC_ID) chosenHomogenization - call homogenization_RGC_partitionDeformation(& + call mech_RGC_partitionDeformation(& crystallite_partionedF(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & materialpoint_subF(1:3,1:3,ip,el),& ip, & @@ -643,8 +684,6 @@ end subroutine partitionDeformation !> "happy" with result !-------------------------------------------------------------------------------------------------- function updateState(ip,el) - use homogenization_mech_RGC, only: & - homogenization_RGC_updateState use thermal_adiabatic, only: & thermal_adiabatic_updateState use damage_local, only: & @@ -660,14 +699,14 @@ function updateState(ip,el) case (HOMOGENIZATION_RGC_ID) chosenHomogenization updateState = & updateState .and. & - homogenization_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & - crystallite_partionedF(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & - crystallite_partionedF0(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el),& - materialpoint_subF(1:3,1:3,ip,el),& - materialpoint_subdt(ip,el), & - crystallite_dPdF(1:3,1:3,1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & - ip, & - el) + mech_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & + crystallite_partionedF(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & + crystallite_partionedF0(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el),& + materialpoint_subF(1:3,1:3,ip,el),& + materialpoint_subdt(ip,el), & + crystallite_dPdF(1:3,1:3,1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & + ip, & + el) end select chosenHomogenization chosenThermal: select case (thermal_type(mesh_element(3,el))) @@ -695,8 +734,6 @@ end function updateState !> @brief derive average stress and stiffness from constituent quantities !-------------------------------------------------------------------------------------------------- subroutine averageStressAndItsTangent(ip,el) - use homogenization_mech_RGC, only: & - homogenization_RGC_averageStressAndItsTangent integer, intent(in) :: & ip, & !< integration point @@ -716,7 +753,7 @@ subroutine averageStressAndItsTangent(ip,el) homogenization_typeInstance(mesh_element(3,el))) case (HOMOGENIZATION_RGC_ID) chosenHomogenization - call homogenization_RGC_averageStressAndItsTangent(& + call mech_RGC_averageStressAndItsTangent(& materialpoint_P(1:3,1:3,ip,el), & materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el),& crystallite_P(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index f27733eb5..4e115a498 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -6,20 +6,7 @@ !> @brief Relaxed grain cluster (RGC) homogenization scheme !> Nconstituents is defined as p x q x r (cluster) !-------------------------------------------------------------------------------------------------- -module homogenization_mech_RGC - use prec - use IO - use config - use debug - use math - use material - use numerics - use constitutive -#if defined(PETSc) || defined(DAMASK_HDF5) - use results -#endif - - implicit none +submodule(homogenization) homogenization_mech_RGC enum, bind(c) enumerator :: & @@ -79,7 +66,7 @@ contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all necessary fields, reads information from material configuration file !-------------------------------------------------------------------------------------------------- -subroutine homogenization_RGC_init +module subroutine mech_RGC_init integer :: & Ninstance, & @@ -197,17 +184,17 @@ subroutine homogenization_RGC_init enddo -end subroutine homogenization_RGC_init +end subroutine mech_RGC_init !-------------------------------------------------------------------------------------------------- !> @brief partitions the deformation gradient onto the constituents !-------------------------------------------------------------------------------------------------- -subroutine homogenization_RGC_partitionDeformation(F,avgF,instance,of) +module subroutine mech_RGC_partitionDeformation(F,avgF,instance,of) real(pReal), dimension (:,:,:), intent(out) :: F !< partioned F per grain - real(pReal), dimension (:,:), intent(in) :: avgF !< averaged F + real(pReal), dimension (3,3), intent(in) :: avgF !< averaged F integer, intent(in) :: & instance, & of @@ -247,14 +234,14 @@ subroutine homogenization_RGC_partitionDeformation(F,avgF,instance,of) end associate -end subroutine homogenization_RGC_partitionDeformation +end subroutine mech_RGC_partitionDeformation !-------------------------------------------------------------------------------------------------- !> @brief update the internal state of the homogenization scheme and tell whether "done" and ! "happy" with result !-------------------------------------------------------------------------------------------------- -function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) +module function mech_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) real(pReal), dimension(:,:,:), intent(in) :: & P,& !< array of P @@ -267,8 +254,6 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) ip, & !< integration point number el !< element number - logical, dimension(2) :: homogenization_RGC_updateState - integer, dimension(4) :: intFaceN,intFaceP,faceID integer, dimension(3) :: nGDim,iGr3N,iGr3P integer :: instance,iNum,i,j,nIntFaceTot,iGrN,iGrP,iMun,iFace,k,l,ipert,iGrain,nGrain, of @@ -285,7 +270,7 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) #endif zeroTimeStep: if(dEq0(dt)) then - homogenization_RGC_updateState = .true. ! pretend everything is fine and return + mech_RGC_updateState = .true. ! pretend everything is fine and return return endif zeroTimeStep @@ -404,12 +389,12 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) endif #endif - homogenization_RGC_updateState = .false. + mech_RGC_updateState = .false. !-------------------------------------------------------------------------------------------------- ! If convergence reached => done and happy if (residMax < relTol_RGC*stresMax .or. residMax < absTol_RGC) then - homogenization_RGC_updateState = .true. + mech_RGC_updateState = .true. #ifdef DEBUG if (iand(debug_level(debug_homogenization),debug_levelExtensive) /= 0 & .and. prm%of_debug == of) write(6,'(1x,a55,/)')'... done and happy' @@ -451,7 +436,7 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) !-------------------------------------------------------------------------------------------------- ! if residual blows-up => done but unhappy elseif (residMax > relMax_RGC*stresMax .or. residMax > absMax_RGC) then ! try to restart when residual blows up exceeding maximum bound - homogenization_RGC_updateState = [.true.,.false.] ! with direct cut-back + mech_RGC_updateState = [.true.,.false.] ! with direct cut-back #ifdef DEBUG if (iand(debug_level(debug_homogenization),debug_levelExtensive) /= 0 & @@ -648,7 +633,7 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) enddo; enddo stt%relaxationVector(:,of) = relax + drelax ! Updateing the state variable for the next iteration if (any(abs(drelax) > maxdRelax_RGC)) then ! Forcing cutback when the incremental change of relaxation vector becomes too large - homogenization_RGC_updateState = [.true.,.false.] + mech_RGC_updateState = [.true.,.false.] !$OMP CRITICAL (write2out) write(6,'(1x,a,1x,i3,1x,a,1x,i3,1x,a)')'RGC_updateState: ip',ip,'| el',el,'enforces cutback' write(6,'(1x,a,1x,e15.8)')'due to large relaxation change =',maxval(abs(drelax)) @@ -935,13 +920,13 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) end subroutine grainDeformation -end function homogenization_RGC_updateState +end function mech_RGC_updateState !-------------------------------------------------------------------------------------------------- !> @brief derive average stress and stiffness from constituent quantities !-------------------------------------------------------------------------------------------------- -subroutine homogenization_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) +module subroutine mech_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) real(pReal), dimension (3,3), intent(out) :: avgP !< average stress at material point real(pReal), dimension (3,3,3,3), intent(out) :: dAvgPdAvgF !< average stiffness at material point @@ -953,7 +938,7 @@ subroutine homogenization_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF, avgP = sum(P,3) /real(product(param(instance)%Nconstituents),pReal) dAvgPdAvgF = sum(dPdF,5)/real(product(param(instance)%Nconstituents),pReal) -end subroutine homogenization_RGC_averageStressAndItsTangent +end subroutine mech_RGC_averageStressAndItsTangent !-------------------------------------------------------------------------------------------------- @@ -963,8 +948,9 @@ end subroutine homogenization_RGC_averageStressAndItsTangent subroutine mech_RGC_results(instance,group) #if defined(PETSc) || defined(DAMASK_HDF5) - integer, intent(in) :: instance - character(len=*) :: group + integer, intent(in) :: instance + character(len=*), intent(in) :: group + integer :: o associate(stt => state(instance), dst => dependentState(instance), prm => param(instance)) @@ -1135,7 +1121,7 @@ integer pure function interface4to1(iFace4D, nGDim) else interface4to1 = iFace4D(4) + nGDim(3)*(iFace4D(2)-1) & + nGDim(3)*nGDim(1)*(iFace4D(3)-1) & - + (nGDim(1)-1)*nGDim(2)*nGDim(3) ! total number of interfaces normal //e1 + + (nGDim(1)-1)*nGDim(2)*nGDim(3) ! total # of interfaces normal || e1 endif case(3) @@ -1144,8 +1130,8 @@ integer pure function interface4to1(iFace4D, nGDim) else interface4to1 = iFace4D(2) + nGDim(1)*(iFace4D(3)-1) & + nGDim(1)*nGDim(2)*(iFace4D(4)-1) & - + (nGDim(1)-1)*nGDim(2)*nGDim(3) & ! total number of interfaces normal //e1 - + nGDim(1)*(nGDim(2)-1)*nGDim(3) ! total number of interfaces normal //e2 + + (nGDim(1)-1)*nGDim(2)*nGDim(3) & ! total # of interfaces normal || e1 + + nGDim(1)*(nGDim(2)-1)*nGDim(3) ! total # of interfaces normal || e2 endif case default @@ -1169,23 +1155,23 @@ pure function interface1to4(iFace1D, nGDim) !-------------------------------------------------------------------------------------------------- ! compute the total number of interfaces, which ... - nIntFace = [(nGDim(1)-1)*nGDim(2)*nGDim(3), & ! ... normal //e1 - nGDim(1)*(nGDim(2)-1)*nGDim(3), & ! ... normal //e2 - nGDim(1)*nGDim(2)*(nGDim(3)-1)] ! ... normal //e3 + nIntFace = [(nGDim(1)-1)*nGDim(2)*nGDim(3), & ! ... normal || e1 + nGDim(1)*(nGDim(2)-1)*nGDim(3), & ! ... normal || e2 + nGDim(1)*nGDim(2)*(nGDim(3)-1)] ! ... normal || e3 !-------------------------------------------------------------------------------------------------- ! get the corresponding interface ID in 4D (normal and local position) - if (iFace1D > 0 .and. iFace1D <= nIntFace(1)) then ! interface with normal //e1 + if (iFace1D > 0 .and. iFace1D <= nIntFace(1)) then ! interface with normal || e1 interface1to4(1) = 1 interface1to4(3) = mod((iFace1D-1),nGDim(2))+1 interface1to4(4) = mod(int(real(iFace1D-1,pReal)/real(nGDim(2),pReal)),nGDim(3))+1 interface1to4(2) = int(real(iFace1D-1,pReal)/real(nGDim(2),pReal)/real(nGDim(3),pReal))+1 - elseif (iFace1D > nIntFace(1) .and. iFace1D <= (nIntFace(2) + nIntFace(1))) then ! interface with normal //e2 + elseif (iFace1D > nIntFace(1) .and. iFace1D <= (nIntFace(2) + nIntFace(1))) then ! interface with normal || e2 interface1to4(1) = 2 interface1to4(4) = mod((iFace1D-nIntFace(1)-1),nGDim(3))+1 interface1to4(2) = mod(int(real(iFace1D-nIntFace(1)-1,pReal)/real(nGDim(3),pReal)),nGDim(1))+1 interface1to4(3) = int(real(iFace1D-nIntFace(1)-1,pReal)/real(nGDim(3),pReal)/real(nGDim(1),pReal))+1 - elseif (iFace1D > nIntFace(2) + nIntFace(1) .and. iFace1D <= (nIntFace(3) + nIntFace(2) + nIntFace(1))) then ! interface with normal //e3 + elseif (iFace1D > nIntFace(2) + nIntFace(1) .and. iFace1D <= (nIntFace(3) + nIntFace(2) + nIntFace(1))) then ! interface with normal || e3 interface1to4(1) = 3 interface1to4(2) = mod((iFace1D-nIntFace(2)-nIntFace(1)-1),nGDim(1))+1 interface1to4(3) = mod(int(real(iFace1D-nIntFace(2)-nIntFace(1)-1,pReal)/real(nGDim(1),pReal)),nGDim(2))+1 @@ -1195,4 +1181,4 @@ pure function interface1to4(iFace1D, nGDim) end function interface1to4 -end module homogenization_mech_RGC +end submodule homogenization_mech_RGC diff --git a/src/homogenization_mech_isostrain.f90 b/src/homogenization_mech_isostrain.f90 index b17085b57..c1eba821c 100644 --- a/src/homogenization_mech_isostrain.f90 +++ b/src/homogenization_mech_isostrain.f90 @@ -5,11 +5,6 @@ !> @brief Isostrain (full constraint Taylor assuption) homogenization scheme !-------------------------------------------------------------------------------------------------- submodule(homogenization) homogenization_mech_isostrain - use config - use debug - use IO - - implicit none enum, bind(c) enumerator :: & diff --git a/src/homogenization_mech_none.f90 b/src/homogenization_mech_none.f90 index b8b74c267..d5b24e5d1 100644 --- a/src/homogenization_mech_none.f90 +++ b/src/homogenization_mech_none.f90 @@ -5,10 +5,6 @@ !> @brief dummy homogenization homogenization scheme for 1 constituent per material point !-------------------------------------------------------------------------------------------------- submodule(homogenization) homogenization_mech_none - use config - use debug - - implicit none contains @@ -23,7 +19,7 @@ module subroutine mech_none_init NofMyHomog write(6,'(/,a)') ' <<<+- homogenization_'//HOMOGENIZATION_NONE_label//' init -+>>>' - + Ninstance = count(homogenization_type == HOMOGENIZATION_NONE_ID) if (iand(debug_level(debug_HOMOGENIZATION),debug_levelBasic) /= 0) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance From 34bcd382403d1dff1c7e136d53db3e9672cb269f Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 18 May 2019 07:39:55 +0200 Subject: [PATCH 17/29] cleanup --- src/commercialFEM_fileList.f90 | 2 +- src/homogenization.f90 | 53 ++++++++++++--------------------- src/homogenization_mech_RGC.f90 | 8 ++--- 3 files changed, 24 insertions(+), 39 deletions(-) diff --git a/src/commercialFEM_fileList.f90 b/src/commercialFEM_fileList.f90 index 0ae1323f0..5e1cd71eb 100644 --- a/src/commercialFEM_fileList.f90 +++ b/src/commercialFEM_fileList.f90 @@ -46,7 +46,6 @@ #include "plastic_nonlocal.f90" #include "constitutive.f90" #include "crystallite.f90" -#include "homogenization_mech_RGC.f90" #include "thermal_isothermal.f90" #include "thermal_adiabatic.f90" #include "thermal_conduction.f90" @@ -56,4 +55,5 @@ #include "homogenization.f90" #include "homogenization_mech_none.f90" #include "homogenization_mech_isostrain.f90" +#include "homogenization_mech_RGC.f90" #include "CPFEM.f90" diff --git a/src/homogenization.f90 b/src/homogenization.f90 index cbc6471b4..be07318d7 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -35,7 +35,6 @@ module homogenization materialpoint_results !< results array of material point integer, public, protected :: & materialpoint_sizeResults, & - homogenization_maxSizePostResults, & thermal_maxSizePostResults, & damage_maxSizePostResults @@ -98,25 +97,25 @@ module homogenization module function mech_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) -logical, dimension(2) :: mech_RGC_updateState - real(pReal), dimension(:,:,:), intent(in) :: & - P,& !< array of P - F,& !< array of F - F0 !< array of initial F - real(pReal), dimension(:,:,:,:,:), intent(in) :: dPdF !< array of current grain stiffness - real(pReal), dimension(3,3), intent(in) :: avgF !< average F - real(pReal), intent(in) :: dt !< time increment - integer, intent(in) :: & - ip, & !< integration point number - el !< element number - end function - - module subroutine mech_RGC_results(instance,group) + logical, dimension(2) :: mech_RGC_updateState + real(pReal), dimension(:,:,:), intent(in) :: & + P,& !< partitioned stresses + F,& !< partitioned deformation gradients + F0 !< partitioned initial deformation gradients + real(pReal), dimension(:,:,:,:,:), intent(in) :: dPdF !< partitioned stiffnesses + real(pReal), dimension(3,3), intent(in) :: avgF !< average F + real(pReal), intent(in) :: dt !< time increment + integer, intent(in) :: & + ip, & !< integration point number + el !< element number + end function mech_RGC_updateState - integer, intent(in) :: instance - character(len=*), intent(in) :: group - end subroutine + module subroutine mech_RGC_results(instance,group) + integer, intent(in) :: instance !< homogenization instance + character(len=*), intent(in) :: group !< group name in HDF5 file + end subroutine mech_RGC_results + end interface public :: & @@ -124,11 +123,6 @@ logical, dimension(2) :: mech_RGC_updateState materialpoint_stressAndItsTangent, & materialpoint_postResults, & homogenization_results - private :: & - partitionDeformation, & - updateState, & - averageStressAndItsTangent, & - postResults contains @@ -137,7 +131,6 @@ contains !> @brief module initialization !-------------------------------------------------------------------------------------------------- subroutine homogenization_init - use homogenization_mech_RGC use thermal_isothermal use thermal_adiabatic use thermal_conduction @@ -245,7 +238,7 @@ subroutine homogenization_init ! allocate and initialize global variables allocate(materialpoint_dPdF(3,3,3,3,theMesh%elem%nIPs,theMesh%nElems), source=0.0_pReal) allocate(materialpoint_F0(3,3,theMesh%elem%nIPs,theMesh%nElems), source=0.0_pReal) - materialpoint_F0 = spread(spread(math_I3,3,theMesh%elem%nIPs),4,theMesh%nElems) ! initialize to identity + materialpoint_F0 = spread(spread(math_I3,3,theMesh%elem%nIPs),4,theMesh%nElems) ! initialize to identity allocate(materialpoint_F(3,3,theMesh%elem%nIPs,theMesh%nElems), source=0.0_pReal) materialpoint_F = materialpoint_F0 ! initialize to identity allocate(materialpoint_subF0(3,3,theMesh%elem%nIPs,theMesh%nElems), source=0.0_pReal) @@ -260,18 +253,15 @@ subroutine homogenization_init !-------------------------------------------------------------------------------------------------- ! allocate and initialize global state and postresutls variables - homogenization_maxSizePostResults = 0 thermal_maxSizePostResults = 0 damage_maxSizePostResults = 0 do p = 1,size(config_homogenization) - homogenization_maxSizePostResults = max(homogenization_maxSizePostResults,homogState (p)%sizePostResults) thermal_maxSizePostResults = max(thermal_maxSizePostResults, thermalState (p)%sizePostResults) damage_maxSizePostResults = max(damage_maxSizePostResults ,damageState (p)%sizePostResults) enddo materialpoint_sizeResults = 1 & ! grain count - + 1 + homogenization_maxSizePostResults & ! homogSize & homogResult - + thermal_maxSizePostResults & + + 1 + thermal_maxSizePostResults & + damage_maxSizePostResults & + homogenization_maxNgrains * (1 + crystallite_maxSizePostResults & ! crystallite size & crystallite results + 1 + constitutive_plasticity_maxSizePostResults & ! constitutive size & constitutive results @@ -281,11 +271,6 @@ subroutine homogenization_init write(6,'(/,a)') ' <<<+- homogenization init -+>>>' if (iand(debug_level(debug_homogenization), debug_levelBasic) /= 0) then -#ifdef TODO - write(6,'(a32,1x,7(i8,1x))') 'homogenization_state0: ', shape(homogenization_state0) - write(6,'(a32,1x,7(i8,1x))') 'homogenization_subState0: ', shape(homogenization_subState0) - write(6,'(a32,1x,7(i8,1x))') 'homogenization_state: ', shape(homogenization_state) -#endif write(6,'(a32,1x,7(i8,1x))') 'materialpoint_dPdF: ', shape(materialpoint_dPdF) write(6,'(a32,1x,7(i8,1x))') 'materialpoint_F0: ', shape(materialpoint_F0) write(6,'(a32,1x,7(i8,1x))') 'materialpoint_F: ', shape(materialpoint_F) diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index 4e115a498..a1e92d575 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -30,7 +30,7 @@ submodule(homogenization) homogenization_mech_RGC angles integer :: & of_debug = 0 - integer(kind(undefined_ID)), dimension(:), allocatable :: & + integer(kind(undefined_ID)), dimension(:), allocatable :: & outputID end type tParameters @@ -945,7 +945,7 @@ end subroutine mech_RGC_averageStressAndItsTangent !> @brief writes results to HDF5 output file ! ToDo: check wheter units are correct !-------------------------------------------------------------------------------------------------- -subroutine mech_RGC_results(instance,group) +module subroutine mech_RGC_results(instance,group) #if defined(PETSc) || defined(DAMASK_HDF5) integer, intent(in) :: instance @@ -981,8 +981,8 @@ subroutine mech_RGC_results(instance,group) end associate #else - integer, intent(in) :: instance - character(len=*) :: group + integer, intent(in) :: instance + character(len=*), intent(in) :: group #endif end subroutine mech_RGC_results From 86205f508149e90c0e00f85af0b235c16aa0bdb7 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 18 May 2019 09:47:20 +0200 Subject: [PATCH 18/29] gfortran complaints about repeated dimension attribute we need to decide whether we want to repeat the declaration of the interface or not --- src/homogenization_mech_RGC.f90 | 18 ++++-------------- 1 file changed, 4 insertions(+), 14 deletions(-) diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index a1e92d575..1cbe837d2 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -241,18 +241,7 @@ end subroutine mech_RGC_partitionDeformation !> @brief update the internal state of the homogenization scheme and tell whether "done" and ! "happy" with result !-------------------------------------------------------------------------------------------------- -module function mech_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) - - real(pReal), dimension(:,:,:), intent(in) :: & - P,& !< array of P - F,& !< array of F - F0 !< array of initial F - real(pReal), dimension(:,:,:,:,:), intent(in) :: dPdF !< array of current grain stiffness - real(pReal), dimension(3,3), intent(in) :: avgF !< average F - real(pReal), intent(in) :: dt !< time increment - integer, intent(in) :: & - ip, & !< integration point number - el !< element number +module procedure mech_RGC_updateState integer, dimension(4) :: intFaceN,intFaceP,faceID integer, dimension(3) :: nGDim,iGr3N,iGr3P @@ -579,9 +568,10 @@ module function mech_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) !-------------------------------------------------------------------------------------------------- ! ... of the numerical viscosity traction "rmatrix" allocate(rmatrix(3*nIntFaceTot,3*nIntFaceTot),source=0.0_pReal) - forall (i=1:3*nIntFaceTot) & + do i=1,3*nIntFaceTot rmatrix(i,i) = viscModus_RGC*viscPower_RGC/(refRelaxRate_RGC*dt)* & ! tangent due to numerical viscosity traction appears (abs(drelax(i))/(refRelaxRate_RGC*dt))**(viscPower_RGC - 1.0_pReal) ! only in the main diagonal term + enddo #ifdef DEBUG if (iand(debug_level(debug_homogenization), debug_levelExtensive) /= 0) then @@ -920,7 +910,7 @@ module function mech_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) end subroutine grainDeformation -end function mech_RGC_updateState +end procedure mech_RGC_updateState !-------------------------------------------------------------------------------------------------- From 60ed514e0cdc013e764256773332ad7729e807ae Mon Sep 17 00:00:00 2001 From: Test User Date: Wed, 15 May 2019 23:51:29 +0200 Subject: [PATCH 19/29] [skip ci] updated version information after successful test of v2.0.3-301-g789420c9 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index acf2aa8dc..69faed072 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v2.0.3-297-gae084bb2 +v2.0.3-301-g789420c9 From 60c2a5fc068393033a778bdf2a823a07bdf9dcd8 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 16 May 2019 18:24:54 +0200 Subject: [PATCH 20/29] loop (forall) over integration points wrong this was done for each integration point, but this was not detected for the forall loop --- src/homogenization.f90 | 119 ++++++++++++++++++++--------------------- 1 file changed, 59 insertions(+), 60 deletions(-) diff --git a/src/homogenization.f90 b/src/homogenization.f90 index be07318d7..3210f02d4 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -323,43 +323,46 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) ! initialize restoration points of ... do e = FEsolving_execElem(1),FEsolving_execElem(2) myNgrains = homogenization_Ngrains(mesh_element(3,e)) - do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do g = 1,myNgrains + do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); + do g = 1,myNgrains + + plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = & + plasticState (phaseAt(g,i,e))%state0( :,phasememberAt(g,i,e)) + do mySource = 1, phase_Nsources(phaseAt(g,i,e)) + sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) = & + sourceState(phaseAt(g,i,e))%p(mySource)%state0( :,phasememberAt(g,i,e)) + enddo + + crystallite_partionedFp0(1:3,1:3,g,i,e) = crystallite_Fp0(1:3,1:3,g,i,e) + crystallite_partionedLp0(1:3,1:3,g,i,e) = crystallite_Lp0(1:3,1:3,g,i,e) + crystallite_partionedFi0(1:3,1:3,g,i,e) = crystallite_Fi0(1:3,1:3,g,i,e) + crystallite_partionedLi0(1:3,1:3,g,i,e) = crystallite_Li0(1:3,1:3,g,i,e) + crystallite_partionedF0(1:3,1:3,g,i,e) = crystallite_F0(1:3,1:3,g,i,e) + crystallite_partionedS0(1:3,1:3,g,i,e) = crystallite_S0(1:3,1:3,g,i,e) - plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = & - plasticState (phaseAt(g,i,e))%state0( :,phasememberAt(g,i,e)) - do mySource = 1, phase_Nsources(phaseAt(g,i,e)) - sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) = & - sourceState(phaseAt(g,i,e))%p(mySource)%state0( :,phasememberAt(g,i,e)) enddo - crystallite_partionedFp0(1:3,1:3,g,i,e) = crystallite_Fp0(1:3,1:3,g,i,e) ! ...plastic def grads - crystallite_partionedLp0(1:3,1:3,g,i,e) = crystallite_Lp0(1:3,1:3,g,i,e) ! ...plastic velocity grads - crystallite_partionedFi0(1:3,1:3,g,i,e) = crystallite_Fi0(1:3,1:3,g,i,e) ! ...intermediate def grads - crystallite_partionedLi0(1:3,1:3,g,i,e) = crystallite_Li0(1:3,1:3,g,i,e) ! ...intermediate velocity grads - crystallite_partionedF0(1:3,1:3,g,i,e) = crystallite_F0(1:3,1:3,g,i,e) ! ...def grads - crystallite_partionedS0(1:3,1:3,g,i,e) = crystallite_S0(1:3,1:3,g,i,e) ! ...2nd PK stress - enddo; enddo - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e)) - materialpoint_subF0(1:3,1:3,i,e) = materialpoint_F0(1:3,1:3,i,e) ! ...def grad + materialpoint_subF0(1:3,1:3,i,e) = materialpoint_F0(1:3,1:3,i,e) materialpoint_subFrac(i,e) = 0.0_pReal materialpoint_subStep(i,e) = 1.0_pReal/subStepSizeHomog ! <> materialpoint_converged(i,e) = .false. ! pretend failed step of twice the required size materialpoint_requested(i,e) = .true. ! everybody requires calculation - endforall - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - homogState(material_homogenizationAt(e))%sizeState > 0) & - homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - homogState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal homogenization state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - thermalState(material_homogenizationAt(e))%sizeState > 0) & - thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - thermalState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal thermal state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - damageState(material_homogenizationAt(e))%sizeState > 0) & - damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - damageState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal damage state + + if (homogState(material_homogenizationAt(e))%sizeState > 0) & + homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & + homogState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal homogenization state + + if (thermalState(material_homogenizationAt(e))%sizeState > 0) & + thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & + thermalState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal thermal state + + if (damageState(material_homogenizationAt(e))%sizeState > 0) & + damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & + damageState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal damage state + enddo enddo + NiterationHomog = 0 cutBackLooping: do while (.not. terminallyIll .and. & @@ -370,7 +373,7 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) myNgrains = homogenization_Ngrains(mesh_element(3,e)) IpLooping1: do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) - converged: if ( materialpoint_converged(i,e) ) then + converged: if (materialpoint_converged(i,e)) then #ifdef DEBUG if (iand(debug_level(debug_homogenization), debug_levelExtensive) /= 0 & .and. ((e == debug_e .and. i == debug_i) & @@ -391,22 +394,22 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) ! wind forward grain starting point of... crystallite_partionedF0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedF(1:3,1:3,1:myNgrains,i,e) ! ...def grads + crystallite_partionedF(1:3,1:3,1:myNgrains,i,e) crystallite_partionedFp0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_Fp (1:3,1:3,1:myNgrains,i,e) ! ...plastic def grads + crystallite_Fp (1:3,1:3,1:myNgrains,i,e) crystallite_partionedLp0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_Lp (1:3,1:3,1:myNgrains,i,e) ! ...plastic velocity grads + crystallite_Lp (1:3,1:3,1:myNgrains,i,e) crystallite_partionedFi0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_Fi (1:3,1:3,1:myNgrains,i,e) ! ...intermediate def grads + crystallite_Fi (1:3,1:3,1:myNgrains,i,e) crystallite_partionedLi0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_Li (1:3,1:3,1:myNgrains,i,e) ! ...intermediate velocity grads + crystallite_Li (1:3,1:3,1:myNgrains,i,e) crystallite_partionedS0 (1:3,1:3,1:myNgrains,i,e) = & - crystallite_S (1:3,1:3,1:myNgrains,i,e) ! ...2nd PK stress + crystallite_S (1:3,1:3,1:myNgrains,i,e) do g = 1,myNgrains plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = & @@ -417,23 +420,22 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) enddo enddo - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - homogState(material_homogenizationAt(e))%sizeState > 0) & + if(homogState(material_homogenizationAt(e))%sizeState > 0) & homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - homogState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) ! ...internal homogenization state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - thermalState(material_homogenizationAt(e))%sizeState > 0) & + homogState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) + if(thermalState(material_homogenizationAt(e))%sizeState > 0) & thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - thermalState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) ! ...internal thermal state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - damageState(material_homogenizationAt(e))%sizeState > 0) & + thermalState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) + if(damageState(material_homogenizationAt(e))%sizeState > 0) & damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - damageState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) ! ...internal damage state - materialpoint_subF0(1:3,1:3,i,e) = materialpoint_subF(1:3,1:3,i,e) ! ...def grad + damageState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) + + materialpoint_subF0(1:3,1:3,i,e) = materialpoint_subF(1:3,1:3,i,e) + endif steppingNeeded else converged - if ( (myNgrains == 1 .and. materialpoint_subStep(i,e) <= 1.0 ) .or. & ! single grain already tried internal subStepping in crystallite + if ( (myNgrains == 1 .and. materialpoint_subStep(i,e) <= 1.0 ) .or. & ! single grain already tried internal subStepping in crystallite subStepSizeHomog * materialpoint_subStep(i,e) <= subStepMinHomog ) then ! would require too small subStep ! cutback makes no sense !$OMP FLUSH(terminallyIll) @@ -462,16 +464,16 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) ! restore... if (materialpoint_subStep(i,e) < 1.0_pReal) then ! protect against fake cutback from \Delta t = 2 to 1. Maybe that "trick" is not necessary anymore at all? I.e. start with \Delta t = 1 crystallite_Lp(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedLp0(1:3,1:3,1:myNgrains,i,e) ! ...plastic velocity grads + crystallite_partionedLp0(1:3,1:3,1:myNgrains,i,e) crystallite_Li(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedLi0(1:3,1:3,1:myNgrains,i,e) ! ...intermediate velocity grads + crystallite_partionedLi0(1:3,1:3,1:myNgrains,i,e) endif ! maybe protecting everything from overwriting (not only L) makes even more sense crystallite_Fp(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedFp0(1:3,1:3,1:myNgrains,i,e) ! ...plastic def grads + crystallite_partionedFp0(1:3,1:3,1:myNgrains,i,e) crystallite_Fi(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedFi0(1:3,1:3,1:myNgrains,i,e) ! ...intermediate def grads + crystallite_partionedFi0(1:3,1:3,1:myNgrains,i,e) crystallite_S(1:3,1:3,1:myNgrains,i,e) = & - crystallite_partionedS0(1:3,1:3,1:myNgrains,i,e) ! ...2nd PK stress + crystallite_partionedS0(1:3,1:3,1:myNgrains,i,e) do g = 1, myNgrains plasticState (phaseAt(g,i,e))%state( :,phasememberAt(g,i,e)) = & plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) @@ -480,18 +482,15 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) enddo enddo - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - homogState(material_homogenizationAt(e))%sizeState > 0) & + if(homogState(material_homogenizationAt(e))%sizeState > 0) & homogState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = & - homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) ! ...internal homogenization state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - thermalState(material_homogenizationAt(e))%sizeState > 0) & + homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) + if(thermalState(material_homogenizationAt(e))%sizeState > 0) & thermalState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = & - thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) ! ...internal thermal state - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & - damageState(material_homogenizationAt(e))%sizeState > 0) & + thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) + if(damageState(material_homogenizationAt(e))%sizeState > 0) & damageState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = & - damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) ! ...internal damage state + damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) endif endif converged From c7036e39701c3e7b0fcbf4cfc80db47cc76ac92c Mon Sep 17 00:00:00 2001 From: Test User Date: Fri, 17 May 2019 09:06:28 +0200 Subject: [PATCH 21/29] [skip ci] updated version information after successful test of v2.0.3-304-g7b14263c --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index 69faed072..2d8c83361 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v2.0.3-301-g789420c9 +v2.0.3-304-g7b14263c From 1de0c7e652fcc6572859976e645374b6a36bd1fb Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 18 May 2019 14:19:31 +0200 Subject: [PATCH 22/29] wrong name --- src/plastic_dislotwin.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/plastic_dislotwin.f90 b/src/plastic_dislotwin.f90 index 7db1f5f7f..05e9a49fd 100644 --- a/src/plastic_dislotwin.f90 +++ b/src/plastic_dislotwin.f90 @@ -1068,7 +1068,7 @@ subroutine plastic_dislotwin_results(instance,group) call results_writeDataset(group,stt%rho_dip,'rho_dip',& 'dislocation dipole density''1/m²') case (dot_gamma_sl_ID) - call results_writeDataset(group,stt%gamma_sl,'dot_gamma_sl',& + call results_writeDataset(group,stt%gamma_sl,'gamma_sl',& 'plastic shear','1') case (Lambda_sl_ID) call results_writeDataset(group,dst%Lambda_sl,'Lambda_sl',& From c4f07a9ad97102f9e8065dcdab06a1283298fbf3 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 18 May 2019 21:25:05 +0200 Subject: [PATCH 23/29] need to correct tensor order --- src/results.f90 | 20 ++++++++++++++------ 1 file changed, 14 insertions(+), 6 deletions(-) diff --git a/src/results.f90 b/src/results.f90 index c5582b927..05db831f7 100644 --- a/src/results.f90 +++ b/src/results.f90 @@ -299,18 +299,26 @@ end subroutine results_writeVectorDataset_real !-------------------------------------------------------------------------------------------------- subroutine results_writeTensorDataset_real(group,dataset,label,description,SIunit) - character(len=*), intent(in) :: label,group,description - character(len=*), intent(in), optional :: SIunit - real(pReal), intent(inout), dimension(:,:,:) :: dataset + character(len=*), intent(in) :: label,group,description + character(len=*), intent(in), optional :: SIunit + real(pReal), intent(in), dimension(:,:,:) :: dataset + integer :: i integer(HID_T) :: groupHandle - + real(pReal), dimension(:,:,:), allocatable :: dataset_transposed + + + allocate(dataset_transposed,mold=dataset) + do i=1,size(dataset,3) + dataset_transposed(1:3,1:3,i) = transpose(dataset(1:3,1:3,i)) + enddo + groupHandle = results_openGroup(group) #ifdef PETSc - call HDF5_write(groupHandle,dataset,label,.true.) + call HDF5_write(groupHandle,dataset_transposed,label,.true.) #else - call HDF5_write(groupHandle,dataset,label,.false.) + call HDF5_write(groupHandle,dataset_transposed,label,.false.) #endif if (HDF5_objectExists(groupHandle,label)) & From d2e64df6a1813b4b169fdf3fbeb3d401ab7d69a1 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 18 May 2019 21:27:27 +0200 Subject: [PATCH 24/29] some tests based on HDF5 output --- PRIVATE | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PRIVATE b/PRIVATE index 639c6f4a5..3a2f89547 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 639c6f4a5eafc893c83c740c57f417eaaabc45ae +Subproject commit 3a2f89547c264044a7bfab9d33aee78eec495a76 From b35465b591f4d4866f421045d03d1ac7da9a3a99 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 20 May 2019 19:08:56 +0200 Subject: [PATCH 25/29] gamma_slip_ID should be used to write result --- src/plastic_dislotwin.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/plastic_dislotwin.f90 b/src/plastic_dislotwin.f90 index 05e9a49fd..cb4b3cbae 100644 --- a/src/plastic_dislotwin.f90 +++ b/src/plastic_dislotwin.f90 @@ -1067,7 +1067,7 @@ subroutine plastic_dislotwin_results(instance,group) case (rho_dip_ID) call results_writeDataset(group,stt%rho_dip,'rho_dip',& 'dislocation dipole density''1/m²') - case (dot_gamma_sl_ID) + case (gamma_sl_ID) call results_writeDataset(group,stt%gamma_sl,'gamma_sl',& 'plastic shear','1') case (Lambda_sl_ID) From e6cec6ecbe32d601f03579eabee202531f87d4c0 Mon Sep 17 00:00:00 2001 From: Philip Eisenlohr Date: Thu, 23 May 2019 12:03:54 -0400 Subject: [PATCH 26/29] added option to reverse inside/outside of primitive body --- processing/pre/geom_addPrimitive.py | 43 +++++++++++++++++------------ 1 file changed, 26 insertions(+), 17 deletions(-) diff --git a/processing/pre/geom_addPrimitive.py b/processing/pre/geom_addPrimitive.py index 7fcfdbc5c..0dfd06732 100755 --- a/processing/pre/geom_addPrimitive.py +++ b/processing/pre/geom_addPrimitive.py @@ -43,7 +43,7 @@ parser.add_option('-e', '--exponent', dest='exponent', 1 gives a sphere (|x|^(2^1) + |y|^(2^1) + |z|^(2^1) < 1), \ large values produce boxes, negative turns concave.') parser.add_option('-f', '--fill', dest='fill', - type='int', metavar = 'int', + type='float', metavar = 'float', help='grain index to fill primitive. "0" selects maximum microstructure index + 1 [%default]') parser.add_option('-q', '--quaternion', dest='quaternion', type='float', nargs = 4, metavar=' '.join(['float']*4), @@ -60,15 +60,24 @@ parser.add_option( '--nonperiodic', dest='periodic', parser.add_option( '--realspace', dest='realspace', action='store_true', help = '-c and -d span [origin,origin+size] instead of [0,grid] coordinates') +parser.add_option( '--invert', dest='inside', + action='store_false', + help = 'invert the volume filled by the primitive (inside/outside)') +parser.add_option('--float', dest = 'float', + action = 'store_true', + help = 'use float input') parser.set_defaults(center = (.0,.0,.0), - fill = 0, + fill = 0.0, degrees = False, exponent = (20,20,20), # box shape by default periodic = True, realspace = False, + inside = True, + float = False, ) (options, filenames) = parser.parse_args() + if options.dimension is None: parser.error('no dimension specified.') if options.angleaxis is not None: @@ -78,6 +87,8 @@ elif options.quaternion is not None: else: rotation = damask.Rotation() +datatype = 'f' if options.float else 'i' + options.center = np.array(options.center) options.dimension = np.array(options.dimension) # undo logarithmic sense of exponent and generate ellipsoids for negative dimensions (backward compatibility) @@ -97,13 +108,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: %s'%(' x '.join(map(str,info['grid']))), - 'size x y z: %s'%(' x '.join(map(str,info['size']))), - 'origin x y z: %s'%(' : '.join(map(str,info['origin']))), - 'homogenization: %i'%info['homogenization'], - 'microstructures: %i'%info['microstructures'], - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') @@ -115,7 +120,7 @@ for name in filenames: #--- read data ------------------------------------------------------------------------------------ - microstructure = table.microstructure_read(info['grid']) # read microstructure + microstructure = table.microstructure_read(info['grid'],datatype) # read microstructure # --- do work ------------------------------------------------------------------------------------ @@ -123,7 +128,7 @@ for name in filenames: 'microstructures': 0, } - options.fill = microstructure.max()+1 if options.fill == 0 else options.fill + options.fill = np.nanmax(microstructure)+1 if options.fill == 0 else options.fill microstructure = microstructure.reshape(info['grid'],order='F') @@ -193,19 +198,23 @@ for name in filenames: grid[1] * j : grid[1] * (j+1), grid[2] * k : grid[2] * (k+1)])**options.exponent[2] <= 1.0) - microstructure = np.where(inside, options.fill, microstructure) + microstructure = np.where(inside, + options.fill if options.inside else microstructure, + microstructure if options.inside else options.fill) else: # nonperiodic, much lighter on resources microstructure = np.where(np.abs(X)**options.exponent[0] + np.abs(Y)**options.exponent[1] + - np.abs(Z)**options.exponent[2] <= 1.0, options.fill, microstructure) + np.abs(Z)**options.exponent[2] <= 1.0, + options.fill if options.inside else microstructure, + microstructure if options.inside else options.fill) np.seterr(**old_settings) # Reset warnings to old state - newInfo['microstructures'] = microstructure.max() + newInfo['microstructures'] = len(np.unique(microstructure)) # --- report --------------------------------------------------------------------------------------- if (newInfo['microstructures'] != info['microstructures']): - damask.util.croak('--> microstructures: %i'%newInfo['microstructures']) + damask.util.croak('--> microstructures: {}'.format(newInfo['microstructures'])) #--- write header --------------------------------------------------------------------------------- @@ -225,9 +234,9 @@ for name in filenames: # --- write microstructure information ------------------------------------------------------------ - formatwidth = int(math.floor(math.log10(microstructure.max())+1)) + format = '%g' if options.float else '%{}i'.format(int(math.floor(math.log10(np.nanmax(microstructure))+1))) table.data = microstructure.reshape((info['grid'][0],info['grid'][1]*info['grid'][2]),order='F').transpose() - table.data_writeArray('%%%ii'%(formatwidth),delimiter = ' ') + table.data_writeArray(format,delimiter = ' ') #--- output finalization -------------------------------------------------------------------------- From eb13fbc0cee296b9448a8f2eca5ec8c87bda3c9f Mon Sep 17 00:00:00 2001 From: Philip Eisenlohr Date: Thu, 23 May 2019 13:03:24 -0400 Subject: [PATCH 27/29] streamlined geom-info reporting; added --float option to some scripts; hardened against NaN; "microstructures" now reports uniques not max --- processing/pre/geom_canvas.py | 14 ++++----- processing/pre/geom_clean.py | 14 +++------ processing/pre/geom_fromMinimalSurface.py | 7 +---- processing/pre/geom_fromTable.py | 9 ++---- processing/pre/geom_grainGrowth.py | 13 ++------- processing/pre/geom_mirror.py | 26 ++++++++++------- processing/pre/geom_pack.py | 10 ++----- processing/pre/geom_renumber.py | 10 ++----- processing/pre/geom_rescale.py | 26 +++++++++-------- processing/pre/geom_rotate.py | 35 ++++++++++++----------- processing/pre/geom_toTable.py | 14 +++------ processing/pre/geom_translate.py | 22 ++++++-------- processing/pre/geom_unpack.py | 8 +----- processing/pre/geom_vicinityOffset.py | 14 +++------ 14 files changed, 85 insertions(+), 137 deletions(-) diff --git a/processing/pre/geom_canvas.py b/processing/pre/geom_canvas.py index d7fd1614a..01682deb8 100755 --- a/processing/pre/geom_canvas.py +++ b/processing/pre/geom_canvas.py @@ -35,7 +35,7 @@ parser.add_option('-f', type = 'float', metavar = 'float', help = '(background) canvas grain index. "0" selects maximum microstructure index + 1 [%default]') parser.add_option('--float', - dest = 'real', + dest = 'float', action = 'store_true', help = 'use float input') parser.add_option('--blank', @@ -45,13 +45,13 @@ parser.add_option('--blank', parser.set_defaults(grid = ['0','0','0'], offset = (0,0,0), - fill = 0, - real = False, + fill = 0.0, + float = False, ) (options, filenames) = parser.parse_args() -datatype = 'f' if options.real else 'i' +datatype = 'f' if options.float else 'i' options.grid = ['1','1','1'] if options.blank and options.grid == ['0','0','0'] else options.grid options.fill = 1 if options.blank and options.fill == 0 else options.fill @@ -107,7 +107,7 @@ for name in filenames: newInfo['grid'] = np.where(newInfo['grid'] > 0, newInfo['grid'],info['grid']) microstructure_cropped = np.zeros(newInfo['grid'],datatype) - microstructure_cropped.fill(options.fill if options.real or options.fill > 0 else microstructure.max()+1) + microstructure_cropped.fill(options.fill if options.float or options.fill > 0 else np.nanmax(microstructure)+1) if not options.blank: xindex = np.arange(max(options.offset[0],0),min(options.offset[0]+newInfo['grid'][0],info['grid'][0])) @@ -130,7 +130,7 @@ for name in filenames: newInfo['size'] = info['size']/info['grid']*newInfo['grid'] if np.all(info['grid'] > 0) else newInfo['grid'] newInfo['origin'] = info['origin']+(info['size']/info['grid'] if np.all(info['grid'] > 0) \ else newInfo['size']/newInfo['grid'])*options.offset - newInfo['microstructures'] = microstructure_cropped.max() + newInfo['microstructures'] = len(np.unique(microstructure_cropped)) # --- report --------------------------------------------------------------------------------------- @@ -172,7 +172,7 @@ for name in filenames: # --- write microstructure information ------------------------------------------------------------ - format = '%g' if options.real else '%{}i'.format(int(math.floor(math.log10(microstructure_cropped.max())+1))) + format = '%g' if options.float else '%{}i'.format(int(math.floor(math.log10(np.nanmax(microstructure_cropped))+1))) table.data = microstructure_cropped.reshape((newInfo['grid'][0],newInfo['grid'][1]*newInfo['grid'][2]),order='F').transpose() table.data_writeArray(format,delimiter=' ') diff --git a/processing/pre/geom_clean.py b/processing/pre/geom_clean.py index 907431146..1d0769ab3 100755 --- a/processing/pre/geom_clean.py +++ b/processing/pre/geom_clean.py @@ -50,13 +50,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: {}'.format(' x '.join(map(str,info['grid']))), - 'size x y z: {}'.format(' x '.join(map(str,info['size']))), - 'origin x y z: {}'.format(' : '.join(map(str,info['origin']))), - 'homogenization: {}'.format(info['homogenization']), - 'microstructures: {}'.format(info['microstructures']), - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') @@ -73,7 +67,7 @@ for name in filenames: # --- do work ------------------------------------------------------------------------------------ microstructure = ndimage.filters.generic_filter(microstructure,mostFrequent,size=(options.stencil,)*3).astype('int_') - newInfo = {'microstructures': microstructure.max()} + newInfo = {'microstructures': len(np.unique(microstructure))} # --- report --------------------------------------------------------------------------------------- if ( newInfo['microstructures'] != info['microstructures']): @@ -91,9 +85,9 @@ for name in filenames: # --- write microstructure information ------------------------------------------------------------ - formatwidth = int(math.floor(math.log10(microstructure.max())+1)) + formatwidth = int(math.floor(math.log10(np.nanmax(microstructure))+1)) table.data = microstructure.reshape((info['grid'][0],np.prod(info['grid'][1:])),order='F').transpose() - table.data_writeArray('%%%ii'%(formatwidth),delimiter = ' ') + table.data_writeArray('%{}i'.format(formatwidth),delimiter = ' ') # --- output finalization -------------------------------------------------------------------------- diff --git a/processing/pre/geom_fromMinimalSurface.py b/processing/pre/geom_fromMinimalSurface.py index 002b4800b..e0023e7ec 100755 --- a/processing/pre/geom_fromMinimalSurface.py +++ b/processing/pre/geom_fromMinimalSurface.py @@ -90,12 +90,7 @@ for name in filenames: #--- report --------------------------------------------------------------------------------------- - damask.util.croak(['grid a b c: %s'%(' x '.join(map(str,info['grid']))), - 'size x y z: %s'%(' x '.join(map(str,info['size']))), - 'origin x y z: %s'%(' : '.join(map(str,info['origin']))), - 'homogenization: %i'%info['homogenization'], - 'microstructures: %i'%info['microstructures'], - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') diff --git a/processing/pre/geom_fromTable.py b/processing/pre/geom_fromTable.py index 8eb1ed8bf..7a905cd26 100755 --- a/processing/pre/geom_fromTable.py +++ b/processing/pre/geom_fromTable.py @@ -192,12 +192,7 @@ for name in filenames: 'homogenization': options.homogenization, } - damask.util.croak(['grid a b c: {}'.format(' x '.join(map(str,info['grid']))), - 'size x y z: {}'.format(' x '.join(map(str,info['size']))), - 'origin x y z: {}'.format(' : '.join(map(str,info['origin']))), - 'homogenization: {}'.format(info['homogenization']), - 'microstructures: {}'.format(info['microstructures']), - ]) + damask.util.report_geom(info) # --- write header --------------------------------------------------------------------------------- @@ -230,7 +225,7 @@ for name in filenames: # --- write microstructure information ------------------------------------------------------------ table.data = grain.reshape(info['grid'][1]*info['grid'][2],info['grid'][0]) - table.data_writeArray('%%%ii'%(formatwidth),delimiter=' ') + table.data_writeArray('%{}i'.format(formatwidth),delimiter=' ') #--- output finalization -------------------------------------------------------------------------- diff --git a/processing/pre/geom_grainGrowth.py b/processing/pre/geom_grainGrowth.py index 1afb02715..f7c50c2e5 100755 --- a/processing/pre/geom_grainGrowth.py +++ b/processing/pre/geom_grainGrowth.py @@ -69,13 +69,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: {}'.format(' x '.join(list(map(str,info['grid'])))), - 'size x y z: {}'.format(' x '.join(list(map(str,info['size'])))), - 'origin x y z: {}'.format(' : '.join(list(map(str,info['origin'])))), - 'homogenization: {}'.format(info['homogenization']), - 'microstructures: {}'.format(info['microstructures']), - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') @@ -200,8 +194,7 @@ for name in filenames: newID += 1 microstructure = np.where(microstructure == microstructureID, newID, microstructure) - newInfo = {'microstructures': 0,} - newInfo['microstructures'] = microstructure.max() + newInfo = {'microstructures': len(np.unique(microstructure)),} # --- report -------------------------------------------------------------------------------------- @@ -226,7 +219,7 @@ for name in filenames: # --- write microstructure information ------------------------------------------------------------ - formatwidth = int(math.floor(math.log10(microstructure.max())+1)) + formatwidth = int(math.floor(math.log10(np.nanmax(microstructure))+1)) table.data = microstructure[::1 if info['grid'][0]>1 else 2, ::1 if info['grid'][1]>1 else 2, ::1 if info['grid'][2]>1 else 2,].\ diff --git a/processing/pre/geom_mirror.py b/processing/pre/geom_mirror.py index 951fb0842..853b99632 100755 --- a/processing/pre/geom_mirror.py +++ b/processing/pre/geom_mirror.py @@ -23,6 +23,13 @@ parser.add_option('-d','--direction', dest = 'directions', action = 'extend', metavar = '', help = "directions in which to mirror {'x','y','z'}") +parser.add_option('--float', + dest = 'float', + action = 'store_true', + help = 'use float input') + +parser.set_defaults(float = False, + ) (options, filenames) = parser.parse_args() @@ -32,6 +39,8 @@ if not set(options.directions).issubset(validDirections): invalidDirections = [str(e) for e in set(options.directions).difference(validDirections)] parser.error('invalid directions {}. '.format(*invalidDirections)) +datatype = 'f' if options.float else 'i' + # --- loop over input files ------------------------------------------------------------------------- if filenames == []: filenames = [None] @@ -39,7 +48,8 @@ if filenames == []: filenames = [None] for name in filenames: try: table = damask.ASCIItable(name = name, - buffered = False, labeled = False) + buffered = False, + labeled = False) except: continue damask.util.report(scriptName,name) @@ -47,13 +57,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: %s'%(' x '.join(map(str,info['grid']))), - 'size x y z: %s'%(' x '.join(map(str,info['size']))), - 'origin x y z: %s'%(' : '.join(map(str,info['origin']))), - 'homogenization: %i'%info['homogenization'], - 'microstructures: %i'%info['microstructures'], - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') @@ -65,7 +69,7 @@ for name in filenames: # --- read data ------------------------------------------------------------------------------------ - microstructure = table.microstructure_read(info['grid']).reshape(info['grid'],order='F') # read microstructure + microstructure = table.microstructure_read(info['grid'],datatype).reshape(info['grid'],order='F') # read microstructure if 'z' in options.directions: microstructure = np.concatenate([microstructure,microstructure[:,:,::-1]],2) @@ -107,9 +111,9 @@ for name in filenames: # --- write microstructure information ------------------------------------------------------------ - formatwidth = int(math.floor(math.log10(microstructure.max())+1)) + formatwidth = int(math.floor(math.log10(np.nanmax(microstructure))+1)) table.data = microstructure.reshape((newInfo['grid'][0],np.prod(newInfo['grid'][1:])),order='F').transpose() - table.data_writeArray('%%%ii'%(formatwidth),delimiter = ' ') + table.data_writeArray('%{}i'.format(formatwidth),delimiter = ' ') # --- output finalization -------------------------------------------------------------------------- diff --git a/processing/pre/geom_pack.py b/processing/pre/geom_pack.py index 0d864bbf5..2e6080a6b 100755 --- a/processing/pre/geom_pack.py +++ b/processing/pre/geom_pack.py @@ -35,14 +35,8 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: {}'.format(' x '.join(map(str,info['grid']))), - 'size x y z: {}'.format(' x '.join(map(str,info['size']))), - 'origin x y z: {}'.format(' : '.join(map(str,info['origin']))), - 'homogenization: {}'.format(info['homogenization']), - 'microstructures: {}'.format(info['microstructures']), - ]) - + damask.util.report_geom(info) + errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') if np.any(info['size'] <= 0.0): errors.append('invalid size x y z.') diff --git a/processing/pre/geom_renumber.py b/processing/pre/geom_renumber.py index 033b4a566..3faa7f449 100755 --- a/processing/pre/geom_renumber.py +++ b/processing/pre/geom_renumber.py @@ -35,13 +35,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: %s'%(' x '.join(map(str,info['grid']))), - 'size x y z: %s'%(' x '.join(map(str,info['size']))), - 'origin x y z: %s'%(' : '.join(map(str,info['origin']))), - 'homogenization: %i'%info['homogenization'], - 'microstructures: %i'%info['microstructures'], - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') @@ -93,7 +87,7 @@ for name in filenames: # --- write microstructure information ----------------------------------------------------------- - format = '%{}i'.format(int(math.floor(math.log10(newInfo['microstructures'])+1))) + format = '%{}i'.format(int(math.floor(math.log10(np.nanmax(renumbered))+1))) table.data = renumbered.reshape((info['grid'][0],info['grid'][1]*info['grid'][2]),order='F').transpose() table.data_writeArray(format,delimiter = ' ') diff --git a/processing/pre/geom_rescale.py b/processing/pre/geom_rescale.py index b3716bd62..4a14c0050 100755 --- a/processing/pre/geom_rescale.py +++ b/processing/pre/geom_rescale.py @@ -31,14 +31,21 @@ parser.add_option('-r', '--renumber', dest = 'renumber', action = 'store_true', help = 'renumber microstructure indices from 1..N [%default]') +parser.add_option('--float', + dest = 'float', + action = 'store_true', + help = 'use float input') parser.set_defaults(renumber = False, grid = ['0','0','0'], size = ['0.0','0.0','0.0'], + float = False, ) (options, filenames) = parser.parse_args() +datatype = 'f' if options.float else 'i' + # --- loop over input files ------------------------------------------------------------------------- if filenames == []: filenames = [None] @@ -46,7 +53,8 @@ if filenames == []: filenames = [None] for name in filenames: try: table = damask.ASCIItable(name = name, - buffered = False, labeled = False) + buffered = False, + labeled = False) except: continue damask.util.report(scriptName,name) @@ -54,13 +62,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: %s'%(' x '.join(map(str,info['grid']))), - 'size x y z: %s'%(' x '.join(map(str,info['size']))), - 'origin x y z: %s'%(' : '.join(map(str,info['origin']))), - 'homogenization: %i'%info['homogenization'], - 'microstructures: %i'%info['microstructures'], - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') @@ -72,7 +74,7 @@ for name in filenames: # --- read data ------------------------------------------------------------------------------------ - microstructure = table.microstructure_read(info['grid']) # read microstructure + microstructure = table.microstructure_read(info['grid'],datatype) # read microstructure # --- do work ------------------------------------------------------------------------------------ @@ -113,7 +115,7 @@ for name in filenames: newID += 1 microstructure = np.where(microstructure == microstructureID, newID,microstructure).reshape(microstructure.shape) - newInfo['microstructures'] = microstructure.max() + newInfo['microstructures'] = len(np.unique(microstructure)) # --- report --------------------------------------------------------------------------------------- @@ -152,9 +154,9 @@ for name in filenames: # --- write microstructure information ------------------------------------------------------------ - formatwidth = int(math.floor(math.log10(microstructure.max())+1)) + format = '%g' if options.float else '%{}i'.format(int(math.floor(math.log10(np.nanmax(microstructure))+1))) table.data = microstructure.reshape((newInfo['grid'][0],newInfo['grid'][1]*newInfo['grid'][2]),order='F').transpose() - table.data_writeArray('%%%ii'%(formatwidth),delimiter = ' ') + table.data_writeArray(format,delimiter=' ') # --- output finalization -------------------------------------------------------------------------- diff --git a/processing/pre/geom_rotate.py b/processing/pre/geom_rotate.py index 4da59cddf..7cce5800d 100755 --- a/processing/pre/geom_rotate.py +++ b/processing/pre/geom_rotate.py @@ -43,9 +43,15 @@ parser.add_option('-f', '--fill', dest = 'fill', type = 'int', metavar = 'int', help = 'background grain index. "0" selects maximum microstructure index + 1 [%default]') +parser.add_option('--float', + dest = 'float', + action = 'store_true', + help = 'use float input') parser.set_defaults(degrees = False, - fill = 0) + fill = 0, + float = False, + ) (options, filenames) = parser.parse_args() @@ -61,6 +67,8 @@ if options.matrix is not None: if options.eulers is not None: eulers = damask.Rotation.fromEulers(np.array(options.eulers),degrees=True).asEulers(degrees=True) +datatype = 'f' if options.float else 'i' + # --- loop over input files ------------------------------------------------------------------------- if filenames == []: filenames = [None] @@ -77,13 +85,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: {}'.format(' x '.join(map(str,info['grid']))), - 'size x y z: {}'.format(' x '.join(map(str,info['size']))), - 'origin x y z: {}'.format(' : '.join(map(str,info['origin']))), - 'homogenization: {}'.format(info['homogenization']), - 'microstructures: {}'.format(info['microstructures']), - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') @@ -95,9 +97,9 @@ for name in filenames: # --- read data ------------------------------------------------------------------------------------ - microstructure = table.microstructure_read(info['grid']).reshape(info['grid'],order='F') # read microstructure + microstructure = table.microstructure_read(info['grid'],datatype).reshape(info['grid'],order='F') # read microstructure - newGrainID = options.fill if options.fill != 0 else microstructure.max()+1 + newGrainID = options.fill if options.fill != 0 else np.nanmax(microstructure)+1 microstructure = ndimage.rotate(microstructure,eulers[2],(0,1),order=0,prefilter=False,output=int,cval=newGrainID) # rotation around Z microstructure = ndimage.rotate(microstructure,eulers[1],(1,2),order=0,prefilter=False,output=int,cval=newGrainID) # rotation around X microstructure = ndimage.rotate(microstructure,eulers[0],(0,1),order=0,prefilter=False,output=int,cval=newGrainID) # rotation around Z @@ -107,19 +109,18 @@ for name in filenames: newInfo = { 'size': microstructure.shape*info['size']/info['grid'], 'grid': microstructure.shape, - 'microstructures': microstructure.max(), + 'microstructures': len(np.unique(microstructure)), } - # --- report --------------------------------------------------------------------------------------- remarks = [] if (any(newInfo['grid'] != info['grid'])): - remarks.append('--> grid a b c: %s'%(' x '.join(map(str,newInfo['grid'])))) + remarks.append('--> grid a b c: {}'.format(' x '.join(map(str,newInfo['grid'])))) if (any(newInfo['size'] != info['size'])): - remarks.append('--> size x y z: %s'%(' x '.join(map(str,newInfo['size'])))) + remarks.append('--> size x y z: {}'.format(' x '.join(map(str,newInfo['size'])))) if ( newInfo['microstructures'] != info['microstructures']): - remarks.append('--> microstructures: %i'%newInfo['microstructures']) + remarks.append('--> microstructures: {}'.format(newInfo['microstructures'])) if remarks != []: damask.util.croak(remarks) # --- write header --------------------------------------------------------------------------------- @@ -138,9 +139,9 @@ for name in filenames: # --- write microstructure information ------------------------------------------------------------ - formatwidth = int(math.floor(math.log10(microstructure.max())+1)) + format = '%g' if options.float else '%{}i'.format(int(math.floor(math.log10(np.nanmax(microstructure))+1))) table.data = microstructure.reshape((newInfo['grid'][0],np.prod(newInfo['grid'][1:])),order='F').transpose() - table.data_writeArray('%%%ii'%(formatwidth),delimiter = ' ') + table.data_writeArray(format,delimiter=' ') # --- output finalization -------------------------------------------------------------------------- diff --git a/processing/pre/geom_toTable.py b/processing/pre/geom_toTable.py index 73e4888d1..0a71b335e 100755 --- a/processing/pre/geom_toTable.py +++ b/processing/pre/geom_toTable.py @@ -20,15 +20,15 @@ Translate geom description into ASCIItable containing position and microstructur """, version = scriptID) parser.add_option('--float', - dest = 'real', + dest = 'float', action = 'store_true', help = 'use float input') -parser.set_defaults(real = False, +parser.set_defaults(float = False, ) (options, filenames) = parser.parse_args() -datatype = 'f' if options.real else 'i' +datatype = 'f' if options.float else 'i' # --- loop over input files ------------------------------------------------------------------------- @@ -47,13 +47,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: {}'.format(' x '.join(list(map(str,info['grid'])))), - 'size x y z: {}'.format(' x '.join(list(map(str,info['size'])))), - 'origin x y z: {}'.format(' : '.join(list(map(str,info['origin'])))), - 'homogenization: {}'.format(info['homogenization']), - 'microstructures: {}'.format(info['microstructures']), - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') diff --git a/processing/pre/geom_translate.py b/processing/pre/geom_translate.py index 59aaac5d5..072c270ea 100755 --- a/processing/pre/geom_translate.py +++ b/processing/pre/geom_translate.py @@ -31,19 +31,19 @@ parser.add_option('-s', '--substitute', action = 'extend', metavar = '', help = 'substitutions of microstructure indices from,to,from,to,...') parser.add_option('--float', - dest = 'real', + dest = 'float', action = 'store_true', help = 'use float input') parser.set_defaults(origin = (0.0,0.0,0.0), microstructure = 0, substitute = [], - real = False, + float = False, ) (options, filenames) = parser.parse_args() -datatype = 'f' if options.real else 'i' +datatype = 'f' if options.float else 'i' sub = {} for i in range(len(options.substitute)//2): # split substitution list into "from" -> "to" @@ -64,13 +64,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: %s'%(' x '.join(map(str,info['grid']))), - 'size x y z: %s'%(' x '.join(map(str,info['size']))), - 'origin x y z: %s'%(' : '.join(map(str,info['origin']))), - 'homogenization: %i'%info['homogenization'], - 'microstructures: %i'%info['microstructures'], - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') @@ -92,7 +86,7 @@ for name in filenames: } substituted = np.copy(microstructure) - for k, v in sub.items(): substituted[microstructure==k] = v # substitute microstructure indices + for k, v in sub.items(): substituted[microstructure==k] = v # substitute microstructure indices substituted += options.microstructure # shift microstructure indices @@ -103,9 +97,9 @@ for name in filenames: remarks = [] if (any(newInfo['origin'] != info['origin'])): - remarks.append('--> origin x y z: %s'%(' : '.join(map(str,newInfo['origin'])))) + remarks.append('--> origin x y z: {}'.format(' : '.join(map(str,newInfo['origin'])))) if ( newInfo['microstructures'] != info['microstructures']): - remarks.append('--> microstructures: %i'%newInfo['microstructures']) + remarks.append('--> microstructures: {}'.format(newInfo['microstructures'])) if remarks != []: damask.util.croak(remarks) # --- write header ------------------------------------------------------------------------------- @@ -124,7 +118,7 @@ for name in filenames: # --- write microstructure information ----------------------------------------------------------- - format = '%g' if options.real else '%{}i'.format(int(math.floor(math.log10(microstructure.max())+1))) + format = '%g' if options.float else '%{}i'.format(int(math.floor(math.log10(np.nanmax(substituted))+1))) table.data = substituted.reshape((info['grid'][0],info['grid'][1]*info['grid'][2]),order='F').transpose() table.data_writeArray(format,delimiter = ' ') diff --git a/processing/pre/geom_unpack.py b/processing/pre/geom_unpack.py index 726e4ef04..4cac76c5f 100755 --- a/processing/pre/geom_unpack.py +++ b/processing/pre/geom_unpack.py @@ -43,13 +43,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: %s'%(' x '.join(map(str,info['grid']))), - 'size x y z: %s'%(' x '.join(map(str,info['size']))), - 'origin x y z: %s'%(' : '.join(map(str,info['origin']))), - 'homogenization: %i'%info['homogenization'], - 'microstructures: %i'%info['microstructures'], - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') diff --git a/processing/pre/geom_vicinityOffset.py b/processing/pre/geom_vicinityOffset.py index 9fce7201a..733276d01 100755 --- a/processing/pre/geom_vicinityOffset.py +++ b/processing/pre/geom_vicinityOffset.py @@ -73,13 +73,7 @@ for name in filenames: table.head_read() info,extra_header = table.head_getGeom() - - damask.util.croak(['grid a b c: %s'%(' x '.join(map(str,info['grid']))), - 'size x y z: %s'%(' x '.join(map(str,info['size']))), - 'origin x y z: %s'%(' : '.join(map(str,info['origin']))), - 'homogenization: %i'%info['homogenization'], - 'microstructures: %i'%info['microstructures'], - ]) + damask.util.report_geom(info) errors = [] if np.any(info['grid'] < 1): errors.append('invalid grid a b c.') @@ -108,7 +102,7 @@ for name in filenames: extra_keywords={"trigger":options.trigger,"size":1+2*options.vicinity}), microstructure + options.offset,microstructure) - newInfo['microstructures'] = microstructure.max() + newInfo['microstructures'] = len(np.unique(microstructure)) # --- report --------------------------------------------------------------------------------------- @@ -131,9 +125,9 @@ for name in filenames: # --- write microstructure information ------------------------------------------------------------ - formatwidth = int(math.floor(math.log10(microstructure.max())+1)) + formatwidth = int(math.floor(math.log10(np.nanmax(microstructure))+1)) table.data = microstructure.reshape((info['grid'][0],info['grid'][1]*info['grid'][2]),order='F').transpose() - table.data_writeArray('%%%ii'%(formatwidth),delimiter = ' ') + table.data_writeArray('%{}i'.format(formatwidth),delimiter = ' ') # --- output finalization -------------------------------------------------------------------------- From e743a55ff2254a0b870ed65d60bbba0e6bba9f0c Mon Sep 17 00:00:00 2001 From: Test User Date: Thu, 23 May 2019 23:01:28 +0200 Subject: [PATCH 28/29] [skip ci] updated version information after successful test of v2.0.3-307-geb13fbc0 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index 2d8c83361..98f27ab76 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v2.0.3-304-g7b14263c +v2.0.3-307-geb13fbc0 From 144361295d362a7179adfaf76184f74b176e2208 Mon Sep 17 00:00:00 2001 From: Test User Date: Mon, 27 May 2019 15:39:37 +0200 Subject: [PATCH 29/29] [skip ci] updated version information after successful test of v2.0.3-332-g5abcca50 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index 98f27ab76..d961b0f73 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v2.0.3-307-geb13fbc0 +v2.0.3-332-g5abcca50