diff --git a/src/constitutive.f90 b/src/constitutive.f90 index e2c9dbc05..41baab7ab 100644 --- a/src/constitutive.f90 +++ b/src/constitutive.f90 @@ -35,29 +35,44 @@ module constitutive interface - module subroutine plastic_none_init + module subroutine plastic_none_init(debug_constitutive) + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_none_init - module subroutine plastic_isotropic_init + module subroutine plastic_isotropic_init(debug_constitutive) + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_isotropic_init - - module subroutine plastic_phenopowerlaw_init + + module subroutine plastic_phenopowerlaw_init(debug_constitutive) + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_phenopowerlaw_init - module subroutine plastic_kinehardening_init + module subroutine plastic_kinehardening_init(debug_constitutive) + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_kinehardening_init - module subroutine plastic_dislotwin_init + module subroutine plastic_dislotwin_init(debug_constitutive) + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_dislotwin_init - module subroutine plastic_disloUCLA_init + module subroutine plastic_disloUCLA_init(debug_constitutive) + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_disloUCLA_init - module subroutine plastic_nonlocal_init + module subroutine plastic_nonlocal_init(debug_constitutive) + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_nonlocal_init - module subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of) + module subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of, & + debug_constitutive) real(pReal), dimension(3,3), intent(out) :: & Lp !< plastic velocity gradient real(pReal), dimension(3,3,3,3), intent(out) :: & @@ -68,6 +83,8 @@ module constitutive integer, intent(in) :: & instance, & of + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_isotropic_LpAndItsTangent pure module subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of) @@ -145,7 +162,8 @@ module constitutive end subroutine plastic_nonlocal_LpAndItsTangent - module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,of) + module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,of, & + debug_constitutive) real(pReal), dimension(3,3), intent(out) :: & Li !< inleastic velocity gradient real(pReal), dimension(3,3,3,3), intent(out) :: & @@ -156,6 +174,8 @@ module constitutive integer, intent(in) :: & instance, & of + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_isotropic_LiAndItsTangent @@ -204,7 +224,7 @@ module constitutive end subroutine plastic_disloUCLA_dotState module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, & - instance,of,ip,el) + instance,of,ip,el,debug_constitutive) real(pReal), dimension(3,3), intent(in) ::& Mp !< MandelStress real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: & @@ -218,6 +238,8 @@ module constitutive of, & ip, & !< current integration point el !< current element number + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_nonlocal_dotState @@ -235,7 +257,8 @@ module constitutive of end subroutine plastic_disloUCLA_dependentState - module subroutine plastic_nonlocal_dependentState(F, Fp, instance, of, ip, el) + module subroutine plastic_nonlocal_dependentState(F, Fp, instance, of, ip, el, & + debug_constitutive) real(pReal), dimension(3,3), intent(in) :: & F, & Fp @@ -244,18 +267,22 @@ module constitutive of, & ip, & el + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_nonlocal_dependentState - module subroutine plastic_kinehardening_deltaState(Mp,instance,of) + module subroutine plastic_kinehardening_deltaState(Mp,instance,of,debug_constitutive) real(pReal), dimension(3,3), intent(in) :: & Mp !< Mandel stress integer, intent(in) :: & instance, & of + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_kinehardening_deltaState - module subroutine plastic_nonlocal_deltaState(Mp,instance,of,ip,el) + module subroutine plastic_nonlocal_deltaState(Mp,instance,of,ip,el,debug_constitutive) real(pReal), dimension(3,3), intent(in) :: & Mp integer, intent(in) :: & @@ -263,6 +290,8 @@ module constitutive of, & ip, & el + class(tNode), pointer , intent(in) :: & + debug_constitutive end subroutine plastic_nonlocal_deltaState @@ -341,34 +370,37 @@ subroutine constitutive_init integer :: & ph, & !< counter in phase loop s !< counter in source loop + class(tNode), pointer :: & + debug_constitutive + debug_constitutive => debug_root%get('constitutuve',defaultVal=emptyList) !-------------------------------------------------------------------------------------------------- ! initialized plasticity - if (any(phase_plasticity == PLASTICITY_NONE_ID)) call plastic_none_init - if (any(phase_plasticity == PLASTICITY_ISOTROPIC_ID)) call plastic_isotropic_init - if (any(phase_plasticity == PLASTICITY_PHENOPOWERLAW_ID)) call plastic_phenopowerlaw_init - if (any(phase_plasticity == PLASTICITY_KINEHARDENING_ID)) call plastic_kinehardening_init - if (any(phase_plasticity == PLASTICITY_DISLOTWIN_ID)) call plastic_dislotwin_init - if (any(phase_plasticity == PLASTICITY_DISLOUCLA_ID)) call plastic_disloucla_init + if (any(phase_plasticity == PLASTICITY_NONE_ID)) call plastic_none_init(debug_constitutive) + if (any(phase_plasticity == PLASTICITY_ISOTROPIC_ID)) call plastic_isotropic_init(debug_constitutive) + if (any(phase_plasticity == PLASTICITY_PHENOPOWERLAW_ID)) call plastic_phenopowerlaw_init(debug_constitutive) + if (any(phase_plasticity == PLASTICITY_KINEHARDENING_ID)) call plastic_kinehardening_init(debug_constitutive) + if (any(phase_plasticity == PLASTICITY_DISLOTWIN_ID)) call plastic_dislotwin_init(debug_constitutive) + if (any(phase_plasticity == PLASTICITY_DISLOUCLA_ID)) call plastic_disloucla_init(debug_constitutive) if (any(phase_plasticity == PLASTICITY_NONLOCAL_ID)) then - call plastic_nonlocal_init + call plastic_nonlocal_init(debug_constitutive) else call geometry_plastic_nonlocal_disable endif !-------------------------------------------------------------------------------------------------- ! initialize source mechanisms - if (any(phase_source == SOURCE_thermal_dissipation_ID)) call source_thermal_dissipation_init - if (any(phase_source == SOURCE_thermal_externalheat_ID)) call source_thermal_externalheat_init - if (any(phase_source == SOURCE_damage_isoBrittle_ID)) call source_damage_isoBrittle_init - if (any(phase_source == SOURCE_damage_isoDuctile_ID)) call source_damage_isoDuctile_init - if (any(phase_source == SOURCE_damage_anisoBrittle_ID)) call source_damage_anisoBrittle_init - if (any(phase_source == SOURCE_damage_anisoDuctile_ID)) call source_damage_anisoDuctile_init + if (any(phase_source == SOURCE_thermal_dissipation_ID)) call source_thermal_dissipation_init(debug_constitutive) + if (any(phase_source == SOURCE_thermal_externalheat_ID)) call source_thermal_externalheat_init(debug_constitutive) + if (any(phase_source == SOURCE_damage_isoBrittle_ID)) call source_damage_isoBrittle_init(debug_constitutive) + if (any(phase_source == SOURCE_damage_isoDuctile_ID)) call source_damage_isoDuctile_init(debug_constitutive) + if (any(phase_source == SOURCE_damage_anisoBrittle_ID)) call source_damage_anisoBrittle_init(debug_constitutive) + if (any(phase_source == SOURCE_damage_anisoDuctile_ID)) call source_damage_anisoDuctile_init(debug_constitutive) !-------------------------------------------------------------------------------------------------- ! initialize kinematic mechanisms - if (any(phase_kinematics == KINEMATICS_cleavage_opening_ID)) call kinematics_cleavage_opening_init - if (any(phase_kinematics == KINEMATICS_slipplane_opening_ID)) call kinematics_slipplane_opening_init - if (any(phase_kinematics == KINEMATICS_thermal_expansion_ID)) call kinematics_thermal_expansion_init + if (any(phase_kinematics == KINEMATICS_cleavage_opening_ID)) call kinematics_cleavage_opening_init(debug_constitutive) + if (any(phase_kinematics == KINEMATICS_slipplane_opening_ID)) call kinematics_slipplane_opening_init(debug_constitutive) + if (any(phase_kinematics == KINEMATICS_thermal_expansion_ID)) call kinematics_thermal_expansion_init(debug_constitutive) write(6,'(/,a)') ' <<<+- constitutive init -+>>>'; flush(6) @@ -430,7 +462,10 @@ subroutine constitutive_dependentState(F, Fp, ipc, ip, el) ho, & !< homogenization tme, & !< thermal member position instance, of + class(tNode), pointer :: & + debug_constitutive + debug_constitutive => debug_root%get('constitutive',defaultVal=emptyList) ho = material_homogenizationAt(el) tme = thermalMapping(ho)%p(ip,el) of = material_phasememberAt(ipc,ip,el) @@ -442,7 +477,7 @@ subroutine constitutive_dependentState(F, Fp, ipc, ip, el) case (PLASTICITY_DISLOUCLA_ID) plasticityType call plastic_disloUCLA_dependentState(instance,of) case (PLASTICITY_NONLOCAL_ID) plasticityType - call plastic_nonlocal_dependentState (F,Fp,instance,of,ip,el) + call plastic_nonlocal_dependentState (F,Fp,instance,of,ip,el,debug_constitutive) end select plasticityType end subroutine constitutive_dependentState @@ -476,7 +511,11 @@ subroutine constitutive_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, & tme !< thermal member position integer :: & i, j, instance, of - + class(tNode), pointer :: & + debug_constitutive + + debug_constitutive => debug_root%get('constitutive',defaultVal=emptyList) + ho = material_homogenizationAt(el) tme = thermalMapping(ho)%p(ip,el) @@ -491,7 +530,7 @@ subroutine constitutive_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, & dLp_dMp = 0.0_pReal case (PLASTICITY_ISOTROPIC_ID) plasticityType - call plastic_isotropic_LpAndItsTangent (Lp,dLp_dMp,Mp,instance,of) + call plastic_isotropic_LpAndItsTangent (Lp,dLp_dMp,Mp,instance,of,debug_constitutive) case (PLASTICITY_PHENOPOWERLAW_ID) plasticityType call plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of) @@ -551,6 +590,10 @@ subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, & integer :: & k, i, j, & instance, of + class(tNode), pointer :: & + debug_constitutive + + debug_constitutive => debug_root%get('constitutive',defaultVal=emptyList) Li = 0.0_pReal dLi_dS = 0.0_pReal @@ -560,7 +603,7 @@ subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, & case (PLASTICITY_isotropic_ID) plasticityType of = material_phasememberAt(ipc,ip,el) instance = phase_plasticityInstance(material_phaseAt(ipc,el)) - call plastic_isotropic_LiAndItsTangent(my_Li, my_dLi_dS, S ,instance,of) + call plastic_isotropic_LiAndItsTangent(my_Li, my_dLi_dS, S ,instance,of,debug_constitutive) case default plasticityType my_Li = 0.0_pReal my_dLi_dS = 0.0_pReal @@ -733,8 +776,12 @@ function constitutive_collectDotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el tme, & !< thermal member position i, & !< counter in source loop instance + class(tNode), pointer :: & + debug_constitutive logical :: broken + debug_constitutive => debug_root%get('constitutive',defaultVal=emptyList) + ho = material_homogenizationAt(el) tme = thermalMapping(ho)%p(ip,el) instance = phase_plasticityInstance(phase) @@ -760,7 +807,7 @@ function constitutive_collectDotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el case (PLASTICITY_NONLOCAL_ID) plasticityType call plastic_nonlocal_dotState (Mp,FArray,FpArray,temperature(ho)%p(tme),subdt, & - instance,of,ip,el) + instance,of,ip,el,debug_constitutive) end select plasticityType broken = any(IEEE_is_NaN(plasticState(phase)%dotState(:,of))) @@ -812,20 +859,24 @@ function constitutive_deltaState(S, Fe, Fi, ipc, ip, el, phase, of) result(broke instance, & myOffset, & mySize + class(tNode), pointer :: & + debug_constitutive logical :: & broken + debug_constitutive => debug_root%get('constitutive',defaultVal=emptyList) + Mp = matmul(matmul(transpose(Fi),Fi),S) instance = phase_plasticityInstance(phase) plasticityType: select case (phase_plasticity(phase)) case (PLASTICITY_KINEHARDENING_ID) plasticityType - call plastic_kinehardening_deltaState(Mp,instance,of) + call plastic_kinehardening_deltaState(Mp,instance,of,debug_constitutive) broken = any(IEEE_is_NaN(plasticState(phase)%deltaState(:,of))) case (PLASTICITY_NONLOCAL_ID) plasticityType - call plastic_nonlocal_deltaState(Mp,instance,of,ip,el) + call plastic_nonlocal_deltaState(Mp,instance,of,ip,el,debug_constitutive) broken = any(IEEE_is_NaN(plasticState(phase)%deltaState(:,of))) case default diff --git a/src/constitutive_plastic_disloUCLA.f90 b/src/constitutive_plastic_disloUCLA.f90 index 90a933910..3a2946435 100644 --- a/src/constitutive_plastic_disloUCLA.f90 +++ b/src/constitutive_plastic_disloUCLA.f90 @@ -74,7 +74,10 @@ contains !> @brief Perform module initialization. !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -module subroutine plastic_disloUCLA_init +module subroutine plastic_disloUCLA_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & Ninstance, & @@ -97,7 +100,7 @@ module subroutine plastic_disloUCLA_init write(6,'(a)') ' https://dx.doi.org/10.1016/j.ijplas.2015.09.002' Ninstance = count(phase_plasticity == PLASTICITY_DISLOUCLA_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(param(Ninstance)) diff --git a/src/constitutive_plastic_dislotwin.f90 b/src/constitutive_plastic_dislotwin.f90 index 7c7d24ab8..c643376c7 100644 --- a/src/constitutive_plastic_dislotwin.f90 +++ b/src/constitutive_plastic_dislotwin.f90 @@ -122,7 +122,10 @@ contains !> @brief Perform module initialization. !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -module subroutine plastic_dislotwin_init +module subroutine plastic_dislotwin_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & Ninstance, & @@ -151,7 +154,7 @@ module subroutine plastic_dislotwin_init Ninstance = count(phase_plasticity == PLASTICITY_DISLOTWIN_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(param(Ninstance)) diff --git a/src/constitutive_plastic_isotropic.f90 b/src/constitutive_plastic_isotropic.f90 index ecf029124..7e5736d31 100644 --- a/src/constitutive_plastic_isotropic.f90 +++ b/src/constitutive_plastic_isotropic.f90 @@ -49,13 +49,19 @@ contains !> @brief Perform module initialization. !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -module subroutine plastic_isotropic_init +module subroutine plastic_isotropic_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & Ninstance, & p, & NipcMyPhase, & - sizeState, sizeDotState + sizeState, sizeDotState, & + debug_g, & + debug_e, & + debug_i real(pReal) :: & xi_0 !< initial critical stress character(len=pStringLen) :: & @@ -67,7 +73,7 @@ module subroutine plastic_isotropic_init write(6,'(a)') ' https://doi.org/10.1016/j.scriptamat.2017.09.047' Ninstance = count(phase_plasticity == PLASTICITY_ISOTROPIC_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(param(Ninstance)) @@ -84,6 +90,10 @@ module subroutine plastic_isotropic_init prm%output = config%getStrings('(output)',defaultVal=emptyStringArray) #ifdef DEBUG + debug_g = debug_root%get_asInt('grain',defaultVal=1) + debug_e = debug_root%get_asInt('element',defaultVal=1) + debug_i = debug_root%get_asInt('integrationpoint',defaultVal=1) + if (p==material_phaseAt(debug_g,debug_e)) & prm%of_debug = material_phasememberAt(debug_g,debug_i,debug_e) #endif @@ -150,7 +160,7 @@ end subroutine plastic_isotropic_init !-------------------------------------------------------------------------------------------------- !> @brief Calculate plastic velocity gradient and its tangent. !-------------------------------------------------------------------------------------------------- -module subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of) +module subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of,debug_constitutive) real(pReal), dimension(3,3), intent(out) :: & Lp !< plastic velocity gradient @@ -162,7 +172,9 @@ module subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of) integer, intent(in) :: & instance, & of - + class(tNode), pointer, intent(in) :: & + debug_constitutive + real(pReal), dimension(3,3) :: & Mp_dev !< deviatoric part of the Mandel stress real(pReal) :: & @@ -183,8 +195,8 @@ module subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of) Lp = dot_gamma/prm%M * Mp_dev/norm_Mp_dev #ifdef DEBUG - if (iand(debug_level(debug_constitutive), debug_levelExtensive) /= 0 & - .and. (of == prm%of_debug .or. .not. iand(debug_level(debug_constitutive),debug_levelSelective) /= 0)) then + if (debug_constitutive%contains('extensive') & + .and. (of == prm%of_debug .or. .not. debug_constitutive%contains('selective'))) then write(6,'(/,a,/,3(12x,3(f12.4,1x)/))') '<< CONST isotropic >> Tstar (dev) / MPa', & transpose(Mp_dev)*1.0e-6_pReal write(6,'(/,a,/,f12.5)') '<< CONST isotropic >> norm Tstar / MPa', norm_Mp_dev*1.0e-6_pReal @@ -211,7 +223,7 @@ end subroutine plastic_isotropic_LpAndItsTangent !-------------------------------------------------------------------------------------------------- !> @brief Calculate inelastic velocity gradient and its tangent. !-------------------------------------------------------------------------------------------------- -module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,of) +module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,of,debug_constitutive) real(pReal), dimension(3,3), intent(out) :: & Li !< inleastic velocity gradient @@ -223,7 +235,9 @@ module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,of) integer, intent(in) :: & instance, & of - + class(tNode), pointer, intent(in) :: & + debug_constitutive + real(pReal) :: & tr !< trace of spherical part of Mandel stress (= 3 x pressure) integer :: & @@ -239,8 +253,8 @@ module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,of) * tr * abs(tr)**(prm%n-1.0_pReal) #ifdef DEBUG - if (iand(debug_level(debug_constitutive), debug_levelExtensive) /= 0 & - .and. (of == prm%of_debug .or. .not. iand(debug_level(debug_constitutive),debug_levelSelective) /= 0)) then + if (debug_constitutive%contains('extensive') & + .and. (of == prm%of_debug .or. .not. debug_constitutive%contains('selective'))) then write(6,'(/,a,/,f12.5)') '<< CONST isotropic >> pressure / MPa', tr/3.0_pReal*1.0e-6_pReal write(6,'(/,a,/,f12.5)') '<< CONST isotropic >> gdot', prm%dot_gamma_0 * (3.0_pReal*prm%M*stt%xi(of))**(-prm%n) & * tr * abs(tr)**(prm%n-1.0_pReal) diff --git a/src/constitutive_plastic_kinehardening.f90 b/src/constitutive_plastic_kinehardening.f90 index 36b1eedf9..6ff250057 100644 --- a/src/constitutive_plastic_kinehardening.f90 +++ b/src/constitutive_plastic_kinehardening.f90 @@ -58,14 +58,18 @@ contains !> @brief Perform module initialization. !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -module subroutine plastic_kinehardening_init +module subroutine plastic_kinehardening_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & Ninstance, & p, o, & NipcMyPhase, & sizeState, sizeDeltaState, sizeDotState, & - startIndex, endIndex + startIndex, endIndex, & + debug_e, debug_i, debug_g integer, dimension(:), allocatable :: & N_sl real(pReal), dimension(:), allocatable :: & @@ -77,7 +81,7 @@ module subroutine plastic_kinehardening_init write(6,'(/,a)') ' <<<+- plastic_'//PLASTICITY_KINEHARDENING_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_plasticity == PLASTICITY_KINEHARDENING_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(param(Ninstance)) @@ -96,6 +100,10 @@ module subroutine plastic_kinehardening_init prm%output = config%getStrings('(output)',defaultVal=emptyStringArray) #ifdef DEBUG + debug_g = debug_root%get_asInt('grain',defaultVal=1) + debug_e = debug_root%get_asInt('element',defaultVal=1) + debug_i = debug_root%get_asInt('integrationpoint',defaultVal=1) + if (p==material_phaseAt(debug_g,debug_e)) then prm%of_debug = material_phasememberAt(debug_g,debug_i,debug_e) endif @@ -308,13 +316,15 @@ end subroutine plastic_kinehardening_dotState !-------------------------------------------------------------------------------------------------- !> @brief Calculate (instantaneous) incremental change of microstructure. !-------------------------------------------------------------------------------------------------- -module subroutine plastic_kinehardening_deltaState(Mp,instance,of) +module subroutine plastic_kinehardening_deltaState(Mp,instance,of,debug_constitutive) real(pReal), dimension(3,3), intent(in) :: & Mp !< Mandel stress integer, intent(in) :: & instance, & of + class(tNode), pointer , intent(in) :: & + debug_constitutive real(pReal), dimension(param(instance)%sum_N_sl) :: & gdot_pos,gdot_neg, & @@ -328,9 +338,9 @@ module subroutine plastic_kinehardening_deltaState(Mp,instance,of) dEq0(gdot_pos+gdot_neg,1e-10_pReal)) ! current sense of shear direction #ifdef DEBUG - if (iand(debug_level(debug_constitutive), debug_levelExtensive) /= 0 & + if (debug_constitutive%contains('extensive') & .and. (of == prm%of_debug & - .or. .not. iand(debug_level(debug_constitutive),debug_levelSelective) /= 0)) then + .or. .not. debug_constitutive%contains('selective'))) then write(6,'(a)') '======= kinehardening delta state =======' write(6,*) sense,state(instance)%sense(:,of) endif diff --git a/src/constitutive_plastic_none.f90 b/src/constitutive_plastic_none.f90 index 667fe5638..d79b08d93 100644 --- a/src/constitutive_plastic_none.f90 +++ b/src/constitutive_plastic_none.f90 @@ -12,7 +12,10 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -module subroutine plastic_none_init +module subroutine plastic_none_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & Ninstance, & @@ -22,7 +25,7 @@ module subroutine plastic_none_init write(6,'(/,a)') ' <<<+- plastic_'//PLASTICITY_NONE_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_plasticity == PLASTICITY_NONE_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance do p = 1, size(phase_plasticity) diff --git a/src/constitutive_plastic_nonlocal.f90 b/src/constitutive_plastic_nonlocal.f90 index 7f21d0194..596c62c86 100644 --- a/src/constitutive_plastic_nonlocal.f90 +++ b/src/constitutive_plastic_nonlocal.f90 @@ -163,7 +163,10 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -module subroutine plastic_nonlocal_init +module subroutine plastic_nonlocal_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & Ninstance, & @@ -188,7 +191,7 @@ module subroutine plastic_nonlocal_init write(6,'(a)') ' http://publications.rwth-aachen.de/record/229993' Ninstance = count(phase_plasticity == PLASTICITY_NONLOCAL_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(param(Ninstance)) @@ -522,7 +525,7 @@ end subroutine plastic_nonlocal_init !-------------------------------------------------------------------------------------------------- !> @brief calculates quantities characterizing the microstructure !-------------------------------------------------------------------------------------------------- -module subroutine plastic_nonlocal_dependentState(F, Fp, instance, of, ip, el) +module subroutine plastic_nonlocal_dependentState(F, Fp, instance, of, ip, el,debug_constitutive) real(pReal), dimension(3,3), intent(in) :: & F, & @@ -532,6 +535,8 @@ module subroutine plastic_nonlocal_dependentState(F, Fp, instance, of, ip, el) of, & ip, & el + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & no, & !< neighbor offset @@ -541,7 +546,8 @@ module subroutine plastic_nonlocal_dependentState(F, Fp, instance, of, ip, el) c, & ! index of dilsocation character (edge, screw) s, & ! slip system index dir, & - n + n, & + debug_e, debug_i real(pReal) :: & FVsize, & nRealNeighbors ! number of really existing neighbors @@ -710,9 +716,11 @@ module subroutine plastic_nonlocal_dependentState(F, Fp, instance, of, ip, el) endif #ifdef DEBUG - if (iand(debug_level(debug_constitutive),debug_levelExtensive) /= 0 & + debug_e = debug_root%get_asInt('element',defaultVal=1) + debug_i = debug_root%get_asInt('integrationpoint',defaultVal=1) + if (debug_constitutive%contains('extensive') & .and. ((debug_e == el .and. debug_i == ip)& - .or. .not. iand(debug_level(debug_constitutive),debug_levelSelective) /= 0)) then + .or. .not. debug_constitutive%contains('selective'))) then write(6,'(/,a,i8,1x,i2,1x,i1,/)') '<< CONST >> nonlocal_microstructure at el ip ',el,ip write(6,'(a,/,12x,12(e10.3,1x))') '<< CONST >> rhoForest', stt%rho_forest(:,of) write(6,'(a,/,12x,12(f10.5,1x))') '<< CONST >> tauThreshold / MPa', dst%tau_pass(:,of)*1e-6 @@ -836,7 +844,7 @@ end subroutine plastic_nonlocal_LpAndItsTangent !-------------------------------------------------------------------------------------------------- !> @brief (instantaneous) incremental change of microstructure !-------------------------------------------------------------------------------------------------- -module subroutine plastic_nonlocal_deltaState(Mp,instance,of,ip,el) +module subroutine plastic_nonlocal_deltaState(Mp,instance,of,ip,el,debug_constitutive) real(pReal), dimension(3,3), intent(in) :: & Mp !< MandelStress @@ -845,13 +853,16 @@ module subroutine plastic_nonlocal_deltaState(Mp,instance,of,ip,el) of, & !< offset ip, & el + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & ph, & !< phase ns, & ! short notation for the total number of active slip systems c, & ! character of dislocation t, & ! type of dislocation - s ! index of my current slip system + s, & ! index of my current slip system + debug_e, debug_i real(pReal), dimension(param(instance)%sum_N_sl,10) :: & deltaRhoRemobilization, & ! density increment by remobilization deltaRhoDipole2SingleStress ! density increment by dipole dissociation (by stress change) @@ -927,9 +938,11 @@ module subroutine plastic_nonlocal_deltaState(Mp,instance,of,ip,el) del%rho(:,of) = reshape(deltaRhoRemobilization + deltaRhoDipole2SingleStress, [10*ns]) #ifdef DEBUG - if (iand(debug_level(debug_constitutive),debug_levelExtensive) /= 0 & + debug_e = debug_root%get_asInt('element',defaultVal=1) + debug_i = debug_root%get_asInt('integrationpoint',defaultVal=1) + if (debug_constitutive%contains('extensive') & .and. ((debug_e == el .and. debug_i == ip)& - .or. .not. iand(debug_level(debug_constitutive),debug_levelSelective) /= 0 )) then + .or. .not. debug_constitutive%contains('selective'))) then write(6,'(a,/,8(12x,12(e12.5,1x),/))') '<< CONST >> dislocation remobilization', deltaRhoRemobilization(:,1:8) write(6,'(a,/,10(12x,12(e12.5,1x),/),/)') '<< CONST >> dipole dissociation by stress increase', deltaRhoDipole2SingleStress endif @@ -944,7 +957,7 @@ end subroutine plastic_nonlocal_deltaState !> @brief calculates the rate of change of microstructure !--------------------------------------------------------------------------------------------------- module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, & - instance,of,ip,el) + instance,of,ip,el,debug_constitutive) real(pReal), dimension(3,3), intent(in) :: & Mp !< MandelStress @@ -959,13 +972,16 @@ module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, & of, & ip, & !< current integration point el !< current element number + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & ph, & ns, & !< short notation for the total number of active slip systems c, & !< character of dislocation t, & !< type of dislocation - s !< index of my current slip system + s, & !< index of my current slip system + debug_e, debug_i real(pReal), dimension(param(instance)%sum_N_sl,10) :: & rho, & rho0, & !< dislocation density at beginning of time step @@ -1016,9 +1032,11 @@ module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, & gdot = rhoSgl(:,1:4) * v * spread(prm%burgers,2,4) #ifdef DEBUG - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0 & - .and. ((debug_e == el .and. debug_i == ip)& - .or. .not. iand(debug_level(debug_constitutive),debug_levelSelective) /= 0 )) then + debug_e = debug_root%get_asInt('element',defaultVal=1) + debug_i = debug_root%get_asInt('integrationpoint',defaultVal=1) + if (debug_constitutive%contains('basic') & + .and. ((debug_e == el .and. debug_i == ip) & + .or. .not. debug_constitutive%contains('selective') )) then write(6,'(a,/,10(12x,12(e12.5,1x),/))') '<< CONST >> rho / 1/m^2', rhoSgl, rhoDip write(6,'(a,/,4(12x,12(e12.5,1x),/))') '<< CONST >> gdot / 1/s',gdot endif @@ -1117,7 +1135,7 @@ module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, & - rhoDip(s,1) / timestep - rhoDotAthermalAnnihilation(s,9) & - rhoDotSingle2DipoleGlide(s,9)) ! make sure that we do not annihilate more dipoles than we have - rhoDot = rhoDotFlux(F,Fp,timestep, instance,of,ip,el) & + rhoDot = rhoDotFlux(F,Fp,timestep, instance,of,ip,el,debug_constitutive) & + rhoDotMultiplication & + rhoDotSingle2DipoleGlide & + rhoDotAthermalAnnihilation & @@ -1127,7 +1145,7 @@ module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, & if ( any(rho(:,mob) + rhoDot(:,1:4) * timestep < -prm%atol_rho) & .or. any(rho(:,dip) + rhoDot(:,9:10) * timestep < -prm%atol_rho)) then #ifdef DEBUG - if (iand(debug_level(debug_constitutive),debug_levelExtensive) /= 0) then + if (debug_constitutive%contains('extensive')) then write(6,'(a,i5,a,i2)') '<< CONST >> evolution rate leads to negative density at el ',el,' ip ',ip write(6,'(a)') '<< CONST >> enforcing cutback !!!' endif @@ -1146,7 +1164,7 @@ end subroutine plastic_nonlocal_dotState !--------------------------------------------------------------------------------------------------- !> @brief calculates the rate of change of microstructure !--------------------------------------------------------------------------------------------------- -function rhoDotFlux(F,Fp,timestep, instance,of,ip,el) +function rhoDotFlux(F,Fp,timestep, instance,of,ip,el,debug_constitutive) real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: & F, & !< elastic deformation gradient @@ -1158,7 +1176,9 @@ function rhoDotFlux(F,Fp,timestep, instance,of,ip,el) of, & ip, & !< current integration point el !< current element number - + class(tNode), pointer, intent(in) :: & + debug_constitutive + integer :: & ph, & neighbor_instance, & !< instance of my neighbor's plasticity @@ -1239,7 +1259,7 @@ function rhoDotFlux(F,Fp,timestep, instance,of,ip,el) .and. prm%CFLfactor * abs(v0) * timestep & > IPvolume(ip,el) / maxval(IParea(:,ip,el)))) then ! ...with velocity above critical value (we use the reference volume and area for simplicity here) #ifdef DEBUG - if (iand(debug_level(debug_constitutive),debug_levelExtensive) /= 0) then + if (debug_constitutive%contains('extensive')) then write(6,'(a,i5,a,i2)') '<< CONST >> CFL condition not fullfilled at el ',el,' ip ',ip write(6,'(a,e10.3,a,e10.3)') '<< CONST >> velocity is at ', & maxval(abs(v0), abs(gdot) > 0.0_pReal & diff --git a/src/constitutive_plastic_phenopowerlaw.f90 b/src/constitutive_plastic_phenopowerlaw.f90 index fa273cbd3..76e2606f8 100644 --- a/src/constitutive_plastic_phenopowerlaw.f90 +++ b/src/constitutive_plastic_phenopowerlaw.f90 @@ -66,7 +66,10 @@ contains !> @brief Perform module initialization. !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -module subroutine plastic_phenopowerlaw_init +module subroutine plastic_phenopowerlaw_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: & Ninstance, & @@ -86,7 +89,7 @@ module subroutine plastic_phenopowerlaw_init write(6,'(/,a)') ' <<<+- plastic_'//PLASTICITY_PHENOPOWERLAW_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_plasticity == PLASTICITY_PHENOPOWERLAW_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(param(Ninstance)) diff --git a/src/kinematics_cleavage_opening.f90 b/src/kinematics_cleavage_opening.f90 index 3366a7b1e..2ffb27f2e 100644 --- a/src/kinematics_cleavage_opening.f90 +++ b/src/kinematics_cleavage_opening.f90 @@ -43,7 +43,10 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine kinematics_cleavage_opening_init +subroutine kinematics_cleavage_opening_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: Ninstance,p integer, dimension(:), allocatable :: N_cl !< active number of cleavage systems per family @@ -52,7 +55,7 @@ subroutine kinematics_cleavage_opening_init write(6,'(/,a)') ' <<<+- kinematics_'//KINEMATICS_CLEAVAGE_OPENING_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_kinematics == KINEMATICS_CLEAVAGE_OPENING_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(kinematics_cleavage_opening_instance(size(config_phase)), source=0) diff --git a/src/kinematics_slipplane_opening.f90 b/src/kinematics_slipplane_opening.f90 index 833fa8759..278754b79 100644 --- a/src/kinematics_slipplane_opening.f90 +++ b/src/kinematics_slipplane_opening.f90 @@ -45,7 +45,10 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine kinematics_slipplane_opening_init +subroutine kinematics_slipplane_opening_init(debug_constitutive) + + class(tNode), pointer , intent(in) :: & + debug_constitutive integer :: Ninstance,p,i character(len=pStringLen) :: extmsg = '' @@ -55,7 +58,7 @@ subroutine kinematics_slipplane_opening_init write(6,'(/,a)') ' <<<+- kinematics_'//KINEMATICS_SLIPPLANE_OPENING_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_kinematics == KINEMATICS_SLIPPLANE_OPENING_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(kinematics_slipplane_opening_instance(size(config_phase)), source=0) diff --git a/src/kinematics_thermal_expansion.f90 b/src/kinematics_thermal_expansion.f90 index acf3a5067..1a84fb7b8 100644 --- a/src/kinematics_thermal_expansion.f90 +++ b/src/kinematics_thermal_expansion.f90 @@ -38,7 +38,10 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine kinematics_thermal_expansion_init +subroutine kinematics_thermal_expansion_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: Ninstance,p,i real(pReal), dimension(:), allocatable :: temp @@ -46,7 +49,7 @@ subroutine kinematics_thermal_expansion_init write(6,'(/,a)') ' <<<+- kinematics_'//KINEMATICS_thermal_expansion_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_kinematics == KINEMATICS_thermal_expansion_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(kinematics_thermal_expansion_instance(size(config_phase)), source=0) diff --git a/src/source_damage_anisoBrittle.f90 b/src/source_damage_anisoBrittle.f90 index b3af24f38..5826e7160 100644 --- a/src/source_damage_anisoBrittle.f90 +++ b/src/source_damage_anisoBrittle.f90 @@ -53,7 +53,10 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine source_damage_anisoBrittle_init +subroutine source_damage_anisoBrittle_init(debug_constitutive) + + class(tNode), pointer , intent(in) :: & + debug_constitutive integer :: Ninstance,sourceOffset,NipcMyPhase,p integer, dimension(:), allocatable :: N_cl @@ -62,7 +65,7 @@ subroutine source_damage_anisoBrittle_init write(6,'(/,a)') ' <<<+- source_'//SOURCE_DAMAGE_ANISOBRITTLE_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_source == SOURCE_DAMAGE_ANISOBRITTLE_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(source_damage_anisoBrittle_offset (size(config_phase)), source=0) diff --git a/src/source_damage_anisoDuctile.f90 b/src/source_damage_anisoDuctile.f90 index 79cc0c2f7..73e68d021 100644 --- a/src/source_damage_anisoDuctile.f90 +++ b/src/source_damage_anisoDuctile.f90 @@ -46,7 +46,10 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine source_damage_anisoDuctile_init +subroutine source_damage_anisoDuctile_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: Ninstance,sourceOffset,NipcMyPhase,p integer, dimension(:), allocatable :: N_sl @@ -55,7 +58,7 @@ subroutine source_damage_anisoDuctile_init write(6,'(/,a)') ' <<<+- source_'//SOURCE_DAMAGE_ANISODUCTILE_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_source == SOURCE_DAMAGE_ANISODUCTILE_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(source_damage_anisoDuctile_offset (size(config_phase)), source=0) diff --git a/src/source_damage_isoBrittle.f90 b/src/source_damage_isoBrittle.f90 index 9eacb4516..e1c0625e4 100644 --- a/src/source_damage_isoBrittle.f90 +++ b/src/source_damage_isoBrittle.f90 @@ -45,7 +45,10 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine source_damage_isoBrittle_init +subroutine source_damage_isoBrittle_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: Ninstance,sourceOffset,NipcMyPhase,p character(len=pStringLen) :: extmsg = '' @@ -53,7 +56,7 @@ subroutine source_damage_isoBrittle_init write(6,'(/,a)') ' <<<+- source_'//SOURCE_DAMAGE_ISOBRITTLE_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_source == SOURCE_DAMAGE_ISOBRITTLE_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(source_damage_isoBrittle_offset (size(config_phase)), source=0) diff --git a/src/source_damage_isoDuctile.f90 b/src/source_damage_isoDuctile.f90 index 96754725d..8ba28ee1b 100644 --- a/src/source_damage_isoDuctile.f90 +++ b/src/source_damage_isoDuctile.f90 @@ -44,7 +44,10 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine source_damage_isoDuctile_init +subroutine source_damage_isoDuctile_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: Ninstance,sourceOffset,NipcMyPhase,p character(len=pStringLen) :: extmsg = '' @@ -52,7 +55,7 @@ subroutine source_damage_isoDuctile_init write(6,'(/,a)') ' <<<+- source_'//SOURCE_DAMAGE_ISODUCTILE_LABEL//' init -+>>>'; flush(6) Ninstance = count(phase_source == SOURCE_DAMAGE_ISODUCTILE_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(source_damage_isoDuctile_offset (size(config_phase)), source=0) diff --git a/src/source_thermal_dissipation.f90 b/src/source_thermal_dissipation.f90 index c323e68b5..3b08f7e25 100644 --- a/src/source_thermal_dissipation.f90 +++ b/src/source_thermal_dissipation.f90 @@ -37,14 +37,17 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine source_thermal_dissipation_init +subroutine source_thermal_dissipation_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: Ninstance,sourceOffset,NipcMyPhase,p write(6,'(/,a)') ' <<<+- source_'//SOURCE_thermal_dissipation_label//' init -+>>>'; flush(6) Ninstance = count(phase_source == SOURCE_THERMAL_DISSIPATION_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(source_thermal_dissipation_offset (size(config_phase)), source=0) diff --git a/src/source_thermal_externalheat.f90 b/src/source_thermal_externalheat.f90 index 06b8a5197..482de87de 100644 --- a/src/source_thermal_externalheat.f90 +++ b/src/source_thermal_externalheat.f90 @@ -41,14 +41,17 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine source_thermal_externalheat_init +subroutine source_thermal_externalheat_init(debug_constitutive) + + class(tNode), pointer, intent(in) :: & + debug_constitutive integer :: Ninstance,sourceOffset,NipcMyPhase,p write(6,'(/,a)') ' <<<+- source_'//SOURCE_thermal_externalheat_label//' init -+>>>'; flush(6) Ninstance = count(phase_source == SOURCE_thermal_externalheat_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (debug_constitutive%contains('basic')) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance allocate(source_thermal_externalheat_offset (size(config_phase)), source=0)