From f9edeb40a5a3e3fb2a15584365ab16c90022c91f Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 17 Jul 2021 11:50:21 +0200 Subject: [PATCH] descriptive names --- src/CPFEM.f90 | 4 +- src/grid/grid_mech_FEM.f90 | 2 +- src/grid/grid_mech_spectral_basic.f90 | 2 +- src/grid/grid_mech_spectral_polarisation.f90 | 2 +- src/grid/spectral_utilities.f90 | 6 +- src/homogenization.f90 | 74 ++++++++------------ src/homogenization_mechanical.f90 | 8 +-- src/mesh/FEM_utilities.f90 | 4 +- src/phase.f90 | 14 ++-- src/phase_damage.f90 | 14 ++-- src/phase_damage_anisobrittle.f90 | 18 ++--- src/phase_mechanical.f90 | 27 +++---- src/phase_mechanical_eigen.f90 | 4 +- 13 files changed, 83 insertions(+), 96 deletions(-) diff --git a/src/CPFEM.f90 b/src/CPFEM.f90 index 06d1ca6a8..aa532859a 100644 --- a/src/CPFEM.f90 +++ b/src/CPFEM.f90 @@ -190,9 +190,9 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature_inp, dt, elFE, ip, cauchyS else validCalculation if (debugCPFEM%extensive) print'(a,i8,1x,i2)', '<< CPFEM >> calculation for elFE ip ',elFE,ip - call materialpoint_stressAndItsTangent(dt,[ip,ip],[elCP,elCP]) + call homogenization_mechanical_response(dt,[ip,ip],[elCP,elCP]) if (.not. terminallyIll) & - call materialpoint_stressAndItsTangent2(dt,[ip,ip],[elCP,elCP]) + call homogenization_mechanical_response2(dt,[ip,ip],[elCP,elCP]) terminalIllness: if (terminallyIll) then diff --git a/src/grid/grid_mech_FEM.f90 b/src/grid/grid_mech_FEM.f90 index d2fd89b2b..111e9fece 100644 --- a/src/grid/grid_mech_FEM.f90 +++ b/src/grid/grid_mech_FEM.f90 @@ -261,7 +261,7 @@ subroutine grid_mechanical_FEM_init F = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid3) endif restartRead - homogenization_F0 = reshape(F_lastInc, [3,3,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent + homogenization_F0 = reshape(F_lastInc, [3,3,product(grid(1:2))*grid3]) ! set starting condition for homogenization_mechanical_response call utilities_updateCoords(F) call utilities_constitutiveResponse(P_current,P_av,C_volAvg,devNull, & ! stress field, stress avg, global average of stiffness and (min+max)/2 F, & ! target F diff --git a/src/grid/grid_mech_spectral_basic.f90 b/src/grid/grid_mech_spectral_basic.f90 index 3ec14d00e..c4fce61a7 100644 --- a/src/grid/grid_mech_spectral_basic.f90 +++ b/src/grid/grid_mech_spectral_basic.f90 @@ -211,7 +211,7 @@ subroutine grid_mechanical_spectral_basic_init F = reshape(F_lastInc,[9,grid(1),grid(2),grid3]) endif restartRead - homogenization_F0 = reshape(F_lastInc, [3,3,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent + homogenization_F0 = reshape(F_lastInc, [3,3,product(grid(1:2))*grid3]) ! set starting condition for homogenization_mechanical_response call utilities_updateCoords(reshape(F,shape(F_lastInc))) call utilities_constitutiveResponse(P,P_av,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 reshape(F,shape(F_lastInc)), & ! target F diff --git a/src/grid/grid_mech_spectral_polarisation.f90 b/src/grid/grid_mech_spectral_polarisation.f90 index 7c9c85199..a52554e3c 100644 --- a/src/grid/grid_mech_spectral_polarisation.f90 +++ b/src/grid/grid_mech_spectral_polarisation.f90 @@ -237,7 +237,7 @@ subroutine grid_mechanical_spectral_polarisation_init F_tau_lastInc = 2.0_pReal*F_lastInc endif restartRead - homogenization_F0 = reshape(F_lastInc, [3,3,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent + homogenization_F0 = reshape(F_lastInc, [3,3,product(grid(1:2))*grid3]) ! set starting condition for homogenization_mechanical_response call utilities_updateCoords(reshape(F,shape(F_lastInc))) call utilities_constitutiveResponse(P,P_av,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 reshape(F,shape(F_lastInc)), & ! target F diff --git a/src/grid/spectral_utilities.f90 b/src/grid/spectral_utilities.f90 index 607928c1b..bd37e2654 100644 --- a/src/grid/spectral_utilities.f90 +++ b/src/grid/spectral_utilities.f90 @@ -815,11 +815,11 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& homogenization_F = reshape(F,[3,3,product(grid(1:2))*grid3]) ! set materialpoint target F to estimated field - call materialpoint_stressAndItsTangent(timeinc,[1,1],[1,product(grid(1:2))*grid3]) ! calculate P field + call homogenization_mechanical_response(timeinc,[1,1],[1,product(grid(1:2))*grid3]) ! calculate P field if (.not. terminallyIll) & - call materialpoint_stressAndItsTangent3(timeinc,[1,1],[1,product(grid(1:2))*grid3]) + call homogenization_thermal_response(timeinc,[1,1],[1,product(grid(1:2))*grid3]) if (.not. terminallyIll) & - call materialpoint_stressAndItsTangent2(timeinc,[1,1],[1,product(grid(1:2))*grid3]) + call homogenization_mechanical_response2(timeinc,[1,1],[1,product(grid(1:2))*grid3]) P = reshape(homogenization_P, [3,3,grid(1),grid(2),grid3]) P_av = sum(sum(sum(P,dim=5),dim=4),dim=3) * wgt diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 57b3dae4b..03a7c620d 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -101,8 +101,8 @@ module homogenization integer, intent(in) :: ce end subroutine damage_partition - module subroutine mechanical_homogenize(dt,ce) - real(pReal), intent(in) :: dt + module subroutine mechanical_homogenize(Delta_t,ce) + real(pReal), intent(in) :: Delta_t integer, intent(in) :: & ce !< cell end subroutine mechanical_homogenize @@ -178,9 +178,9 @@ module homogenization public :: & homogenization_init, & - materialpoint_stressAndItsTangent, & - materialpoint_stressAndItsTangent2, & - materialpoint_stressAndItsTangent3, & + homogenization_mechanical_response, & + homogenization_mechanical_response2, & + homogenization_thermal_response, & homogenization_mu_T, & homogenization_K_T, & homogenization_f_T, & @@ -231,15 +231,15 @@ end subroutine homogenization_init !-------------------------------------------------------------------------------------------------- !> @brief !-------------------------------------------------------------------------------------------------- -subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execElem) +subroutine homogenization_mechanical_response(Delta_t,FEsolving_execIP,FEsolving_execElem) - real(pReal), intent(in) :: dt !< time increment + real(pReal), intent(in) :: Delta_t !< time increment integer, dimension(2), intent(in) :: FEsolving_execElem, FEsolving_execIP integer :: & NiterationMPstate, & ip, & !< integration point number el, & !< element number - co, ce, ho, en, ph + co, ce, ho, en logical :: & converged logical, dimension(2) :: & @@ -268,13 +268,10 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE NiterationMPstate = NiterationMPstate + 1 call mechanical_partition(homogenization_F(1:3,1:3,ce),ce) - converged = .true. - do co = 1, homogenization_Nconstituents(ho) - converged = converged .and. phase_mechanical_constitutive(dt,co,ip,el) - enddo + converged = all([(phase_mechanical_constitutive(Delta_t,co,ip,el),co=1,homogenization_Nconstituents(ho))]) if (converged) then - doneAndHappy = mechanical_updateState(dt,homogenization_F(1:3,1:3,ce),ce) + doneAndHappy = mechanical_updateState(Delta_t,homogenization_F(1:3,1:3,ce),ce) converged = all(doneAndHappy) else doneAndHappy = [.true.,.false.] @@ -282,34 +279,30 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE enddo convergenceLooping if (.not. converged) then - if (.not. terminallyIll) print*, ' Integration point ', ip,' at element ', el, ' terminally ill' + if (.not. terminallyIll) print*, ' Cell ', ce, ' terminally ill' terminallyIll = .true. endif enddo enddo !$OMP END PARALLEL DO -end subroutine materialpoint_stressAndItsTangent +end subroutine homogenization_mechanical_response !-------------------------------------------------------------------------------------------------- !> @brief !-------------------------------------------------------------------------------------------------- -subroutine materialpoint_stressAndItsTangent3(dt,FEsolving_execIP,FEsolving_execElem) +subroutine homogenization_thermal_response(Delta_t,FEsolving_execIP,FEsolving_execElem) - real(pReal), intent(in) :: dt !< time increment + real(pReal), intent(in) :: Delta_t !< time increment integer, dimension(2), intent(in) :: FEsolving_execElem, FEsolving_execIP integer :: & - NiterationMPstate, & - ip, & !< integration point number - el, & !< element number - co, ce, ho, en, ph - logical :: & - converged - logical, dimension(2) :: & - doneAndHappy + ip, & !< integration point number + el, & !< element number + co, ce, ho - !$OMP PARALLEL DO PRIVATE(ho,ph,ce) + + !$OMP PARALLEL DO PRIVATE(ho,ce) do el = FEsolving_execElem(1),FEsolving_execElem(2) if (terminallyIll) continue do ip = FEsolving_execIP(1),FEsolving_execIP(2) @@ -317,36 +310,29 @@ subroutine materialpoint_stressAndItsTangent3(dt,FEsolving_execIP,FEsolving_exec ho = material_homogenizationID(ce) call thermal_partition(ce) do co = 1, homogenization_Nconstituents(ho) - ph = material_phaseID(co,ce) - if (.not. phase_thermal_constitutive(dt,ph,material_phaseMemberAt(co,ip,el))) then - if (.not. terminallyIll) & ! so first signals terminally ill... - print*, ' Integration point ', ip,' at element ', el, ' terminally ill' - terminallyIll = .true. ! ...and kills all others + if (.not. phase_thermal_constitutive(Delta_t,material_phaseID(co,ce),material_phaseEntry(co,ce))) then + if (.not. terminallyIll) print*, ' Cell ', ce, ' terminally ill' + terminallyIll = .true. endif enddo enddo enddo !$OMP END PARALLEL DO -end subroutine materialpoint_stressAndItsTangent3 +end subroutine homogenization_thermal_response !-------------------------------------------------------------------------------------------------- !> @brief !-------------------------------------------------------------------------------------------------- -subroutine materialpoint_stressAndItsTangent2(dt,FEsolving_execIP,FEsolving_execElem) +subroutine homogenization_mechanical_response2(Delta_t,FEsolving_execIP,FEsolving_execElem) - real(pReal), intent(in) :: dt !< time increment + real(pReal), intent(in) :: Delta_t !< time increment integer, dimension(2), intent(in) :: FEsolving_execElem, FEsolving_execIP integer :: & - NiterationMPstate, & - ip, & !< integration point number - el, & !< element number - co, ce, ho, en, ph - logical :: & - converged - logical, dimension(2) :: & - doneAndHappy + ip, & !< integration point number + el, & !< element number + co, ce, ho !$OMP PARALLEL DO PRIVATE(ho,ce) @@ -357,13 +343,13 @@ subroutine materialpoint_stressAndItsTangent2(dt,FEsolving_execIP,FEsolving_exec do co = 1, homogenization_Nconstituents(ho) call crystallite_orientations(co,ip,el) enddo - call mechanical_homogenize(dt,ce) + call mechanical_homogenize(Delta_t,ce) enddo IpLooping3 enddo elementLooping3 !$OMP END PARALLEL DO -end subroutine materialpoint_stressAndItsTangent2 +end subroutine homogenization_mechanical_response2 !-------------------------------------------------------------------------------------------------- diff --git a/src/homogenization_mechanical.f90 b/src/homogenization_mechanical.f90 index 7ea1f74e9..0138d2468 100644 --- a/src/homogenization_mechanical.f90 +++ b/src/homogenization_mechanical.f90 @@ -123,21 +123,21 @@ end subroutine mechanical_partition !-------------------------------------------------------------------------------------------------- !> @brief Average P and dPdF from the individual constituents. !-------------------------------------------------------------------------------------------------- -module subroutine mechanical_homogenize(dt,ce) +module subroutine mechanical_homogenize(Delta_t,ce) - real(pReal), intent(in) :: dt + real(pReal), intent(in) :: Delta_t integer, intent(in) :: ce integer :: co homogenization_P(1:3,1:3,ce) = phase_P(1,ce) - homogenization_dPdF(1:3,1:3,1:3,1:3,ce) = phase_mechanical_dPdF(dt,1,ce) + homogenization_dPdF(1:3,1:3,1:3,1:3,ce) = phase_mechanical_dPdF(Delta_t,1,ce) do co = 2, homogenization_Nconstituents(material_homogenizationID(ce)) homogenization_P(1:3,1:3,ce) = homogenization_P(1:3,1:3,ce) & + phase_P(co,ce) homogenization_dPdF(1:3,1:3,1:3,1:3,ce) = homogenization_dPdF(1:3,1:3,1:3,1:3,ce) & - + phase_mechanical_dPdF(dt,co,ce) + + phase_mechanical_dPdF(Delta_t,co,ce) enddo homogenization_P(1:3,1:3,ce) = homogenization_P(1:3,1:3,ce) & diff --git a/src/mesh/FEM_utilities.f90 b/src/mesh/FEM_utilities.f90 index 402c6ddf7..5bb6d5f44 100644 --- a/src/mesh/FEM_utilities.f90 +++ b/src/mesh/FEM_utilities.f90 @@ -162,9 +162,9 @@ subroutine utilities_constitutiveResponse(timeinc,P_av,forwardData) print'(/,a)', ' ... evaluating constitutive response ......................................' - call materialpoint_stressAndItsTangent(timeinc,[1,mesh_maxNips],[1,mesh_NcpElems]) ! calculate P field + call homogenization_mechanical_response(timeinc,[1,mesh_maxNips],[1,mesh_NcpElems]) ! calculate P field if (.not. terminallyIll) & - call materialpoint_stressAndItsTangent2(timeinc,[1,mesh_maxNips],[1,mesh_NcpElems]) + call homogenization_mechanical_response2(timeinc,[1,mesh_maxNips],[1,mesh_NcpElems]) cutBack = .false. P_av = sum(homogenization_P,dim=3) * wgt diff --git a/src/phase.f90 b/src/phase.f90 index 9324225a2..b551c4170 100644 --- a/src/phase.f90 +++ b/src/phase.f90 @@ -113,8 +113,8 @@ module phase end subroutine damage_restore - module function phase_mechanical_dPdF(dt,co,ce) result(dPdF) - real(pReal), intent(in) :: dt + module function phase_mechanical_dPdF(Delta_t,co,ce) result(dPdF) + real(pReal), intent(in) :: Delta_t integer, intent(in) :: & co, & !< counter in constituent loop ce @@ -221,21 +221,21 @@ module phase end function phase_thermal_constitutive - module function integrateDamageState(dt,co,ce) result(broken) - real(pReal), intent(in) :: dt + module function integrateDamageState(Delta_t,co,ce) result(broken) + real(pReal), intent(in) :: Delta_t integer, intent(in) :: & ce, & co logical :: broken end function integrateDamageState - module function phase_mechanical_constitutive(dt,co,ip,el) result(converged_) - real(pReal), intent(in) :: dt + module function phase_mechanical_constitutive(Delta_t,co,ip,el) result(converged_) + real(pReal), intent(in) :: Delta_t integer, intent(in) :: co, ip, el logical :: converged_ end function phase_mechanical_constitutive - !ToDo: Try to merge the all stiffness functions + !ToDo: Merge all the stiffness functions module function phase_homogenizedC(ph,en) result(C) integer, intent(in) :: ph, en real(pReal), dimension(6,6) :: C diff --git a/src/phase_damage.f90 b/src/phase_damage.f90 index 700bb2142..89d3b1653 100644 --- a/src/phase_damage.f90 +++ b/src/phase_damage.f90 @@ -197,9 +197,9 @@ end function phase_f_phi !> @brief integrate stress, state with adaptive 1st order explicit Euler method !> using Fixed Point Iteration to adapt the stepsize !-------------------------------------------------------------------------------------------------- -module function integrateDamageState(dt,co,ce) result(broken) +module function integrateDamageState(Delta_t,co,ce) result(broken) - real(pReal), intent(in) :: dt + real(pReal), intent(in) :: Delta_t integer, intent(in) :: & ce, & co @@ -230,10 +230,10 @@ module function integrateDamageState(dt,co,ce) result(broken) broken = phase_damage_collectDotState(ph,me) if(broken) return - size_so = damageState(ph)%sizeDotState - damageState(ph)%state(1:size_so,me) = damageState(ph)%state0 (1:size_so,me) & - + damageState(ph)%dotState(1:size_so,me) * dt - source_dotState(1:size_so,2) = 0.0_pReal + size_so = damageState(ph)%sizeDotState + damageState(ph)%state(1:size_so,me) = damageState(ph)%state0 (1:size_so,me) & + + damageState(ph)%dotState(1:size_so,me) * Delta_t + source_dotState(1:size_so,2) = 0.0_pReal iteration: do NiterationState = 1, num%nState @@ -249,7 +249,7 @@ module function integrateDamageState(dt,co,ce) result(broken) + source_dotState(1:size_so,1)* (1.0_pReal - zeta) r(1:size_so) = damageState(ph)%state (1:size_so,me) & - damageState(ph)%State0 (1:size_so,me) & - - damageState(ph)%dotState(1:size_so,me) * dt + - damageState(ph)%dotState(1:size_so,me) * Delta_t damageState(ph)%state(1:size_so,me) = damageState(ph)%state(1:size_so,me) - r(1:size_so) converged_ = converged_ .and. converged(r(1:size_so), & damageState(ph)%state(1:size_so,me), & diff --git a/src/phase_damage_anisobrittle.f90 b/src/phase_damage_anisobrittle.f90 index e71c6d3cc..cc3628323 100644 --- a/src/phase_damage_anisobrittle.f90 +++ b/src/phase_damage_anisobrittle.f90 @@ -40,7 +40,7 @@ module function anisobrittle_init() result(mySources) phase, & sources, & src - integer :: Nmembers,p + integer :: Nmembers,ph integer, dimension(:), allocatable :: N_cl character(len=pStringLen) :: extmsg = '' @@ -56,12 +56,12 @@ module function anisobrittle_init() result(mySources) allocate(param(phases%length)) - do p = 1, phases%length - if(mySources(p)) then - phase => phases%get(p) + do ph = 1, phases%length + if(mySources(ph)) then + phase => phases%get(ph) sources => phase%get('damage') - associate(prm => param(p)) + associate(prm => param(ph)) src => sources%get(1) N_cl = src%get_as1dInt('N_cl',defaultVal=emptyIntArray) @@ -92,10 +92,10 @@ module function anisobrittle_init() result(mySources) if (any(prm%g_crit < 0.0_pReal)) extmsg = trim(extmsg)//' g_crit' if (any(prm%s_crit < 0.0_pReal)) extmsg = trim(extmsg)//' s_crit' - Nmembers = count(material_phaseID==p) - call phase_allocateState(damageState(p),Nmembers,1,1,0) - damageState(p)%atol = src%get_asFloat('atol_phi',defaultVal=1.0e-9_pReal) - if(any(damageState(p)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' atol_phi' + Nmembers = count(material_phaseID==ph) + call phase_allocateState(damageState(ph),Nmembers,1,1,0) + damageState(ph)%atol = src%get_asFloat('atol_phi',defaultVal=1.0e-9_pReal) + if(any(damageState(ph)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' atol_phi' end associate diff --git a/src/phase_mechanical.f90 b/src/phase_mechanical.f90 index 65dd2ceac..e03998989 100644 --- a/src/phase_mechanical.f90 +++ b/src/phase_mechanical.f90 @@ -44,9 +44,9 @@ submodule(phase) mechanical interface - module subroutine eigendeformation_init(phases) + module subroutine eigen_init(phases) class(tNode), pointer :: phases - end subroutine eigendeformation_init + end subroutine eigen_init module subroutine elastic_init(phases) class(tNode), pointer :: phases @@ -302,7 +302,7 @@ module subroutine mechanical_init(materials,phases) end select - call eigendeformation_init(phases) + call eigen_init(phases) end subroutine mechanical_init @@ -976,9 +976,9 @@ end subroutine mechanical_forward !-------------------------------------------------------------------------------------------------- !> @brief calculate stress (P) !-------------------------------------------------------------------------------------------------- -module function phase_mechanical_constitutive(dt,co,ip,el) result(converged_) +module function phase_mechanical_constitutive(Delta_t,co,ip,el) result(converged_) - real(pReal), intent(in) :: dt + real(pReal), intent(in) :: Delta_t integer, intent(in) :: & co, & ip, & @@ -1054,12 +1054,12 @@ module function phase_mechanical_constitutive(dt,co,ip,el) result(converged_) if (todo) then subF = subF0 & + subStep * (phase_mechanical_F(ph)%data(1:3,1:3,en) - phase_mechanical_F0(ph)%data(1:3,1:3,en)) - converged_ = .not. integrateState(subF0,subF,subFp0,subFi0,subState0(1:sizeDotState),subStep * dt,co,ip,el) + converged_ = .not. integrateState(subF0,subF,subFp0,subFi0,subState0(1:sizeDotState),subStep * Delta_t,co,ip,el) endif enddo cutbackLooping - converged_ = converged_ .and. .not. integrateDamageState(dt,co,(el-1)*discretization_nIPs + ip) + converged_ = converged_ .and. .not. integrateDamageState(Delta_t,co,(el-1)*discretization_nIPs + ip) end function phase_mechanical_constitutive @@ -1094,12 +1094,13 @@ module subroutine mechanical_restore(ce,includeL) end subroutine mechanical_restore + !-------------------------------------------------------------------------------------------------- !> @brief Calculate tangent (dPdF). !-------------------------------------------------------------------------------------------------- -module function phase_mechanical_dPdF(dt,co,ce) result(dPdF) +module function phase_mechanical_dPdF(Delta_t,co,ce) result(dPdF) - real(pReal), intent(in) :: dt + real(pReal), intent(in) :: Delta_t integer, intent(in) :: & co, & !< counter in constituent loop ce @@ -1149,11 +1150,11 @@ module function phase_mechanical_dPdF(dt,co,ce) result(dPdF) lhs_3333 = 0.0_pReal; rhs_3333 = 0.0_pReal do o=1,3; do p=1,3 lhs_3333(1:3,1:3,o,p) = lhs_3333(1:3,1:3,o,p) & - + matmul(invSubFi0,dLidFi(1:3,1:3,o,p)) * dt + + matmul(invSubFi0,dLidFi(1:3,1:3,o,p)) * Delta_t lhs_3333(1:3,o,1:3,p) = lhs_3333(1:3,o,1:3,p) & + invFi*invFi(p,o) rhs_3333(1:3,1:3,o,p) = rhs_3333(1:3,1:3,o,p) & - - matmul(invSubFi0,dLidS(1:3,1:3,o,p)) * dt + - matmul(invSubFi0,dLidS(1:3,1:3,o,p)) * Delta_t enddo; enddo call math_invert(temp_99,error,math_3333to99(lhs_3333)) if (error) then @@ -1182,7 +1183,7 @@ module function phase_mechanical_dPdF(dt,co,ce) result(dPdF) temp_3333(1:3,1:3,p,o) = matmul(matmul(temp_33_2,dLpdS(1:3,1:3,p,o)), invFi) & + matmul(temp_33_3,dLidS(1:3,1:3,p,o)) enddo; enddo - lhs_3333 = math_mul3333xx3333(dSdFe,temp_3333) * dt & + lhs_3333 = math_mul3333xx3333(dSdFe,temp_3333) * Delta_t & + math_mul3333xx3333(dSdFi,dFidS) call math_invert(temp_99,error,math_eye(9)+math_3333to99(lhs_3333)) @@ -1198,7 +1199,7 @@ module function phase_mechanical_dPdF(dt,co,ce) result(dPdF) ! calculate dFpinvdF temp_3333 = math_mul3333xx3333(dLpdS,dSdF) do o=1,3; do p=1,3 - dFpinvdF(1:3,1:3,p,o) = - matmul(invSubFp0, matmul(temp_3333(1:3,1:3,p,o),invFi)) * dt + dFpinvdF(1:3,1:3,p,o) = - matmul(invSubFp0, matmul(temp_3333(1:3,1:3,p,o),invFi)) * Delta_t enddo; enddo !-------------------------------------------------------------------------------------------------- diff --git a/src/phase_mechanical_eigen.f90 b/src/phase_mechanical_eigen.f90 index 019838689..14d86661b 100644 --- a/src/phase_mechanical_eigen.f90 +++ b/src/phase_mechanical_eigen.f90 @@ -32,7 +32,7 @@ submodule(phase:mechanical) eigen contains -module subroutine eigendeformation_init(phases) +module subroutine eigen_init(phases) class(tNode), pointer :: & phases @@ -68,7 +68,7 @@ module subroutine eigendeformation_init(phases) where(damage_anisobrittle_init()) model_damage = KINEMATICS_cleavage_opening_ID -end subroutine eigendeformation_init +end subroutine eigen_init !--------------------------------------------------------------------------------------------------