diff --git a/code/CommercialFEM_fileList.f90 b/code/CommercialFEM_fileList.f90 index c762718e8..0142c07f0 100644 --- a/code/CommercialFEM_fileList.f90 +++ b/code/CommercialFEM_fileList.f90 @@ -15,10 +15,10 @@ #include "material.f90" #include "lattice.f90" #include "damage_none.f90" -#include "damage_brittle.f90" -#include "damage_ductile.f90" +#include "damage_isoBrittle.f90" +#include "damage_isoDuctile.f90" #include "damage_gurson.f90" -#include "damage_anisotropic.f90" +#include "damage_anisoBrittle.f90" #include "thermal_isothermal.f90" #include "thermal_adiabatic.f90" #include "vacancy_constant.f90" diff --git a/code/Makefile b/code/Makefile index 1d1a49013..f8a22dd8d 100644 --- a/code/Makefile +++ b/code/Makefile @@ -326,7 +326,7 @@ COMPILE =$(COMPILE_OPTIONS_$(F90)) $(STANDARD_CHECK_$(F90)) $(OPTIMIZATI COMPILE_MAXOPTI =$(COMPILE_OPTIONS_$(F90)) $(STANDARD_CHECK_$(F90)) $(OPTIMIZATION_$(MAXOPTI)_$(F90)) $(INCLUDE_DIRS) $(PRECISION_$(F90)) ################################################################################################### DAMAGE_FILES = \ - damage_none.o damage_brittle.o damage_ductile.o damage_gurson.o damage_anisotropic.o + damage_none.o damage_isoBrittle.o damage_isoDuctile.o damage_gurson.o damage_anisoBrittle.o THERMAL_FILES = \ thermal_isothermal.o thermal_adiabatic.o @@ -481,17 +481,17 @@ constitutive_none.o: constitutive_none.f90 \ damage_none.o: damage_none.f90 \ lattice.o -damage_brittle.o: damage_brittle.f90 \ +damage_isoBrittle.o: damage_isoBrittle.f90 \ lattice.o -damage_ductile.o: damage_ductile.f90 \ +damage_isoDuctile.o: damage_isoDuctile.f90 \ + lattice.o + +damage_anisoBrittle.o: damage_anisoBrittle.f90 \ lattice.o damage_gurson.o: damage_gurson.f90 \ lattice.o - -damage_anisotropic.o: damage_anisotropic.f90 \ - lattice.o thermal_isothermal.o: thermal_isothermal.f90 \ lattice.o diff --git a/code/constitutive.f90 b/code/constitutive.f90 index 26cd5cffc..1112ce77f 100644 --- a/code/constitutive.f90 +++ b/code/constitutive.f90 @@ -123,19 +123,19 @@ subroutine constitutive_init(temperature_init) PLASTICITY_TITANMOD_label, & PLASTICITY_NONLOCAL_label, & LOCAL_DAMAGE_none_ID, & - LOCAL_DAMAGE_brittle_ID, & - LOCAL_DAMAGE_ductile_ID, & + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID, & LOCAL_THERMAL_isothermal_ID, & LOCAL_THERMAL_adiabatic_ID, & LOCAL_VACANCY_constant_ID, & LOCAL_VACANCY_generation_ID, & - LOCAL_DAMAGE_NONE_label, & - LOCAL_DAMAGE_brittle_label, & - LOCAL_DAMAGE_ductile_label, & - LOCAL_DAMAGE_gurson_label, & - LOCAL_DAMAGE_anisotropic_label, & + LOCAL_DAMAGE_none_LABEL, & + LOCAL_DAMAGE_isoBrittle_LABEL, & + LOCAL_DAMAGE_isoDuctile_LABEL, & + LOCAL_DAMAGE_anisoBrittle_LABEL, & + LOCAL_DAMAGE_gurson_LABEL, & LOCAL_THERMAL_isothermal_label, & LOCAL_THERMAL_adiabatic_label, & LOCAL_VACANCY_constant_label, & @@ -155,10 +155,10 @@ subroutine constitutive_init(temperature_init) use constitutive_titanmod use constitutive_nonlocal use damage_none - use damage_brittle - use damage_ductile + use damage_isoBrittle + use damage_isoDuctile + use damage_anisoBrittle use damage_gurson - use damage_anisotropic use thermal_isothermal use thermal_adiabatic use vacancy_constant @@ -199,11 +199,11 @@ subroutine constitutive_init(temperature_init) ! parse damage from config file if (.not. IO_open_jobFile_stat(FILEUNIT,material_localFileExt)) & ! no local material configuration present... call IO_open_file(FILEUNIT,material_configFile) ! ... open material.config file - if (any(phase_damage == LOCAL_DAMAGE_none_ID)) call damage_none_init(FILEUNIT) - if (any(phase_damage == LOCAL_DAMAGE_brittle_ID)) call damage_brittle_init(FILEUNIT) - if (any(phase_damage == LOCAL_DAMAGE_ductile_ID)) call damage_ductile_init(FILEUNIT) - if (any(phase_damage == LOCAL_DAMAGE_gurson_ID)) call damage_gurson_init(FILEUNIT) - if (any(phase_damage == LOCAL_DAMAGE_anisotropic_ID)) call damage_anisotropic_init(FILEUNIT) + if (any(phase_damage == LOCAL_DAMAGE_none_ID)) call damage_none_init(FILEUNIT) + if (any(phase_damage == LOCAL_DAMAGE_isoBrittle_ID)) call damage_isoBrittle_init(FILEUNIT) + if (any(phase_damage == LOCAL_DAMAGE_isoductile_ID)) call damage_isoDuctile_init(FILEUNIT) + if (any(phase_damage == LOCAL_DAMAGE_gurson_ID)) call damage_gurson_init(FILEUNIT) + if (any(phase_damage == LOCAL_DAMAGE_anisoBrittle_ID)) call damage_anisoBrittle_init(FILEUNIT) close(FILEUNIT) !-------------------------------------------------------------------------------------------------- @@ -292,26 +292,26 @@ subroutine constitutive_init(temperature_init) thisNoutput => null() thisOutput => null() thisSize => null() - case (LOCAL_DAMAGE_brittle_ID) - outputName = LOCAL_DAMAGE_BRITTLE_label - thisNoutput => damage_brittle_Noutput - thisOutput => damage_brittle_output - thisSize => damage_brittle_sizePostResult - case (LOCAL_DAMAGE_ductile_ID) - outputName = LOCAL_DAMAGE_DUCTILE_label - thisNoutput => damage_ductile_Noutput - thisOutput => damage_ductile_output - thisSize => damage_ductile_sizePostResult + case (LOCAL_DAMAGE_isoBrittle_ID) + outputName = LOCAL_DAMAGE_isoBrittle_LABEL + thisNoutput => damage_isoBrittle_Noutput + thisOutput => damage_isoBrittle_output + thisSize => damage_isoBrittle_sizePostResult + case (LOCAL_DAMAGE_isoDuctile_ID) + outputName = LOCAL_DAMAGE_isoDuctile_LABEL + thisNoutput => damage_isoDuctile_Noutput + thisOutput => damage_isoDuctile_output + thisSize => damage_isoDuctile_sizePostResult case (LOCAL_DAMAGE_gurson_ID) outputName = LOCAL_DAMAGE_gurson_label thisNoutput => damage_gurson_Noutput thisOutput => damage_gurson_output thisSize => damage_gurson_sizePostResult - case (LOCAL_DAMAGE_anisotropic_ID) - outputName = LOCAL_DAMAGE_anisotropic_label - thisNoutput => damage_anisotropic_Noutput - thisOutput => damage_anisotropic_output - thisSize => damage_anisotropic_sizePostResult + case (LOCAL_DAMAGE_anisoBrittle_ID) + outputName = LOCAL_DAMAGE_anisoBrittle_label + thisNoutput => damage_anisoBrittle_Noutput + thisOutput => damage_anisoBrittle_output + thisSize => damage_anisoBrittle_sizePostResult case default knownDamage = .false. end select @@ -496,10 +496,10 @@ function constitutive_damagedC(ipc,ip,el) pReal use material, only: & material_phase, & - LOCAL_DAMAGE_brittle_ID, & + LOCAL_DAMAGE_isoBrittle_ID, & phase_damage - use damage_brittle, only: & - damage_brittle_getDamage + use damage_isoBrittle, only: & + damage_isoBrittle_getDamage implicit none real(pReal), dimension(6,6) :: constitutive_damagedC @@ -510,8 +510,8 @@ function constitutive_damagedC(ipc,ip,el) real(pReal) :: damage !< element number select case (phase_damage(material_phase(ipc,ip,el))) - case (LOCAL_DAMAGE_brittle_ID) - damage = damage_brittle_getDamage(ipc, ip, el) + case (LOCAL_DAMAGE_isoBrittle_ID) + damage = damage_isoBrittle_getDamage(ipc, ip, el) constitutive_damagedC = damage*damage* & constitutive_homogenizedC(ipc,ip,el) case default @@ -535,10 +535,10 @@ subroutine constitutive_microstructure(Tstar_v, Fe, Fp, ipc, ip, el) PLASTICITY_dislokmc_ID, & PLASTICITY_titanmod_ID, & PLASTICITY_nonlocal_ID, & - LOCAL_DAMAGE_brittle_ID, & - LOCAL_DAMAGE_ductile_ID, & - LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & + LOCAL_DAMAGE_gurson_ID use constitutive_titanmod, only: & constitutive_titanmod_microstructure @@ -548,11 +548,11 @@ subroutine constitutive_microstructure(Tstar_v, Fe, Fp, ipc, ip, el) constitutive_dislotwin_microstructure use constitutive_dislokmc, only: & constitutive_dislokmc_microstructure - use damage_brittle, only: & - damage_brittle_microstructure, & - damage_brittle_getDamage - use damage_ductile, only: & - damage_ductile_microstructure + use damage_isoBrittle, only: & + damage_isoBrittle_microstructure, & + damage_isoBrittle_getDamage + use damage_isoDuctile, only: & + damage_isoDuctile_microstructure use damage_gurson, only: & damage_gurson_microstructure @@ -588,13 +588,13 @@ subroutine constitutive_microstructure(Tstar_v, Fe, Fp, ipc, ip, el) end select select case (phase_damage(material_phase(ipc,ip,el))) - case (LOCAL_DAMAGE_brittle_ID) - damage = damage_brittle_getDamage(ipc,ip,el) + case (LOCAL_DAMAGE_isoBrittle_ID) + damage = damage_isoBrittle_getDamage(ipc,ip,el) Tstar_v_effective = Tstar_v/(damage*damage) - call damage_brittle_microstructure(Tstar_v_effective, Fe, ipc, ip, el) - case (LOCAL_DAMAGE_ductile_ID) + call damage_isoBrittle_microstructure(Tstar_v_effective, Fe, ipc, ip, el) + case (LOCAL_DAMAGE_isoDuctile_ID) call constitutive_getAccumulatedSlip(nSlip,accumulatedSlip,Fp,ipc, ip, el) - call damage_ductile_microstructure(nSlip,accumulatedSlip,ipc, ip, el) + call damage_isoDuctile_microstructure(nSlip,accumulatedSlip,ipc, ip, el) case (LOCAL_DAMAGE_gurson_ID) call damage_gurson_microstructure(ipc, ip, el) @@ -806,10 +806,10 @@ subroutine constitutive_collectDotState(Tstar_v, Lp, FeArray, FpArray, subdt, su PLASTICITY_dislokmc_ID, & PLASTICITY_titanmod_ID, & PLASTICITY_nonlocal_ID, & - LOCAL_DAMAGE_brittle_ID, & - LOCAL_DAMAGE_ductile_ID, & + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID, & LOCAL_THERMAL_adiabatic_ID, & LOCAL_VACANCY_generation_ID use constitutive_j2, only: & @@ -824,14 +824,14 @@ subroutine constitutive_collectDotState(Tstar_v, Lp, FeArray, FpArray, subdt, su constitutive_titanmod_dotState use constitutive_nonlocal, only: & constitutive_nonlocal_dotState - use damage_brittle, only: & - damage_brittle_dotState - use damage_ductile, only: & - damage_ductile_dotState + use damage_isoBrittle, only: & + damage_isoBrittle_dotState + use damage_isoDuctile, only: & + damage_isoDuctile_dotState use damage_gurson, only: & damage_gurson_dotState - use damage_anisotropic, only: & - damage_anisotropic_dotState + use damage_anisoBrittle, only: & + damage_anisoBrittle_dotState use thermal_adiabatic, only: & thermal_adiabatic_dotState use vacancy_generation, only: & @@ -882,14 +882,14 @@ subroutine constitutive_collectDotState(Tstar_v, Lp, FeArray, FpArray, subdt, su end select select case (phase_damage(material_phase(ipc,ip,el))) - case (LOCAL_DAMAGE_brittle_ID) - call damage_brittle_dotState(ipc, ip, el) - case (LOCAL_DAMAGE_ductile_ID) - call damage_ductile_dotState(ipc, ip, el) + case (LOCAL_DAMAGE_isoBrittle_ID) + call damage_isoBrittle_dotState(ipc, ip, el) + case (LOCAL_DAMAGE_isoDuctile_ID) + call damage_isoDuctile_dotState(ipc, ip, el) case (LOCAL_DAMAGE_gurson_ID) call damage_gurson_dotState(Tstar_v, Lp, ipc, ip, el) - case (LOCAL_DAMAGE_anisotropic_ID) - call damage_anisotropic_dotState(Tstar_v, ipc, ip, el) + case (LOCAL_DAMAGE_anisoBrittle_ID) + call damage_anisoBrittle_dotState(Tstar_v, ipc, ip, el) end select select case (phase_thermal(material_phase(ipc,ip,el))) @@ -985,19 +985,19 @@ function constitutive_getLocalDamage(ipc, ip, el) use material, only: & material_phase, & LOCAL_DAMAGE_none_ID, & - LOCAL_DAMAGE_brittle_ID, & - LOCAL_DAMAGE_ductile_ID, & + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID, & phase_damage - use damage_brittle, only: & - damage_brittle_getLocalDamage - use damage_ductile, only: & - damage_ductile_getLocalDamage + use damage_isoBrittle, only: & + damage_isoBrittle_getLocalDamage + use damage_isoDuctile, only: & + damage_isoDuctile_getLocalDamage use damage_gurson, only: & damage_gurson_getLocalDamage - use damage_anisotropic, only: & - damage_anisotropic_getLocalDamage + use damage_anisoBrittle, only: & + damage_anisoBrittle_getLocalDamage implicit none integer(pInt), intent(in) :: & @@ -1010,17 +1010,17 @@ function constitutive_getLocalDamage(ipc, ip, el) case (LOCAL_DAMAGE_none_ID) constitutive_getLocalDamage = 1.0_pReal - case (LOCAL_DAMAGE_brittle_ID) - constitutive_getLocalDamage = damage_brittle_getLocalDamage(ipc, ip, el) + case (LOCAL_DAMAGE_isoBrittle_ID) + constitutive_getLocalDamage = damage_isoBrittle_getLocalDamage(ipc, ip, el) - case (LOCAL_DAMAGE_ductile_ID) - constitutive_getLocalDamage = damage_ductile_getLocalDamage(ipc, ip, el) + case (LOCAL_DAMAGE_isoDuctile_ID) + constitutive_getLocalDamage = damage_isoDuctile_getLocalDamage(ipc, ip, el) case (LOCAL_DAMAGE_gurson_ID) constitutive_getLocalDamage = damage_gurson_getLocalDamage(ipc, ip, el) - case (LOCAL_DAMAGE_anisotropic_ID) - constitutive_getLocalDamage = damage_anisotropic_getLocalDamage(ipc, ip, el) + case (LOCAL_DAMAGE_anisoBrittle_ID) + constitutive_getLocalDamage = damage_anisoBrittle_getLocalDamage(ipc, ip, el) end select end function constitutive_getLocalDamage @@ -1033,19 +1033,19 @@ subroutine constitutive_putLocalDamage(ipc, ip, el, localDamage) pReal use material, only: & material_phase, & - LOCAL_DAMAGE_BRITTLE_ID, & - LOCAL_DAMAGE_DUCTILE_ID, & + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID, & phase_damage - use damage_brittle, only: & - damage_brittle_putLocalDamage - use damage_ductile, only: & - damage_ductile_putLocalDamage + use damage_isoBrittle, only: & + damage_isoBrittle_putLocalDamage + use damage_isoDuctile, only: & + damage_isoDuctile_putLocalDamage use damage_gurson, only: & damage_gurson_putLocalDamage - use damage_anisotropic, only: & - damage_anisotropic_putLocalDamage + use damage_anisoBrittle, only: & + damage_anisoBrittle_putLocalDamage implicit none integer(pInt), intent(in) :: & @@ -1056,17 +1056,17 @@ subroutine constitutive_putLocalDamage(ipc, ip, el, localDamage) localDamage select case (phase_damage(material_phase(ipc,ip,el))) - case (LOCAL_DAMAGE_BRITTLE_ID) - call damage_brittle_putLocalDamage(ipc, ip, el, localDamage) + case (LOCAL_DAMAGE_isoBrittle_ID) + call damage_isoBrittle_putLocalDamage(ipc, ip, el, localDamage) - case (LOCAL_DAMAGE_DUCTILE_ID) - call damage_ductile_putLocalDamage(ipc, ip, el, localDamage) + case (LOCAL_DAMAGE_isoDuctile_ID) + call damage_isoDuctile_putLocalDamage(ipc, ip, el, localDamage) case (LOCAL_DAMAGE_gurson_ID) call damage_gurson_putLocalDamage(ipc, ip, el, localDamage) - case (LOCAL_DAMAGE_anisotropic_ID) - call damage_anisotropic_putLocalDamage(ipc, ip, el, localDamage) + case (LOCAL_DAMAGE_anisoBrittle_ID) + call damage_anisoBrittle_putLocalDamage(ipc, ip, el, localDamage) end select @@ -1081,19 +1081,19 @@ function constitutive_getDamage(ipc, ip, el) use material, only: & material_phase, & LOCAL_DAMAGE_none_ID, & - LOCAL_DAMAGE_BRITTLE_ID, & - LOCAL_DAMAGE_DUCTILE_ID, & + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID, & phase_damage - use damage_brittle, only: & - damage_brittle_getDamage - use damage_ductile, only: & - damage_ductile_getDamage + use damage_isoBrittle, only: & + damage_isoBrittle_getDamage + use damage_isoDuctile, only: & + damage_isoDuctile_getDamage use damage_gurson, only: & damage_gurson_getDamage - use damage_anisotropic, only: & - damage_anisotropic_getDamage + use damage_anisoBrittle, only: & + damage_anisoBrittle_getDamage implicit none integer(pInt), intent(in) :: & @@ -1106,17 +1106,17 @@ function constitutive_getDamage(ipc, ip, el) case (LOCAL_DAMAGE_none_ID) constitutive_getDamage = 1.0_pReal - case (LOCAL_DAMAGE_brittle_ID) - constitutive_getDamage = damage_brittle_getDamage(ipc, ip, el) + case (LOCAL_DAMAGE_isoBrittle_ID) + constitutive_getDamage = damage_isoBrittle_getDamage(ipc, ip, el) - case (LOCAL_DAMAGE_ductile_ID) - constitutive_getDamage = damage_ductile_getDamage(ipc, ip, el) + case (LOCAL_DAMAGE_isoDuctile_ID) + constitutive_getDamage = damage_isoDuctile_getDamage(ipc, ip, el) case (LOCAL_DAMAGE_gurson_ID) constitutive_getDamage = damage_gurson_getDamage(ipc, ip, el) - case (LOCAL_DAMAGE_anisotropic_ID) - constitutive_getDamage = damage_anisotropic_getDamage(ipc, ip, el) + case (LOCAL_DAMAGE_anisoBrittle_ID) + constitutive_getDamage = damage_anisoBrittle_getDamage(ipc, ip, el) end select end function constitutive_getDamage @@ -1129,11 +1129,11 @@ function constitutive_getSlipDamage(nSlip, Tstar_v, ipc, ip, el) pReal use material, only: & material_phase, & - LOCAL_DAMAGE_ductile_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & LOCAL_DAMAGE_gurson_ID, & phase_damage - use damage_ductile, only: & - damage_ductile_getSlipDamage + use damage_isoDuctile, only: & + damage_isoDuctile_getSlipDamage use damage_gurson, only: & damage_gurson_getSlipDamage @@ -1148,8 +1148,8 @@ function constitutive_getSlipDamage(nSlip, Tstar_v, ipc, ip, el) real(pReal) :: constitutive_getSlipDamage(nSlip) select case (phase_damage(material_phase(ipc,ip,el))) - case (LOCAL_DAMAGE_ductile_ID) - constitutive_getSlipDamage = damage_ductile_getSlipDamage(ipc, ip, el) + case (LOCAL_DAMAGE_isoDuctile_ID) + constitutive_getSlipDamage = damage_isoDuctile_getSlipDamage(ipc, ip, el) case (LOCAL_DAMAGE_gurson_ID) constitutive_getSlipDamage = damage_gurson_getSlipDamage(Tstar_v, ipc, ip, el) @@ -1170,10 +1170,10 @@ function constitutive_getDamageStrain(ipc, ip, el) math_I3 use material, only: & material_phase, & - LOCAL_DAMAGE_anisotropic_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & phase_damage - use damage_anisotropic, only: & - damage_anisotropic_getDamageStrain + use damage_anisoBrittle, only: & + damage_anisoBrittle_getDamageStrain implicit none @@ -1186,8 +1186,8 @@ function constitutive_getDamageStrain(ipc, ip, el) select case (phase_damage(material_phase(ipc,ip,el))) - case (LOCAL_DAMAGE_anisotropic_ID) - constitutive_getDamageStrain = damage_anisotropic_getDamageStrain(ipc, ip, el) + case (LOCAL_DAMAGE_anisoBrittle_ID) + constitutive_getDamageStrain = damage_anisoBrittle_getDamageStrain(ipc, ip, el) case default constitutive_getDamageStrain = math_I3 @@ -1206,13 +1206,13 @@ function constitutive_getDamageDiffusion33(ipc, ip, el) use material, only: & material_phase, & LOCAL_DAMAGE_none_ID, & - LOCAL_DAMAGE_brittle_ID, & - LOCAL_DAMAGE_ductile_ID, & + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID, & phase_damage - use damage_brittle, only: & - damage_brittle_getDamageDiffusion33 + use damage_isoBrittle, only: & + damage_isoBrittle_getDamageDiffusion33 implicit none integer(pInt), intent(in) :: & @@ -1224,8 +1224,8 @@ function constitutive_getDamageDiffusion33(ipc, ip, el) constitutive_getDamageDiffusion33 = lattice_DamageDiffusion33(1:3,1:3,material_phase(ipc,ip,el)) select case(phase_damage(material_phase(ipc,ip,el))) - case (LOCAL_DAMAGE_brittle_ID) - constitutive_getDamageDiffusion33 = damage_brittle_getDamageDiffusion33(ipc, ip, el) + case (LOCAL_DAMAGE_isoBrittle_ID) + constitutive_getDamageDiffusion33 = damage_isoBrittle_getDamageDiffusion33(ipc, ip, el) end select @@ -1267,7 +1267,7 @@ function constitutive_getAdiabaticTemperature(ipc, ip, el) end function constitutive_getAdiabaticTemperature !-------------------------------------------------------------------------------------------------- -!> @brief Returns the local(unregularised) damage +!> @brief assigns the local/nonlocal value of temperature to local thermal state !-------------------------------------------------------------------------------------------------- subroutine constitutive_putAdiabaticTemperature(ipc, ip, el, localTemperature) use prec, only: & @@ -1656,10 +1656,10 @@ function constitutive_postResults(Tstar_v, FeArray, ipc, ip, el) PLASTICITY_DISLOKMC_ID, & PLASTICITY_TITANMOD_ID, & PLASTICITY_NONLOCAL_ID, & - LOCAL_DAMAGE_BRITTLE_ID, & - LOCAL_DAMAGE_DUCTILE_ID, & + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID, & LOCAL_THERMAL_ADIABATIC_ID, & LOCAL_VACANCY_generation_ID use constitutive_j2, only: & @@ -1678,14 +1678,14 @@ function constitutive_postResults(Tstar_v, FeArray, ipc, ip, el) use constitutive_nonlocal, only: & constitutive_nonlocal_postResults #ifdef multiphysicsOut - use damage_brittle, only: & - damage_brittle_postResults - use damage_ductile, only: & - damage_ductile_postResults + use damage_isoBrittle, only: & + damage_isoBrittle_postResults + use damage_isoDuctile, only: & + damage_isoDuctile_postResults use damage_gurson, only: & damage_gurson_postResults - use damage_anisotropic, only: & - damage_anisotropic_postResults + use damage_anisoBrittle, only: & + damage_anisoBrittle_postResults use thermal_adiabatic, only: & thermal_adiabatic_postResults use vacancy_generation, only: & @@ -1741,14 +1741,14 @@ function constitutive_postResults(Tstar_v, FeArray, ipc, ip, el) startPos = endPos + 1_pInt endPos = endPos + damageState(material_phase(ipc,ip,el))%sizePostResults select case (phase_damage(material_phase(ipc,ip,el))) - case (LOCAL_DAMAGE_BRITTLE_ID) - constitutive_postResults(startPos:endPos) = damage_brittle_postResults(ipc, ip, el) - case (LOCAL_DAMAGE_DUCTILE_ID) - constitutive_postResults(startPos:endPos) = damage_ductile_postResults(ipc, ip, el) + case (LOCAL_DAMAGE_isoBrittle_ID) + constitutive_postResults(startPos:endPos) = damage_isoBrittle_postResults(ipc, ip, el) + case (LOCAL_DAMAGE_isoDuctile_ID) + constitutive_postResults(startPos:endPos) = damage_isoDuctile_postResults(ipc, ip, el) case (LOCAL_DAMAGE_gurson_ID) constitutive_postResults(startPos:endPos) = damage_gurson_postResults(ipc, ip, el) - case (LOCAL_DAMAGE_anisotropic_ID) - constitutive_postResults(startPos:endPos) = damage_anisotropic_postResults(ipc, ip, el) + case (LOCAL_DAMAGE_anisoBrittle_ID) + constitutive_postResults(startPos:endPos) = damage_anisoBrittle_postResults(ipc, ip, el) end select startPos = endPos + 1_pInt diff --git a/code/damage_anisoBrittle.f90 b/code/damage_anisoBrittle.f90 new file mode 100644 index 000000000..16df31125 --- /dev/null +++ b/code/damage_anisoBrittle.f90 @@ -0,0 +1,517 @@ +!-------------------------------------------------------------------------------------------------- +! $Id: damage_anisoBrittle.f90 3210 2014-06-17 15:24:44Z MPIE\m.diehl $ +!-------------------------------------------------------------------------------------------------- +!> @author Luv Sharma, Max-Planck-Institut für Eisenforschung GmbH +!> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH +!> @brief material subroutine incorporating anisotropic ductile damage +!> @details to be done +!-------------------------------------------------------------------------------------------------- +module damage_anisoBrittle + use prec, only: & + pReal, & + pInt + + implicit none + private + integer(pInt), dimension(:), allocatable, public, protected :: & + damage_anisoBrittle_sizePostResults !< cumulative size of post results + + integer(pInt), dimension(:,:), allocatable, target, public :: & + damage_anisoBrittle_sizePostResult !< size of each post result output + + character(len=64), dimension(:,:), allocatable, target, public :: & + damage_anisoBrittle_output !< name of each post result output + + integer(pInt), dimension(:), allocatable, target, public :: & + damage_anisoBrittle_Noutput !< number of outputs per instance of this damage + + integer(pInt), dimension(:), allocatable, private :: & + damage_anisoBrittle_totalNslip !< Todo + + integer(pInt), dimension(:,:), allocatable, private :: & + damage_anisoBrittle_Nslip !< Todo + + real(pReal), dimension(:), allocatable, private :: & + damage_anisoBrittle_aTol, & + damage_anisoBrittle_sdot_0, & + damage_anisoBrittle_N + + real(pReal), dimension(:,:), allocatable, private :: & + damage_anisoBrittle_critDisp, & + damage_anisoBrittle_critLoad + + enum, bind(c) + enumerator :: undefined_ID, & + local_damage_ID + end enum !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!11 ToDo + + integer(kind(undefined_ID)), dimension(:,:), allocatable, private :: & + damage_anisoBrittle_outputID !< ID of each post result output + + + public :: & + damage_anisoBrittle_init, & + damage_anisoBrittle_stateInit, & + damage_anisoBrittle_aTolState, & + damage_anisoBrittle_dotState, & + damage_anisoBrittle_getDamage, & + damage_anisoBrittle_putLocalDamage, & + damage_anisoBrittle_getLocalDamage, & + damage_anisoBrittle_getDamageStrain, & + damage_anisoBrittle_postResults + +contains + + +!-------------------------------------------------------------------------------------------------- +!> @brief module initialization +!> @details reads in material parameters, allocates arrays, and does sanity checks +!-------------------------------------------------------------------------------------------------- +subroutine damage_anisoBrittle_init(fileUnit) + use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment) + use debug, only: & + debug_level,& + debug_constitutive,& + debug_levelBasic + use mesh, only: & + mesh_maxNips, & + mesh_NcpElems + use IO, only: & + IO_read, & + IO_lc, & + IO_getTag, & + IO_isBlank, & + IO_stringPos, & + IO_stringValue, & + IO_floatValue, & + IO_intValue, & + IO_warning, & + IO_error, & + IO_timeStamp, & + IO_EOF + use material, only: & + homogenization_maxNgrains, & + phase_damage, & + phase_damageInstance, & + phase_Noutput, & + LOCAL_damage_anisoBrittle_label, & + LOCAL_damage_anisoBrittle_ID, & + material_phase, & + damageState, & + MATERIAL_partPhase + use numerics,only: & + worldrank, & + numerics_integrator + use lattice, only: & + lattice_maxNslipFamily + + implicit none + integer(pInt), intent(in) :: fileUnit + + integer(pInt), parameter :: MAXNCHUNKS = 7_pInt + integer(pInt), dimension(1+2*MAXNCHUNKS) :: positions + integer(pInt) :: maxNinstance,mySize=0_pInt,phase,instance,o + integer(pInt) :: sizeState, sizeDotState + integer(pInt) :: NofMyPhase + integer(pInt) :: Nchunks_SlipFamilies, j + character(len=65536) :: & + tag = '', & + line = '' + + mainProcess: if (worldrank == 0) then + write(6,'(/,a)') ' <<<+- damage_'//LOCAL_damage_anisoBrittle_LABEL//' init -+>>>' + write(6,'(a)') ' $Id: damage_anisoBrittle.f90 3210 2014-06-17 15:24:44Z MPIE\m.diehl $' + write(6,'(a15,a)') ' Current time: ',IO_timeStamp() +#include "compilation_info.f90" + endif mainProcess + + maxNinstance = int(count(phase_damage == LOCAL_damage_anisoBrittle_ID),pInt) + if (maxNinstance == 0_pInt) return + + if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) & + write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance + + allocate(damage_anisoBrittle_sizePostResults(maxNinstance), source=0_pInt) + allocate(damage_anisoBrittle_sizePostResult(maxval(phase_Noutput),maxNinstance),source=0_pInt) + allocate(damage_anisoBrittle_output(maxval(phase_Noutput),maxNinstance)) + damage_anisoBrittle_output = '' + allocate(damage_anisoBrittle_outputID(maxval(phase_Noutput),maxNinstance), source=undefined_ID) + allocate(damage_anisoBrittle_Noutput(maxNinstance), source=0_pInt) + allocate(damage_anisoBrittle_critDisp(lattice_maxNslipFamily,maxNinstance), source=0.0_pReal) + allocate(damage_anisoBrittle_critLoad (lattice_maxNslipFamily,maxNinstance), source=0.0_pReal) + allocate(damage_anisoBrittle_Nslip(lattice_maxNslipFamily,maxNinstance), source=0_pInt) + allocate(damage_anisoBrittle_totalNslip(maxNinstance), source=0_pInt) + allocate(damage_anisoBrittle_aTol(maxNinstance), source=0.0_pReal) + allocate(damage_anisoBrittle_sdot_0(maxNinstance), source=0.0_pReal) + allocate(damage_anisoBrittle_N(maxNinstance), source=0.0_pReal) + + rewind(fileUnit) + phase = 0_pInt + do while (trim(line) /= IO_EOF .and. IO_lc(IO_getTag(line,'<','>')) /= MATERIAL_partPhase) ! wind forward to + line = IO_read(fileUnit) + enddo + + parsingFile: do while (trim(line) /= IO_EOF) ! read through sections of phase part + line = IO_read(fileUnit) + if (IO_isBlank(line)) cycle ! skip empty lines + if (IO_getTag(line,'<','>') /= '') then ! stop at next part + line = IO_read(fileUnit, .true.) ! reset IO_read + exit + endif + if (IO_getTag(line,'[',']') /= '') then ! next phase section + phase = phase + 1_pInt ! advance phase section counter + cycle ! skip to next line + endif + if (phase > 0_pInt ) then; if (phase_damage(phase) == LOCAL_damage_anisoBrittle_ID) then ! do not short-circuit here (.and. with next if statemen). It's not safe in Fortran + instance = phase_damageInstance(phase) ! which instance of my damage is present phase + positions = IO_stringPos(line,MAXNCHUNKS) + tag = IO_lc(IO_stringValue(line,positions,1_pInt)) ! extract key + select case(tag) + case ('(output)') + select case(IO_lc(IO_stringValue(line,positions,2_pInt))) + case ('local_damage') + damage_anisoBrittle_Noutput(instance) = damage_anisoBrittle_Noutput(instance) + 1_pInt + damage_anisoBrittle_outputID(damage_anisoBrittle_Noutput(instance),instance) = local_damage_ID + damage_anisoBrittle_output(damage_anisoBrittle_Noutput(instance),instance) = & + IO_lc(IO_stringValue(line,positions,2_pInt)) + end select + + case ('atol_damage') + damage_anisoBrittle_aTol(instance) = IO_floatValue(line,positions,2_pInt) + + case ('sdot_0') + damage_anisoBrittle_sdot_0(instance) = IO_floatValue(line,positions,2_pInt) + + case ('n_damage') + damage_anisoBrittle_N(instance) = IO_floatValue(line,positions,2_pInt) + + case ('Nslip') ! + Nchunks_SlipFamilies = positions(1) - 1_pInt + do j = 1_pInt, Nchunks_SlipFamilies + damage_anisoBrittle_Nslip(j,instance) = IO_intValue(line,positions,1_pInt+j) + enddo + damage_anisoBrittle_totalNslip(instance) = sum(damage_anisoBrittle_Nslip(:,instance)) + + case ('critical_displacement') + do j = 1_pInt, Nchunks_SlipFamilies + damage_anisoBrittle_critDisp(j,instance) = IO_floatValue(line,positions,1_pInt+j) + enddo + + case ('critical_load') + do j = 1_pInt, Nchunks_SlipFamilies + damage_anisoBrittle_critLoad(j,instance) = IO_floatValue(line,positions,1_pInt+j) + enddo + + end select + endif; endif + enddo parsingFile + + initializeInstances: do phase = 1_pInt, size(phase_damage) + if (phase_damage(phase) == LOCAL_damage_anisoBrittle_ID) then + NofMyPhase=count(material_phase==phase) + instance = phase_damageInstance(phase) + +!-------------------------------------------------------------------------------------------------- +! Determine size of postResults array + outputsLoop: do o = 1_pInt,damage_anisoBrittle_Noutput(instance) + select case(damage_anisoBrittle_outputID(o,instance)) + case(local_damage_ID) + mySize = damage_anisoBrittle_totalNslip(instance) + end select + + if (mySize > 0_pInt) then ! any meaningful output found + damage_anisoBrittle_sizePostResult(o,instance) = mySize + damage_anisoBrittle_sizePostResults(instance) = damage_anisoBrittle_sizePostResults(instance) + mySize + endif + enddo outputsLoop +! Determine size of state array + sizeDotState = 2_pInt + & ! viscous and non-viscous damage values + 9_pInt + & ! damage deformation gradient + damage_anisoBrittle_totalNslip(instance) ! opening on each damage system + sizeState = sizeDotState + + damageState(phase)%sizeState = sizeState + damageState(phase)%sizeDotState = sizeDotState + damageState(phase)%sizePostResults = damage_anisoBrittle_sizePostResults(instance) + allocate(damageState(phase)%aTolState (sizeState), source=0.0_pReal) + allocate(damageState(phase)%state0 (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%partionedState0 (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%subState0 (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%state (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%state_backup (sizeState,NofMyPhase), source=0.0_pReal) + + allocate(damageState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%deltaState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal) + if (any(numerics_integrator == 1_pInt)) then + allocate(damageState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal) + endif + if (any(numerics_integrator == 4_pInt)) & + allocate(damageState(phase)%RK4dotState (sizeDotState,NofMyPhase), source=0.0_pReal) + if (any(numerics_integrator == 5_pInt)) & + allocate(damageState(phase)%RKCK45dotState (6,sizeDotState,NofMyPhase),source=0.0_pReal) + + call damage_anisoBrittle_stateInit(phase) + call damage_anisoBrittle_aTolState(phase,instance) + endif + + enddo initializeInstances +end subroutine damage_anisoBrittle_init + +!-------------------------------------------------------------------------------------------------- +!> @brief sets the relevant state values for a given instance of this damage +!-------------------------------------------------------------------------------------------------- +subroutine damage_anisoBrittle_stateInit(phase) + use material, only: & + damageState + + implicit none + integer(pInt), intent(in) :: phase !< number specifying the phase of the damage + + real(pReal), dimension(damageState(phase)%sizeState) :: tempState + + tempState = 0.0_pReal + tempState(1) = 1.0_pReal + tempState(2) = 1.0_pReal + tempState(3) = 1.0_pReal + tempState(4) = 0.0_pReal + tempState(5) = 0.0_pReal + tempState(6) = 0.0_pReal + tempState(7) = 1.0_pReal + tempState(8) = 0.0_pReal + tempState(9) = 0.0_pReal + tempState(10) = 0.0_pReal + tempState(11) = 1.0_pReal + damageState(phase)%state = spread(tempState,2,size(damageState(phase)%state(1,:))) + damageState(phase)%state0 = damageState(phase)%state + damageState(phase)%partionedState0 = damageState(phase)%state +end subroutine damage_anisoBrittle_stateInit + +!-------------------------------------------------------------------------------------------------- +!> @brief sets the relevant state values for a given instance of this damage +!-------------------------------------------------------------------------------------------------- +subroutine damage_anisoBrittle_aTolState(phase,instance) + use material, only: & + damageState + + implicit none + integer(pInt), intent(in) :: & + phase, & + instance ! number specifying the current instance of the damage + real(pReal), dimension(damageState(phase)%sizeState) :: tempTol + + tempTol = damage_anisoBrittle_aTol(instance) + damageState(phase)%aTolState = tempTol +end subroutine damage_anisoBrittle_aTolState + +!-------------------------------------------------------------------------------------------------- +!> @brief calculates derived quantities from state +!-------------------------------------------------------------------------------------------------- +subroutine damage_anisoBrittle_dotState(Tstar_v,ipc, ip, el) + use material, only: & + mappingConstitutive, & + phase_damageInstance, & + damageState + use math, only: & + math_mul33x33 + use lattice, only: & + lattice_Sslip, & + lattice_Sslip_v, & + lattice_maxNslipFamily, & + lattice_NslipSystem, & + lattice_DamageMobility + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< component-ID of integration point + ip, & !< integration point + el !< element + real(pReal), intent(in), dimension(6) :: & + Tstar_v !< 2nd Piola Kirchhoff stress tensor (Mandel) + integer(pInt) :: & + phase, & + constituent, & + instance, & + j, f, i, index_myFamily + real(pReal), dimension(3,3) :: & + Ld + real(pReal) :: & + tau, & + tau_critical, & + nonLocalFactor + + phase = mappingConstitutive(2,ipc,ip,el) + constituent = mappingConstitutive(1,ipc,ip,el) + instance = phase_damageInstance(phase) + + damageState(phase)%dotState(1,constituent) = & + (1.0_pReal/lattice_DamageMobility(phase))* & + (damageState(phase)%state(2,constituent) - & + damageState(phase)%state(1,constituent)) + + nonLocalFactor = 1.0_pReal + & + (damageState(phase)%state(1,constituent) - & + damage_anisoBrittle_getDamage(ipc, ip, el)) + Ld = 0.0_pReal + j = 0_pInt + slipFamiliesLoop: do f = 1_pInt,lattice_maxNslipFamily + index_myFamily = sum(lattice_NslipSystem(1:f-1_pInt,phase)) ! at which index starts my family + do i = 1_pInt,damage_anisoBrittle_Nslip(f,instance) ! process each (active) slip system in family + j = j+1_pInt + tau = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,phase)) + tau_critical = (1.0_pReal - damageState(phase)%state(11+j,constituent)/& + damage_anisoBrittle_critDisp(f,instance))* & + damage_anisoBrittle_critLoad(f,instance)*nonLocalFactor + damageState(phase)%dotState(11+j,constituent) = & + damage_anisoBrittle_sdot_0(instance)*(tau/tau_critical)**damage_anisoBrittle_N(instance) + damageState(phase)%dotState(2,constituent) = & + damageState(phase)%dotState(2,constituent) - & + 2.0_pReal*tau*damageState(phase)%dotState(11+j,constituent)/ & + (damage_anisoBrittle_critDisp(f,instance)*damage_anisoBrittle_critLoad(f,instance)) + Ld = Ld + damageState(phase)%dotState(11+j,constituent)* & + lattice_Sslip(1:3,1:3,1,index_myFamily+i,phase) + enddo + enddo slipFamiliesLoop + damageState(phase)%dotState(3:11,constituent) = & + reshape(math_mul33x33(Ld,reshape(damageState(phase)%state(3:11,constituent),shape=[3,3])),shape=[9]) + +end subroutine damage_anisoBrittle_dotState + +!-------------------------------------------------------------------------------------------------- +!> @brief returns damage +!-------------------------------------------------------------------------------------------------- +function damage_anisoBrittle_getDamage(ipc, ip, el) + use material, only: & + material_homog, & + mappingHomogenization, & + fieldDamage, & + field_damage_type, & + FIELD_DAMAGE_LOCAL_ID, & + FIELD_DAMAGE_NONLOCAL_ID + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal) :: damage_anisoBrittle_getDamage + + select case(field_damage_type(material_homog(ip,el))) + case (FIELD_DAMAGE_LOCAL_ID) + damage_anisoBrittle_getDamage = damage_anisoBrittle_getLocalDamage(ipc, ip, el) + + case (FIELD_DAMAGE_NONLOCAL_ID) + damage_anisoBrittle_getDamage = fieldDamage(material_homog(ip,el))% & + field(1,mappingHomogenization(1,ip,el)) ! Taylor type + + end select + +end function damage_anisoBrittle_getDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief returns damage value based on local damage +!-------------------------------------------------------------------------------------------------- +subroutine damage_anisoBrittle_putLocalDamage(ipc, ip, el, localDamage) + use material, only: & + mappingConstitutive, & + phase_damageInstance, & + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal), intent(in) :: & + localDamage + + damageState(mappingConstitutive(2,ipc,ip,el))%state(1,mappingConstitutive(1,ipc,ip,el)) = & + localDamage + +end subroutine damage_anisoBrittle_putLocalDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief returns local damage +!-------------------------------------------------------------------------------------------------- +function damage_anisoBrittle_getLocalDamage(ipc, ip, el) + use material, only: & + mappingConstitutive, & + phase_damageInstance, & + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal) :: & + damage_anisoBrittle_getLocalDamage + + damage_anisoBrittle_getLocalDamage = & + damageState(mappingConstitutive(2,ipc,ip,el))%state(1,mappingConstitutive(1,ipc,ip,el)) + +end function damage_anisoBrittle_getLocalDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief returns local damage deformation gradient +!-------------------------------------------------------------------------------------------------- +function damage_anisoBrittle_getDamageStrain(ipc, ip, el) + use material, only: & + mappingConstitutive, & + phase_damageInstance, & + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal), dimension(3,3) :: & + damage_anisoBrittle_getDamageStrain + + damage_anisoBrittle_getDamageStrain = & + reshape(damageState(mappingConstitutive(2,ipc,ip,el))%state(3:11,mappingConstitutive(1,ipc,ip,el)), & + shape=[3,3]) + +end function damage_anisoBrittle_getDamageStrain + +!-------------------------------------------------------------------------------------------------- +!> @brief return array of constitutive results +!-------------------------------------------------------------------------------------------------- +function damage_anisoBrittle_postResults(ipc,ip,el) + use material, only: & + mappingConstitutive, & + phase_damageInstance,& + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< component-ID of integration point + ip, & !< integration point + el !< element + real(pReal), dimension(damage_anisoBrittle_sizePostResults(phase_damageInstance(mappingConstitutive(2,ipc,ip,el)))) :: & + damage_anisoBrittle_postResults + + integer(pInt) :: & + instance, phase, constituent, o, c + + phase = mappingConstitutive(2,ipc,ip,el) + constituent = mappingConstitutive(1,ipc,ip,el) + instance = phase_damageInstance(phase) + + c = 0_pInt + damage_anisoBrittle_postResults = 0.0_pReal + + do o = 1_pInt,damage_anisoBrittle_Noutput(instance) + select case(damage_anisoBrittle_outputID(o,instance)) + case (local_damage_ID) + damage_anisoBrittle_postResults(c+1_pInt:c+damage_anisoBrittle_totalNslip(instance)) = & + damageState(phase)%state(1,constituent) + c = c + damage_anisoBrittle_totalNslip(instance) + + end select + enddo +end function damage_anisoBrittle_postResults + +end module damage_anisoBrittle diff --git a/code/damage_isoBrittle.f90 b/code/damage_isoBrittle.f90 new file mode 100644 index 000000000..b88f37561 --- /dev/null +++ b/code/damage_isoBrittle.f90 @@ -0,0 +1,462 @@ +!-------------------------------------------------------------------------------------------------- +! $Id$ +!-------------------------------------------------------------------------------------------------- +!> @author Franz Roters, Max-Planck-Institut für Eisenforschung GmbH +!> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH +!> @brief material subroutine incoprorating dislocation and twinning physics +!> @details to be done +!-------------------------------------------------------------------------------------------------- +module damage_isoBrittle + use prec, only: & + pReal, & + pInt + + implicit none + private + integer(pInt), dimension(:), allocatable, public, protected :: & + damage_isoBrittle_sizePostResults !< cumulative size of post results + + integer(pInt), dimension(:,:), allocatable, target, public :: & + damage_isoBrittle_sizePostResult !< size of each post result output + + character(len=64), dimension(:,:), allocatable, target, public :: & + damage_isoBrittle_output !< name of each post result output + + integer(pInt), dimension(:), allocatable, target, public :: & + damage_isoBrittle_Noutput !< number of outputs per instance of this damage + + real(pReal), dimension(:), allocatable, private :: & + damage_isoBrittle_aTol, & + damage_isoBrittle_critStrainEnergy + + enum, bind(c) + enumerator :: undefined_ID, & + local_damage_ID + end enum !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!11 ToDo + + integer(kind(undefined_ID)), dimension(:,:), allocatable, private :: & + damage_isoBrittle_outputID !< ID of each post result output + + + public :: & + damage_isoBrittle_init, & + damage_isoBrittle_stateInit, & + damage_isoBrittle_aTolState, & + damage_isoBrittle_dotState, & + damage_isoBrittle_microstructure, & + damage_isoBrittle_getDamage, & + damage_isoBrittle_putLocalDamage, & + damage_isoBrittle_getLocalDamage, & + damage_isoBrittle_getDamageDiffusion33, & + damage_isoBrittle_postResults + +contains + + +!-------------------------------------------------------------------------------------------------- +!> @brief module initialization +!> @details reads in material parameters, allocates arrays, and does sanity checks +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoBrittle_init(fileUnit) + use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment) + use debug, only: & + debug_level,& + debug_constitutive,& + debug_levelBasic + use mesh, only: & + mesh_maxNips, & + mesh_NcpElems + use IO, only: & + IO_read, & + IO_lc, & + IO_getTag, & + IO_isBlank, & + IO_stringPos, & + IO_stringValue, & + IO_floatValue, & + IO_intValue, & + IO_warning, & + IO_error, & + IO_timeStamp, & + IO_EOF + use material, only: & + homogenization_maxNgrains, & + phase_damage, & + phase_damageInstance, & + phase_Noutput, & + LOCAL_damage_isoBrittle_label, & + LOCAL_damage_isoBrittle_ID, & + material_phase, & + damageState, & + MATERIAL_partPhase + use numerics,only: & + worldrank, & + numerics_integrator + + implicit none + integer(pInt), intent(in) :: fileUnit + + integer(pInt), parameter :: MAXNCHUNKS = 7_pInt + integer(pInt), dimension(1+2*MAXNCHUNKS) :: positions + integer(pInt) :: maxNinstance,mySize=0_pInt,phase,instance,o + integer(pInt) :: sizeState, sizeDotState + integer(pInt) :: NofMyPhase + character(len=65536) :: & + tag = '', & + line = '' + + mainProcess: if (worldrank == 0) then + write(6,'(/,a)') ' <<<+- damage_'//LOCAL_damage_isoBrittle_label//' init -+>>>' + write(6,'(a)') ' $Id$' + write(6,'(a15,a)') ' Current time: ',IO_timeStamp() +#include "compilation_info.f90" + endif mainProcess + + maxNinstance = int(count(phase_damage == LOCAL_damage_isoBrittle_ID),pInt) + if (maxNinstance == 0_pInt) return + + if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) & + write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance + + allocate(damage_isoBrittle_sizePostResults(maxNinstance), source=0_pInt) + allocate(damage_isoBrittle_sizePostResult(maxval(phase_Noutput),maxNinstance),source=0_pInt) + allocate(damage_isoBrittle_output(maxval(phase_Noutput),maxNinstance)) + damage_isoBrittle_output = '' + allocate(damage_isoBrittle_outputID(maxval(phase_Noutput),maxNinstance), source=undefined_ID) + allocate(damage_isoBrittle_Noutput(maxNinstance), source=0_pInt) + allocate(damage_isoBrittle_critStrainEnergy(maxNinstance), source=0.0_pReal) + allocate(damage_isoBrittle_aTol(maxNinstance), source=0.0_pReal) + + rewind(fileUnit) + phase = 0_pInt + do while (trim(line) /= IO_EOF .and. IO_lc(IO_getTag(line,'<','>')) /= MATERIAL_partPhase) ! wind forward to + line = IO_read(fileUnit) + enddo + + parsingFile: do while (trim(line) /= IO_EOF) ! read through sections of phase part + line = IO_read(fileUnit) + if (IO_isBlank(line)) cycle ! skip empty lines + if (IO_getTag(line,'<','>') /= '') then ! stop at next part + line = IO_read(fileUnit, .true.) ! reset IO_read + exit + endif + if (IO_getTag(line,'[',']') /= '') then ! next phase section + phase = phase + 1_pInt ! advance phase section counter + cycle ! skip to next line + endif + if (phase > 0_pInt ) then; if (phase_damage(phase) == LOCAL_damage_isoBrittle_ID) then ! do not short-circuit here (.and. with next if statemen). It's not safe in Fortran + instance = phase_damageInstance(phase) ! which instance of my damage is present phase + positions = IO_stringPos(line,MAXNCHUNKS) + tag = IO_lc(IO_stringValue(line,positions,1_pInt)) ! extract key + select case(tag) + case ('(output)') + select case(IO_lc(IO_stringValue(line,positions,2_pInt))) + case ('local_damage') + damage_isoBrittle_Noutput(instance) = damage_isoBrittle_Noutput(instance) + 1_pInt + damage_isoBrittle_outputID(damage_isoBrittle_Noutput(instance),instance) = local_damage_ID + damage_isoBrittle_output(damage_isoBrittle_Noutput(instance),instance) = & + IO_lc(IO_stringValue(line,positions,2_pInt)) + end select + + case ('critical_strain_energy') + damage_isoBrittle_critStrainEnergy(instance) = IO_floatValue(line,positions,2_pInt) + + case ('atol_damage') + damage_isoBrittle_aTol(instance) = IO_floatValue(line,positions,2_pInt) + + end select + endif; endif + enddo parsingFile + + initializeInstances: do phase = 1_pInt, size(phase_damage) + if (phase_damage(phase) == LOCAL_damage_isoBrittle_ID) then + NofMyPhase=count(material_phase==phase) + instance = phase_damageInstance(phase) + +!-------------------------------------------------------------------------------------------------- +! Determine size of postResults array + outputsLoop: do o = 1_pInt,damage_isoBrittle_Noutput(instance) + select case(damage_isoBrittle_outputID(o,instance)) + case(local_damage_ID) + mySize = 1_pInt + end select + + if (mySize > 0_pInt) then ! any meaningful output found + damage_isoBrittle_sizePostResult(o,instance) = mySize + damage_isoBrittle_sizePostResults(instance) = damage_isoBrittle_sizePostResults(instance) + mySize + endif + enddo outputsLoop +! Determine size of state array + sizeDotState = 1_pInt + sizeState = 2_pInt + + damageState(phase)%sizeState = sizeState + damageState(phase)%sizeDotState = sizeDotState + damageState(phase)%sizePostResults = damage_isoBrittle_sizePostResults(instance) + allocate(damageState(phase)%aTolState (sizeState), source=0.0_pReal) + allocate(damageState(phase)%state0 (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%partionedState0 (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%subState0 (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%state (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%state_backup (sizeState,NofMyPhase), source=0.0_pReal) + + allocate(damageState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%deltaState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal) + if (any(numerics_integrator == 1_pInt)) then + allocate(damageState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal) + endif + if (any(numerics_integrator == 4_pInt)) & + allocate(damageState(phase)%RK4dotState (sizeDotState,NofMyPhase), source=0.0_pReal) + if (any(numerics_integrator == 5_pInt)) & + allocate(damageState(phase)%RKCK45dotState (6,sizeDotState,NofMyPhase),source=0.0_pReal) + + call damage_isoBrittle_stateInit(phase) + call damage_isoBrittle_aTolState(phase,instance) + endif + + enddo initializeInstances +end subroutine damage_isoBrittle_init + +!-------------------------------------------------------------------------------------------------- +!> @brief sets the relevant NEW state values for a given instance of this damage +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoBrittle_stateInit(phase) + use material, only: & + damageState + + implicit none + integer(pInt), intent(in) :: phase !< number specifying the phase of the damage + + real(pReal), dimension(damageState(phase)%sizeState) :: tempState + + tempState(1) = 1.0_pReal + tempState(2) = 1.0_pReal + damageState(phase)%state = spread(tempState,2,size(damageState(phase)%state(1,:))) + damageState(phase)%state0 = damageState(phase)%state + damageState(phase)%partionedState0 = damageState(phase)%state +end subroutine damage_isoBrittle_stateInit + +!-------------------------------------------------------------------------------------------------- +!> @brief sets the relevant state values for a given instance of this damage +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoBrittle_aTolState(phase,instance) + use material, only: & + damageState + + implicit none + integer(pInt), intent(in) :: & + phase, & + instance ! number specifying the current instance of the damage + real(pReal), dimension(damageState(phase)%sizeState) :: tempTol + + tempTol = damage_isoBrittle_aTol(instance) + damageState(phase)%aTolState = tempTol +end subroutine damage_isoBrittle_aTolState + +!-------------------------------------------------------------------------------------------------- +!> @brief calculates derived quantities from state +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoBrittle_dotState(ipc, ip, el) + use material, only: & + mappingConstitutive, & + phase_damageInstance, & + damageState + use lattice, only: & + lattice_DamageMobility + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< component-ID of integration point + ip, & !< integration point + el !< element + integer(pInt) :: & + phase, constituent + + phase = mappingConstitutive(2,ipc,ip,el) + constituent = mappingConstitutive(1,ipc,ip,el) + + damageState(phase)%dotState(1,constituent) = & + (1.0_pReal/lattice_DamageMobility(phase))* & + (damageState(phase)%state(2,constituent) - & + damageState(phase)%state(1,constituent)) + +end subroutine damage_isoBrittle_dotState + +!-------------------------------------------------------------------------------------------------- +!> @brief calculates derived quantities from state +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoBrittle_microstructure(Tstar_v, Fe, ipc, ip, el) + use material, only: & + mappingConstitutive, & + phase_damageInstance, & + damageState + use math, only: & + math_Mandel6to33, & + math_mul33x33, & + math_transpose33, & + math_I3 + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< component-ID of integration point + ip, & !< integration point + el !< element + real(pReal), intent(in), dimension(6) :: & + Tstar_v !< 2nd Piola Kirchhoff stress tensor (Mandel) + real(pReal), intent(in), dimension(3,3) :: & + Fe + integer(pInt) :: & + phase, constituent, instance + real(pReal) :: & + strain(3,3) + + phase = mappingConstitutive(2,ipc,ip,el) + constituent = mappingConstitutive(1,ipc,ip,el) + instance = phase_damageInstance(phase) + + strain = 0.5_pReal*(math_mul33x33(math_transpose33(Fe),Fe)-math_I3) + damageState(phase)%state(2,constituent) = min(damageState(phase)%state0(2,constituent), & + damage_isoBrittle_critStrainEnergy(instance)/ & + sum(abs(math_Mandel6to33(Tstar_v)*strain))) + +end subroutine damage_isoBrittle_microstructure + +!-------------------------------------------------------------------------------------------------- +!> @brief returns damage +!-------------------------------------------------------------------------------------------------- +function damage_isoBrittle_getDamage(ipc, ip, el) + use material, only: & + material_homog, & + mappingHomogenization, & + fieldDamage, & + field_damage_type, & + FIELD_DAMAGE_LOCAL_ID, & + FIELD_DAMAGE_NONLOCAL_ID + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal) :: damage_isoBrittle_getDamage + + select case(field_damage_type(material_homog(ip,el))) + case (FIELD_DAMAGE_LOCAL_ID) + damage_isoBrittle_getDamage = damage_isoBrittle_getLocalDamage(ipc, ip, el) + + case (FIELD_DAMAGE_NONLOCAL_ID) + damage_isoBrittle_getDamage = fieldDamage(material_homog(ip,el))% & + field(1,mappingHomogenization(1,ip,el)) ! Taylor type + + end select + +end function damage_isoBrittle_getDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief returns temperature based on local damage model state layout +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoBrittle_putLocalDamage(ipc, ip, el, localDamage) + use material, only: & + mappingConstitutive, & + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal), intent(in) :: localDamage + + damageState(mappingConstitutive(2,ipc,ip,el))%state(1,mappingConstitutive(1,ipc,ip,el)) = & + localDamage + +end subroutine damage_isoBrittle_putLocalDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief returns local damage +!-------------------------------------------------------------------------------------------------- +function damage_isoBrittle_getLocalDamage(ipc, ip, el) + use material, only: & + mappingConstitutive, & + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal) :: damage_isoBrittle_getLocalDamage + + damage_isoBrittle_getLocalDamage = & + damageState(mappingConstitutive(2,ipc,ip,el))%state(1,mappingConstitutive(1,ipc,ip,el)) + +end function damage_isoBrittle_getLocalDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief returns brittle damage diffusion tensor +!-------------------------------------------------------------------------------------------------- +function damage_isoBrittle_getDamageDiffusion33(ipc, ip, el) + use lattice, only: & + lattice_DamageDiffusion33 + use material, only: & + mappingConstitutive, & + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal), dimension(3,3) :: & + damage_isoBrittle_getDamageDiffusion33 + integer(pInt) :: & + phase, constituent + + phase = mappingConstitutive(2,ipc,ip,el) + constituent = mappingConstitutive(1,ipc,ip,el) + damage_isoBrittle_getDamageDiffusion33 = & + damageState(phase)%state(2,constituent)* & + lattice_DamageDiffusion33(1:3,1:3,phase) + +end function damage_isoBrittle_getDamageDiffusion33 + +!-------------------------------------------------------------------------------------------------- +!> @brief return array of constitutive results +!-------------------------------------------------------------------------------------------------- +function damage_isoBrittle_postResults(ipc,ip,el) + use material, only: & + mappingConstitutive, & + phase_damageInstance,& + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< component-ID of integration point + ip, & !< integration point + el !< element + real(pReal), dimension(damage_isoBrittle_sizePostResults(phase_damageInstance(mappingConstitutive(2,ipc,ip,el)))) :: & + damage_isoBrittle_postResults + + integer(pInt) :: & + instance, phase, constituent, o, c + + phase = mappingConstitutive(2,ipc,ip,el) + constituent = mappingConstitutive(1,ipc,ip,el) + instance = phase_damageInstance(phase) + + c = 0_pInt + damage_isoBrittle_postResults = 0.0_pReal + + do o = 1_pInt,damage_isoBrittle_Noutput(instance) + select case(damage_isoBrittle_outputID(o,instance)) + case (local_damage_ID) + damage_isoBrittle_postResults(c+1_pInt) = damageState(phase)%state(1,constituent) + c = c + 1 + + end select + enddo +end function damage_isoBrittle_postResults + +end module damage_isoBrittle diff --git a/code/damage_isoDuctile.f90 b/code/damage_isoDuctile.f90 new file mode 100644 index 000000000..1b50f43c6 --- /dev/null +++ b/code/damage_isoDuctile.f90 @@ -0,0 +1,450 @@ +!-------------------------------------------------------------------------------------------------- +! $Id$ +!-------------------------------------------------------------------------------------------------- +!> @author Luv Sharma, Max-Planck-Institut für Eisenforschung GmbH +!> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH +!> @brief material subroutine incoprorating ductile damage +!> @details to be done +!-------------------------------------------------------------------------------------------------- +module damage_isoDuctile + use prec, only: & + pReal, & + pInt + + implicit none + private + integer(pInt), dimension(:), allocatable, public, protected :: & + damage_isoDuctile_sizePostResults !< cumulative size of post results + + integer(pInt), dimension(:,:), allocatable, target, public :: & + damage_isoDuctile_sizePostResult !< size of each post result output + + character(len=64), dimension(:,:), allocatable, target, public :: & + damage_isoDuctile_output !< name of each post result output + + integer(pInt), dimension(:), allocatable, target, public :: & + damage_isoDuctile_Noutput !< number of outputs per instance of this damage + + real(pReal), dimension(:), allocatable, private :: & + damage_isoDuctile_aTol, & + damage_isoDuctile_critpStrain + + enum, bind(c) + enumerator :: undefined_ID, & + local_damage_ID + end enum !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!11 ToDo + + integer(kind(undefined_ID)), dimension(:,:), allocatable, private :: & + damage_isoDuctile_outputID !< ID of each post result output + + + public :: & + damage_isoDuctile_init, & + damage_isoDuctile_stateInit, & + damage_isoDuctile_aTolState, & + damage_isoDuctile_dotState, & + damage_isoDuctile_microstructure, & + damage_isoDuctile_getDamage, & + damage_isoDuctile_getSlipDamage, & + damage_isoDuctile_putLocalDamage, & + damage_isoDuctile_getLocalDamage, & + damage_isoDuctile_postResults + +contains + + +!-------------------------------------------------------------------------------------------------- +!> @brief module initialization +!> @details reads in material parameters, allocates arrays, and does sanity checks +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoDuctile_init(fileUnit) + use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment) + use debug, only: & + debug_level,& + debug_constitutive,& + debug_levelBasic + use mesh, only: & + mesh_maxNips, & + mesh_NcpElems + use IO, only: & + IO_read, & + IO_lc, & + IO_getTag, & + IO_isBlank, & + IO_stringPos, & + IO_stringValue, & + IO_floatValue, & + IO_intValue, & + IO_warning, & + IO_error, & + IO_timeStamp, & + IO_EOF + use material, only: & + homogenization_maxNgrains, & + phase_damage, & + phase_damageInstance, & + phase_Noutput, & + LOCAL_damage_isoDuctile_label, & + LOCAL_damage_isoDuctile_ID, & + material_phase, & + damageState, & + MATERIAL_partPhase + use numerics,only: & + worldrank, & + numerics_integrator + + implicit none + integer(pInt), intent(in) :: fileUnit + + integer(pInt), parameter :: MAXNCHUNKS = 7_pInt + integer(pInt), dimension(1+2*MAXNCHUNKS) :: positions + integer(pInt) :: maxNinstance,mySize=0_pInt,phase,instance,o + integer(pInt) :: sizeState, sizeDotState + integer(pInt) :: NofMyPhase + character(len=65536) :: & + tag = '', & + line = '' + + mainProcess: if (worldrank == 0) then + write(6,'(/,a)') ' <<<+- damage_'//LOCAL_damage_isoDuctile_LABEL//' init -+>>>' + write(6,'(a)') ' $Id$' + write(6,'(a15,a)') ' Current time: ',IO_timeStamp() +#include "compilation_info.f90" + endif mainProcess + + maxNinstance = int(count(phase_damage == LOCAL_damage_isoDuctile_ID),pInt) + if (maxNinstance == 0_pInt) return + + if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) & + write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance + + allocate(damage_isoDuctile_sizePostResults(maxNinstance), source=0_pInt) + allocate(damage_isoDuctile_sizePostResult(maxval(phase_Noutput),maxNinstance),source=0_pInt) + allocate(damage_isoDuctile_output(maxval(phase_Noutput),maxNinstance)) + damage_isoDuctile_output = '' + allocate(damage_isoDuctile_outputID(maxval(phase_Noutput),maxNinstance), source=undefined_ID) + allocate(damage_isoDuctile_Noutput(maxNinstance), source=0_pInt) + allocate(damage_isoDuctile_critpStrain(maxNinstance), source=0.0_pReal) + allocate(damage_isoDuctile_aTol(maxNinstance), source=0.0_pReal) + + rewind(fileUnit) + phase = 0_pInt + do while (trim(line) /= IO_EOF .and. IO_lc(IO_getTag(line,'<','>')) /= MATERIAL_partPhase) ! wind forward to + line = IO_read(fileUnit) + enddo + + parsingFile: do while (trim(line) /= IO_EOF) ! read through sections of phase part + line = IO_read(fileUnit) + if (IO_isBlank(line)) cycle ! skip empty lines + if (IO_getTag(line,'<','>') /= '') then ! stop at next part + line = IO_read(fileUnit, .true.) ! reset IO_read + exit + endif + if (IO_getTag(line,'[',']') /= '') then ! next phase section + phase = phase + 1_pInt ! advance phase section counter + cycle ! skip to next line + endif + if (phase > 0_pInt ) then; if (phase_damage(phase) == LOCAL_damage_isoDuctile_ID) then ! do not short-circuit here (.and. with next if statemen). It's not safe in Fortran + instance = phase_damageInstance(phase) ! which instance of my damage is present phase + positions = IO_stringPos(line,MAXNCHUNKS) + tag = IO_lc(IO_stringValue(line,positions,1_pInt)) ! extract key + select case(tag) + case ('(output)') + select case(IO_lc(IO_stringValue(line,positions,2_pInt))) + case ('local_damage') + damage_isoDuctile_Noutput(instance) = damage_isoDuctile_Noutput(instance) + 1_pInt + damage_isoDuctile_outputID(damage_isoDuctile_Noutput(instance),instance) = local_damage_ID + damage_isoDuctile_output(damage_isoDuctile_Noutput(instance),instance) = & + IO_lc(IO_stringValue(line,positions,2_pInt)) + end select + + case ('critical_plastic_strain') + damage_isoDuctile_critpStrain(instance) = IO_floatValue(line,positions,2_pInt) + + case ('atol_damage') + damage_isoDuctile_aTol(instance) = IO_floatValue(line,positions,2_pInt) + + end select + endif; endif + enddo parsingFile + + initializeInstances: do phase = 1_pInt, size(phase_damage) + if (phase_damage(phase) == LOCAL_damage_isoDuctile_ID) then + NofMyPhase=count(material_phase==phase) + instance = phase_damageInstance(phase) + +!-------------------------------------------------------------------------------------------------- +! Determine size of postResults array + outputsLoop: do o = 1_pInt,damage_isoDuctile_Noutput(instance) + select case(damage_isoDuctile_outputID(o,instance)) + case(local_damage_ID) + mySize = 1_pInt + end select + + if (mySize > 0_pInt) then ! any meaningful output found + damage_isoDuctile_sizePostResult(o,instance) = mySize + damage_isoDuctile_sizePostResults(instance) = damage_isoDuctile_sizePostResults(instance) + mySize + endif + enddo outputsLoop +! Determine size of state array + sizeDotState = 1_pInt + sizeState = 2_pInt + + damageState(phase)%sizeState = sizeState + damageState(phase)%sizeDotState = sizeDotState + damageState(phase)%sizePostResults = damage_isoDuctile_sizePostResults(instance) + allocate(damageState(phase)%aTolState (sizeState), source=0.0_pReal) + allocate(damageState(phase)%state0 (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%partionedState0 (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%subState0 (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%state (sizeState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%state_backup (sizeState,NofMyPhase), source=0.0_pReal) + + allocate(damageState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%deltaState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal) + if (any(numerics_integrator == 1_pInt)) then + allocate(damageState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(damageState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal) + endif + if (any(numerics_integrator == 4_pInt)) & + allocate(damageState(phase)%RK4dotState (sizeDotState,NofMyPhase), source=0.0_pReal) + if (any(numerics_integrator == 5_pInt)) & + allocate(damageState(phase)%RKCK45dotState (6,sizeDotState,NofMyPhase),source=0.0_pReal) + + call damage_isoDuctile_stateInit(phase) + call damage_isoDuctile_aTolState(phase,instance) + endif + + enddo initializeInstances +end subroutine damage_isoDuctile_init + +!-------------------------------------------------------------------------------------------------- +!> @brief sets the relevant state values for a given instance of this damage +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoDuctile_stateInit(phase) + use material, only: & + damageState + + implicit none + integer(pInt), intent(in) :: phase !< number specifying the phase of the damage + + real(pReal), dimension(damageState(phase)%sizeState) :: tempState + + tempState(1) = 1.0_pReal + tempState(2) = 1.0_pReal + + + damageState(phase)%state = spread(tempState,2,size(damageState(phase)%state(1,:))) + damageState(phase)%state0 = damageState(phase)%state + damageState(phase)%partionedState0 = damageState(phase)%state +end subroutine damage_isoDuctile_stateInit + +!-------------------------------------------------------------------------------------------------- +!> @brief sets the relevant state values for a given instance of this damage +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoDuctile_aTolState(phase,instance) + use material, only: & + damageState + + implicit none + integer(pInt), intent(in) :: & + phase, & + instance ! number specifying the current instance of the damage + real(pReal), dimension(damageState(phase)%sizeState) :: tempTol + + tempTol = damage_isoDuctile_aTol(instance) + damageState(phase)%aTolState = tempTol +end subroutine damage_isoDuctile_aTolState + +!-------------------------------------------------------------------------------------------------- +!> @brief calculates derived quantities from state +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoDuctile_dotState(ipc, ip, el) + use material, only: & + mappingConstitutive, & + damageState + use math, only: & + math_norm33 + use lattice, only: & + lattice_DamageMobility + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< component-ID of integration point + ip, & !< integration point + el !< element + integer(pInt) :: & + phase, constituent + + phase = mappingConstitutive(2,ipc,ip,el) + constituent = mappingConstitutive(1,ipc,ip,el) + + damageState(phase)%dotState(1,constituent) = & + (1.0_pReal/lattice_DamageMobility(phase))* & + (damageState(phase)%state(2,constituent) - & + damageState(phase)%state(1,constituent)) + +end subroutine damage_isoDuctile_dotState + +!-------------------------------------------------------------------------------------------------- +!> @brief calculates derived quantities from state +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoDuctile_microstructure(nSlip,accumulatedSlip,ipc, ip, el) + use material, only: & + mappingConstitutive, & + phase_damageInstance, & + damageState + use math, only: & + math_Mandel6to33, & + math_mul33x33, & + math_transpose33, & + math_I3, & + math_norm33 + + implicit none + integer(pInt), intent(in) :: & + nSlip, & + ipc, & !< component-ID of integration point + ip, & !< integration point + el !< element + real(pReal), dimension(nSlip), intent(in) :: & + accumulatedSlip + integer(pInt) :: & + phase, constituent, instance + + phase = mappingConstitutive(2,ipc,ip,el) + constituent = mappingConstitutive(1,ipc,ip,el) + instance = phase_damageInstance(phase) + damageState(phase)%state(2,constituent) = min(damageState(phase)%state(2,constituent), & + damage_isoDuctile_critpStrain(instance)/ & + sum(accumulatedSlip)) !< akin to damage surface + +end subroutine damage_isoDuctile_microstructure + +!-------------------------------------------------------------------------------------------------- +!> @brief returns damage +!-------------------------------------------------------------------------------------------------- +function damage_isoDuctile_getDamage(ipc, ip, el) + use material, only: & + material_homog, & + mappingHomogenization, & + fieldDamage, & + field_damage_type, & + FIELD_DAMAGE_LOCAL_ID, & + FIELD_DAMAGE_NONLOCAL_ID + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal) :: damage_isoDuctile_getDamage + + select case(field_damage_type(material_homog(ip,el))) + case (FIELD_DAMAGE_LOCAL_ID) + damage_isoDuctile_getDamage = damage_isoDuctile_getLocalDamage(ipc, ip, el) + + case (FIELD_DAMAGE_NONLOCAL_ID) + damage_isoDuctile_getDamage = fieldDamage(material_homog(ip,el))% & + field(1,mappingHomogenization(1,ip,el)) ! Taylor type + + end select + +end function damage_isoDuctile_getDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief returns slip damage +!-------------------------------------------------------------------------------------------------- +function damage_isoDuctile_getSlipDamage(ipc, ip, el) + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal) :: damage_isoDuctile_getSlipDamage, damage + + damage = damage_isoDuctile_getDamage(ipc, ip, el) + damage_isoDuctile_getSlipDamage = damage*damage + +end function damage_isoDuctile_getSlipDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief puts local damage +!-------------------------------------------------------------------------------------------------- +subroutine damage_isoDuctile_putLocalDamage(ipc, ip, el, localDamage) + use material, only: & + mappingConstitutive, & + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal), intent(in) :: localDamage + + damageState(mappingConstitutive(2,ipc,ip,el))%state(1,mappingConstitutive(1,ipc,ip,el)) = & + localDamage + +end subroutine damage_isoDuctile_putLocalDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief returns local damage +!-------------------------------------------------------------------------------------------------- +function damage_isoDuctile_getLocalDamage(ipc, ip, el) + use material, only: & + mappingConstitutive, & + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal) :: damage_isoDuctile_getLocalDamage + + damage_isoDuctile_getLocalDamage = & + damageState(mappingConstitutive(2,ipc,ip,el))%state(1,mappingConstitutive(1,ipc,ip,el)) + +end function damage_isoDuctile_getLocalDamage + +!-------------------------------------------------------------------------------------------------- +!> @brief return array of constitutive results +!-------------------------------------------------------------------------------------------------- +function damage_isoDuctile_postResults(ipc,ip,el) + use material, only: & + mappingConstitutive, & + phase_damageInstance,& + damageState + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< component-ID of integration point + ip, & !< integration point + el !< element + real(pReal), dimension(damage_isoDuctile_sizePostResults(phase_damageInstance(mappingConstitutive(2,ipc,ip,el)))) :: & + damage_isoDuctile_postResults + + integer(pInt) :: & + instance, phase, constituent, o, c + + phase = mappingConstitutive(2,ipc,ip,el) + constituent = mappingConstitutive(1,ipc,ip,el) + instance = phase_damageInstance(phase) + + c = 0_pInt + damage_isoDuctile_postResults = 0.0_pReal + + do o = 1_pInt,damage_isoDuctile_Noutput(instance) + select case(damage_isoDuctile_outputID(o,instance)) + case (local_damage_ID) + damage_isoDuctile_postResults(c+1_pInt) = damageState(phase)%state(1,constituent) + c = c + 1 + + end select + enddo +end function damage_isoDuctile_postResults + +end module damage_isoDuctile diff --git a/code/material.f90 b/code/material.f90 index 4c5d803e7..94dadb864 100644 --- a/code/material.f90 +++ b/code/material.f90 @@ -27,11 +27,11 @@ module material PLASTICITY_dislokmc_label = 'dislokmc', & PLASTICITY_titanmod_label = 'titanmod', & PLASTICITY_nonlocal_label = 'nonlocal', & - LOCAL_DAMAGE_none_label = 'none', & - LOCAL_DAMAGE_brittle_label = 'brittle', & - LOCAL_DAMAGE_ductile_label = 'ductile', & - LOCAL_DAMAGE_gurson_label = 'gurson', & - LOCAL_DAMAGE_anisotropic_label = 'anisotropic', & + LOCAL_DAMAGE_none_LABEL = 'none', & + LOCAL_DAMAGE_isoBrittle_LABEL = 'isoBrittle', & + LOCAL_DAMAGE_isoDuctile_LABEL = 'isoDuctile', & + LOCAL_DAMAGE_anisoBrittle_LABEL= 'anisoBrittle', & + LOCAL_DAMAGE_gurson_LABEL = 'gurson', & LOCAL_THERMAL_isothermal_label = 'isothermal', & LOCAL_THERMAL_adiabatic_label = 'adiabatic', & LOCAL_VACANCY_constant_label = 'constant', & @@ -64,10 +64,10 @@ module material end enum enum, bind(c) enumerator :: LOCAL_DAMAGE_none_ID, & - LOCAL_DAMAGE_brittle_ID, & - LOCAL_DAMAGE_ductile_ID, & - LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & + LOCAL_DAMAGE_gurson_ID end enum enum, bind(c) enumerator :: LOCAL_THERMAL_isothermal_ID, & @@ -237,10 +237,10 @@ module material PLASTICITY_titanmod_ID, & PLASTICITY_nonlocal_ID, & LOCAL_DAMAGE_none_ID, & - LOCAL_DAMAGE_brittle_ID, & - LOCAL_DAMAGE_ductile_ID, & + LOCAL_DAMAGE_isoBrittle_ID, & + LOCAL_DAMAGE_isoDuctile_ID, & + LOCAL_DAMAGE_anisoBrittle_ID, & LOCAL_DAMAGE_gurson_ID, & - LOCAL_DAMAGE_anisotropic_ID, & LOCAL_THERMAL_isothermal_ID, & LOCAL_THERMAL_adiabatic_ID, & LOCAL_VACANCY_constant_ID, & @@ -846,14 +846,14 @@ subroutine material_parsePhase(fileUnit,myPart) select case (IO_lc(IO_stringValue(line,positions,2_pInt))) case (LOCAL_DAMAGE_none_label) phase_damage(section) = LOCAL_DAMAGE_none_ID - case (LOCAL_DAMAGE_brittle_label) - phase_damage(section) = LOCAL_DAMAGE_BRITTLE_ID - case (LOCAL_DAMAGE_ductile_label) - phase_damage(section) = LOCAL_DAMAGE_DUCTILE_ID + case (LOCAL_DAMAGE_isoBrittle_label) + phase_damage(section) = LOCAL_DAMAGE_isoBrittle_ID + case (LOCAL_DAMAGE_isoDuctile_label) + phase_damage(section) = LOCAL_DAMAGE_isoDuctile_ID + case (LOCAL_DAMAGE_anisoBrittle_label) + phase_damage(section) = LOCAL_DAMAGE_anisoBrittle_ID case (LOCAL_DAMAGE_gurson_label) phase_damage(section) = LOCAL_DAMAGE_gurson_ID - case (LOCAL_DAMAGE_anisotropic_label) - phase_damage(section) = LOCAL_DAMAGE_anisotropic_ID case default call IO_error(200_pInt,ext_msg=trim(IO_stringValue(line,positions,2_pInt))) end select diff --git a/code/vacancy_generation.f90 b/code/vacancy_generation.f90 index e040ffead..682db1f6f 100644 --- a/code/vacancy_generation.f90 +++ b/code/vacancy_generation.f90 @@ -215,7 +215,7 @@ subroutine vacancy_generation_init(fileUnit) allocate(vacancyState(phase)%state_backup (sizeState,NofMyPhase), source=0.0_pReal) allocate(vacancyState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal) - allocate(vacancyState(phase)%deltaState (sizeDotState,NofMyPhase), source=0.0_pReal) + allocate(vacancyState(phase)%deltaState (sizeDotState,NofMyPhase), source=0.0_pReal) allocate(vacancyState(phase)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal) if (any(numerics_integrator == 1_pInt)) then allocate(vacancyState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)