diff --git a/PRIVATE b/PRIVATE index 50eb21714..67fcdfd00 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 50eb21714e2f501b111bb62096ebb6a5bfc6708a +Subproject commit 67fcdfd001a9e4f3da56d72ebeb0ea16d6d613e6 diff --git a/VERSION b/VERSION index fca0385fb..6ab575d5c 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v2.0.2-261-gbc3f6ae9 +v2.0.2-282-gc05337a2 diff --git a/lib/damask/config/material.py b/lib/damask/config/material.py index bb184b4d2..02658019d 100644 --- a/lib/damask/config/material.py +++ b/lib/damask/config/material.py @@ -277,5 +277,16 @@ class Material(): self.data[part.lower()][section.lower()][key.lower()] = value if newlen is not oldlen: print('Length of value was changed from %i to %i!'%(oldlen,newlen)) - + + def add_value(self, part=None, + section=None, + key=None, + value=None): + if not isinstance(value,list): + if not isinstance(value,str): + value = '%s'%value + value = [value] + print('adding %s:%s:%s with value %s '%(part.lower(),section.lower(),key.lower(),value)) + self.data[part.lower()][section.lower()][key.lower()] = value + self.data[part.lower()][section.lower()]['__order__'] += [key.lower()] diff --git a/src/DAMASK_abaqus.f b/src/DAMASK_abaqus.f index e91cbb0bb..69f6fba4b 100644 --- a/src/DAMASK_abaqus.f +++ b/src/DAMASK_abaqus.f @@ -3,38 +3,42 @@ !> @author Franz Roters, Max-Planck-Institut für Eisenforschung GmbH !> @author Koen Janssens, Paul Scherrer Institut !> @author Arun Prakash, Fraunhofer IWM +!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH !> @brief interfaces DAMASK with Abaqus/Standard !> @details put the included file abaqus_v6.env in either your home or model directory, !> it is a minimum Abaqus environment file containing all changes necessary to use the !> DAMASK subroutine (see Abaqus documentation for more information on the use of abaqus_v6.env) !-------------------------------------------------------------------------------------------------- - -#ifndef INT -#define INT 4 -#endif - -#ifndef FLOAT -#define FLOAT 8 -#endif - #define Abaqus #include "prec.f90" module DAMASK_interface -implicit none -character(len=4), dimension(2), parameter :: INPUTFILEEXTENSION = ['.pes','.inp'] -character(len=4), parameter :: LOGFILEEXTENSION = '.log' + implicit none + private + character(len=4), dimension(2), parameter, public :: INPUTFILEEXTENSION = ['.pes','.inp'] + character(len=4), parameter, public :: LOGFILEEXTENSION = '.log' + + public :: & + DAMASK_interface_init, & + getSolverJobName contains !-------------------------------------------------------------------------------------------------- -!> @brief just reporting +!> @brief reports and sets working directory !-------------------------------------------------------------------------------------------------- subroutine DAMASK_interface_init + use ifport, only: & + CHDIR + + implicit none integer, dimension(8) :: & dateAndTime ! type default integer + integer :: lenOutDir,ierr + character(len=256) :: wd + call date_and_time(values = dateAndTime) write(6,'(/,a)') ' <<<+- DAMASK_abaqus_std -+>>>' write(6,'(/,a)') ' Roters et al., Computational Materials Science, 2018' @@ -46,26 +50,19 @@ subroutine DAMASK_interface_init dateAndTime(6),':',& dateAndTime(7) write(6,'(/,a)') ' <<<+- DAMASK_interface init -+>>>' + + call getoutdir(wd, lenOutDir) + ierr = CHDIR(wd) + if (ierr /= 0) then + write(6,'(a20,a,a16)') ' working directory "',trim(wd),'" does not exist' + call quit(1) + endif + #include "compilation_info.f90" end subroutine DAMASK_interface_init -!-------------------------------------------------------------------------------------------------- -!> @brief using Abaqus/Standard function to get working directory name -!-------------------------------------------------------------------------------------------------- -character(1024) function getSolverWorkingDirectoryName() - - implicit none - integer :: lenOutDir - - getSolverWorkingDirectoryName='' - call getoutdir(getSolverWorkingDirectoryName, lenOutDir) - getSolverWorkingDirectoryName=trim(getSolverWorkingDirectoryName)//'/' - -end function getSolverWorkingDirectoryName - - !-------------------------------------------------------------------------------------------------- !> @brief using Abaqus/Standard function to get solver job name !-------------------------------------------------------------------------------------------------- @@ -79,10 +76,17 @@ character(1024) function getSolverJobName() end function getSolverJobName + end module DAMASK_interface -#include "commercialFEM_fileList.f90" + + +#include "commercialFEM_fileList.f90" + +!-------------------------------------------------------------------------------------------------- +!> @brief This is the Abaqus std user subroutine for defining material behavior +!-------------------------------------------------------------------------------------------------- subroutine UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,& RPL,DDSDDT,DRPLDE,DRPLDT,STRAN,DSTRAN,& TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMNAME,NDI,NSHR,NTENS,& diff --git a/src/DAMASK_marc.f90 b/src/DAMASK_marc.f90 index 600db3c2e..f3130c5cd 100644 --- a/src/DAMASK_marc.f90 +++ b/src/DAMASK_marc.f90 @@ -1,15 +1,3 @@ -#define QUOTE(x) #x -#define PASTE(x,y) x ## y - -#ifndef INT -#define INT 4 -#endif - -#ifndef FLOAT -#define FLOAT 8 -#endif - -#include "prec.f90" !-------------------------------------------------------------------------------------------------- !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH !> @author Franz Roters, Max-Planck-Institut für Eisenforschung GmbH @@ -17,13 +5,12 @@ !> @author W.A. Counts !> @author Denny Tjahjanto, Max-Planck-Institut für Eisenforschung GmbH !> @author Christoph Kords, Max-Planck-Institut für Eisenforschung GmbH -!> @brief Material subroutine for MSC.Marc +!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH +!> @brief Interfaces DAMASK with MSC.Marc !> @details Usage: !> @details - choose material as hypela2 !> @details - set statevariable 2 to index of homogenization !> @details - set statevariable 3 to index of microstructure -!> @details - make sure the file "material.config" exists in the working directory -!> @details - make sure the file "numerics.config" exists in the working directory !> @details - use nonsymmetric option for solver (e.g. direct profile or multifrontal sparse, the latter seems to be faster!) !> @details - in case of ddm (domain decomposition) a SYMMETRIC solver has to be used, i.e uncheck "non-symmetric" !> @details Marc subroutines used: @@ -34,23 +21,36 @@ !> @details - concom: lovl, inc !> @details - creeps: timinc !-------------------------------------------------------------------------------------------------- +#define QUOTE(x) #x +#define PASTE(x,y) x ## y + +#include "prec.f90" + module DAMASK_interface implicit none - character(len=4), parameter :: InputFileExtension = '.dat' - character(len=4), parameter :: LogFileExtension = '.log' + private + character(len=4), parameter, public :: InputFileExtension = '.dat' + character(len=4), parameter, public :: LogFileExtension = '.log' + + public :: & + DAMASK_interface_init, & + getSolverJobName contains - !-------------------------------------------------------------------------------------------------- -!> @brief only output of current version +!> @brief reports and sets working directory !-------------------------------------------------------------------------------------------------- subroutine DAMASK_interface_init + use ifport, only: & + CHDIR implicit none integer, dimension(8) :: & dateAndTime ! type default integer + integer :: ierr + character(len=1024) :: wd call date_and_time(values = dateAndTime) write(6,'(/,a)') ' <<<+- DAMASK_Marc -+>>>' @@ -64,27 +64,17 @@ subroutine DAMASK_interface_init dateAndTime(7) write(6,'(/,a)') ' <<<+- DAMASK_interface init -+>>>' #include "compilation_info.f90" + inquire(5, name=wd) ! determine inputputfile + wd = wd(1:scan(wd,'/',back=.true.)) + ierr = CHDIR(wd) + if (ierr /= 0) then + write(6,'(a20,a,a16)') ' working directory "',trim(wd),'" does not exist' + call quit(1) + endif end subroutine DAMASK_interface_init -!-------------------------------------------------------------------------------------------------- -!> @brief returns the current workingDir -!-------------------------------------------------------------------------------------------------- -function getSolverWorkingDirectoryName() - - implicit none - character(1024) getSolverWorkingDirectoryName, inputName - character(len=*), parameter :: pathSep = achar(47)//achar(92) ! forward and backward slash - - getSolverWorkingDirectoryName='' - inputName='' - inquire(5, name=inputName) ! determine inputputfile - getSolverWorkingDirectoryName=inputName(1:scan(inputName,pathSep,back=.true.)) - -end function getSolverWorkingDirectoryName - - !-------------------------------------------------------------------------------------------------- !> @brief solver job name (no extension) as combination of geometry and load case name !-------------------------------------------------------------------------------------------------- @@ -109,6 +99,9 @@ end function getSolverJobName end module DAMASK_interface + + + #include "commercialFEM_fileList.f90" !-------------------------------------------------------------------------------------------------- @@ -118,17 +111,6 @@ end module DAMASK_interface !> @details !> @details (2) Use the -> 'Plasticity,3' card(=update+finite+large disp+constant d) !> @details in the parameter section of input deck (updated Lagrangian formulation). -!> @details -!> @details The following operation obtains U (stretch tensor) at t=n+1 : -!> @details -!> @details call scla(un1,0.d0,itel,itel,1) -!> @details do k=1,3 -!> @details do i=1,3 -!> @details do j=1,3 -!> @details un1(i,j)=un1(i,j)+dsqrt(strechn1(k))*eigvn1(i,k)*eigvn1(j,k) -!> @details enddo -!> @details enddo -!> @details enddo !-------------------------------------------------------------------------------------------------- subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,nshear,disp, & dispt,coord,ffn,frotn,strechn,eigvn,ffn1,frotn1, & diff --git a/src/DAMASK_spectral.f90 b/src/DAMASK_spectral.f90 index 0d77c57f5..2ed94d06a 100644 --- a/src/DAMASK_spectral.f90 +++ b/src/DAMASK_spectral.f90 @@ -20,11 +20,12 @@ program DAMASK_spectral pReal, & tol_math_check, & dNeq + use system_routines, only: & + getCWD use DAMASK_interface, only: & DAMASK_interface_init, & loadCaseFile, & geometryFile, & - getSolverWorkingDirectoryName, & getSolverJobName, & appendToOutFile use IO, only: & @@ -133,7 +134,9 @@ program DAMASK_spectral lastRestartWritten = 0_pInt, & !< total increment # at which last restart information was written stagIter character(len=6) :: loadcase_string - character(len=1024) :: incInfo !< string parsed to solution with information about current load case + character(len=1024) :: & + incInfo, & !< string parsed to solution with information about current load case + workingDir type(tLoadCase), allocatable, dimension(:) :: loadCases !< array of all load cases type(tSolutionState), allocatable, dimension(:) :: solres integer(MPI_OFFSET_KIND) :: fileOffset @@ -381,10 +384,11 @@ program DAMASK_spectral ! write header of output file if (worldrank == 0) then if (.not. appendToOutFile) then ! after restart, append to existing results file - open(newunit=resUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//& + if (getCWD(workingDir)) call IO_error(106_pInt,ext_msg=trim(workingDir)) + open(newunit=resUnit,file=trim(getSolverJobName())//& '.spectralOut',form='UNFORMATTED',status='REPLACE') write(resUnit) 'load:', trim(loadCaseFile) ! ... and write header - write(resUnit) 'workingdir:', trim(getSolverWorkingDirectoryName()) + write(resUnit) 'workingdir:', trim(workingDir) write(resUnit) 'geometry:', trim(geometryFile) write(resUnit) 'grid:', grid write(resUnit) 'size:', geomSize @@ -397,14 +401,14 @@ program DAMASK_spectral write(resUnit) 'startingIncrement:', restartInc ! start with writing out the previous inc write(resUnit) 'eoh' close(resUnit) ! end of header - open(newunit=statUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//& + open(newunit=statUnit,file=trim(getSolverJobName())//& '.sta',form='FORMATTED',status='REPLACE') write(statUnit,'(a)') 'Increment Time CutbackLevel Converged IterationsNeeded' ! statistics file if (iand(debug_level(debug_spectral),debug_levelBasic) /= 0) & write(6,'(/,a)') ' header of result and statistics file written out' flush(6) else ! open new files ... - open(newunit=statUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//& + open(newunit=statUnit,file=trim(getSolverJobName())//& '.sta',form='FORMATTED', position='APPEND', status='OLD') endif endif @@ -415,8 +419,7 @@ program DAMASK_spectral outputSize(worldrank+1) = size(materialpoint_results,kind=MPI_OFFSET_KIND)*int(pReal,MPI_OFFSET_KIND) call MPI_allreduce(MPI_IN_PLACE,outputSize,worldsize,MPI_LONG,MPI_SUM,PETSC_COMM_WORLD,ierr) ! get total output size over each process if (ierr /= 0_pInt) call IO_error(error_ID=894_pInt, ext_msg='MPI_allreduce') - call MPI_file_open(PETSC_COMM_WORLD, & - trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//'.spectralOut', & + call MPI_file_open(PETSC_COMM_WORLD, trim(getSolverJobName())//'.spectralOut', & MPI_MODE_WRONLY + MPI_MODE_APPEND, & MPI_INFO_NULL, & resUnit, & diff --git a/src/IO.f90 b/src/IO.f90 index a7e77f0f4..807686e86 100644 --- a/src/IO.f90 +++ b/src/IO.f90 @@ -195,18 +195,14 @@ end subroutine IO_checkAndRewind !> @details like IO_open_file_stat, but error is handled via call to IO_error and not via return !! value !-------------------------------------------------------------------------------------------------- -subroutine IO_open_file(fileUnit,relPath) - use DAMASK_interface, only: & - getSolverWorkingDirectoryName +subroutine IO_open_file(fileUnit,path) implicit none integer(pInt), intent(in) :: fileUnit !< file unit - character(len=*), intent(in) :: relPath !< relative path from working directory + character(len=*), intent(in) :: path !< relative path from working directory integer(pInt) :: myStat - character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//relPath open(fileUnit,status='old',iostat=myStat,file=path) if (myStat /= 0_pInt) call IO_error(100_pInt,el=myStat,ext_msg=path) @@ -218,18 +214,14 @@ end subroutine IO_open_file !! directory !> @details Like IO_open_file, but error is handled via return value and not via call to IO_error !-------------------------------------------------------------------------------------------------- -logical function IO_open_file_stat(fileUnit,relPath) - use DAMASK_interface, only: & - getSolverWorkingDirectoryName +logical function IO_open_file_stat(fileUnit,path) implicit none integer(pInt), intent(in) :: fileUnit !< file unit - character(len=*), intent(in) :: relPath !< relative path from working directory + character(len=*), intent(in) :: path !< relative path from working directory integer(pInt) :: myStat - character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//relPath open(fileUnit,status='old',iostat=myStat,file=path) IO_open_file_stat = (myStat == 0_pInt) @@ -244,7 +236,6 @@ end function IO_open_file_stat !-------------------------------------------------------------------------------------------------- subroutine IO_open_jobFile(fileUnit,ext) use DAMASK_interface, only: & - getSolverWorkingDirectoryName, & getSolverJobName implicit none @@ -254,7 +245,7 @@ subroutine IO_open_jobFile(fileUnit,ext) integer(pInt) :: myStat character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//'.'//ext + path = trim(getSolverJobName())//'.'//ext open(fileUnit,status='old',iostat=myStat,file=path) if (myStat /= 0_pInt) call IO_error(100_pInt,el=myStat,ext_msg=path) @@ -269,7 +260,6 @@ end subroutine IO_open_jobFile !-------------------------------------------------------------------------------------------------- logical function IO_open_jobFile_stat(fileUnit,ext) use DAMASK_interface, only: & - getSolverWorkingDirectoryName, & getSolverJobName implicit none @@ -279,7 +269,7 @@ logical function IO_open_jobFile_stat(fileUnit,ext) integer(pInt) :: myStat character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//'.'//ext + path = trim(getSolverJobName())//'.'//ext open(fileUnit,status='old',iostat=myStat,file=path) IO_open_jobFile_stat = (myStat == 0_pInt) @@ -292,7 +282,6 @@ end function IO_open_JobFile_stat !-------------------------------------------------------------------------------------------------- subroutine IO_open_inputFile(fileUnit,modelName) use DAMASK_interface, only: & - getSolverWorkingDirectoryName,& getSolverJobName, & inputFileExtension @@ -306,23 +295,23 @@ subroutine IO_open_inputFile(fileUnit,modelName) integer(pInt) :: fileType fileType = 1_pInt ! assume .pes - path = trim(getSolverWorkingDirectoryName())//trim(modelName)//inputFileExtension(fileType) ! attempt .pes, if it exists: it should be used + path = trim(modelName)//inputFileExtension(fileType) ! attempt .pes, if it exists: it should be used open(fileUnit+1,status='old',iostat=myStat,file=path) if(myStat /= 0_pInt) then ! if .pes does not work / exist; use conventional extension, i.e.".inp" fileType = 2_pInt - path = trim(getSolverWorkingDirectoryName())//trim(modelName)//inputFileExtension(fileType) + path = trim(modelName)//inputFileExtension(fileType) open(fileUnit+1,status='old',iostat=myStat,file=path) endif if (myStat /= 0_pInt) call IO_error(100_pInt,el=myStat,ext_msg=path) - path = trim(getSolverWorkingDirectoryName())//trim(modelName)//inputFileExtension(fileType)//'_assembly' + path = trim(modelName)//inputFileExtension(fileType)//'_assembly' open(fileUnit,iostat=myStat,file=path) if (myStat /= 0_pInt) call IO_error(100_pInt,el=myStat,ext_msg=path) if (.not.abaqus_assembleInputFile(fileUnit,fileUnit+1_pInt)) call IO_error(103_pInt) ! strip comments and concatenate any "include"s close(fileUnit+1_pInt) #endif #ifdef Marc4DAMASK - path = trim(getSolverWorkingDirectoryName())//trim(modelName)//inputFileExtension + path = trim(modelName)//inputFileExtension open(fileUnit,status='old',iostat=myStat,file=path) if (myStat /= 0_pInt) call IO_error(100_pInt,el=myStat,ext_msg=path) #endif @@ -336,7 +325,6 @@ end subroutine IO_open_inputFile !-------------------------------------------------------------------------------------------------- subroutine IO_open_logFile(fileUnit) use DAMASK_interface, only: & - getSolverWorkingDirectoryName, & getSolverJobName, & LogFileExtension @@ -346,7 +334,7 @@ subroutine IO_open_logFile(fileUnit) integer(pInt) :: myStat character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//LogFileExtension + path = trim(getSolverJobName())//LogFileExtension open(fileUnit,status='old',iostat=myStat,file=path) if (myStat /= 0_pInt) call IO_error(100_pInt,el=myStat,ext_msg=path) @@ -360,7 +348,6 @@ end subroutine IO_open_logFile !-------------------------------------------------------------------------------------------------- subroutine IO_write_jobFile(fileUnit,ext) use DAMASK_interface, only: & - getSolverWorkingDirectoryName, & getSolverJobName implicit none @@ -370,7 +357,7 @@ subroutine IO_write_jobFile(fileUnit,ext) integer(pInt) :: myStat character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//'.'//ext + path = trim(getSolverJobName())//'.'//ext open(fileUnit,status='replace',iostat=myStat,file=path) if (myStat /= 0_pInt) call IO_error(100_pInt,el=myStat,ext_msg=path) @@ -383,7 +370,6 @@ end subroutine IO_write_jobFile !-------------------------------------------------------------------------------------------------- subroutine IO_write_jobRealFile(fileUnit,ext,recMultiplier) use DAMASK_interface, only: & - getSolverWorkingDirectoryName, & getSolverJobName implicit none @@ -394,7 +380,7 @@ subroutine IO_write_jobRealFile(fileUnit,ext,recMultiplier) integer(pInt) :: myStat character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//'.'//ext + path = trim(getSolverJobName())//'.'//ext if (present(recMultiplier)) then open(fileUnit,status='replace',form='unformatted',access='direct', & recl=pReal*recMultiplier,iostat=myStat,file=path) @@ -414,7 +400,6 @@ end subroutine IO_write_jobRealFile !-------------------------------------------------------------------------------------------------- subroutine IO_write_jobIntFile(fileUnit,ext,recMultiplier) use DAMASK_interface, only: & - getSolverWorkingDirectoryName, & getSolverJobName implicit none @@ -425,7 +410,7 @@ subroutine IO_write_jobIntFile(fileUnit,ext,recMultiplier) integer(pInt) :: myStat character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//'.'//ext + path = trim(getSolverJobName())//'.'//ext if (present(recMultiplier)) then open(fileUnit,status='replace',form='unformatted',access='direct', & recl=pInt*recMultiplier,iostat=myStat,file=path) @@ -444,8 +429,6 @@ end subroutine IO_write_jobIntFile !! located in current working directory !-------------------------------------------------------------------------------------------------- subroutine IO_read_realFile(fileUnit,ext,modelName,recMultiplier) - use DAMASK_interface, only: & - getSolverWorkingDirectoryName implicit none integer(pInt), intent(in) :: fileUnit !< file unit @@ -456,7 +439,7 @@ subroutine IO_read_realFile(fileUnit,ext,modelName,recMultiplier) integer(pInt) :: myStat character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//trim(modelName)//'.'//ext + path = trim(modelName)//'.'//ext if (present(recMultiplier)) then open(fileUnit,status='old',form='unformatted',access='direct', & recl=pReal*recMultiplier,iostat=myStat,file=path) @@ -474,8 +457,6 @@ end subroutine IO_read_realFile !! located in current working directory !-------------------------------------------------------------------------------------------------- subroutine IO_read_intFile(fileUnit,ext,modelName,recMultiplier) - use DAMASK_interface, only: & - getSolverWorkingDirectoryName implicit none integer(pInt), intent(in) :: fileUnit !< file unit @@ -486,7 +467,7 @@ subroutine IO_read_intFile(fileUnit,ext,modelName,recMultiplier) integer(pInt) :: myStat character(len=1024) :: path - path = trim(getSolverWorkingDirectoryName())//trim(modelName)//'.'//ext + path = trim(modelName)//'.'//ext if (present(recMultiplier)) then open(fileUnit,status='old',form='unformatted',access='direct', & recl=pInt*recMultiplier,iostat=myStat,file=path) @@ -1007,11 +988,7 @@ function IO_stringValue(string,chunkPos,myChunk,silent) logical :: warn - if (.not. present(silent)) then - warn = .false. - else - warn = silent - endif + warn = merge(silent,.false.,present(silent)) IO_stringValue = '' valuePresent: if (myChunk > chunkPos(1) .or. myChunk < 1_pInt) then @@ -1534,6 +1511,8 @@ subroutine IO_error(error_ID,el,ip,g,instance,ext_msg) msg = '{input} recursion limit reached' case (105_pInt) msg = 'unknown output:' + case (106_pInt) + msg = 'working directory does not exist:' !-------------------------------------------------------------------------------------------------- ! lattice error messages @@ -1905,8 +1884,6 @@ end function IO_verifyFloatValue !> including "include"s !-------------------------------------------------------------------------------------------------- recursive function abaqus_assembleInputFile(unit1,unit2) result(createSuccess) - use DAMASK_interface, only: & - getSolverWorkingDirectoryName implicit none integer(pInt), intent(in) :: unit1, & @@ -1923,7 +1900,7 @@ recursive function abaqus_assembleInputFile(unit1,unit2) result(createSuccess) chunkPos = IO_stringPos(line) if (IO_lc(IO_StringValue(line,chunkPos,1_pInt))=='*include') then - fname = trim(getSolverWorkingDirectoryName())//trim(line(9+scan(line(9:),'='):)) + fname = trim(line(9+scan(line(9:),'='):)) inquire(file=fname, exist=fexist) if (.not.(fexist)) then !$OMP CRITICAL (write2out) diff --git a/src/constitutive.f90 b/src/constitutive.f90 index 88d75dec1..9c3989a9c 100644 --- a/src/constitutive.f90 +++ b/src/constitutive.f90 @@ -160,7 +160,7 @@ subroutine constitutive_init() ! parse plasticities from config file if (any(phase_plasticity == PLASTICITY_NONE_ID)) call plastic_none_init if (any(phase_plasticity == PLASTICITY_ISOTROPIC_ID)) call plastic_isotropic_init - if (any(phase_plasticity == PLASTICITY_PHENOPOWERLAW_ID)) call plastic_phenopowerlaw_init(FILEUNIT) + if (any(phase_plasticity == PLASTICITY_PHENOPOWERLAW_ID)) call plastic_phenopowerlaw_init if (any(phase_plasticity == PLASTICITY_KINEHARDENING_ID)) call plastic_kinehardening_init(FILEUNIT) if (any(phase_plasticity == PLASTICITY_DISLOTWIN_ID)) call plastic_dislotwin_init(FILEUNIT) if (any(phase_plasticity == PLASTICITY_DISLOUCLA_ID)) call plastic_disloucla_init(FILEUNIT) diff --git a/src/math.f90 b/src/math.f90 index edf2ff5a6..955be4457 100644 --- a/src/math.f90 +++ b/src/math.f90 @@ -379,6 +379,9 @@ pure function math_expand(what,how) real(pReal), dimension(sum(how)) :: math_expand integer(pInt) :: i + if (sum(how) == 0_pInt) & + return + do i = 1_pInt, size(how) math_expand(sum(how(1:i-1))+1:sum(how(1:i))) = what(mod(i-1_pInt,size(what))+1_pInt) enddo diff --git a/src/plastic_none.f90 b/src/plastic_none.f90 index 140754556..5470c4a43 100644 --- a/src/plastic_none.f90 +++ b/src/plastic_none.f90 @@ -4,16 +4,9 @@ !> @brief material subroutine for purely elastic material !-------------------------------------------------------------------------------------------------- module plastic_none - use prec, only: & - pInt implicit none private - integer(pInt), dimension(:), allocatable, public, protected :: & - plastic_none_sizePostResults - - integer(pInt), dimension(:,:), allocatable, target, public :: & - plastic_none_sizePostResult !< size of each post result output public :: & plastic_none_init @@ -31,6 +24,8 @@ subroutine plastic_none_init compiler_version, & compiler_options #endif + use prec, only: & + pInt use debug, only: & debug_level, & debug_constitutive, & @@ -51,18 +46,13 @@ subroutine plastic_none_init integer(pInt) :: & maxNinstance, & phase, & - NofMyPhase, & - sizeState, & - sizeDotState, & - sizeDeltaState + NofMyPhase write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_NONE_label//' init -+>>>' write(6,'(a15,a)') ' Current time: ',IO_timeStamp() #include "compilation_info.f90" maxNinstance = int(count(phase_plasticity == PLASTICITY_none_ID),pInt) - if (maxNinstance == 0_pInt) return - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) & write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance @@ -70,37 +60,25 @@ subroutine plastic_none_init if (phase_plasticity(phase) == PLASTICITY_none_ID) then NofMyPhase=count(material_phase==phase) - sizeState = 0_pInt - plasticState(phase)%sizeState = sizeState - sizeDotState = sizeState - plasticState(phase)%sizeDotState = sizeDotState - sizeDeltaState = 0_pInt - plasticState(phase)%sizeDeltaState = sizeDeltaState - plasticState(phase)%sizePostResults = 0_pInt - plasticState(phase)%nSlip = 0_pInt - plasticState(phase)%nTwin = 0_pInt - plasticState(phase)%nTrans = 0_pInt - allocate(plasticState(phase)%aTolState (sizeState)) - allocate(plasticState(phase)%state0 (sizeState,NofMyPhase)) - allocate(plasticState(phase)%partionedState0 (sizeState,NofMyPhase)) - allocate(plasticState(phase)%subState0 (sizeState,NofMyPhase)) - allocate(plasticState(phase)%state (sizeState,NofMyPhase)) + allocate(plasticState(phase)%aTolState (0_pInt)) + allocate(plasticState(phase)%state0 (0_pInt,NofMyPhase)) + allocate(plasticState(phase)%partionedState0 (0_pInt,NofMyPhase)) + allocate(plasticState(phase)%subState0 (0_pInt,NofMyPhase)) + allocate(plasticState(phase)%state (0_pInt,NofMyPhase)) - allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase)) - allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase)) + allocate(plasticState(phase)%dotState (0_pInt,NofMyPhase)) + allocate(plasticState(phase)%deltaState (0_pInt,NofMyPhase)) if (any(numerics_integrator == 1_pInt)) then - allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase)) - allocate(plasticState(phase)%previousDotState2(sizeDotState,NofMyPhase)) + allocate(plasticState(phase)%previousDotState (0_pInt,NofMyPhase)) + allocate(plasticState(phase)%previousDotState2(0_pInt,NofMyPhase)) endif if (any(numerics_integrator == 4_pInt)) & - allocate(plasticState(phase)%RK4dotState (sizeDotState,NofMyPhase)) + allocate(plasticState(phase)%RK4dotState (0_pInt,NofMyPhase)) if (any(numerics_integrator == 5_pInt)) & - allocate(plasticState(phase)%RKCK45dotState (6,sizeDotState,NofMyPhase)) + allocate(plasticState(phase)%RKCK45dotState (6,0_pInt,NofMyPhase)) endif enddo initializeInstances - allocate(plastic_none_sizePostResults(maxNinstance), source=0_pInt) - end subroutine plastic_none_init end module plastic_none diff --git a/src/plastic_phenopowerlaw.f90 b/src/plastic_phenopowerlaw.f90 index 8a6d8b145..bdc6e12a6 100644 --- a/src/plastic_phenopowerlaw.f90 +++ b/src/plastic_phenopowerlaw.f90 @@ -1,5 +1,7 @@ +!-------------------------------------------------------------------------------------------------- !> @author Franz Roters, Max-Planck-Institut für Eisenforschung GmbH !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH +!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH !> @brief material subroutine for phenomenological crystal plasticity formulation using a powerlaw !! fitting !-------------------------------------------------------------------------------------------------- @@ -10,91 +12,75 @@ module plastic_phenopowerlaw implicit none private - integer(pInt), dimension(:), allocatable, public, protected :: & - plastic_phenopowerlaw_sizePostResults !< cumulative size of post results - integer(pInt), dimension(:,:), allocatable, target, public :: & - plastic_phenopowerlaw_sizePostResult !< size of each post result output - + plastic_phenopowerlaw_sizePostResult !< size of each post result output character(len=64), dimension(:,:), allocatable, target, public :: & - plastic_phenopowerlaw_output !< name of each post result output - - integer(pInt), dimension(:), allocatable, target, public :: & - plastic_phenopowerlaw_Noutput !< number of outputs per instance of this constitution - - integer(pInt), dimension(:), allocatable, public, protected :: & - plastic_phenopowerlaw_totalNslip, & !< no. of slip system used in simulation - plastic_phenopowerlaw_totalNtwin, & !< no. of twin system used in simulation - plastic_phenopowerlaw_totalNtrans !< no. of trans system used in simulation - - integer(pInt), dimension(:,:), allocatable, private :: & - plastic_phenopowerlaw_Nslip, & !< active number of slip systems per family (input parameter, per family) - plastic_phenopowerlaw_Ntwin, & !< active number of twin systems per family (input parameter, per family) - plastic_phenopowerlaw_Ntrans !< active number of trans systems per family (input parameter, per family) - - real(pReal), dimension(:), allocatable, private :: & - plastic_phenopowerlaw_gdot0_slip, & !< reference shear strain rate for slip (input parameter) - plastic_phenopowerlaw_gdot0_twin, & !< reference shear strain rate for twin (input parameter) - plastic_phenopowerlaw_n_slip, & !< stress exponent for slip (input parameter) - plastic_phenopowerlaw_n_twin, & !< stress exponent for twin (input parameter) - plastic_phenopowerlaw_spr, & !< push-up factor for slip saturation due to twinning - plastic_phenopowerlaw_twinB, & - plastic_phenopowerlaw_twinC, & - plastic_phenopowerlaw_twinD, & - plastic_phenopowerlaw_twinE, & - plastic_phenopowerlaw_h0_SlipSlip, & !< reference hardening slip - slip (input parameter) - plastic_phenopowerlaw_h0_TwinSlip, & !< reference hardening twin - slip (input parameter) - plastic_phenopowerlaw_h0_TwinTwin, & !< reference hardening twin - twin (input parameter) - plastic_phenopowerlaw_a_slip, & - plastic_phenopowerlaw_aTolResistance, & - plastic_phenopowerlaw_aTolShear, & - plastic_phenopowerlaw_aTolTwinfrac, & - plastic_phenopowerlaw_aTolTransfrac, & - plastic_phenopowerlaw_Cnuc, & !< coefficient for strain-induced martensite nucleation - plastic_phenopowerlaw_Cdwp, & !< coefficient for double well potential - plastic_phenopowerlaw_Cgro, & !< coefficient for stress-assisted martensite growth - plastic_phenopowerlaw_deltaG !< free energy difference between austensite and martensite [MPa] - - real(pReal), dimension(:,:), allocatable, private :: & - plastic_phenopowerlaw_tau0_slip, & !< initial critical shear stress for slip (input parameter, per family) - plastic_phenopowerlaw_tau0_twin, & !< initial critical shear stress for twin (input parameter, per family) - plastic_phenopowerlaw_tausat_slip, & !< maximum critical shear stress for slip (input parameter, per family) - plastic_phenopowerlaw_H_int, & !< per family hardening activity(input parameter(optional), per family) - plastic_phenopowerlaw_nonSchmidCoeff, & - - plastic_phenopowerlaw_interaction_SlipSlip, & !< interaction factors slip - slip (input parameter) - plastic_phenopowerlaw_interaction_SlipTwin, & !< interaction factors slip - twin (input parameter) - plastic_phenopowerlaw_interaction_TwinSlip, & !< interaction factors twin - slip (input parameter) - plastic_phenopowerlaw_interaction_TwinTwin !< interaction factors twin - twin (input parameter) - - real(pReal), dimension(:,:,:), allocatable, private :: & - plastic_phenopowerlaw_hardeningMatrix_SlipSlip, & - plastic_phenopowerlaw_hardeningMatrix_SlipTwin, & - plastic_phenopowerlaw_hardeningMatrix_TwinSlip, & - plastic_phenopowerlaw_hardeningMatrix_TwinTwin + plastic_phenopowerlaw_output !< name of each post result output enum, bind(c) - enumerator :: undefined_ID, & - resistance_slip_ID, & - accumulatedshear_slip_ID, & - shearrate_slip_ID, & - resolvedstress_slip_ID, & - totalshear_ID, & - resistance_twin_ID, & - accumulatedshear_twin_ID, & - shearrate_twin_ID, & - resolvedstress_twin_ID, & - totalvolfrac_twin_ID + enumerator :: & + undefined_ID, & + resistance_slip_ID, & + accumulatedshear_slip_ID, & + shearrate_slip_ID, & + resolvedstress_slip_ID, & + totalshear_ID, & + resistance_twin_ID, & + accumulatedshear_twin_ID, & + shearrate_twin_ID, & + resolvedstress_twin_ID, & + totalvolfrac_twin_ID end enum - integer(kind(undefined_ID)), dimension(:,:), allocatable, private :: & - plastic_phenopowerlaw_outputID !< ID of each post result output + + type, private :: tParameters !< container type for internal constitutive parameters + integer(pInt) :: & + totalNslip, & + totalNtwin + real(pReal) :: & + gdot0_slip, & !< reference shear strain rate for slip + gdot0_twin, & !< reference shear strain rate for twin + n_slip, & !< stress exponent for slip + n_twin, & !< stress exponent for twin + spr, & !< push-up factor for slip saturation due to twinning + twinB, & + twinC, & + twinD, & + twinE, & + h0_SlipSlip, & !< reference hardening slip - slip + h0_TwinSlip, & !< reference hardening twin - slip + h0_TwinTwin, & !< reference hardening twin - twin + a_slip, & + aTolResistance, & ! default absolute tolerance 1 Pa + aTolShear, & ! default absolute tolerance 1e-6 + aTolTwinfrac ! default absolute tolerance 1e-6 + integer(pInt), dimension(:), allocatable :: & + Nslip, & !< active number of slip systems per family + Ntwin !< active number of twin systems per family + real(pReal), dimension(:), allocatable :: & + tau0_slip, & !< initial critical shear stress for slip + tau0_twin, & !< initial critical shear stress for twin + tausat_slip, & !< maximum critical shear stress for slip + nonSchmidCoeff, & + H_int !< per family hardening activity (optional) + real(pReal), dimension(:,:), allocatable :: & + interaction_SlipSlip, & !< slip resistance from slip activity + interaction_SlipTwin, & !< slip resistance from twin activity + interaction_TwinSlip, & !< twin resistance from slip activity + interaction_TwinTwin !< twin resistance from twin activity + + integer(kind(undefined_ID)), dimension(:), allocatable :: & + outputID !< ID of each post result output + end type + + type(tParameters), dimension(:), allocatable, private :: param !< containers of constitutive parameters (len Ninstance) type, private :: tPhenopowerlawState real(pReal), pointer, dimension(:,:) :: & s_slip, & s_twin, & accshear_slip, & - accshear_twin + accshear_twin, & + whole real(pReal), pointer, dimension(:) :: & sumGamma, & sumF @@ -102,18 +88,13 @@ module plastic_phenopowerlaw type(tPhenopowerlawState), allocatable, dimension(:), private :: & dotState, & - state, & - state0 + state public :: & plastic_phenopowerlaw_init, & plastic_phenopowerlaw_LpAndItsTangent, & plastic_phenopowerlaw_dotState, & plastic_phenopowerlaw_postResults - private :: & - plastic_phenopowerlaw_aTolState, & - plastic_phenopowerlaw_stateInit - contains @@ -122,7 +103,7 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine plastic_phenopowerlaw_init(fileUnit) +subroutine plastic_phenopowerlaw_init #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800 use, intrinsic :: iso_fortran_env, only: & compiler_version, & @@ -136,20 +117,12 @@ subroutine plastic_phenopowerlaw_init(fileUnit) debug_levelBasic use math, only: & math_Mandel3333to66, & - math_Voigt66to3333 + math_Voigt66to3333, & + math_expand use IO, only: & - IO_read, & - IO_lc, & - IO_getTag, & - IO_isBlank, & - IO_stringPos, & - IO_stringValue, & - IO_floatValue, & - IO_intValue, & IO_warning, & IO_error, & - IO_timeStamp, & - IO_EOF + IO_timeStamp use material, only: & phase_plasticity, & phase_plasticityInstance, & @@ -159,615 +132,338 @@ subroutine plastic_phenopowerlaw_init(fileUnit) material_phase, & plasticState use config, only: & - MATERIAL_partPhase + MATERIAL_partPhase, & + config_phase use lattice use numerics,only: & numerics_integrator implicit none - integer(pInt), intent(in) :: fileUnit - integer(pInt), allocatable, dimension(:) :: chunkPos integer(pInt) :: & maxNinstance, & - instance,phase,j,k, f,o, & - Nchunks_SlipSlip = 0_pInt, Nchunks_SlipTwin = 0_pInt, & - Nchunks_TwinSlip = 0_pInt, Nchunks_TwinTwin = 0_pInt, & - Nchunks_SlipFamilies = 0_pInt, Nchunks_TwinFamilies = 0_pInt, & - Nchunks_TransFamilies = 0_pInt, Nchunks_nonSchmid = 0_pInt, & - NipcMyPhase, & - offset_slip, index_myFamily, index_otherFamily, & - mySize=0_pInt,sizeState,sizeDotState, sizeDeltaState, & + instance,p,j,k, f,o, i,& + NipcMyPhase, outputSize, & + index_myFamily, index_otherFamily, & + sizeState,sizeDotState, & startIndex, endIndex - character(len=65536) :: & - tag = '', & - line = '' - real(pReal), dimension(:), allocatable :: tempPerSlip + + real(pReal), dimension(:,:), allocatable :: temp1, temp2 + + integer(pInt), dimension(0), parameter :: emptyIntArray = [integer(pInt)::] + real(pReal), dimension(0), parameter :: emptyRealArray = [real(pReal)::] + character(len=65536), dimension(0), parameter :: emptyStringArray = [character(len=65536)::] + + type(tParameters) :: prm + + integer(kind(undefined_ID)) :: & + outputID !< ID of each post result output + + character(len=512) :: & + extmsg = '' + character(len=65536), dimension(:), allocatable :: outputs write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_PHENOPOWERLAW_label//' init -+>>>' write(6,'(a15,a)') ' Current time: ',IO_timeStamp() #include "compilation_info.f90" maxNinstance = int(count(phase_plasticity == PLASTICITY_PHENOPOWERLAW_ID),pInt) - if (maxNinstance == 0_pInt) return - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) & write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance - - allocate(plastic_phenopowerlaw_sizePostResults(maxNinstance), source=0_pInt) - allocate(plastic_phenopowerlaw_sizePostResult(maxval(phase_Noutput),maxNinstance), & - source=0_pInt) + allocate(plastic_phenopowerlaw_sizePostResult(maxval(phase_Noutput),maxNinstance),source=0_pInt) allocate(plastic_phenopowerlaw_output(maxval(phase_Noutput),maxNinstance)) - plastic_phenopowerlaw_output = '' - allocate(plastic_phenopowerlaw_outputID(maxval(phase_Noutput),maxNinstance),source=undefined_ID) - allocate(plastic_phenopowerlaw_Noutput(maxNinstance), source=0_pInt) - allocate(plastic_phenopowerlaw_Nslip(lattice_maxNslipFamily,maxNinstance), source=0_pInt) - allocate(plastic_phenopowerlaw_Ntwin(lattice_maxNtwinFamily,maxNinstance), source=0_pInt) - allocate(plastic_phenopowerlaw_Ntrans(lattice_maxNtransFamily,maxNinstance),source=0_pInt) - allocate(plastic_phenopowerlaw_totalNslip(maxNinstance), source=0_pInt) - allocate(plastic_phenopowerlaw_totalNtwin(maxNinstance), source=0_pInt) - allocate(plastic_phenopowerlaw_totalNtrans(maxNinstance), source=0_pInt) - allocate(plastic_phenopowerlaw_gdot0_slip(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_n_slip(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_tau0_slip(lattice_maxNslipFamily,maxNinstance),source=0.0_pReal) - allocate(plastic_phenopowerlaw_tausat_slip(lattice_maxNslipFamily,maxNinstance),source=0.0_pReal) - allocate(plastic_phenopowerlaw_H_int(lattice_maxNslipFamily,maxNinstance),source=0.0_pReal) - allocate(plastic_phenopowerlaw_gdot0_twin(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_n_twin(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_tau0_twin(lattice_maxNtwinFamily,maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_spr(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_twinB(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_twinC(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_twinD(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_twinE(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_h0_SlipSlip(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_h0_TwinSlip(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_h0_TwinTwin(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_interaction_SlipSlip(lattice_maxNinteraction,maxNinstance), & - source=0.0_pReal) - allocate(plastic_phenopowerlaw_interaction_SlipTwin(lattice_maxNinteraction,maxNinstance), & - source=0.0_pReal) - allocate(plastic_phenopowerlaw_interaction_TwinSlip(lattice_maxNinteraction,maxNinstance), & - source=0.0_pReal) - allocate(plastic_phenopowerlaw_interaction_TwinTwin(lattice_maxNinteraction,maxNinstance), & - source=0.0_pReal) - allocate(plastic_phenopowerlaw_a_slip(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_aTolResistance(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_aTolShear(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_aTolTwinfrac(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_aTolTransfrac(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_nonSchmidCoeff(lattice_maxNnonSchmid,maxNinstance), & - source=0.0_pReal) - allocate(plastic_phenopowerlaw_Cnuc(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_Cdwp(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_Cgro(maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_deltaG(maxNinstance), source=0.0_pReal) + plastic_phenopowerlaw_output = '' - rewind(fileUnit) - phase = 0_pInt - do while (trim(line) /= IO_EOF .and. IO_lc(IO_getTag(line,'<','>')) /= material_partPhase) ! wind forward to - line = IO_read(fileUnit) - enddo - - parsingFile: do while (trim(line) /= IO_EOF) ! read through sections of phase part - line = IO_read(fileUnit) - if (IO_isBlank(line)) cycle ! skip empty lines - if (IO_getTag(line,'<','>') /= '') then ! stop at next part - line = IO_read(fileUnit, .true.) ! reset IO_read - exit - endif - if (IO_getTag(line,'[',']') /= '') then ! next phase - phase = phase + 1_pInt ! advance phase section counter - if (phase_plasticity(phase) == PLASTICITY_PHENOPOWERLAW_ID) then - Nchunks_SlipFamilies = count(lattice_NslipSystem(:,phase) > 0_pInt) ! maximum number of slip families according to lattice type of current phase - Nchunks_TwinFamilies = count(lattice_NtwinSystem(:,phase) > 0_pInt) ! maximum number of twin families according to lattice type of current phase - Nchunks_TransFamilies = count(lattice_NtransSystem(:,phase) > 0_pInt) ! maximum number of trans families according to lattice type of current phase - Nchunks_SlipSlip = maxval(lattice_interactionSlipSlip(:,:,phase)) - Nchunks_SlipTwin = maxval(lattice_interactionSlipTwin(:,:,phase)) - Nchunks_TwinSlip = maxval(lattice_interactionTwinSlip(:,:,phase)) - Nchunks_TwinTwin = maxval(lattice_interactionTwinTwin(:,:,phase)) - Nchunks_nonSchmid = lattice_NnonSchmid(phase) - if(allocated(tempPerSlip)) deallocate(tempPerSlip) - allocate(tempPerSlip(Nchunks_SlipFamilies)) - endif - cycle ! skip to next line - endif - if (phase > 0_pInt ) then; if (phase_plasticity(phase) == PLASTICITY_PHENOPOWERLAW_ID) then ! one of my phases. Do not short-circuit here (.and. between if-statements), it's not safe in Fortran - instance = phase_plasticityInstance(phase) ! which instance of my plasticity is present phase - chunkPos = IO_stringPos(line) - tag = IO_lc(IO_stringValue(line,chunkPos,1_pInt)) ! extract key - select case(tag) - case ('(output)') - select case(IO_lc(IO_stringValue(line,chunkPos,2_pInt))) - case ('resistance_slip') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = resistance_slip_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('accumulatedshear_slip','accumulated_shear_slip') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = accumulatedshear_slip_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('shearrate_slip') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = shearrate_slip_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('resolvedstress_slip') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = resolvedstress_slip_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('totalshear') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = totalshear_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('resistance_twin') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = resistance_twin_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('accumulatedshear_twin','accumulated_shear_twin') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = accumulatedshear_twin_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('shearrate_twin') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = shearrate_twin_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('resolvedstress_twin') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = resolvedstress_twin_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('totalvolfrac_twin') - plastic_phenopowerlaw_Noutput(instance) = plastic_phenopowerlaw_Noutput(instance) + 1_pInt - plastic_phenopowerlaw_outputID(plastic_phenopowerlaw_Noutput(instance),instance) = totalvolfrac_twin_ID - plastic_phenopowerlaw_output(plastic_phenopowerlaw_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case default - - end select -!-------------------------------------------------------------------------------------------------- -! parameters depending on number of slip families - case ('nslip') - if (chunkPos(1) < Nchunks_SlipFamilies + 1_pInt) & - call IO_warning(50_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - if (chunkPos(1) > Nchunks_SlipFamilies + 1_pInt) & - call IO_error(150_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - Nchunks_SlipFamilies = chunkPos(1) - 1_pInt ! user specified number of (possibly) active slip families (e.g. 6 0 6 --> 3) - do j = 1_pInt, Nchunks_SlipFamilies - plastic_phenopowerlaw_Nslip(j,instance) = IO_intValue(line,chunkPos,1_pInt+j) - enddo - case ('tausat_slip','tau0_slip','H_int') - tempPerSlip = 0.0_pReal - do j = 1_pInt, Nchunks_SlipFamilies - if (plastic_phenopowerlaw_Nslip(j,instance) > 0_pInt) & - tempPerSlip(j) = IO_floatValue(line,chunkPos,1_pInt+j) - enddo - select case(tag) - case ('tausat_slip') - plastic_phenopowerlaw_tausat_slip(1:Nchunks_SlipFamilies,instance) = tempPerSlip(1:Nchunks_SlipFamilies) - case ('tau0_slip') - plastic_phenopowerlaw_tau0_slip(1:Nchunks_SlipFamilies,instance) = tempPerSlip(1:Nchunks_SlipFamilies) - case ('H_int') - plastic_phenopowerlaw_H_int(1:Nchunks_SlipFamilies,instance) = tempPerSlip(1:Nchunks_SlipFamilies) - end select -!-------------------------------------------------------------------------------------------------- -! parameters depending on number of twin families - case ('ntwin') - if (chunkPos(1) < Nchunks_TwinFamilies + 1_pInt) & - call IO_warning(51_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - if (chunkPos(1) > Nchunks_TwinFamilies + 1_pInt) & - call IO_error(150_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - Nchunks_TwinFamilies = chunkPos(1) - 1_pInt - do j = 1_pInt, Nchunks_TwinFamilies - plastic_phenopowerlaw_Ntwin(j,instance) = IO_intValue(line,chunkPos,1_pInt+j) - enddo - case ('tau0_twin') - do j = 1_pInt, Nchunks_TwinFamilies - if (plastic_phenopowerlaw_Ntwin(j,instance) > 0_pInt) & - plastic_phenopowerlaw_tau0_twin(j,instance) = IO_floatValue(line,chunkPos,1_pInt+j) - enddo -!-------------------------------------------------------------------------------------------------- -! parameters depending on number of transformation families - case ('ntrans') - if (chunkPos(1) < Nchunks_TransFamilies + 1_pInt) & - call IO_warning(53_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - if (chunkPos(1) > Nchunks_TransFamilies + 1_pInt) & - call IO_error(150_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - Nchunks_TransFamilies = chunkPos(1) - 1_pInt - do j = 1_pInt, Nchunks_TransFamilies - plastic_phenopowerlaw_Ntrans(j,instance) = IO_intValue(line,chunkPos,1_pInt+j) - enddo -!-------------------------------------------------------------------------------------------------- -! parameters depending on number of interactions - case ('interaction_slipslip') - if (chunkPos(1) < 1_pInt + Nchunks_SlipSlip) & - call IO_warning(52_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - do j = 1_pInt, Nchunks_SlipSlip - plastic_phenopowerlaw_interaction_SlipSlip(j,instance) = IO_floatValue(line,chunkPos,1_pInt+j) - enddo - case ('interaction_sliptwin') - if (chunkPos(1) < 1_pInt + Nchunks_SlipTwin) & - call IO_warning(52_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - do j = 1_pInt, Nchunks_SlipTwin - plastic_phenopowerlaw_interaction_SlipTwin(j,instance) = IO_floatValue(line,chunkPos,1_pInt+j) - enddo - case ('interaction_twinslip') - if (chunkPos(1) < 1_pInt + Nchunks_TwinSlip) & - call IO_warning(52_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - do j = 1_pInt, Nchunks_TwinSlip - plastic_phenopowerlaw_interaction_TwinSlip(j,instance) = IO_floatValue(line,chunkPos,1_pInt+j) - enddo - case ('interaction_twintwin') - if (chunkPos(1) < 1_pInt + Nchunks_TwinTwin) & - call IO_warning(52_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - do j = 1_pInt, Nchunks_TwinTwin - plastic_phenopowerlaw_interaction_TwinTwin(j,instance) = IO_floatValue(line,chunkPos,1_pInt+j) - enddo - case ('nonschmid_coefficients') - if (chunkPos(1) < 1_pInt + Nchunks_nonSchmid) & - call IO_warning(52_pInt,ext_msg=trim(tag)//' ('//PLASTICITY_PHENOPOWERLAW_label//')') - do j = 1_pInt,Nchunks_nonSchmid - plastic_phenopowerlaw_nonSchmidCoeff(j,instance) = IO_floatValue(line,chunkPos,1_pInt+j) - enddo -!-------------------------------------------------------------------------------------------------- -! parameters independent of number of slip/twin systems - case ('gdot0_slip') - plastic_phenopowerlaw_gdot0_slip(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('n_slip') - plastic_phenopowerlaw_n_slip(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('a_slip', 'w0_slip') - plastic_phenopowerlaw_a_slip(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('gdot0_twin') - plastic_phenopowerlaw_gdot0_twin(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('n_twin') - plastic_phenopowerlaw_n_twin(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('s_pr') - plastic_phenopowerlaw_spr(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('twin_b') - plastic_phenopowerlaw_twinB(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('twin_c') - plastic_phenopowerlaw_twinC(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('twin_d') - plastic_phenopowerlaw_twinD(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('twin_e') - plastic_phenopowerlaw_twinE(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('h0_slipslip') - plastic_phenopowerlaw_h0_SlipSlip(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('h0_twinslip') - plastic_phenopowerlaw_h0_TwinSlip(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('h0_twintwin') - plastic_phenopowerlaw_h0_TwinTwin(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('atol_resistance') - plastic_phenopowerlaw_aTolResistance(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('atol_shear') - plastic_phenopowerlaw_aTolShear(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('atol_twinfrac') - plastic_phenopowerlaw_aTolTwinfrac(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('atol_transfrac') - plastic_phenopowerlaw_aTolTransfrac(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('cnuc') - plastic_phenopowerlaw_Cnuc(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('cdwp') - plastic_phenopowerlaw_Cdwp(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('cgro') - plastic_phenopowerlaw_Cgro(instance) = IO_floatValue(line,chunkPos,2_pInt) - case ('deltag') - plastic_phenopowerlaw_deltaG(instance) = IO_floatValue(line,chunkPos,2_pInt) - case default - - end select - endif; endif - enddo parsingFile - - sanityChecks: do phase = 1_pInt, size(phase_plasticity) - myPhase: if (phase_plasticity(phase) == PLASTICITY_phenopowerlaw_ID) then - instance = phase_plasticityInstance(phase) - plastic_phenopowerlaw_Nslip(1:lattice_maxNslipFamily,instance) = & - min(lattice_NslipSystem(1:lattice_maxNslipFamily,phase),& ! limit active slip systems per family to min of available and requested - plastic_phenopowerlaw_Nslip(1:lattice_maxNslipFamily,instance)) - plastic_phenopowerlaw_Ntwin(1:lattice_maxNtwinFamily,instance) = & - min(lattice_NtwinSystem(1:lattice_maxNtwinFamily,phase),& ! limit active twin systems per family to min of available and requested - plastic_phenopowerlaw_Ntwin(:,instance)) - plastic_phenopowerlaw_totalNslip(instance) = sum(plastic_phenopowerlaw_Nslip(:,instance)) ! how many slip systems altogether - plastic_phenopowerlaw_totalNtwin(instance) = sum(plastic_phenopowerlaw_Ntwin(:,instance)) ! how many twin systems altogether - plastic_phenopowerlaw_totalNtrans(instance) = sum(plastic_phenopowerlaw_Ntrans(:,instance)) ! how many trans systems altogether - - if (any(plastic_phenopowerlaw_tau0_slip(:,instance) < 0.0_pReal .and. & - plastic_phenopowerlaw_Nslip(:,instance) > 0)) & - call IO_error(211_pInt,el=instance,ext_msg='tau0_slip ('//PLASTICITY_PHENOPOWERLAW_label//')') - if (plastic_phenopowerlaw_gdot0_slip(instance) <= 0.0_pReal) & - call IO_error(211_pInt,el=instance,ext_msg='gdot0_slip ('//PLASTICITY_PHENOPOWERLAW_label//')') - if (plastic_phenopowerlaw_n_slip(instance) <= 0.0_pReal) & - call IO_error(211_pInt,el=instance,ext_msg='n_slip ('//PLASTICITY_PHENOPOWERLAW_label//')') - if (any(plastic_phenopowerlaw_tausat_slip(:,instance) <= 0.0_pReal .and. & - plastic_phenopowerlaw_Nslip(:,instance) > 0)) & - call IO_error(211_pInt,el=instance,ext_msg='tausat_slip ('//PLASTICITY_PHENOPOWERLAW_label//')') - if (any(dEq0(plastic_phenopowerlaw_a_slip(instance)) .and. plastic_phenopowerlaw_Nslip(:,instance) > 0)) & - call IO_error(211_pInt,el=instance,ext_msg='a_slip ('//PLASTICITY_PHENOPOWERLAW_label//')') - if (any(plastic_phenopowerlaw_tau0_twin(:,instance) < 0.0_pReal .and. & - plastic_phenopowerlaw_Ntwin(:,instance) > 0)) & - call IO_error(211_pInt,el=instance,ext_msg='tau0_twin ('//PLASTICITY_PHENOPOWERLAW_label//')') - if ( plastic_phenopowerlaw_gdot0_twin(instance) <= 0.0_pReal .and. & - any(plastic_phenopowerlaw_Ntwin(:,instance) > 0)) & - call IO_error(211_pInt,el=instance,ext_msg='gdot0_twin ('//PLASTICITY_PHENOPOWERLAW_label//')') - if ( plastic_phenopowerlaw_n_twin(instance) <= 0.0_pReal .and. & - any(plastic_phenopowerlaw_Ntwin(:,instance) > 0)) & - call IO_error(211_pInt,el=instance,ext_msg='n_twin ('//PLASTICITY_PHENOPOWERLAW_label//')') - if (plastic_phenopowerlaw_aTolResistance(instance) <= 0.0_pReal) & - plastic_phenopowerlaw_aTolResistance(instance) = 1.0_pReal ! default absolute tolerance 1 Pa - if (plastic_phenopowerlaw_aTolShear(instance) <= 0.0_pReal) & - plastic_phenopowerlaw_aTolShear(instance) = 1.0e-6_pReal ! default absolute tolerance 1e-6 - if (plastic_phenopowerlaw_aTolTwinfrac(instance) <= 0.0_pReal) & - plastic_phenopowerlaw_aTolTwinfrac(instance) = 1.0e-6_pReal ! default absolute tolerance 1e-6 - if (plastic_phenopowerlaw_aTolTransfrac(instance) <= 0.0_pReal) & - plastic_phenopowerlaw_aTolTransfrac(instance) = 1.0e-6_pReal ! default absolute tolerance 1e-6 - endif myPhase - enddo sanityChecks - -!-------------------------------------------------------------------------------------------------- -! allocation of variables whose size depends on the total number of active slip systems - allocate(plastic_phenopowerlaw_hardeningMatrix_SlipSlip(maxval(plastic_phenopowerlaw_totalNslip),& ! slip resistance from slip activity - maxval(plastic_phenopowerlaw_totalNslip),& - maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_hardeningMatrix_SlipTwin(maxval(plastic_phenopowerlaw_totalNslip),& ! slip resistance from twin activity - maxval(plastic_phenopowerlaw_totalNtwin),& - maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_hardeningMatrix_TwinSlip(maxval(plastic_phenopowerlaw_totalNtwin),& ! twin resistance from slip activity - maxval(plastic_phenopowerlaw_totalNslip),& - maxNinstance), source=0.0_pReal) - allocate(plastic_phenopowerlaw_hardeningMatrix_TwinTwin(maxval(plastic_phenopowerlaw_totalNtwin),& ! twin resistance from twin activity - maxval(plastic_phenopowerlaw_totalNtwin),& - maxNinstance), source=0.0_pReal) + allocate(param(maxNinstance)) ! one container of parameters per instance allocate(state(maxNinstance)) - allocate(state0(maxNinstance)) allocate(dotState(maxNinstance)) - initializeInstances: do phase = 1_pInt, size(phase_plasticity) ! loop through all phases in material.config - myPhase2: if (phase_plasticity(phase) == PLASTICITY_phenopowerlaw_ID) then ! only consider my phase - NipcMyPhase = count(material_phase == phase) ! number of IPCs containing my phase - instance = phase_plasticityInstance(phase) ! which instance of my phase + do p = 1_pInt, size(phase_plasticityInstance) + if (phase_plasticity(p) /= PLASTICITY_PHENOPOWERLAW_ID) cycle + instance = phase_plasticityInstance(p) + associate(prm => param(instance)) -!-------------------------------------------------------------------------------------------------- -! Determine size of postResults array - outputsLoop: do o = 1_pInt,plastic_phenopowerlaw_Noutput(instance) - select case(plastic_phenopowerlaw_outputID(o,instance)) - case(resistance_slip_ID, & - shearrate_slip_ID, & - accumulatedshear_slip_ID, & - resolvedstress_slip_ID & - ) - mySize = plastic_phenopowerlaw_totalNslip(instance) - case(resistance_twin_ID, & - shearrate_twin_ID, & - accumulatedshear_twin_ID, & - resolvedstress_twin_ID & - ) - mySize = plastic_phenopowerlaw_totalNtwin(instance) - case(totalshear_ID, & - totalvolfrac_twin_ID & - ) - mySize = 1_pInt - case default - end select + prm%Nslip = config_phase(p)%getInts('nslip',defaultVal=emptyIntArray) + if (size(prm%Nslip) > count(lattice_NslipSystem(:,p) > 0_pInt)) call IO_error(150_pInt,ext_msg='Nslip') + if (any(lattice_NslipSystem(1:size(prm%Nslip),p)-prm%Nslip < 0_pInt)) call IO_error(150_pInt,ext_msg='Nslip') + prm%totalNslip = sum(prm%Nslip) + + if (prm%totalNslip > 0_pInt) then + prm%tau0_slip = config_phase(p)%getFloats('tau0_slip') + prm%tausat_slip = config_phase(p)%getFloats('tausat_slip') + prm%interaction_SlipSlip = spread(config_phase(p)%getFloats('interaction_slipslip'),2,1) + prm%H_int = config_phase(p)%getFloats('h_int',& + defaultVal=[(0.0_pReal,i=1_pInt,size(prm%Nslip))]) + prm%nonSchmidCoeff = config_phase(p)%getFloats('nonschmid_coefficients',& + defaultVal = emptyRealArray ) + + prm%gdot0_slip = config_phase(p)%getFloat('gdot0_slip') + prm%n_slip = config_phase(p)%getFloat('n_slip') + prm%a_slip = config_phase(p)%getFloat('a_slip') + prm%h0_SlipSlip = config_phase(p)%getFloat('h0_slipslip') + endif + + prm%Ntwin = config_phase(p)%getInts('ntwin', defaultVal=emptyIntArray) + if (size(prm%Ntwin) > count(lattice_NtwinSystem(:,p) > 0_pInt)) call IO_error(150_pInt,ext_msg='Ntwin') + if (any(lattice_NtwinSystem(1:size(prm%Ntwin),p)-prm%Ntwin < 0_pInt)) call IO_error(150_pInt,ext_msg='Ntwin') + prm%totalNtwin = sum(prm%Ntwin) + + if (prm%totalNtwin > 0_pInt) then + prm%tau0_twin = config_phase(p)%getFloats('tau0_twin') + prm%interaction_TwinTwin = spread(config_phase(p)%getFloats('interaction_twintwin'),2,1) + + prm%gdot0_twin = config_phase(p)%getFloat('gdot0_twin') + prm%n_twin = config_phase(p)%getFloat('n_twin') + prm%spr = config_phase(p)%getFloat('s_pr') + prm%twinB = config_phase(p)%getFloat('twin_b') + prm%twinC = config_phase(p)%getFloat('twin_c') + prm%twinD = config_phase(p)%getFloat('twin_d') + prm%twinE = config_phase(p)%getFloat('twin_e') + prm%h0_TwinTwin = config_phase(p)%getFloat('h0_twintwin') + endif + + if (prm%totalNslip > 0_pInt .and. prm%totalNtwin > 0_pInt) then + prm%interaction_SlipTwin = spread(config_phase(p)%getFloats('interaction_sliptwin'),2,1) + prm%interaction_TwinSlip = spread(config_phase(p)%getFloats('interaction_twinslip'),2,1) + prm%h0_TwinSlip = config_phase(p)%getFloat('h0_twinslip') + endif + + + prm%aTolResistance = config_phase(p)%getFloat('atol_resistance',defaultVal=1.0_pReal) + prm%aTolShear = config_phase(p)%getFloat('atol_shear',defaultVal=1.0e-6_pReal) + prm%aTolTwinfrac = config_phase(p)%getFloat('atol_twinfrac',defaultVal=1.0e-6_pReal) + + outputs = config_phase(p)%getStrings('(output)',defaultVal=emptyStringArray) + allocate(prm%outputID(0)) + do i=1_pInt, size(outputs) + outputID = undefined_ID + select case(outputs(i)) + case ('resistance_slip') + outputID = resistance_slip_ID + outputSize = sum(prm%Nslip) + case ('accumulatedshear_slip') + outputID = accumulatedshear_slip_ID + outputSize = sum(prm%Nslip) + case ('shearrate_slip') + outputID = shearrate_slip_ID + outputSize = sum(prm%Nslip) + case ('resolvedstress_slip') + outputID = resolvedstress_slip_ID + outputSize = sum(prm%Nslip) + + case ('resistance_twin') + outputID = resistance_twin_ID + outputSize = sum(prm%Ntwin) + case ('accumulatedshear_twin') + outputID = accumulatedshear_twin_ID + outputSize = sum(prm%Ntwin) + case ('shearrate_twin') + outputID = shearrate_twin_ID + outputSize = sum(prm%Ntwin) + case ('resolvedstress_twin') + outputID = resolvedstress_twin_ID + outputSize = sum(prm%Ntwin) + + case ('totalvolfrac_twin') + outputID = totalvolfrac_twin_ID + outputSize = 1_pInt + case ('totalshear') + outputID = totalshear_ID + outputSize = 1_pInt + end select + + if (outputID /= undefined_ID) then + plastic_phenopowerlaw_output(i,instance) = outputs(i) + plastic_phenopowerlaw_sizePostResult(i,instance) = outputSize + prm%outputID = [prm%outputID , outputID] + endif + + end do + + extmsg = '' + if (sum(prm%Nslip) > 0_pInt) then + if (size(prm%tau0_slip) /= size(prm%Nslip)) call IO_error(211_pInt,ip=instance, & + ext_msg='shape(tau0_slip) ('//PLASTICITY_PHENOPOWERLAW_label//')') + if (size(prm%tausat_slip) /= size(prm%Nslip)) call IO_error(211_pInt,ip=instance, & + ext_msg='shape(tausat_slip) ('//PLASTICITY_PHENOPOWERLAW_label//')') + if (size(prm%H_int) /= size(prm%Nslip)) call IO_error(211_pInt,ip=instance, & + ext_msg='shape(H_int) ('//PLASTICITY_PHENOPOWERLAW_label//')') + + if (any(prm%tau0_slip < 0.0_pReal .and. prm%Nslip > 0_pInt)) & + extmsg = trim(extmsg)//"tau0_slip " + if (any(prm%tausat_slip < prm%tau0_slip .and. prm%Nslip > 0_pInt)) & + extmsg = trim(extmsg)//"tausat_slip " + + if (prm%gdot0_slip <= 0.0_pReal) extmsg = trim(extmsg)//" gdot0_slip " + if (dEq0(prm%a_slip)) extmsg = trim(extmsg)//" a_slip " ! ToDo: negative values ok? + if (dEq0(prm%n_slip)) extmsg = trim(extmsg)//" n_slip " ! ToDo: negative values ok? + endif + + if (sum(prm%Ntwin) > 0_pInt) then + if (size(prm%tau0_twin) /= size(prm%ntwin)) call IO_error(211_pInt,ip=instance,& + ext_msg='shape(tau0_twin) ('//PLASTICITY_PHENOPOWERLAW_label//')') + + if (any(prm%tau0_twin < 0.0_pReal .and. prm%Ntwin > 0_pInt)) & + extmsg = trim(extmsg)//"tau0_twin " + + if (prm%gdot0_twin <= 0.0_pReal) extmsg = trim(extmsg)//"gdot0_twin " + if (dEq0(prm%n_twin)) extmsg = trim(extmsg)//"n_twin " ! ToDo: negative values ok? + endif + + if (prm%aTolResistance <= 0.0_pReal) extmsg = trim(extmsg)//"aTolresistance " + if (prm%aTolShear <= 0.0_pReal) extmsg = trim(extmsg)//"aTolShear " + if (prm%aTolTwinfrac <= 0.0_pReal) extmsg = trim(extmsg)//"atoltwinfrac " + + if (extmsg /= '') call IO_error(211_pInt,ip=instance,& + ext_msg=trim(extmsg)//'('//PLASTICITY_PHENOPOWERLAW_label//')') - outputFound: if (mySize > 0_pInt) then - plastic_phenopowerlaw_sizePostResult(o,instance) = mySize - plastic_phenopowerlaw_sizePostResults(instance) = plastic_phenopowerlaw_sizePostResults(instance) + mySize - endif outputFound - enddo outputsLoop !-------------------------------------------------------------------------------------------------- ! allocate state arrays - sizeState = plastic_phenopowerlaw_totalNslip(instance) & ! s_slip - + plastic_phenopowerlaw_totalNtwin(instance) & ! s_twin - + 2_pInt & ! sum(gamma) + sum(f) - + plastic_phenopowerlaw_totalNslip(instance) & ! accshear_slip - + plastic_phenopowerlaw_totalNtwin(instance) ! accshear_twin + NipcMyPhase = count(material_phase == p) ! number of IPCs containing my phase + sizeState = size(['tau_slip ','accshear_slip']) * prm%TotalNslip & + + size(['tau_twin ','accshear_twin']) * prm%TotalNtwin & + + size(['sum(gamma)', 'sum(f) ']) - sizeDotState = sizeState - sizeDeltaState = 0_pInt - plasticState(phase)%sizeState = sizeState - plasticState(phase)%sizeDotState = sizeDotState - plasticState(phase)%sizeDeltaState = sizeDeltaState - plasticState(phase)%sizePostResults = plastic_phenopowerlaw_sizePostResults(instance) - plasticState(phase)%nSlip =plastic_phenopowerlaw_totalNslip(instance) - plasticState(phase)%nTwin =plastic_phenopowerlaw_totalNtwin(instance) - plasticState(phase)%nTrans=plastic_phenopowerlaw_totalNtrans(instance) - allocate(plasticState(phase)%aTolState ( sizeState), source=0.0_pReal) - allocate(plasticState(phase)%state0 ( sizeState,NipcMyPhase), source=0.0_pReal) - allocate(plasticState(phase)%partionedState0 ( sizeState,NipcMyPhase), source=0.0_pReal) - allocate(plasticState(phase)%subState0 ( sizeState,NipcMyPhase), source=0.0_pReal) - allocate(plasticState(phase)%state ( sizeState,NipcMyPhase), source=0.0_pReal) - allocate(plasticState(phase)%dotState (sizeDotState,NipcMyPhase), source=0.0_pReal) - allocate(plasticState(phase)%deltaState (sizeDeltaState,NipcMyPhase), source=0.0_pReal) - if (any(numerics_integrator == 1_pInt)) then - allocate(plasticState(phase)%previousDotState (sizeDotState,NipcMyPhase),source=0.0_pReal) - allocate(plasticState(phase)%previousDotState2(sizeDotState,NipcMyPhase),source=0.0_pReal) - endif - if (any(numerics_integrator == 4_pInt)) & - allocate(plasticState(phase)%RK4dotState (sizeDotState,NipcMyPhase), source=0.0_pReal) - if (any(numerics_integrator == 5_pInt)) & - allocate(plasticState(phase)%RKCK45dotState (6,sizeDotState,NipcMyPhase), source=0.0_pReal) + sizeDotState = sizeState + plasticState(p)%sizeState = sizeState + plasticState(p)%sizeDotState = sizeDotState + plasticState(p)%sizePostResults = sum(plastic_phenopowerlaw_sizePostResult(:,instance)) + plasticState(p)%nSlip = sum(prm%Nslip) + plasticState(p)%nTwin = sum(prm%Ntwin) + allocate(plasticState(p)%aTolState ( sizeState), source=0.0_pReal) + allocate(plasticState(p)%state0 ( sizeState,NipcMyPhase), source=0.0_pReal) + allocate(plasticState(p)%partionedState0 ( sizeState,NipcMyPhase), source=0.0_pReal) + allocate(plasticState(p)%subState0 ( sizeState,NipcMyPhase), source=0.0_pReal) + allocate(plasticState(p)%state ( sizeState,NipcMyPhase), source=0.0_pReal) - offset_slip = plasticState(phase)%nSlip+plasticState(phase)%nTwin+2_pInt - plasticState(phase)%slipRate => & - plasticState(phase)%dotState(offset_slip+1:offset_slip+plasticState(phase)%nSlip,1:NipcMyPhase) - plasticState(phase)%accumulatedSlip => & - plasticState(phase)%state(offset_slip+1:offset_slip+plasticState(phase)%nSlip,1:NipcMyPhase) - - do f = 1_pInt,lattice_maxNslipFamily ! >>> interaction slip -- X - index_myFamily = sum(plastic_phenopowerlaw_Nslip(1:f-1_pInt,instance)) - do j = 1_pInt,plastic_phenopowerlaw_Nslip(f,instance) ! loop over (active) systems in my family (slip) - do o = 1_pInt,lattice_maxNslipFamily - index_otherFamily = sum(plastic_phenopowerlaw_Nslip(1:o-1_pInt,instance)) - do k = 1_pInt,plastic_phenopowerlaw_Nslip(o,instance) ! loop over (active) systems in other family (slip) - plastic_phenopowerlaw_hardeningMatrix_SlipSlip(index_myFamily+j,index_otherFamily+k,instance) = & - plastic_phenopowerlaw_interaction_SlipSlip(lattice_interactionSlipSlip( & - sum(lattice_NslipSystem(1:f-1,phase))+j, & - sum(lattice_NslipSystem(1:o-1,phase))+k, & - phase), instance ) - enddo; enddo - - do o = 1_pInt,lattice_maxNtwinFamily - index_otherFamily = sum(plastic_phenopowerlaw_Ntwin(1:o-1_pInt,instance)) - do k = 1_pInt,plastic_phenopowerlaw_Ntwin(o,instance) ! loop over (active) systems in other family (twin) - plastic_phenopowerlaw_hardeningMatrix_SlipTwin(index_myFamily+j,index_otherFamily+k,instance) = & - plastic_phenopowerlaw_interaction_SlipTwin(lattice_interactionSlipTwin( & - sum(lattice_NslipSystem(1:f-1_pInt,phase))+j, & - sum(lattice_NtwinSystem(1:o-1_pInt,phase))+k, & - phase), instance ) - enddo; enddo - - enddo; enddo - - do f = 1_pInt,lattice_maxNtwinFamily ! >>> interaction twin -- X - index_myFamily = sum(plastic_phenopowerlaw_Ntwin(1:f-1_pInt,instance)) - do j = 1_pInt,plastic_phenopowerlaw_Ntwin(f,instance) ! loop over (active) systems in my family (twin) - - do o = 1_pInt,lattice_maxNslipFamily - index_otherFamily = sum(plastic_phenopowerlaw_Nslip(1:o-1_pInt,instance)) - do k = 1_pInt,plastic_phenopowerlaw_Nslip(o,instance) ! loop over (active) systems in other family (slip) - plastic_phenopowerlaw_hardeningMatrix_TwinSlip(index_myFamily+j,index_otherFamily+k,instance) = & - plastic_phenopowerlaw_interaction_TwinSlip(lattice_interactionTwinSlip( & - sum(lattice_NtwinSystem(1:f-1_pInt,phase))+j, & - sum(lattice_NslipSystem(1:o-1_pInt,phase))+k, & - phase), instance ) - enddo; enddo - - do o = 1_pInt,lattice_maxNtwinFamily - index_otherFamily = sum(plastic_phenopowerlaw_Ntwin(1:o-1_pInt,instance)) - do k = 1_pInt,plastic_phenopowerlaw_Ntwin(o,instance) ! loop over (active) systems in other family (twin) - plastic_phenopowerlaw_hardeningMatrix_TwinTwin(index_myFamily+j,index_otherFamily+k,instance) = & - plastic_phenopowerlaw_interaction_TwinTwin(lattice_interactionTwinTwin( & - sum(lattice_NtwinSystem(1:f-1_pInt,phase))+j, & - sum(lattice_NtwinSystem(1:o-1_pInt,phase))+k, & - phase), instance ) - enddo; enddo - - enddo; enddo - startIndex = 1_pInt - endIndex = plastic_phenopowerlaw_totalNslip(instance) - state (instance)%s_slip=>plasticState(phase)%state (startIndex:endIndex,:) - state0 (instance)%s_slip=>plasticState(phase)%state0 (startIndex:endIndex,:) - dotState(instance)%s_slip=>plasticState(phase)%dotState(startIndex:endIndex,:) - - startIndex = endIndex + 1_pInt - endIndex = endIndex + plastic_phenopowerlaw_totalNtwin(instance) - state (instance)%s_twin=>plasticState(phase)%state (startIndex:endIndex,:) - state0 (instance)%s_twin=>plasticState(phase)%state0 (startIndex:endIndex,:) - dotState(instance)%s_twin=>plasticState(phase)%dotState(startIndex:endIndex,:) - - startIndex = endIndex + 1_pInt - endIndex = endIndex + 1_pInt - state (instance)%sumGamma=>plasticState(phase)%state (startIndex,:) - state0 (instance)%sumGamma=>plasticState(phase)%state0 (startIndex,:) - dotState(instance)%sumGamma=>plasticState(phase)%dotState(startIndex,:) - - startIndex = endIndex + 1_pInt - endIndex = endIndex + 1_pInt - state (instance)%sumF=>plasticState(phase)%state (startIndex,:) - state0 (instance)%sumF=>plasticState(phase)%state0 (startIndex,:) - dotState(instance)%sumF=>plasticState(phase)%dotState(startIndex,:) - - startIndex = endIndex + 1_pInt - endIndex = endIndex +plastic_phenopowerlaw_totalNslip(instance) - state (instance)%accshear_slip=>plasticState(phase)%state (startIndex:endIndex,:) - state0 (instance)%accshear_slip=>plasticState(phase)%state0 (startIndex:endIndex,:) - dotState(instance)%accshear_slip=>plasticState(phase)%dotState(startIndex:endIndex,:) - - startIndex = endIndex + 1_pInt - endIndex = endIndex +plastic_phenopowerlaw_totalNtwin(instance) - state (instance)%accshear_twin=>plasticState(phase)%state (startIndex:endIndex,:) - state0 (instance)%accshear_twin=>plasticState(phase)%state0 (startIndex:endIndex,:) - dotState(instance)%accshear_twin=>plasticState(phase)%dotState(startIndex:endIndex,:) + allocate(plasticState(p)%dotState (sizeDotState,NipcMyPhase), source=0.0_pReal) + allocate(plasticState(p)%deltaState (0_pInt,NipcMyPhase), source=0.0_pReal) + if (any(numerics_integrator == 1_pInt)) then + allocate(plasticState(p)%previousDotState (sizeDotState,NipcMyPhase),source=0.0_pReal) + allocate(plasticState(p)%previousDotState2(sizeDotState,NipcMyPhase),source=0.0_pReal) + endif + if (any(numerics_integrator == 4_pInt)) & + allocate(plasticState(p)%RK4dotState (sizeDotState,NipcMyPhase), source=0.0_pReal) + if (any(numerics_integrator == 5_pInt)) & + allocate(plasticState(p)%RKCK45dotState (6,sizeDotState,NipcMyPhase), source=0.0_pReal) - call plastic_phenopowerlaw_stateInit(phase,instance) - call plastic_phenopowerlaw_aTolState(phase,instance) - endif myPhase2 - enddo initializeInstances +!-------------------------------------------------------------------------------------------------- +! calculate hardening matrices + allocate(temp1(sum(prm%Nslip),sum(prm%Nslip)),source =0.0_pReal) + allocate(temp2(sum(prm%Nslip),sum(prm%Ntwin)),source =0.0_pReal) + mySlipFamilies: do f = 1_pInt,size(prm%Nslip,1) ! >>> interaction slip -- X + index_myFamily = sum(prm%Nslip(1:f-1_pInt)) + + mySlipSystems: do j = 1_pInt,prm%Nslip(f) + otherSlipFamilies: do o = 1_pInt,size(prm%Nslip,1) + index_otherFamily = sum(prm%Nslip(1:o-1_pInt)) + otherSlipSystems: do k = 1_pInt,prm%Nslip(o) + temp1(index_myFamily+j,index_otherFamily+k) = & + prm%interaction_SlipSlip(lattice_interactionSlipSlip( & + sum(lattice_NslipSystem(1:f-1,p))+j, & + sum(lattice_NslipSystem(1:o-1,p))+k, & + p),1) + enddo otherSlipSystems; enddo otherSlipFamilies + + twinFamilies: do o = 1_pInt,size(prm%Ntwin,1) + index_otherFamily = sum(prm%Ntwin(1:o-1_pInt)) + twinSystems: do k = 1_pInt,prm%Ntwin(o) + temp2(index_myFamily+j,index_otherFamily+k) = & + prm%interaction_SlipTwin(lattice_interactionSlipTwin( & + sum(lattice_NslipSystem(1:f-1_pInt,p))+j, & + sum(lattice_NtwinSystem(1:o-1_pInt,p))+k, & + p),1) + enddo twinSystems; enddo twinFamilies + enddo mySlipSystems + enddo mySlipFamilies + prm%interaction_SlipSlip = temp1; deallocate(temp1) + prm%interaction_SlipTwin = temp2; deallocate(temp2) + + + allocate(temp1(sum(prm%Ntwin),sum(prm%Nslip)),source =0.0_pReal) + allocate(temp2(sum(prm%Ntwin),sum(prm%Ntwin)),source =0.0_pReal) + myTwinFamilies: do f = 1_pInt,size(prm%Ntwin,1) ! >>> interaction twin -- X + index_myFamily = sum(prm%Ntwin(1:f-1_pInt)) + myTwinSystems: do j = 1_pInt,prm%Ntwin(f) + slipFamilies: do o = 1_pInt,size(prm%Nslip,1) + index_otherFamily = sum(prm%Nslip(1:o-1_pInt)) + slipSystems: do k = 1_pInt,prm%Nslip(o) + temp1(index_myFamily+j,index_otherFamily+k) = & + prm%interaction_TwinSlip(lattice_interactionTwinSlip( & + sum(lattice_NtwinSystem(1:f-1_pInt,p))+j, & + sum(lattice_NslipSystem(1:o-1_pInt,p))+k, & + p),1) + enddo slipSystems; enddo slipFamilies + + otherTwinFamilies: do o = 1_pInt,size(prm%Ntwin,1) + index_otherFamily = sum(prm%Ntwin(1:o-1_pInt)) + otherTwinSystems: do k = 1_pInt,prm%Ntwin(o) + temp2(index_myFamily+j,index_otherFamily+k) = & + prm%interaction_TwinTwin(lattice_interactionTwinTwin( & + sum(lattice_NtwinSystem(1:f-1_pInt,p))+j, & + sum(lattice_NtwinSystem(1:o-1_pInt,p))+k, & + p),1) + enddo otherTwinSystems; enddo otherTwinFamilies + enddo myTwinSystems + enddo myTwinFamilies + prm%interaction_TwinSlip = temp1; deallocate(temp1) + prm%interaction_TwinTwin = temp2; deallocate(temp2) + +!-------------------------------------------------------------------------------------------------- +! locally defined state aliases and initialization of state0 and aTolState + startIndex = 1_pInt + endIndex = plasticState(p)%nSlip + state (instance)%s_slip=>plasticState(p)%state (startIndex:endIndex,:) + dotState(instance)%s_slip=>plasticState(p)%dotState(startIndex:endIndex,:) + plasticState(p)%state0(startIndex:endIndex,:) = & + spread(math_expand(prm%tau0_slip, prm%Nslip), 2, NipcMyPhase) + plasticState(p)%aTolState(startIndex:endIndex) = prm%aTolResistance + + startIndex = endIndex + 1_pInt + endIndex = endIndex + plasticState(p)%nTwin + state (instance)%s_twin=>plasticState(p)%state (startIndex:endIndex,:) + dotState(instance)%s_twin=>plasticState(p)%dotState(startIndex:endIndex,:) + plasticState(p)%state0(startIndex:endIndex,:) = & + spread(math_expand(prm%tau0_twin, prm%Ntwin), 2, NipcMyPhase) + plasticState(p)%aTolState(startIndex:endIndex) = prm%aTolResistance + + startIndex = endIndex + 1_pInt + endIndex = endIndex + 1_pInt + state (instance)%sumGamma=>plasticState(p)%state (startIndex,:) + dotState(instance)%sumGamma=>plasticState(p)%dotState(startIndex,:) + plasticState(p)%aTolState(startIndex:endIndex) = prm%aTolShear + + startIndex = endIndex + 1_pInt + endIndex = endIndex + 1_pInt + state (instance)%sumF=>plasticState(p)%state (startIndex,:) + dotState(instance)%sumF=>plasticState(p)%dotState(startIndex,:) + plasticState(p)%aTolState(startIndex:endIndex) = prm%aTolTwinFrac + + startIndex = endIndex + 1_pInt + endIndex = endIndex + plasticState(p)%nSlip + state (instance)%accshear_slip=>plasticState(p)%state (startIndex:endIndex,:) + dotState(instance)%accshear_slip=>plasticState(p)%dotState(startIndex:endIndex,:) + plasticState(p)%aTolState(startIndex:endIndex) = prm%aTolShear + ! global alias + plasticState(p)%slipRate =>plasticState(p)%dotState(startIndex:endIndex,:) + plasticState(p)%accumulatedSlip =>plasticState(p)%state(startIndex:endIndex,:) + + startIndex = endIndex + 1_pInt + endIndex = endIndex + plasticState(p)%nTwin + state (instance)%accshear_twin=>plasticState(p)%state (startIndex:endIndex,:) + dotState(instance)%accshear_twin=>plasticState(p)%dotState(startIndex:endIndex,:) + plasticState(p)%aTolState(startIndex:endIndex) = prm%aTolShear + + dotState(instance)%whole =>plasticState(p)%dotState + + end associate + enddo end subroutine plastic_phenopowerlaw_init -!-------------------------------------------------------------------------------------------------- -!> @brief sets the initial microstructural state for a given instance of this plasticity -!-------------------------------------------------------------------------------------------------- -subroutine plastic_phenopowerlaw_stateInit(ph,instance) - use lattice, only: & - lattice_maxNslipFamily, & - lattice_maxNtwinFamily - use material, only: & - plasticState - - implicit none - integer(pInt), intent(in) :: & - instance, & !< number specifying the instance of the plasticity - ph - integer(pInt) :: & - i - real(pReal), dimension(plasticState(ph)%sizeState) :: & - tempState - - tempState = 0.0_pReal - do i = 1_pInt,lattice_maxNslipFamily - tempState(1+sum(plastic_phenopowerlaw_Nslip(1:i-1,instance)) : & - sum(plastic_phenopowerlaw_Nslip(1:i ,instance))) = & - plastic_phenopowerlaw_tau0_slip(i,instance) - enddo - - do i = 1_pInt,lattice_maxNtwinFamily - tempState(1+sum(plastic_phenopowerlaw_Nslip(:,instance))+& - sum(plastic_phenopowerlaw_Ntwin(1:i-1,instance)) : & - sum(plastic_phenopowerlaw_Nslip(:,instance))+& - sum(plastic_phenopowerlaw_Ntwin(1:i ,instance))) = & - plastic_phenopowerlaw_tau0_twin(i,instance) - enddo - - plasticState(ph)%state0(:,:) = spread(tempState, & ! spread single tempstate array - 2, & ! along dimension 2 - size(plasticState(ph)%state0(1,:))) ! number of copies (number of IPCs) - -end subroutine plastic_phenopowerlaw_stateInit - - -!-------------------------------------------------------------------------------------------------- -!> @brief sets the relevant state values for a given instance of this plasticity -!-------------------------------------------------------------------------------------------------- -subroutine plastic_phenopowerlaw_aTolState(ph,instance) - use material, only: & - plasticState - - implicit none - integer(pInt), intent(in) :: & - instance, & !< number specifying the instance of the plasticity - ph - - plasticState(ph)%aTolState(1:plastic_phenopowerlaw_totalNslip(instance)+ & - plastic_phenopowerlaw_totalNtwin(instance)) = & - plastic_phenopowerlaw_aTolResistance(instance) - plasticState(ph)%aTolState(1+plastic_phenopowerlaw_totalNslip(instance)+ & - plastic_phenopowerlaw_totalNtwin(instance)) = & - plastic_phenopowerlaw_aTolShear(instance) - plasticState(ph)%aTolState(2+plastic_phenopowerlaw_totalNslip(instance)+ & - plastic_phenopowerlaw_totalNtwin(instance)) = & - plastic_phenopowerlaw_aTolTwinFrac(instance) - plasticState(ph)%aTolState(3+plastic_phenopowerlaw_totalNslip(instance)+ & - plastic_phenopowerlaw_totalNtwin(instance): & - 2+2*(plastic_phenopowerlaw_totalNslip(instance)+ & - plastic_phenopowerlaw_totalNtwin(instance))) = & - plastic_phenopowerlaw_aTolShear(instance) - -end subroutine plastic_phenopowerlaw_aTolState - - !-------------------------------------------------------------------------------------------------- !> @brief calculates plastic velocity gradient and its tangent !-------------------------------------------------------------------------------------------------- @@ -782,13 +478,11 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip, lattice_Sslip_v, & lattice_Stwin, & lattice_Stwin_v, & - lattice_maxNslipFamily, & - lattice_maxNtwinFamily, & lattice_NslipSystem, & - lattice_NtwinSystem, & - lattice_NnonSchmid + lattice_NtwinSystem use material, only: & - phaseAt, phasememberAt, & + phasememberAt, & + material_phase, & phase_plasticityInstance implicit none @@ -805,7 +499,6 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip, Tstar_v !< 2nd Piola Kirchhoff stress tensor in Mandel notation integer(pInt) :: & - instance, & index_myFamily, & f,i,j,k,l,m,n, & of, & @@ -819,10 +512,13 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip, dLp_dTstar3333 !< derivative of Lp with respect to Tstar as 4th order tensor real(pReal), dimension(3,3,2) :: & nonSchmid_tensor + type(tParameters) :: prm + type(tPhenopowerlawState) :: stt of = phasememberAt(ipc,ip,el) - ph = phaseAt(ipc,ip,el) - instance = phase_plasticityInstance(ph) + ph = material_phase(ipc,ip,el) + + associate(prm => param(phase_plasticityInstance(ph)), stt => state(phase_plasticityInstance(ph))) Lp = 0.0_pReal dLp_dTstar3333 = 0.0_pReal @@ -831,9 +527,9 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip, !-------------------------------------------------------------------------------------------------- ! Slip part j = 0_pInt - slipFamilies: do f = 1_pInt,lattice_maxNslipFamily + slipFamilies: do f = 1_pInt,size(prm%Nslip,1) index_myFamily = sum(lattice_NslipSystem(1:f-1_pInt,ph)) ! at which index starts my family - slipSystems: do i = 1_pInt,plastic_phenopowerlaw_Nslip(f,instance) + slipSystems: do i = 1_pInt,prm%Nslip(f) j = j+1_pInt ! Calculation of Lp @@ -841,38 +537,36 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip, tau_slip_neg = tau_slip_pos nonSchmid_tensor(1:3,1:3,1) = lattice_Sslip(1:3,1:3,1,index_myFamily+i,ph) nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,1) - do k = 1,lattice_NnonSchmid(ph) - tau_slip_pos = tau_slip_pos + plastic_phenopowerlaw_nonSchmidCoeff(k,instance)* & + do k = 1,size(prm%nonSchmidCoeff) + tau_slip_pos = tau_slip_pos + prm%nonSchmidCoeff(k)* & dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k,index_myFamily+i,ph)) - tau_slip_neg = tau_slip_neg + plastic_phenopowerlaw_nonSchmidCoeff(k,instance)* & + tau_slip_neg = tau_slip_neg + prm%nonSchmidCoeff(k)* & dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k+1,index_myFamily+i,ph)) - nonSchmid_tensor(1:3,1:3,1) = nonSchmid_tensor(1:3,1:3,1) + plastic_phenopowerlaw_nonSchmidCoeff(k,instance)*& + nonSchmid_tensor(1:3,1:3,1) = nonSchmid_tensor(1:3,1:3,1) + prm%nonSchmidCoeff(k)*& lattice_Sslip(1:3,1:3,2*k,index_myFamily+i,ph) - nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,2) + plastic_phenopowerlaw_nonSchmidCoeff(k,instance)*& + nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,2) + prm%nonSchmidCoeff(k)*& lattice_Sslip(1:3,1:3,2*k+1,index_myFamily+i,ph) enddo - gdot_slip_pos = 0.5_pReal*plastic_phenopowerlaw_gdot0_slip(instance)* & - ((abs(tau_slip_pos)/(state(instance)%s_slip(j,of))) & - **plastic_phenopowerlaw_n_slip(instance))*sign(1.0_pReal,tau_slip_pos) + gdot_slip_pos = 0.5_pReal*prm%gdot0_slip* & + ((abs(tau_slip_pos)/(stt%s_slip(j,of)))**prm%n_slip)*sign(1.0_pReal,tau_slip_pos) - gdot_slip_neg = 0.5_pReal*plastic_phenopowerlaw_gdot0_slip(instance)* & - ((abs(tau_slip_neg)/(state(instance)%s_slip(j,of))) & - **plastic_phenopowerlaw_n_slip(instance))*sign(1.0_pReal,tau_slip_neg) + gdot_slip_neg = 0.5_pReal*prm%gdot0_slip* & + ((abs(tau_slip_neg)/(stt%s_slip(j,of)))**prm%n_slip)*sign(1.0_pReal,tau_slip_neg) - Lp = Lp + (1.0_pReal-state(instance)%sumF(of))*& ! 1-F + Lp = Lp + (1.0_pReal-stt%sumF(of))*& (gdot_slip_pos+gdot_slip_neg)*lattice_Sslip(1:3,1:3,1,index_myFamily+i,ph) ! Calculation of the tangent of Lp - if (dNeq0(gdot_slip_pos)) then - dgdot_dtauslip_pos = gdot_slip_pos*plastic_phenopowerlaw_n_slip(instance)/tau_slip_pos + if (dNeq0(tau_slip_pos)) then + dgdot_dtauslip_pos = gdot_slip_pos*prm%n_slip/tau_slip_pos forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) & dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n) + & dgdot_dtauslip_pos*lattice_Sslip(k,l,1,index_myFamily+i,ph)* & nonSchmid_tensor(m,n,1) endif - if (dNeq0(gdot_slip_neg)) then - dgdot_dtauslip_neg = gdot_slip_neg*plastic_phenopowerlaw_n_slip(instance)/tau_slip_neg + if (dNeq0(tau_slip_neg)) then + dgdot_dtauslip_neg = gdot_slip_neg*prm%n_slip/tau_slip_neg forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) & dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n) + & dgdot_dtauslip_neg*lattice_Sslip(k,l,1,index_myFamily+i,ph)* & @@ -884,22 +578,21 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip, !-------------------------------------------------------------------------------------------------- ! Twinning part j = 0_pInt - twinFamilies: do f = 1_pInt,lattice_maxNtwinFamily + twinFamilies: do f = 1_pInt,size(prm%Ntwin,1) index_myFamily = sum(lattice_NtwinSystem(1:f-1_pInt,ph)) ! at which index starts my family - twinSystems: do i = 1_pInt,plastic_phenopowerlaw_Ntwin(f,instance) + twinSystems: do i = 1_pInt,prm%Ntwin(f) j = j+1_pInt ! Calculation of Lp tau_twin = dot_product(Tstar_v,lattice_Stwin_v(1:6,index_myFamily+i,ph)) - gdot_twin = (1.0_pReal-state(instance)%sumF(of))*& ! 1-F - plastic_phenopowerlaw_gdot0_twin(instance)*& - (abs(tau_twin)/state(instance)%s_twin(j,of))**& - plastic_phenopowerlaw_n_twin(instance)*max(0.0_pReal,sign(1.0_pReal,tau_twin)) + gdot_twin = (1.0_pReal-stt%sumF(of))*prm%gdot0_twin*& + (abs(tau_twin)/stt%s_twin(j,of))**& + prm%n_twin*max(0.0_pReal,sign(1.0_pReal,tau_twin)) Lp = Lp + gdot_twin*lattice_Stwin(1:3,1:3,index_myFamily+i,ph) ! Calculation of the tangent of Lp if (dNeq0(gdot_twin)) then - dgdot_dtautwin = gdot_twin*plastic_phenopowerlaw_n_twin(instance)/tau_twin + dgdot_dtautwin = gdot_twin*prm%n_twin/tau_twin forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) & dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n) + & dgdot_dtautwin*lattice_Stwin(k,l,index_myFamily+i,ph)* & @@ -910,9 +603,10 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip, dLp_dTstar99 = math_Plain3333to99(dLp_dTstar3333) - + end associate end subroutine plastic_phenopowerlaw_LpAndItsTangent + !-------------------------------------------------------------------------------------------------- !> @brief calculates the rate of change of microstructure !-------------------------------------------------------------------------------------------------- @@ -920,16 +614,12 @@ subroutine plastic_phenopowerlaw_dotState(Tstar_v,ipc,ip,el) use lattice, only: & lattice_Sslip_v, & lattice_Stwin_v, & - lattice_maxNslipFamily, & - lattice_maxNtwinFamily, & lattice_NslipSystem, & lattice_NtwinSystem, & - lattice_shearTwin, & - lattice_NnonSchmid + lattice_shearTwin use material, only: & material_phase, & - phaseAt, phasememberAt, & - plasticState, & + phasememberAt, & phase_plasticityInstance implicit none @@ -941,142 +631,108 @@ subroutine plastic_phenopowerlaw_dotState(Tstar_v,ipc,ip,el) el !< element !< microstructure state integer(pInt) :: & - instance,ph, & - nSlip,nTwin, & + ph, & f,i,j,k, & - index_Gamma,index_F,index_myFamily, & - offset_accshear_slip,offset_accshear_twin, & + index_myFamily, & of real(pReal) :: & c_SlipSlip,c_TwinSlip,c_TwinTwin, & ssat_offset, & tau_slip_pos,tau_slip_neg,tau_twin - real(pReal), dimension(plastic_phenopowerlaw_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))) :: & - gdot_slip,left_SlipSlip,left_SlipTwin,right_SlipSlip,right_TwinSlip - real(pReal), dimension(plastic_phenopowerlaw_totalNtwin(phase_plasticityInstance(material_phase(ipc,ip,el)))) :: & - gdot_twin,left_TwinSlip,left_TwinTwin,right_SlipTwin,right_TwinTwin + real(pReal), dimension(param(phase_plasticityInstance(material_phase(ipc,ip,el)))%totalNslip) :: & + gdot_slip,left_SlipSlip,right_SlipSlip + real(pReal), dimension(param(phase_plasticityInstance(material_phase(ipc,ip,el)))%totalNtwin) :: & + gdot_twin + + type(tParameters) :: prm + type(tPhenopowerlawState) :: dst,stt of = phasememberAt(ipc,ip,el) - ph = phaseAt(ipc,ip,el) - instance = phase_plasticityInstance(ph) + ph = material_phase(ipc,ip,el) + associate( prm => param(phase_plasticityInstance(ph)), & + stt => state(phase_plasticityInstance(ph)), & + dst => dotState(phase_plasticityInstance(ph))) - nSlip = plastic_phenopowerlaw_totalNslip(instance) - nTwin = plastic_phenopowerlaw_totalNtwin(instance) - - index_Gamma = nSlip + nTwin + 1_pInt - index_F = nSlip + nTwin + 2_pInt - offset_accshear_slip = nSlip + nTwin + 2_pInt - offset_accshear_twin = nSlip + nTwin + 2_pInt + nSlip - plasticState(ph)%dotState(:,of) = 0.0_pReal + dst%whole(:,of) = 0.0_pReal !-------------------------------------------------------------------------------------------------- ! system-independent (nonlinear) prefactors to M_Xx (X influenced by x) matrices - c_SlipSlip = plastic_phenopowerlaw_h0_SlipSlip(instance)*& - (1.0_pReal + plastic_phenopowerlaw_twinC(instance)*plasticState(ph)%state(index_F,of)**& - plastic_phenopowerlaw_twinB(instance)) - c_TwinSlip = plastic_phenopowerlaw_h0_TwinSlip(instance)*& - plasticState(ph)%state(index_Gamma,of)**plastic_phenopowerlaw_twinE(instance) - c_TwinTwin = plastic_phenopowerlaw_h0_TwinTwin(instance)*& - plasticState(ph)%state(index_F,of)**plastic_phenopowerlaw_twinD(instance) + c_SlipSlip = prm%h0_slipslip * (1.0_pReal + prm%twinC*stt%sumF(of)** prm%twinB) + c_TwinSlip = prm%h0_TwinSlip * stt%sumGamma(of)**prm%twinE + c_TwinTwin = prm%h0_TwinTwin * stt%sumF(of)**prm%twinD !-------------------------------------------------------------------------------------------------- ! calculate left and right vectors and calculate dot gammas - ssat_offset = plastic_phenopowerlaw_spr(instance)*sqrt(plasticState(ph)%state(index_F,of)) + ssat_offset = prm%spr*sqrt(stt%sumF(of)) j = 0_pInt - slipFamilies1: do f = 1_pInt,lattice_maxNslipFamily + slipFamilies1: do f =1_pInt,size(prm%Nslip,1) index_myFamily = sum(lattice_NslipSystem(1:f-1_pInt,ph)) ! at which index starts my family - slipSystems1: do i = 1_pInt,plastic_phenopowerlaw_Nslip(f,instance) + slipSystems1: do i = 1_pInt,prm%Nslip(f) j = j+1_pInt - left_SlipSlip(j) = 1.0_pReal + plastic_phenopowerlaw_H_int(f,instance) ! modified no system-dependent left part - left_SlipTwin(j) = 1.0_pReal ! no system-dependent left part - right_SlipSlip(j) = abs(1.0_pReal-plasticState(ph)%state(j,of) / & - (plastic_phenopowerlaw_tausat_slip(f,instance)+ssat_offset)) & - **plastic_phenopowerlaw_a_slip(instance)& - *sign(1.0_pReal,1.0_pReal-plasticState(ph)%state(j,of) / & - (plastic_phenopowerlaw_tausat_slip(f,instance)+ssat_offset)) - right_TwinSlip(j) = 1.0_pReal ! no system-dependent part + left_SlipSlip(j) = 1.0_pReal + prm%H_int(f) ! modified no system-dependent left part + right_SlipSlip(j) = abs(1.0_pReal-stt%s_slip(j,of) / (prm%tausat_slip(f)+ssat_offset)) **prm%a_slip & + * sign(1.0_pReal,1.0_pReal-stt%s_slip(j,of) / (prm%tausat_slip(f)+ssat_offset)) !-------------------------------------------------------------------------------------------------- ! Calculation of dot gamma tau_slip_pos = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,ph)) tau_slip_neg = tau_slip_pos - nonSchmidSystems: do k = 1,lattice_NnonSchmid(ph) - tau_slip_pos = tau_slip_pos + plastic_phenopowerlaw_nonSchmidCoeff(k,instance)* & + nonSchmidSystems: do k = 1,size(prm%nonSchmidCoeff) + tau_slip_pos = tau_slip_pos + prm%nonSchmidCoeff(k)* & dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k, index_myFamily+i,ph)) - tau_slip_neg = tau_slip_neg + plastic_phenopowerlaw_nonSchmidCoeff(k,instance)* & + tau_slip_neg = tau_slip_neg +prm%nonSchmidCoeff(k)* & dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k+1,index_myFamily+i,ph)) enddo nonSchmidSystems - gdot_slip(j) = plastic_phenopowerlaw_gdot0_slip(instance)*0.5_pReal* & - ((abs(tau_slip_pos)/(plasticState(ph)%state(j,of)))**plastic_phenopowerlaw_n_slip(instance) & - *sign(1.0_pReal,tau_slip_pos) & - +(abs(tau_slip_neg)/(plasticState(ph)%state(j,of)))**plastic_phenopowerlaw_n_slip(instance) & - *sign(1.0_pReal,tau_slip_neg)) + gdot_slip(j) = prm%gdot0_slip*0.5_pReal* & + ( (abs(tau_slip_pos)/(stt%s_slip(j,of)))**prm%n_slip*sign(1.0_pReal,tau_slip_pos) & + +(abs(tau_slip_neg)/(stt%s_slip(j,of)))**prm%n_slip*sign(1.0_pReal,tau_slip_neg)) enddo slipSystems1 enddo slipFamilies1 - - j = 0_pInt - twinFamilies1: do f = 1_pInt,lattice_maxNtwinFamily + twinFamilies1: do f = 1_pInt,size(prm%Ntwin,1) index_myFamily = sum(lattice_NtwinSystem(1:f-1_pInt,ph)) ! at which index starts my family - twinSystems1: do i = 1_pInt,plastic_phenopowerlaw_Ntwin(f,instance) + twinSystems1: do i = 1_pInt,prm%Ntwin(f) j = j+1_pInt - left_TwinSlip(j) = 1.0_pReal ! no system-dependent left part - left_TwinTwin(j) = 1.0_pReal ! no system-dependent left part - right_SlipTwin(j) = 1.0_pReal ! no system-dependent right part - right_TwinTwin(j) = 1.0_pReal ! no system-dependent right part !-------------------------------------------------------------------------------------------------- ! Calculation of dot vol frac tau_twin = dot_product(Tstar_v,lattice_Stwin_v(1:6,index_myFamily+i,ph)) - gdot_twin(j) = (1.0_pReal-plasticState(ph)%state(index_F,of))*& ! 1-F - plastic_phenopowerlaw_gdot0_twin(instance)*& - (abs(tau_twin)/plasticState(ph)%state(nslip+j,of))**& - plastic_phenopowerlaw_n_twin(instance)*max(0.0_pReal,sign(1.0_pReal,tau_twin)) + gdot_twin(j) = (1.0_pReal-stt%sumF(of))*& ! 1-F + prm%gdot0_twin*& + (abs(tau_twin)/stt%s_twin(j,of))**& + prm%n_twin*max(0.0_pReal,sign(1.0_pReal,tau_twin)) enddo twinSystems1 enddo twinFamilies1 !-------------------------------------------------------------------------------------------------- ! calculate the overall hardening based on above - j = 0_pInt - slipFamilies2: do f = 1_pInt,lattice_maxNslipFamily - slipSystems2: do i = 1_pInt,plastic_phenopowerlaw_Nslip(f,instance) - j = j+1_pInt - plasticState(ph)%dotState(j,of) = & ! evolution of slip resistance j - c_SlipSlip * left_SlipSlip(j) * & - dot_product(plastic_phenopowerlaw_hardeningMatrix_SlipSlip(j,1:nSlip,instance), & - right_SlipSlip*abs(gdot_slip)) + & ! dot gamma_slip modulated by right-side slip factor - dot_product(plastic_phenopowerlaw_hardeningMatrix_SlipTwin(j,1:nTwin,instance), & - right_SlipTwin*gdot_twin) ! dot gamma_twin modulated by right-side twin factor - plasticState(ph)%dotState(index_Gamma,of) = plasticState(ph)%dotState(index_Gamma,of) + & - abs(gdot_slip(j)) - plasticState(ph)%dotState(offset_accshear_slip+j,of) = abs(gdot_slip(j)) - enddo slipSystems2 - enddo slipFamilies2 + do j = 1_pInt,prm%totalNslip + dst%s_slip(j,of) = c_SlipSlip * left_SlipSlip(j) * & ! evolution of slip resistance j + dot_product(prm%interaction_SlipSlip(j,1:prm%totalNslip),right_SlipSlip*abs(gdot_slip)) + & ! dot gamma_slip modulated by right-side slip factor + dot_product(prm%interaction_SlipTwin(j,1:prm%totalNtwin),gdot_twin) ! dot gamma_twin modulated by right-side twin factor + enddo + dst%sumGamma(of) = dst%sumGamma(of) + sum(abs(gdot_slip)) + dst%accshear_slip(1:prm%totalNslip,of) = abs(gdot_slip) j = 0_pInt - twinFamilies2: do f = 1_pInt,lattice_maxNtwinFamily + twinFamilies2: do f = 1_pInt,size(prm%Ntwin,1) index_myFamily = sum(lattice_NtwinSystem(1:f-1_pInt,ph)) ! at which index starts my family - twinSystems2: do i = 1_pInt,plastic_phenopowerlaw_Ntwin(f,instance) + twinSystems2: do i = 1_pInt,prm%Ntwin(f) j = j+1_pInt - plasticState(ph)%dotState(j+nSlip,of) = & ! evolution of twin resistance j - c_TwinSlip * left_TwinSlip(j) * & - dot_product(plastic_phenopowerlaw_hardeningMatrix_TwinSlip(j,1:nSlip,instance), & - right_TwinSlip*abs(gdot_slip)) + & ! dot gamma_slip modulated by right-side slip factor - c_TwinTwin * left_TwinTwin(j) * & - dot_product(plastic_phenopowerlaw_hardeningMatrix_TwinTwin(j,1:nTwin,instance), & - right_TwinTwin*gdot_twin) ! dot gamma_twin modulated by right-side twin factor - if (plasticState(ph)%state(index_F,of) < 0.98_pReal) & ! ensure twin volume fractions stays below 1.0 - plasticState(ph)%dotState(index_F,of) = plasticState(ph)%dotState(index_F,of) + & - gdot_twin(j)/lattice_shearTwin(index_myFamily+i,ph) - plasticState(ph)%dotState(offset_accshear_twin+j,of) = abs(gdot_twin(j)) + dst%s_twin(j,of) = & ! evolution of twin resistance j + c_TwinSlip * dot_product(prm%interaction_TwinSlip(j,1:prm%totalNslip),abs(gdot_slip)) + & ! dot gamma_slip modulated by right-side slip factor + c_TwinTwin * dot_product(prm%interaction_TwinTwin(j,1:prm%totalNtwin),gdot_twin) ! dot gamma_twin modulated by right-side twin factor + if (stt%sumF(of) < 0.98_pReal) & ! ensure twin volume fractions stays below 1.0 + dst%sumF(of) = dst%sumF(of) + gdot_twin(j)/lattice_shearTwin(index_myFamily+i,ph) + dst%accshear_twin(j,of) = abs(gdot_twin(j)) enddo twinSystems2 enddo twinFamilies2 - - + end associate end subroutine plastic_phenopowerlaw_dotState + !-------------------------------------------------------------------------------------------------- !> @brief return array of constitutive results !-------------------------------------------------------------------------------------------------- @@ -1084,13 +740,11 @@ function plastic_phenopowerlaw_postResults(Tstar_v,ipc,ip,el) use material, only: & material_phase, & plasticState, & - phaseAt, phasememberAt, & + phasememberAt, & phase_plasticityInstance use lattice, only: & lattice_Sslip_v, & lattice_Stwin_v, & - lattice_maxNslipFamily, & - lattice_maxNtwinFamily, & lattice_NslipSystem, & lattice_NtwinSystem, & lattice_NnonSchmid @@ -1103,126 +757,122 @@ function plastic_phenopowerlaw_postResults(Tstar_v,ipc,ip,el) ip, & !< integration point el !< element !< microstructure state - real(pReal), dimension(plastic_phenopowerlaw_sizePostResults(phase_plasticityInstance(material_phase(ipc,ip,el)))) :: & + real(pReal), dimension(plasticState(material_phase(ipc,ip,el))%sizePostResults) :: & plastic_phenopowerlaw_postResults integer(pInt) :: & - instance,ph, of, & - nSlip,nTwin, & + ph, of, & o,f,i,c,j,k, & - index_Gamma,index_F,index_accshear_slip,index_accshear_twin,index_myFamily + index_myFamily real(pReal) :: & tau_slip_pos,tau_slip_neg,tau + type(tParameters) :: prm + type(tPhenopowerlawState) :: stt, dst + of = phasememberAt(ipc,ip,el) - ph = phaseAt(ipc,ip,el) - instance = phase_plasticityInstance(ph) + ph = material_phase(ipc,ip,el) - nSlip = plastic_phenopowerlaw_totalNslip(instance) - nTwin = plastic_phenopowerlaw_totalNtwin(instance) - - index_Gamma = nSlip + nTwin + 1_pInt - index_F = nSlip + nTwin + 2_pInt - index_accshear_slip = nSlip + nTwin + 3_pInt - index_accshear_twin = nSlip + nTwin + 3_pInt + nSlip + associate( prm => param(phase_plasticityInstance(ph)), & + stt => state(phase_plasticityInstance(ph)), & + dst => dotState(phase_plasticityInstance(ph))) plastic_phenopowerlaw_postResults = 0.0_pReal c = 0_pInt - outputsLoop: do o = 1_pInt,plastic_phenopowerlaw_Noutput(instance) - select case(plastic_phenopowerlaw_outputID(o,instance)) + outputsLoop: do o = 1_pInt,size(prm%outputID) + select case(prm%outputID(o)) case (resistance_slip_ID) - plastic_phenopowerlaw_postResults(c+1_pInt:c+nSlip) = plasticState(ph)%state(1:nSlip,of) - c = c + nSlip + plastic_phenopowerlaw_postResults(c+1_pInt:c+prm%totalNslip) = stt%s_slip(1:prm%totalNslip,of) + c = c + prm%totalNslip case (accumulatedshear_slip_ID) - plastic_phenopowerlaw_postResults(c+1_pInt:c+nSlip) = plasticState(ph)%state(index_accshear_slip:& - index_accshear_slip+nSlip-1_pInt,of) - c = c + nSlip + plastic_phenopowerlaw_postResults(c+1_pInt:c+prm%totalNslip) = stt%accshear_slip(1:prm%totalNslip,of) + c = c + prm%totalNslip case (shearrate_slip_ID) j = 0_pInt - slipFamilies1: do f = 1_pInt,lattice_maxNslipFamily + slipFamilies1: do f = 1_pInt,size(prm%Nslip,1) index_myFamily = sum(lattice_NslipSystem(1:f-1_pInt,ph)) ! at which index starts my family - slipSystems1: do i = 1_pInt,plastic_phenopowerlaw_Nslip(f,instance) + slipSystems1: do i = 1_pInt,prm%Nslip(f) j = j + 1_pInt tau_slip_pos = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,ph)) tau_slip_neg = tau_slip_pos do k = 1,lattice_NnonSchmid(ph) - tau_slip_pos = tau_slip_pos + plastic_phenopowerlaw_nonSchmidCoeff(k,instance)* & + tau_slip_pos = tau_slip_pos +prm%nonSchmidCoeff(k)* & dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k,index_myFamily+i,ph)) - tau_slip_neg = tau_slip_neg + plastic_phenopowerlaw_nonSchmidCoeff(k,instance)* & + tau_slip_neg = tau_slip_neg +prm%nonSchmidCoeff(k)* & dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k+1,index_myFamily+i,ph)) enddo - plastic_phenopowerlaw_postResults(c+j) = plastic_phenopowerlaw_gdot0_slip(instance)*0.5_pReal* & - ((abs(tau_slip_pos)/plasticState(ph)%state(j,of))**plastic_phenopowerlaw_n_slip(instance) & + plastic_phenopowerlaw_postResults(c+j) = prm%gdot0_slip*0.5_pReal* & + ((abs(tau_slip_pos)/stt%s_slip(j,of))**prm%n_slip & *sign(1.0_pReal,tau_slip_pos) & - +(abs(tau_slip_neg)/(plasticState(ph)%state(j,of)))**plastic_phenopowerlaw_n_slip(instance) & + +(abs(tau_slip_neg)/(stt%s_slip(j,of)))**prm%n_slip & *sign(1.0_pReal,tau_slip_neg)) enddo slipSystems1 enddo slipFamilies1 - c = c + nSlip + c = c + prm%totalNslip case (resolvedstress_slip_ID) j = 0_pInt - slipFamilies2: do f = 1_pInt,lattice_maxNslipFamily + slipFamilies2: do f = 1_pInt,size(prm%Nslip,1) index_myFamily = sum(lattice_NslipSystem(1:f-1_pInt,ph)) ! at which index starts my family - slipSystems2: do i = 1_pInt,plastic_phenopowerlaw_Nslip(f,instance) + slipSystems2: do i = 1_pInt,prm%Nslip(f) j = j + 1_pInt plastic_phenopowerlaw_postResults(c+j) = & dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,ph)) enddo slipSystems2 enddo slipFamilies2 - c = c + nSlip + c = c + prm%totalNslip case (totalshear_ID) - plastic_phenopowerlaw_postResults(c+1_pInt) = & - plasticState(ph)%state(index_Gamma,of) + plastic_phenopowerlaw_postResults(c+1_pInt) = stt%sumGamma(of) c = c + 1_pInt case (resistance_twin_ID) - plastic_phenopowerlaw_postResults(c+1_pInt:c+nTwin) = & - plasticState(ph)%state(1_pInt+nSlip:1_pInt+nSlip+nTwin-1_pInt,of) - c = c + nTwin + plastic_phenopowerlaw_postResults(c+1_pInt:c+prm%totalNtwin) = & + stt%s_twin(1:prm%totalNtwin,of) + c = c + prm%totalNtwin case (accumulatedshear_twin_ID) - plastic_phenopowerlaw_postResults(c+1_pInt:c+nTwin) = & - plasticState(ph)%state(index_accshear_twin:index_accshear_twin+nTwin-1_pInt,of) - c = c + nTwin + plastic_phenopowerlaw_postResults(c+1_pInt:c+prm%totalNtwin) = & + stt%accshear_twin(1:prm%totalNtwin,of) + c = c + prm%totalNtwin + case (shearrate_twin_ID) j = 0_pInt - twinFamilies1: do f = 1_pInt,lattice_maxNtwinFamily + twinFamilies1: do f = 1_pInt,size(prm%Ntwin,1) index_myFamily = sum(lattice_NtwinSystem(1:f-1_pInt,ph)) ! at which index starts my family - twinSystems1: do i = 1_pInt,plastic_phenopowerlaw_Ntwin(f,instance) + twinSystems1: do i = 1_pInt,prm%Ntwin(f) j = j + 1_pInt tau = dot_product(Tstar_v,lattice_Stwin_v(1:6,index_myFamily+i,ph)) - plastic_phenopowerlaw_postResults(c+j) = (1.0_pReal-plasticState(ph)%state(index_F,of))*& ! 1-F - plastic_phenopowerlaw_gdot0_twin(instance)*& - (abs(tau)/plasticState(ph)%state(j+nSlip,of))**& - plastic_phenopowerlaw_n_twin(instance)*max(0.0_pReal,sign(1.0_pReal,tau)) + plastic_phenopowerlaw_postResults(c+j) = (1.0_pReal-stt%sumF(of))*& ! 1-F + prm%gdot0_twin*& + (abs(tau)/stt%s_twin(j,of))**& + prm%n_twin*max(0.0_pReal,sign(1.0_pReal,tau)) enddo twinSystems1 enddo twinFamilies1 - c = c + nTwin + c = c + prm%totalNtwin case (resolvedstress_twin_ID) j = 0_pInt - twinFamilies2: do f = 1_pInt,lattice_maxNtwinFamily + twinFamilies2: do f = 1_pInt,size(prm%Ntwin,1) index_myFamily = sum(lattice_NtwinSystem(1:f-1_pInt,ph)) ! at which index starts my family - twinSystems2: do i = 1_pInt,plastic_phenopowerlaw_Ntwin(f,instance) + twinSystems2: do i = 1_pInt,prm%Ntwin(f) j = j + 1_pInt plastic_phenopowerlaw_postResults(c+j) = & dot_product(Tstar_v,lattice_Stwin_v(1:6,index_myFamily+i,ph)) enddo twinSystems2 enddo twinFamilies2 - c = c + nTwin + c = c + prm%totalNtwin case (totalvolfrac_twin_ID) - plastic_phenopowerlaw_postResults(c+1_pInt) = plasticState(ph)%state(index_F,of) + plastic_phenopowerlaw_postResults(c+1_pInt) = stt%sumF(of) c = c + 1_pInt end select enddo outputsLoop - + end associate end function plastic_phenopowerlaw_postResults end module plastic_phenopowerlaw diff --git a/src/spectral_interface.f90 b/src/spectral_interface.f90 index d2adcf9ba..1ab92a178 100644 --- a/src/spectral_interface.f90 +++ b/src/spectral_interface.f90 @@ -19,21 +19,19 @@ module DAMASK_interface character(len=1024), public, protected :: & geometryFile = '', & !< parameter given for geometry file loadCaseFile = '' !< parameter given for load case file - character(len=1024), private :: workingDirectory !< accessed by getSolverWorkingDirectoryName for compatibility reasons + character(len=1024), private :: workingDirectory public :: & - getSolverWorkingDirectoryName, & getSolverJobName, & DAMASK_interface_init private :: & - storeWorkingDirectory, & + setWorkingDirectory, & getGeometryFile, & getLoadCaseFile, & rectifyPath, & makeRelativePath, & IIO_stringValue, & IIO_intValue, & - IIO_lc, & IIO_stringPos contains @@ -57,9 +55,9 @@ subroutine DAMASK_interface_init() implicit none character(len=1024) :: & commandLine, & !< command line call as string - loadCaseArg ='', & !< -l argument given to DAMASK_spectral.exe - geometryArg ='', & !< -g argument given to DAMASK_spectral.exe - workingDirArg ='', & !< -w argument given to DAMASK_spectral.exe + loadcaseArg = '', & !< -l argument given to DAMASK_spectral.exe + geometryArg = '', & !< -g argument given to DAMASK_spectral.exe + workingDirArg = '', & !< -w argument given to DAMASK_spectral.exe hostName, & !< name of machine on which DAMASK_spectral.exe is execute (might require export HOSTNAME) userName, & !< name of user calling DAMASK_spectral.exe tag @@ -114,7 +112,7 @@ subroutine DAMASK_interface_init() call date_and_time(values = dateAndTime) write(6,'(/,a)') ' <<<+- DAMASK_spectral -+>>>' - write(6,'(/,a)') ' Roters et al., Computational Materials Science, 2018' + write(6,'(a,/)') ' Roters et al., Computational Materials Science, 2018' write(6,'(/,a)') ' Version: '//DAMASKVERSION write(6,'(a,2(i2.2,a),i4.4)') ' Date: ',dateAndTime(3),'/',& dateAndTime(2),'/',& @@ -128,9 +126,8 @@ subroutine DAMASK_interface_init() call get_command(commandLine) chunkPos = IIO_stringPos(commandLine) - do i = 1, chunkPos(1) - tag = IIO_lc(IIO_stringValue(commandLine,chunkPos,i)) ! extract key - select case(tag) + do i = 2_pInt, chunkPos(1) + select case(IIO_stringValue(commandLine,chunkPos,i)) ! extract key case ('-h','--help') write(6,'(a)') ' #######################################################################' write(6,'(a)') ' DAMASK_spectral:' @@ -179,23 +176,25 @@ subroutine DAMASK_interface_init() write(6,'(a,/)')' Prints this message and exits' call quit(0_pInt) ! normal Termination case ('-l', '--load', '--loadcase') - loadcaseArg = IIO_stringValue(commandLine,chunkPos,i+1_pInt) + if ( i < chunkPos(1)) loadcaseArg = trim(IIO_stringValue(commandLine,chunkPos,i+1_pInt)) case ('-g', '--geom', '--geometry') - geometryArg = IIO_stringValue(commandLine,chunkPos,i+1_pInt) + if (i < chunkPos(1)) geometryArg = trim(IIO_stringValue(commandLine,chunkPos,i+1_pInt)) case ('-w', '-d', '--wd', '--directory', '--workingdir', '--workingdirectory') - workingDirArg = IIO_stringValue(commandLine,chunkPos,i+1_pInt) + if (i < chunkPos(1)) workingDirArg = trim(IIO_stringValue(commandLine,chunkPos,i+1_pInt)) case ('-r', '--rs', '--restart') - spectralRestartInc = IIO_IntValue(commandLine,chunkPos,i+1_pInt) - appendToOutFile = .true. + if (i < chunkPos(1)) then + spectralRestartInc = IIO_IntValue(commandLine,chunkPos,i+1_pInt) + appendToOutFile = .true. + endif end select enddo - - if (len(trim(loadcaseArg)) == 0 .or. len(trim(geometryArg)) == 0) then + + if (len_trim(loadcaseArg) == 0 .or. len_trim(geometryArg) == 0) then write(6,'(a)') ' Please specify geometry AND load case (-h for help)' call quit(1_pInt) endif - workingDirectory = trim(storeWorkingDirectory(trim(workingDirArg),trim(geometryArg))) + workingDirectory = trim(setWorkingDirectory(trim(workingDirArg))) geometryFile = getGeometryFile(geometryArg) loadCaseFile = getLoadCaseFile(loadCaseArg) @@ -208,7 +207,7 @@ subroutine DAMASK_interface_init() write(6,'(a,a)') ' Working dir argument: ', trim(workingDirArg) write(6,'(a,a)') ' Geometry argument: ', trim(geometryArg) write(6,'(a,a)') ' Loadcase argument: ', trim(loadcaseArg) - write(6,'(a,a)') ' Working directory: ', trim(getSolverWorkingDirectoryName()) + write(6,'(a,a)') ' Working directory: ', trim(workingDirectory) write(6,'(a,a)') ' Geometry file: ', trim(geometryFile) write(6,'(a,a)') ' Loadcase file: ', trim(loadCaseFile) write(6,'(a,a)') ' Solver job name: ', trim(getSolverJobName()) @@ -222,59 +221,39 @@ end subroutine DAMASK_interface_init !-------------------------------------------------------------------------------------------------- !> @brief extract working directory from given argument or from location of geometry file, !! possibly converting relative arguments to absolut path -!> @todo change working directory with call chdir(storeWorkingDirectory)? !-------------------------------------------------------------------------------------------------- -character(len=1024) function storeWorkingDirectory(workingDirectoryArg,geometryArg) +character(len=1024) function setWorkingDirectory(workingDirectoryArg) use system_routines, only: & - isDirectory, & - getCWD + getCWD, & + setCWD implicit none character(len=*), intent(in) :: workingDirectoryArg !< working directory argument - character(len=*), intent(in) :: geometryArg !< geometry argument - character(len=1024) :: cwd logical :: error external :: quit wdGiven: if (len(workingDirectoryArg)>0) then absolutePath: if (workingDirectoryArg(1:1) == '/') then - storeWorkingDirectory = workingDirectoryArg + setWorkingDirectory = workingDirectoryArg else absolutePath - error = getCWD(cwd) + error = getCWD(setWorkingDirectory) if (error) call quit(1_pInt) - storeWorkingDirectory = trim(cwd)//'/'//workingDirectoryArg + setWorkingDirectory = trim(setWorkingDirectory)//'/'//workingDirectoryArg endif absolutePath - if (storeWorkingDirectory(len(trim(storeWorkingDirectory)):len(trim(storeWorkingDirectory))) /= '/') & - storeWorkingDirectory = trim(storeWorkingDirectory)//'/' ! if path seperator is not given, append it else wdGiven - if (geometryArg(1:1) == '/') then ! absolute path given as command line argument - storeWorkingDirectory = geometryArg(1:scan(geometryArg,'/',back=.true.)) - else - error = getCWD(cwd) ! relative path given as command line argument - if (error) call quit(1_pInt) - storeWorkingDirectory = trim(cwd)//'/'//geometryArg(1:scan(geometryArg,'/',back=.true.)) - endif + error = getCWD(setWorkingDirectory) ! relative path given as command line argument + if (error) call quit(1_pInt) endif wdGiven - storeWorkingDirectory = trim(rectifyPath(storeWorkingDirectory)) - if(.not. isDirectory(trim(storeWorkingDirectory))) then ! check if the directory exists - write(6,'(a20,a,a16)') ' working directory "',trim(storeWorkingDirectory),'" does not exist' - call quit(1_pInt) + setWorkingDirectory = trim(rectifyPath(setWorkingDirectory)) + + error = setCWD(trim(setWorkingDirectory)) + if(error) then + write(6,'(a20,a,a16)') ' working directory "',trim(setWorkingDirectory),'" does not exist' + call quit(1_pInt) endif -end function storeWorkingDirectory - - -!-------------------------------------------------------------------------------------------------- -!> @brief simply returns the private string workingDir -!-------------------------------------------------------------------------------------------------- -character(len=1024) function getSolverWorkingDirectoryName() - - implicit none - - getSolverWorkingDirectoryName = workingDirectory - -end function getSolverWorkingDirectoryName +end function setWorkingDirectory !-------------------------------------------------------------------------------------------------- @@ -306,32 +285,23 @@ end function getSolverJobName !> @brief basename of geometry file with extension from command line arguments !-------------------------------------------------------------------------------------------------- character(len=1024) function getGeometryFile(geometryParameter) - use system_routines, only: & - getCWD implicit none character(len=1024), intent(in) :: & geometryParameter - character(len=1024) :: & - cwd integer :: posExt, posSep - logical :: error external :: quit - getGeometryFile = geometryParameter + getGeometryFile = trim(geometryParameter) posExt = scan(getGeometryFile,'.',back=.true.) posSep = scan(getGeometryFile,'/',back=.true.) - if (posExt <= posSep) getGeometryFile = trim(getGeometryFile)//('.geom') ! no extension present - if (scan(getGeometryFile,'/') /= 1) then ! relative path given as command line argument - error = getcwd(cwd) - if (error) call quit(1_pInt) - getGeometryFile = rectifyPath(trim(cwd)//'/'//getGeometryFile) - else - getGeometryFile = rectifyPath(getGeometryFile) - endif + if (posExt <= posSep) getGeometryFile = trim(getGeometryFile)//('.geom') + if (scan(getGeometryFile,'/') /= 1) & + getGeometryFile = trim(workingDirectory)//'/'//trim(getGeometryFile) + + getGeometryFile = makeRelativePath(workingDirectory, getGeometryFile) - getGeometryFile = makeRelativePath(getSolverWorkingDirectoryName(), getGeometryFile) end function getGeometryFile @@ -340,38 +310,29 @@ end function getGeometryFile !> @brief relative path of loadcase from command line arguments !-------------------------------------------------------------------------------------------------- character(len=1024) function getLoadCaseFile(loadCaseParameter) - use system_routines, only: & - getCWD implicit none character(len=1024), intent(in) :: & loadCaseParameter - character(len=1024) :: & - cwd integer :: posExt, posSep - logical :: error external :: quit - getLoadCaseFile = loadcaseParameter + getLoadCaseFile = trim(loadCaseParameter) posExt = scan(getLoadCaseFile,'.',back=.true.) posSep = scan(getLoadCaseFile,'/',back=.true.) - if (posExt <= posSep) getLoadCaseFile = trim(getLoadCaseFile)//('.load') ! no extension present - if (scan(getLoadCaseFile,'/') /= 1) then ! relative path given as command line argument - error = getcwd(cwd) - if (error) call quit(1_pInt) - getLoadCaseFile = rectifyPath(trim(cwd)//'/'//getLoadCaseFile) - else - getLoadCaseFile = rectifyPath(getLoadCaseFile) - endif + if (posExt <= posSep) getLoadCaseFile = trim(getLoadCaseFile)//('.load') + if (scan(getLoadCaseFile,'/') /= 1) & + getLoadCaseFile = trim(workingDirectory)//'/'//trim(getLoadCaseFile) - getLoadCaseFile = makeRelativePath(getSolverWorkingDirectoryName(), getLoadCaseFile) + getLoadCaseFile = makeRelativePath(workingDirectory, getLoadCaseFile) end function getLoadCaseFile !-------------------------------------------------------------------------------------------------- -!> @brief remove ../ and /./ from path +!> @brief remove ../, /./, and // from path. +!> @details works only if absolute path is given !-------------------------------------------------------------------------------------------------- function rectifyPath(path) @@ -385,14 +346,21 @@ function rectifyPath(path) l = len_trim(path) rectifyPath = path do i = l,3,-1 - if (rectifyPath(i-2:i) == '/'//'.'//'/') & - rectifyPath(i-1:l) = rectifyPath(i+1:l)//' ' + if (rectifyPath(i-2:i) == '/./') rectifyPath(i-1:l) = rectifyPath(i+1:l)//' ' + enddo + +!-------------------------------------------------------------------------------------------------- +! remove // from path + l = len_trim(path) + rectifyPath = path + do i = l,2,-1 + if (rectifyPath(i-1:i) == '//') rectifyPath(i-1:l) = rectifyPath(i:l)//' ' enddo !-------------------------------------------------------------------------------------------------- ! remove ../ and corresponding directory from rectifyPath l = len_trim(rectifyPath) - i = index(rectifyPath(i:l),'..'//'/') + i = index(rectifyPath(i:l),'../') j = 0 do while (i > j) j = scan(rectifyPath(1:i-2),'/',back=.true.) @@ -402,7 +370,7 @@ function rectifyPath(path) rectifyPath(j+1:k-1) = rectifyPath(j+2:k) rectifyPath(k:k) = ' ' endif - i = j+index(rectifyPath(j+1:l),'..'//'/') + i = j+index(rectifyPath(j+1:l),'../') enddo if(len_trim(rectifyPath) == 0) rectifyPath = '/' @@ -415,20 +383,24 @@ end function rectifyPath character(len=1024) function makeRelativePath(a,b) implicit none - character (len=*) :: a,b + character (len=*), intent(in) :: a,b + character (len=1024) :: a_cleaned,b_cleaned integer :: i,posLastCommonSlash,remainingSlashes !no pInt posLastCommonSlash = 0 remainingSlashes = 0 + a_cleaned = rectifyPath(trim(a)//'/') + b_cleaned = rectifyPath(b) - do i = 1, min(1024,len_trim(a),len_trim(b)) - if (a(i:i) /= b(i:i)) exit - if (a(i:i) == '/') posLastCommonSlash = i + do i = 1, min(1024,len_trim(a_cleaned),len_trim(rectifyPath(b_cleaned))) + if (a_cleaned(i:i) /= b_cleaned(i:i)) exit + if (a_cleaned(i:i) == '/') posLastCommonSlash = i enddo - do i = posLastCommonSlash+1,len_trim(a) - if (a(i:i) == '/') remainingSlashes = remainingSlashes + 1 + do i = posLastCommonSlash+1,len_trim(a_cleaned) + if (a_cleaned(i:i) == '/') remainingSlashes = remainingSlashes + 1 enddo - makeRelativePath = repeat('..'//'/',remainingSlashes)//b(posLastCommonSlash+1:len_trim(b)) + + makeRelativePath = repeat('..'//'/',remainingSlashes)//b_cleaned(posLastCommonSlash+1:len_trim(b_cleaned)) end function makeRelativePath @@ -439,17 +411,12 @@ end function makeRelativePath pure function IIO_stringValue(string,chunkPos,myChunk) implicit none - integer(pInt), dimension(:), intent(in) :: chunkPos !< positions of start and end of each tag/chunk in given string - integer(pInt), intent(in) :: myChunk !< position number of desired chunk - character(len=1+chunkPos(myChunk*2+1)-chunkPos(myChunk*2)) :: IIO_stringValue - character(len=*), intent(in) :: string !< raw input with known start and end of each chunk + integer(pInt), dimension(:), intent(in) :: chunkPos !< positions of start and end of each tag/chunk in given string + integer(pInt), intent(in) :: myChunk !< position number of desired chunk + character(len=chunkPos(myChunk*2+1)-chunkPos(myChunk*2)+1) :: IIO_stringValue + character(len=*), intent(in) :: string !< raw input with known start and end of each chunk - - valuePresent: if (myChunk > chunkPos(1) .or. myChunk < 1_pInt) then - IIO_stringValue = '' - else valuePresent - IIO_stringValue = string(chunkPos(myChunk*2):chunkPos(myChunk*2+1)) - endif valuePresent + IIO_stringValue = string(chunkPos(myChunk*2):chunkPos(myChunk*2+1)) end function IIO_stringValue @@ -476,29 +443,6 @@ integer(pInt) pure function IIO_intValue(string,chunkPos,myChunk) end function IIO_intValue -!-------------------------------------------------------------------------------------------------- -!> @brief taken from IO, check IO_lc for documentation -!-------------------------------------------------------------------------------------------------- -pure function IIO_lc(string) - - implicit none - character(len=*), intent(in) :: string !< string to convert - character(len=len(string)) :: IIO_lc - - character(26), parameter :: LOWER = 'abcdefghijklmnopqrstuvwxyz' - character(26), parameter :: UPPER = 'ABCDEFGHIJKLMNOPQRSTUVWXYZ' - - integer :: i,n ! no pInt (len returns default integer) - - IIO_lc = string - do i=1,len(string) - n = index(UPPER,IIO_lc(i:i)) - if (n/=0) IIO_lc(i:i) = LOWER(n:n) - enddo - -end function IIO_lc - - !-------------------------------------------------------------------------------------------------- !> @brief taken from IO, check IO_stringPos for documentation !--------------------------------------------------------------------------------------------------