From 869307c5ecda89b712267de5103b8c2af16cf220 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 7 Feb 2021 15:05:38 +0100 Subject: [PATCH 01/10] temperature not needed for pure mechanics would also require to define thermal homogenization --- PRIVATE | 2 +- examples/SpectralMethod/Polycrystal/shearXY.yaml | 3 --- examples/SpectralMethod/Polycrystal/shearZX.yaml | 3 --- examples/SpectralMethod/Polycrystal/tensionX.yaml | 3 --- python/damask/_configmaterial.py | 14 +++++++------- 5 files changed, 8 insertions(+), 17 deletions(-) diff --git a/PRIVATE b/PRIVATE index d90babadf..91a4329a6 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit d90babadfb0a33afa1b793044cc4efd4d7430731 +Subproject commit 91a4329a6fe073fc2ef17e5c2c8e2f796e3b897b diff --git a/examples/SpectralMethod/Polycrystal/shearXY.yaml b/examples/SpectralMethod/Polycrystal/shearXY.yaml index fe746e7ac..58471e8f1 100644 --- a/examples/SpectralMethod/Polycrystal/shearXY.yaml +++ b/examples/SpectralMethod/Polycrystal/shearXY.yaml @@ -3,9 +3,6 @@ solver: mechanical: spectral_basic -initial_conditions: - T: 300 #in Kelvin - loadstep: - boundary_conditions: mechanical: diff --git a/examples/SpectralMethod/Polycrystal/shearZX.yaml b/examples/SpectralMethod/Polycrystal/shearZX.yaml index 5d977b4ba..a32fafb85 100644 --- a/examples/SpectralMethod/Polycrystal/shearZX.yaml +++ b/examples/SpectralMethod/Polycrystal/shearZX.yaml @@ -3,9 +3,6 @@ solver: mechanical: spectral_basic -initial_conditions: - T: 300 #in Kelvin - loadstep: - boundary_conditions: mechanical: diff --git a/examples/SpectralMethod/Polycrystal/tensionX.yaml b/examples/SpectralMethod/Polycrystal/tensionX.yaml index 0809fd53d..870755d58 100644 --- a/examples/SpectralMethod/Polycrystal/tensionX.yaml +++ b/examples/SpectralMethod/Polycrystal/tensionX.yaml @@ -3,9 +3,6 @@ solver: mechanical: spectral_basic -initial_conditions: - T: 300 #in Kelvin - loadstep: - boundary_conditions: mechanical: diff --git a/python/damask/_configmaterial.py b/python/damask/_configmaterial.py index d3ed972e3..0adc494b8 100644 --- a/python/damask/_configmaterial.py +++ b/python/damask/_configmaterial.py @@ -168,11 +168,11 @@ class ConfigMaterial(Config): ok = False if 'material' in self: - for i,v in enumerate(self['material']): - if 'constituents' in v: - f = 0.0 - for c in v['constituents']: - f+= float(c['v']) + for i,m in enumerate(self['material']): + if 'constituents' in m: + v = 0.0 + for c in m['constituents']: + v+= float(c['v']) if 'O' in c: try: Rotation.from_quaternion(c['O']) @@ -180,8 +180,8 @@ class ConfigMaterial(Config): o = c['O'] print(f"Invalid orientation: '{o}' in material '{i}'") ok = False - if not np.isclose(f,1.0): - print(f"Invalid total fraction '{f}' in material '{i}'") + if not np.isclose(v,1.0): + print(f"Invalid total fraction (v) '{v}' in material '{i}'") ok = False return ok From 5b8e1996273572de5c36e8d1fcaec22ac89a6820 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 9 Feb 2021 23:09:41 +0100 Subject: [PATCH 02/10] avoid errors related to CRLF (windows) file endings --- src/IO.f90 | 21 +++++++++++++++++++-- 1 file changed, 19 insertions(+), 2 deletions(-) diff --git a/src/IO.f90 b/src/IO.f90 index fd87907aa..33ab7ee97 100644 --- a/src/IO.f90 +++ b/src/IO.f90 @@ -65,8 +65,8 @@ end subroutine IO_init function IO_readlines(fileName) result(fileContent) character(len=*), intent(in) :: fileName - character(len=pStringLen), dimension(:), allocatable :: fileContent !< file content, separated per lines + character(len=pStringLen) :: line character(len=:), allocatable :: rawData integer :: & @@ -75,6 +75,7 @@ function IO_readlines(fileName) result(fileContent) l logical :: warned + rawData = IO_read(fileName) !-------------------------------------------------------------------------------------------------- @@ -112,17 +113,20 @@ end function IO_readlines !-------------------------------------------------------------------------------------------------- !> @brief Read whole file. -!> @details ensures that the string ends with a new line (expected UNIX behavior) +!> @details ensures that the string ends with a new line (expected UNIX behavior) and rejects +! windows (CRLF) line endings !-------------------------------------------------------------------------------------------------- function IO_read(fileName) result(fileContent) character(len=*), intent(in) :: fileName character(len=:), allocatable :: fileContent + integer :: & fileLength, & fileUnit, & myStat + inquire(file = fileName, size=fileLength) open(newunit=fileUnit, file=fileName, access='stream',& status='old', position='rewind', action='read',iostat=myStat) @@ -137,6 +141,8 @@ function IO_read(fileName) result(fileContent) if(myStat /= 0) call IO_error(102,ext_msg=trim(fileName)) close(fileUnit) + if(scan(fileContent,achar(13)) /= 0) call IO_error(115) + if(fileContent(fileLength:fileLength) /= IO_EOL) fileContent = fileContent//IO_EOL ! ensure EOL@EOF end function IO_read @@ -151,6 +157,7 @@ logical pure function IO_isBlank(string) integer :: posNonBlank + posNonBlank = verify(string,IO_WHITESPACE) IO_isBlank = posNonBlank == 0 .or. posNonBlank == scan(string,IO_COMMENT) @@ -170,6 +177,7 @@ pure function IO_stringPos(string) integer :: left, right + allocate(IO_stringPos(1), source=0) right = 0 @@ -249,6 +257,7 @@ pure function IO_lc(string) integer :: i,n + do i=1,len(string) n = index(UPPER,string(i:i)) if(n/=0) then @@ -271,6 +280,7 @@ function IO_rmComment(line) character(len=:), allocatable :: IO_rmComment integer :: split + split = index(line,IO_COMMENT) if (split == 0) then @@ -292,6 +302,7 @@ integer function IO_stringAsInt(string) integer :: readStatus character(len=*), parameter :: VALIDCHARS = '0123456789+- ' + valid: if (verify(string,VALIDCHARS) == 0) then read(string,*,iostat=readStatus) IO_stringAsInt if (readStatus /= 0) call IO_error(111,ext_msg=string) @@ -313,6 +324,7 @@ real(pReal) function IO_stringAsFloat(string) integer :: readStatus character(len=*), parameter :: VALIDCHARS = '0123456789eE.+- ' + valid: if (verify(string,VALIDCHARS) == 0) then read(string,*,iostat=readStatus) IO_stringAsFloat if (readStatus /= 0) call IO_error(112,ext_msg=string) @@ -331,6 +343,7 @@ logical function IO_stringAsBool(string) character(len=*), intent(in) :: string !< string for conversion to int value + if (trim(adjustl(string)) == 'True' .or. trim(adjustl(string)) == 'true') then IO_stringAsBool = .true. elseif (trim(adjustl(string)) == 'False' .or. trim(adjustl(string)) == 'false') then @@ -356,6 +369,7 @@ subroutine IO_error(error_ID,el,ip,g,instance,ext_msg) character(len=:), allocatable :: msg character(len=pStringLen) :: formatString + select case (error_ID) !-------------------------------------------------------------------------------------------------- @@ -382,6 +396,9 @@ subroutine IO_error(error_ID,el,ip,g,instance,ext_msg) msg = 'invalid character for logical:' case (114) msg = 'cannot decode base64 string:' + case (115) + msg = 'found CR. Windows file endings (CRLF) are not supported.' + !-------------------------------------------------------------------------------------------------- ! lattice error messages From ef45e856a15cb1f265b08f5c8934a8fc28ad6489 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 10 Feb 2021 09:07:42 +0100 Subject: [PATCH 03/10] don't scan the whole file in case of proper line endings might lead to strange behavior if people randomly distribute CRs in their file. But that actually deserves to get strange behavior + Test --- PRIVATE | 2 +- src/IO.f90 | 8 ++++++-- 2 files changed, 7 insertions(+), 3 deletions(-) diff --git a/PRIVATE b/PRIVATE index d90babadf..5d27d879a 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit d90babadfb0a33afa1b793044cc4efd4d7430731 +Subproject commit 5d27d879abaeb1542e3f9f3065172be740af2899 diff --git a/src/IO.f90 b/src/IO.f90 index 33ab7ee97..36b774191 100644 --- a/src/IO.f90 +++ b/src/IO.f90 @@ -124,7 +124,9 @@ function IO_read(fileName) result(fileContent) integer :: & fileLength, & fileUnit, & - myStat + myStat, & + firstEOL + character, parameter :: CR = achar(13) inquire(file = fileName, size=fileLength) @@ -141,10 +143,12 @@ function IO_read(fileName) result(fileContent) if(myStat /= 0) call IO_error(102,ext_msg=trim(fileName)) close(fileUnit) - if(scan(fileContent,achar(13)) /= 0) call IO_error(115) if(fileContent(fileLength:fileLength) /= IO_EOL) fileContent = fileContent//IO_EOL ! ensure EOL@EOF + firstEOL = index(fileContent,IO_EOL) + if(scan(fileContent(firstEOL:firstEOL),CR) /= 0) call IO_error(115) + end function IO_read From 6895ef6b18b2cc3b3c26d0ab09c81e0115016b21 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 10 Feb 2021 10:03:35 +0100 Subject: [PATCH 04/10] always write LF line endings small pitfall: I windows users use a filehandle that results from a call to open() without the newline option, they get still CRLF line endings --- python/damask/_colormap.py | 16 ++++++++-------- python/damask/_config.py | 2 +- python/damask/_result.py | 2 +- python/damask/_table.py | 2 +- python/damask/solver/_marc.py | 2 +- 5 files changed, 12 insertions(+), 12 deletions(-) diff --git a/python/damask/_colormap.py b/python/damask/_colormap.py index 5a22f049b..8af14e081 100644 --- a/python/damask/_colormap.py +++ b/python/damask/_colormap.py @@ -238,7 +238,7 @@ class Colormap(mpl.colors.ListedColormap): fhandle = None else: try: - fhandle = open(fname,'w') + fhandle = open(fname,'w',newline='\n') except TypeError: fhandle = fname @@ -254,7 +254,7 @@ class Colormap(mpl.colors.ListedColormap): 'RGBPoints':colors }] - with open(self.name.replace(' ','_')+'.json', 'w') if fhandle is None else fhandle as f: + with open(self.name.replace(' ','_')+'.json','w',newline='\n') if fhandle is None else fhandle as f: json.dump(out, f,indent=4) @@ -273,14 +273,14 @@ class Colormap(mpl.colors.ListedColormap): fhandle = None else: try: - fhandle = open(fname,'w') + fhandle = open(fname,'w',newline='\n') except TypeError: fhandle = fname labels = {'RGBA':4} if self.colors.shape[1] == 4 else {'RGB': 3} t = Table(self.colors,labels,f'Creator: {util.execution_stamp("Colormap")}') - with open(self.name.replace(' ','_')+'.txt', 'w') if fhandle is None else fhandle as f: + with open(self.name.replace(' ','_')+'.txt','w',newline='\n') if fhandle is None else fhandle as f: t.save(f) @@ -299,7 +299,7 @@ class Colormap(mpl.colors.ListedColormap): fhandle = None else: try: - fhandle = open(fname,'w') + fhandle = open(fname,'w',newline='\n') except TypeError: fhandle = fname # ToDo: test in GOM @@ -308,7 +308,7 @@ class Colormap(mpl.colors.ListedColormap): + f'30 NO_UNIT 1 1 64 64 64 255 1 0 0 0 0 0 0 3 0 {len(self.colors)}' \ + ' '.join([f' 0 {c[0]} {c[1]} {c[2]} 255 1' for c in reversed((self.colors*255).astype(int))]) \ + '\n' - with open(self.name.replace(' ','_')+'.legend', 'w') if fhandle is None else fhandle as f: + with open(self.name.replace(' ','_')+'.legend','w',newline='\n') if fhandle is None else fhandle as f: f.write(GOM_str) @@ -327,14 +327,14 @@ class Colormap(mpl.colors.ListedColormap): fhandle = None else: try: - fhandle = open(fname,'w') + fhandle = open(fname,'w',newline='\n') except TypeError: fhandle = fname # ToDo: test in gmsh gmsh_str = 'View.ColorTable = {\n' \ +'\n'.join([f'{c[0]},{c[1]},{c[2]},' for c in self.colors[:,:3]*255]) \ +'\n}\n' - with open(self.name.replace(' ','_')+'.msh', 'w') if fhandle is None else fhandle as f: + with open(self.name.replace(' ','_')+'.msh','w',newline='\n') if fhandle is None else fhandle as f: f.write(gmsh_str) diff --git a/python/damask/_config.py b/python/damask/_config.py index 9aa031ff0..edf3a5676 100644 --- a/python/damask/_config.py +++ b/python/damask/_config.py @@ -75,7 +75,7 @@ class Config(dict): """ try: - fhandle = open(fname,'w') + fhandle = open(fname,'w',newline='\n') except TypeError: fhandle = fname diff --git a/python/damask/_result.py b/python/damask/_result.py index 3d8368911..6f2494299 100644 --- a/python/damask/_result.py +++ b/python/damask/_result.py @@ -1287,7 +1287,7 @@ class Result: np.prod(shape))} data_items[-1].text=f'{os.path.split(self.fname)[1]}:{name}' - with open(self.fname.with_suffix('.xdmf').name,'w') as f: + with open(self.fname.with_suffix('.xdmf').name,'w',newline='\n') as f: f.write(xml.dom.minidom.parseString(ET.tostring(xdmf).decode()).toprettyxml()) diff --git a/python/damask/_table.py b/python/damask/_table.py index 78a8a276e..e9acb371d 100644 --- a/python/damask/_table.py +++ b/python/damask/_table.py @@ -380,7 +380,7 @@ class Table: [f'# {comment}' for comment in self.comments] try: - fhandle = open(fname,'w') + fhandle = open(fname,'w',newline='\n') except TypeError: fhandle = fname diff --git a/python/damask/solver/_marc.py b/python/damask/solver/_marc.py index d4aadb7ff..f6c81da9f 100644 --- a/python/damask/solver/_marc.py +++ b/python/damask/solver/_marc.py @@ -66,7 +66,7 @@ class Marc: if logfile is not None: try: - f = open(logfile,'w+') + f = open(logfile,'w+',newline='\n') except TypeError: f = logfile else: From 4e31862f0ff9cae6c89687bdad26ff660dd216d4 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 10 Feb 2021 18:35:13 +0100 Subject: [PATCH 05/10] avoid repetition --- python/damask/_colormap.py | 75 ++++++++++++++++++-------------------- python/damask/_table.py | 2 +- 2 files changed, 36 insertions(+), 41 deletions(-) diff --git a/python/damask/_colormap.py b/python/damask/_colormap.py index 8af14e081..7e8860dae 100644 --- a/python/damask/_colormap.py +++ b/python/damask/_colormap.py @@ -223,25 +223,46 @@ class Colormap(mpl.colors.ListedColormap): return Colormap(np.array(rev.colors),rev.name[:-4] if rev.name.endswith('_r_r') else rev.name) - def save_paraview(self,fname=None): + def _get_file_handle(self,fname,extension): """ - Save as JSON file for use in Paraview. + Provide file handle. Parameters ---------- - fname : file, str, or pathlib.Path, optional. - Filename to store results. If not given, the filename will - consist of the name of the colormap and extension '.json'. + fname : file, str, pathlib.Path, or None + Filename or filehandle, will be name of the colormap+extension if None. + + extension: str + Extension of the filename. + + Returns + ------- + f + File handle """ if fname is None: - fhandle = None + fhandle = open(self.name.replace(' ','_')+'.'+extension,'w',newline='\n') else: try: fhandle = open(fname,'w',newline='\n') except TypeError: fhandle = fname + return fhandle + + + def save_paraview(self,fname=None): + """ + Save as JSON file for use in Paraview. + + Parameters + ---------- + fname : file, str, or pathlib.Path, optional + Filename to store results. If not given, the filename will + consist of the name of the colormap and extension '.json'. + + """ colors = [] for i,c in enumerate(np.round(self.colors,6).tolist()): colors+=[i]+c @@ -254,8 +275,7 @@ class Colormap(mpl.colors.ListedColormap): 'RGBPoints':colors }] - with open(self.name.replace(' ','_')+'.json','w',newline='\n') if fhandle is None else fhandle as f: - json.dump(out, f,indent=4) + json.dump(out,self._get_file_handle(fname,'json'),indent=4) def save_ASCII(self,fname=None): @@ -264,24 +284,14 @@ class Colormap(mpl.colors.ListedColormap): Parameters ---------- - fname : file, str, or pathlib.Path, optional. + fname : file, str, or pathlib.Path, optional Filename to store results. If not given, the filename will consist of the name of the colormap and extension '.txt'. """ - if fname is None: - fhandle = None - else: - try: - fhandle = open(fname,'w',newline='\n') - except TypeError: - fhandle = fname - labels = {'RGBA':4} if self.colors.shape[1] == 4 else {'RGB': 3} t = Table(self.colors,labels,f'Creator: {util.execution_stamp("Colormap")}') - - with open(self.name.replace(' ','_')+'.txt','w',newline='\n') if fhandle is None else fhandle as f: - t.save(f) + t.save(self._get_file_handle(fname,'txt')) def save_GOM(self,fname=None): @@ -290,26 +300,19 @@ class Colormap(mpl.colors.ListedColormap): Parameters ---------- - fname : file, str, or pathlib.Path, optional. + fname : file, str, or pathlib.Path, optional Filename to store results. If not given, the filename will consist of the name of the colormap and extension '.legend'. """ - if fname is None: - fhandle = None - else: - try: - fhandle = open(fname,'w',newline='\n') - except TypeError: - fhandle = fname # ToDo: test in GOM GOM_str = '1 1 {name} 9 {name} '.format(name=self.name.replace(" ","_")) \ + '0 1 0 3 0 0 -1 9 \\ 0 0 0 255 255 255 0 0 255 ' \ + f'30 NO_UNIT 1 1 64 64 64 255 1 0 0 0 0 0 0 3 0 {len(self.colors)}' \ + ' '.join([f' 0 {c[0]} {c[1]} {c[2]} 255 1' for c in reversed((self.colors*255).astype(int))]) \ + '\n' - with open(self.name.replace(' ','_')+'.legend','w',newline='\n') if fhandle is None else fhandle as f: - f.write(GOM_str) + + self._get_file_handle(fname,'legend').write(GOM_str) def save_gmsh(self,fname=None): @@ -318,24 +321,16 @@ class Colormap(mpl.colors.ListedColormap): Parameters ---------- - fname : file, str, or pathlib.Path, optional. + fname : file, str, or pathlib.Path, optional Filename to store results. If not given, the filename will consist of the name of the colormap and extension '.msh'. """ - if fname is None: - fhandle = None - else: - try: - fhandle = open(fname,'w',newline='\n') - except TypeError: - fhandle = fname # ToDo: test in gmsh gmsh_str = 'View.ColorTable = {\n' \ +'\n'.join([f'{c[0]},{c[1]},{c[2]},' for c in self.colors[:,:3]*255]) \ +'\n}\n' - with open(self.name.replace(' ','_')+'.msh','w',newline='\n') if fhandle is None else fhandle as f: - f.write(gmsh_str) + self._get_file_handle(fname,'msh').write(gmsh_str) @staticmethod diff --git a/python/damask/_table.py b/python/damask/_table.py index e9acb371d..68d6f94bf 100644 --- a/python/damask/_table.py +++ b/python/damask/_table.py @@ -26,7 +26,7 @@ class Table: comments_ = [comments] if isinstance(comments,str) else comments self.comments = [] if comments_ is None else [c for c in comments_] self.data = pd.DataFrame(data=data) - self.shapes = { k:(v,) if isinstance(v,(np.int,int)) else v for k,v in shapes.items() } + self.shapes = { k:(v,) if isinstance(v,(np.int64,np.int32,int)) else v for k,v in shapes.items() } self._label_uniform() def __repr__(self): From 992b4a7e6dab945b6d5dc3bb848e7ab310ffa0ab Mon Sep 17 00:00:00 2001 From: Test User Date: Thu, 11 Feb 2021 11:42:33 +0100 Subject: [PATCH 06/10] [skip ci] updated version information after successful test of v3.0.0-alpha2-421-ge96352b0e --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index b91ec55e4..a1f129287 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v3.0.0-alpha2-415-g09c2d7f3f +v3.0.0-alpha2-421-ge96352b0e From e855083964794f44aa0793597c945091b5b34ff2 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 8 Feb 2021 23:21:53 +0100 Subject: [PATCH 07/10] systematic names --- src/CPFEM.f90 | 6 +- src/CPFEM2.f90 | 10 +- src/commercialFEM_fileList.f90 | 34 +- src/grid/DAMASK_grid.f90 | 66 +-- src/grid/grid_mech_FEM.f90 | 135 +++--- src/grid/grid_mech_spectral_basic.f90 | 40 +- src/grid/grid_mech_spectral_polarisation.f90 | 40 +- src/homogenization.f90 | 32 +- src/homogenization_damage.f90 | 4 +- ...nics.f90 => homogenization_mechanical.f90} | 94 ++--- ....f90 => homogenization_mechanical_RGC.f90} | 24 +- ...> homogenization_mechanical_isostrain.f90} | 14 +- ...f90 => homogenization_mechanical_pass.f90} | 8 +- src/homogenization_thermal.f90 | 6 +- src/mesh/DAMASK_mesh.f90 | 8 +- src/mesh/FEM_utilities.f90 | 12 +- src/mesh/mesh_mech_FEM.f90 | 110 ++--- src/phase.f90 | 180 ++++---- src/phase_damage.f90 | 34 +- src/phase_damage_anisobrittle.f90 | 2 +- src/phase_damage_anisoductile.f90 | 2 +- src/phase_damage_isobrittle.f90 | 2 +- src/phase_damage_isoductile.f90 | 2 +- ...ase_mechanics.f90 => phase_mechanical.f90} | 384 +++++++++--------- ...rmation.f90 => phase_mechanical_eigen.f90} | 14 +- ...hase_mechanical_eigen_cleavageopening.f90} | 0 ...ase_mechanical_eigen_slipplaneopening.f90} | 0 ...ase_mechanical_eigen_thermalexpansion.f90} | 0 ...astic.f90 => phase_mechanical_plastic.f90} | 64 +-- ...hase_mechanical_plastic_dislotungsten.f90} | 16 +- ...=> phase_mechanical_plastic_dislotwin.f90} | 36 +- ...=> phase_mechanical_plastic_isotropic.f90} | 8 +- ...hase_mechanical_plastic_kinehardening.f90} | 16 +- ....f90 => phase_mechanical_plastic_none.f90} | 2 +- ... => phase_mechanical_plastic_nonlocal.f90} | 64 +-- ...hase_mechanical_plastic_phenopowerlaw.f90} | 22 +- src/phase_thermal.f90 | 14 +- src/phase_thermal_dissipation.f90 | 4 +- src/phase_thermal_externalheat.f90 | 2 +- 39 files changed, 769 insertions(+), 742 deletions(-) rename src/{homogenization_mechanics.f90 => homogenization_mechanical.f90} (79%) rename src/{homogenization_mechanics_RGC.f90 => homogenization_mechanical_RGC.f90} (98%) rename src/{homogenization_mechanics_isostrain.f90 => homogenization_mechanical_isostrain.f90} (91%) rename src/{homogenization_mechanics_none.f90 => homogenization_mechanical_pass.f90} (87%) rename src/{phase_mechanics.f90 => phase_mechanical.f90} (82%) rename src/{phase_mechanics_eigendeformation.f90 => phase_mechanical_eigen.f90} (96%) rename src/{phase_mechanics_eigendeformation_cleavageopening.f90 => phase_mechanical_eigen_cleavageopening.f90} (100%) rename src/{phase_mechanics_eigendeformation_slipplaneopening.f90 => phase_mechanical_eigen_slipplaneopening.f90} (100%) rename src/{phase_mechanics_eigendeformation_thermalexpansion.f90 => phase_mechanical_eigen_thermalexpansion.f90} (100%) rename src/{phase_mechanics_plastic.f90 => phase_mechanical_plastic.f90} (91%) rename src/{phase_mechanics_plastic_dislotungsten.f90 => phase_mechanical_plastic_dislotungsten.f90} (97%) rename src/{phase_mechanics_plastic_dislotwin.f90 => phase_mechanical_plastic_dislotwin.f90} (97%) rename src/{phase_mechanics_plastic_isotropic.f90 => phase_mechanical_plastic_isotropic.f90} (97%) rename src/{phase_mechanics_plastic_kinehardening.f90 => phase_mechanical_plastic_kinehardening.f90} (97%) rename src/{phase_mechanics_plastic_none.f90 => phase_mechanical_plastic_none.f90} (96%) rename src/{phase_mechanics_plastic_nonlocal.f90 => phase_mechanical_plastic_nonlocal.f90} (97%) rename src/{phase_mechanics_plastic_phenopowerlaw.f90 => phase_mechanical_plastic_phenopowerlaw.f90} (96%) diff --git a/src/CPFEM.f90 b/src/CPFEM.f90 index e441339cb..f0e510433 100644 --- a/src/CPFEM.f90 +++ b/src/CPFEM.f90 @@ -85,7 +85,7 @@ subroutine CPFEM_initAll call discretization_marc_init call lattice_init call material_init(.false.) - call constitutive_init + call phase_init call homogenization_init call crystallite_init call CPFEM_init @@ -257,7 +257,7 @@ end subroutine CPFEM_general subroutine CPFEM_forward call homogenization_forward - call constitutive_forward + call phase_forward end subroutine CPFEM_forward @@ -272,7 +272,7 @@ subroutine CPFEM_results(inc,time) call results_openJobFile call results_addIncrement(inc,time) - call constitutive_results + call phase_results call homogenization_results call discretization_results call results_finalizeIncrement diff --git a/src/CPFEM2.f90 b/src/CPFEM2.f90 index 2b32a0cbb..e57de7f67 100644 --- a/src/CPFEM2.f90 +++ b/src/CPFEM2.f90 @@ -60,7 +60,7 @@ subroutine CPFEM_initAll call discretization_grid_init(restart=interface_restartInc>0) #endif call material_init(restart=interface_restartInc>0) - call constitutive_init + call phase_init call homogenization_init call crystallite_init call CPFEM_init @@ -87,7 +87,7 @@ subroutine CPFEM_init fileHandle = HDF5_openFile(fileName) call homogenization_restartRead(fileHandle) - call constitutive_restartRead(fileHandle) + call phase_restartRead(fileHandle) call HDF5_closeFile(fileHandle) endif @@ -110,7 +110,7 @@ subroutine CPFEM_restartWrite fileHandle = HDF5_openFile(fileName,'a') call homogenization_restartWrite(fileHandle) - call constitutive_restartWrite(fileHandle) + call phase_restartWrite(fileHandle) call HDF5_closeFile(fileHandle) @@ -123,7 +123,7 @@ end subroutine CPFEM_restartWrite subroutine CPFEM_forward call homogenization_forward - call constitutive_forward + call phase_forward end subroutine CPFEM_forward @@ -138,7 +138,7 @@ subroutine CPFEM_results(inc,time) call results_openJobFile call results_addIncrement(inc,time) - call constitutive_results + call phase_results call homogenization_results call discretization_results call results_finalizeIncrement diff --git a/src/commercialFEM_fileList.f90 b/src/commercialFEM_fileList.f90 index 6480f7382..a191298a3 100644 --- a/src/commercialFEM_fileList.f90 +++ b/src/commercialFEM_fileList.f90 @@ -20,19 +20,19 @@ #include "material.f90" #include "lattice.f90" #include "phase.f90" -#include "phase_mechanics.f90" -#include "phase_mechanics_plastic.f90" -#include "phase_mechanics_plastic_none.f90" -#include "phase_mechanics_plastic_isotropic.f90" -#include "phase_mechanics_plastic_phenopowerlaw.f90" -#include "phase_mechanics_plastic_kinehardening.f90" -#include "phase_mechanics_plastic_dislotwin.f90" -#include "phase_mechanics_plastic_dislotungsten.f90" -#include "phase_mechanics_plastic_nonlocal.f90" -#include "phase_mechanics_eigendeformation.f90" -#include "phase_mechanics_eigendeformation_cleavageopening.f90" -#include "phase_mechanics_eigendeformation_slipplaneopening.f90" -#include "phase_mechanics_eigendeformation_thermalexpansion.f90" +#include "phase_mechanical.f90" +#include "phase_mechanical_plastic.f90" +#include "phase_mechanical_plastic_none.f90" +#include "phase_mechanical_plastic_isotropic.f90" +#include "phase_mechanical_plastic_phenopowerlaw.f90" +#include "phase_mechanical_plastic_kinehardening.f90" +#include "phase_mechanical_plastic_dislotwin.f90" +#include "phase_mechanical_plastic_dislotungsten.f90" +#include "phase_mechanical_plastic_nonlocal.f90" +#include "phase_mechanical_eigen.f90" +#include "phase_mechanical_eigen_cleavageopening.f90" +#include "phase_mechanical_eigen_slipplaneopening.f90" +#include "phase_mechanical_eigen_thermalexpansion.f90" #include "phase_thermal.f90" #include "phase_thermal_dissipation.f90" #include "phase_thermal_externalheat.f90" @@ -44,10 +44,10 @@ #include "damage_none.f90" #include "damage_nonlocal.f90" #include "homogenization.f90" -#include "homogenization_mechanics.f90" -#include "homogenization_mechanics_none.f90" -#include "homogenization_mechanics_isostrain.f90" -#include "homogenization_mechanics_RGC.f90" +#include "homogenization_mechanical.f90" +#include "homogenization_mechanical_pass.f90" +#include "homogenization_mechanical_isostrain.f90" +#include "homogenization_mechanical_RGC.f90" #include "homogenization_thermal.f90" #include "homogenization_damage.f90" #include "CPFEM.f90" diff --git a/src/grid/DAMASK_grid.f90 b/src/grid/DAMASK_grid.f90 index 7e649ff38..ed1ada171 100644 --- a/src/grid/DAMASK_grid.f90 +++ b/src/grid/DAMASK_grid.f90 @@ -18,9 +18,9 @@ program DAMASK_grid use CPFEM2 use material use spectral_utilities - use grid_mech_spectral_basic - use grid_mech_spectral_polarisation - use grid_mech_FEM + use grid_mechanical_spectral_basic + use grid_mechanical_spectral_polarisation + use grid_mechanical_FEM use grid_damage_spectral use grid_thermal_spectral use results @@ -83,16 +83,16 @@ program DAMASK_grid type(tLoadCase), allocatable, dimension(:) :: loadCases !< array of all load cases type(tSolutionState), allocatable, dimension(:) :: solres - procedure(grid_mech_spectral_basic_init), pointer :: & - mech_init - procedure(grid_mech_spectral_basic_forward), pointer :: & - mech_forward - procedure(grid_mech_spectral_basic_solution), pointer :: & - mech_solution - procedure(grid_mech_spectral_basic_updateCoords), pointer :: & - mech_updateCoords - procedure(grid_mech_spectral_basic_restartWrite), pointer :: & - mech_restartWrite + procedure(grid_mechanical_spectral_basic_init), pointer :: & + mechanical_init + procedure(grid_mechanical_spectral_basic_forward), pointer :: & + mechanical_forward + procedure(grid_mechanical_spectral_basic_solution), pointer :: & + mechanical_solution + procedure(grid_mechanical_spectral_basic_updateCoords), pointer :: & + mechanical_updateCoords + procedure(grid_mechanical_spectral_basic_restartWrite), pointer :: & + mechanical_restartWrite external :: & quit @@ -139,25 +139,25 @@ program DAMASK_grid debug_grid => config_debug%get('grid',defaultVal=emptyList) select case (trim(num_grid%get_asString('solver', defaultVal = 'Basic'))) case ('Basic') - mech_init => grid_mech_spectral_basic_init - mech_forward => grid_mech_spectral_basic_forward - mech_solution => grid_mech_spectral_basic_solution - mech_updateCoords => grid_mech_spectral_basic_updateCoords - mech_restartWrite => grid_mech_spectral_basic_restartWrite + mechanical_init => grid_mechanical_spectral_basic_init + mechanical_forward => grid_mechanical_spectral_basic_forward + mechanical_solution => grid_mechanical_spectral_basic_solution + mechanical_updateCoords => grid_mechanical_spectral_basic_updateCoords + mechanical_restartWrite => grid_mechanical_spectral_basic_restartWrite case ('Polarisation') - mech_init => grid_mech_spectral_polarisation_init - mech_forward => grid_mech_spectral_polarisation_forward - mech_solution => grid_mech_spectral_polarisation_solution - mech_updateCoords => grid_mech_spectral_polarisation_updateCoords - mech_restartWrite => grid_mech_spectral_polarisation_restartWrite + mechanical_init => grid_mechanical_spectral_polarisation_init + mechanical_forward => grid_mechanical_spectral_polarisation_forward + mechanical_solution => grid_mechanical_spectral_polarisation_solution + mechanical_updateCoords => grid_mechanical_spectral_polarisation_updateCoords + mechanical_restartWrite => grid_mechanical_spectral_polarisation_restartWrite case ('FEM') - mech_init => grid_mech_FEM_init - mech_forward => grid_mech_FEM_forward - mech_solution => grid_mech_FEM_solution - mech_updateCoords => grid_mech_FEM_updateCoords - mech_restartWrite => grid_mech_FEM_restartWrite + mechanical_init => grid_mechanical_FEM_init + mechanical_forward => grid_mechanical_FEM_forward + mechanical_solution => grid_mechanical_FEM_solution + mechanical_updateCoords => grid_mechanical_FEM_updateCoords + mechanical_restartWrite => grid_mechanical_FEM_restartWrite case default call IO_error(error_ID = 891, ext_msg = trim(num_grid%get_asString('solver'))) @@ -303,7 +303,7 @@ program DAMASK_grid do field = 1, nActiveFields select case (loadCases(1)%ID(field)) case(FIELD_MECH_ID) - call mech_init + call mechanical_init case(FIELD_THERMAL_ID) call grid_thermal_spectral_init @@ -379,7 +379,7 @@ program DAMASK_grid do field = 1, nActiveFields select case(loadCases(l)%ID(field)) case(FIELD_MECH_ID) - call mech_forward (& + call mechanical_forward (& cutBack,guess,timeinc,timeIncOld,remainingLoadCaseTime, & deformation_BC = loadCases(l)%deformation, & stress_BC = loadCases(l)%stress, & @@ -399,7 +399,7 @@ program DAMASK_grid do field = 1, nActiveFields select case(loadCases(l)%ID(field)) case(FIELD_MECH_ID) - solres(field) = mech_solution(incInfo) + solres(field) = mechanical_solution(incInfo) case(FIELD_THERMAL_ID) solres(field) = grid_thermal_spectral_solution(timeinc) case(FIELD_DAMAGE_ID) @@ -420,7 +420,7 @@ program DAMASK_grid if ( (all(solres(:)%converged .and. solres(:)%stagConverged)) & ! converged .and. .not. solres(1)%termIll) then ! and acceptable solution found - call mech_updateCoords + call mechanical_updateCoords timeIncOld = timeinc cutBack = .false. guess = .true. ! start guessing after first converged (sub)inc @@ -463,7 +463,7 @@ program DAMASK_grid call MPI_Allreduce(interface_SIGUSR2,signal,1,MPI_LOGICAL,MPI_LOR,PETSC_COMM_WORLD,ierr) if (ierr /= 0) error stop 'MPI error' if (mod(inc,loadCases(l)%f_restart) == 0 .or. signal) then - call mech_restartWrite + call mechanical_restartWrite call CPFEM_restartWrite endif if(signal) call interface_setSIGUSR2(.false.) diff --git a/src/grid/grid_mech_FEM.f90 b/src/grid/grid_mech_FEM.f90 index 003f568c6..806973a4d 100644 --- a/src/grid/grid_mech_FEM.f90 +++ b/src/grid/grid_mech_FEM.f90 @@ -4,7 +4,7 @@ !> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH !> @brief Grid solver for mechanics: FEM !-------------------------------------------------------------------------------------------------- -module grid_mech_FEM +module grid_mechanical_FEM #include #include use PETScdmda @@ -45,8 +45,8 @@ module grid_mech_FEM !-------------------------------------------------------------------------------------------------- ! PETSc data - DM :: mech_grid - SNES :: mech_snes + DM :: mechanical_grid + SNES :: mechanical_snes Vec :: solution_current, solution_lastInc, solution_rate !-------------------------------------------------------------------------------------------------- @@ -79,18 +79,18 @@ module grid_mech_FEM totalIter = 0 !< total iteration in current increment public :: & - grid_mech_FEM_init, & - grid_mech_FEM_solution, & - grid_mech_FEM_forward, & - grid_mech_FEM_updateCoords, & - grid_mech_FEM_restartWrite + grid_mechanical_FEM_init, & + grid_mechanical_FEM_solution, & + grid_mechanical_FEM_forward, & + grid_mechanical_FEM_updateCoords, & + grid_mechanical_FEM_restartWrite contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all necessary fields and fills them with data, potentially from restart info !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_FEM_init +subroutine grid_mechanical_FEM_init real(pReal), parameter :: HGCoeff = 0.0e-2_pReal real(pReal), parameter, dimension(4,8) :: & @@ -114,7 +114,7 @@ subroutine grid_mech_FEM_init num_grid, & debug_grid - print'(/,a)', ' <<<+- grid_mech_FEM init -+>>>'; flush(IO_STDOUT) + print'(/,a)', ' <<<+- grid_mechanical_FEM init -+>>>'; flush(IO_STDOUT) !------------------------------------------------------------------------------------------------- ! debugging options @@ -141,8 +141,11 @@ subroutine grid_mech_FEM_init !-------------------------------------------------------------------------------------------------- ! set default and user defined options for PETSc - call PetscOptionsInsertString(PETSC_NULL_OPTIONS,'-mech_snes_type newtonls -mech_ksp_type fgmres & - &-mech_ksp_max_it 25 -mech_pc_type ml -mech_mg_levels_ksp_type chebyshev',ierr) + call PetscOptionsInsertString(PETSC_NULL_OPTIONS, & + '-mechanical_snes_type newtonls -mechanical_ksp_type fgmres & + &-mechanical_ksp_max_it 25 -mechanical_pc_type ml & + &-mechanical_mg_levels_ksp_type chebyshev', & + ierr) CHKERRQ(ierr) call PetscOptionsInsertString(PETSC_NULL_OPTIONS,num_grid%get_asString('petsc_options',defaultVal=''),ierr) CHKERRQ(ierr) @@ -155,8 +158,10 @@ subroutine grid_mech_FEM_init !-------------------------------------------------------------------------------------------------- ! initialize solver specific parts of PETSc - call SNESCreate(PETSC_COMM_WORLD,mech_snes,ierr); CHKERRQ(ierr) - call SNESSetOptionsPrefix(mech_snes,'mech_',ierr);CHKERRQ(ierr) + call SNESCreate(PETSC_COMM_WORLD,mechanical_snes,ierr) + CHKERRQ(ierr) + call SNESSetOptionsPrefix(mechanical_snes,'mechanical_',ierr) + CHKERRQ(ierr) localK = 0 localK(worldrank) = grid3 call MPI_Allreduce(MPI_IN_PLACE,localK,worldsize,MPI_INTEGER,MPI_SUM,PETSC_COMM_WORLD,ierr) @@ -167,34 +172,44 @@ subroutine grid_mech_FEM_init 1, 1, worldsize, & 3, 1, & [grid(1)],[grid(2)],localK, & - mech_grid,ierr) + mechanical_grid,ierr) CHKERRQ(ierr) - call SNESSetDM(mech_snes,mech_grid,ierr); CHKERRQ(ierr) - call DMsetFromOptions(mech_grid,ierr); CHKERRQ(ierr) - call DMsetUp(mech_grid,ierr); CHKERRQ(ierr) - call DMDASetUniformCoordinates(mech_grid,0.0_pReal,geomSize(1),0.0_pReal,geomSize(2),0.0_pReal,geomSize(3),ierr) + call SNESSetDM(mechanical_snes,mechanical_grid,ierr) CHKERRQ(ierr) - call DMCreateGlobalVector(mech_grid,solution_current,ierr); CHKERRQ(ierr) - call DMCreateGlobalVector(mech_grid,solution_lastInc,ierr); CHKERRQ(ierr) - call DMCreateGlobalVector(mech_grid,solution_rate ,ierr); CHKERRQ(ierr) - call DMSNESSetFunctionLocal(mech_grid,formResidual,PETSC_NULL_SNES,ierr) + call DMsetFromOptions(mechanical_grid,ierr) CHKERRQ(ierr) - call DMSNESSetJacobianLocal(mech_grid,formJacobian,PETSC_NULL_SNES,ierr) + call DMsetUp(mechanical_grid,ierr) CHKERRQ(ierr) - call SNESSetConvergenceTest(mech_snes,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,ierr) ! specify custom convergence check function "_converged" + call DMDASetUniformCoordinates(mechanical_grid,0.0_pReal,geomSize(1),0.0_pReal,geomSize(2),0.0_pReal,geomSize(3),ierr) + CHKERRQ(ierr) + call DMCreateGlobalVector(mechanical_grid,solution_current,ierr) + CHKERRQ(ierr) + call DMCreateGlobalVector(mechanical_grid,solution_lastInc,ierr) + CHKERRQ(ierr) + call DMCreateGlobalVector(mechanical_grid,solution_rate ,ierr) + CHKERRQ(ierr) + call DMSNESSetFunctionLocal(mechanical_grid,formResidual,PETSC_NULL_SNES,ierr) + CHKERRQ(ierr) + call DMSNESSetJacobianLocal(mechanical_grid,formJacobian,PETSC_NULL_SNES,ierr) + CHKERRQ(ierr) + call SNESSetConvergenceTest(mechanical_snes,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,ierr) ! specify custom convergence check function "_converged" + CHKERRQ(ierr) + call SNESSetMaxLinearSolveFailures(mechanical_snes, huge(1), ierr) ! ignore linear solve failures + CHKERRQ(ierr) + call SNESSetFromOptions(mechanical_snes,ierr) ! pull it all together with additional cli arguments CHKERRQ(ierr) - call SNESSetMaxLinearSolveFailures(mech_snes, huge(1), ierr); CHKERRQ(ierr) ! ignore linear solve failures - call SNESSetFromOptions(mech_snes,ierr); CHKERRQ(ierr) ! pull it all together with additional cli arguments !-------------------------------------------------------------------------------------------------- ! init fields call VecSet(solution_current,0.0_pReal,ierr);CHKERRQ(ierr) call VecSet(solution_lastInc,0.0_pReal,ierr);CHKERRQ(ierr) call VecSet(solution_rate ,0.0_pReal,ierr);CHKERRQ(ierr) - call DMDAVecGetArrayF90(mech_grid,solution_current,u_current,ierr); CHKERRQ(ierr) - call DMDAVecGetArrayF90(mech_grid,solution_lastInc,u_lastInc,ierr); CHKERRQ(ierr) + call DMDAVecGetArrayF90(mechanical_grid,solution_current,u_current,ierr) + CHKERRQ(ierr) + call DMDAVecGetArrayF90(mechanical_grid,solution_lastInc,u_lastInc,ierr) + CHKERRQ(ierr) - call DMDAGetCorners(mech_grid,xstart,ystart,zstart,xend,yend,zend,ierr) ! local grid extent + call DMDAGetCorners(mechanical_grid,xstart,ystart,zstart,xend,yend,zend,ierr) ! local grid extent CHKERRQ(ierr) xend = xstart+xend-1 yend = ystart+yend-1 @@ -242,9 +257,9 @@ subroutine grid_mech_FEM_init call utilities_constitutiveResponse(P_current,P_av,C_volAvg,devNull, & ! stress field, stress avg, global average of stiffness and (min+max)/2 F, & ! target F 0.0_pReal) ! time increment - call DMDAVecRestoreArrayF90(mech_grid,solution_current,u_current,ierr) + call DMDAVecRestoreArrayF90(mechanical_grid,solution_current,u_current,ierr) CHKERRQ(ierr) - call DMDAVecRestoreArrayF90(mech_grid,solution_lastInc,u_lastInc,ierr) + call DMDAVecRestoreArrayF90(mechanical_grid,solution_lastInc,u_lastInc,ierr) CHKERRQ(ierr) restartRead2: if (interface_restartInc > 0) then @@ -257,13 +272,13 @@ subroutine grid_mech_FEM_init endif restartRead2 -end subroutine grid_mech_FEM_init +end subroutine grid_mechanical_FEM_init !-------------------------------------------------------------------------------------------------- !> @brief solution for the FEM scheme with internal iterations !-------------------------------------------------------------------------------------------------- -function grid_mech_FEM_solution(incInfoIn) result(solution) +function grid_mechanical_FEM_solution(incInfoIn) result(solution) !-------------------------------------------------------------------------------------------------- ! input data for solution @@ -284,11 +299,13 @@ function grid_mech_FEM_solution(incInfoIn) result(solution) !-------------------------------------------------------------------------------------------------- ! solve BVP - call SNESsolve(mech_snes,PETSC_NULL_VEC,solution_current,ierr); CHKERRQ(ierr) + call SNESsolve(mechanical_snes,PETSC_NULL_VEC,solution_current,ierr) + CHKERRQ(ierr) !-------------------------------------------------------------------------------------------------- ! check convergence - call SNESGetConvergedReason(mech_snes,reason,ierr); CHKERRQ(ierr) + call SNESGetConvergedReason(mechanical_snes,reason,ierr) + CHKERRQ(ierr) solution%converged = reason > 0 solution%iterationsNeeded = totalIter @@ -296,14 +313,14 @@ function grid_mech_FEM_solution(incInfoIn) result(solution) terminallyIll = .false. P_aim = merge(P_aim,P_av,params%stress_mask) -end function grid_mech_FEM_solution +end function grid_mechanical_FEM_solution !-------------------------------------------------------------------------------------------------- !> @brief forwarding routine !> @details find new boundary conditions and best F estimate for end of current timestep !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_FEM_forward(cutBack,guess,Delta_t,Delta_t_old,t_remaining,& +subroutine grid_mechanical_FEM_forward(cutBack,guess,Delta_t,Delta_t_old,t_remaining,& deformation_BC,stress_BC,rotation_BC) logical, intent(in) :: & @@ -323,8 +340,10 @@ subroutine grid_mech_FEM_forward(cutBack,guess,Delta_t,Delta_t_old,t_remaining,& u_current,u_lastInc - call DMDAVecGetArrayF90(mech_grid,solution_current,u_current,ierr); CHKERRQ(ierr) - call DMDAVecGetArrayF90(mech_grid,solution_lastInc,u_lastInc,ierr); CHKERRQ(ierr) + call DMDAVecGetArrayF90(mechanical_grid,solution_current,u_current,ierr) + CHKERRQ(ierr) + call DMDAVecGetArrayF90(mechanical_grid,solution_lastInc,u_lastInc,ierr) + CHKERRQ(ierr) if (cutBack) then C_volAvg = C_volAvgLastInc @@ -371,8 +390,10 @@ subroutine grid_mech_FEM_forward(cutBack,guess,Delta_t,Delta_t_old,t_remaining,& call VecAXPY(solution_current,Delta_t,solution_rate,ierr); CHKERRQ(ierr) - call DMDAVecRestoreArrayF90(mech_grid,solution_current,u_current,ierr); CHKERRQ(ierr) - call DMDAVecRestoreArrayF90(mech_grid,solution_lastInc,u_lastInc,ierr); CHKERRQ(ierr) + call DMDAVecRestoreArrayF90(mechanical_grid,solution_current,u_current,ierr) + CHKERRQ(ierr) + call DMDAVecRestoreArrayF90(mechanical_grid,solution_lastInc,u_lastInc,ierr) + CHKERRQ(ierr) !-------------------------------------------------------------------------------------------------- ! set module wide available data @@ -380,31 +401,33 @@ subroutine grid_mech_FEM_forward(cutBack,guess,Delta_t,Delta_t_old,t_remaining,& params%rotation_BC = rotation_BC params%timeinc = Delta_t -end subroutine grid_mech_FEM_forward +end subroutine grid_mechanical_FEM_forward !-------------------------------------------------------------------------------------------------- !> @brief Update coordinates !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_FEM_updateCoords +subroutine grid_mechanical_FEM_updateCoords call utilities_updateCoords(F) -end subroutine grid_mech_FEM_updateCoords +end subroutine grid_mechanical_FEM_updateCoords !-------------------------------------------------------------------------------------------------- !> @brief Write current solver and constitutive data for restart to file !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_FEM_restartWrite +subroutine grid_mechanical_FEM_restartWrite PetscErrorCode :: ierr integer(HID_T) :: fileHandle, groupHandle PetscScalar, dimension(:,:,:,:), pointer :: u_current,u_lastInc character(len=pStringLen) :: fileName - call DMDAVecGetArrayF90(mech_grid,solution_current,u_current,ierr); CHKERRQ(ierr) - call DMDAVecGetArrayF90(mech_grid,solution_lastInc,u_lastInc,ierr); CHKERRQ(ierr) + call DMDAVecGetArrayF90(mechanical_grid,solution_current,u_current,ierr) + CHKERRQ(ierr) + call DMDAVecGetArrayF90(mechanical_grid,solution_lastInc,u_lastInc,ierr) + CHKERRQ(ierr) print*, 'writing solver data required for restart to file'; flush(IO_STDOUT) @@ -427,10 +450,12 @@ subroutine grid_mech_FEM_restartWrite call HDF5_closeGroup(groupHandle) call HDF5_closeFile(fileHandle) - call DMDAVecRestoreArrayF90(mech_grid,solution_current,u_current,ierr);CHKERRQ(ierr) - call DMDAVecRestoreArrayF90(mech_grid,solution_lastInc,u_lastInc,ierr);CHKERRQ(ierr) + call DMDAVecRestoreArrayF90(mechanical_grid,solution_current,u_current,ierr) + CHKERRQ(ierr) + call DMDAVecRestoreArrayF90(mechanical_grid,solution_lastInc,u_lastInc,ierr) + CHKERRQ(ierr) -end subroutine grid_mech_FEM_restartWrite +end subroutine grid_mechanical_FEM_restartWrite !-------------------------------------------------------------------------------------------------- @@ -498,8 +523,10 @@ subroutine formResidual(da_local,x_local, & PetscErrorCode :: ierr real(pReal), dimension(3,3,3,3) :: devNull - call SNESGetNumberFunctionEvals(mech_snes,nfuncs,ierr); CHKERRQ(ierr) - call SNESGetIterationNumber(mech_snes,PETScIter,ierr); CHKERRQ(ierr) + call SNESGetNumberFunctionEvals(mechanical_snes,nfuncs,ierr) + CHKERRQ(ierr) + call SNESGetIterationNumber(mechanical_snes,PETScIter,ierr) + CHKERRQ(ierr) if (nfuncs == 0 .and. PETScIter == 0) totalIter = -1 ! new increment @@ -679,4 +706,4 @@ subroutine formJacobian(da_local,x_local,Jac_pre,Jac,dummy,ierr) end subroutine formJacobian -end module grid_mech_FEM +end module grid_mechanical_FEM diff --git a/src/grid/grid_mech_spectral_basic.f90 b/src/grid/grid_mech_spectral_basic.f90 index 9bc36165f..f3f30c0af 100644 --- a/src/grid/grid_mech_spectral_basic.f90 +++ b/src/grid/grid_mech_spectral_basic.f90 @@ -4,7 +4,7 @@ !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH !> @brief Grid solver for mechanics: Spectral basic !-------------------------------------------------------------------------------------------------- -module grid_mech_spectral_basic +module grid_mechanical_spectral_basic #include #include use PETScdmda @@ -79,18 +79,18 @@ module grid_mech_spectral_basic totalIter = 0 !< total iteration in current increment public :: & - grid_mech_spectral_basic_init, & - grid_mech_spectral_basic_solution, & - grid_mech_spectral_basic_forward, & - grid_mech_spectral_basic_updateCoords, & - grid_mech_spectral_basic_restartWrite + grid_mechanical_spectral_basic_init, & + grid_mechanical_spectral_basic_solution, & + grid_mechanical_spectral_basic_forward, & + grid_mechanical_spectral_basic_updateCoords, & + grid_mechanical_spectral_basic_restartWrite contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all necessary fields and fills them with data, potentially from restart info !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_spectral_basic_init +subroutine grid_mechanical_spectral_basic_init real(pReal), dimension(3,3,grid(1),grid(2),grid3) :: P PetscErrorCode :: ierr @@ -105,7 +105,7 @@ subroutine grid_mech_spectral_basic_init num_grid, & debug_grid - print'(/,a)', ' <<<+- grid_mech_spectral_basic init -+>>>'; flush(IO_STDOUT) + print'(/,a)', ' <<<+- grid_mechanical_spectral_basic init -+>>>'; flush(IO_STDOUT) print*, 'Eisenlohr et al., International Journal of Plasticity 46:37–53, 2013' print*, 'https://doi.org/10.1016/j.ijplas.2012.09.012'//IO_EOL @@ -139,7 +139,7 @@ subroutine grid_mech_spectral_basic_init !-------------------------------------------------------------------------------------------------- ! set default and user defined options for PETSc - call PetscOptionsInsertString(PETSC_NULL_OPTIONS,'-mech_snes_type ngmres',ierr) + call PetscOptionsInsertString(PETSC_NULL_OPTIONS,'-mechanical_snes_type ngmres',ierr) CHKERRQ(ierr) call PetscOptionsInsertString(PETSC_NULL_OPTIONS,num_grid%get_asString('petsc_options',defaultVal=''),ierr) CHKERRQ(ierr) @@ -152,7 +152,7 @@ subroutine grid_mech_spectral_basic_init !-------------------------------------------------------------------------------------------------- ! initialize solver specific parts of PETSc call SNESCreate(PETSC_COMM_WORLD,snes,ierr); CHKERRQ(ierr) - call SNESSetOptionsPrefix(snes,'mech_',ierr);CHKERRQ(ierr) + call SNESSetOptionsPrefix(snes,'mechanical_',ierr);CHKERRQ(ierr) localK = 0 localK(worldrank) = grid3 call MPI_Allreduce(MPI_IN_PLACE,localK,worldsize,MPI_INTEGER,MPI_SUM,PETSC_COMM_WORLD,ierr) @@ -222,13 +222,13 @@ subroutine grid_mech_spectral_basic_init call utilities_updateGamma(C_minMaxAvg) call utilities_saveReferenceStiffness -end subroutine grid_mech_spectral_basic_init +end subroutine grid_mechanical_spectral_basic_init !-------------------------------------------------------------------------------------------------- !> @brief solution for the basic scheme with internal iterations !-------------------------------------------------------------------------------------------------- -function grid_mech_spectral_basic_solution(incInfoIn) result(solution) +function grid_mechanical_spectral_basic_solution(incInfoIn) result(solution) !-------------------------------------------------------------------------------------------------- ! input data for solution @@ -262,14 +262,14 @@ function grid_mech_spectral_basic_solution(incInfoIn) result(solution) terminallyIll = .false. P_aim = merge(P_aim,P_av,params%stress_mask) -end function grid_mech_spectral_basic_solution +end function grid_mechanical_spectral_basic_solution !-------------------------------------------------------------------------------------------------- !> @brief forwarding routine !> @details find new boundary conditions and best F estimate for end of current timestep !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_spectral_basic_forward(cutBack,guess,Delta_t,Delta_t_old,t_remaining,& +subroutine grid_mechanical_spectral_basic_forward(cutBack,guess,Delta_t,Delta_t_old,t_remaining,& deformation_BC,stress_BC,rotation_BC) logical, intent(in) :: & @@ -339,13 +339,13 @@ subroutine grid_mech_spectral_basic_forward(cutBack,guess,Delta_t,Delta_t_old,t_ params%rotation_BC = rotation_BC params%timeinc = Delta_t -end subroutine grid_mech_spectral_basic_forward +end subroutine grid_mechanical_spectral_basic_forward !-------------------------------------------------------------------------------------------------- !> @brief Update coordinates !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_spectral_basic_updateCoords +subroutine grid_mechanical_spectral_basic_updateCoords PetscErrorCode :: ierr PetscScalar, dimension(:,:,:,:), pointer :: F @@ -354,13 +354,13 @@ subroutine grid_mech_spectral_basic_updateCoords call utilities_updateCoords(F) call DMDAVecRestoreArrayF90(da,solution_vec,F,ierr); CHKERRQ(ierr) -end subroutine grid_mech_spectral_basic_updateCoords +end subroutine grid_mechanical_spectral_basic_updateCoords !-------------------------------------------------------------------------------------------------- !> @brief Write current solver and constitutive data for restart to file !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_spectral_basic_restartWrite +subroutine grid_mechanical_spectral_basic_restartWrite PetscErrorCode :: ierr integer(HID_T) :: fileHandle, groupHandle @@ -393,7 +393,7 @@ subroutine grid_mech_spectral_basic_restartWrite call DMDAVecRestoreArrayF90(da,solution_vec,F,ierr); CHKERRQ(ierr) -end subroutine grid_mech_spectral_basic_restartWrite +end subroutine grid_mechanical_spectral_basic_restartWrite !-------------------------------------------------------------------------------------------------- @@ -506,4 +506,4 @@ subroutine formResidual(in, F, & end subroutine formResidual -end module grid_mech_spectral_basic +end module grid_mechanical_spectral_basic diff --git a/src/grid/grid_mech_spectral_polarisation.f90 b/src/grid/grid_mech_spectral_polarisation.f90 index 7160c1adc..633fced7f 100644 --- a/src/grid/grid_mech_spectral_polarisation.f90 +++ b/src/grid/grid_mech_spectral_polarisation.f90 @@ -4,7 +4,7 @@ !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH !> @brief Grid solver for mechanics: Spectral Polarisation !-------------------------------------------------------------------------------------------------- -module grid_mech_spectral_polarisation +module grid_mechanical_spectral_polarisation #include #include use PETScdmda @@ -90,18 +90,18 @@ module grid_mech_spectral_polarisation totalIter = 0 !< total iteration in current increment public :: & - grid_mech_spectral_polarisation_init, & - grid_mech_spectral_polarisation_solution, & - grid_mech_spectral_polarisation_forward, & - grid_mech_spectral_polarisation_updateCoords, & - grid_mech_spectral_polarisation_restartWrite + grid_mechanical_spectral_polarisation_init, & + grid_mechanical_spectral_polarisation_solution, & + grid_mechanical_spectral_polarisation_forward, & + grid_mechanical_spectral_polarisation_updateCoords, & + grid_mechanical_spectral_polarisation_restartWrite contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all necessary fields and fills them with data, potentially from restart info !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_spectral_polarisation_init +subroutine grid_mechanical_spectral_polarisation_init real(pReal), dimension(3,3,grid(1),grid(2),grid3) :: P PetscErrorCode :: ierr @@ -118,7 +118,7 @@ subroutine grid_mech_spectral_polarisation_init num_grid, & debug_grid - print'(/,a)', ' <<<+- grid_mech_spectral_polarisation init -+>>>'; flush(IO_STDOUT) + print'(/,a)', ' <<<+- grid_mechanical_spectral_polarisation init -+>>>'; flush(IO_STDOUT) print*, 'Shanthraj et al., International Journal of Plasticity 66:31–45, 2015' print*, 'https://doi.org/10.1016/j.ijplas.2014.02.006' @@ -157,7 +157,7 @@ subroutine grid_mech_spectral_polarisation_init !-------------------------------------------------------------------------------------------------- ! set default and user defined options for PETSc - call PetscOptionsInsertString(PETSC_NULL_OPTIONS,'-mech_snes_type ngmres',ierr) + call PetscOptionsInsertString(PETSC_NULL_OPTIONS,'-mechanical_snes_type ngmres',ierr) CHKERRQ(ierr) call PetscOptionsInsertString(PETSC_NULL_OPTIONS,num_grid%get_asString('petsc_options',defaultVal=''),ierr) CHKERRQ(ierr) @@ -172,7 +172,7 @@ subroutine grid_mech_spectral_polarisation_init !-------------------------------------------------------------------------------------------------- ! initialize solver specific parts of PETSc call SNESCreate(PETSC_COMM_WORLD,snes,ierr); CHKERRQ(ierr) - call SNESSetOptionsPrefix(snes,'mech_',ierr);CHKERRQ(ierr) + call SNESSetOptionsPrefix(snes,'mechanical_',ierr);CHKERRQ(ierr) localK = 0 localK(worldrank) = grid3 call MPI_Allreduce(MPI_IN_PLACE,localK,worldsize,MPI_INTEGER,MPI_SUM,PETSC_COMM_WORLD,ierr) @@ -250,13 +250,13 @@ subroutine grid_mech_spectral_polarisation_init C_scale = C_minMaxAvg S_scale = math_invSym3333(C_minMaxAvg) -end subroutine grid_mech_spectral_polarisation_init +end subroutine grid_mechanical_spectral_polarisation_init !-------------------------------------------------------------------------------------------------- !> @brief solution for the Polarisation scheme with internal iterations !-------------------------------------------------------------------------------------------------- -function grid_mech_spectral_polarisation_solution(incInfoIn) result(solution) +function grid_mechanical_spectral_polarisation_solution(incInfoIn) result(solution) !-------------------------------------------------------------------------------------------------- ! input data for solution @@ -294,14 +294,14 @@ function grid_mech_spectral_polarisation_solution(incInfoIn) result(solution) terminallyIll = .false. P_aim = merge(P_aim,P_av,params%stress_mask) -end function grid_mech_spectral_polarisation_solution +end function grid_mechanical_spectral_polarisation_solution !-------------------------------------------------------------------------------------------------- !> @brief forwarding routine !> @details find new boundary conditions and best F estimate for end of current timestep !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_spectral_polarisation_forward(cutBack,guess,Delta_t,Delta_t_old,t_remaining,& +subroutine grid_mechanical_spectral_polarisation_forward(cutBack,guess,Delta_t,Delta_t_old,t_remaining,& deformation_BC,stress_BC,rotation_BC) logical, intent(in) :: & @@ -393,13 +393,13 @@ subroutine grid_mech_spectral_polarisation_forward(cutBack,guess,Delta_t,Delta_t params%rotation_BC = rotation_BC params%timeinc = Delta_t -end subroutine grid_mech_spectral_polarisation_forward +end subroutine grid_mechanical_spectral_polarisation_forward !-------------------------------------------------------------------------------------------------- !> @brief Update coordinates !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_spectral_polarisation_updateCoords +subroutine grid_mechanical_spectral_polarisation_updateCoords PetscErrorCode :: ierr PetscScalar, dimension(:,:,:,:), pointer :: FandF_tau @@ -408,13 +408,13 @@ subroutine grid_mech_spectral_polarisation_updateCoords call utilities_updateCoords(FandF_tau(0:8,:,:,:)) call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr) -end subroutine grid_mech_spectral_polarisation_updateCoords +end subroutine grid_mechanical_spectral_polarisation_updateCoords !-------------------------------------------------------------------------------------------------- !> @brief Write current solver and constitutive data for restart to file !-------------------------------------------------------------------------------------------------- -subroutine grid_mech_spectral_polarisation_restartWrite +subroutine grid_mechanical_spectral_polarisation_restartWrite PetscErrorCode :: ierr integer(HID_T) :: fileHandle, groupHandle @@ -450,7 +450,7 @@ subroutine grid_mech_spectral_polarisation_restartWrite call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr) -end subroutine grid_mech_spectral_polarisation_restartWrite +end subroutine grid_mechanical_spectral_polarisation_restartWrite !-------------------------------------------------------------------------------------------------- @@ -618,4 +618,4 @@ subroutine formResidual(in, FandF_tau, & end subroutine formResidual -end module grid_mech_spectral_polarisation +end module grid_mechanical_spectral_polarisation diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 2255e83b0..62ac52f06 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -45,10 +45,10 @@ module homogenization !-------------------------------------------------------------------------------------------------- interface - module subroutine mech_init(num_homog) + module subroutine mechanical_init(num_homog) class(tNode), pointer, intent(in) :: & num_homog !< pointer to mechanical homogenization numerics data - end subroutine mech_init + end subroutine mechanical_init module subroutine thermal_init end subroutine thermal_init @@ -56,13 +56,13 @@ module homogenization module subroutine damage_init end subroutine damage_init - module subroutine mech_partition(subF,ip,el) + module subroutine mechanical_partition(subF,ip,el) real(pReal), intent(in), dimension(3,3) :: & subF integer, intent(in) :: & ip, & !< integration point el !< element number - end subroutine mech_partition + end subroutine mechanical_partition module subroutine thermal_partition(ce) integer, intent(in) :: ce @@ -76,19 +76,19 @@ module homogenization integer, intent(in) :: ip,el end subroutine thermal_homogenize - module subroutine mech_homogenize(dt,ip,el) + module subroutine mechanical_homogenize(dt,ip,el) real(pReal), intent(in) :: dt integer, intent(in) :: & ip, & !< integration point el !< element number - end subroutine mech_homogenize + end subroutine mechanical_homogenize - module subroutine mech_results(group_base,h) + module subroutine mechanical_results(group_base,h) character(len=*), intent(in) :: group_base integer, intent(in) :: h - end subroutine mech_results + end subroutine mechanical_results - module function mech_updateState(subdt,subF,ip,el) result(doneAndHappy) + module function mechanical_updateState(subdt,subF,ip,el) result(doneAndHappy) real(pReal), intent(in) :: & subdt !< current time step real(pReal), intent(in), dimension(3,3) :: & @@ -97,7 +97,7 @@ module homogenization ip, & !< integration point el !< element number logical, dimension(2) :: doneAndHappy - end function mech_updateState + end function mechanical_updateState module function thermal_conduction_getConductivity(ip,el) result(K) @@ -216,7 +216,7 @@ subroutine homogenization_init() if (num%nMPstate < 1) call IO_error(301,ext_msg='nMPstate') - call mech_init(num_homog) + call mechanical_init(num_homog) call thermal_init() call damage_init() @@ -253,7 +253,7 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE ce = (el-1)*discretization_nIPs + ip me = material_homogenizationMemberAt2(ce) - call constitutive_restore(ce,.false.) ! wrong name (is more a forward function) + call phase_restore(ce,.false.) ! wrong name (is more a forward function) if(homogState(ho)%sizeState > 0) homogState(ho)%State(:,me) = homogState(ho)%State0(:,me) if(damageState_h(ho)%sizeState > 0) damageState_h(ho)%State(:,me) = damageState_h(ho)%State0(:,me) @@ -267,7 +267,7 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE if (.not. doneAndHappy(1)) then - call mech_partition(homogenization_F(1:3,1:3,ce),ip,el) + call mechanical_partition(homogenization_F(1:3,1:3,ce),ip,el) converged = .true. do co = 1, myNgrains converged = converged .and. crystallite_stress(dt,co,ip,el) @@ -276,7 +276,7 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE if (.not. converged) then doneAndHappy = [.true.,.false.] else - doneAndHappy = mech_updateState(dt,homogenization_F(1:3,1:3,ce),ip,el) + doneAndHappy = mechanical_updateState(dt,homogenization_F(1:3,1:3,ce),ip,el) converged = all(doneAndHappy) endif endif @@ -338,7 +338,7 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE do co = 1, homogenization_Nconstituents(ho) call crystallite_orientations(co,ip,el) enddo - call mech_homogenize(dt,ip,el) + call mechanical_homogenize(dt,ip,el) enddo IpLooping3 enddo elementLooping3 !$OMP END DO @@ -365,7 +365,7 @@ subroutine homogenization_results group_base = 'current/homogenization/'//trim(material_name_homogenization(ho)) call results_closeGroup(results_addGroup(group_base)) - call mech_results(group_base,ho) + call mechanical_results(group_base,ho) group = trim(group_base)//'/damage' call results_closeGroup(results_addGroup(group)) diff --git a/src/homogenization_damage.f90 b/src/homogenization_damage.f90 index ba021a0e0..8b73ed54a 100644 --- a/src/homogenization_damage.f90 +++ b/src/homogenization_damage.f90 @@ -74,7 +74,7 @@ module subroutine damage_partition(ce) phi = current(material_homogenizationAt2(ce))%phi(material_homogenizationMemberAt2(ce)) do co = 1, homogenization_Nconstituents(material_homogenizationAt2(ce)) - call constitutive_damage_set_phi(phi,co,ce) + call phase_damage_set_phi(phi,co,ce) enddo end subroutine damage_partition @@ -120,7 +120,7 @@ module subroutine damage_nonlocal_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, p phiDot = 0.0_pReal dPhiDot_dPhi = 0.0_pReal - call constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) + call phase_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) phiDot = phiDot/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) diff --git a/src/homogenization_mechanics.f90 b/src/homogenization_mechanical.f90 similarity index 79% rename from src/homogenization_mechanics.f90 rename to src/homogenization_mechanical.f90 index 4eae859ee..91e1d3d84 100644 --- a/src/homogenization_mechanics.f90 +++ b/src/homogenization_mechanical.f90 @@ -7,52 +7,52 @@ submodule(homogenization) mechanics interface - module subroutine mech_none_init - end subroutine mech_none_init + module subroutine mechanical_pass_init + end subroutine mechanical_pass_init - module subroutine mech_isostrain_init - end subroutine mech_isostrain_init + module subroutine mechanical_isostrain_init + end subroutine mechanical_isostrain_init - module subroutine mech_RGC_init(num_homogMech) + module subroutine mechanical_RGC_init(num_homogMech) class(tNode), pointer, intent(in) :: & num_homogMech !< pointer to mechanical homogenization numerics data - end subroutine mech_RGC_init + end subroutine mechanical_RGC_init - module subroutine mech_isostrain_partitionDeformation(F,avgF) + module subroutine mechanical_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 + end subroutine mechanical_isostrain_partitionDeformation - module subroutine mech_RGC_partitionDeformation(F,avgF,instance,of) + module subroutine mechanical_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 + end subroutine mechanical_RGC_partitionDeformation - module subroutine mech_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) + module subroutine mechanical_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 real(pReal), dimension (:,:,:), intent(in) :: P !< partitioned stresses real(pReal), dimension (:,:,:,:,:), intent(in) :: dPdF !< partitioned stiffnesses integer, intent(in) :: instance - end subroutine mech_isostrain_averageStressAndItsTangent + end subroutine mechanical_isostrain_averageStressAndItsTangent - module subroutine mech_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) + module subroutine mechanical_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 + end subroutine mechanical_RGC_averageStressAndItsTangent - module function mech_RGC_updateState(P,F,avgF,dt,dPdF,ip,el) result(doneAndHappy) + module function mechanical_RGC_updateState(P,F,avgF,dt,dPdF,ip,el) result(doneAndHappy) logical, dimension(2) :: doneAndHappy real(pReal), dimension(:,:,:), intent(in) :: & P,& !< partitioned stresses @@ -63,13 +63,13 @@ submodule(homogenization) mechanics integer, intent(in) :: & ip, & !< integration point number el !< element number - end function mech_RGC_updateState + end function mechanical_RGC_updateState - module subroutine mech_RGC_results(instance,group) + module subroutine mechanical_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 subroutine mechanical_RGC_results end interface @@ -78,7 +78,7 @@ contains !-------------------------------------------------------------------------------------------------- !> @brief Allocate variables and set parameters. !-------------------------------------------------------------------------------------------------- -module subroutine mech_init(num_homog) +module subroutine mechanical_init(num_homog) class(tNode), pointer, intent(in) :: & num_homog @@ -94,17 +94,17 @@ module subroutine mech_init(num_homog) allocate(homogenization_P(3,3,discretization_nIPs*discretization_Nelems), source=0.0_pReal) num_homogMech => num_homog%get('mech',defaultVal=emptyDict) - 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 mech_RGC_init(num_homogMech) + if (any(homogenization_type == HOMOGENIZATION_NONE_ID)) call mechanical_pass_init + if (any(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID)) call mechanical_isostrain_init + if (any(homogenization_type == HOMOGENIZATION_RGC_ID)) call mechanical_RGC_init(num_homogMech) -end subroutine mech_init +end subroutine mechanical_init !-------------------------------------------------------------------------------------------------- !> @brief Partition F onto the individual constituents. !-------------------------------------------------------------------------------------------------- -module subroutine mech_partition(subF,ip,el) +module subroutine mechanical_partition(subF,ip,el) real(pReal), intent(in), dimension(3,3) :: & subF @@ -122,25 +122,25 @@ module subroutine mech_partition(subF,ip,el) Fs(1:3,1:3,1) = subF case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization - call mech_isostrain_partitionDeformation(Fs,subF) + call mechanical_isostrain_partitionDeformation(Fs,subF) case (HOMOGENIZATION_RGC_ID) chosenHomogenization - call mech_RGC_partitionDeformation(Fs,subF,ip,el) + call mechanical_RGC_partitionDeformation(Fs,subF,ip,el) end select chosenHomogenization do co = 1,homogenization_Nconstituents(material_homogenizationAt(el)) - call constitutive_mech_setF(Fs(1:3,1:3,co),co,ip,el) + call phase_mechanical_setF(Fs(1:3,1:3,co),co,ip,el) enddo -end subroutine mech_partition +end subroutine mechanical_partition !-------------------------------------------------------------------------------------------------- !> @brief Average P and dPdF from the individual constituents. !-------------------------------------------------------------------------------------------------- -module subroutine mech_homogenize(dt,ip,el) +module subroutine mechanical_homogenize(dt,ip,el) real(pReal), intent(in) :: dt integer, intent(in) :: & @@ -156,15 +156,15 @@ module subroutine mech_homogenize(dt,ip,el) chosenHomogenization: select case(homogenization_type(material_homogenizationAt(el))) case (HOMOGENIZATION_NONE_ID) chosenHomogenization - homogenization_P(1:3,1:3,ce) = constitutive_mech_getP(1,ip,el) - homogenization_dPdF(1:3,1:3,1:3,1:3,ce) = constitutive_mech_dPdF(dt,1,ip,el) + homogenization_P(1:3,1:3,ce) = phase_mechanical_getP(1,ip,el) + homogenization_dPdF(1:3,1:3,1:3,1:3,ce) = phase_mechanical_dPdF(dt,1,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization do co = 1, homogenization_Nconstituents(material_homogenizationAt(el)) - dPdFs(:,:,:,:,co) = constitutive_mech_dPdF(dt,co,ip,el) - Ps(:,:,co) = constitutive_mech_getP(co,ip,el) + dPdFs(:,:,:,:,co) = phase_mechanical_dPdF(dt,co,ip,el) + Ps(:,:,co) = phase_mechanical_getP(co,ip,el) enddo - call mech_isostrain_averageStressAndItsTangent(& + call mechanical_isostrain_averageStressAndItsTangent(& homogenization_P(1:3,1:3,ce), & homogenization_dPdF(1:3,1:3,1:3,1:3,ce),& Ps,dPdFs, & @@ -172,10 +172,10 @@ module subroutine mech_homogenize(dt,ip,el) case (HOMOGENIZATION_RGC_ID) chosenHomogenization do co = 1, homogenization_Nconstituents(material_homogenizationAt(el)) - dPdFs(:,:,:,:,co) = constitutive_mech_dPdF(dt,co,ip,el) - Ps(:,:,co) = constitutive_mech_getP(co,ip,el) + dPdFs(:,:,:,:,co) = phase_mechanical_dPdF(dt,co,ip,el) + Ps(:,:,co) = phase_mechanical_getP(co,ip,el) enddo - call mech_RGC_averageStressAndItsTangent(& + call mechanical_RGC_averageStressAndItsTangent(& homogenization_P(1:3,1:3,ce), & homogenization_dPdF(1:3,1:3,1:3,1:3,ce),& Ps,dPdFs, & @@ -183,14 +183,14 @@ module subroutine mech_homogenize(dt,ip,el) end select chosenHomogenization -end subroutine mech_homogenize +end subroutine mechanical_homogenize !-------------------------------------------------------------------------------------------------- !> @brief update the internal state of the homogenization scheme and tell whether "done" and !> "happy" with result !-------------------------------------------------------------------------------------------------- -module function mech_updateState(subdt,subF,ip,el) result(doneAndHappy) +module function mechanical_updateState(subdt,subF,ip,el) result(doneAndHappy) real(pReal), intent(in) :: & subdt !< current time step @@ -209,22 +209,22 @@ module function mech_updateState(subdt,subF,ip,el) result(doneAndHappy) if (homogenization_type(material_homogenizationAt(el)) == HOMOGENIZATION_RGC_ID) then do co = 1, homogenization_Nconstituents(material_homogenizationAt(el)) - dPdFs(:,:,:,:,co) = constitutive_mech_dPdF(subdt,co,ip,el) - Fs(:,:,co) = constitutive_mech_getF(co,ip,el) - Ps(:,:,co) = constitutive_mech_getP(co,ip,el) + dPdFs(:,:,:,:,co) = phase_mechanical_dPdF(subdt,co,ip,el) + Fs(:,:,co) = phase_mechanical_getF(co,ip,el) + Ps(:,:,co) = phase_mechanical_getP(co,ip,el) enddo - doneAndHappy = mech_RGC_updateState(Ps,Fs,subF,subdt,dPdFs,ip,el) + doneAndHappy = mechanical_RGC_updateState(Ps,Fs,subF,subdt,dPdFs,ip,el) else doneAndHappy = .true. endif -end function mech_updateState +end function mechanical_updateState !-------------------------------------------------------------------------------------------------- !> @brief Write results to file. !-------------------------------------------------------------------------------------------------- -module subroutine mech_results(group_base,h) +module subroutine mechanical_results(group_base,h) use material, only: & material_homogenization_type => homogenization_type @@ -239,7 +239,7 @@ module subroutine mech_results(group_base,h) select case(material_homogenization_type(h)) case(HOMOGENIZATION_rgc_ID) - call mech_RGC_results(homogenization_typeInstance(h),group) + call mechanical_RGC_results(homogenization_typeInstance(h),group) end select @@ -250,7 +250,7 @@ module subroutine mech_results(group_base,h) !call results_writeDataset(group,temp,'P',& ! '1st Piola-Kirchhoff stress','Pa') -end subroutine mech_results +end subroutine mechanical_results end submodule mechanics diff --git a/src/homogenization_mechanics_RGC.f90 b/src/homogenization_mechanical_RGC.f90 similarity index 98% rename from src/homogenization_mechanics_RGC.f90 rename to src/homogenization_mechanical_RGC.f90 index 89a84314b..546c2b65c 100644 --- a/src/homogenization_mechanics_RGC.f90 +++ b/src/homogenization_mechanical_RGC.f90 @@ -71,7 +71,7 @@ contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all necessary fields, reads information from material configuration file !-------------------------------------------------------------------------------------------------- -module subroutine mech_RGC_init(num_homogMech) +module subroutine mechanical_RGC_init(num_homogMech) class(tNode), pointer, intent(in) :: & num_homogMech !< pointer to mechanical homogenization numerics data @@ -155,7 +155,7 @@ module subroutine mech_RGC_init(num_homogMech) prm%N_constituents = homogMech%get_asInts('cluster_size',requiredSize=3) if (homogenization_Nconstituents(h) /= product(prm%N_constituents)) & - call IO_error(211,ext_msg='N_constituents (mech_RGC)') + call IO_error(211,ext_msg='N_constituents (mechanical_RGC)') prm%xi_alpha = homogMech%get_asFloat('xi_alpha') prm%c_alpha = homogMech%get_asFloat('c_alpha') @@ -190,13 +190,13 @@ module subroutine mech_RGC_init(num_homogMech) enddo -end subroutine mech_RGC_init +end subroutine mechanical_RGC_init !-------------------------------------------------------------------------------------------------- !> @brief partitions the deformation gradient onto the constituents !-------------------------------------------------------------------------------------------------- -module subroutine mech_RGC_partitionDeformation(F,avgF,instance,of) +module subroutine mechanical_RGC_partitionDeformation(F,avgF,instance,of) real(pReal), dimension (:,:,:), intent(out) :: F !< partitioned F per grain @@ -229,14 +229,14 @@ module subroutine mech_RGC_partitionDeformation(F,avgF,instance,of) end associate -end subroutine mech_RGC_partitionDeformation +end subroutine mechanical_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,avgF,dt,dPdF,ip,el) result(doneAndHappy) +module function mechanical_RGC_updateState(P,F,avgF,dt,dPdF,ip,el) result(doneAndHappy) logical, dimension(2) :: doneAndHappy real(pReal), dimension(:,:,:), intent(in) :: & P,& !< partitioned stresses @@ -658,7 +658,7 @@ module function mech_RGC_updateState(P,F,avgF,dt,dPdF,ip,el) result(doneAndHappy real(pReal), dimension(6,6) :: C - C = constitutive_homogenizedC(material_phaseAt(grainID,el),material_phaseMemberAt(grainID,ip,el)) + C = phase_homogenizedC(material_phaseAt(grainID,el),material_phaseMemberAt(grainID,ip,el)) equivalentMu = lattice_equivalent_mu(C,'voigt') end function equivalentMu @@ -704,13 +704,13 @@ module function mech_RGC_updateState(P,F,avgF,dt,dPdF,ip,el) result(doneAndHappy end subroutine grainDeformation -end function mech_RGC_updateState +end function mechanical_RGC_updateState !-------------------------------------------------------------------------------------------------- !> @brief derive average stress and stiffness from constituent quantities !-------------------------------------------------------------------------------------------------- -module subroutine mech_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) +module subroutine mechanical_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 @@ -722,13 +722,13 @@ module subroutine mech_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,ins avgP = sum(P,3) /real(product(param(instance)%N_constituents),pReal) dAvgPdAvgF = sum(dPdF,5)/real(product(param(instance)%N_constituents),pReal) -end subroutine mech_RGC_averageStressAndItsTangent +end subroutine mechanical_RGC_averageStressAndItsTangent !-------------------------------------------------------------------------------------------------- !> @brief writes results to HDF5 output file !-------------------------------------------------------------------------------------------------- -module subroutine mech_RGC_results(instance,group) +module subroutine mechanical_RGC_results(instance,group) integer, intent(in) :: instance character(len=*), intent(in) :: group @@ -754,7 +754,7 @@ module subroutine mech_RGC_results(instance,group) enddo outputsLoop end associate -end subroutine mech_RGC_results +end subroutine mechanical_RGC_results !-------------------------------------------------------------------------------------------------- diff --git a/src/homogenization_mechanics_isostrain.f90 b/src/homogenization_mechanical_isostrain.f90 similarity index 91% rename from src/homogenization_mechanics_isostrain.f90 rename to src/homogenization_mechanical_isostrain.f90 index 9634b2a38..f9616de42 100644 --- a/src/homogenization_mechanics_isostrain.f90 +++ b/src/homogenization_mechanical_isostrain.f90 @@ -26,7 +26,7 @@ contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all neccessary fields, reads information from material configuration file !-------------------------------------------------------------------------------------------------- -module subroutine mech_isostrain_init +module subroutine mechanical_isostrain_init integer :: & Ninstances, & @@ -58,7 +58,7 @@ module subroutine mech_isostrain_init case ('avg') prm%mapping = average_ID case default - call IO_error(211,ext_msg='sum'//' (mech_isostrain)') + call IO_error(211,ext_msg='sum'//' (mechanical_isostrain)') end select Nmaterialpoints = count(material_homogenizationAt == h) @@ -70,13 +70,13 @@ module subroutine mech_isostrain_init enddo -end subroutine mech_isostrain_init +end subroutine mechanical_isostrain_init !-------------------------------------------------------------------------------------------------- !> @brief partitions the deformation gradient onto the constituents !-------------------------------------------------------------------------------------------------- -module subroutine mech_isostrain_partitionDeformation(F,avgF) +module subroutine mechanical_isostrain_partitionDeformation(F,avgF) real(pReal), dimension (:,:,:), intent(out) :: F !< partitioned deformation gradient @@ -84,13 +84,13 @@ module subroutine mech_isostrain_partitionDeformation(F,avgF) F = spread(avgF,3,size(F,3)) -end subroutine mech_isostrain_partitionDeformation +end subroutine mechanical_isostrain_partitionDeformation !-------------------------------------------------------------------------------------------------- !> @brief derive average stress and stiffness from constituent quantities !-------------------------------------------------------------------------------------------------- -module subroutine mech_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) +module subroutine mechanical_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 @@ -112,6 +112,6 @@ module subroutine mech_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dP end associate -end subroutine mech_isostrain_averageStressAndItsTangent +end subroutine mechanical_isostrain_averageStressAndItsTangent end submodule isostrain diff --git a/src/homogenization_mechanics_none.f90 b/src/homogenization_mechanical_pass.f90 similarity index 87% rename from src/homogenization_mechanics_none.f90 rename to src/homogenization_mechanical_pass.f90 index ebe2ea8f1..9e8f3e44c 100644 --- a/src/homogenization_mechanics_none.f90 +++ b/src/homogenization_mechanical_pass.f90 @@ -11,14 +11,14 @@ contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all necessary fields, reads information from material configuration file !-------------------------------------------------------------------------------------------------- -module subroutine mech_none_init +module subroutine mechanical_pass_init integer :: & Ninstances, & h, & Nmaterialpoints - print'(/,a)', ' <<<+- homogenization:mechanics:none init -+>>>' + print'(/,a)', ' <<<+- homogenization:mechanics:pass init -+>>>' Ninstances = count(homogenization_type == HOMOGENIZATION_NONE_ID) print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) @@ -27,7 +27,7 @@ module subroutine mech_none_init if(homogenization_type(h) /= HOMOGENIZATION_NONE_ID) cycle if(homogenization_Nconstituents(h) /= 1) & - call IO_error(211,ext_msg='N_constituents (mech_none)') + call IO_error(211,ext_msg='N_constituents (mechanical_pass)') Nmaterialpoints = count(material_homogenizationAt == h) homogState(h)%sizeState = 0 @@ -36,6 +36,6 @@ module subroutine mech_none_init enddo -end subroutine mech_none_init +end subroutine mechanical_pass_init end submodule none diff --git a/src/homogenization_thermal.f90 b/src/homogenization_thermal.f90 index b7170e6ec..67027f52a 100644 --- a/src/homogenization_thermal.f90 +++ b/src/homogenization_thermal.f90 @@ -78,7 +78,7 @@ module subroutine thermal_partition(ce) T = current(material_homogenizationAt2(ce))%T(material_homogenizationMemberAt2(ce)) dot_T = current(material_homogenizationAt2(ce))%dot_T(material_homogenizationMemberAt2(ce)) do co = 1, homogenization_Nconstituents(material_homogenizationAt2(ce)) - call constitutive_thermal_setField(T,dot_T,co,ce) + call phase_thermal_setField(T,dot_T,co,ce) enddo end subroutine thermal_partition @@ -91,7 +91,7 @@ module subroutine thermal_homogenize(ip,el) integer, intent(in) :: ip,el - !call constitutive_thermal_getRate(homogenization_dot_T((el-1)*discretization_nIPs+ip), ip,el) + !call phase_thermal_getRate(homogenization_dot_T((el-1)*discretization_nIPs+ip), ip,el) end subroutine thermal_homogenize @@ -235,7 +235,7 @@ module subroutine thermal_conduction_getSource(Tdot, ip,el) do co = 1, homogenization_Nconstituents(ho) ph = material_phaseAt(co,el) me = material_phasememberAt(co,ip,el) - call constitutive_thermal_getRate(dot_T_temp, ph,me) + call phase_thermal_getRate(dot_T_temp, ph,me) Tdot = Tdot + dot_T_temp enddo diff --git a/src/mesh/DAMASK_mesh.f90 b/src/mesh/DAMASK_mesh.f90 index 7369520c1..5ef0f7a36 100644 --- a/src/mesh/DAMASK_mesh.f90 +++ b/src/mesh/DAMASK_mesh.f90 @@ -18,7 +18,7 @@ program DAMASK_mesh use config use discretization_mesh use FEM_Utilities - use mesh_mech_FEM + use mesh_mechanical_FEM implicit none @@ -242,7 +242,7 @@ program DAMASK_mesh do field = 1, nActiveFields select case (loadCases(1)%fieldBC(field)%ID) case(FIELD_MECH_ID) - call FEM_mech_init(loadCases(1)%fieldBC(field)) + call FEM_mechanical_init(loadCases(1)%fieldBC(field)) end select enddo @@ -306,7 +306,7 @@ program DAMASK_mesh do field = 1, nActiveFields select case (loadCases(currentLoadCase)%fieldBC(field)%ID) case(FIELD_MECH_ID) - call FEM_mech_forward (& + call FEM_mechanical_forward (& guess,timeinc,timeIncOld,loadCases(currentLoadCase)%fieldBC(field)) end select @@ -320,7 +320,7 @@ program DAMASK_mesh do field = 1, nActiveFields select case (loadCases(currentLoadCase)%fieldBC(field)%ID) case(FIELD_MECH_ID) - solres(field) = FEM_mech_solution (& + solres(field) = FEM_mechanical_solution (& incInfo,timeinc,timeIncOld,loadCases(currentLoadCase)%fieldBC(field)) end select diff --git a/src/mesh/FEM_utilities.f90 b/src/mesh/FEM_utilities.f90 index 2f3633e11..4b3be8a42 100644 --- a/src/mesh/FEM_utilities.f90 +++ b/src/mesh/FEM_utilities.f90 @@ -127,12 +127,12 @@ subroutine FEM_utilities_init CHKERRQ(ierr) if(debugPETSc) call PetscOptionsInsertString(PETSC_NULL_OPTIONS,trim(PETSCDEBUG),ierr) CHKERRQ(ierr) - call PetscOptionsInsertString(PETSC_NULL_OPTIONS,'-mech_snes_type newtonls & - &-mech_snes_linesearch_type cp -mech_snes_ksp_ew & - &-mech_snes_ksp_ew_rtol0 0.01 -mech_snes_ksp_ew_rtolmax 0.01 & - &-mech_ksp_type fgmres -mech_ksp_max_it 25 & - &-mech_pc_type ml -mech_mg_levels_ksp_type chebyshev & - &-mech_mg_levels_pc_type sor -mech_pc_ml_nullspace user',ierr) + call PetscOptionsInsertString(PETSC_NULL_OPTIONS,'-mechanical_snes_type newtonls & + &-mechanical_snes_linesearch_type cp -mechanical_snes_ksp_ew & + &-mechanical_snes_ksp_ew_rtol0 0.01 -mechanical_snes_ksp_ew_rtolmax 0.01 & + &-mechanical_ksp_type fgmres -mechanical_ksp_max_it 25 & + &-mechanical_pc_type ml -mechanical_mg_levels_ksp_type chebyshev & + &-mechanical_mg_levels_pc_type sor -mechanical_pc_ml_nullspace user',ierr) CHKERRQ(ierr) call PetscOptionsInsertString(PETSC_NULL_OPTIONS,num_mesh%get_asString('petsc_options',defaultVal=''),ierr) CHKERRQ(ierr) diff --git a/src/mesh/mesh_mech_FEM.f90 b/src/mesh/mesh_mech_FEM.f90 index e19c35998..c811d842b 100644 --- a/src/mesh/mesh_mech_FEM.f90 +++ b/src/mesh/mesh_mech_FEM.f90 @@ -4,7 +4,7 @@ !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH !> @brief FEM PETSc solver !-------------------------------------------------------------------------------------------------- -module mesh_mech_FEM +module mesh_mechanical_FEM #include #include #include @@ -50,7 +50,7 @@ module mesh_mech_FEM type(tNumerics), private :: num !-------------------------------------------------------------------------------------------------- ! PETSc data - SNES :: mech_snes + SNES :: mechanical_snes Vec :: solution, solution_rate, solution_local PetscInt :: dimPlex, cellDof, nQuadrature, nBasis PetscReal, allocatable, target :: qPoints(:), qWeights(:) @@ -65,20 +65,20 @@ module mesh_mech_FEM real(pReal), parameter :: eps = 1.0e-18_pReal public :: & - FEM_mech_init, & - FEM_mech_solution, & - FEM_mech_forward + FEM_mechanical_init, & + FEM_mechanical_solution, & + FEM_mechanical_forward contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all neccessary fields and fills them with data !-------------------------------------------------------------------------------------------------- -subroutine FEM_mech_init(fieldBC) +subroutine FEM_mechanical_init(fieldBC) type(tFieldBC), intent(in) :: fieldBC - DM :: mech_mesh + DM :: mechanical_mesh PetscFE :: mechFE PetscQuadrature :: mechQuad, functional PetscDS :: mechDS @@ -126,8 +126,8 @@ subroutine FEM_mech_init(fieldBC) !-------------------------------------------------------------------------------------------------- ! Setup FEM mech mesh - call DMClone(geomMesh,mech_mesh,ierr); CHKERRQ(ierr) - call DMGetDimension(mech_mesh,dimPlex,ierr); CHKERRQ(ierr) + call DMClone(geomMesh,mechanical_mesh,ierr); CHKERRQ(ierr) + call DMGetDimension(mechanical_mesh,dimPlex,ierr); CHKERRQ(ierr) !-------------------------------------------------------------------------------------------------- ! Setup FEM mech discretization @@ -146,22 +146,22 @@ subroutine FEM_mech_init(fieldBC) call PetscFESetQuadrature(mechFE,mechQuad,ierr); CHKERRQ(ierr) call PetscFEGetDimension(mechFE,nBasis,ierr); CHKERRQ(ierr) nBasis = nBasis/nc - call DMAddField(mech_mesh,PETSC_NULL_DMLABEL,mechFE,ierr); CHKERRQ(ierr) - call DMCreateDS(mech_mesh,ierr); CHKERRQ(ierr) - call DMGetDS(mech_mesh,mechDS,ierr); CHKERRQ(ierr) + call DMAddField(mechanical_mesh,PETSC_NULL_DMLABEL,mechFE,ierr); CHKERRQ(ierr) + call DMCreateDS(mechanical_mesh,ierr); CHKERRQ(ierr) + call DMGetDS(mechanical_mesh,mechDS,ierr); CHKERRQ(ierr) call PetscDSGetTotalDimension(mechDS,cellDof,ierr); CHKERRQ(ierr) call PetscFEDestroy(mechFE,ierr); CHKERRQ(ierr) call PetscQuadratureDestroy(mechQuad,ierr); CHKERRQ(ierr) !-------------------------------------------------------------------------------------------------- ! Setup FEM mech boundary conditions - call DMGetLabel(mech_mesh,'Face Sets',BCLabel,ierr); CHKERRQ(ierr) - call DMPlexLabelComplete(mech_mesh,BCLabel,ierr); CHKERRQ(ierr) - call DMGetLocalSection(mech_mesh,section,ierr); CHKERRQ(ierr) + call DMGetLabel(mechanical_mesh,'Face Sets',BCLabel,ierr); CHKERRQ(ierr) + call DMPlexLabelComplete(mechanical_mesh,BCLabel,ierr); CHKERRQ(ierr) + call DMGetLocalSection(mechanical_mesh,section,ierr); CHKERRQ(ierr) allocate(pnumComp(1), source=dimPlex) allocate(pnumDof(0:dimPlex), source = 0) do topologDim = 0, dimPlex - call DMPlexGetDepthStratum(mech_mesh,topologDim,cellStart,cellEnd,ierr) + call DMPlexGetDepthStratum(mechanical_mesh,topologDim,cellStart,cellEnd,ierr) CHKERRQ(ierr) call PetscSectionGetDof(section,cellStart,pnumDof(topologDim),ierr) CHKERRQ(ierr) @@ -179,10 +179,10 @@ subroutine FEM_mech_init(fieldBC) numBC = numBC + 1 call ISCreateGeneral(PETSC_COMM_WORLD,1,[field-1],PETSC_COPY_VALUES,pbcComps(numBC),ierr) CHKERRQ(ierr) - call DMGetStratumSize(mech_mesh,'Face Sets',mesh_boundaries(faceSet),bcSize,ierr) + call DMGetStratumSize(mechanical_mesh,'Face Sets',mesh_boundaries(faceSet),bcSize,ierr) CHKERRQ(ierr) if (bcSize > 0) then - call DMGetStratumIS(mech_mesh,'Face Sets',mesh_boundaries(faceSet),bcPoint,ierr) + call DMGetStratumIS(mechanical_mesh,'Face Sets',mesh_boundaries(faceSet),bcPoint,ierr) CHKERRQ(ierr) call ISGetIndicesF90(bcPoint,pBcPoint,ierr); CHKERRQ(ierr) call ISCreateGeneral(PETSC_COMM_WORLD,bcSize,pBcPoint,PETSC_COPY_VALUES,pbcPoints(numBC),ierr) @@ -195,32 +195,32 @@ subroutine FEM_mech_init(fieldBC) endif endif enddo; enddo - call DMPlexCreateSection(mech_mesh,nolabel,pNumComp,pNumDof, & + call DMPlexCreateSection(mechanical_mesh,nolabel,pNumComp,pNumDof, & numBC,pBcField,pBcComps,pBcPoints,PETSC_NULL_IS,section,ierr) CHKERRQ(ierr) - call DMSetSection(mech_mesh,section,ierr); CHKERRQ(ierr) + call DMSetSection(mechanical_mesh,section,ierr); CHKERRQ(ierr) do faceSet = 1, numBC call ISDestroy(pbcPoints(faceSet),ierr); CHKERRQ(ierr) enddo !-------------------------------------------------------------------------------------------------- ! initialize solver specific parts of PETSc - call SNESCreate(PETSC_COMM_WORLD,mech_snes,ierr);CHKERRQ(ierr) - call SNESSetOptionsPrefix(mech_snes,'mech_',ierr);CHKERRQ(ierr) - call SNESSetDM(mech_snes,mech_mesh,ierr); CHKERRQ(ierr) !< set the mesh for non-linear solver - call DMCreateGlobalVector(mech_mesh,solution ,ierr); CHKERRQ(ierr) !< locally owned displacement Dofs - call DMCreateGlobalVector(mech_mesh,solution_rate ,ierr); CHKERRQ(ierr) !< locally owned velocity Dofs to guess solution at next load step - call DMCreateLocalVector (mech_mesh,solution_local ,ierr); CHKERRQ(ierr) !< locally owned velocity Dofs to guess solution at next load step - call DMSNESSetFunctionLocal(mech_mesh,FEM_mech_formResidual,PETSC_NULL_VEC,ierr) !< function to evaluate residual forces + call SNESCreate(PETSC_COMM_WORLD,mechanical_snes,ierr);CHKERRQ(ierr) + call SNESSetOptionsPrefix(mechanical_snes,'mechanical_',ierr);CHKERRQ(ierr) + call SNESSetDM(mechanical_snes,mechanical_mesh,ierr); CHKERRQ(ierr) !< set the mesh for non-linear solver + call DMCreateGlobalVector(mechanical_mesh,solution ,ierr); CHKERRQ(ierr) !< locally owned displacement Dofs + call DMCreateGlobalVector(mechanical_mesh,solution_rate ,ierr); CHKERRQ(ierr) !< locally owned velocity Dofs to guess solution at next load step + call DMCreateLocalVector (mechanical_mesh,solution_local ,ierr); CHKERRQ(ierr) !< locally owned velocity Dofs to guess solution at next load step + call DMSNESSetFunctionLocal(mechanical_mesh,FEM_mechanical_formResidual,PETSC_NULL_VEC,ierr) !< function to evaluate residual forces CHKERRQ(ierr) - call DMSNESSetJacobianLocal(mech_mesh,FEM_mech_formJacobian,PETSC_NULL_VEC,ierr) !< function to evaluate stiffness matrix + call DMSNESSetJacobianLocal(mechanical_mesh,FEM_mechanical_formJacobian,PETSC_NULL_VEC,ierr) !< function to evaluate stiffness matrix CHKERRQ(ierr) - call SNESSetMaxLinearSolveFailures(mech_snes, huge(1), ierr); CHKERRQ(ierr) !< ignore linear solve failures - call SNESSetConvergenceTest(mech_snes,FEM_mech_converged,PETSC_NULL_VEC,PETSC_NULL_FUNCTION,ierr) + call SNESSetMaxLinearSolveFailures(mechanical_snes, huge(1), ierr); CHKERRQ(ierr) !< ignore linear solve failures + call SNESSetConvergenceTest(mechanical_snes,FEM_mechanical_converged,PETSC_NULL_VEC,PETSC_NULL_FUNCTION,ierr) CHKERRQ(ierr) - call SNESSetTolerances(mech_snes,1.0,0.0,0.0,num%itmax,num%itmax,ierr) + call SNESSetTolerances(mechanical_snes,1.0,0.0,0.0,num%itmax,num%itmax,ierr) CHKERRQ(ierr) - call SNESSetFromOptions(mech_snes,ierr); CHKERRQ(ierr) + call SNESSetFromOptions(mechanical_snes,ierr); CHKERRQ(ierr) !-------------------------------------------------------------------------------------------------- ! init fields @@ -236,11 +236,11 @@ subroutine FEM_mech_init(fieldBC) call PetscDSGetDiscretization(mechDS,0,mechFE,ierr) CHKERRQ(ierr) call PetscFEGetDualSpace(mechFE,mechDualSpace,ierr); CHKERRQ(ierr) - call DMPlexGetHeightStratum(mech_mesh,0,cellStart,cellEnd,ierr) + call DMPlexGetHeightStratum(mechanical_mesh,0,cellStart,cellEnd,ierr) CHKERRQ(ierr) do cell = cellStart, cellEnd-1 !< loop over all elements x_scal = 0.0_pReal - call DMPlexComputeCellGeometryAffineFEM(mech_mesh,cell,pV0,pCellJ,pInvcellJ,detJ,ierr) + call DMPlexComputeCellGeometryAffineFEM(mechanical_mesh,cell,pV0,pCellJ,pInvcellJ,detJ,ierr) CHKERRQ(ierr) cellJMat = reshape(pCellJ,shape=[dimPlex,dimPlex]) do basis = 0, nBasis*dimPlex-1, dimPlex @@ -251,17 +251,17 @@ subroutine FEM_mech_init(fieldBC) x_scal(basis+1:basis+dimPlex) = pV0 + matmul(transpose(cellJMat),nodalPointsP + 1.0_pReal) enddo px_scal => x_scal - call DMPlexVecSetClosure(mech_mesh,section,solution_local,cell,px_scal,5,ierr) + call DMPlexVecSetClosure(mechanical_mesh,section,solution_local,cell,px_scal,5,ierr) CHKERRQ(ierr) enddo -end subroutine FEM_mech_init +end subroutine FEM_mechanical_init !-------------------------------------------------------------------------------------------------- !> @brief solution for the FEM load step !-------------------------------------------------------------------------------------------------- -type(tSolutionState) function FEM_mech_solution( & +type(tSolutionState) function FEM_mechanical_solution( & incInfoIn,timeinc,timeinc_old,fieldBC) !-------------------------------------------------------------------------------------------------- @@ -278,35 +278,35 @@ type(tSolutionState) function FEM_mech_solution( & SNESConvergedReason :: reason incInfo = incInfoIn - FEM_mech_solution%converged =.false. + FEM_mechanical_solution%converged =.false. !-------------------------------------------------------------------------------------------------- ! set module wide availabe data params%timeinc = timeinc params%fieldBC = fieldBC - call SNESSolve(mech_snes,PETSC_NULL_VEC,solution,ierr); CHKERRQ(ierr) ! solve mech_snes based on solution guess (result in solution) - call SNESGetConvergedReason(mech_snes,reason,ierr); CHKERRQ(ierr) ! solution converged? + call SNESSolve(mechanical_snes,PETSC_NULL_VEC,solution,ierr); CHKERRQ(ierr) ! solve mechanical_snes based on solution guess (result in solution) + call SNESGetConvergedReason(mechanical_snes,reason,ierr); CHKERRQ(ierr) ! solution converged? terminallyIll = .false. if (reason < 1) then ! 0: still iterating (will not occur), negative -> convergence error - FEM_mech_solution%converged = .false. - FEM_mech_solution%iterationsNeeded = num%itmax + FEM_mechanical_solution%converged = .false. + FEM_mechanical_solution%iterationsNeeded = num%itmax else ! >= 1 proper convergence (or terminally ill) - FEM_mech_solution%converged = .true. - call SNESGetIterationNumber(mech_snes,FEM_mech_solution%iterationsNeeded,ierr) + FEM_mechanical_solution%converged = .true. + call SNESGetIterationNumber(mechanical_snes,FEM_mechanical_solution%iterationsNeeded,ierr) CHKERRQ(ierr) endif print'(/,a)', ' ===========================================================================' flush(IO_STDOUT) -end function FEM_mech_solution +end function FEM_mechanical_solution !-------------------------------------------------------------------------------------------------- !> @brief forms the FEM residual vector !-------------------------------------------------------------------------------------------------- -subroutine FEM_mech_formResidual(dm_local,xx_local,f_local,dummy,ierr) +subroutine FEM_mechanical_formResidual(dm_local,xx_local,f_local,dummy,ierr) DM :: dm_local PetscObject,intent(in) :: dummy @@ -431,13 +431,13 @@ subroutine FEM_mech_formResidual(dm_local,xx_local,f_local,dummy,ierr) enddo call DMRestoreLocalVector(dm_local,x_local,ierr); CHKERRQ(ierr) -end subroutine FEM_mech_formResidual +end subroutine FEM_mechanical_formResidual !-------------------------------------------------------------------------------------------------- !> @brief forms the FEM stiffness matrix !-------------------------------------------------------------------------------------------------- -subroutine FEM_mech_formJacobian(dm_local,xx_local,Jac_pre,Jac,dummy,ierr) +subroutine FEM_mechanical_formJacobian(dm_local,xx_local,Jac_pre,Jac,dummy,ierr) DM :: dm_local @@ -575,13 +575,13 @@ subroutine FEM_mech_formJacobian(dm_local,xx_local,Jac_pre,Jac,dummy,ierr) call MatSetNearNullSpace(Jac,matnull,ierr); CHKERRQ(ierr) call MatNullSpaceDestroy(matnull,ierr); CHKERRQ(ierr) -end subroutine FEM_mech_formJacobian +end subroutine FEM_mechanical_formJacobian !-------------------------------------------------------------------------------------------------- !> @brief forwarding routine !-------------------------------------------------------------------------------------------------- -subroutine FEM_mech_forward(guess,timeinc,timeinc_old,fieldBC) +subroutine FEM_mechanical_forward(guess,timeinc,timeinc_old,fieldBC) type(tFieldBC), intent(in) :: & fieldBC @@ -603,7 +603,7 @@ subroutine FEM_mech_forward(guess,timeinc,timeinc_old,fieldBC) if (guess .and. .not. cutBack) then ForwardData = .True. homogenization_F0 = homogenization_F - call SNESGetDM(mech_snes,dm_local,ierr); CHKERRQ(ierr) !< retrieve mesh info from mech_snes into dm_local + call SNESGetDM(mechanical_snes,dm_local,ierr); CHKERRQ(ierr) !< retrieve mesh info from mechanical_snes into dm_local call DMGetSection(dm_local,section,ierr); CHKERRQ(ierr) call DMGetLocalVector(dm_local,x_local,ierr); CHKERRQ(ierr) call VecSet(x_local,0.0_pReal,ierr); CHKERRQ(ierr) @@ -634,13 +634,13 @@ subroutine FEM_mech_forward(guess,timeinc,timeinc_old,fieldBC) call VecCopy(solution_rate,solution,ierr); CHKERRQ(ierr) call VecScale(solution,timeinc,ierr); CHKERRQ(ierr) -end subroutine FEM_mech_forward +end subroutine FEM_mechanical_forward !-------------------------------------------------------------------------------------------------- !> @brief reporting !-------------------------------------------------------------------------------------------------- -subroutine FEM_mech_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,dummy,ierr) +subroutine FEM_mechanical_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,dummy,ierr) SNES :: snes_local PetscInt :: PETScIter @@ -662,6 +662,6 @@ subroutine FEM_mech_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,dumm ' Piola--Kirchhoff stress / MPa =',transpose(P_av)*1.e-6_pReal flush(IO_STDOUT) -end subroutine FEM_mech_converged +end subroutine FEM_mechanical_converged -end module mesh_mech_FEM +end module mesh_mechanical_FEM diff --git a/src/phase.f90 b/src/phase.f90 index dd14e0831..628bfc443 100644 --- a/src/phase.f90 +++ b/src/phase.f90 @@ -62,7 +62,7 @@ module phase phase_Nsources, & !< number of source mechanisms active in each phase phase_Nkinematics, & !< number of kinematic mechanisms active in each phase phase_NstiffnessDegradations, & !< number of stiffness degradation mechanisms active in each phase - phase_plasticityInstance, & !< instance of particular plasticity of each phase + phase_plasticInstance, & !< instance of particular plasticity of each phase phase_elasticityInstance !< instance of particular elasticity of each phase logical, dimension(:), allocatable, public :: & ! ToDo: should be protected (bug in Intel Compiler) @@ -75,15 +75,15 @@ module phase integer, public, protected :: & - constitutive_plasticity_maxSizeDotState, & - constitutive_source_maxSizeDotState + phase_plasticity_maxSizeDotState, & + phase_source_maxSizeDotState interface ! == cleaned:begin ================================================================================= - module subroutine mech_init(phases) + module subroutine mechanical_init(phases) class(tNode), pointer :: phases - end subroutine mech_init + end subroutine mechanical_init module subroutine damage_init end subroutine damage_init @@ -93,83 +93,83 @@ module phase end subroutine thermal_init - module subroutine mech_results(group,ph) + module subroutine mechanical_results(group,ph) character(len=*), intent(in) :: group integer, intent(in) :: ph - end subroutine mech_results + end subroutine mechanical_results module subroutine damage_results(group,ph) character(len=*), intent(in) :: group integer, intent(in) :: ph end subroutine damage_results - module subroutine mech_windForward(ph,me) + module subroutine mechanical_windForward(ph,me) integer, intent(in) :: ph, me - end subroutine mech_windForward + end subroutine mechanical_windForward - module subroutine mech_forward() - end subroutine mech_forward + module subroutine mechanical_forward() + end subroutine mechanical_forward module subroutine thermal_forward() end subroutine thermal_forward - module subroutine mech_restore(ce,includeL) + module subroutine mechanical_restore(ce,includeL) integer, intent(in) :: ce logical, intent(in) :: includeL - end subroutine mech_restore + end subroutine mechanical_restore - module function constitutive_mech_dPdF(dt,co,ip,el) result(dPdF) + module function phase_mechanical_dPdF(dt,co,ip,el) result(dPdF) real(pReal), intent(in) :: dt integer, intent(in) :: & co, & !< counter in constituent loop ip, & !< counter in integration point loop el !< counter in element loop real(pReal), dimension(3,3,3,3) :: dPdF - end function constitutive_mech_dPdF + end function phase_mechanical_dPdF - module subroutine mech_restartWrite(groupHandle,ph) + module subroutine mechanical_restartWrite(groupHandle,ph) integer(HID_T), intent(in) :: groupHandle integer, intent(in) :: ph - end subroutine mech_restartWrite + end subroutine mechanical_restartWrite - module subroutine mech_restartRead(groupHandle,ph) + module subroutine mechanical_restartRead(groupHandle,ph) integer(HID_T), intent(in) :: groupHandle integer, intent(in) :: ph - end subroutine mech_restartRead + end subroutine mechanical_restartRead - module function mech_S(ph,me) result(S) + module function mechanical_S(ph,me) result(S) integer, intent(in) :: ph,me real(pReal), dimension(3,3) :: S - end function mech_S + end function mechanical_S - module function mech_L_p(ph,me) result(L_p) + module function mechanical_L_p(ph,me) result(L_p) integer, intent(in) :: ph,me real(pReal), dimension(3,3) :: L_p - end function mech_L_p + end function mechanical_L_p - module function constitutive_mech_getF(co,ip,el) result(F) + module function phase_mechanical_getF(co,ip,el) result(F) integer, intent(in) :: co, ip, el real(pReal), dimension(3,3) :: F - end function constitutive_mech_getF + end function phase_mechanical_getF - module function mech_F_e(ph,me) result(F_e) + module function mechanical_F_e(ph,me) result(F_e) integer, intent(in) :: ph,me real(pReal), dimension(3,3) :: F_e - end function mech_F_e + end function mechanical_F_e - module function constitutive_mech_getP(co,ip,el) result(P) + module function phase_mechanical_getP(co,ip,el) result(P) integer, intent(in) :: co, ip, el real(pReal), dimension(3,3) :: P - end function constitutive_mech_getP + end function phase_mechanical_getP - module function constitutive_damage_get_phi(co,ip,el) result(phi) + module function phase_damage_get_phi(co,ip,el) result(phi) integer, intent(in) :: co, ip, el real(pReal) :: phi - end function constitutive_damage_get_phi + end function phase_damage_get_phi module function thermal_T(ph,me) result(T) integer, intent(in) :: ph,me @@ -182,20 +182,20 @@ module phase end function thermal_dot_T - module subroutine constitutive_mech_setF(F,co,ip,el) + module subroutine phase_mechanical_setF(F,co,ip,el) real(pReal), dimension(3,3), intent(in) :: F integer, intent(in) :: co, ip, el - end subroutine constitutive_mech_setF + end subroutine phase_mechanical_setF - module subroutine constitutive_thermal_setField(T,dot_T, co,ce) + module subroutine phase_thermal_setField(T,dot_T, co,ce) real(pReal), intent(in) :: T, dot_T integer, intent(in) :: ce, co - end subroutine constitutive_thermal_setField + end subroutine phase_thermal_setField - module subroutine constitutive_damage_set_phi(phi,co,ce) + module subroutine phase_damage_set_phi(phi,co,ce) real(pReal), intent(in) :: phi integer, intent(in) :: co, ce - end subroutine constitutive_damage_set_phi + end subroutine phase_damage_set_phi ! == cleaned:end =================================================================================== @@ -222,13 +222,13 @@ module phase logical :: converged_ end function crystallite_stress - module function constitutive_homogenizedC(ph,me) result(C) + module function phase_homogenizedC(ph,me) result(C) integer, intent(in) :: ph, me real(pReal), dimension(6,6) :: C - end function constitutive_homogenizedC + end function phase_homogenizedC - module subroutine constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) + module subroutine phase_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) integer, intent(in) :: & ip, & !< integration point number el !< element number @@ -237,13 +237,13 @@ module phase real(pReal), intent(inout) :: & phiDot, & dPhiDot_dPhi - end subroutine constitutive_damage_getRateAndItsTangents + end subroutine phase_damage_getRateAndItsTangents - module subroutine constitutive_thermal_getRate(TDot, ph,me) + module subroutine phase_thermal_getRate(TDot, ph,me) integer, intent(in) :: ph, me real(pReal), intent(out) :: & TDot - end subroutine constitutive_thermal_getRate + end subroutine phase_thermal_getRate module subroutine plastic_nonlocal_updateCompatibility(orientation,instance,i,e) integer, intent(in) :: & @@ -281,39 +281,39 @@ module phase #endif public :: & - constitutive_init, & - constitutive_homogenizedC, & - constitutive_damage_getRateAndItsTangents, & - constitutive_thermal_getRate, & - constitutive_results, & - constitutive_allocateState, & - constitutive_forward, & - constitutive_restore, & + phase_init, & + phase_homogenizedC, & + phase_damage_getRateAndItsTangents, & + phase_thermal_getRate, & + phase_results, & + phase_allocateState, & + phase_forward, & + phase_restore, & plastic_nonlocal_updateCompatibility, & converged, & crystallite_init, & crystallite_stress, & thermal_stress, & - constitutive_mech_dPdF, & + phase_mechanical_dPdF, & crystallite_orientations, & crystallite_push33ToRef, & - constitutive_restartWrite, & - constitutive_restartRead, & + phase_restartWrite, & + phase_restartRead, & integrateDamageState, & - constitutive_thermal_setField, & - constitutive_damage_set_phi, & - constitutive_damage_get_phi, & - constitutive_mech_getP, & - constitutive_mech_setF, & - constitutive_mech_getF, & - constitutive_windForward + phase_thermal_setField, & + phase_damage_set_phi, & + phase_damage_get_phi, & + phase_mechanical_getP, & + phase_mechanical_setF, & + phase_mechanical_getF, & + phase_windForward contains !-------------------------------------------------------------------------------------------------- !> @brief Initialze constitutive models for individual physics !-------------------------------------------------------------------------------------------------- -subroutine constitutive_init +subroutine phase_init integer :: & ph, & !< counter in phase loop @@ -336,12 +336,12 @@ subroutine constitutive_init phases => config_material%get('phase') - call mech_init(phases) + call mechanical_init(phases) call damage_init call thermal_init(phases) - constitutive_source_maxSizeDotState = 0 + phase_source_maxSizeDotState = 0 PhaseLoop2:do ph = 1,phases%length !-------------------------------------------------------------------------------------------------- ! partition and initialize state @@ -350,18 +350,18 @@ subroutine constitutive_init damageState(ph)%p(so)%state = damageState(ph)%p(so)%state0 end forall - constitutive_source_maxSizeDotState = max(constitutive_source_maxSizeDotState, & + phase_source_maxSizeDotState = max(phase_source_maxSizeDotState, & maxval(damageState(ph)%p%sizeDotState)) enddo PhaseLoop2 - constitutive_plasticity_maxSizeDotState = maxval(plasticState%sizeDotState) + phase_plasticity_maxSizeDotState = maxval(plasticState%sizeDotState) -end subroutine constitutive_init +end subroutine phase_init !-------------------------------------------------------------------------------------------------- !> @brief Allocate the components of the state structure for a given phase !-------------------------------------------------------------------------------------------------- -subroutine constitutive_allocateState(state, & +subroutine phase_allocateState(state, & Nconstituents,sizeState,sizeDotState,sizeDeltaState) class(tState), intent(out) :: & @@ -387,13 +387,13 @@ subroutine constitutive_allocateState(state, & allocate(state%deltaState (sizeDeltaState,Nconstituents), source=0.0_pReal) -end subroutine constitutive_allocateState +end subroutine phase_allocateState !-------------------------------------------------------------------------------------------------- !> @brief Restore data after homog cutback. !-------------------------------------------------------------------------------------------------- -subroutine constitutive_restore(ce,includeL) +subroutine phase_restore(ce,includeL) logical, intent(in) :: includeL integer, intent(in) :: ce @@ -410,21 +410,21 @@ subroutine constitutive_restore(ce,includeL) enddo enddo - call mech_restore(ce,includeL) + call mechanical_restore(ce,includeL) -end subroutine constitutive_restore +end subroutine phase_restore !-------------------------------------------------------------------------------------------------- !> @brief Forward data after successful increment. ! ToDo: Any guessing for the current states possible? !-------------------------------------------------------------------------------------------------- -subroutine constitutive_forward() +subroutine phase_forward() integer :: ph, so - call mech_forward() + call mechanical_forward() call thermal_forward() do ph = 1, size(damageState) @@ -432,13 +432,13 @@ subroutine constitutive_forward() damageState(ph)%p(so)%state0 = damageState(ph)%p(so)%state enddo; enddo -end subroutine constitutive_forward +end subroutine phase_forward !-------------------------------------------------------------------------------------------------- !> @brief writes constitutive results to HDF5 output file !-------------------------------------------------------------------------------------------------- -subroutine constitutive_results() +subroutine phase_results() integer :: ph character(len=:), allocatable :: group @@ -451,12 +451,12 @@ subroutine constitutive_results() group = '/current/phase/'//trim(material_name_phase(ph))//'/' call results_closeGroup(results_addGroup(group)) - call mech_results(group,ph) + call mechanical_results(group,ph) call damage_results(group,ph) enddo -end subroutine constitutive_results +end subroutine phase_results !-------------------------------------------------------------------------------------------------- @@ -557,7 +557,7 @@ end subroutine crystallite_init !-------------------------------------------------------------------------------------------------- !> @brief Wind homog inc forward. !-------------------------------------------------------------------------------------------------- -subroutine constitutive_windForward(ip,el) +subroutine phase_windForward(ip,el) integer, intent(in) :: & ip, & !< integration point number @@ -572,7 +572,7 @@ subroutine constitutive_windForward(ip,el) ph = material_phaseAt(co,el) me = material_phaseMemberAt(co,ip,el) - call mech_windForward(ph,me) + call mechanical_windForward(ph,me) do so = 1, phase_Nsources(material_phaseAt(co,el)) damageState(ph)%p(so)%state0(:,me) = damageState(ph)%p(so)%state(:,me) @@ -580,7 +580,7 @@ subroutine constitutive_windForward(ip,el) enddo -end subroutine constitutive_windForward +end subroutine phase_windForward !-------------------------------------------------------------------------------------------------- @@ -595,11 +595,11 @@ subroutine crystallite_orientations(co,ip,el) call crystallite_orientation(co,ip,el)%fromMatrix(transpose(math_rotationalPart(& - mech_F_e(material_phaseAt(co,el),material_phaseMemberAt(co,ip,el))))) + mechanical_F_e(material_phaseAt(co,el),material_phaseMemberAt(co,ip,el))))) if (plasticState(material_phaseAt(1,el))%nonlocal) & call plastic_nonlocal_updateCompatibility(crystallite_orientation, & - phase_plasticityInstance(material_phaseAt(1,el)),ip,el) + phase_plasticInstance(material_phaseAt(1,el)),ip,el) end subroutine crystallite_orientations @@ -620,7 +620,7 @@ function crystallite_push33ToRef(co,ip,el, tensor33) real(pReal), dimension(3,3) :: T - T = matmul(material_orientation0(co,ip,el)%asMatrix(),transpose(math_inv33(constitutive_mech_getF(co,ip,el)))) ! ToDo: initial orientation correct? + T = matmul(material_orientation0(co,ip,el)%asMatrix(),transpose(math_inv33(phase_mechanical_getF(co,ip,el)))) ! ToDo: initial orientation correct? crystallite_push33ToRef = matmul(transpose(T),matmul(tensor33,T)) @@ -648,7 +648,7 @@ end function converged !> @brief Write current restart information (Field and constitutive data) to file. ! ToDo: Merge data into one file for MPI !-------------------------------------------------------------------------------------------------- -subroutine constitutive_restartWrite(fileHandle) +subroutine phase_restartWrite(fileHandle) integer(HID_T), intent(in) :: fileHandle @@ -662,7 +662,7 @@ subroutine constitutive_restartWrite(fileHandle) groupHandle(2) = HDF5_addGroup(groupHandle(1),material_name_phase(ph)) - call mech_restartWrite(groupHandle(2),ph) + call mechanical_restartWrite(groupHandle(2),ph) call HDF5_closeGroup(groupHandle(2)) @@ -670,14 +670,14 @@ subroutine constitutive_restartWrite(fileHandle) call HDF5_closeGroup(groupHandle(1)) -end subroutine constitutive_restartWrite +end subroutine phase_restartWrite !-------------------------------------------------------------------------------------------------- !> @brief Read data for restart ! ToDo: Merge data into one file for MPI !-------------------------------------------------------------------------------------------------- -subroutine constitutive_restartRead(fileHandle) +subroutine phase_restartRead(fileHandle) integer(HID_T), intent(in) :: fileHandle @@ -691,7 +691,7 @@ subroutine constitutive_restartRead(fileHandle) groupHandle(2) = HDF5_openGroup(groupHandle(1),material_name_phase(ph)) - call mech_restartRead(groupHandle(2),ph) + call mechanical_restartRead(groupHandle(2),ph) call HDF5_closeGroup(groupHandle(2)) @@ -699,7 +699,7 @@ subroutine constitutive_restartRead(fileHandle) call HDF5_closeGroup(groupHandle(1)) -end subroutine constitutive_restartRead +end subroutine phase_restartRead end module phase diff --git a/src/phase_damage.f90 b/src/phase_damage.f90 index 014347c1a..57145c550 100644 --- a/src/phase_damage.f90 +++ b/src/phase_damage.f90 @@ -196,7 +196,7 @@ end subroutine damage_init !---------------------------------------------------------------------------------------------- !< @brief returns local part of nonlocal damage driving force !---------------------------------------------------------------------------------------------- -module subroutine constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) +module subroutine phase_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) integer, intent(in) :: & ip, & !< integration point number @@ -246,7 +246,7 @@ module subroutine constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi enddo enddo -end subroutine constitutive_damage_getRateAndItsTangents +end subroutine phase_damage_getRateAndItsTangents @@ -272,9 +272,9 @@ module function integrateDamageState(dt,co,ip,el) result(broken) size_so real(pReal) :: & zeta - real(pReal), dimension(constitutive_source_maxSizeDotState) :: & + real(pReal), dimension(phase_source_maxSizeDotState) :: & r ! state residuum - real(pReal), dimension(constitutive_source_maxSizeDotState,2,maxval(phase_Nsources)) :: source_dotState + real(pReal), dimension(phase_source_maxSizeDotState,2,maxval(phase_Nsources)) :: source_dotState logical :: & converged_ @@ -283,7 +283,7 @@ module function integrateDamageState(dt,co,ip,el) result(broken) me = material_phaseMemberAt(co,ip,el) converged_ = .true. - broken = constitutive_damage_collectDotState(co,ip,el,ph,me) + broken = phase_damage_collectDotState(co,ip,el,ph,me) if(broken) return do so = 1, phase_Nsources(ph) @@ -300,7 +300,7 @@ module function integrateDamageState(dt,co,ip,el) result(broken) source_dotState(1:size_so(so),1,so) = damageState(ph)%p(so)%dotState(:,me) enddo - broken = constitutive_damage_collectDotState(co,ip,el,ph,me) + broken = phase_damage_collectDotState(co,ip,el,ph,me) if(broken) exit iteration do so = 1, phase_Nsources(ph) @@ -320,7 +320,7 @@ module function integrateDamageState(dt,co,ip,el) result(broken) enddo if(converged_) then - broken = constitutive_damage_deltaState(mech_F_e(ph,me),ph,me) + broken = phase_damage_deltaState(mechanical_F_e(ph,me),ph,me) exit iteration endif @@ -393,7 +393,7 @@ end subroutine damage_results !-------------------------------------------------------------------------------------------------- !> @brief contains the constitutive equation for calculating the rate of change of microstructure !-------------------------------------------------------------------------------------------------- -function constitutive_damage_collectDotState(co,ip,el,ph,me) result(broken) +function phase_damage_collectDotState(co,ip,el,ph,me) result(broken) integer, intent(in) :: & co, & !< component-ID me integration point @@ -419,7 +419,7 @@ function constitutive_damage_collectDotState(co,ip,el,ph,me) result(broken) call anisoductile_dotState(co, ip, el) case (DAMAGE_ANISOBRITTLE_ID) sourceType - call anisobrittle_dotState(mech_S(ph,me),co, ip, el) ! correct stress? + call anisobrittle_dotState(mechanical_S(ph,me),co, ip, el) ! correct stress? end select sourceType @@ -427,7 +427,7 @@ function constitutive_damage_collectDotState(co,ip,el,ph,me) result(broken) enddo SourceLoop -end function constitutive_damage_collectDotState +end function phase_damage_collectDotState @@ -435,7 +435,7 @@ end function constitutive_damage_collectDotState !> @brief for constitutive models having an instantaneous change of state !> will return false if delta state is not needed/supported by the constitutive model !-------------------------------------------------------------------------------------------------- -function constitutive_damage_deltaState(Fe, ph, me) result(broken) +function phase_damage_deltaState(Fe, ph, me) result(broken) integer, intent(in) :: & ph, & @@ -457,7 +457,7 @@ function constitutive_damage_deltaState(Fe, ph, me) result(broken) sourceType: select case (phase_source(so,ph)) case (DAMAGE_ISOBRITTLE_ID) sourceType - call source_damage_isoBrittle_deltaState(constitutive_homogenizedC(ph,me), Fe, ph,me) + call source_damage_isoBrittle_deltaState(phase_homogenizedC(ph,me), Fe, ph,me) broken = any(IEEE_is_NaN(damageState(ph)%p(so)%deltaState(:,me))) if(.not. broken) then myOffset = damageState(ph)%p(so)%offsetDeltaState @@ -470,7 +470,7 @@ function constitutive_damage_deltaState(Fe, ph, me) result(broken) enddo SourceLoop -end function constitutive_damage_deltaState +end function phase_damage_deltaState !-------------------------------------------------------------------------------------------------- @@ -507,7 +507,7 @@ end function source_active !---------------------------------------------------------------------------------------------- !< @brief Set damage parameter !---------------------------------------------------------------------------------------------- -module subroutine constitutive_damage_set_phi(phi,co,ce) +module subroutine phase_damage_set_phi(phi,co,ce) real(pReal), intent(in) :: phi integer, intent(in) :: ce, co @@ -515,17 +515,17 @@ module subroutine constitutive_damage_set_phi(phi,co,ce) current(material_phaseAt2(co,ce))%phi(material_phaseMemberAt2(co,ce)) = phi -end subroutine constitutive_damage_set_phi +end subroutine phase_damage_set_phi -module function constitutive_damage_get_phi(co,ip,el) result(phi) +module function phase_damage_get_phi(co,ip,el) result(phi) integer, intent(in) :: co, ip, el real(pReal) :: phi phi = current(material_phaseAt(co,el))%phi(material_phaseMemberAt(co,ip,el)) -end function constitutive_damage_get_phi +end function phase_damage_get_phi end submodule damagee diff --git a/src/phase_damage_anisobrittle.f90 b/src/phase_damage_anisobrittle.f90 index 2ae5ca951..51ca7786a 100644 --- a/src/phase_damage_anisobrittle.f90 +++ b/src/phase_damage_anisobrittle.f90 @@ -101,7 +101,7 @@ module function anisobrittle_init(source_length) result(mySources) if (any(prm%s_crit < 0.0_pReal)) extmsg = trim(extmsg)//' s_crit' Nconstituents = count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(damageState(p)%p(sourceOffset),Nconstituents,1,1,0) + call phase_allocateState(damageState(p)%p(sourceOffset),Nconstituents,1,1,0) damageState(p)%p(sourceOffset)%atol = src%get_asFloat('anisobrittle_atol',defaultVal=1.0e-3_pReal) if(any(damageState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisobrittle_atol' diff --git a/src/phase_damage_anisoductile.f90 b/src/phase_damage_anisoductile.f90 index 8d5904e6b..e54eff201 100644 --- a/src/phase_damage_anisoductile.f90 +++ b/src/phase_damage_anisoductile.f90 @@ -87,7 +87,7 @@ module function anisoductile_init(source_length) result(mySources) if (any(prm%gamma_crit < 0.0_pReal)) extmsg = trim(extmsg)//' gamma_crit' Nconstituents=count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(damageState(p)%p(sourceOffset),Nconstituents,1,1,0) + call phase_allocateState(damageState(p)%p(sourceOffset),Nconstituents,1,1,0) damageState(p)%p(sourceOffset)%atol = src%get_asFloat('anisoDuctile_atol',defaultVal=1.0e-3_pReal) if(any(damageState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisoductile_atol' diff --git a/src/phase_damage_isobrittle.f90 b/src/phase_damage_isobrittle.f90 index 091377171..529ecd404 100644 --- a/src/phase_damage_isobrittle.f90 +++ b/src/phase_damage_isobrittle.f90 @@ -74,7 +74,7 @@ module function isobrittle_init(source_length) result(mySources) if (prm%W_crit <= 0.0_pReal) extmsg = trim(extmsg)//' W_crit' Nconstituents = count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(damageState(p)%p(sourceOffset),Nconstituents,1,1,1) + call phase_allocateState(damageState(p)%p(sourceOffset),Nconstituents,1,1,1) damageState(p)%p(sourceOffset)%atol = src%get_asFloat('isoBrittle_atol',defaultVal=1.0e-3_pReal) if(any(damageState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isobrittle_atol' diff --git a/src/phase_damage_isoductile.f90 b/src/phase_damage_isoductile.f90 index 6f0a2d0fb..71f99078b 100644 --- a/src/phase_damage_isoductile.f90 +++ b/src/phase_damage_isoductile.f90 @@ -78,7 +78,7 @@ module function isoductile_init(source_length) result(mySources) if (prm%gamma_crit <= 0.0_pReal) extmsg = trim(extmsg)//' gamma_crit' Nconstituents=count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(damageState(p)%p(sourceOffset),Nconstituents,1,1,0) + call phase_allocateState(damageState(p)%p(sourceOffset),Nconstituents,1,1,0) damageState(p)%p(sourceOffset)%atol = src%get_asFloat('isoDuctile_atol',defaultVal=1.0e-3_pReal) if(any(damageState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isoductile_atol' diff --git a/src/phase_mechanics.f90 b/src/phase_mechanical.f90 similarity index 82% rename from src/phase_mechanics.f90 rename to src/phase_mechanical.f90 index a84c1a385..efd4d51d2 100644 --- a/src/phase_mechanics.f90 +++ b/src/phase_mechanical.f90 @@ -31,21 +31,21 @@ submodule(phase) mechanics type(tTensorContainer), dimension(:), allocatable :: & ! current value - constitutive_mech_Fe, & - constitutive_mech_Fi, & - constitutive_mech_Fp, & - constitutive_mech_F, & - constitutive_mech_Li, & - constitutive_mech_Lp, & - constitutive_mech_S, & - constitutive_mech_P, & + phase_mechanical_Fe, & + phase_mechanical_Fi, & + phase_mechanical_Fp, & + phase_mechanical_F, & + phase_mechanical_Li, & + phase_mechanical_Lp, & + phase_mechanical_S, & + phase_mechanical_P, & ! converged value at end of last solver increment - constitutive_mech_Fi0, & - constitutive_mech_Fp0, & - constitutive_mech_F0, & - constitutive_mech_Li0, & - constitutive_mech_Lp0, & - constitutive_mech_S0 + phase_mechanical_Fi0, & + phase_mechanical_Fp0, & + phase_mechanical_F0, & + phase_mechanical_Li0, & + phase_mechanical_Lp0, & + phase_mechanical_S0 integer(kind(PLASTICITY_undefined_ID)), dimension(:), allocatable :: & @@ -97,7 +97,7 @@ submodule(phase) mechanics broken end function plastic_deltaState - module subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, & + module subroutine phase_LiAndItsTangents(Li, dLi_dS, dLi_dFi, & S, Fi, co, ip, el) integer, intent(in) :: & @@ -114,7 +114,7 @@ submodule(phase) mechanics dLi_dS, & !< derivative of Li with respect to S dLi_dFi - end subroutine constitutive_LiAndItsTangents + end subroutine phase_LiAndItsTangents module subroutine plastic_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, & @@ -186,7 +186,7 @@ contains !> @brief Initialize mechanical field related constitutive models !> @details Initialize elasticity, plasticity and stiffness degradation models. !-------------------------------------------------------------------------------------------------- -module subroutine mech_init(phases) +module subroutine mechanical_init(phases) class(tNode), pointer :: & phases @@ -215,38 +215,38 @@ module subroutine mech_init(phases) allocate(phase_NstiffnessDegradations(phases%length),source=0) allocate(output_constituent(phases%length)) - allocate(constitutive_mech_Fe(phases%length)) - allocate(constitutive_mech_Fi(phases%length)) - allocate(constitutive_mech_Fi0(phases%length)) - allocate(constitutive_mech_Fp(phases%length)) - allocate(constitutive_mech_Fp0(phases%length)) - allocate(constitutive_mech_F(phases%length)) - allocate(constitutive_mech_F0(phases%length)) - allocate(constitutive_mech_Li(phases%length)) - allocate(constitutive_mech_Li0(phases%length)) - allocate(constitutive_mech_Lp0(phases%length)) - allocate(constitutive_mech_Lp(phases%length)) - allocate(constitutive_mech_S(phases%length)) - allocate(constitutive_mech_P(phases%length)) - allocate(constitutive_mech_S0(phases%length)) + allocate(phase_mechanical_Fe(phases%length)) + allocate(phase_mechanical_Fi(phases%length)) + allocate(phase_mechanical_Fi0(phases%length)) + allocate(phase_mechanical_Fp(phases%length)) + allocate(phase_mechanical_Fp0(phases%length)) + allocate(phase_mechanical_F(phases%length)) + allocate(phase_mechanical_F0(phases%length)) + allocate(phase_mechanical_Li(phases%length)) + allocate(phase_mechanical_Li0(phases%length)) + allocate(phase_mechanical_Lp0(phases%length)) + allocate(phase_mechanical_Lp(phases%length)) + allocate(phase_mechanical_S(phases%length)) + allocate(phase_mechanical_P(phases%length)) + allocate(phase_mechanical_S0(phases%length)) do ph = 1, phases%length Nconstituents = count(material_phaseAt == ph) * discretization_nIPs - allocate(constitutive_mech_Fi(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_Fe(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_Fi0(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_Fp(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_Fp0(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_Li(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_Li0(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_Lp0(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_Lp(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_S(ph)%data(3,3,Nconstituents),source=0.0_pReal) - allocate(constitutive_mech_P(ph)%data(3,3,Nconstituents),source=0.0_pReal) - allocate(constitutive_mech_S0(ph)%data(3,3,Nconstituents),source=0.0_pReal) - allocate(constitutive_mech_F(ph)%data(3,3,Nconstituents)) - allocate(constitutive_mech_F0(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_Fi(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_Fe(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_Fi0(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_Fp(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_Fp0(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_Li(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_Li0(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_Lp0(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_Lp(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_S(ph)%data(3,3,Nconstituents),source=0.0_pReal) + allocate(phase_mechanical_P(ph)%data(3,3,Nconstituents),source=0.0_pReal) + allocate(phase_mechanical_S0(ph)%data(3,3,Nconstituents),source=0.0_pReal) + allocate(phase_mechanical_F(ph)%data(3,3,Nconstituents)) + allocate(phase_mechanical_F0(ph)%data(3,3,Nconstituents)) phase => phases%get(ph) mech => phase%get('mechanics') @@ -287,17 +287,17 @@ module subroutine mech_init(phases) ph = material_phaseAt(co,el) me = material_phaseMemberAt(co,ip,el) - constitutive_mech_Fp0(ph)%data(1:3,1:3,me) = material_orientation0(co,ip,el)%asMatrix() ! Fp reflects initial orientation (see 10.1016/j.actamat.2006.01.005) - constitutive_mech_Fp0(ph)%data(1:3,1:3,me) = constitutive_mech_Fp0(ph)%data(1:3,1:3,me) & - / math_det33(constitutive_mech_Fp0(ph)%data(1:3,1:3,me))**(1.0_pReal/3.0_pReal) - constitutive_mech_Fi0(ph)%data(1:3,1:3,me) = math_I3 - constitutive_mech_F0(ph)%data(1:3,1:3,me) = math_I3 + phase_mechanical_Fp0(ph)%data(1:3,1:3,me) = material_orientation0(co,ip,el)%asMatrix() ! Fp reflects initial orientation (see 10.1016/j.actamat.2006.01.005) + phase_mechanical_Fp0(ph)%data(1:3,1:3,me) = phase_mechanical_Fp0(ph)%data(1:3,1:3,me) & + / math_det33(phase_mechanical_Fp0(ph)%data(1:3,1:3,me))**(1.0_pReal/3.0_pReal) + phase_mechanical_Fi0(ph)%data(1:3,1:3,me) = math_I3 + phase_mechanical_F0(ph)%data(1:3,1:3,me) = math_I3 - constitutive_mech_Fe(ph)%data(1:3,1:3,me) = math_inv33(matmul(constitutive_mech_Fi0(ph)%data(1:3,1:3,me), & - constitutive_mech_Fp0(ph)%data(1:3,1:3,me))) ! assuming that euler angles are given in internal strain free configuration - constitutive_mech_Fp(ph)%data(1:3,1:3,me) = constitutive_mech_Fp0(ph)%data(1:3,1:3,me) - constitutive_mech_Fi(ph)%data(1:3,1:3,me) = constitutive_mech_Fi0(ph)%data(1:3,1:3,me) - constitutive_mech_F(ph)%data(1:3,1:3,me) = constitutive_mech_F0(ph)%data(1:3,1:3,me) + phase_mechanical_Fe(ph)%data(1:3,1:3,me) = math_inv33(matmul(phase_mechanical_Fi0(ph)%data(1:3,1:3,me), & + phase_mechanical_Fp0(ph)%data(1:3,1:3,me))) ! assuming that euler angles are given in internal strain free configuration + phase_mechanical_Fp(ph)%data(1:3,1:3,me) = phase_mechanical_Fp0(ph)%data(1:3,1:3,me) + phase_mechanical_Fi(ph)%data(1:3,1:3,me) = phase_mechanical_Fi0(ph)%data(1:3,1:3,me) + phase_mechanical_F(ph)%data(1:3,1:3,me) = phase_mechanical_F0(ph)%data(1:3,1:3,me) enddo enddo; enddo @@ -307,14 +307,14 @@ module subroutine mech_init(phases) ! initialize plasticity allocate(plasticState(phases%length)) allocate(phase_plasticity(phases%length),source = PLASTICITY_undefined_ID) - allocate(phase_plasticityInstance(phases%length),source = 0) + allocate(phase_plasticInstance(phases%length),source = 0) allocate(phase_localPlasticity(phases%length), source=.true.) call plastic_init() do ph = 1, phases%length phase_elasticityInstance(ph) = count(phase_elasticity(1:ph) == phase_elasticity(ph)) - phase_plasticityInstance(ph) = count(phase_plasticity(1:ph) == phase_plasticity(ph)) + phase_plasticInstance(ph) = count(phase_plasticity(1:ph) == phase_plasticity(ph)) enddo num_crystallite => config_numerics%get('crystallite',defaultVal=emptyDict) @@ -345,14 +345,14 @@ module subroutine mech_init(phases) call eigendeformation_init(phases) -end subroutine mech_init +end subroutine mechanical_init !-------------------------------------------------------------------------------------------------- !> @brief returns the 2nd Piola-Kirchhoff stress tensor and its tangent with respect to !> the elastic and intermediate deformation gradients using Hooke's law !-------------------------------------------------------------------------------------------------- -subroutine constitutive_hooke_SandItsTangents(S, dS_dFe, dS_dFi, & +subroutine phase_hooke_SandItsTangents(S, dS_dFe, dS_dFi, & Fe, Fi, co, ip, el) integer, intent(in) :: & @@ -376,7 +376,7 @@ subroutine constitutive_hooke_SandItsTangents(S, dS_dFe, dS_dFi, & i, j, ph, me ho = material_homogenizationAt(el) - C = math_66toSym3333(constitutive_homogenizedC(material_phaseAt(co,el),material_phaseMemberAt(co,ip,el))) + C = math_66toSym3333(phase_homogenizedC(material_phaseAt(co,el),material_phaseMemberAt(co,ip,el))) DegradationLoop: do d = 1, phase_NstiffnessDegradations(material_phaseAt(co,el)) degradationType: select case(phase_stiffnessDegradation(d,material_phaseAt(co,el))) @@ -393,10 +393,10 @@ subroutine constitutive_hooke_SandItsTangents(S, dS_dFe, dS_dFi, & dS_dFi(i,j,1:3,1:3) = 2.0_pReal*matmul(matmul(E,Fi),C(i,j,1:3,1:3)) !< dS_ij/dFi_kl = C_ijln * E_km * Fe_mn enddo; enddo -end subroutine constitutive_hooke_SandItsTangents +end subroutine phase_hooke_SandItsTangents -module subroutine mech_results(group,ph) +module subroutine mechanical_results(group,ph) character(len=*), intent(in) :: group integer, intent(in) :: ph @@ -407,28 +407,28 @@ module subroutine mech_results(group,ph) select case(phase_plasticity(ph)) case(PLASTICITY_ISOTROPIC_ID) - call plastic_isotropic_results(phase_plasticityInstance(ph),group//'plastic/') + call plastic_isotropic_results(phase_plasticInstance(ph),group//'plastic/') case(PLASTICITY_PHENOPOWERLAW_ID) - call plastic_phenopowerlaw_results(phase_plasticityInstance(ph),group//'plastic/') + call plastic_phenopowerlaw_results(phase_plasticInstance(ph),group//'plastic/') case(PLASTICITY_KINEHARDENING_ID) - call plastic_kinehardening_results(phase_plasticityInstance(ph),group//'plastic/') + call plastic_kinehardening_results(phase_plasticInstance(ph),group//'plastic/') case(PLASTICITY_DISLOTWIN_ID) - call plastic_dislotwin_results(phase_plasticityInstance(ph),group//'plastic/') + call plastic_dislotwin_results(phase_plasticInstance(ph),group//'plastic/') case(PLASTICITY_DISLOTUNGSTEN_ID) - call plastic_dislotungsten_results(phase_plasticityInstance(ph),group//'plastic/') + call plastic_dislotungsten_results(phase_plasticInstance(ph),group//'plastic/') case(PLASTICITY_NONLOCAL_ID) - call plastic_nonlocal_results(phase_plasticityInstance(ph),group//'plastic/') + call plastic_nonlocal_results(phase_plasticInstance(ph),group//'plastic/') end select call crystallite_results(group,ph) -end subroutine mech_results +end subroutine mechanical_results !-------------------------------------------------------------------------------------------------- @@ -503,8 +503,8 @@ function integrateStress(F,subFp0,subFi0,Delta_t,co,ip,el) result(broken) call plastic_dependentState(co,ip,el) - Lpguess = constitutive_mech_Lp(ph)%data(1:3,1:3,me) ! take as first guess - Liguess = constitutive_mech_Li(ph)%data(1:3,1:3,me) ! take as first guess + Lpguess = phase_mechanical_Lp(ph)%data(1:3,1:3,me) ! take as first guess + Liguess = phase_mechanical_Li(ph)%data(1:3,1:3,me) ! take as first guess call math_invert33(invFp_current,devNull,error,subFp0) if (error) return ! error @@ -538,7 +538,7 @@ function integrateStress(F,subFp0,subFi0,Delta_t,co,ip,el) result(broken) B = math_I3 - Delta_t*Lpguess Fe = matmul(matmul(A,B), invFi_new) - call constitutive_hooke_SandItsTangents(S, dS_dFe, dS_dFi, & + call phase_hooke_SandItsTangents(S, dS_dFe, dS_dFi, & Fe, Fi_new, co, ip, el) call plastic_LpAndItsTangents(Lp_constitutive, dLp_dS, dLp_dFi, & @@ -582,7 +582,7 @@ function integrateStress(F,subFp0,subFi0,Delta_t,co,ip,el) result(broken) + deltaLp * steplengthLp enddo LpLoop - call constitutive_LiAndItsTangents(Li_constitutive, dLi_dS, dLi_dFi, & + call phase_LiAndItsTangents(Li_constitutive, dLi_dS, dLi_dFi, & S, Fi_new, co, ip, el) !* update current residuum and check for convergence of loop @@ -633,13 +633,13 @@ function integrateStress(F,subFp0,subFi0,Delta_t,co,ip,el) result(broken) call math_invert33(Fp_new,devNull,error,invFp_new) if (error) return ! error - constitutive_mech_P(ph)%data(1:3,1:3,me) = matmul(matmul(F,invFp_new),matmul(S,transpose(invFp_new))) - constitutive_mech_S(ph)%data(1:3,1:3,me) = S - constitutive_mech_Lp(ph)%data(1:3,1:3,me) = Lpguess - constitutive_mech_Li(ph)%data(1:3,1:3,me) = Liguess - constitutive_mech_Fp(ph)%data(1:3,1:3,me) = Fp_new / math_det33(Fp_new)**(1.0_pReal/3.0_pReal) ! regularize - constitutive_mech_Fi(ph)%data(1:3,1:3,me) = Fi_new - constitutive_mech_Fe(ph)%data(1:3,1:3,me) = matmul(matmul(F,invFp_new),invFi_new) + phase_mechanical_P(ph)%data(1:3,1:3,me) = matmul(matmul(F,invFp_new),matmul(S,transpose(invFp_new))) + phase_mechanical_S(ph)%data(1:3,1:3,me) = S + phase_mechanical_Lp(ph)%data(1:3,1:3,me) = Lpguess + phase_mechanical_Li(ph)%data(1:3,1:3,me) = Liguess + phase_mechanical_Fp(ph)%data(1:3,1:3,me) = Fp_new / math_det33(Fp_new)**(1.0_pReal/3.0_pReal) ! regularize + phase_mechanical_Fi(ph)%data(1:3,1:3,me) = Fi_new + phase_mechanical_Fe(ph)%data(1:3,1:3,me) = matmul(matmul(F,invFp_new),invFi_new) broken = .false. end function integrateStress @@ -668,9 +668,9 @@ function integrateStateFPI(F_0,F,subFp0,subFi0,subState0,Delta_t,co,ip,el) resul sizeDotState real(pReal) :: & zeta - real(pReal), dimension(constitutive_plasticity_maxSizeDotState) :: & + real(pReal), dimension(phase_plasticity_maxSizeDotState) :: & r ! state residuum - real(pReal), dimension(constitutive_plasticity_maxSizeDotState,2) :: & + real(pReal), dimension(phase_plasticity_maxSizeDotState,2) :: & dotState @@ -796,7 +796,7 @@ function integrateStateAdaptiveEuler(F_0,F,subFp0,subFi0,subState0,Delta_t,co,ip ph, & me, & sizeDotState - real(pReal), dimension(constitutive_plasticity_maxSizeDotState) :: residuum_plastic + real(pReal), dimension(phase_plasticity_maxSizeDotState) :: residuum_plastic ph = material_phaseAt(co,el) @@ -914,7 +914,7 @@ function integrateStateRK(F_0,F,subFp0,subFi0,subState0,Delta_t,co,ip,el,A,B,C,D ph, & me, & sizeDotState - real(pReal), dimension(constitutive_plasticity_maxSizeDotState,size(B)) :: plastic_RKdotState + real(pReal), dimension(phase_plasticity_maxSizeDotState,size(B)) :: plastic_RKdotState ph = material_phaseAt(co,el) @@ -987,28 +987,28 @@ subroutine crystallite_results(group,ph) select case (output_constituent(ph)%label(ou)) case('F') - call results_writeDataset(group//'/mechanics/',constitutive_mech_F(ph)%data,'F',& + call results_writeDataset(group//'/mechanics/',phase_mechanical_F(ph)%data,'F',& 'deformation gradient','1') case('F_e') - call results_writeDataset(group//'/mechanics/',constitutive_mech_Fe(ph)%data,'F_e',& + call results_writeDataset(group//'/mechanics/',phase_mechanical_Fe(ph)%data,'F_e',& 'elastic deformation gradient','1') case('F_p') - call results_writeDataset(group//'/mechanics/',constitutive_mech_Fp(ph)%data,'F_p', & + call results_writeDataset(group//'/mechanics/',phase_mechanical_Fp(ph)%data,'F_p', & 'plastic deformation gradient','1') case('F_i') - call results_writeDataset(group//'/mechanics/',constitutive_mech_Fi(ph)%data,'F_i', & + call results_writeDataset(group//'/mechanics/',phase_mechanical_Fi(ph)%data,'F_i', & 'inelastic deformation gradient','1') case('L_p') - call results_writeDataset(group//'/mechanics/',constitutive_mech_Lp(ph)%data,'L_p', & + call results_writeDataset(group//'/mechanics/',phase_mechanical_Lp(ph)%data,'L_p', & 'plastic velocity gradient','1/s') case('L_i') - call results_writeDataset(group//'/mechanics/',constitutive_mech_Li(ph)%data,'L_i', & + call results_writeDataset(group//'/mechanics/',phase_mechanical_Li(ph)%data,'L_i', & 'inelastic velocity gradient','1/s') case('P') - call results_writeDataset(group//'/mechanics/',constitutive_mech_P(ph)%data,'P', & + call results_writeDataset(group//'/mechanics/',phase_mechanical_P(ph)%data,'P', & 'First Piola-Kirchhoff stress','Pa') case('S') - call results_writeDataset(group//'/mechanics/',constitutive_mech_S(ph)%data,'S', & + call results_writeDataset(group//'/mechanics/',phase_mechanical_S(ph)%data,'S', & 'Second Piola-Kirchhoff stress','Pa') case('O') select case(lattice_structure(ph)) @@ -1067,43 +1067,43 @@ end subroutine crystallite_results !-------------------------------------------------------------------------------------------------- !> @brief Wind homog inc forward. !-------------------------------------------------------------------------------------------------- -module subroutine mech_windForward(ph,me) +module subroutine mechanical_windForward(ph,me) integer, intent(in) :: ph, me - constitutive_mech_Fp0(ph)%data(1:3,1:3,me) = constitutive_mech_Fp(ph)%data(1:3,1:3,me) - constitutive_mech_Fi0(ph)%data(1:3,1:3,me) = constitutive_mech_Fi(ph)%data(1:3,1:3,me) - constitutive_mech_F0(ph)%data(1:3,1:3,me) = constitutive_mech_F(ph)%data(1:3,1:3,me) - constitutive_mech_Li0(ph)%data(1:3,1:3,me) = constitutive_mech_Li(ph)%data(1:3,1:3,me) - constitutive_mech_Lp0(ph)%data(1:3,1:3,me) = constitutive_mech_Lp(ph)%data(1:3,1:3,me) - constitutive_mech_S0(ph)%data(1:3,1:3,me) = constitutive_mech_S(ph)%data(1:3,1:3,me) + phase_mechanical_Fp0(ph)%data(1:3,1:3,me) = phase_mechanical_Fp(ph)%data(1:3,1:3,me) + phase_mechanical_Fi0(ph)%data(1:3,1:3,me) = phase_mechanical_Fi(ph)%data(1:3,1:3,me) + phase_mechanical_F0(ph)%data(1:3,1:3,me) = phase_mechanical_F(ph)%data(1:3,1:3,me) + phase_mechanical_Li0(ph)%data(1:3,1:3,me) = phase_mechanical_Li(ph)%data(1:3,1:3,me) + phase_mechanical_Lp0(ph)%data(1:3,1:3,me) = phase_mechanical_Lp(ph)%data(1:3,1:3,me) + phase_mechanical_S0(ph)%data(1:3,1:3,me) = phase_mechanical_S(ph)%data(1:3,1:3,me) plasticState(ph)%State0(:,me) = plasticState(ph)%state(:,me) -end subroutine mech_windForward +end subroutine mechanical_windForward !-------------------------------------------------------------------------------------------------- !> @brief Forward data after successful increment. ! ToDo: Any guessing for the current states possible? !-------------------------------------------------------------------------------------------------- -module subroutine mech_forward() +module subroutine mechanical_forward() integer :: ph do ph = 1, size(plasticState) - constitutive_mech_Fi0(ph) = constitutive_mech_Fi(ph) - constitutive_mech_Fp0(ph) = constitutive_mech_Fp(ph) - constitutive_mech_F0(ph) = constitutive_mech_F(ph) - constitutive_mech_Li0(ph) = constitutive_mech_Li(ph) - constitutive_mech_Lp0(ph) = constitutive_mech_Lp(ph) - constitutive_mech_S0(ph) = constitutive_mech_S(ph) + phase_mechanical_Fi0(ph) = phase_mechanical_Fi(ph) + phase_mechanical_Fp0(ph) = phase_mechanical_Fp(ph) + phase_mechanical_F0(ph) = phase_mechanical_F(ph) + phase_mechanical_Li0(ph) = phase_mechanical_Li(ph) + phase_mechanical_Lp0(ph) = phase_mechanical_Lp(ph) + phase_mechanical_S0(ph) = phase_mechanical_S(ph) plasticState(ph)%state0 = plasticState(ph)%state enddo -end subroutine mech_forward +end subroutine mechanical_forward @@ -1111,19 +1111,19 @@ end subroutine mech_forward !> @brief returns the homogenize elasticity matrix !> ToDo: homogenizedC66 would be more consistent !-------------------------------------------------------------------------------------------------- -module function constitutive_homogenizedC(ph,me) result(C) +module function phase_homogenizedC(ph,me) result(C) real(pReal), dimension(6,6) :: C integer, intent(in) :: ph, me - plasticityType: select case (phase_plasticity(ph)) - case (PLASTICITY_DISLOTWIN_ID) plasticityType + plasticType: select case (phase_plasticity(ph)) + case (PLASTICITY_DISLOTWIN_ID) plasticType C = plastic_dislotwin_homogenizedC(ph,me) - case default plasticityType + case default plasticType C = lattice_C66(1:6,1:6,ph) - end select plasticityType + end select plasticType -end function constitutive_homogenizedC +end function phase_homogenizedC !-------------------------------------------------------------------------------------------------- @@ -1158,17 +1158,17 @@ module function crystallite_stress(dt,co,ip,el) result(converged_) me = material_phaseMemberAt(co,ip,el) sizeDotState = plasticState(ph)%sizeDotState - subLi0 = constitutive_mech_Li0(ph)%data(1:3,1:3,me) - subLp0 = constitutive_mech_Lp0(ph)%data(1:3,1:3,me) + subLi0 = phase_mechanical_Li0(ph)%data(1:3,1:3,me) + subLp0 = phase_mechanical_Lp0(ph)%data(1:3,1:3,me) subState0 = plasticState(ph)%State0(:,me) do so = 1, phase_Nsources(ph) damageState(ph)%p(so)%subState0(:,me) = damageState(ph)%p(so)%state0(:,me) enddo - subFp0 = constitutive_mech_Fp0(ph)%data(1:3,1:3,me) - subFi0 = constitutive_mech_Fi0(ph)%data(1:3,1:3,me) - subF0 = constitutive_mech_F0(ph)%data(1:3,1:3,me) + subFp0 = phase_mechanical_Fp0(ph)%data(1:3,1:3,me) + subFi0 = phase_mechanical_Fi0(ph)%data(1:3,1:3,me) + subF0 = phase_mechanical_F0(ph)%data(1:3,1:3,me) subFrac = 0.0_pReal subStep = 1.0_pReal/num%subStepSizeCryst todo = .true. @@ -1186,10 +1186,10 @@ module function crystallite_stress(dt,co,ip,el) result(converged_) if (todo) then subF0 = subF - subLp0 = constitutive_mech_Lp(ph)%data(1:3,1:3,me) - subLi0 = constitutive_mech_Li(ph)%data(1:3,1:3,me) - subFp0 = constitutive_mech_Fp(ph)%data(1:3,1:3,me) - subFi0 = constitutive_mech_Fi(ph)%data(1:3,1:3,me) + subLp0 = phase_mechanical_Lp(ph)%data(1:3,1:3,me) + subLi0 = phase_mechanical_Li(ph)%data(1:3,1:3,me) + subFp0 = phase_mechanical_Fp(ph)%data(1:3,1:3,me) + subFi0 = phase_mechanical_Fi(ph)%data(1:3,1:3,me) subState0 = plasticState(ph)%state(:,me) do so = 1, phase_Nsources(ph) damageState(ph)%p(so)%subState0(:,me) = damageState(ph)%p(so)%state(:,me) @@ -1199,12 +1199,12 @@ module function crystallite_stress(dt,co,ip,el) result(converged_) ! cut back (reduced time and restore) else subStep = num%subStepSizeCryst * subStep - constitutive_mech_Fp(ph)%data(1:3,1:3,me) = subFp0 - constitutive_mech_Fi(ph)%data(1:3,1:3,me) = subFi0 - constitutive_mech_S(ph)%data(1:3,1:3,me) = constitutive_mech_S0(ph)%data(1:3,1:3,me) ! why no subS0 ? is S0 of any use? + phase_mechanical_Fp(ph)%data(1:3,1:3,me) = subFp0 + phase_mechanical_Fi(ph)%data(1:3,1:3,me) = subFi0 + phase_mechanical_S(ph)%data(1:3,1:3,me) = phase_mechanical_S0(ph)%data(1:3,1:3,me) ! why no subS0 ? is S0 of any use? if (subStep < 1.0_pReal) then ! actual (not initial) cutback - constitutive_mech_Lp(ph)%data(1:3,1:3,me) = subLp0 - constitutive_mech_Li(ph)%data(1:3,1:3,me) = subLi0 + phase_mechanical_Lp(ph)%data(1:3,1:3,me) = subLp0 + phase_mechanical_Li(ph)%data(1:3,1:3,me) = subLi0 endif plasticState(ph)%state(:,me) = subState0 do so = 1, phase_Nsources(ph) @@ -1218,9 +1218,9 @@ module function crystallite_stress(dt,co,ip,el) result(converged_) ! prepare for integration if (todo) then subF = subF0 & - + subStep * (constitutive_mech_F(ph)%data(1:3,1:3,me) - constitutive_mech_F0(ph)%data(1:3,1:3,me)) - constitutive_mech_Fe(ph)%data(1:3,1:3,me) = matmul(subF,math_inv33(matmul(constitutive_mech_Fi(ph)%data(1:3,1:3,me), & - constitutive_mech_Fp(ph)%data(1:3,1:3,me)))) + + subStep * (phase_mechanical_F(ph)%data(1:3,1:3,me) - phase_mechanical_F0(ph)%data(1:3,1:3,me)) + phase_mechanical_Fe(ph)%data(1:3,1:3,me) = matmul(subF,math_inv33(matmul(phase_mechanical_Fi(ph)%data(1:3,1:3,me), & + phase_mechanical_Fp(ph)%data(1:3,1:3,me)))) converged_ = .not. integrateState(subF0,subF,subFp0,subFi0,subState0(1:sizeDotState),subStep * dt,co,ip,el) converged_ = converged_ .and. .not. integrateDamageState(subStep * dt,co,ip,el) endif @@ -1233,7 +1233,7 @@ end function crystallite_stress !-------------------------------------------------------------------------------------------------- !> @brief Restore data after homog cutback. !-------------------------------------------------------------------------------------------------- -module subroutine mech_restore(ce,includeL) +module subroutine mechanical_restore(ce,includeL) integer, intent(in) :: ce logical, intent(in) :: & @@ -1247,23 +1247,23 @@ module subroutine mech_restore(ce,includeL) ph = material_phaseAt2(co,ce) me = material_phaseMemberAt2(co,ce) if (includeL) then - constitutive_mech_Lp(ph)%data(1:3,1:3,me) = constitutive_mech_Lp0(ph)%data(1:3,1:3,me) - constitutive_mech_Li(ph)%data(1:3,1:3,me) = constitutive_mech_Li0(ph)%data(1:3,1:3,me) + phase_mechanical_Lp(ph)%data(1:3,1:3,me) = phase_mechanical_Lp0(ph)%data(1:3,1:3,me) + phase_mechanical_Li(ph)%data(1:3,1:3,me) = phase_mechanical_Li0(ph)%data(1:3,1:3,me) endif ! maybe protecting everything from overwriting makes more sense - constitutive_mech_Fp(ph)%data(1:3,1:3,me) = constitutive_mech_Fp0(ph)%data(1:3,1:3,me) - constitutive_mech_Fi(ph)%data(1:3,1:3,me) = constitutive_mech_Fi0(ph)%data(1:3,1:3,me) - constitutive_mech_S(ph)%data(1:3,1:3,me) = constitutive_mech_S0(ph)%data(1:3,1:3,me) + phase_mechanical_Fp(ph)%data(1:3,1:3,me) = phase_mechanical_Fp0(ph)%data(1:3,1:3,me) + phase_mechanical_Fi(ph)%data(1:3,1:3,me) = phase_mechanical_Fi0(ph)%data(1:3,1:3,me) + phase_mechanical_S(ph)%data(1:3,1:3,me) = phase_mechanical_S0(ph)%data(1:3,1:3,me) plasticState(ph)%state(:,me) = plasticState(ph)%State0(:,me) enddo -end subroutine mech_restore +end subroutine mechanical_restore !-------------------------------------------------------------------------------------------------- !> @brief Calculate tangent (dPdF). !-------------------------------------------------------------------------------------------------- -module function constitutive_mech_dPdF(dt,co,ip,el) result(dPdF) +module function phase_mechanical_dPdF(dt,co,ip,el) result(dPdF) real(pReal), intent(in) :: dt integer, intent(in) :: & @@ -1297,18 +1297,18 @@ module function constitutive_mech_dPdF(dt,co,ip,el) result(dPdF) ph = material_phaseAt(co,el) me = material_phaseMemberAt(co,ip,el) - call constitutive_hooke_SandItsTangents(devNull,dSdFe,dSdFi, & - constitutive_mech_Fe(ph)%data(1:3,1:3,me), & - constitutive_mech_Fi(ph)%data(1:3,1:3,me),co,ip,el) - call constitutive_LiAndItsTangents(devNull,dLidS,dLidFi, & - constitutive_mech_S(ph)%data(1:3,1:3,me), & - constitutive_mech_Fi(ph)%data(1:3,1:3,me), & + call phase_hooke_SandItsTangents(devNull,dSdFe,dSdFi, & + phase_mechanical_Fe(ph)%data(1:3,1:3,me), & + phase_mechanical_Fi(ph)%data(1:3,1:3,me),co,ip,el) + call phase_LiAndItsTangents(devNull,dLidS,dLidFi, & + phase_mechanical_S(ph)%data(1:3,1:3,me), & + phase_mechanical_Fi(ph)%data(1:3,1:3,me), & co,ip,el) - invFp = math_inv33(constitutive_mech_Fp(ph)%data(1:3,1:3,me)) - invFi = math_inv33(constitutive_mech_Fi(ph)%data(1:3,1:3,me)) - invSubFp0 = math_inv33(constitutive_mech_Fp0(ph)%data(1:3,1:3,me)) - invSubFi0 = math_inv33(constitutive_mech_Fi0(ph)%data(1:3,1:3,me)) + invFp = math_inv33(phase_mechanical_Fp(ph)%data(1:3,1:3,me)) + invFi = math_inv33(phase_mechanical_Fi(ph)%data(1:3,1:3,me)) + invSubFp0 = math_inv33(phase_mechanical_Fp0(ph)%data(1:3,1:3,me)) + invSubFi0 = math_inv33(phase_mechanical_Fi0(ph)%data(1:3,1:3,me)) if (sum(abs(dLidS)) < tol_math_check) then dFidS = 0.0_pReal @@ -1334,15 +1334,15 @@ module function constitutive_mech_dPdF(dt,co,ip,el) result(dPdF) endif call plastic_LpAndItsTangents(devNull,dLpdS,dLpdFi, & - constitutive_mech_S(ph)%data(1:3,1:3,me), & - constitutive_mech_Fi(ph)%data(1:3,1:3,me),co,ip,el) + phase_mechanical_S(ph)%data(1:3,1:3,me), & + phase_mechanical_Fi(ph)%data(1:3,1:3,me),co,ip,el) dLpdS = math_mul3333xx3333(dLpdFi,dFidS) + dLpdS !-------------------------------------------------------------------------------------------------- ! calculate dSdF temp_33_1 = transpose(matmul(invFp,invFi)) - temp_33_2 = matmul(constitutive_mech_F(ph)%data(1:3,1:3,me),invSubFp0) - temp_33_3 = matmul(matmul(constitutive_mech_F(ph)%data(1:3,1:3,me),invFp), invSubFi0) + temp_33_2 = matmul(phase_mechanical_F(ph)%data(1:3,1:3,me),invSubFp0) + temp_33_3 = matmul(matmul(phase_mechanical_F(ph)%data(1:3,1:3,me),invFp), invSubFi0) do o=1,3; do p=1,3 rhs_3333(p,o,1:3,1:3) = matmul(dSdFe(p,o,1:3,1:3),temp_33_1) @@ -1370,9 +1370,9 @@ module function constitutive_mech_dPdF(dt,co,ip,el) result(dPdF) !-------------------------------------------------------------------------------------------------- ! assemble dPdF - temp_33_1 = matmul(constitutive_mech_S(ph)%data(1:3,1:3,me),transpose(invFp)) - temp_33_2 = matmul(constitutive_mech_F(ph)%data(1:3,1:3,me),invFp) - temp_33_3 = matmul(temp_33_2,constitutive_mech_S(ph)%data(1:3,1:3,me)) + temp_33_1 = matmul(phase_mechanical_S(ph)%data(1:3,1:3,me),transpose(invFp)) + temp_33_2 = matmul(phase_mechanical_F(ph)%data(1:3,1:3,me),invFp) + temp_33_3 = matmul(temp_33_2,phase_mechanical_S(ph)%data(1:3,1:3,me)) dPdF = 0.0_pReal do p=1,3 @@ -1380,129 +1380,129 @@ module function constitutive_mech_dPdF(dt,co,ip,el) result(dPdF) enddo do o=1,3; do p=1,3 dPdF(1:3,1:3,p,o) = dPdF(1:3,1:3,p,o) & - + matmul(matmul(constitutive_mech_F(ph)%data(1:3,1:3,me),dFpinvdF(1:3,1:3,p,o)),temp_33_1) & + + matmul(matmul(phase_mechanical_F(ph)%data(1:3,1:3,me),dFpinvdF(1:3,1:3,p,o)),temp_33_1) & + matmul(matmul(temp_33_2,dSdF(1:3,1:3,p,o)),transpose(invFp)) & + matmul(temp_33_3,transpose(dFpinvdF(1:3,1:3,p,o))) enddo; enddo -end function constitutive_mech_dPdF +end function phase_mechanical_dPdF -module subroutine mech_restartWrite(groupHandle,ph) +module subroutine mechanical_restartWrite(groupHandle,ph) integer(HID_T), intent(in) :: groupHandle integer, intent(in) :: ph call HDF5_write(groupHandle,plasticState(ph)%state,'omega') - call HDF5_write(groupHandle,constitutive_mech_Fi(ph)%data,'F_i') - call HDF5_write(groupHandle,constitutive_mech_Li(ph)%data,'L_i') - call HDF5_write(groupHandle,constitutive_mech_Lp(ph)%data,'L_p') - call HDF5_write(groupHandle,constitutive_mech_Fp(ph)%data,'F_p') - call HDF5_write(groupHandle,constitutive_mech_S(ph)%data,'S') - call HDF5_write(groupHandle,constitutive_mech_F(ph)%data,'F') + call HDF5_write(groupHandle,phase_mechanical_Fi(ph)%data,'F_i') + call HDF5_write(groupHandle,phase_mechanical_Li(ph)%data,'L_i') + call HDF5_write(groupHandle,phase_mechanical_Lp(ph)%data,'L_p') + call HDF5_write(groupHandle,phase_mechanical_Fp(ph)%data,'F_p') + call HDF5_write(groupHandle,phase_mechanical_S(ph)%data,'S') + call HDF5_write(groupHandle,phase_mechanical_F(ph)%data,'F') -end subroutine mech_restartWrite +end subroutine mechanical_restartWrite -module subroutine mech_restartRead(groupHandle,ph) +module subroutine mechanical_restartRead(groupHandle,ph) integer(HID_T), intent(in) :: groupHandle integer, intent(in) :: ph call HDF5_read(groupHandle,plasticState(ph)%state0,'omega') - call HDF5_read(groupHandle,constitutive_mech_Fi0(ph)%data,'F_i') - call HDF5_read(groupHandle,constitutive_mech_Li0(ph)%data,'L_i') - call HDF5_read(groupHandle,constitutive_mech_Lp0(ph)%data,'L_p') - call HDF5_read(groupHandle,constitutive_mech_Fp0(ph)%data,'F_p') - call HDF5_read(groupHandle,constitutive_mech_S0(ph)%data,'S') - call HDF5_read(groupHandle,constitutive_mech_F0(ph)%data,'F') + call HDF5_read(groupHandle,phase_mechanical_Fi0(ph)%data,'F_i') + call HDF5_read(groupHandle,phase_mechanical_Li0(ph)%data,'L_i') + call HDF5_read(groupHandle,phase_mechanical_Lp0(ph)%data,'L_p') + call HDF5_read(groupHandle,phase_mechanical_Fp0(ph)%data,'F_p') + call HDF5_read(groupHandle,phase_mechanical_S0(ph)%data,'S') + call HDF5_read(groupHandle,phase_mechanical_F0(ph)%data,'F') -end subroutine mech_restartRead +end subroutine mechanical_restartRead !---------------------------------------------------------------------------------------------- !< @brief Get first Piola-Kichhoff stress (for use by non-mech physics) !---------------------------------------------------------------------------------------------- -module function mech_S(ph,me) result(S) +module function mechanical_S(ph,me) result(S) integer, intent(in) :: ph,me real(pReal), dimension(3,3) :: S - S = constitutive_mech_S(ph)%data(1:3,1:3,me) + S = phase_mechanical_S(ph)%data(1:3,1:3,me) -end function mech_S +end function mechanical_S !---------------------------------------------------------------------------------------------- !< @brief Get plastic velocity gradient (for use by non-mech physics) !---------------------------------------------------------------------------------------------- -module function mech_L_p(ph,me) result(L_p) +module function mechanical_L_p(ph,me) result(L_p) integer, intent(in) :: ph,me real(pReal), dimension(3,3) :: L_p - L_p = constitutive_mech_Lp(ph)%data(1:3,1:3,me) + L_p = phase_mechanical_Lp(ph)%data(1:3,1:3,me) -end function mech_L_p +end function mechanical_L_p !---------------------------------------------------------------------------------------------- !< @brief Get deformation gradient (for use by homogenization) !---------------------------------------------------------------------------------------------- -module function constitutive_mech_getF(co,ip,el) result(F) +module function phase_mechanical_getF(co,ip,el) result(F) integer, intent(in) :: co, ip, el real(pReal), dimension(3,3) :: F - F = constitutive_mech_F(material_phaseAt(co,el))%data(1:3,1:3,material_phaseMemberAt(co,ip,el)) + F = phase_mechanical_F(material_phaseAt(co,el))%data(1:3,1:3,material_phaseMemberAt(co,ip,el)) -end function constitutive_mech_getF +end function phase_mechanical_getF !---------------------------------------------------------------------------------------------- !< @brief Get elastic deformation gradient (for use by non-mech physics) !---------------------------------------------------------------------------------------------- -module function mech_F_e(ph,me) result(F_e) +module function mechanical_F_e(ph,me) result(F_e) integer, intent(in) :: ph,me real(pReal), dimension(3,3) :: F_e - F_e = constitutive_mech_Fe(ph)%data(1:3,1:3,me) + F_e = phase_mechanical_Fe(ph)%data(1:3,1:3,me) -end function mech_F_e +end function mechanical_F_e !---------------------------------------------------------------------------------------------- !< @brief Get second Piola-Kichhoff stress (for use by homogenization) !---------------------------------------------------------------------------------------------- -module function constitutive_mech_getP(co,ip,el) result(P) +module function phase_mechanical_getP(co,ip,el) result(P) integer, intent(in) :: co, ip, el real(pReal), dimension(3,3) :: P - P = constitutive_mech_P(material_phaseAt(co,el))%data(1:3,1:3,material_phaseMemberAt(co,ip,el)) + P = phase_mechanical_P(material_phaseAt(co,el))%data(1:3,1:3,material_phaseMemberAt(co,ip,el)) -end function constitutive_mech_getP +end function phase_mechanical_getP ! setter for homogenization -module subroutine constitutive_mech_setF(F,co,ip,el) +module subroutine phase_mechanical_setF(F,co,ip,el) real(pReal), dimension(3,3), intent(in) :: F integer, intent(in) :: co, ip, el - constitutive_mech_F(material_phaseAt(co,el))%data(1:3,1:3,material_phaseMemberAt(co,ip,el)) = F + phase_mechanical_F(material_phaseAt(co,el))%data(1:3,1:3,material_phaseMemberAt(co,ip,el)) = F -end subroutine constitutive_mech_setF +end subroutine phase_mechanical_setF end submodule mechanics diff --git a/src/phase_mechanics_eigendeformation.f90 b/src/phase_mechanical_eigen.f90 similarity index 96% rename from src/phase_mechanics_eigendeformation.f90 rename to src/phase_mechanical_eigen.f90 index 45cfd82d3..6df993565 100644 --- a/src/phase_mechanics_eigendeformation.f90 +++ b/src/phase_mechanical_eigen.f90 @@ -126,7 +126,7 @@ end function kinematics_active !> @brief contains the constitutive equation for calculating the velocity gradient ! ToDo: MD: S is Mi? !-------------------------------------------------------------------------------------------------- -module subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, & +module subroutine phase_LiAndItsTangents(Li, dLi_dS, dLi_dFi, & S, Fi, co, ip, el) integer, intent(in) :: & @@ -159,15 +159,15 @@ module subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, & dLi_dS = 0.0_pReal dLi_dFi = 0.0_pReal - plasticityType: select case (phase_plasticity(material_phaseAt(co,el))) - case (PLASTICITY_isotropic_ID) plasticityType + plasticType: select case (phase_plasticity(material_phaseAt(co,el))) + case (PLASTICITY_isotropic_ID) plasticType of = material_phasememberAt(co,ip,el) - instance = phase_plasticityInstance(material_phaseAt(co,el)) + instance = phase_plasticInstance(material_phaseAt(co,el)) call plastic_isotropic_LiAndItsTangent(my_Li, my_dLi_dS, S ,instance,of) - case default plasticityType + case default plasticType my_Li = 0.0_pReal my_dLi_dS = 0.0_pReal - end select plasticityType + end select plasticType Li = Li + my_Li dLi_dS = dLi_dS + my_dLi_dS @@ -201,7 +201,7 @@ module subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, & dLi_dFi(1:3,i,1:3,j) = dLi_dFi(1:3,i,1:3,j) + math_I3*temp_33(j,i) + Li*FiInv(j,i) enddo; enddo -end subroutine constitutive_LiAndItsTangents +end subroutine phase_LiAndItsTangents end submodule eigendeformation diff --git a/src/phase_mechanics_eigendeformation_cleavageopening.f90 b/src/phase_mechanical_eigen_cleavageopening.f90 similarity index 100% rename from src/phase_mechanics_eigendeformation_cleavageopening.f90 rename to src/phase_mechanical_eigen_cleavageopening.f90 diff --git a/src/phase_mechanics_eigendeformation_slipplaneopening.f90 b/src/phase_mechanical_eigen_slipplaneopening.f90 similarity index 100% rename from src/phase_mechanics_eigendeformation_slipplaneopening.f90 rename to src/phase_mechanical_eigen_slipplaneopening.f90 diff --git a/src/phase_mechanics_eigendeformation_thermalexpansion.f90 b/src/phase_mechanical_eigen_thermalexpansion.f90 similarity index 100% rename from src/phase_mechanics_eigendeformation_thermalexpansion.f90 rename to src/phase_mechanical_eigen_thermalexpansion.f90 diff --git a/src/phase_mechanics_plastic.f90 b/src/phase_mechanical_plastic.f90 similarity index 91% rename from src/phase_mechanics_plastic.f90 rename to src/phase_mechanical_plastic.f90 index fff1fbdcf..ed1dc64fb 100644 --- a/src/phase_mechanics_plastic.f90 +++ b/src/phase_mechanical_plastic.f90 @@ -270,31 +270,31 @@ module subroutine plastic_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, & me = material_phasememberAt(co,ip,el) ph = material_phaseAt(co,el) - plasticityType: select case (phase_plasticity(material_phaseAt(co,el))) + plasticType: select case (phase_plasticity(material_phaseAt(co,el))) - case (PLASTICITY_NONE_ID) plasticityType + case (PLASTICITY_NONE_ID) plasticType Lp = 0.0_pReal dLp_dMp = 0.0_pReal - case (PLASTICITY_ISOTROPIC_ID) plasticityType + case (PLASTICITY_ISOTROPIC_ID) plasticType call isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,ph,me) - case (PLASTICITY_PHENOPOWERLAW_ID) plasticityType + case (PLASTICITY_PHENOPOWERLAW_ID) plasticType call phenopowerlaw_LpAndItsTangent(Lp,dLp_dMp,Mp,ph,me) - case (PLASTICITY_KINEHARDENING_ID) plasticityType + case (PLASTICITY_KINEHARDENING_ID) plasticType call kinehardening_LpAndItsTangent(Lp,dLp_dMp,Mp,ph,me) - case (PLASTICITY_NONLOCAL_ID) plasticityType + case (PLASTICITY_NONLOCAL_ID) plasticType call nonlocal_LpAndItsTangent(Lp,dLp_dMp,Mp, thermal_T(ph,me),ph,me,ip,el) - case (PLASTICITY_DISLOTWIN_ID) plasticityType + case (PLASTICITY_DISLOTWIN_ID) plasticType call dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp, thermal_T(ph,me),ph,me) - case (PLASTICITY_DISLOTUNGSTEN_ID) plasticityType + case (PLASTICITY_DISLOTUNGSTEN_ID) plasticType call dislotungsten_LpAndItsTangent(Lp,dLp_dMp,Mp, thermal_T(ph,me),ph,me) - end select plasticityType + end select plasticType do i=1,3; do j=1,3 dLp_dFi(i,j,1:3,1:3) = matmul(matmul(Fi,S),transpose(dLp_dMp(i,j,1:3,1:3))) + & @@ -323,29 +323,29 @@ module function plastic_dotState(subdt,co,ip,el,ph,me) result(broken) logical :: broken - Mp = matmul(matmul(transpose(constitutive_mech_Fi(ph)%data(1:3,1:3,me)),& - constitutive_mech_Fi(ph)%data(1:3,1:3,me)),constitutive_mech_S(ph)%data(1:3,1:3,me)) + Mp = matmul(matmul(transpose(phase_mechanical_Fi(ph)%data(1:3,1:3,me)),& + phase_mechanical_Fi(ph)%data(1:3,1:3,me)),phase_mechanical_S(ph)%data(1:3,1:3,me)) - plasticityType: select case (phase_plasticity(ph)) + plasticType: select case (phase_plasticity(ph)) - case (PLASTICITY_ISOTROPIC_ID) plasticityType + case (PLASTICITY_ISOTROPIC_ID) plasticType call isotropic_dotState(Mp,ph,me) - case (PLASTICITY_PHENOPOWERLAW_ID) plasticityType + case (PLASTICITY_PHENOPOWERLAW_ID) plasticType call phenopowerlaw_dotState(Mp,ph,me) - case (PLASTICITY_KINEHARDENING_ID) plasticityType + case (PLASTICITY_KINEHARDENING_ID) plasticType call plastic_kinehardening_dotState(Mp,ph,me) - case (PLASTICITY_DISLOTWIN_ID) plasticityType + case (PLASTICITY_DISLOTWIN_ID) plasticType call dislotwin_dotState(Mp,thermal_T(ph,me),ph,me) - case (PLASTICITY_DISLOTUNGSTEN_ID) plasticityType + case (PLASTICITY_DISLOTUNGSTEN_ID) plasticType call dislotungsten_dotState(Mp,thermal_T(ph,me),ph,me) - case (PLASTICITY_NONLOCAL_ID) plasticityType + case (PLASTICITY_NONLOCAL_ID) plasticType call nonlocal_dotState(Mp,thermal_T(ph,me),subdt,ph,me,ip,el) - end select plasticityType + end select plasticType broken = any(IEEE_is_NaN(plasticState(ph)%dotState(:,me))) @@ -369,20 +369,20 @@ module subroutine plastic_dependentState(co, ip, el) ph = material_phaseAt(co,el) me = material_phasememberAt(co,ip,el) - instance = phase_plasticityInstance(ph) + instance = phase_plasticInstance(ph) - plasticityType: select case (phase_plasticity(material_phaseAt(co,el))) + plasticType: select case (phase_plasticity(material_phaseAt(co,el))) - case (PLASTICITY_DISLOTWIN_ID) plasticityType + case (PLASTICITY_DISLOTWIN_ID) plasticType call dislotwin_dependentState(thermal_T(ph,me),instance,me) - case (PLASTICITY_DISLOTUNGSTEN_ID) plasticityType + case (PLASTICITY_DISLOTUNGSTEN_ID) plasticType call dislotungsten_dependentState(instance,me) - case (PLASTICITY_NONLOCAL_ID) plasticityType + case (PLASTICITY_NONLOCAL_ID) plasticType call nonlocal_dependentState(instance,me,ip,el) - end select plasticityType + end select plasticType end subroutine plastic_dependentState @@ -410,24 +410,24 @@ module function plastic_deltaState(co, ip, el, ph, me) result(broken) mySize - Mp = matmul(matmul(transpose(constitutive_mech_Fi(ph)%data(1:3,1:3,me)),& - constitutive_mech_Fi(ph)%data(1:3,1:3,me)),constitutive_mech_S(ph)%data(1:3,1:3,me)) - instance = phase_plasticityInstance(ph) + Mp = matmul(matmul(transpose(phase_mechanical_Fi(ph)%data(1:3,1:3,me)),& + phase_mechanical_Fi(ph)%data(1:3,1:3,me)),phase_mechanical_S(ph)%data(1:3,1:3,me)) + instance = phase_plasticInstance(ph) - plasticityType: select case (phase_plasticity(ph)) + plasticType: select case (phase_plasticity(ph)) - case (PLASTICITY_KINEHARDENING_ID) plasticityType + case (PLASTICITY_KINEHARDENING_ID) plasticType call plastic_kinehardening_deltaState(Mp,instance,me) broken = any(IEEE_is_NaN(plasticState(ph)%deltaState(:,me))) - case (PLASTICITY_NONLOCAL_ID) plasticityType + case (PLASTICITY_NONLOCAL_ID) plasticType call plastic_nonlocal_deltaState(Mp,instance,me,ip,el) broken = any(IEEE_is_NaN(plasticState(ph)%deltaState(:,me))) case default broken = .false. - end select plasticityType + end select plasticType if(.not. broken) then select case(phase_plasticity(ph)) diff --git a/src/phase_mechanics_plastic_dislotungsten.f90 b/src/phase_mechanical_plastic_dislotungsten.f90 similarity index 97% rename from src/phase_mechanics_plastic_dislotungsten.f90 rename to src/phase_mechanical_plastic_dislotungsten.f90 index a47132c63..a9946b6e7 100644 --- a/src/phase_mechanics_plastic_dislotungsten.f90 +++ b/src/phase_mechanical_plastic_dislotungsten.f90 @@ -226,7 +226,7 @@ module function plastic_dislotungsten_init() result(myPlasticity) sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) + call phase_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization @@ -289,16 +289,16 @@ pure module subroutine dislotungsten_LpAndItsTangent(Lp,dLp_dMp, & integer :: & i,k,l,m,n - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & dot_gamma_pos,dot_gamma_neg, & ddot_gamma_dtau_pos,ddot_gamma_dtau_neg Lp = 0.0_pReal dLp_dMp = 0.0_pReal - associate(prm => param(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph))) - call kinetics(Mp,T,phase_plasticityInstance(ph),me,dot_gamma_pos,dot_gamma_neg,ddot_gamma_dtau_pos,ddot_gamma_dtau_neg) + call kinetics(Mp,T,phase_plasticInstance(ph),me,dot_gamma_pos,dot_gamma_neg,ddot_gamma_dtau_pos,ddot_gamma_dtau_neg) do i = 1, prm%sum_N_sl Lp = Lp + (dot_gamma_pos(i)+dot_gamma_neg(i))*prm%P_sl(1:3,1:3,i) forall (k=1:3,l=1:3,m=1:3,n=1:3) & @@ -327,7 +327,7 @@ module subroutine dislotungsten_dotState(Mp,T,ph,me) real(pReal) :: & VacancyDiffusion - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & gdot_pos, gdot_neg,& tau_pos,& tau_neg, & @@ -336,10 +336,10 @@ module subroutine dislotungsten_dotState(Mp,T,ph,me) dot_rho_dip_climb, & dip_distance - associate(prm => param(phase_plasticityInstance(ph)), stt => state(phase_plasticityInstance(ph)),& - dot => dotState(phase_plasticityInstance(ph)), dst => dependentState(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph)),& + dot => dotState(phase_plasticInstance(ph)), dst => dependentState(phase_plasticInstance(ph))) - call kinetics(Mp,T,phase_plasticityInstance(ph),me,& + call kinetics(Mp,T,phase_plasticInstance(ph),me,& gdot_pos,gdot_neg, & tau_pos_out = tau_pos,tau_neg_out = tau_neg) diff --git a/src/phase_mechanics_plastic_dislotwin.f90 b/src/phase_mechanical_plastic_dislotwin.f90 similarity index 97% rename from src/phase_mechanics_plastic_dislotwin.f90 rename to src/phase_mechanical_plastic_dislotwin.f90 index 9f7464323..e1259fbd2 100644 --- a/src/phase_mechanics_plastic_dislotwin.f90 +++ b/src/phase_mechanical_plastic_dislotwin.f90 @@ -415,7 +415,7 @@ module function plastic_dislotwin_init() result(myPlasticity) sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) + call phase_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! locally defined state aliases and initialization of state0 and atol @@ -496,8 +496,8 @@ module function plastic_dislotwin_homogenizedC(ph,me) result(homogenizedC) real(pReal) :: f_unrotated - associate(prm => param(phase_plasticityInstance(ph)),& - stt => state(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)),& + stt => state(phase_plasticInstance(ph))) f_unrotated = 1.0_pReal & - sum(stt%f_tw(1:prm%sum_N_tw,me)) & @@ -535,11 +535,11 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me) BoltzmannRatio, & ddot_gamma_dtau, & tau - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & dot_gamma_sl,ddot_gamma_dtau_slip - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_tw) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_tw) :: & dot_gamma_twin,ddot_gamma_dtau_twin - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_tr) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_tr) :: & dot_gamma_tr,ddot_gamma_dtau_trans real(pReal):: dot_gamma_sb real(pReal), dimension(3,3) :: eigVectors, P_sb @@ -564,7 +564,7 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me) 0, 1, 1 & ],pReal),[ 3,6]) - associate(prm => param(phase_plasticityInstance(ph)), stt => state(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph))) f_unrotated = 1.0_pReal & - sum(stt%f_tw(1:prm%sum_N_tw,me)) & @@ -573,7 +573,7 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me) Lp = 0.0_pReal dLp_dMp = 0.0_pReal - call kinetics_slip(Mp,T,phase_plasticityInstance(ph),me,dot_gamma_sl,ddot_gamma_dtau_slip) + call kinetics_slip(Mp,T,phase_plasticInstance(ph),me,dot_gamma_sl,ddot_gamma_dtau_slip) slipContribution: do i = 1, prm%sum_N_sl Lp = Lp + dot_gamma_sl(i)*prm%P_sl(1:3,1:3,i) forall (k=1:3,l=1:3,m=1:3,n=1:3) & @@ -581,7 +581,7 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me) + ddot_gamma_dtau_slip(i) * prm%P_sl(k,l,i) * prm%P_sl(m,n,i) enddo slipContribution - call kinetics_twin(Mp,T,dot_gamma_sl,phase_plasticityInstance(ph),me,dot_gamma_twin,ddot_gamma_dtau_twin) + call kinetics_twin(Mp,T,dot_gamma_sl,phase_plasticInstance(ph),me,dot_gamma_twin,ddot_gamma_dtau_twin) twinContibution: do i = 1, prm%sum_N_tw Lp = Lp + dot_gamma_twin(i)*prm%P_tw(1:3,1:3,i) forall (k=1:3,l=1:3,m=1:3,n=1:3) & @@ -589,7 +589,7 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me) + ddot_gamma_dtau_twin(i)* prm%P_tw(k,l,i)*prm%P_tw(m,n,i) enddo twinContibution - call kinetics_trans(Mp,T,dot_gamma_sl,phase_plasticityInstance(ph),me,dot_gamma_tr,ddot_gamma_dtau_trans) + call kinetics_trans(Mp,T,dot_gamma_sl,phase_plasticInstance(ph),me,dot_gamma_tr,ddot_gamma_dtau_trans) transContibution: do i = 1, prm%sum_N_tr Lp = Lp + dot_gamma_tr(i)*prm%P_tr(1:3,1:3,i) forall (k=1:3,l=1:3,m=1:3,n=1:3) & @@ -653,24 +653,24 @@ module subroutine dislotwin_dotState(Mp,T,ph,me) tau, & sigma_cl, & !< climb stress b_d !< ratio of Burgers vector to stacking fault width - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & dot_rho_dip_formation, & dot_rho_dip_climb, & rho_dip_distance_min, & dot_gamma_sl - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_tw) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_tw) :: & dot_gamma_twin - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_tr) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_tr) :: & dot_gamma_tr - associate(prm => param(phase_plasticityInstance(ph)), stt => state(phase_plasticityInstance(ph)), & - dot => dotState(phase_plasticityInstance(ph)), dst => dependentState(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph)), & + dot => dotState(phase_plasticInstance(ph)), dst => dependentState(phase_plasticInstance(ph))) f_unrotated = 1.0_pReal & - sum(stt%f_tw(1:prm%sum_N_tw,me)) & - sum(stt%f_tr(1:prm%sum_N_tr,me)) - call kinetics_slip(Mp,T,phase_plasticityInstance(ph),me,dot_gamma_sl) + call kinetics_slip(Mp,T,phase_plasticInstance(ph),me,dot_gamma_sl) dot%gamma_sl(:,me) = abs(dot_gamma_sl) rho_dip_distance_min = prm%D_a*prm%b_sl @@ -721,10 +721,10 @@ module subroutine dislotwin_dotState(Mp,T,ph,me) - 2.0_pReal*rho_dip_distance_min/prm%b_sl * stt%rho_dip(:,me)*abs(dot_gamma_sl) & - dot_rho_dip_climb - call kinetics_twin(Mp,T,dot_gamma_sl,phase_plasticityInstance(ph),me,dot_gamma_twin) + call kinetics_twin(Mp,T,dot_gamma_sl,phase_plasticInstance(ph),me,dot_gamma_twin) dot%f_tw(:,me) = f_unrotated*dot_gamma_twin/prm%gamma_char - call kinetics_trans(Mp,T,dot_gamma_sl,phase_plasticityInstance(ph),me,dot_gamma_tr) + call kinetics_trans(Mp,T,dot_gamma_sl,phase_plasticInstance(ph),me,dot_gamma_tr) dot%f_tr(:,me) = f_unrotated*dot_gamma_tr end associate diff --git a/src/phase_mechanics_plastic_isotropic.f90 b/src/phase_mechanical_plastic_isotropic.f90 similarity index 97% rename from src/phase_mechanics_plastic_isotropic.f90 rename to src/phase_mechanical_plastic_isotropic.f90 index 6789e74b4..1a4e1449a 100644 --- a/src/phase_mechanics_plastic_isotropic.f90 +++ b/src/phase_mechanical_plastic_isotropic.f90 @@ -135,7 +135,7 @@ module function plastic_isotropic_init() result(myPlasticity) sizeDotState = size(['xi ','gamma']) sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) + call phase_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization @@ -190,7 +190,7 @@ module subroutine isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,ph,me) integer :: & k, l, m, n - associate(prm => param(phase_plasticityInstance(ph)), stt => state(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph))) Mp_dev = math_deviatoric33(Mp) squarenorm_Mp_dev = math_tensordot(Mp_dev,Mp_dev) @@ -275,8 +275,8 @@ module subroutine isotropic_dotState(Mp,ph,me) xi_inf_star, & !< saturation xi norm_Mp !< norm of the (deviatoric) Mandel stress - associate(prm => param(phase_plasticityInstance(ph)), stt => state(phase_plasticityInstance(ph)), & - dot => dotState(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph)), & + dot => dotState(phase_plasticInstance(ph))) if (prm%dilatation) then norm_Mp = sqrt(math_tensordot(Mp,Mp)) diff --git a/src/phase_mechanics_plastic_kinehardening.f90 b/src/phase_mechanical_plastic_kinehardening.f90 similarity index 97% rename from src/phase_mechanics_plastic_kinehardening.f90 rename to src/phase_mechanical_plastic_kinehardening.f90 index 919302bc1..f00e4171b 100644 --- a/src/phase_mechanics_plastic_kinehardening.f90 +++ b/src/phase_mechanical_plastic_kinehardening.f90 @@ -180,7 +180,7 @@ module function plastic_kinehardening_init() result(myPlasticity) sizeDeltaState = size(['sense ', 'chi0 ', 'gamma0' ]) * prm%sum_N_sl !ToDo: adjust names sizeState = sizeDotState + sizeDeltaState - call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,sizeDeltaState) + call phase_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,sizeDeltaState) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization @@ -255,16 +255,16 @@ pure module subroutine kinehardening_LpAndItsTangent(Lp,dLp_dMp,Mp,ph,me) integer :: & i,k,l,m,n - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & gdot_pos,gdot_neg, & dgdot_dtau_pos,dgdot_dtau_neg Lp = 0.0_pReal dLp_dMp = 0.0_pReal - associate(prm => param(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph))) - call kinetics(Mp,phase_plasticityInstance(ph),me,gdot_pos,gdot_neg,dgdot_dtau_pos,dgdot_dtau_neg) + call kinetics(Mp,phase_plasticInstance(ph),me,gdot_pos,gdot_neg,dgdot_dtau_pos,dgdot_dtau_neg) do i = 1, prm%sum_N_sl Lp = Lp + (gdot_pos(i)+gdot_neg(i))*prm%P(1:3,1:3,i) @@ -292,14 +292,14 @@ module subroutine plastic_kinehardening_dotState(Mp,ph,me) real(pReal) :: & sumGamma - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & gdot_pos,gdot_neg - associate(prm => param(phase_plasticityInstance(ph)), stt => state(phase_plasticityInstance(ph)),& - dot => dotState(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph)),& + dot => dotState(phase_plasticInstance(ph))) - call kinetics(Mp,phase_plasticityInstance(ph),me,gdot_pos,gdot_neg) + call kinetics(Mp,phase_plasticInstance(ph),me,gdot_pos,gdot_neg) dot%accshear(:,me) = abs(gdot_pos+gdot_neg) sumGamma = sum(stt%accshear(:,me)) diff --git a/src/phase_mechanics_plastic_none.f90 b/src/phase_mechanical_plastic_none.f90 similarity index 96% rename from src/phase_mechanics_plastic_none.f90 rename to src/phase_mechanical_plastic_none.f90 index b8d1678eb..b95eaa039 100644 --- a/src/phase_mechanics_plastic_none.f90 +++ b/src/phase_mechanical_plastic_none.f90 @@ -44,7 +44,7 @@ module function plastic_none_init() result(myPlasticity) phase => phases%get(p) if(.not. myPlasticity(p)) cycle Nconstituents = count(material_phaseAt2 == p) - call constitutive_allocateState(plasticState(p),Nconstituents,0,0,0) + call phase_allocateState(plasticState(p),Nconstituents,0,0,0) enddo end function plastic_none_init diff --git a/src/phase_mechanics_plastic_nonlocal.f90 b/src/phase_mechanical_plastic_nonlocal.f90 similarity index 97% rename from src/phase_mechanics_plastic_nonlocal.f90 rename to src/phase_mechanical_plastic_nonlocal.f90 index 693ffcf93..724086916 100644 --- a/src/phase_mechanics_plastic_nonlocal.f90 +++ b/src/phase_mechanical_plastic_nonlocal.f90 @@ -407,7 +407,7 @@ module function plastic_nonlocal_init() result(myPlasticity) 'maxDipoleHeightEdge ','maxDipoleHeightScrew' ]) * prm%sum_N_sl !< other dependent state variables that are not updated by microstructure sizeDeltaState = sizeDotState - call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,sizeDeltaState) + call phase_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,sizeDeltaState) plasticState(p)%nonlocal = pl%get_asBool('nonlocal') if(plasticState(p)%nonlocal .and. .not. allocated(IPneighborhood)) & @@ -642,8 +642,8 @@ module subroutine nonlocal_dependentState(instance, me, ip, el) rho0 = getRho0(instance,me,ip,el) if (.not. phase_localPlasticity(material_phaseAt(1,el)) .and. prm%shortRangeStressCorrection) then ph = material_phaseAt(1,el) - invFp = math_inv33(constitutive_mech_Fp(ph)%data(1:3,1:3,me)) - invFe = math_inv33(constitutive_mech_Fe(ph)%data(1:3,1:3,me)) + invFp = math_inv33(phase_mechanical_Fp(ph)%data(1:3,1:3,me)) + invFe = math_inv33(phase_mechanical_Fe(ph)%data(1:3,1:3,me)) rho_edg_delta = rho0(:,mob_edg_pos) - rho0(:,mob_edg_neg) rho_scr_delta = rho0(:,mob_scr_pos) - rho0(:,mob_scr_neg) @@ -662,7 +662,7 @@ module subroutine nonlocal_dependentState(instance, me, ip, el) neighbor_ip = IPneighborhood(2,n,ip,el) no = material_phasememberAt(1,neighbor_ip,neighbor_el) if (neighbor_el > 0 .and. neighbor_ip > 0) then - neighbor_instance = phase_plasticityInstance(material_phaseAt(1,neighbor_el)) + neighbor_instance = phase_plasticInstance(material_phaseAt(1,neighbor_el)) if (neighbor_instance == instance) then nRealNeighbors = nRealNeighbors + 1.0_pReal @@ -782,25 +782,25 @@ module subroutine nonlocal_LpAndItsTangent(Lp,dLp_dMp, & l, & t, & !< dislocation type s !< index of my current slip system - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl,8) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl,8) :: & rhoSgl !< single dislocation densities (including blocked) - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl,10) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl,10) :: & rho - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl,4) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl,4) :: & v, & !< velocity tauNS, & !< resolved shear stress including non Schmid and backstress terms dv_dtau, & !< velocity derivative with respect to the shear stress dv_dtauNS !< velocity derivative with respect to the shear stress - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & tau, & !< resolved shear stress including backstress terms gdotTotal !< shear rate - associate(prm => param(phase_plasticityInstance(ph)),dst=>microstructure(phase_plasticityInstance(ph)),& - stt=>state(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)),dst=>microstructure(phase_plasticInstance(ph)),& + stt=>state(phase_plasticInstance(ph))) ns = prm%sum_N_sl !*** shortcut to state variables - rho = getRho(phase_plasticityInstance(ph),me,ip,el) + rho = getRho(phase_plasticInstance(ph),me,ip,el) rhoSgl = rho(:,sgl) do s = 1,ns @@ -820,7 +820,7 @@ module subroutine nonlocal_LpAndItsTangent(Lp,dLp_dMp, & ! edges call kinetics(v(:,1), dv_dtau(:,1), dv_dtauNS(:,1), & - tau, tauNS(:,1), dst%tau_pass(:,me),1,Temperature, phase_plasticityInstance(ph)) + tau, tauNS(:,1), dst%tau_pass(:,me),1,Temperature, phase_plasticInstance(ph)) v(:,2) = v(:,1) dv_dtau(:,2) = dv_dtau(:,1) dv_dtauNS(:,2) = dv_dtauNS(:,1) @@ -833,7 +833,7 @@ module subroutine nonlocal_LpAndItsTangent(Lp,dLp_dMp, & else do t = 3,4 call kinetics(v(:,t), dv_dtau(:,t), dv_dtauNS(:,t), & - tau, tauNS(:,t), dst%tau_pass(:,me),2,Temperature, phase_plasticityInstance(ph)) + tau, tauNS(:,t), dst%tau_pass(:,me),2,Temperature, phase_plasticInstance(ph)) enddo endif @@ -992,7 +992,7 @@ module subroutine nonlocal_dotState(Mp, Temperature,timestep, & c, & !< character of dislocation t, & !< type of dislocation s !< index of my current slip system - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl,10) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl,10) :: & rho, & rho0, & !< dislocation density at beginning of time step rhoDot, & !< density evolution @@ -1000,17 +1000,17 @@ module subroutine nonlocal_dotState(Mp, Temperature,timestep, & rhoDotSingle2DipoleGlide, & !< density evolution by dipole formation (by glide) rhoDotAthermalAnnihilation, & !< density evolution by athermal annihilation rhoDotThermalAnnihilation !< density evolution by thermal annihilation - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl,8) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl,8) :: & rhoSgl, & !< current single dislocation densities (positive/negative screw and edge without dipoles) my_rhoSgl0 !< single dislocation densities of central ip (positive/negative screw and edge without dipoles) - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl,4) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl,4) :: & v, & !< current dislocation glide velocity v0, & gdot !< shear rates - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & tau, & !< current resolved shear stress vClimb !< climb velocity of edge dipoles - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl,2) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl,2) :: & rhoDip, & !< current dipole dislocation densities (screw and edge dipoles) dLower, & !< minimum stable dipole distance for edges and screws dUpper !< current maximum stable dipole distance for edges and screws @@ -1022,22 +1022,22 @@ module subroutine nonlocal_dotState(Mp, Temperature,timestep, & return endif - associate(prm => param(phase_plasticityInstance(ph)), & - dst => microstructure(phase_plasticityInstance(ph)), & - dot => dotState(phase_plasticityInstance(ph)), & - stt => state(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)), & + dst => microstructure(phase_plasticInstance(ph)), & + dot => dotState(phase_plasticInstance(ph)), & + stt => state(phase_plasticInstance(ph))) ns = prm%sum_N_sl tau = 0.0_pReal gdot = 0.0_pReal - rho = getRho(phase_plasticityInstance(ph),me,ip,el) + rho = getRho(phase_plasticInstance(ph),me,ip,el) rhoSgl = rho(:,sgl) rhoDip = rho(:,dip) - rho0 = getRho0(phase_plasticityInstance(ph),me,ip,el) + rho0 = getRho0(phase_plasticInstance(ph),me,ip,el) my_rhoSgl0 = rho0(:,sgl) - forall (s = 1:ns, t = 1:4) v(s,t) = plasticState(ph)%state(iV(s,t,phase_plasticityInstance(ph)),me) + forall (s = 1:ns, t = 1:4) v(s,t) = plasticState(ph)%state(iV(s,t,phase_plasticInstance(ph)),me) gdot = rhoSgl(:,1:4) * v * spread(prm%b_sl,2,4) #ifdef DEBUG @@ -1086,7 +1086,7 @@ module subroutine nonlocal_dotState(Mp, Temperature,timestep, & * sqrt(stt%rho_forest(:,me)) / prm%i_sl / prm%b_sl, 2, 4) endif isBCC - forall (s = 1:ns, t = 1:4) v0(s,t) = plasticState(ph)%state0(iV(s,t,phase_plasticityInstance(ph)),me) + forall (s = 1:ns, t = 1:4) v0(s,t) = plasticState(ph)%state0(iV(s,t,phase_plasticInstance(ph)),me) !**************************************************************************** @@ -1142,7 +1142,7 @@ module subroutine nonlocal_dotState(Mp, Temperature,timestep, & - rhoDip(s,1) / timestep - rhoDotAthermalAnnihilation(s,9) & - rhoDotSingle2DipoleGlide(s,9)) ! make sure that we do not annihilate more dipoles than we have - rhoDot = rhoDotFlux(timestep, phase_plasticityInstance(ph),me,ip,el) & + rhoDot = rhoDotFlux(timestep, phase_plasticInstance(ph),me,ip,el) & + rhoDotMultiplication & + rhoDotSingle2DipoleGlide & + rhoDotAthermalAnnihilation & @@ -1284,8 +1284,8 @@ function rhoDotFlux(timestep,instance,me,ip,el) m(1:3,:,3) = -prm%slip_transverse m(1:3,:,4) = prm%slip_transverse - my_F = constitutive_mech_F(ph)%data(1:3,1:3,me) - my_Fe = matmul(my_F, math_inv33(constitutive_mech_Fp(ph)%data(1:3,1:3,me))) + my_F = phase_mechanical_F(ph)%data(1:3,1:3,me) + my_Fe = matmul(my_F, math_inv33(phase_mechanical_Fp(ph)%data(1:3,1:3,me))) neighbors: do n = 1,nIPneighbors @@ -1301,9 +1301,9 @@ function rhoDotFlux(timestep,instance,me,ip,el) opposite_n = IPneighborhood(3,opposite_neighbor,ip,el) if (neighbor_n > 0) then ! if neighbor exists, average deformation gradient - neighbor_instance = phase_plasticityInstance(material_phaseAt(1,neighbor_el)) - neighbor_F = constitutive_mech_F(np)%data(1:3,1:3,no) - neighbor_Fe = matmul(neighbor_F, math_inv33(constitutive_mech_Fp(np)%data(1:3,1:3,no))) + neighbor_instance = phase_plasticInstance(material_phaseAt(1,neighbor_el)) + neighbor_F = phase_mechanical_F(np)%data(1:3,1:3,no) + neighbor_Fe = matmul(neighbor_F, math_inv33(phase_mechanical_Fp(np)%data(1:3,1:3,no))) Favg = 0.5_pReal * (my_F + neighbor_F) else ! if no neighbor, take my value as average Favg = my_F diff --git a/src/phase_mechanics_plastic_phenopowerlaw.f90 b/src/phase_mechanical_plastic_phenopowerlaw.f90 similarity index 96% rename from src/phase_mechanics_plastic_phenopowerlaw.f90 rename to src/phase_mechanical_plastic_phenopowerlaw.f90 index ea2f53100..f9791faa6 100644 --- a/src/phase_mechanics_plastic_phenopowerlaw.f90 +++ b/src/phase_mechanical_plastic_phenopowerlaw.f90 @@ -231,7 +231,7 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) + call phase_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization @@ -300,18 +300,18 @@ pure module subroutine phenopowerlaw_LpAndItsTangent(Lp,dLp_dMp,Mp,ph,me) integer :: & i,k,l,m,n - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & gdot_slip_pos,gdot_slip_neg, & dgdot_dtauslip_pos,dgdot_dtauslip_neg - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_tw) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_tw) :: & gdot_twin,dgdot_dtautwin Lp = 0.0_pReal dLp_dMp = 0.0_pReal - associate(prm => param(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph))) - call kinetics_slip(Mp,phase_plasticityInstance(ph),me,gdot_slip_pos,gdot_slip_neg,dgdot_dtauslip_pos,dgdot_dtauslip_neg) + call kinetics_slip(Mp,phase_plasticInstance(ph),me,gdot_slip_pos,gdot_slip_neg,dgdot_dtauslip_pos,dgdot_dtauslip_neg) slipSystems: do i = 1, prm%sum_N_sl Lp = Lp + (gdot_slip_pos(i)+gdot_slip_neg(i))*prm%P_sl(1:3,1:3,i) forall (k=1:3,l=1:3,m=1:3,n=1:3) & @@ -320,7 +320,7 @@ pure module subroutine phenopowerlaw_LpAndItsTangent(Lp,dLp_dMp,Mp,ph,me) + dgdot_dtauslip_neg(i) * prm%P_sl(k,l,i) * prm%nonSchmid_neg(m,n,i) enddo slipSystems - call kinetics_twin(Mp,phase_plasticityInstance(ph),me,gdot_twin,dgdot_dtautwin) + call kinetics_twin(Mp,phase_plasticInstance(ph),me,gdot_twin,dgdot_dtautwin) twinSystems: do i = 1, prm%sum_N_tw Lp = Lp + gdot_twin(i)*prm%P_tw(1:3,1:3,i) forall (k=1:3,l=1:3,m=1:3,n=1:3) & @@ -348,12 +348,12 @@ module subroutine phenopowerlaw_dotState(Mp,ph,me) c_SlipSlip,c_TwinSlip,c_TwinTwin, & xi_slip_sat_offset,& sumGamma,sumF - real(pReal), dimension(param(phase_plasticityInstance(ph))%sum_N_sl) :: & + real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: & left_SlipSlip,right_SlipSlip, & gdot_slip_pos,gdot_slip_neg - associate(prm => param(phase_plasticityInstance(ph)), stt => state(phase_plasticityInstance(ph)), & - dot => dotState(phase_plasticityInstance(ph))) + associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph)), & + dot => dotState(phase_plasticInstance(ph))) sumGamma = sum(stt%gamma_slip(:,me)) sumF = sum(stt%gamma_twin(:,me)/prm%gamma_tw_char) @@ -373,9 +373,9 @@ module subroutine phenopowerlaw_dotState(Mp,ph,me) !-------------------------------------------------------------------------------------------------- ! shear rates - call kinetics_slip(Mp,phase_plasticityInstance(ph),me,gdot_slip_pos,gdot_slip_neg) + call kinetics_slip(Mp,phase_plasticInstance(ph),me,gdot_slip_pos,gdot_slip_neg) dot%gamma_slip(:,me) = abs(gdot_slip_pos+gdot_slip_neg) - call kinetics_twin(Mp,phase_plasticityInstance(ph),me,dot%gamma_twin(:,me)) + call kinetics_twin(Mp,phase_plasticInstance(ph),me,dot%gamma_twin(:,me)) !-------------------------------------------------------------------------------------------------- ! hardening diff --git a/src/phase_thermal.f90 b/src/phase_thermal.f90 index d5d52b010..85ca2f7a5 100644 --- a/src/phase_thermal.f90 +++ b/src/phase_thermal.f90 @@ -124,7 +124,7 @@ end subroutine thermal_init !---------------------------------------------------------------------------------------------- !< @brief calculates thermal dissipation rate !---------------------------------------------------------------------------------------------- -module subroutine constitutive_thermal_getRate(TDot, ph,me) +module subroutine phase_thermal_getRate(TDot, ph,me) integer, intent(in) :: ph, me real(pReal), intent(out) :: & @@ -153,13 +153,13 @@ module subroutine constitutive_thermal_getRate(TDot, ph,me) enddo -end subroutine constitutive_thermal_getRate +end subroutine phase_thermal_getRate !-------------------------------------------------------------------------------------------------- !> @brief contains the constitutive equation for calculating the rate of change of microstructure !-------------------------------------------------------------------------------------------------- -function constitutive_thermal_collectDotState(ph,me) result(broken) +function phase_thermal_collectDotState(ph,me) result(broken) integer, intent(in) :: ph, me logical :: broken @@ -178,7 +178,7 @@ function constitutive_thermal_collectDotState(ph,me) result(broken) enddo SourceLoop -end function constitutive_thermal_collectDotState +end function phase_thermal_collectDotState module function thermal_stress(Delta_t,ph,me) result(converged_) @@ -207,7 +207,7 @@ function integrateThermalState(Delta_t, ph,me) result(broken) so, & sizeDotState - broken = constitutive_thermal_collectDotState(ph,me) + broken = phase_thermal_collectDotState(ph,me) if(broken) return do so = 1, thermal_Nsources(ph) @@ -264,7 +264,7 @@ end function thermal_dot_T !---------------------------------------------------------------------------------------------- !< @brief Set temperature !---------------------------------------------------------------------------------------------- -module subroutine constitutive_thermal_setField(T,dot_T, co,ce) +module subroutine phase_thermal_setField(T,dot_T, co,ce) real(pReal), intent(in) :: T, dot_T integer, intent(in) :: ce, co @@ -273,7 +273,7 @@ module subroutine constitutive_thermal_setField(T,dot_T, co,ce) current(material_phaseAt2(co,ce))%T(material_phaseMemberAt2(co,ce)) = T current(material_phaseAt2(co,ce))%dot_T(material_phaseMemberAt2(co,ce)) = dot_T -end subroutine constitutive_thermal_setField +end subroutine phase_thermal_setField diff --git a/src/phase_thermal_dissipation.f90 b/src/phase_thermal_dissipation.f90 index 7a043308d..3c0095401 100644 --- a/src/phase_thermal_dissipation.f90 +++ b/src/phase_thermal_dissipation.f90 @@ -56,7 +56,7 @@ module function dissipation_init(source_length) result(mySources) prm%kappa = src%get_asFloat('kappa') Nconstituents = count(material_phaseAt2 == ph) - call constitutive_allocateState(thermalState(ph)%p(so),Nconstituents,0,0,0) + call phase_allocateState(thermalState(ph)%p(so),Nconstituents,0,0,0) end associate endif @@ -78,7 +78,7 @@ module subroutine dissipation_getRate(TDot, ph,me) associate(prm => param(ph)) - TDot = prm%kappa*sum(abs(mech_S(ph,me)*mech_L_p(ph,me))) + TDot = prm%kappa*sum(abs(mechanical_S(ph,me)*mechanical_L_p(ph,me))) end associate end subroutine dissipation_getRate diff --git a/src/phase_thermal_externalheat.f90 b/src/phase_thermal_externalheat.f90 index ddae6763b..d2874ded0 100644 --- a/src/phase_thermal_externalheat.f90 +++ b/src/phase_thermal_externalheat.f90 @@ -69,7 +69,7 @@ module function externalheat_init(source_length) result(mySources) prm%f_T = src%get_asFloats('f_T',requiredSize = size(prm%t_n)) Nconstituents = count(material_phaseAt2 == ph) - call constitutive_allocateState(thermalState(ph)%p(so),Nconstituents,1,1,0) + call phase_allocateState(thermalState(ph)%p(so),Nconstituents,1,1,0) end associate endif enddo From e8fae6b2a719eab09d2cdc21386bf11e28e7c320 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 11 Feb 2021 21:56:53 +0100 Subject: [PATCH 08/10] simplified --- python/damask/_vtk.py | 30 ++++++++++++++---------------- 1 file changed, 14 insertions(+), 16 deletions(-) diff --git a/python/damask/_vtk.py b/python/damask/_vtk.py index ad017ee1e..00a07efa5 100644 --- a/python/damask/_vtk.py +++ b/python/damask/_vtk.py @@ -347,23 +347,21 @@ class VTK: See http://compilatrix.com/article/vtk-1 for further ideas. """ - def screen_size(): + try: + import wx + _ = wx.App(False) # noqa + width, height = wx.GetDisplaySize() + except ImportError: try: - import wx - _ = wx.App(False) # noqa - width, height = wx.GetDisplaySize() - except ImportError: - try: - import tkinter - tk = tkinter.Tk() - width = tk.winfo_screenwidth() - height = tk.winfo_screenheight() - tk.destroy() - except Exception as e: - width = 1024 - height = 768 + import tkinter + tk = tkinter.Tk() + width = tk.winfo_screenwidth() + height = tk.winfo_screenheight() + tk.destroy() + except Exception as e: + width = 1024 + height = 768 - return (width,height) mapper = vtk.vtkDataSetMapper() mapper.SetInputData(self.vtk_data) actor = vtk.vtkActor() @@ -377,7 +375,7 @@ class VTK: ren.AddActor(actor) ren.SetBackground(0.2,0.2,0.2) - window.SetSize(screen_size[0],screen_size[1]) + window.SetSize(width,height) iren = vtk.vtkRenderWindowInteractor() iren.SetRenderWindow(window) From 09c330b8d2932003d5e4942386cef80913a13cef Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 11 Feb 2021 21:57:44 +0100 Subject: [PATCH 09/10] test for CRLF error (got lost) --- PRIVATE | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PRIVATE b/PRIVATE index 91a4329a6..e74cf0062 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 91a4329a6fe073fc2ef17e5c2c8e2f796e3b897b +Subproject commit e74cf00628285a587ced1e551cc9837c1011ca1c From 26768375781c781db96ab3e0c1f425becf835411 Mon Sep 17 00:00:00 2001 From: Test User Date: Fri, 12 Feb 2021 22:30:53 +0100 Subject: [PATCH 10/10] [skip ci] updated version information after successful test of v3.0.0-alpha2-428-gb7b764ab5 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index a1f129287..ce989b710 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v3.0.0-alpha2-421-ge96352b0e +v3.0.0-alpha2-428-gb7b764ab5