diff --git a/src/CPFEM.f90 b/src/CPFEM.f90 index 1f76e6c25..172b79122 100644 --- a/src/CPFEM.f90 +++ b/src/CPFEM.f90 @@ -26,7 +26,7 @@ module CPFEM implicit none private - + real(pReal), parameter, private :: & CPFEM_odd_stress = 1e15_pReal, & !< return value for stress in case of ping pong dummy cycle CPFEM_odd_jacobian = 1e50_pReal !< return value for jacobian in case of ping pong dummy cycle @@ -42,8 +42,8 @@ module CPFEM lastLovl = 0_pInt !< lovl in previous call to marc hypela2 real(pReal), public :: & theTime = 0.0_pReal, & !< needs description - theDelta = 0.0_pReal - logical, public :: & + theDelta = 0.0_pReal + logical, public :: & outdatedFFN1 = .false., & !< needs description lastIncConverged = .false., & !< needs description outdatedByNewInc = .false. !< needs description @@ -72,31 +72,32 @@ contains !> @brief call (thread safe) all module initializations !-------------------------------------------------------------------------------------------------- subroutine CPFEM_initAll(el,ip) - integer(pInt), intent(in) :: el, & !< FE el number - ip !< FE integration point number + + integer(pInt), intent(in) :: el, & !< FE el number + ip !< FE integration point number - !$OMP CRITICAL (init) - if (.not. CPFEM_init_done) then - call DAMASK_interface_init ! Spectral and FEM interface to commandline - call prec_init - call IO_init - call numerics_init - call debug_init - call config_init - call math_init - call rotations_init - call HDF5_utilities_init - call results_init - call mesh_init(ip, el) - call lattice_init - call material_init - call constitutive_init - call crystallite_init - call homogenization_init - call CPFEM_init - CPFEM_init_done = .true. - endif - !$OMP END CRITICAL (init) + !$OMP CRITICAL(init) + if (.not. CPFEM_init_done) then + call DAMASK_interface_init + call prec_init + call IO_init + call numerics_init + call debug_init + call config_init + call math_init + call rotations_init + call HDF5_utilities_init + call results_init + call mesh_init(ip, el) + call lattice_init + call material_init + call constitutive_init + call crystallite_init + call homogenization_init + call CPFEM_init + CPFEM_init_done = .true. + endif + !$OMP END CRITICAL(init) end subroutine CPFEM_initAll @@ -174,35 +175,7 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt CPFEM_dcsde = CPFEM_dcsde_knownGood !*** age results - if (iand(mode, CPFEM_AGERESULTS) /= 0_pInt) then - crystallite_F0 = crystallite_partionedF ! crystallite deformation - crystallite_Fp0 = crystallite_Fp ! crystallite plastic deformation - crystallite_Lp0 = crystallite_Lp ! crystallite plastic velocity - crystallite_Fi0 = crystallite_Fi ! crystallite intermediate deformation - crystallite_Li0 = crystallite_Li ! crystallite intermediate velocity - crystallite_S0 = crystallite_S ! crystallite 2nd Piola Kirchhoff stress - - forall (i = 1:size(plasticState)) plasticState(i)%state0 = plasticState(i)%state - do i = 1, size(sourceState) - do mySource = 1,phase_Nsources(i) - sourceState(i)%p(mySource)%state0 = sourceState(i)%p(mySource)%state - enddo; enddo - if (iand(debug_level(debug_CPFEM), debug_levelBasic) /= 0_pInt) then - write(6,'(a)') '<< CPFEM >> aging states' - if (debug_e <= discretization_nElem .and. debug_i <=discretization_nIP) then - write(6,'(a,1x,i8,1x,i2,1x,i4,/,(12x,6(e20.8,1x)),/)') & - '<< CPFEM >> aged state of elFE ip grain',debug_e, debug_i, 1, & - plasticState(material_phaseAt(1,debug_e))%state(:,material_phasememberAt(1,debug_i,debug_e)) - endif - endif - - do homog = 1_pInt, material_Nhomogenization - homogState (homog)%state0 = homogState (homog)%state - thermalState (homog)%state0 = thermalState (homog)%state - damageState (homog)%state0 = damageState (homog)%state - enddo - endif - + if (iand(mode, CPFEM_AGERESULTS) /= 0_pInt) call CPFEM_forward !*** collection of FEM input with returning of randomize odd stress and jacobian @@ -358,9 +331,18 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt end subroutine CPFEM_general +!-------------------------------------------------------------------------------------------------- +!> @brief Forward data for new time increment. +!-------------------------------------------------------------------------------------------------- +subroutine CPFEM_forward + + call crystallite_forward + +end subroutine CPFEM_forward + !-------------------------------------------------------------------------------------------------- -!> @brief triggers writing of the results +!> @brief Trigger writing of results. !-------------------------------------------------------------------------------------------------- subroutine CPFEM_results(inc,time) diff --git a/src/CPFEM2.f90 b/src/CPFEM2.f90 index 68296697d..882fe8ae3 100644 --- a/src/CPFEM2.f90 +++ b/src/CPFEM2.f90 @@ -64,131 +64,45 @@ end subroutine CPFEM_initAll !-------------------------------------------------------------------------------------------------- -!> @brief allocate the arrays defined in module CPFEM and initialize them +!> @brief Read restart information if needed. !-------------------------------------------------------------------------------------------------- subroutine CPFEM_init - - integer :: i - integer(HID_T) :: fileHandle, groupHandle - character(len=pStringLen) :: fileName, datasetName write(6,'(/,a)') ' <<<+- CPFEM init -+>>>'; flush(6) - if (interface_restartInc > 0) then - write(6,'(/,a,i0,a)') ' reading restart information of increment ', interface_restartInc, ' from file' - - write(fileName,'(a,i0,a)') trim(getSolverJobName())//'_',worldrank,'.hdf5' - fileHandle = HDF5_openFile(fileName) - - call HDF5_read(fileHandle,crystallite_F0, 'F') - call HDF5_read(fileHandle,crystallite_Fp0,'Fp') - call HDF5_read(fileHandle,crystallite_Fi0,'Fi') - call HDF5_read(fileHandle,crystallite_Lp0,'Lp') - call HDF5_read(fileHandle,crystallite_Li0,'Li') - call HDF5_read(fileHandle,crystallite_S0, 'S') - - groupHandle = HDF5_openGroup(fileHandle,'constituent') - do i = 1,size(phase_plasticity) - write(datasetName,'(i0,a)') i,'_omega_plastic' - call HDF5_read(groupHandle,plasticState(i)%state0,datasetName) - enddo - call HDF5_closeGroup(groupHandle) - - groupHandle = HDF5_openGroup(fileHandle,'materialpoint') - do i = 1, material_Nhomogenization - write(datasetName,'(i0,a)') i,'_omega_homogenization' - call HDF5_read(groupHandle,homogState(i)%state0,datasetName) - enddo - call HDF5_closeGroup(groupHandle) - - call HDF5_closeFile(fileHandle) - endif + if (interface_restartInc > 0) call crystallite_restartRead end subroutine CPFEM_init !-------------------------------------------------------------------------------------------------- -!> @brief Forward data after successful increment. -! ToDo: Any guessing for the current states possible? -!-------------------------------------------------------------------------------------------------- -subroutine CPFEM_forward - - integer :: i, j - - if (iand(debug_level(debug_CPFEM), debug_levelBasic) /= 0) & - write(6,'(a)') '<< CPFEM >> aging states' - - crystallite_F0 = crystallite_partionedF - crystallite_Fp0 = crystallite_Fp - crystallite_Lp0 = crystallite_Lp - crystallite_Fi0 = crystallite_Fi - crystallite_Li0 = crystallite_Li - crystallite_S0 = crystallite_S - - do i = 1, size(plasticState) - plasticState(i)%state0 = plasticState(i)%state - enddo - do i = 1, size(sourceState) - do j = 1,phase_Nsources(i) - sourceState(i)%p(j)%state0 = sourceState(i)%p(j)%state - enddo; enddo - do i = 1, material_Nhomogenization - homogState (i)%state0 = homogState (i)%state - thermalState(i)%state0 = thermalState(i)%state - damageState (i)%state0 = damageState (i)%state - enddo - -end subroutine CPFEM_forward - - -!-------------------------------------------------------------------------------------------------- -!> @brief Write current restart information (Field and constitutive data) to file. +!> @brief Write restart information. !-------------------------------------------------------------------------------------------------- subroutine CPFEM_restartWrite - integer :: i - integer(HID_T) :: fileHandle, groupHandle - character(len=pStringLen) :: fileName, datasetName - - write(6,'(a)') ' writing field and constitutive data required for restart to file';flush(6) - - write(fileName,'(a,i0,a)') trim(getSolverJobName())//'_',worldrank,'.hdf5' - fileHandle = HDF5_openFile(fileName,'a') - - call HDF5_write(fileHandle,crystallite_partionedF,'F') - call HDF5_write(fileHandle,crystallite_Fp, 'Fp') - call HDF5_write(fileHandle,crystallite_Fi, 'Fi') - call HDF5_write(fileHandle,crystallite_Lp, 'Lp') - call HDF5_write(fileHandle,crystallite_Li, 'Li') - call HDF5_write(fileHandle,crystallite_S, 'S') - - groupHandle = HDF5_addGroup(fileHandle,'constituent') - do i = 1,size(phase_plasticity) - write(datasetName,'(i0,a)') i,'_omega_plastic' - call HDF5_write(groupHandle,plasticState(i)%state,datasetName) - enddo - call HDF5_closeGroup(groupHandle) - - groupHandle = HDF5_addGroup(fileHandle,'materialpoint') - do i = 1, material_Nhomogenization - write(datasetName,'(i0,a)') i,'_omega_homogenization' - call HDF5_write(groupHandle,homogState(i)%state,datasetName) - enddo - call HDF5_closeGroup(groupHandle) - - call HDF5_closeFile(fileHandle) + call crystallite_restartWrite end subroutine CPFEM_restartWrite +!-------------------------------------------------------------------------------------------------- +!> @brief Forward data for new time increment. +!-------------------------------------------------------------------------------------------------- +subroutine CPFEM_forward + + call crystallite_forward + +end subroutine CPFEM_forward + + !-------------------------------------------------------------------------------------------------- !> @brief Trigger writing of results. !-------------------------------------------------------------------------------------------------- subroutine CPFEM_results(inc,time) - + integer, intent(in) :: inc real(pReal), intent(in) :: time - + call results_openJobFile call results_addIncrement(inc,time) call constitutive_results diff --git a/src/DAMASK_marc.f90 b/src/DAMASK_marc.f90 index f5e3c4744..782657da6 100644 --- a/src/DAMASK_marc.f90 +++ b/src/DAMASK_marc.f90 @@ -40,9 +40,9 @@ module DAMASK_interface implicit none private - logical, public :: symmetricSolver + logical, public :: symmetricSolver character(len=*), parameter, public :: INPUTFILEEXTENSION = '.dat' - + public :: & DAMASK_interface_init, & getSolverJobName @@ -57,14 +57,14 @@ subroutine DAMASK_interface_init integer, dimension(8) :: dateAndTime integer :: ierr character(len=pPathLen) :: wd - + write(6,'(/,a)') ' <<<+- DAMASK_marc init -+>>>' - + write(6,'(/,a)') ' Roters et al., Computational Materials Science 158:420–478, 2019' write(6,'(a)') ' https://doi.org/10.1016/j.commatsci.2018.04.030' - + write(6,'(/,a)') ' Version: '//DAMASKVERSION - + ! https://github.com/jeffhammond/HPCInfo/blob/master/docs/Preprocessor-Macros.md #if __INTEL_COMPILER >= 1800 write(6,'(/,a)') ' Compiled with: '//compiler_version() @@ -73,13 +73,13 @@ subroutine DAMASK_interface_init write(6,'(/,a,i4.4,a,i8.8)') ' Compiled with Intel fortran version :', __INTEL_COMPILER,& ', build date :', __INTEL_COMPILER_BUILD_DATE #endif - + write(6,'(/,a)') ' Compiled on: '//__DATE__//' at '//__TIME__ - + call date_and_time(values = dateAndTime) write(6,'(/,a,2(i2.2,a),i4.4)') ' Date: ',dateAndTime(3),'/',dateAndTime(2),'/', dateAndTime(1) write(6,'(a,2(i2.2,a),i2.2)') ' Time: ',dateAndTime(5),':', dateAndTime(6),':', dateAndTime(7) - + inquire(5, name=wd) wd = wd(1:scan(wd,'/',back=.true.)) ierr = CHDIR(wd) diff --git a/src/crystallite.f90 b/src/crystallite.f90 index 1a8252cce..d1f073a91 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -11,6 +11,8 @@ module crystallite use prec use IO + use HDF5_utilities + use DAMASK_interface use config use debug use numerics @@ -36,25 +38,25 @@ module crystallite crystallite_orientation !< current orientation real(pReal), dimension(:,:,:,:,:), allocatable, public, protected :: & crystallite_Fe, & !< current "elastic" def grad (end of converged time step) - crystallite_P !< 1st Piola-Kirchhoff stress per grain + crystallite_P, & !< 1st Piola-Kirchhoff stress per grain + crystallite_S0, & !< 2nd Piola-Kirchhoff stress vector at start of FE inc + crystallite_Fp0, & !< plastic def grad at start of FE inc + crystallite_Fi0, & !< intermediate def grad at start of FE inc + crystallite_F0, & !< def grad at start of FE inc + crystallite_Lp0, & !< plastic velocitiy grad at start of FE inc + crystallite_Li0 !< intermediate velocitiy grad at start of FE inc real(pReal), dimension(:,:,:,:,:), allocatable, public :: & crystallite_S, & !< current 2nd Piola-Kirchhoff stress vector (end of converged time step) - crystallite_S0, & !< 2nd Piola-Kirchhoff stress vector at start of FE inc crystallite_partionedS0, & !< 2nd Piola-Kirchhoff stress vector at start of homog inc crystallite_Fp, & !< current plastic def grad (end of converged time step) - crystallite_Fp0, & !< plastic def grad at start of FE inc crystallite_partionedFp0,& !< plastic def grad at start of homog inc crystallite_Fi, & !< current intermediate def grad (end of converged time step) - crystallite_Fi0, & !< intermediate def grad at start of FE inc crystallite_partionedFi0,& !< intermediate def grad at start of homog inc - crystallite_F0, & !< def grad at start of FE inc crystallite_partionedF, & !< def grad to be reached at end of homog inc crystallite_partionedF0, & !< def grad at start of homog inc crystallite_Lp, & !< current plastic velocitiy grad (end of converged time step) - crystallite_Lp0, & !< plastic velocitiy grad at start of FE inc crystallite_partionedLp0, & !< plastic velocity grad at start of homog inc crystallite_Li, & !< current intermediate velocitiy grad (end of converged time step) - crystallite_Li0, & !< intermediate velocitiy grad at start of FE inc crystallite_partionedLi0 !< intermediate velocity grad at start of homog inc real(pReal), dimension(:,:,:,:,:), allocatable :: & crystallite_subFp0,& !< plastic def grad at start of crystallite inc @@ -104,7 +106,10 @@ module crystallite crystallite_stressTangent, & crystallite_orientations, & crystallite_push33ToRef, & - crystallite_results + crystallite_results, & + crystallite_restartWrite, & + crystallite_restartRead, & + crystallite_forward contains @@ -130,38 +135,30 @@ subroutine crystallite_init iMax = discretization_nIP eMax = discretization_nElem - allocate(crystallite_S0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_partionedS0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_S(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_P(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_F0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_partionedF0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_partionedF(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_subF0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_subF(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_Fp0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_partionedFp0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_subFp0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_Fp(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_Fi0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_partionedFi0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_subFi0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_Fi(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_Fe(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_Lp0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_partionedLp0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_subLp0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_Lp(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_Li0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_partionedLi0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_subLi0(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_Li(3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_dPdF(3,3,3,3,cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_dt(cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_subdt(cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_subFrac(cMax,iMax,eMax), source=0.0_pReal) - allocate(crystallite_subStep(cMax,iMax,eMax), source=0.0_pReal) + allocate(crystallite_partionedF(3,3,cMax,iMax,eMax),source=0.0_pReal) + + allocate(crystallite_S0, & + crystallite_F0, crystallite_Fi0,crystallite_Fp0, & + crystallite_Li0,crystallite_Lp0, & + crystallite_partionedS0, & + crystallite_partionedF0,crystallite_partionedFp0,crystallite_partionedFi0, & + crystallite_partionedLp0,crystallite_partionedLi0, & + crystallite_S,crystallite_P, & + crystallite_Fe,crystallite_Fi,crystallite_Fp, & + crystallite_Li,crystallite_Lp, & + crystallite_subF,crystallite_subF0, & + crystallite_subFp0,crystallite_subFi0, & + crystallite_subLi0,crystallite_subLp0, & + source = crystallite_partionedF) + + allocate(crystallite_dPdF(3,3,3,3,cMax,iMax,eMax),source=0.0_pReal) + + allocate(crystallite_dt(cMax,iMax,eMax),source=0.0_pReal) + allocate(crystallite_subdt,crystallite_subFrac,crystallite_subStep, & + source = crystallite_dt) + allocate(crystallite_orientation(cMax,iMax,eMax)) + allocate(crystallite_localPlasticity(cMax,iMax,eMax), source=.true.) allocate(crystallite_requested(cMax,iMax,eMax), source=.false.) allocate(crystallite_todo(cMax,iMax,eMax), source=.false.) @@ -1844,4 +1841,117 @@ logical function stateJump(ipc,ip,el) end function stateJump + +!-------------------------------------------------------------------------------------------------- +!> @brief Write current restart information (Field and constitutive data) to file. +! ToDo: Merge data into one file for MPI, move state to constitutive and homogenization, respectively +!-------------------------------------------------------------------------------------------------- +subroutine crystallite_restartWrite + + integer :: i + integer(HID_T) :: fileHandle, groupHandle + character(len=pStringLen) :: fileName, datasetName + + write(6,'(a)') ' writing field and constitutive data required for restart to file';flush(6) + + write(fileName,'(a,i0,a)') trim(getSolverJobName())//'_',worldrank,'.hdf5' + fileHandle = HDF5_openFile(fileName,'a') + + call HDF5_write(fileHandle,crystallite_partionedF,'F') + call HDF5_write(fileHandle,crystallite_Fp, 'Fp') + call HDF5_write(fileHandle,crystallite_Fi, 'Fi') + call HDF5_write(fileHandle,crystallite_Lp, 'Lp') + call HDF5_write(fileHandle,crystallite_Li, 'Li') + call HDF5_write(fileHandle,crystallite_S, 'S') + + groupHandle = HDF5_addGroup(fileHandle,'constituent') + do i = 1,size(phase_plasticity) + write(datasetName,'(i0,a)') i,'_omega_plastic' + call HDF5_write(groupHandle,plasticState(i)%state,datasetName) + enddo + call HDF5_closeGroup(groupHandle) + + groupHandle = HDF5_addGroup(fileHandle,'materialpoint') + do i = 1, material_Nhomogenization + write(datasetName,'(i0,a)') i,'_omega_homogenization' + call HDF5_write(groupHandle,homogState(i)%state,datasetName) + enddo + call HDF5_closeGroup(groupHandle) + + call HDF5_closeFile(fileHandle) + +end subroutine crystallite_restartWrite + + +!-------------------------------------------------------------------------------------------------- +!> @brief Read data for restart +! ToDo: Merge data into one file for MPI, move state to constitutive and homogenization, respectively +!-------------------------------------------------------------------------------------------------- +subroutine crystallite_restartRead + + integer :: i + integer(HID_T) :: fileHandle, groupHandle + character(len=pStringLen) :: fileName, datasetName + + write(6,'(/,a,i0,a)') ' reading restart information of increment from file' + + write(fileName,'(a,i0,a)') trim(getSolverJobName())//'_',worldrank,'.hdf5' + fileHandle = HDF5_openFile(fileName) + + call HDF5_read(fileHandle,crystallite_F0, 'F') + call HDF5_read(fileHandle,crystallite_Fp0,'Fp') + call HDF5_read(fileHandle,crystallite_Fi0,'Fi') + call HDF5_read(fileHandle,crystallite_Lp0,'Lp') + call HDF5_read(fileHandle,crystallite_Li0,'Li') + call HDF5_read(fileHandle,crystallite_S0, 'S') + + groupHandle = HDF5_openGroup(fileHandle,'constituent') + do i = 1,size(phase_plasticity) + write(datasetName,'(i0,a)') i,'_omega_plastic' + call HDF5_read(groupHandle,plasticState(i)%state0,datasetName) + enddo + call HDF5_closeGroup(groupHandle) + + groupHandle = HDF5_openGroup(fileHandle,'materialpoint') + do i = 1, material_Nhomogenization + write(datasetName,'(i0,a)') i,'_omega_homogenization' + call HDF5_read(groupHandle,homogState(i)%state0,datasetName) + enddo + call HDF5_closeGroup(groupHandle) + + call HDF5_closeFile(fileHandle) + +end subroutine crystallite_restartRead + + +!-------------------------------------------------------------------------------------------------- +!> @brief Forward data after successful increment. +! ToDo: Any guessing for the current states possible? +!-------------------------------------------------------------------------------------------------- +subroutine crystallite_forward + + integer :: i, j + + crystallite_F0 = crystallite_partionedF + crystallite_Fp0 = crystallite_Fp + crystallite_Lp0 = crystallite_Lp + crystallite_Fi0 = crystallite_Fi + crystallite_Li0 = crystallite_Li + crystallite_S0 = crystallite_S + + do i = 1, size(plasticState) + plasticState(i)%state0 = plasticState(i)%state + enddo + do i = 1, size(sourceState) + do j = 1,phase_Nsources(i) + sourceState(i)%p(j)%state0 = sourceState(i)%p(j)%state + enddo; enddo + do i = 1, material_Nhomogenization + homogState (i)%state0 = homogState (i)%state + thermalState(i)%state0 = thermalState(i)%state + damageState (i)%state0 = damageState (i)%state + enddo + +end subroutine crystallite_forward + end module crystallite