diff --git a/src/CPFEM.f90 b/src/CPFEM.f90 index fa774de66..c49ecbcb6 100644 --- a/src/CPFEM.f90 +++ b/src/CPFEM.f90 @@ -107,9 +107,9 @@ subroutine CPFEM_init print'(/,a)', ' <<<+- CPFEM init -+>>>'; flush(IO_STDOUT) - allocate(CPFEM_cs( 6,discretization_nIP,discretization_nElem), source= 0.0_pReal) - allocate(CPFEM_dcsdE( 6,6,discretization_nIP,discretization_nElem), source= 0.0_pReal) - allocate(CPFEM_dcsdE_knownGood(6,6,discretization_nIP,discretization_nElem), source= 0.0_pReal) + allocate(CPFEM_cs( 6,discretization_nIPs,discretization_Nelems), source= 0.0_pReal) + allocate(CPFEM_dcsdE( 6,6,discretization_nIPs,discretization_Nelems), source= 0.0_pReal) + allocate(CPFEM_dcsdE_knownGood(6,6,discretization_nIPs,discretization_Nelems), source= 0.0_pReal) !------------------------------------------------------------------------------ ! read debug options @@ -184,8 +184,8 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature_inp, dt, elFE, ip, cauchyS temperature(material_homogenizationAt(elCP))%p(thermalMapping(material_homogenizationAt(elCP))%p(ip,elCP)) = & temperature_inp end select chosenThermal1 - materialpoint_F0(1:3,1:3,ip,elCP) = ffn - materialpoint_F(1:3,1:3,ip,elCP) = ffn1 + homogenization_F0(1:3,1:3,ip,elCP) = ffn + homogenization_F(1:3,1:3,ip,elCP) = ffn1 if (iand(mode, CPFEM_CALCRESULTS) /= 0_pInt) then @@ -212,17 +212,17 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature_inp, dt, elFE, ip, cauchyS else terminalIllness ! translate from P to sigma - Kirchhoff = matmul(materialpoint_P(1:3,1:3,ip,elCP), transpose(materialpoint_F(1:3,1:3,ip,elCP))) - J_inverse = 1.0_pReal / math_det33(materialpoint_F(1:3,1:3,ip,elCP)) + Kirchhoff = matmul(homogenization_P(1:3,1:3,ip,elCP), transpose(homogenization_F(1:3,1:3,ip,elCP))) + J_inverse = 1.0_pReal / math_det33(homogenization_F(1:3,1:3,ip,elCP)) CPFEM_cs(1:6,ip,elCP) = math_sym33to6(J_inverse * Kirchhoff,weighted=.false.) ! translate from dP/dF to dCS/dE H = 0.0_pReal do i=1,3; do j=1,3; do k=1,3; do l=1,3; do m=1,3; do n=1,3 H(i,j,k,l) = H(i,j,k,l) & - + materialpoint_F(j,m,ip,elCP) * materialpoint_F(l,n,ip,elCP) & - * materialpoint_dPdF(i,m,k,n,ip,elCP) & - - math_delta(j,l) * materialpoint_F(i,m,ip,elCP) * materialpoint_P(k,m,ip,elCP) & + + homogenization_F(j,m,ip,elCP) * homogenization_F(l,n,ip,elCP) & + * homogenization_dPdF(i,m,k,n,ip,elCP) & + - math_delta(j,l) * homogenization_F(i,m,ip,elCP) * homogenization_P(k,m,ip,elCP) & + 0.5_pReal * ( Kirchhoff(j,l)*math_delta(i,k) + Kirchhoff(i,k)*math_delta(j,l) & + Kirchhoff(j,k)*math_delta(i,l) + Kirchhoff(i,l)*math_delta(j,k)) enddo; enddo; enddo; enddo; enddo; enddo diff --git a/src/YAML_types.f90 b/src/YAML_types.f90 index f11faf54b..b71261d9c 100644 --- a/src/YAML_types.f90 +++ b/src/YAML_types.f90 @@ -583,9 +583,9 @@ function tNode_get_byIndex_asStrings(self,i) result(nodeAsStrings) end function tNode_get_byIndex_asStrings -!------------------------------------------------------------------------------------------------------- +!-------------------------------------------------------------------------------------------------- !> @brief Returns the key in a dictionary as a string -!------------------------------------------------------------------------------------------------------- +!-------------------------------------------------------------------------------------------------- function tNode_getKey_byIndex(self,i) result(key) class(tNode), intent(in), target :: self @@ -1169,17 +1169,21 @@ function tList_asStrings(self) type(tScalar), pointer :: scalar len_max = 0 - allocate(character(len=pStringLen) :: tList_asStrings(self%length)) + item => self%first + do i = 1, self%length + scalar => item%node%asScalar() + len_max = max(len_max, len_trim(scalar%asString())) + item => item%next + enddo + + allocate(character(len=len_max) :: tList_asStrings(self%length)) item => self%first do i = 1, self%length scalar => item%node%asScalar() tList_asStrings(i) = scalar%asString() - len_max = max(len_max, len_trim(tList_asStrings(i))) item => item%next enddo - !ToDo: trim to len_max - end function tList_asStrings diff --git a/src/constitutive.f90 b/src/constitutive.f90 index 64d3f2d31..358937e4b 100644 --- a/src/constitutive.f90 +++ b/src/constitutive.f90 @@ -127,7 +127,7 @@ module constitutive instance,of,ip,el) real(pReal), dimension(3,3), intent(in) :: & Mp !< MandelStress - real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: & + real(pReal), dimension(3,3,homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems), intent(in) :: & F, & !< deformation gradient Fp !< plastic deformation gradient real(pReal), intent(in) :: & @@ -218,7 +218,7 @@ module constitutive instance, & i, & e - type(rotation), dimension(1,discretization_nIP,discretization_nElem), intent(in) :: & + type(rotation), dimension(1,discretization_nIPs,discretization_Nelems), intent(in) :: & orientation !< crystal orientation end subroutine plastic_nonlocal_updateCompatibility @@ -753,7 +753,7 @@ function constitutive_collectDotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el of real(pReal), intent(in) :: & subdt !< timestep - real(pReal), intent(in), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: & + real(pReal), intent(in), dimension(3,3,homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems) :: & FArray, & !< elastic deformation gradient FpArray !< plastic deformation gradient real(pReal), intent(in), dimension(3,3) :: & @@ -905,12 +905,12 @@ end function constitutive_deltaState !> @brief Allocate the components of the state structure for a given phase !-------------------------------------------------------------------------------------------------- subroutine constitutive_allocateState(state, & - NipcMyPhase,sizeState,sizeDotState,sizeDeltaState) + Nconstituents,sizeState,sizeDotState,sizeDeltaState) class(tState), intent(out) :: & state integer, intent(in) :: & - NipcMyPhase, & + Nconstituents, & sizeState, & sizeDotState, & sizeDeltaState @@ -921,14 +921,14 @@ subroutine constitutive_allocateState(state, & state%offsetDeltaState = sizeState-sizeDeltaState ! deltaState occupies latter part of state by definition allocate(state%atol (sizeState), source=0.0_pReal) - allocate(state%state0 (sizeState,NipcMyPhase), source=0.0_pReal) - allocate(state%partitionedState0(sizeState,NipcMyPhase), source=0.0_pReal) - allocate(state%subState0 (sizeState,NipcMyPhase), source=0.0_pReal) - allocate(state%state (sizeState,NipcMyPhase), source=0.0_pReal) + allocate(state%state0 (sizeState,Nconstituents), source=0.0_pReal) + allocate(state%partitionedState0(sizeState,Nconstituents), source=0.0_pReal) + allocate(state%subState0 (sizeState,Nconstituents), source=0.0_pReal) + allocate(state%state (sizeState,Nconstituents), source=0.0_pReal) - allocate(state%dotState (sizeDotState,NipcMyPhase), source=0.0_pReal) + allocate(state%dotState (sizeDotState,Nconstituents), source=0.0_pReal) - allocate(state%deltaState(sizeDeltaState,NipcMyPhase), source=0.0_pReal) + allocate(state%deltaState(sizeDeltaState,Nconstituents), source=0.0_pReal) end subroutine constitutive_allocateState diff --git a/src/constitutive_damage.f90 b/src/constitutive_damage.f90 index 5572192f7..56198987e 100644 --- a/src/constitutive_damage.f90 +++ b/src/constitutive_damage.f90 @@ -184,7 +184,7 @@ module subroutine constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi phiDot = 0.0_pReal dPhiDot_dPhi = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) phase = material_phaseAt(grain,el) constituent = material_phasememberAt(grain,ip,el) do source = 1, phase_Nsources(phase) diff --git a/src/constitutive_plastic_disloTungsten.f90 b/src/constitutive_plastic_disloTungsten.f90 index 54c01b912..aef3bef89 100644 --- a/src/constitutive_plastic_disloTungsten.f90 +++ b/src/constitutive_plastic_disloTungsten.f90 @@ -78,9 +78,9 @@ module function plastic_disloTungsten_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState, & startIndex, endIndex integer, dimension(:), allocatable :: & @@ -99,17 +99,17 @@ module function plastic_disloTungsten_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_dislotungsten init -+>>>' myPlasticity = plastic_active('disloTungsten') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return print*, 'Cereceda et al., International Journal of Plasticity 78:242–256, 2016' print*, 'https://dx.doi.org/10.1016/j.ijplas.2015.09.002' - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) - allocate(dependentState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) + allocate(dependentState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -221,18 +221,18 @@ module function plastic_disloTungsten_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,0) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization startIndex = 1 endIndex = prm%sum_N_sl stt%rho_mob => plasticState(p)%state(startIndex:endIndex,:) - stt%rho_mob = spread(rho_mob_0,2,NipcMyPhase) + stt%rho_mob = spread(rho_mob_0,2,Nconstituents) dot%rho_mob => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_rho',defaultVal=1.0_pReal) if (any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_rho' @@ -240,7 +240,7 @@ module function plastic_disloTungsten_init() result(myPlasticity) startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_sl stt%rho_dip => plasticState(p)%state(startIndex:endIndex,:) - stt%rho_dip = spread(rho_dip_0,2,NipcMyPhase) + stt%rho_dip = spread(rho_dip_0,2,Nconstituents) dot%rho_dip => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_rho',defaultVal=1.0_pReal) @@ -252,8 +252,8 @@ module function plastic_disloTungsten_init() result(myPlasticity) ! global alias plasticState(p)%slipRate => plasticState(p)%dotState(startIndex:endIndex,:) - allocate(dst%Lambda_sl(prm%sum_N_sl,NipcMyPhase), source=0.0_pReal) - allocate(dst%threshold_stress(prm%sum_N_sl,NipcMyPhase), source=0.0_pReal) + allocate(dst%Lambda_sl(prm%sum_N_sl,Nconstituents), source=0.0_pReal) + allocate(dst%threshold_stress(prm%sum_N_sl,Nconstituents), source=0.0_pReal) plasticState(p)%state0 = plasticState(p)%state ! ToDo: this could be done centrally diff --git a/src/constitutive_plastic_dislotwin.f90 b/src/constitutive_plastic_dislotwin.f90 index 62dcdd83e..30527e301 100644 --- a/src/constitutive_plastic_dislotwin.f90 +++ b/src/constitutive_plastic_dislotwin.f90 @@ -126,9 +126,9 @@ module function plastic_dislotwin_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState, & startIndex, endIndex integer, dimension(:), allocatable :: & @@ -146,9 +146,9 @@ module function plastic_dislotwin_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_dislotwin init -+>>>' myPlasticity = plastic_active('dislotwin') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return print*, 'Ma and Roters, Acta Materialia 52(12):3603–3612, 2004' print*, 'https://doi.org/10.1016/j.actamat.2004.04.012'//IO_EOL @@ -159,10 +159,10 @@ module function plastic_dislotwin_init() result(myPlasticity) print*, 'Wong et al., Acta Materialia 118:140–151, 2016' print*, 'https://doi.org/10.1016/j.actamat.2016.07.032' - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) - allocate(dependentState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) + allocate(dependentState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -407,21 +407,21 @@ module function plastic_dislotwin_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl & + size(['f_tw']) * prm%sum_N_tw & + size(['f_tr']) * prm%sum_N_tr sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,0) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! locally defined state aliases and initialization of state0 and atol startIndex = 1 endIndex = prm%sum_N_sl stt%rho_mob=>plasticState(p)%state(startIndex:endIndex,:) - stt%rho_mob= spread(rho_mob_0,2,NipcMyPhase) + stt%rho_mob= spread(rho_mob_0,2,Nconstituents) dot%rho_mob=>plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_rho',defaultVal=1.0_pReal) if (any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_rho' @@ -429,7 +429,7 @@ module function plastic_dislotwin_init() result(myPlasticity) startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_sl stt%rho_dip=>plasticState(p)%state(startIndex:endIndex,:) - stt%rho_dip= spread(rho_dip_0,2,NipcMyPhase) + stt%rho_dip= spread(rho_dip_0,2,Nconstituents) dot%rho_dip=>plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_rho',defaultVal=1.0_pReal) @@ -455,18 +455,18 @@ module function plastic_dislotwin_init() result(myPlasticity) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('f_trans',defaultVal=1.0e-6_pReal) if (any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' f_trans' - allocate(dst%Lambda_sl (prm%sum_N_sl,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_pass (prm%sum_N_sl,NipcMyPhase),source=0.0_pReal) + allocate(dst%Lambda_sl (prm%sum_N_sl,Nconstituents),source=0.0_pReal) + allocate(dst%tau_pass (prm%sum_N_sl,Nconstituents),source=0.0_pReal) - allocate(dst%Lambda_tw (prm%sum_N_tw,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_hat_tw (prm%sum_N_tw,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_r_tw (prm%sum_N_tw,NipcMyPhase),source=0.0_pReal) - allocate(dst%V_tw (prm%sum_N_tw,NipcMyPhase),source=0.0_pReal) + allocate(dst%Lambda_tw (prm%sum_N_tw,Nconstituents),source=0.0_pReal) + allocate(dst%tau_hat_tw (prm%sum_N_tw,Nconstituents),source=0.0_pReal) + allocate(dst%tau_r_tw (prm%sum_N_tw,Nconstituents),source=0.0_pReal) + allocate(dst%V_tw (prm%sum_N_tw,Nconstituents),source=0.0_pReal) - allocate(dst%Lambda_tr (prm%sum_N_tr,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_hat_tr (prm%sum_N_tr,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_r_tr (prm%sum_N_tr,NipcMyPhase),source=0.0_pReal) - allocate(dst%V_tr (prm%sum_N_tr,NipcMyPhase),source=0.0_pReal) + allocate(dst%Lambda_tr (prm%sum_N_tr,Nconstituents),source=0.0_pReal) + allocate(dst%tau_hat_tr (prm%sum_N_tr,Nconstituents),source=0.0_pReal) + allocate(dst%tau_r_tr (prm%sum_N_tr,Nconstituents),source=0.0_pReal) + allocate(dst%V_tr (prm%sum_N_tr,Nconstituents),source=0.0_pReal) plasticState(p)%state0 = plasticState(p)%state ! ToDo: this could be done centrally diff --git a/src/constitutive_plastic_isotropic.f90 b/src/constitutive_plastic_isotropic.f90 index 2c9028671..db418663f 100644 --- a/src/constitutive_plastic_isotropic.f90 +++ b/src/constitutive_plastic_isotropic.f90 @@ -53,10 +53,10 @@ module function plastic_isotropic_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, & i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState real(pReal) :: & xi_0 !< initial critical stress @@ -70,16 +70,16 @@ module function plastic_isotropic_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_isotropic init -+>>>' myPlasticity = plastic_active('isotropic') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return print*, 'Maiti and Eisenlohr, Scripta Materialia 145:37–40, 2018' print*, 'https://doi.org/10.1016/j.scriptamat.2017.09.047' - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -130,11 +130,11 @@ module function plastic_isotropic_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['xi ','gamma']) sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,0) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization diff --git a/src/constitutive_plastic_kinehardening.f90 b/src/constitutive_plastic_kinehardening.f90 index 777691242..3faf7dc41 100644 --- a/src/constitutive_plastic_kinehardening.f90 +++ b/src/constitutive_plastic_kinehardening.f90 @@ -62,9 +62,9 @@ module function plastic_kinehardening_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, o, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDeltaState, sizeDotState, & startIndex, endIndex integer, dimension(:), allocatable :: & @@ -82,14 +82,14 @@ module function plastic_kinehardening_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_kinehardening init -+>>>' myPlasticity = plastic_active('kinehardening') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) - allocate(deltaState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) + allocate(deltaState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -174,19 +174,19 @@ module function plastic_kinehardening_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['crss ','crss_back', 'accshear ']) * prm%sum_N_sl!ToDo: adjust names, ask Philip sizeDeltaState = size(['sense ', 'chi0 ', 'gamma0' ]) * prm%sum_N_sl !ToDo: adjust names sizeState = sizeDotState + sizeDeltaState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,sizeDeltaState) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,sizeDeltaState) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization startIndex = 1 endIndex = prm%sum_N_sl stt%crss => plasticState(p)%state (startIndex:endIndex,:) - stt%crss = spread(xi_0, 2, NipcMyPhase) + stt%crss = spread(xi_0, 2, Nconstituents) dot%crss => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_xi',defaultVal=1.0_pReal) if(any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_xi' diff --git a/src/constitutive_plastic_none.f90 b/src/constitutive_plastic_none.f90 index ab5f32d3f..923ae354a 100644 --- a/src/constitutive_plastic_none.f90 +++ b/src/constitutive_plastic_none.f90 @@ -16,9 +16,9 @@ module function plastic_none_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, & - NipcMyPhase + Nconstituents class(tNode), pointer :: & phases, & phase, & @@ -34,15 +34,15 @@ module function plastic_none_init() result(myPlasticity) if(pl%get_asString('type') == 'none') myPlasticity(p) = .true. enddo - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return do p = 1, phases%length phase => phases%get(p) if(.not. myPlasticity(p)) cycle - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP - call constitutive_allocateState(plasticState(p),NipcMyPhase,0,0,0) + Nconstituents = count(material_phaseAt == p) * discretization_nIPs + call constitutive_allocateState(plasticState(p),Nconstituents,0,0,0) enddo end function plastic_none_init diff --git a/src/constitutive_plastic_nonlocal.f90 b/src/constitutive_plastic_nonlocal.f90 index 8238f17c9..472415ba9 100644 --- a/src/constitutive_plastic_nonlocal.f90 +++ b/src/constitutive_plastic_nonlocal.f90 @@ -153,7 +153,7 @@ submodule(constitutive:constitutive_plastic) plastic_nonlocal state, & state0 - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) type(tNonlocalMicrostructure), dimension(:), allocatable :: microstructure @@ -168,9 +168,9 @@ module function plastic_nonlocal_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState, sizeDependentState, sizeDeltaState, & s1, s2, & s, t, l @@ -188,9 +188,9 @@ module function plastic_nonlocal_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_nonlocal init -+>>>' myPlasticity = plastic_active('nonlocal') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) then + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) then call geometry_plastic_nonlocal_disable return endif @@ -201,12 +201,12 @@ module function plastic_nonlocal_init() result(myPlasticity) print*, 'Kords, Dissertation RWTH Aachen, 2014' print*, 'http://publications.rwth-aachen.de/record/229993' - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(state0(Ninstance)) - allocate(dotState(Ninstance)) - allocate(deltaState(Ninstance)) - allocate(microstructure(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(state0(Ninstances)) + allocate(dotState(Ninstances)) + allocate(deltaState(Ninstances)) + allocate(microstructure(Ninstances)) phases => config_material%get('phase') i = 0 @@ -391,7 +391,7 @@ module function plastic_nonlocal_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP + Nconstituents = count(material_phaseAt==p) * discretization_nIPs sizeDotState = size([ 'rhoSglEdgePosMobile ','rhoSglEdgeNegMobile ', & 'rhoSglScrewPosMobile ','rhoSglScrewNegMobile ', & 'rhoSglEdgePosImmobile ','rhoSglEdgeNegImmobile ', & @@ -405,7 +405,7 @@ module function plastic_nonlocal_init() result(myPlasticity) 'maxDipoleHeightEdge ','maxDipoleHeightScrew' ]) * prm%sum_N_sl !< other dependent state variables that are not updated by microstructure sizeDeltaState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,sizeDeltaState) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,sizeDeltaState) plasticState(p)%nonlocal = pl%get_asBool('nonlocal') if(plasticState(p)%nonlocal .and. .not. allocated(IPneighborhood)) & @@ -476,26 +476,26 @@ module function plastic_nonlocal_init() result(myPlasticity) dot%rho_dip_scr => plasticState(p)%dotState (9*prm%sum_N_sl+1:10*prm%sum_N_sl,:) del%rho_dip_scr => plasticState(p)%deltaState (9*prm%sum_N_sl+1:10*prm%sum_N_sl,:) - stt%gamma => plasticState(p)%state (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:NipcMyPhase) - dot%gamma => plasticState(p)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:NipcMyPhase) - del%gamma => plasticState(p)%deltaState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:NipcMyPhase) + stt%gamma => plasticState(p)%state (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nconstituents) + dot%gamma => plasticState(p)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nconstituents) + del%gamma => plasticState(p)%deltaState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nconstituents) plasticState(p)%atol(10*prm%sum_N_sl+1:11*prm%sum_N_sl ) = pl%get_asFloat('atol_gamma', defaultVal = 1.0e-2_pReal) if(any(plasticState(p)%atol(10*prm%sum_N_sl+1:11*prm%sum_N_sl) < 0.0_pReal)) & extmsg = trim(extmsg)//' atol_gamma' - plasticState(p)%slipRate => plasticState(p)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:NipcMyPhase) + plasticState(p)%slipRate => plasticState(p)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nconstituents) - stt%rho_forest => plasticState(p)%state (11*prm%sum_N_sl + 1:12*prm%sum_N_sl,1:NipcMyPhase) - stt%v => plasticState(p)%state (12*prm%sum_N_sl + 1:16*prm%sum_N_sl,1:NipcMyPhase) - stt%v_edg_pos => plasticState(p)%state (12*prm%sum_N_sl + 1:13*prm%sum_N_sl,1:NipcMyPhase) - stt%v_edg_neg => plasticState(p)%state (13*prm%sum_N_sl + 1:14*prm%sum_N_sl,1:NipcMyPhase) - stt%v_scr_pos => plasticState(p)%state (14*prm%sum_N_sl + 1:15*prm%sum_N_sl,1:NipcMyPhase) - stt%v_scr_neg => plasticState(p)%state (15*prm%sum_N_sl + 1:16*prm%sum_N_sl,1:NipcMyPhase) + stt%rho_forest => plasticState(p)%state (11*prm%sum_N_sl + 1:12*prm%sum_N_sl,1:Nconstituents) + stt%v => plasticState(p)%state (12*prm%sum_N_sl + 1:16*prm%sum_N_sl,1:Nconstituents) + stt%v_edg_pos => plasticState(p)%state (12*prm%sum_N_sl + 1:13*prm%sum_N_sl,1:Nconstituents) + stt%v_edg_neg => plasticState(p)%state (13*prm%sum_N_sl + 1:14*prm%sum_N_sl,1:Nconstituents) + stt%v_scr_pos => plasticState(p)%state (14*prm%sum_N_sl + 1:15*prm%sum_N_sl,1:Nconstituents) + stt%v_scr_neg => plasticState(p)%state (15*prm%sum_N_sl + 1:16*prm%sum_N_sl,1:Nconstituents) - allocate(dst%tau_pass(prm%sum_N_sl,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_back(prm%sum_N_sl,NipcMyPhase),source=0.0_pReal) + allocate(dst%tau_pass(prm%sum_N_sl,Nconstituents),source=0.0_pReal) + allocate(dst%tau_back(prm%sum_N_sl,Nconstituents),source=0.0_pReal) end associate - if (NipcMyPhase > 0) call stateInit(ini,p,NipcMyPhase,i) + if (Nconstituents > 0) call stateInit(ini,p,Nconstituents,i) plasticState(p)%state0 = plasticState(p)%state !-------------------------------------------------------------------------------------------------- @@ -505,12 +505,12 @@ module function plastic_nonlocal_init() result(myPlasticity) enddo allocate(compatibility(2,maxval(param%sum_N_sl),maxval(param%sum_N_sl),nIPneighbors,& - discretization_nIP,discretization_nElem), source=0.0_pReal) + discretization_nIPs,discretization_Nelems), source=0.0_pReal) ! BEGIN DEPRECATED---------------------------------------------------------------------------------- - allocate(iRhoU(maxval(param%sum_N_sl),4,Ninstance), source=0) - allocate(iV(maxval(param%sum_N_sl),4,Ninstance), source=0) - allocate(iD(maxval(param%sum_N_sl),2,Ninstance), source=0) + allocate(iRhoU(maxval(param%sum_N_sl),4,Ninstances), source=0) + allocate(iV(maxval(param%sum_N_sl),4,Ninstances), source=0) + allocate(iD(maxval(param%sum_N_sl),2,Ninstances), source=0) i = 0 do p = 1, phases%length @@ -519,7 +519,7 @@ module function plastic_nonlocal_init() result(myPlasticity) if(.not. myPlasticity(p)) cycle i = i + 1 - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP + Nconstituents = count(material_phaseAt==p) * discretization_nIPs l = 0 do t = 1,4 do s = 1,param(i)%sum_N_sl @@ -976,7 +976,7 @@ module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, & real(pReal), dimension(3,3), intent(in) :: & Mp !< MandelStress - real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: & + real(pReal), dimension(3,3,homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems), intent(in) :: & F, & !< elastic deformation gradient Fp !< plastic deformation gradient real(pReal), intent(in) :: & @@ -1176,7 +1176,7 @@ end subroutine plastic_nonlocal_dotState !--------------------------------------------------------------------------------------------------- function rhoDotFlux(F,Fp,timestep, instance,of,ip,el) - real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: & + real(pReal), dimension(3,3,homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems), intent(in) :: & F, & !< elastic deformation gradient Fp !< plastic deformation gradient real(pReal), intent(in) :: & @@ -1416,7 +1416,7 @@ end function rhoDotFlux !-------------------------------------------------------------------------------------------------- module subroutine plastic_nonlocal_updateCompatibility(orientation,instance,i,e) - type(rotation), dimension(1,discretization_nIP,discretization_nElem), intent(in) :: & + type(rotation), dimension(1,discretization_nIPs,discretization_Nelems), intent(in) :: & orientation ! crystal orientation integer, intent(in) :: & instance, & @@ -1601,13 +1601,13 @@ end subroutine plastic_nonlocal_results !-------------------------------------------------------------------------------------------------- !> @brief populates the initial dislocation density !-------------------------------------------------------------------------------------------------- -subroutine stateInit(ini,phase,NipcMyPhase,instance) +subroutine stateInit(ini,phase,Nconstituents,instance) type(tInitialParameters) :: & ini integer,intent(in) :: & phase, & - NipcMyPhase, & + Nconstituents, & instance integer :: & e, & @@ -1625,15 +1625,15 @@ subroutine stateInit(ini,phase,NipcMyPhase,instance) totalVolume, & densityBinning, & minimumIpVolume - real(pReal), dimension(NipcMyPhase) :: & + real(pReal), dimension(Nconstituents) :: & volume associate(stt => state(instance)) if (ini%random_rho_u > 0.0_pReal) then ! randomly distribute dislocation segments on random slip system and of random type in the volume - do e = 1,discretization_nElem - do i = 1,discretization_nIP + do e = 1,discretization_Nelems + do i = 1,discretization_nIPs if (material_phaseAt(1,e) == phase) volume(material_phasememberAt(1,i,e)) = IPvolume(i,e) enddo enddo @@ -1645,13 +1645,13 @@ subroutine stateInit(ini,phase,NipcMyPhase,instance) meanDensity = 0.0_pReal do while(meanDensity < ini%random_rho_u) call random_number(rnd) - phasemember = nint(rnd(1)*real(NipcMyPhase,pReal) + 0.5_pReal) + phasemember = nint(rnd(1)*real(Nconstituents,pReal) + 0.5_pReal) s = nint(rnd(2)*real(sum(ini%N_sl),pReal)*4.0_pReal + 0.5_pReal) meanDensity = meanDensity + densityBinning * volume(phasemember) / totalVolume stt%rhoSglMobile(s,phasemember) = densityBinning enddo else ! homogeneous distribution with noise - do e = 1, NipcMyPhase + do e = 1, Nconstituents do f = 1,size(ini%N_sl,1) from = 1 + sum(ini%N_sl(1:f-1)) upto = sum(ini%N_sl(1:f)) diff --git a/src/constitutive_plastic_phenopowerlaw.f90 b/src/constitutive_plastic_phenopowerlaw.f90 index fd1b5fbbb..72fa0e9e6 100644 --- a/src/constitutive_plastic_phenopowerlaw.f90 +++ b/src/constitutive_plastic_phenopowerlaw.f90 @@ -70,9 +70,9 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState, & startIndex, endIndex integer, dimension(:), allocatable :: & @@ -91,13 +91,13 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_phenopowerlaw init -+>>>' myPlasticity = plastic_active('phenopowerlaw') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -224,20 +224,20 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['xi_sl ','gamma_sl']) * prm%sum_N_sl & + size(['xi_tw ','gamma_tw']) * prm%sum_N_tw sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,0) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization startIndex = 1 endIndex = prm%sum_N_sl stt%xi_slip => plasticState(p)%state (startIndex:endIndex,:) - stt%xi_slip = spread(xi_0_sl, 2, NipcMyPhase) + stt%xi_slip = spread(xi_0_sl, 2, Nconstituents) dot%xi_slip => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_xi',defaultVal=1.0_pReal) if(any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_xi' @@ -245,7 +245,7 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_tw stt%xi_twin => plasticState(p)%state (startIndex:endIndex,:) - stt%xi_twin = spread(xi_0_tw, 2, NipcMyPhase) + stt%xi_twin = spread(xi_0_tw, 2, Nconstituents) dot%xi_twin => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_xi',defaultVal=1.0_pReal) if(any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_xi' diff --git a/src/constitutive_thermal.f90 b/src/constitutive_thermal.f90 index 3aefb99a7..a7d5d3259 100644 --- a/src/constitutive_thermal.f90 +++ b/src/constitutive_thermal.f90 @@ -95,7 +95,7 @@ module subroutine constitutive_thermal_getRateAndItsTangents(TDot, dTDot_dT, T, homog = material_homogenizationAt(el) instance = thermal_typeInstance(homog) - do grain = 1, homogenization_Ngrains(homog) + do grain = 1, homogenization_Nconstituents(homog) phase = material_phaseAt(grain,el) constituent = material_phasememberAt(grain,ip,el) do source = 1, phase_Nsources(phase) diff --git a/src/crystallite.f90 b/src/crystallite.f90 index dc6b9f7da..c00b076a3 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -135,7 +135,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine crystallite_init - logical, dimension(discretization_nIP,discretization_nElem) :: devNull + logical, dimension(discretization_nIPs,discretization_Nelems) :: devNull integer :: & c, & !< counter in integration point component loop i, & !< counter in integration point loop @@ -162,9 +162,9 @@ subroutine crystallite_init debugCrystallite%ip = config_debug%get_asInt('integrationpoint', defaultVal=1) debugCrystallite%grain = config_debug%get_asInt('grain', defaultVal=1) - cMax = homogenization_maxNgrains - iMax = discretization_nIP - eMax = discretization_nElem + cMax = homogenization_maxNconstituents + iMax = discretization_nIPs + eMax = discretization_Nelems allocate(crystallite_partitionedF(3,3,cMax,iMax,eMax),source=0.0_pReal) @@ -253,7 +253,7 @@ subroutine crystallite_init ! initialize !$OMP PARALLEL DO PRIVATE(i,c) do e = FEsolving_execElem(1),FEsolving_execElem(2) - do i = FEsolving_execIP(1), FEsolving_execIP(2); do c = 1, homogenization_Ngrains(material_homogenizationAt(e)) + do i = FEsolving_execIP(1), FEsolving_execIP(2); do c = 1, homogenization_Nconstituents(material_homogenizationAt(e)) crystallite_Fp0(1:3,1:3,c,i,e) = material_orientation0(c,i,e)%asMatrix() ! Fp reflects initial orientation (see 10.1016/j.actamat.2006.01.005) crystallite_Fp0(1:3,1:3,c,i,e) = crystallite_Fp0(1:3,1:3,c,i,e) & / math_det33(crystallite_Fp0(1:3,1:3,c,i,e))**(1.0_pReal/3.0_pReal) @@ -279,7 +279,7 @@ subroutine crystallite_init !$OMP PARALLEL DO do e = FEsolving_execElem(1),FEsolving_execElem(2) do i = FEsolving_execIP(1),FEsolving_execIP(2) - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) call constitutive_dependentState(crystallite_partitionedF0(1:3,1:3,c,i,e), & crystallite_partitionedFp0(1:3,1:3,c,i,e), & c,i,e) ! update dependent state variables to be consistent with basic states @@ -307,7 +307,7 @@ end subroutine crystallite_init !-------------------------------------------------------------------------------------------------- function crystallite_stress() - logical, dimension(discretization_nIP,discretization_nElem) :: crystallite_stress + logical, dimension(discretization_nIPs,discretization_Nelems) :: crystallite_stress real(pReal) :: & formerSubStep integer :: & @@ -317,7 +317,7 @@ function crystallite_stress() e, & !< counter in element loop startIP, endIP, & s - logical, dimension(homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: todo !ToDo: need to set some values to false for different Ngrains + logical, dimension(homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems) :: todo !ToDo: need to set some values to false for different Ngrains real(pReal), dimension(:,:,:,:,:), allocatable :: & subLp0,& !< plastic velocity grad at start of crystallite inc subLi0 !< intermediate velocity grad at start of crystallite inc @@ -335,7 +335,7 @@ function crystallite_stress() crystallite_subStep = 0.0_pReal !$OMP PARALLEL DO elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2) - do i = FEsolving_execIP(1),FEsolving_execIP(2); do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do i = FEsolving_execIP(1),FEsolving_execIP(2); do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) homogenizationRequestsCalculation: if (crystallite_requested(c,i,e)) then plasticState (material_phaseAt(c,e))%subState0( :,material_phaseMemberAt(c,i,e)) = & plasticState (material_phaseAt(c,e))%partitionedState0(:,material_phaseMemberAt(c,i,e)) @@ -362,7 +362,7 @@ function crystallite_stress() endIP = startIP else singleRun startIP = 1 - endIP = discretization_nIP + endIP = discretization_nIPs endif singleRun NiterationCrystallite = 0 @@ -376,7 +376,7 @@ function crystallite_stress() !$OMP PARALLEL DO PRIVATE(formerSubStep) elementLooping3: do e = FEsolving_execElem(1),FEsolving_execElem(2) do i = FEsolving_execIP(1),FEsolving_execIP(2) - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) !-------------------------------------------------------------------------------------------------- ! wind forward if (crystallite_converged(c,i,e)) then @@ -472,7 +472,7 @@ subroutine crystallite_initializeRestorationPoints(i,e) c, & !< constituent number s - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) crystallite_partitionedFp0(1:3,1:3,c,i,e) = crystallite_Fp0(1:3,1:3,c,i,e) crystallite_partitionedLp0(1:3,1:3,c,i,e) = crystallite_Lp0(1:3,1:3,c,i,e) crystallite_partitionedFi0(1:3,1:3,c,i,e) = crystallite_Fi0(1:3,1:3,c,i,e) @@ -503,7 +503,7 @@ subroutine crystallite_windForward(i,e) c, & !< constituent number s - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) crystallite_partitionedF0 (1:3,1:3,c,i,e) = crystallite_partitionedF(1:3,1:3,c,i,e) crystallite_partitionedFp0(1:3,1:3,c,i,e) = crystallite_Fp (1:3,1:3,c,i,e) crystallite_partitionedLp0(1:3,1:3,c,i,e) = crystallite_Lp (1:3,1:3,c,i,e) @@ -536,7 +536,7 @@ subroutine crystallite_restore(i,e,includeL) c, & !< constituent number s - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) if (includeL) then crystallite_Lp(1:3,1:3,c,i,e) = crystallite_partitionedLp0(1:3,1:3,c,i,e) crystallite_Li(1:3,1:3,c,i,e) = crystallite_partitionedLi0(1:3,1:3,c,i,e) @@ -697,7 +697,7 @@ subroutine crystallite_orientations !$OMP PARALLEL DO do e = FEsolving_execElem(1),FEsolving_execElem(2) do i = FEsolving_execIP(1),FEsolving_execIP(2) - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) call crystallite_orientation(c,i,e)%fromMatrix(transpose(math_rotationalPart(crystallite_Fe(1:3,1:3,c,i,e)))) enddo; enddo; enddo !$OMP END PARALLEL DO @@ -821,11 +821,11 @@ subroutine crystallite_results real(pReal), allocatable, dimension(:,:,:) :: select_tensors integer :: e,i,c,j - allocate(select_tensors(3,3,count(material_phaseAt==instance)*discretization_nIP)) + allocate(select_tensors(3,3,count(material_phaseAt==instance)*discretization_nIPs)) j=0 do e = 1, size(material_phaseAt,2) - do i = 1, discretization_nIP + do i = 1, discretization_nIPs do c = 1, size(material_phaseAt,1) !ToDo: this needs to be changed for varying Ngrains if (material_phaseAt(c,e) == instance) then j = j + 1 @@ -848,11 +848,11 @@ subroutine crystallite_results type(rotation), allocatable, dimension(:) :: select_rotations integer :: e,i,c,j - allocate(select_rotations(count(material_phaseAt==instance)*homogenization_maxNgrains*discretization_nIP)) + allocate(select_rotations(count(material_phaseAt==instance)*homogenization_maxNconstituents*discretization_nIPs)) j=0 do e = 1, size(material_phaseAt,2) - do i = 1, discretization_nIP + do i = 1, discretization_nIPs do c = 1, size(material_phaseAt,1) !ToDo: this needs to be changed for varying Ngrains if (material_phaseAt(c,e) == instance) then j = j + 1 @@ -1571,7 +1571,7 @@ subroutine crystallite_restartWrite call HDF5_closeGroup(groupHandle) groupHandle = HDF5_addGroup(fileHandle,'materialpoint') - do i = 1, material_Nhomogenization + do i = 1, size(material_name_homogenization) write(datasetName,'(i0,a)') i,'_omega_homogenization' call HDF5_write(groupHandle,homogState(i)%state,datasetName) enddo @@ -1612,7 +1612,7 @@ subroutine crystallite_restartRead call HDF5_closeGroup(groupHandle) groupHandle = HDF5_openGroup(fileHandle,'materialpoint') - do i = 1, material_Nhomogenization + do i = 1,size(material_name_homogenization) write(datasetName,'(i0,a)') i,'_omega_homogenization' call HDF5_read(groupHandle,homogState(i)%state0,datasetName) enddo @@ -1645,7 +1645,7 @@ subroutine crystallite_forward do j = 1,phase_Nsources(i) sourceState(i)%p(j)%state0 = sourceState(i)%p(j)%state enddo; enddo - do i = 1, material_Nhomogenization + do i = 1,size(material_name_homogenization) homogState (i)%state0 = homogState (i)%state thermalState(i)%state0 = thermalState(i)%state damageState (i)%state0 = damageState (i)%state diff --git a/src/damage_local.f90 b/src/damage_local.f90 index af2532184..e63db90b0 100644 --- a/src/damage_local.f90 +++ b/src/damage_local.f90 @@ -42,7 +42,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine damage_local_init - integer :: Ninstance,NofMyHomog,h + integer :: Ninstances,Nmaterialpoints,h class(tNode), pointer :: & num_generic, & material_homogenization, & @@ -57,8 +57,8 @@ subroutine damage_local_init num%residualStiffness = num_generic%get_asFloat('residualStiffness', defaultVal=1.0e-6_pReal) if (num%residualStiffness < 0.0_pReal) call IO_error(301,ext_msg='residualStiffness') - Ninstance = count(damage_type == DAMAGE_local_ID) - allocate(param(Ninstance)) + Ninstances = count(damage_type == DAMAGE_local_ID) + allocate(param(Ninstances)) material_homogenization => config_material%get('homogenization') do h = 1, material_homogenization%length @@ -73,11 +73,11 @@ subroutine damage_local_init prm%output = homogDamage%get_asStrings('output',defaultVal=emptyStringArray) #endif - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) damageState(h)%sizeState = 1 - allocate(damageState(h)%state0 (1,NofMyHomog), source=damage_initialPhi(h)) - allocate(damageState(h)%subState0(1,NofMyHomog), source=damage_initialPhi(h)) - allocate(damageState(h)%state (1,NofMyHomog), source=damage_initialPhi(h)) + allocate(damageState(h)%state0 (1,Nmaterialpoints), source=damage_initialPhi(h)) + allocate(damageState(h)%subState0(1,Nmaterialpoints), source=damage_initialPhi(h)) + allocate(damageState(h)%state (1,Nmaterialpoints), source=damage_initialPhi(h)) nullify(damageMapping(h)%p) damageMapping(h)%p => material_homogenizationMemberAt @@ -143,8 +143,8 @@ subroutine damage_local_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, phi, ip, el call constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) - phiDot = phiDot/real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) - dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + phiDot = phiDot/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) + dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end subroutine damage_local_getSourceAndItsTangent diff --git a/src/damage_none.f90 b/src/damage_none.f90 index 4e01998bc..2279bc06b 100644 --- a/src/damage_none.f90 +++ b/src/damage_none.f90 @@ -16,18 +16,18 @@ contains !-------------------------------------------------------------------------------------------------- subroutine damage_none_init - integer :: h,NofMyHomog + integer :: h,Nmaterialpoints print'(/,a)', ' <<<+- damage_none init -+>>>'; flush(6) - do h = 1, material_Nhomogenization + do h = 1, size(material_name_homogenization) if (damage_type(h) /= DAMAGE_NONE_ID) cycle - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) damageState(h)%sizeState = 0 - allocate(damageState(h)%state0 (0,NofMyHomog)) - allocate(damageState(h)%subState0(0,NofMyHomog)) - allocate(damageState(h)%state (0,NofMyHomog)) + allocate(damageState(h)%state0 (0,Nmaterialpoints)) + allocate(damageState(h)%subState0(0,Nmaterialpoints)) + allocate(damageState(h)%state (0,Nmaterialpoints)) deallocate(damage(h)%p) allocate (damage(h)%p(1), source=damage_initialPhi(h)) diff --git a/src/damage_nonlocal.f90 b/src/damage_nonlocal.f90 index f31b10242..24a51cf54 100644 --- a/src/damage_nonlocal.f90 +++ b/src/damage_nonlocal.f90 @@ -46,7 +46,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine damage_nonlocal_init - integer :: Ninstance,NofMyHomog,h + integer :: Ninstances,Nmaterialpoints,h class(tNode), pointer :: & num_generic, & material_homogenization, & @@ -60,8 +60,8 @@ subroutine damage_nonlocal_init num_generic => config_numerics%get('generic',defaultVal= emptyDict) num%charLength = num_generic%get_asFloat('charLength',defaultVal=1.0_pReal) - Ninstance = count(damage_type == DAMAGE_nonlocal_ID) - allocate(param(Ninstance)) + Ninstances = count(damage_type == DAMAGE_nonlocal_ID) + allocate(param(Ninstances)) material_homogenization => config_material%get('homogenization') do h = 1, material_homogenization%length @@ -76,11 +76,11 @@ subroutine damage_nonlocal_init prm%output = homogDamage%get_asStrings('output',defaultVal=emptyStringArray) #endif - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) damageState(h)%sizeState = 1 - allocate(damageState(h)%state0 (1,NofMyHomog), source=damage_initialPhi(h)) - allocate(damageState(h)%subState0(1,NofMyHomog), source=damage_initialPhi(h)) - allocate(damageState(h)%state (1,NofMyHomog), source=damage_initialPhi(h)) + allocate(damageState(h)%state0 (1,Nmaterialpoints), source=damage_initialPhi(h)) + allocate(damageState(h)%subState0(1,Nmaterialpoints), source=damage_initialPhi(h)) + allocate(damageState(h)%state (1,Nmaterialpoints), source=damage_initialPhi(h)) nullify(damageMapping(h)%p) damageMapping(h)%p => material_homogenizationMemberAt @@ -110,8 +110,8 @@ subroutine damage_nonlocal_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, phi, ip, dPhiDot_dPhi = 0.0_pReal call constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) - phiDot = phiDot/real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) - dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + phiDot = phiDot/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) + dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end subroutine damage_nonlocal_getSourceAndItsTangent @@ -132,13 +132,13 @@ function damage_nonlocal_getDiffusion(ip,el) homog = material_homogenizationAt(el) damage_nonlocal_getDiffusion = 0.0_pReal - do grain = 1, homogenization_Ngrains(homog) + do grain = 1, homogenization_Nconstituents(homog) damage_nonlocal_getDiffusion = damage_nonlocal_getDiffusion + & crystallite_push33ToRef(grain,ip,el,lattice_D(1:3,1:3,material_phaseAt(grain,el))) enddo damage_nonlocal_getDiffusion = & - num%charLength**2*damage_nonlocal_getDiffusion/real(homogenization_Ngrains(homog),pReal) + num%charLength**2*damage_nonlocal_getDiffusion/real(homogenization_Nconstituents(homog),pReal) end function damage_nonlocal_getDiffusion @@ -156,12 +156,12 @@ real(pReal) function damage_nonlocal_getMobility(ip,el) damage_nonlocal_getMobility = 0.0_pReal - do ipc = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do ipc = 1, homogenization_Nconstituents(material_homogenizationAt(el)) damage_nonlocal_getMobility = damage_nonlocal_getMobility + lattice_M(material_phaseAt(ipc,el)) enddo damage_nonlocal_getMobility = damage_nonlocal_getMobility/& - real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function damage_nonlocal_getMobility diff --git a/src/discretization.f90 b/src/discretization.f90 index e6e53fcf4..0b8925e4a 100644 --- a/src/discretization.f90 +++ b/src/discretization.f90 @@ -11,8 +11,8 @@ module discretization private integer, public, protected :: & - discretization_nIP, & - discretization_nElem + discretization_nIPs, & + discretization_Nelems integer, public, protected, dimension(:), allocatable :: & discretization_materialAt @@ -51,8 +51,8 @@ subroutine discretization_init(materialAt,& print'(/,a)', ' <<<+- discretization init -+>>>'; flush(6) - discretization_nElem = size(materialAt,1) - discretization_nIP = size(IPcoords0,2)/discretization_nElem + discretization_Nelems = size(materialAt,1) + discretization_nIPs = size(IPcoords0,2)/discretization_Nelems discretization_materialAt = materialAt diff --git a/src/grid/DAMASK_grid.f90 b/src/grid/DAMASK_grid.f90 index 4b8302def..1cb3acf16 100644 --- a/src/grid/DAMASK_grid.f90 +++ b/src/grid/DAMASK_grid.f90 @@ -210,7 +210,7 @@ program DAMASK_grid if(.not. step_discretization%contains('t')) call IO_error(error_ID=837,ext_msg = 't missing') if(.not. step_discretization%contains('N')) call IO_error(error_ID=837,ext_msg = 'N missing') loadCases(l)%time = step_discretization%get_asFloat('t') - loadCases(l)%incs = step_discretization%get_asFloat('N') + loadCases(l)%incs = step_discretization%get_asInt ('N') loadCases(l)%logscale = step_discretization%get_asBool ('log_timestep', defaultVal= .false.) loadCases(l)%outputfrequency = step_discretization%get_asInt ('f_out', defaultVal=1) loadCases(l)%restartfrequency = step_discretization%get_asInt ('f_restart', defaultVal=huge(0)) diff --git a/src/grid/grid_mech_FEM.f90 b/src/grid/grid_mech_FEM.f90 index dfbd3f2f3..bf3d7752d 100644 --- a/src/grid/grid_mech_FEM.f90 +++ b/src/grid/grid_mech_FEM.f90 @@ -236,7 +236,7 @@ subroutine grid_mech_FEM_init F = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid3) endif restartRead - materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent + homogenization_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent call utilities_updateCoords(F) call utilities_constitutiveResponse(P_current,temp33_Real,C_volAvg,devNull, & ! stress field, stress avg, global average of stiffness and (min+max)/2 F, & ! target F @@ -364,7 +364,7 @@ subroutine grid_mech_FEM_forward(cutBack,guess,timeinc,timeinc_old,loadCaseTime, F_lastInc = F - materialpoint_F0 = reshape(F, [3,3,1,product(grid(1:2))*grid3]) + homogenization_F0 = reshape(F, [3,3,1,product(grid(1:2))*grid3]) endif !-------------------------------------------------------------------------------------------------- @@ -557,9 +557,9 @@ subroutine formResidual(da_local,x_local, & ii = i-xstart+1; jj = j-ystart+1; kk = k-zstart+1 ele = ele + 1 f_elem = matmul(transpose(BMat),transpose(P_current(1:3,1:3,ii,jj,kk)))*detJ + & - matmul(HGMat,x_elem)*(materialpoint_dPdF(1,1,1,1,1,ele) + & - materialpoint_dPdF(2,2,2,2,1,ele) + & - materialpoint_dPdF(3,3,3,3,1,ele))/3.0_pReal + matmul(HGMat,x_elem)*(homogenization_dPdF(1,1,1,1,1,ele) + & + homogenization_dPdF(2,2,2,2,1,ele) + & + homogenization_dPdF(3,3,3,3,1,ele))/3.0_pReal ctr = 0 do kk = 0, 1; do jj = 0, 1; do ii = 0, 1 ctr = ctr + 1 @@ -636,18 +636,18 @@ subroutine formJacobian(da_local,x_local,Jac_pre,Jac,dummy,ierr) row = col ele = ele + 1 K_ele = 0.0 - K_ele(1 :8 ,1 :8 ) = HGMat*(materialpoint_dPdF(1,1,1,1,1,ele) + & - materialpoint_dPdF(2,2,2,2,1,ele) + & - materialpoint_dPdF(3,3,3,3,1,ele))/3.0_pReal - K_ele(9 :16,9 :16) = HGMat*(materialpoint_dPdF(1,1,1,1,1,ele) + & - materialpoint_dPdF(2,2,2,2,1,ele) + & - materialpoint_dPdF(3,3,3,3,1,ele))/3.0_pReal - K_ele(17:24,17:24) = HGMat*(materialpoint_dPdF(1,1,1,1,1,ele) + & - materialpoint_dPdF(2,2,2,2,1,ele) + & - materialpoint_dPdF(3,3,3,3,1,ele))/3.0_pReal + K_ele(1 :8 ,1 :8 ) = HGMat*(homogenization_dPdF(1,1,1,1,1,ele) + & + homogenization_dPdF(2,2,2,2,1,ele) + & + homogenization_dPdF(3,3,3,3,1,ele))/3.0_pReal + K_ele(9 :16,9 :16) = HGMat*(homogenization_dPdF(1,1,1,1,1,ele) + & + homogenization_dPdF(2,2,2,2,1,ele) + & + homogenization_dPdF(3,3,3,3,1,ele))/3.0_pReal + K_ele(17:24,17:24) = HGMat*(homogenization_dPdF(1,1,1,1,1,ele) + & + homogenization_dPdF(2,2,2,2,1,ele) + & + homogenization_dPdF(3,3,3,3,1,ele))/3.0_pReal K_ele = K_ele + & matmul(transpose(BMatFull), & - matmul(reshape(reshape(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,ele), & + matmul(reshape(reshape(homogenization_dPdF(1:3,1:3,1:3,1:3,1,ele), & shape=[3,3,3,3], order=[2,1,4,3]),shape=[9,9]),BMatFull))*detJ call MatSetValuesStencil(Jac,24,row,24,col,K_ele,ADD_VALUES,ierr) CHKERRQ(ierr) diff --git a/src/grid/grid_mech_spectral_basic.f90 b/src/grid/grid_mech_spectral_basic.f90 index d2f6b40da..ec6c3a540 100644 --- a/src/grid/grid_mech_spectral_basic.f90 +++ b/src/grid/grid_mech_spectral_basic.f90 @@ -198,7 +198,7 @@ subroutine grid_mech_spectral_basic_init F = reshape(F_lastInc,[9,grid(1),grid(2),grid3]) endif restartRead - materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent + homogenization_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent call utilities_updateCoords(reshape(F,shape(F_lastInc))) call utilities_constitutiveResponse(P,temp33_Real,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 reshape(F,shape(F_lastInc)), & ! target F @@ -324,7 +324,7 @@ subroutine grid_mech_spectral_basic_forward(cutBack,guess,timeinc,timeinc_old,lo rotation_BC%rotate(F_aimDot,active=.true.)) F_lastInc = reshape(F,[3,3,grid(1),grid(2),grid3]) - materialpoint_F0 = reshape(F, [3,3,1,product(grid(1:2))*grid3]) + homogenization_F0 = reshape(F, [3,3,1,product(grid(1:2))*grid3]) endif !-------------------------------------------------------------------------------------------------- diff --git a/src/grid/grid_mech_spectral_polarisation.f90 b/src/grid/grid_mech_spectral_polarisation.f90 index 3d495ddf0..8f9ea81b3 100644 --- a/src/grid/grid_mech_spectral_polarisation.f90 +++ b/src/grid/grid_mech_spectral_polarisation.f90 @@ -224,7 +224,7 @@ subroutine grid_mech_spectral_polarisation_init F_tau_lastInc = 2.0_pReal*F_lastInc endif restartRead - materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent + homogenization_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent call utilities_updateCoords(reshape(F,shape(F_lastInc))) call utilities_constitutiveResponse(P,temp33_Real,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 reshape(F,shape(F_lastInc)), & ! target F @@ -364,7 +364,7 @@ subroutine grid_mech_spectral_polarisation_forward(cutBack,guess,timeinc,timeinc F_lastInc = reshape(F, [3,3,grid(1),grid(2),grid3]) F_tau_lastInc = reshape(F_tau,[3,3,grid(1),grid(2),grid3]) - materialpoint_F0 = reshape(F,[3,3,1,product(grid(1:2))*grid3]) + homogenization_F0 = reshape(F,[3,3,1,product(grid(1:2))*grid3]) endif !-------------------------------------------------------------------------------------------------- @@ -606,7 +606,7 @@ subroutine formResidual(in, FandF_tau, & do k = 1, grid3; do j = 1, grid(2); do i = 1, grid(1) e = e + 1 residual_F(1:3,1:3,i,j,k) = & - math_mul3333xx33(math_invSym3333(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,e) + C_scale), & + math_mul3333xx33(math_invSym3333(homogenization_dPdF(1:3,1:3,1:3,1:3,1,e) + C_scale), & residual_F(1:3,1:3,i,j,k) - matmul(F(1:3,1:3,i,j,k), & math_mul3333xx33(C_scale,F_tau(1:3,1:3,i,j,k) - F(1:3,1:3,i,j,k) - math_I3))) & + residual_F_tau(1:3,1:3,i,j,k) diff --git a/src/grid/spectral_utilities.f90 b/src/grid/spectral_utilities.f90 index 77047a317..fddd1885f 100644 --- a/src/grid/spectral_utilities.f90 +++ b/src/grid/spectral_utilities.f90 @@ -802,7 +802,7 @@ end subroutine utilities_fourierTensorDivergence !-------------------------------------------------------------------------------------------------- -!> @brief calculate constitutive response from materialpoint_F0 to F during timeinc +!> @brief calculate constitutive response from homogenization_F0 to F during timeinc !-------------------------------------------------------------------------------------------------- subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& F,timeinc,rotation_BC) @@ -824,11 +824,11 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& print'(/,a)', ' ... evaluating constitutive response ......................................' flush(IO_STDOUT) - materialpoint_F = reshape(F,[3,3,1,product(grid(1:2))*grid3]) ! set materialpoint target F to estimated field + homogenization_F = reshape(F,[3,3,1,product(grid(1:2))*grid3]) ! set materialpoint target F to estimated field call materialpoint_stressAndItsTangent(timeinc) ! calculate P field - P = reshape(materialpoint_P, [3,3,grid(1),grid(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 ! average of P call MPI_Allreduce(MPI_IN_PLACE,P_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) if (debugRotation) & @@ -845,13 +845,13 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& dPdF_min = huge(1.0_pReal) dPdF_norm_min = huge(1.0_pReal) do i = 1, product(grid(1:2))*grid3 - if (dPdF_norm_max < sum(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal)) then - dPdF_max = materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i) - dPdF_norm_max = sum(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal) + if (dPdF_norm_max < sum(homogenization_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal)) then + dPdF_max = homogenization_dPdF(1:3,1:3,1:3,1:3,1,i) + dPdF_norm_max = sum(homogenization_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal) endif - if (dPdF_norm_min > sum(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal)) then - dPdF_min = materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i) - dPdF_norm_min = sum(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal) + if (dPdF_norm_min > sum(homogenization_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal)) then + dPdF_min = homogenization_dPdF(1:3,1:3,1:3,1:3,1,i) + dPdF_norm_min = sum(homogenization_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal) endif end do @@ -869,7 +869,7 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& C_minmaxAvg = 0.5_pReal*(dPdF_max + dPdF_min) - C_volAvg = sum(sum(materialpoint_dPdF,dim=6),dim=5) + C_volAvg = sum(sum(homogenization_dPdF,dim=6),dim=5) call MPI_Allreduce(MPI_IN_PLACE,C_volAvg,81,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) C_volAvg = C_volAvg * wgt diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 5aa10fcf3..02f1d3bf0 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -31,12 +31,12 @@ module homogenization !-------------------------------------------------------------------------------------------------- ! General variables for the homogenization at a material point real(pReal), dimension(:,:,:,:), allocatable, public :: & - materialpoint_F0, & !< def grad of IP at start of FE increment - materialpoint_F !< def grad of IP to be reached at end of FE increment + homogenization_F0, & !< def grad of IP at start of FE increment + homogenization_F !< def grad of IP to be reached at end of FE increment real(pReal), dimension(:,:,:,:), allocatable, public, protected :: & - materialpoint_P !< first P--K stress of IP + homogenization_P !< first P--K stress of IP real(pReal), dimension(:,:,:,:,:,:), allocatable, public, protected :: & - materialpoint_dPdF !< tangent of first P--K stress at IP + homogenization_dPdF !< tangent of first P--K stress at IP type :: tNumerics integer :: & @@ -158,7 +158,7 @@ subroutine homogenization_init debugHomog%grain = config_debug%get_asInt('grain',defaultVal = 1) if (debugHomog%grain < 1 & - .or. debugHomog%grain > homogenization_Ngrains(material_homogenizationAt(debugHomog%element))) & + .or. debugHomog%grain > homogenization_Nconstituents(material_homogenizationAt(debugHomog%element))) & call IO_error(602,ext_msg='constituent', el=debugHomog%element, g=debugHomog%grain) @@ -181,10 +181,10 @@ subroutine homogenization_init !-------------------------------------------------------------------------------------------------- ! allocate and initialize global variables - allocate(materialpoint_dPdF(3,3,3,3,discretization_nIP,discretization_nElem), source=0.0_pReal) - materialpoint_F0 = spread(spread(math_I3,3,discretization_nIP),4,discretization_nElem) ! initialize to identity - materialpoint_F = materialpoint_F0 ! initialize to identity - allocate(materialpoint_P(3,3,discretization_nIP,discretization_nElem), source=0.0_pReal) + allocate(homogenization_dPdF(3,3,3,3,discretization_nIPs,discretization_Nelems), source=0.0_pReal) + homogenization_F0 = spread(spread(math_I3,3,discretization_nIPs),4,discretization_Nelems) ! initialize to identity + homogenization_F = homogenization_F0 ! initialize to identity + allocate(homogenization_P(3,3,discretization_nIPs,discretization_Nelems), source=0.0_pReal) print'(/,a)', ' <<<+- homogenization init -+>>>'; flush(IO_STDOUT) @@ -213,13 +213,13 @@ subroutine materialpoint_stressAndItsTangent(dt) i, & !< integration point number e, & !< element number myNgrains - real(pReal), dimension(discretization_nIP,discretization_nElem) :: & + real(pReal), dimension(discretization_nIPs,discretization_Nelems) :: & subFrac, & subStep - logical, dimension(discretization_nIP,discretization_nElem) :: & + logical, dimension(discretization_nIPs,discretization_Nelems) :: & requested, & converged - logical, dimension(2,discretization_nIP,discretization_nElem) :: & + logical, dimension(2,discretization_nIPs,discretization_Nelems) :: & doneAndHappy @@ -257,7 +257,7 @@ subroutine materialpoint_stressAndItsTangent(dt) !$OMP PARALLEL DO elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2) - myNgrains = homogenization_Ngrains(material_homogenizationAt(e)) + myNgrains = homogenization_Nconstituents(material_homogenizationAt(e)) IpLooping1: do i = FEsolving_execIP(1),FEsolving_execIP(2) if (converged(i,e)) then @@ -327,11 +327,11 @@ subroutine materialpoint_stressAndItsTangent(dt) ! deformation partitioning !$OMP PARALLEL DO PRIVATE(myNgrains) elementLooping2: do e = FEsolving_execElem(1),FEsolving_execElem(2) - myNgrains = homogenization_Ngrains(material_homogenizationAt(e)) + myNgrains = homogenization_Nconstituents(material_homogenizationAt(e)) IpLooping2: do i = FEsolving_execIP(1),FEsolving_execIP(2) if(requested(i,e) .and. .not. doneAndHappy(1,i,e)) then ! requested but not yet done - call partitionDeformation(materialpoint_F0(1:3,1:3,i,e) & - + (materialpoint_F(1:3,1:3,i,e)-materialpoint_F0(1:3,1:3,i,e))& + call partitionDeformation(homogenization_F0(1:3,1:3,i,e) & + + (homogenization_F(1:3,1:3,i,e)-homogenization_F0(1:3,1:3,i,e))& *(subStep(i,e)+subFrac(i,e)), & i,e) crystallite_dt(1:myNgrains,i,e) = dt*subStep(i,e) ! propagate materialpoint dt to grains @@ -357,8 +357,8 @@ subroutine materialpoint_stressAndItsTangent(dt) doneAndHappy(1:2,i,e) = [.true.,.false.] else doneAndHappy(1:2,i,e) = updateState(dt*subStep(i,e), & - materialpoint_F0(1:3,1:3,i,e) & - + (materialpoint_F(1:3,1:3,i,e)-materialpoint_F0(1:3,1:3,i,e)) & + homogenization_F0(1:3,1:3,i,e) & + + (homogenization_F(1:3,1:3,i,e)-homogenization_F0(1:3,1:3,i,e)) & *(subStep(i,e)+subFrac(i,e)), & i,e) converged(i,e) = all(doneAndHappy(1:2,i,e)) ! converged if done and happy @@ -408,12 +408,12 @@ subroutine partitionDeformation(subF,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization call mech_isostrain_partitionDeformation(& - crystallite_partitionedF(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & subF) case (HOMOGENIZATION_RGC_ID) chosenHomogenization call mech_RGC_partitionDeformation(& - crystallite_partitionedF(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & subF,& ip, & el) @@ -437,19 +437,19 @@ function updateState(subdt,subF,ip,el) el !< element number integer :: c logical, dimension(2) :: updateState - real(pReal) :: dPdFs(3,3,3,3,homogenization_Ngrains(material_homogenizationAt(el))) + real(pReal) :: dPdFs(3,3,3,3,homogenization_Nconstituents(material_homogenizationAt(el))) updateState = .true. chosenHomogenization: select case(homogenization_type(material_homogenizationAt(el))) case (HOMOGENIZATION_RGC_ID) chosenHomogenization - do c=1,homogenization_Ngrains(material_homogenizationAt(el)) + do c=1,homogenization_Nconstituents(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo updateState = & updateState .and. & - mech_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & - crystallite_partitionedF(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & - crystallite_partitionedF0(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el),& + mech_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF0(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el),& subF,& subdt, & dPdFs, & @@ -487,33 +487,33 @@ subroutine averageStressAndItsTangent(ip,el) ip, & !< integration point el !< element number integer :: c - real(pReal) :: dPdFs(3,3,3,3,homogenization_Ngrains(material_homogenizationAt(el))) + real(pReal) :: dPdFs(3,3,3,3,homogenization_Nconstituents(material_homogenizationAt(el))) chosenHomogenization: select case(homogenization_type(material_homogenizationAt(el))) case (HOMOGENIZATION_NONE_ID) chosenHomogenization - materialpoint_P(1:3,1:3,ip,el) = crystallite_P(1:3,1:3,1,ip,el) - materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el) = crystallite_stressTangent(1,ip,el) + homogenization_P(1:3,1:3,ip,el) = crystallite_P(1:3,1:3,1,ip,el) + homogenization_dPdF(1:3,1:3,1:3,1:3,ip,el) = crystallite_stressTangent(1,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization - do c = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do c = 1, homogenization_Nconstituents(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo call mech_isostrain_averageStressAndItsTangent(& - materialpoint_P(1:3,1:3,ip,el), & - materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el),& - crystallite_P(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & + homogenization_P(1:3,1:3,ip,el), & + homogenization_dPdF(1:3,1:3,1:3,1:3,ip,el),& + crystallite_P(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & dPdFs, & homogenization_typeInstance(material_homogenizationAt(el))) case (HOMOGENIZATION_RGC_ID) chosenHomogenization - do c = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do c = 1, homogenization_Nconstituents(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo call mech_RGC_averageStressAndItsTangent(& - materialpoint_P(1:3,1:3,ip,el), & - materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el),& - crystallite_P(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & + homogenization_P(1:3,1:3,ip,el), & + homogenization_dPdF(1:3,1:3,1:3,1:3,ip,el),& + crystallite_P(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & dPdFs, & homogenization_typeInstance(material_homogenizationAt(el))) end select chosenHomogenization @@ -539,10 +539,10 @@ subroutine homogenization_results group = trim(group_base)//'/generic' call results_closeGroup(results_addGroup(group)) - !temp = reshape(materialpoint_F,[3,3,discretization_nIP*discretization_nElem]) + !temp = reshape(homogenization_F,[3,3,discretization_nIPs*discretization_Nelems]) !call results_writeDataset(group,temp,'F',& ! 'deformation gradient','1') - !temp = reshape(materialpoint_P,[3,3,discretization_nIP*discretization_nElem]) + !temp = reshape(homogenization_P,[3,3,discretization_nIPs*discretization_Nelems]) !call results_writeDataset(group,temp,'P',& ! '1st Piola-Kirchhoff stress','Pa') diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index ea5bffd28..5fcea1f8d 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -81,9 +81,9 @@ module subroutine mech_RGC_init(num_homogMech) num_homogMech !< pointer to mechanical homogenization numerics data integer :: & - Ninstance, & + Ninstances, & h, & - NofMyHomog, & + Nmaterialpoints, & sizeState, nIntFaceTot class (tNode), pointer :: & @@ -94,8 +94,8 @@ module subroutine mech_RGC_init(num_homogMech) print'(/,a)', ' <<<+- homogenization_mech_rgc init -+>>>' - Ninstance = count(homogenization_type == HOMOGENIZATION_RGC_ID) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) + Ninstances = count(homogenization_type == HOMOGENIZATION_RGC_ID) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) print*, 'Tjahjanto et al., International Journal of Material Forming 2(1):939–942, 2009' print*, 'https://doi.org/10.1007/s12289-009-0619-1'//IO_EOL @@ -105,10 +105,10 @@ module subroutine mech_RGC_init(num_homogMech) - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(state0(Ninstance)) - allocate(dependentState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(state0(Ninstances)) + allocate(dependentState(Ninstances)) num_RGC => num_homogMech%get('RGC',defaultVal=emptyDict) @@ -164,7 +164,7 @@ module subroutine mech_RGC_init(num_homogMech) #endif prm%N_constituents = homogMech%get_asInts('cluster_size',requiredSize=3) - if (homogenization_Ngrains(h) /= product(prm%N_constituents)) & + if (homogenization_Nconstituents(h) /= product(prm%N_constituents)) & call IO_error(211,ext_msg='N_constituents (mech_RGC)') prm%xi_alpha = homogMech%get_asFloat('xi_alpha') @@ -173,7 +173,7 @@ module subroutine mech_RGC_init(num_homogMech) prm%D_alpha = homogMech%get_asFloats('D_alpha', requiredSize=3) prm%a_g = homogMech%get_asFloats('a_g', requiredSize=3) - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) nIntFaceTot = 3*( (prm%N_constituents(1)-1)*prm%N_constituents(2)*prm%N_constituents(3) & + prm%N_constituents(1)*(prm%N_constituents(2)-1)*prm%N_constituents(3) & + prm%N_constituents(1)*prm%N_constituents(2)*(prm%N_constituents(3)-1)) @@ -181,24 +181,24 @@ module subroutine mech_RGC_init(num_homogMech) + size(['avg constitutive work ','average penalty energy']) homogState(h)%sizeState = sizeState - allocate(homogState(h)%state0 (sizeState,NofMyHomog), source=0.0_pReal) - allocate(homogState(h)%subState0(sizeState,NofMyHomog), source=0.0_pReal) - allocate(homogState(h)%state (sizeState,NofMyHomog), source=0.0_pReal) + allocate(homogState(h)%state0 (sizeState,Nmaterialpoints), source=0.0_pReal) + allocate(homogState(h)%subState0(sizeState,Nmaterialpoints), source=0.0_pReal) + allocate(homogState(h)%state (sizeState,Nmaterialpoints), source=0.0_pReal) stt%relaxationVector => homogState(h)%state(1:nIntFaceTot,:) st0%relaxationVector => homogState(h)%state0(1:nIntFaceTot,:) stt%work => homogState(h)%state(nIntFaceTot+1,:) stt%penaltyEnergy => homogState(h)%state(nIntFaceTot+2,:) - allocate(dst%volumeDiscrepancy( NofMyHomog), source=0.0_pReal) - allocate(dst%relaxationRate_avg( NofMyHomog), source=0.0_pReal) - allocate(dst%relaxationRate_max( NofMyHomog), source=0.0_pReal) - allocate(dst%mismatch( 3,NofMyHomog), source=0.0_pReal) + allocate(dst%volumeDiscrepancy( Nmaterialpoints), source=0.0_pReal) + allocate(dst%relaxationRate_avg( Nmaterialpoints), source=0.0_pReal) + allocate(dst%relaxationRate_max( Nmaterialpoints), source=0.0_pReal) + allocate(dst%mismatch( 3,Nmaterialpoints), source=0.0_pReal) !-------------------------------------------------------------------------------------------------- ! assigning cluster orientations - dependentState(homogenization_typeInstance(h))%orientation = spread(eu2om(prm%a_g*inRad),3,NofMyHomog) - !dst%orientation = spread(eu2om(prm%a_g*inRad),3,NofMyHomog) ifort version 18.0.1 crashes (for whatever reason) + dependentState(homogenization_typeInstance(h))%orientation = spread(eu2om(prm%a_g*inRad),3,Nmaterialpoints) + !dst%orientation = spread(eu2om(prm%a_g*inRad),3,Nmaterialpoints) ifort version 18.0.1 crashes (for whatever reason) end associate diff --git a/src/homogenization_mech_isostrain.f90 b/src/homogenization_mech_isostrain.f90 index f064578c8..994c1b410 100644 --- a/src/homogenization_mech_isostrain.f90 +++ b/src/homogenization_mech_isostrain.f90 @@ -18,7 +18,7 @@ submodule(homogenization) homogenization_mech_isostrain mapping end type - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -29,9 +29,9 @@ contains module subroutine mech_isostrain_init integer :: & - Ninstance, & + Ninstances, & h, & - NofMyHomog + Nmaterialpoints class(tNode), pointer :: & material_homogenization, & homog, & @@ -39,10 +39,10 @@ module subroutine mech_isostrain_init print'(/,a)', ' <<<+- homogenization_mech_isostrain init -+>>>' - Ninstance = count(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) + Ninstances = count(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) - allocate(param(Ninstance)) ! one container of parameters per instance + allocate(param(Ninstances)) ! one container of parameters per instance material_homogenization => config_material%get('homogenization') do h = 1, size(homogenization_type) @@ -51,7 +51,7 @@ module subroutine mech_isostrain_init homogMech => homog%get('mech') associate(prm => param(homogenization_typeInstance(h))) - prm%N_constituents = homogenization_Ngrains(h) + prm%N_constituents = homogenization_Nconstituents(h) select case(homogMech%get_asString('mapping',defaultVal = 'sum')) case ('sum') prm%mapping = parallel_ID @@ -61,11 +61,11 @@ module subroutine mech_isostrain_init call IO_error(211,ext_msg='sum'//' (mech_isostrain)') end select - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) homogState(h)%sizeState = 0 - allocate(homogState(h)%state0 (0,NofMyHomog)) - allocate(homogState(h)%subState0(0,NofMyHomog)) - allocate(homogState(h)%state (0,NofMyHomog)) + allocate(homogState(h)%state0 (0,Nmaterialpoints)) + allocate(homogState(h)%subState0(0,Nmaterialpoints)) + allocate(homogState(h)%state (0,Nmaterialpoints)) end associate diff --git a/src/homogenization_mech_none.f90 b/src/homogenization_mech_none.f90 index a58147c45..5b12247cd 100644 --- a/src/homogenization_mech_none.f90 +++ b/src/homogenization_mech_none.f90 @@ -14,26 +14,26 @@ contains module subroutine mech_none_init integer :: & - Ninstance, & + Ninstances, & h, & - NofMyHomog + Nmaterialpoints print'(/,a)', ' <<<+- homogenization_mech_none init -+>>>' - Ninstance = count(homogenization_type == HOMOGENIZATION_NONE_ID) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) + Ninstances = count(homogenization_type == HOMOGENIZATION_NONE_ID) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) do h = 1, size(homogenization_type) if(homogenization_type(h) /= HOMOGENIZATION_NONE_ID) cycle - if(homogenization_Ngrains(h) /= 1) & + if(homogenization_Nconstituents(h) /= 1) & call IO_error(211,ext_msg='N_constituents (mech_none)') - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) homogState(h)%sizeState = 0 - allocate(homogState(h)%state0 (0,NofMyHomog)) - allocate(homogState(h)%subState0(0,NofMyHomog)) - allocate(homogState(h)%state (0,NofMyHomog)) + allocate(homogState(h)%state0 (0,Nmaterialpoints)) + allocate(homogState(h)%subState0(0,Nmaterialpoints)) + allocate(homogState(h)%state (0,Nmaterialpoints)) enddo diff --git a/src/kinematics_cleavage_opening.f90 b/src/kinematics_cleavage_opening.f90 index d52fdbc1c..44bbcb4f4 100644 --- a/src/kinematics_cleavage_opening.f90 +++ b/src/kinematics_cleavage_opening.f90 @@ -20,7 +20,7 @@ submodule(constitutive:constitutive_damage) kinematics_cleavage_opening cleavage_systems end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -35,7 +35,7 @@ module function kinematics_cleavage_opening_init(kinematics_length) result(myKin integer, intent(in) :: kinematics_length logical, dimension(:,:), allocatable :: myKinematics - integer :: Ninstance,p,k + integer :: Ninstances,p,k integer, dimension(:), allocatable :: N_cl !< active number of cleavage systems per family character(len=pStringLen) :: extmsg = '' class(tNode), pointer :: & @@ -48,12 +48,12 @@ module function kinematics_cleavage_opening_init(kinematics_length) result(myKin print'(/,a)', ' <<<+- kinematics_cleavage_opening init -+>>>' myKinematics = kinematics_active('cleavage_opening',kinematics_length) - Ninstance = count(myKinematics) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myKinematics) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(kinematics_cleavage_opening_instance(phases%length), source=0) do p = 1, phases%length diff --git a/src/kinematics_slipplane_opening.f90 b/src/kinematics_slipplane_opening.f90 index e0de5e181..ea8f51427 100644 --- a/src/kinematics_slipplane_opening.f90 +++ b/src/kinematics_slipplane_opening.f90 @@ -22,7 +22,7 @@ submodule(constitutive:constitutive_damage) kinematics_slipplane_opening P_n end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -37,7 +37,7 @@ module function kinematics_slipplane_opening_init(kinematics_length) result(myKi integer, intent(in) :: kinematics_length logical, dimension(:,:), allocatable :: myKinematics - integer :: Ninstance,p,i,k + integer :: Ninstances,p,i,k character(len=pStringLen) :: extmsg = '' integer, dimension(:), allocatable :: N_sl real(pReal), dimension(:,:), allocatable :: d,n,t @@ -51,13 +51,13 @@ module function kinematics_slipplane_opening_init(kinematics_length) result(myKi print'(/,a)', ' <<<+- kinematics_slipplane init -+>>>' myKinematics = kinematics_active('slipplane_opening',kinematics_length) - Ninstance = count(myKinematics) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myKinematics) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') allocate(kinematics_slipplane_opening_instance(phases%length), source=0) - allocate(param(Ninstance)) + allocate(param(Ninstances)) do p = 1, phases%length if(any(myKinematics(:,p))) kinematics_slipplane_opening_instance(p) = count(myKinematics(:,1:p)) diff --git a/src/kinematics_thermal_expansion.f90 b/src/kinematics_thermal_expansion.f90 index 772f5abbf..4afde82e3 100644 --- a/src/kinematics_thermal_expansion.f90 +++ b/src/kinematics_thermal_expansion.f90 @@ -29,7 +29,7 @@ module function kinematics_thermal_expansion_init(kinematics_length) result(myKi integer, intent(in) :: kinematics_length logical, dimension(:,:), allocatable :: myKinematics - integer :: Ninstance,p,i,k + integer :: Ninstances,p,i,k real(pReal), dimension(:), allocatable :: temp class(tNode), pointer :: & phases, & @@ -41,12 +41,12 @@ module function kinematics_thermal_expansion_init(kinematics_length) result(myKi print'(/,a)', ' <<<+- kinematics_thermal_expansion init -+>>>' myKinematics = kinematics_active('thermal_expansion',kinematics_length) - Ninstance = count(myKinematics) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myKinematics) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(kinematics_thermal_expansion_instance(phases%length), source=0) do p = 1, phases%length diff --git a/src/material.f90 b/src/material.f90 index 11dfeb42e..a0c81a519 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -52,7 +52,7 @@ module material HOMOGENIZATION_RGC_ID end enum - character(len=pStringLen), public, protected, allocatable, dimension(:) :: & + character(len=pStringLen), public, protected, allocatable, dimension(:) :: & material_name_phase, & !< name of each phase material_name_homogenization !< name of each homogenization @@ -64,13 +64,10 @@ module material homogenization_type !< type of each homogenization integer, public, protected :: & - material_Nhomogenization !< number of homogenizations - - integer, public, protected :: & - homogenization_maxNgrains !< max number of grains in any USED homogenization + homogenization_maxNconstituents !< max number of grains in any USED homogenization integer, dimension(:), allocatable, public, protected :: & - homogenization_Ngrains, & !< number of grains in each homogenization + homogenization_Nconstituents, & !< number of grains in each homogenization homogenization_typeInstance, & !< instance of particular type of each homogenization thermal_typeInstance, & !< instance of particular type of each thermal transport damage_typeInstance !< instance of particular type of each nonlocal damage @@ -83,9 +80,9 @@ module material material_homogenizationAt !< homogenization ID of each element integer, dimension(:,:), allocatable, public, target :: & ! (ip,elem) ToDo: ugly target for mapping hack material_homogenizationMemberAt !< position of the element within its homogenization instance - integer, dimension(:,:), allocatable, public, protected :: & ! (constituent,elem) + integer, dimension(:,:), allocatable, public, protected :: & ! (constituent,elem) material_phaseAt !< phase ID of each element - integer, dimension(:,:,:), allocatable, public, protected :: & ! (constituent,elem) + integer, dimension(:,:,:), allocatable, public, protected :: & ! (constituent,IP,elem) material_phaseMemberAt !< position of the element within its phase instance type(tState), allocatable, dimension(:), public :: & @@ -96,11 +93,6 @@ module material type(Rotation), dimension(:,:,:), allocatable, public, protected :: & material_orientation0 !< initial orientation of each grain,IP,element - integer, dimension(:), allocatable, private :: & - material_Nconstituents !< number of constituents in each material - - - ! BEGIN DEPRECATED integer, dimension(:,:), allocatable, private, target :: mappingHomogenizationConst !< mapping from material points to offset in constant state/field ! END DEPRECATED @@ -157,28 +149,10 @@ contains subroutine material_init(restart) logical, intent(in) :: restart - - integer :: ph, myHomog - class(tNode), pointer :: & - phases, & - material_homogenization - character(len=pStringLen) :: sectionName + integer :: myHomog print'(/,a)', ' <<<+- material init -+>>>'; flush(IO_STDOUT) - phases => config_material%get('phase') - allocate(material_name_phase(phases%length)) - do ph = 1, phases%length - write(sectionName,'(i0,a)') ph,'_' - material_name_phase(ph) = trim(adjustl(sectionName))//phases%getKey(ph) !ToDO: No reason to do. Update damage tests - enddo - - material_homogenization => config_material%get('homogenization') - allocate(material_name_homogenization(material_homogenization%length)) - do myHomog = 1, material_homogenization%length - write(sectionName,'(i0,a)') myHomog,'_' - material_name_homogenization(myHomog) = trim(adjustl(sectionName))//material_homogenization%getKey(myHomog) - enddo call material_parseMaterial print*, 'Material parsed' @@ -187,34 +161,32 @@ subroutine material_init(restart) print*, 'Homogenization parsed' - if(homogenization_maxNgrains > size(material_phaseAt,1)) call IO_error(148) + allocate(homogState (size(material_name_homogenization))) + allocate(thermalState (size(material_name_homogenization))) + allocate(damageState (size(material_name_homogenization))) - allocate(homogState (material_Nhomogenization)) - allocate(thermalState (material_Nhomogenization)) - allocate(damageState (material_Nhomogenization)) + allocate(thermalMapping (size(material_name_homogenization))) + allocate(damageMapping (size(material_name_homogenization))) - allocate(thermalMapping (material_Nhomogenization)) - allocate(damageMapping (material_Nhomogenization)) + allocate(temperature (size(material_name_homogenization))) + allocate(damage (size(material_name_homogenization))) - allocate(temperature (material_Nhomogenization)) - allocate(damage (material_Nhomogenization)) - - allocate(temperatureRate (material_Nhomogenization)) + allocate(temperatureRate (size(material_name_homogenization))) if (.not. restart) then call results_openJobFile call results_mapping_constituent(material_phaseAt,material_phaseMemberAt,material_name_phase) - call results_mapping_materialpoint(material_homogenizationAt,material_homogenizationMemberAt,material_name_homogenization) + call results_mapping_homogenization(material_homogenizationAt,material_homogenizationMemberAt,material_name_homogenization) call results_closeJobFile endif !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! BEGIN DEPRECATED - allocate(mappingHomogenizationConst( discretization_nIP,discretization_nElem),source=1) + allocate(mappingHomogenizationConst( discretization_nIPs,discretization_Nelems),source=1) ! hack needed to initialize field values used during constitutive initialization - do myHomog = 1,material_Nhomogenization + do myHomog = 1, size(material_name_homogenization) thermalMapping (myHomog)%p => mappingHomogenizationConst damageMapping (myHomog)%p => mappingHomogenizationConst allocate(temperature (myHomog)%p(1), source=thermal_initialT(myHomog)) @@ -225,6 +197,7 @@ subroutine material_init(restart) end subroutine material_init + !-------------------------------------------------------------------------------------------------- !> @brief parses the homogenization part from the material configuration ! ToDo: This should be done in homogenization @@ -241,22 +214,19 @@ subroutine material_parseHomogenization integer :: h material_homogenization => config_material%get('homogenization') - material_Nhomogenization = material_homogenization%length - allocate(homogenization_type(material_Nhomogenization), source=HOMOGENIZATION_undefined_ID) - allocate(thermal_type(material_Nhomogenization), source=THERMAL_isothermal_ID) - allocate(damage_type (material_Nhomogenization), source=DAMAGE_none_ID) - allocate(homogenization_typeInstance(material_Nhomogenization), source=0) - allocate(thermal_typeInstance(material_Nhomogenization), source=0) - allocate(damage_typeInstance(material_Nhomogenization), source=0) - allocate(homogenization_Ngrains(material_Nhomogenization), source=0) - allocate(thermal_initialT(material_Nhomogenization), source=300.0_pReal) - allocate(damage_initialPhi(material_Nhomogenization), source=1.0_pReal) + allocate(homogenization_type(size(material_name_homogenization)), source=HOMOGENIZATION_undefined_ID) + allocate(thermal_type(size(material_name_homogenization)), source=THERMAL_isothermal_ID) + allocate(damage_type (size(material_name_homogenization)), source=DAMAGE_none_ID) + allocate(homogenization_typeInstance(size(material_name_homogenization)), source=0) + allocate(thermal_typeInstance(size(material_name_homogenization)), source=0) + allocate(damage_typeInstance(size(material_name_homogenization)), source=0) + allocate(thermal_initialT(size(material_name_homogenization)), source=300.0_pReal) + allocate(damage_initialPhi(size(material_name_homogenization)), source=1.0_pReal) - do h=1, material_Nhomogenization + do h=1, size(material_name_homogenization) homog => material_homogenization%get(h) homogMech => homog%get('mech') - homogenization_Ngrains(h) = homog%get_asInt('N_constituents') select case (homogMech%get_asString('type')) case('none') homogenization_type(h) = HOMOGENIZATION_NONE_ID @@ -302,15 +272,12 @@ subroutine material_parseHomogenization endif enddo - do h=1, material_Nhomogenization + do h=1, size(material_name_homogenization) homogenization_typeInstance(h) = count(homogenization_type(1:h) == homogenization_type(h)) thermal_typeInstance(h) = count(thermal_type (1:h) == thermal_type (h)) damage_typeInstance(h) = count(damage_type (1:h) == damage_type (h)) enddo - homogenization_maxNgrains = maxval(homogenization_Ngrains) - - end subroutine material_parseHomogenization @@ -324,7 +291,8 @@ subroutine material_parseMaterial constituents, & !> list of constituents constituent, & !> constituent definition phases, & - homogenizations + homogenizations, & + homogenization integer, dimension(:), allocatable :: & counterPhase, & @@ -333,65 +301,108 @@ subroutine material_parseMaterial real(pReal) :: & frac integer :: & - e, & - i, & - m, & - c, & - maxNconstituents + e, i, c, & + h - materials => config_material%get('material') - if(any(discretization_materialAt > materials%length)) & - call IO_error(155,ext_msg='More materials requested than found in material.yaml') - - allocate(material_Nconstituents(materials%length),source=0) - do m = 1, materials%length - material => materials%get(m) - constituents => material%get('constituents') - material_Nconstituents(m) = constituents%length - enddo - maxNconstituents = maxval(material_Nconstituents) - - allocate(material_homogenizationAt(discretization_nElem),source=0) - allocate(material_homogenizationMemberAt(discretization_nIP,discretization_nElem),source=0) - allocate(material_phaseAt(maxNconstituents,discretization_nElem),source=0) - allocate(material_phaseMemberAt(maxNconstituents,discretization_nIP,discretization_nElem),source=0) - - allocate(material_orientation0(maxNconstituents,discretization_nIP,discretization_nElem)) - - phases => config_material%get('phase') - allocate(counterPhase(phases%length),source=0) + materials => config_material%get('material') + phases => config_material%get('phase') homogenizations => config_material%get('homogenization') + + call sanityCheck(materials, homogenizations) + material_name_phase = getKeys(phases) + material_name_homogenization = getKeys(homogenizations) + + allocate(homogenization_Nconstituents(homogenizations%length)) + do h=1, homogenizations%length + homogenization => homogenizations%get(h) + homogenization_Nconstituents(h) = homogenization%get_asInt('N_constituents') + enddo + homogenization_maxNconstituents = maxval(homogenization_Nconstituents) + + allocate(counterPhase(phases%length),source=0) allocate(counterHomogenization(homogenizations%length),source=0) - do e = 1, discretization_nElem + allocate(material_homogenizationAt(discretization_Nelems),source=0) + allocate(material_homogenizationMemberAt(discretization_nIPs,discretization_Nelems),source=0) + allocate(material_phaseAt(homogenization_maxNconstituents,discretization_Nelems),source=0) + allocate(material_phaseMemberAt(homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems),source=0) + + allocate(material_orientation0(homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems)) + + do e = 1, discretization_Nelems material => materials%get(discretization_materialAt(e)) constituents => material%get('constituents') - + material_homogenizationAt(e) = homogenizations%getIndex(material%get_asString('homogenization')) - do i = 1, discretization_nIP + do i = 1, discretization_nIPs counterHomogenization(material_homogenizationAt(e)) = counterHomogenization(material_homogenizationAt(e)) + 1 material_homogenizationMemberAt(i,e) = counterHomogenization(material_homogenizationAt(e)) enddo - + frac = 0.0_pReal do c = 1, constituents%length constituent => constituents%get(c) frac = frac + constituent%get_asFloat('fraction') - + material_phaseAt(c,e) = phases%getIndex(constituent%get_asString('phase')) - do i = 1, discretization_nIP + do i = 1, discretization_nIPs counterPhase(material_phaseAt(c,e)) = counterPhase(material_phaseAt(c,e)) + 1 material_phaseMemberAt(c,i,e) = counterPhase(material_phaseAt(c,e)) - - call material_orientation0(c,i,e)%fromQuaternion(constituent%get_asFloats('O',requiredSize=4)) + + call material_orientation0(c,i,e)%fromQuaternion(constituent%get_asFloats('O',requiredSize=4)) ! should be done in crystallite enddo - + enddo if (dNeq(frac,1.0_pReal)) call IO_error(153,ext_msg='constituent') - + enddo end subroutine material_parseMaterial +!-------------------------------------------------------------------------------------------------- +!> @brief Check if material.yaml is consistent and contains sufficient # of materials +!-------------------------------------------------------------------------------------------------- +subroutine sanityCheck(materials,homogenizations) + + class(tNode), intent(in) :: materials, & + homogenizations + + class(tNode), pointer :: material, & + homogenization, & + constituents + integer :: m + + if(maxval(discretization_materialAt) > materials%length) & + call IO_error(155,ext_msg='More materials requested than found in material.yaml') + + do m = 1, materials%length + material => materials%get(m) + constituents => material%get('constituents') + homogenization => homogenizations%get(material%get_asString('homogenization')) + if(constituents%length /= homogenization%get_asInt('N_constituents')) call IO_error(148) + enddo + +end subroutine sanityCheck + + +!-------------------------------------------------------------------------------------------------- +!> @brief Get all keys from a dictionary (currently with #_ prefix) +!-------------------------------------------------------------------------------------------------- +function getKeys(dict) + + class(tNode), intent(in) :: dict + character(len=pStringLen), dimension(:), allocatable :: getKeys + + integer :: i + character(len=pStringLen) :: sectionName + + allocate(getKeys(dict%length)) + do i=1, dict%length + write(sectionName,'(i0,a)') i,'_' + getKeys(i) = trim(adjustl(sectionName))//dict%getKey(i) !ToDo: remove prefix + enddo + +end function getKeys + end module material diff --git a/src/mesh/FEM_utilities.f90 b/src/mesh/FEM_utilities.f90 index 4927d0c1c..4c958ee2e 100644 --- a/src/mesh/FEM_utilities.f90 +++ b/src/mesh/FEM_utilities.f90 @@ -164,7 +164,7 @@ subroutine utilities_constitutiveResponse(timeinc,P_av,forwardData) cutBack = .false. ! reset cutBack status - P_av = sum(sum(materialpoint_P,dim=4),dim=3) * wgt ! average of P + P_av = sum(sum(homogenization_P,dim=4),dim=3) * wgt ! average of P call MPI_Allreduce(MPI_IN_PLACE,P_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) end subroutine utilities_constitutiveResponse diff --git a/src/mesh/mesh_mech_FEM.f90 b/src/mesh/mesh_mech_FEM.f90 index f7d33adcf..8aa084ac8 100644 --- a/src/mesh/mesh_mech_FEM.f90 +++ b/src/mesh/mesh_mech_FEM.f90 @@ -400,15 +400,15 @@ subroutine FEM_mech_formResidual(dm_local,xx_local,f_local,dummy,ierr) (((qPt*nBasis + basis)*dimPlex + comp)*dimPlex+comp+1)*dimPlex)) enddo enddo - materialpoint_F(1:dimPlex,1:dimPlex,qPt+1,cell+1) = & + homogenization_F(1:dimPlex,1:dimPlex,qPt+1,cell+1) = & reshape(matmul(BMat,x_scal),shape=[dimPlex,dimPlex], order=[2,1]) enddo if (num%BBarStabilisation) then - detFAvg = math_det33(sum(materialpoint_F(1:3,1:3,1:nQuadrature,cell+1),dim=3)/real(nQuadrature)) + detFAvg = math_det33(sum(homogenization_F(1:3,1:3,1:nQuadrature,cell+1),dim=3)/real(nQuadrature)) do qPt = 1, nQuadrature - materialpoint_F(1:dimPlex,1:dimPlex,qPt,cell+1) = & - materialpoint_F(1:dimPlex,1:dimPlex,qPt,cell+1)* & - (detFAvg/math_det33(materialpoint_F(1:3,1:3,qPt,cell+1)))**(1.0/real(dimPlex)) + homogenization_F(1:dimPlex,1:dimPlex,qPt,cell+1) = & + homogenization_F(1:dimPlex,1:dimPlex,qPt,cell+1)* & + (detFAvg/math_det33(homogenization_F(1:3,1:3,qPt,cell+1)))**(1.0/real(dimPlex)) enddo endif @@ -443,7 +443,7 @@ subroutine FEM_mech_formResidual(dm_local,xx_local,f_local,dummy,ierr) enddo f_scal = f_scal + & matmul(transpose(BMat), & - reshape(transpose(materialpoint_P(1:dimPlex,1:dimPlex,qPt+1,cell+1)), & + reshape(transpose(homogenization_P(1:dimPlex,1:dimPlex,qPt+1,cell+1)), & shape=[dimPlex*dimPlex]))*qWeights(qPt+1) enddo f_scal = f_scal*abs(detJ) @@ -545,7 +545,7 @@ subroutine FEM_mech_formJacobian(dm_local,xx_local,Jac_pre,Jac,dummy,ierr) (((qPt*nBasis + basis)*dimPlex + comp)*dimPlex+comp+1)*dimPlex)) enddo enddo - MatA = matmul(reshape(reshape(materialpoint_dPdF(1:dimPlex,1:dimPlex,1:dimPlex,1:dimPlex,qPt+1,cell+1), & + MatA = matmul(reshape(reshape(homogenization_dPdF(1:dimPlex,1:dimPlex,1:dimPlex,1:dimPlex,qPt+1,cell+1), & shape=[dimPlex,dimPlex,dimPlex,dimPlex], order=[2,1,4,3]), & shape=[dimPlex*dimPlex,dimPlex*dimPlex]),BMat)*qWeights(qPt+1) if (num%BBarStabilisation) then @@ -553,12 +553,12 @@ subroutine FEM_mech_formJacobian(dm_local,xx_local,Jac_pre,Jac,dummy,ierr) FInv = math_inv33(F) K_eA = K_eA + matmul(transpose(BMat),MatA)*math_det33(FInv)**(1.0/real(dimPlex)) K_eB = K_eB - & - matmul(transpose(matmul(reshape(materialpoint_F(1:dimPlex,1:dimPlex,qPt+1,cell+1), & + matmul(transpose(matmul(reshape(homogenization_F(1:dimPlex,1:dimPlex,qPt+1,cell+1), & shape=[dimPlex*dimPlex,1]), & matmul(reshape(FInv(1:dimPlex,1:dimPlex), & shape=[1,dimPlex*dimPlex],order=[2,1]),BMat))),MatA) MatB = MatB + & - matmul(reshape(materialpoint_F(1:dimPlex,1:dimPlex,qPt+1,cell+1),shape=[1,dimPlex*dimPlex]),MatA) + matmul(reshape(homogenization_F(1:dimPlex,1:dimPlex,qPt+1,cell+1),shape=[1,dimPlex*dimPlex]),MatA) FAvg = FAvg + F BMatAvg = BMatAvg + BMat else @@ -630,7 +630,7 @@ subroutine FEM_mech_forward(guess,timeinc,timeinc_old,fieldBC) ! forward last inc if (guess .and. .not. cutBack) then ForwardData = .True. - materialpoint_F0 = materialpoint_F + homogenization_F0 = homogenization_F call SNESGetDM(mech_snes,dm_local,ierr); CHKERRQ(ierr) !< retrieve mesh info from mech_snes into dm_local call DMGetSection(dm_local,section,ierr); CHKERRQ(ierr) call DMGetLocalVector(dm_local,x_local,ierr); CHKERRQ(ierr) diff --git a/src/results.f90 b/src/results.f90 index aec90d7be..1ccc6bfab 100644 --- a/src/results.f90 +++ b/src/results.f90 @@ -56,7 +56,7 @@ module results results_addAttribute, & results_removeLink, & results_mapping_constituent, & - results_mapping_materialpoint + results_mapping_homogenization contains subroutine results_init(restart) @@ -644,7 +644,7 @@ end subroutine results_mapping_constituent !-------------------------------------------------------------------------------------------------- !> @brief adds the unique mapping from spatial position and constituent ID to results !-------------------------------------------------------------------------------------------------- -subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) +subroutine results_mapping_homogenization(homogenizationAt,memberAtLocal,label) integer, dimension(:), intent(in) :: homogenizationAt !< homogenization section at (element) integer, dimension(:,:), intent(in) :: memberAtLocal !< homogenization member at (IP,element) @@ -711,13 +711,13 @@ subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) ! MPI settings and communication #ifdef PETSc call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, ierr) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5pset_dxpl_mpio_f') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5pset_dxpl_mpio_f') call MPI_allreduce(MPI_IN_PLACE,writeSize,worldsize,MPI_INT,MPI_SUM,PETSC_COMM_WORLD,ierr) ! get output at each process - if (ierr /= 0) call IO_error(894,ext_msg='results_mapping_materialpoint: MPI_allreduce/writeSize') + if (ierr /= 0) call IO_error(894,ext_msg='results_mapping_homogenization: MPI_allreduce/writeSize') call MPI_allreduce(MPI_IN_PLACE,memberOffset,size(memberOffset),MPI_INT,MPI_SUM,PETSC_COMM_WORLD,ierr)! get offset at each process - if (ierr /= 0) call IO_error(894,ext_msg='results_mapping_materialpoint: MPI_allreduce/memberOffset') + if (ierr /= 0) call IO_error(894,ext_msg='results_mapping_homogenization: MPI_allreduce/memberOffset') #endif myShape = int([writeSize(worldrank)], HSIZE_T) @@ -727,13 +727,13 @@ subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) !-------------------------------------------------------------------------------------------------- ! create dataspace in memory (local shape = hyperslab) and in file (global shape) call h5screate_simple_f(1,myShape,memspace_id,ierr,myShape) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5screate_simple_f/memspace_id') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5screate_simple_f/memspace_id') call h5screate_simple_f(1,totalShape,filespace_id,ierr,totalShape) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5screate_simple_f/filespace_id') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5screate_simple_f/filespace_id') call h5sselect_hyperslab_f(filespace_id, H5S_SELECT_SET_F, myOffset, myShape, ierr) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5sselect_hyperslab_f') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5sselect_hyperslab_f') !--------------------------------------------------------------------------------------------------- ! expand phaseAt to consider IPs (is not stored per IP) @@ -753,14 +753,14 @@ subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) loc_id = results_openGroup('/mapping') call h5dcreate_f(loc_id, 'homogenization', dtype_id, filespace_id, dset_id, ierr) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5dcreate_f') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5dcreate_f') call h5dwrite_f(dset_id, name_id, reshape(label(pack(homogenizationAtMaterialpoint,.true.)),myShape), & myShape, ierr, file_space_id = filespace_id, mem_space_id = memspace_id, xfer_prp = plist_id) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5dwrite_f/name_id') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5dwrite_f/name_id') call h5dwrite_f(dset_id, position_id, reshape(pack(memberAtGlobal,.true.),myShape), & myShape, ierr, file_space_id = filespace_id, mem_space_id = memspace_id, xfer_prp = plist_id) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5dwrite_f/position_id') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5dwrite_f/position_id') !-------------------------------------------------------------------------------------------------- ! close all @@ -776,7 +776,7 @@ subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) ! for backward compatibility call results_setLink('/mapping/homogenization','/mapping/cellResults/materialpoint') -end subroutine results_mapping_materialpoint +end subroutine results_mapping_homogenization !-------------------------------------------------------------------------------------------------- diff --git a/src/source_damage_anisoBrittle.f90 b/src/source_damage_anisoBrittle.f90 index 1954e6ea3..ca8d6ec2b 100644 --- a/src/source_damage_anisoBrittle.f90 +++ b/src/source_damage_anisoBrittle.f90 @@ -25,7 +25,7 @@ submodule (constitutive:constitutive_damage) source_damage_anisoBrittle output end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -45,19 +45,19 @@ module function source_damage_anisoBrittle_init(source_length) result(mySources) phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p integer, dimension(:), allocatable :: N_cl character(len=pStringLen) :: extmsg = '' print'(/,a)', ' <<<+- source_damage_anisoBrittle init -+>>>' mySources = source_active('damage_anisoBrittle',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_damage_anisoBrittle_offset (phases%length), source=0) allocate(source_damage_anisoBrittle_instance(phases%length), source=0) @@ -100,8 +100,8 @@ module function source_damage_anisoBrittle_init(source_length) 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' - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) + Nconstituents = count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,0) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('anisobrittle_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisobrittle_atol' diff --git a/src/source_damage_anisoDuctile.f90 b/src/source_damage_anisoDuctile.f90 index 3b8eb2428..2fdd46bb4 100644 --- a/src/source_damage_anisoDuctile.f90 +++ b/src/source_damage_anisoDuctile.f90 @@ -19,7 +19,7 @@ submodule(constitutive:constitutive_damage) source_damage_anisoDuctile output end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -39,19 +39,19 @@ module function source_damage_anisoDuctile_init(source_length) result(mySources) pl, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p integer, dimension(:), allocatable :: N_sl character(len=pStringLen) :: extmsg = '' print'(/,a)', ' <<<+- source_damage_anisoDuctile init -+>>>' mySources = source_active('damage_anisoDuctile',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_damage_anisoDuctile_offset (phases%length), source=0) allocate(source_damage_anisoDuctile_instance(phases%length), source=0) @@ -84,8 +84,8 @@ module function source_damage_anisoDuctile_init(source_length) result(mySources) if (prm%q <= 0.0_pReal) extmsg = trim(extmsg)//' q' if (any(prm%gamma_crit < 0.0_pReal)) extmsg = trim(extmsg)//' gamma_crit' - NipcMyPhase=count(material_phaseAt==p) * discretization_nIP - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) + Nconstituents=count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,0) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('anisoDuctile_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisoductile_atol' diff --git a/src/source_damage_isoBrittle.f90 b/src/source_damage_isoBrittle.f90 index 714e71ef1..7fcf17ee0 100644 --- a/src/source_damage_isoBrittle.f90 +++ b/src/source_damage_isoBrittle.f90 @@ -17,7 +17,7 @@ submodule(constitutive:constitutive_damage) source_damage_isoBrittle output end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -36,18 +36,18 @@ module function source_damage_isoBrittle_init(source_length) result(mySources) phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p character(len=pStringLen) :: extmsg = '' print'(/,a)', ' <<<+- source_damage_isoBrittle init -+>>>' mySources = source_active('damage_isoBrittle',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_damage_isoBrittle_offset (phases%length), source=0) allocate(source_damage_isoBrittle_instance(phases%length), source=0) @@ -73,8 +73,8 @@ module function source_damage_isoBrittle_init(source_length) result(mySources) ! sanity checks if (prm%W_crit <= 0.0_pReal) extmsg = trim(extmsg)//' W_crit' - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,1) + Nconstituents = count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,1) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('isoBrittle_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isobrittle_atol' diff --git a/src/source_damage_isoDuctile.f90 b/src/source_damage_isoDuctile.f90 index d958aed6a..1bff20570 100644 --- a/src/source_damage_isoDuctile.f90 +++ b/src/source_damage_isoDuctile.f90 @@ -18,7 +18,7 @@ submodule (constitutive:constitutive_damage) source_damage_isoDuctile output end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -38,18 +38,18 @@ module function source_damage_isoDuctile_init(source_length) result(mySources) phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p character(len=pStringLen) :: extmsg = '' print'(/,a)', ' <<<+- source_damage_isoDuctile init -+>>>' mySources = source_active('damage_isoDuctile',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_damage_isoDuctile_offset (phases%length), source=0) allocate(source_damage_isoDuctile_instance(phases%length), source=0) @@ -77,8 +77,8 @@ module function source_damage_isoDuctile_init(source_length) result(mySources) if (prm%q <= 0.0_pReal) extmsg = trim(extmsg)//' q' if (prm%gamma_crit <= 0.0_pReal) extmsg = trim(extmsg)//' gamma_crit' - NipcMyPhase=count(material_phaseAt==p) * discretization_nIP - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) + Nconstituents=count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,0) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('isoDuctile_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isoductile_atol' diff --git a/src/source_thermal_dissipation.f90 b/src/source_thermal_dissipation.f90 index 5cc740424..f28567aa7 100644 --- a/src/source_thermal_dissipation.f90 +++ b/src/source_thermal_dissipation.f90 @@ -15,7 +15,7 @@ submodule(constitutive:constitutive_thermal) source_thermal_dissipation kappa !< TAYLOR-QUINNEY factor end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -35,17 +35,17 @@ module function source_thermal_dissipation_init(source_length) result(mySources) phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p print'(/,a)', ' <<<+- source_thermal_dissipation init -+>>>' mySources = source_active('thermal_dissipation',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_thermal_dissipation_offset (phases%length), source=0) allocate(source_thermal_dissipation_instance(phases%length), source=0) @@ -61,8 +61,8 @@ module function source_thermal_dissipation_init(source_length) result(mySources) src => sources%get(sourceOffset) prm%kappa = src%get_asFloat('kappa') - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,0,0,0) + Nconstituents = count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,0,0,0) end associate endif @@ -74,7 +74,7 @@ end function source_thermal_dissipation_init !-------------------------------------------------------------------------------------------------- -!> @brief Ninstances dissipation rate +!> @brief Ninstancess dissipation rate !-------------------------------------------------------------------------------------------------- module subroutine source_thermal_dissipation_getRateAndItsTangent(TDot, dTDot_dT, Tstar, Lp, phase) diff --git a/src/source_thermal_externalheat.f90 b/src/source_thermal_externalheat.f90 index 2eeeb47df..9ba4a051b 100644 --- a/src/source_thermal_externalheat.f90 +++ b/src/source_thermal_externalheat.f90 @@ -19,7 +19,7 @@ submodule(constitutive:constitutive_thermal) source_thermal_externalheat nIntervals end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -39,17 +39,17 @@ module function source_thermal_externalheat_init(source_length) result(mySources phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p print'(/,a)', ' <<<+- source_thermal_externalHeat init -+>>>' mySources = source_active('thermal_externalheat',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_thermal_externalheat_offset (phases%length), source=0) allocate(source_thermal_externalheat_instance(phases%length), source=0) @@ -69,8 +69,8 @@ module function source_thermal_externalheat_init(source_length) result(mySources prm%f_T = src%get_asFloats('f_T',requiredSize = size(prm%t_n)) - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) + Nconstituents = count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,0) end associate endif diff --git a/src/thermal_adiabatic.f90 b/src/thermal_adiabatic.f90 index 63deb3cd5..aa807924c 100644 --- a/src/thermal_adiabatic.f90 +++ b/src/thermal_adiabatic.f90 @@ -40,7 +40,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine thermal_adiabatic_init - integer :: maxNinstance,h,NofMyHomog + integer :: maxNinstances,h,Nmaterialpoints class(tNode), pointer :: & material_homogenization, & homog, & @@ -48,13 +48,13 @@ subroutine thermal_adiabatic_init print'(/,a)', ' <<<+- thermal_adiabatic init -+>>>'; flush(6) - maxNinstance = count(thermal_type == THERMAL_adiabatic_ID) - if (maxNinstance == 0) return + maxNinstances = count(thermal_type == THERMAL_adiabatic_ID) + if (maxNinstances == 0) return - allocate(param(maxNinstance)) + allocate(param(maxNinstances)) material_homogenization => config_material%get('homogenization') - do h = 1, material_Nhomogenization + do h = 1, size(material_name_homogenization) if (thermal_type(h) /= THERMAL_adiabatic_ID) cycle homog => material_homogenization%get(h) homogThermal => homog%get('thermal') @@ -67,17 +67,17 @@ subroutine thermal_adiabatic_init prm%output = homogThermal%get_asStrings('output',defaultVal=emptyStringArray) #endif - NofMyHomog=count(material_homogenizationAt==h) + Nmaterialpoints=count(material_homogenizationAt==h) thermalState(h)%sizeState = 1 - allocate(thermalState(h)%state0 (1,NofMyHomog), source=thermal_initialT(h)) - allocate(thermalState(h)%subState0(1,NofMyHomog), source=thermal_initialT(h)) - allocate(thermalState(h)%state (1,NofMyHomog), source=thermal_initialT(h)) + allocate(thermalState(h)%state0 (1,Nmaterialpoints), source=thermal_initialT(h)) + allocate(thermalState(h)%subState0(1,Nmaterialpoints), source=thermal_initialT(h)) + allocate(thermalState(h)%state (1,Nmaterialpoints), source=thermal_initialT(h)) thermalMapping(h)%p => material_homogenizationMemberAt deallocate(temperature(h)%p) temperature(h)%p => thermalState(h)%state(1,:) deallocate(temperatureRate(h)%p) - allocate (temperatureRate(h)%p(NofMyHomog), source=0.0_pReal) + allocate (temperatureRate(h)%p(Nmaterialpoints), source=0.0_pReal) end associate enddo @@ -145,8 +145,8 @@ subroutine thermal_adiabatic_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el) homog = material_homogenizationAt(el) call constitutive_thermal_getRateAndItsTangents(TDot, dTDot_dT, T, crystallite_S, crystallite_Lp, ip, el) - Tdot = Tdot/real(homogenization_Ngrains(homog),pReal) - dTdot_dT = dTdot_dT/real(homogenization_Ngrains(homog),pReal) + Tdot = Tdot/real(homogenization_Nconstituents(homog),pReal) + dTdot_dT = dTdot_dT/real(homogenization_Nconstituents(homog),pReal) end subroutine thermal_adiabatic_getSourceAndItsTangent @@ -167,13 +167,13 @@ function thermal_adiabatic_getSpecificHeat(ip,el) thermal_adiabatic_getSpecificHeat = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_adiabatic_getSpecificHeat = thermal_adiabatic_getSpecificHeat & + lattice_c_p(material_phaseAt(grain,el)) enddo thermal_adiabatic_getSpecificHeat = thermal_adiabatic_getSpecificHeat & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_adiabatic_getSpecificHeat @@ -193,13 +193,13 @@ function thermal_adiabatic_getMassDensity(ip,el) thermal_adiabatic_getMassDensity = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_adiabatic_getMassDensity = thermal_adiabatic_getMassDensity & + lattice_rho(material_phaseAt(grain,el)) enddo thermal_adiabatic_getMassDensity = thermal_adiabatic_getMassDensity & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_adiabatic_getMassDensity diff --git a/src/thermal_conduction.f90 b/src/thermal_conduction.f90 index 69ce8025e..daa7391a9 100644 --- a/src/thermal_conduction.f90 +++ b/src/thermal_conduction.f90 @@ -41,7 +41,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine thermal_conduction_init - integer :: Ninstance,NofMyHomog,h + integer :: Ninstances,Nmaterialpoints,h class(tNode), pointer :: & material_homogenization, & homog, & @@ -49,11 +49,11 @@ subroutine thermal_conduction_init print'(/,a)', ' <<<+- thermal_conduction init -+>>>'; flush(6) - Ninstance = count(thermal_type == THERMAL_conduction_ID) - allocate(param(Ninstance)) + Ninstances = count(thermal_type == THERMAL_conduction_ID) + allocate(param(Ninstances)) material_homogenization => config_material%get('homogenization') - do h = 1, material_Nhomogenization + do h = 1, size(material_name_homogenization) if (thermal_type(h) /= THERMAL_conduction_ID) cycle homog => material_homogenization%get(h) homogThermal => homog%get('thermal') @@ -65,17 +65,17 @@ subroutine thermal_conduction_init prm%output = homogThermal%get_asStrings('output',defaultVal=emptyStringArray) #endif - NofMyHomog=count(material_homogenizationAt==h) + Nmaterialpoints=count(material_homogenizationAt==h) thermalState(h)%sizeState = 0 - allocate(thermalState(h)%state0 (0,NofMyHomog)) - allocate(thermalState(h)%subState0(0,NofMyHomog)) - allocate(thermalState(h)%state (0,NofMyHomog)) + allocate(thermalState(h)%state0 (0,Nmaterialpoints)) + allocate(thermalState(h)%subState0(0,Nmaterialpoints)) + allocate(thermalState(h)%state (0,Nmaterialpoints)) thermalMapping(h)%p => material_homogenizationMemberAt deallocate(temperature (h)%p) - allocate (temperature (h)%p(NofMyHomog), source=thermal_initialT(h)) + allocate (temperature (h)%p(Nmaterialpoints), source=thermal_initialT(h)) deallocate(temperatureRate(h)%p) - allocate (temperatureRate(h)%p(NofMyHomog), source=0.0_pReal) + allocate (temperatureRate(h)%p(Nmaterialpoints), source=0.0_pReal) end associate enddo @@ -104,8 +104,8 @@ subroutine thermal_conduction_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el) homog = material_homogenizationAt(el) call constitutive_thermal_getRateAndItsTangents(TDot, dTDot_dT, T, crystallite_S,crystallite_Lp ,ip, el) - Tdot = Tdot/real(homogenization_Ngrains(homog),pReal) - dTdot_dT = dTdot_dT/real(homogenization_Ngrains(homog),pReal) + Tdot = Tdot/real(homogenization_Nconstituents(homog),pReal) + dTdot_dT = dTdot_dT/real(homogenization_Nconstituents(homog),pReal) end subroutine thermal_conduction_getSourceAndItsTangent @@ -125,13 +125,13 @@ function thermal_conduction_getConductivity(ip,el) thermal_conduction_getConductivity = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_conduction_getConductivity = thermal_conduction_getConductivity + & crystallite_push33ToRef(grain,ip,el,lattice_K(:,:,material_phaseAt(grain,el))) enddo thermal_conduction_getConductivity = thermal_conduction_getConductivity & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_conduction_getConductivity @@ -151,13 +151,13 @@ function thermal_conduction_getSpecificHeat(ip,el) thermal_conduction_getSpecificHeat = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_conduction_getSpecificHeat = thermal_conduction_getSpecificHeat & + lattice_c_p(material_phaseAt(grain,el)) enddo thermal_conduction_getSpecificHeat = thermal_conduction_getSpecificHeat & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_conduction_getSpecificHeat @@ -178,13 +178,13 @@ function thermal_conduction_getMassDensity(ip,el) thermal_conduction_getMassDensity = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_conduction_getMassDensity = thermal_conduction_getMassDensity & + lattice_rho(material_phaseAt(grain,el)) enddo thermal_conduction_getMassDensity = thermal_conduction_getMassDensity & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_conduction_getMassDensity diff --git a/src/thermal_isothermal.f90 b/src/thermal_isothermal.f90 index 2dda358ac..39c8efe91 100644 --- a/src/thermal_isothermal.f90 +++ b/src/thermal_isothermal.f90 @@ -16,18 +16,18 @@ contains !-------------------------------------------------------------------------------------------------- subroutine thermal_isothermal_init - integer :: h,NofMyHomog + integer :: h,Nmaterialpoints print'(/,a)', ' <<<+- thermal_isothermal init -+>>>'; flush(6) - do h = 1, material_Nhomogenization + do h = 1, size(material_name_homogenization) if (thermal_type(h) /= THERMAL_isothermal_ID) cycle - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) thermalState(h)%sizeState = 0 - allocate(thermalState(h)%state0 (0,NofMyHomog)) - allocate(thermalState(h)%subState0(0,NofMyHomog)) - allocate(thermalState(h)%state (0,NofMyHomog)) + allocate(thermalState(h)%state0 (0,Nmaterialpoints)) + allocate(thermalState(h)%subState0(0,Nmaterialpoints)) + allocate(thermalState(h)%state (0,Nmaterialpoints)) deallocate(temperature (h)%p) allocate (temperature (h)%p(1), source=thermal_initialT(h))