From 731057f7e037dcf64d6ff33d8aea55c3d6a1e3fd Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 4 Dec 2023 16:25:31 +0100 Subject: [PATCH 1/5] consistent with other plasticity laws --- ...phase_mechanical_plastic_kinehardening.f90 | 20 +++++----- ...phase_mechanical_plastic_phenopowerlaw.f90 | 39 +++++++++---------- 2 files changed, 27 insertions(+), 32 deletions(-) diff --git a/src/phase_mechanical_plastic_kinehardening.f90 b/src/phase_mechanical_plastic_kinehardening.f90 index 36cc8a37a..47dce0ce5 100644 --- a/src/phase_mechanical_plastic_kinehardening.f90 +++ b/src/phase_mechanical_plastic_kinehardening.f90 @@ -74,7 +74,6 @@ module function plastic_kinehardening_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & ph, o, & - N_fam, & Nmembers, & sizeState, sizeDeltaState, sizeDotState, & startIndex, endIndex @@ -138,7 +137,6 @@ module function plastic_kinehardening_init() result(myPlasticity) N_sl = pl%get_as1dInt('N_sl',defaultVal=emptyIntArray) prm%sum_N_sl = sum(abs(N_sl)) slipActive: if (prm%sum_N_sl > 0) then - N_fam = size(N_sl) prm%systems_sl = crystal_labels_slip(N_sl,phase_lattice(ph)) prm%P = crystal_SchmidMatrix_slip(N_sl,phase_lattice(ph),phase_cOverA(ph)) @@ -155,15 +153,15 @@ module function plastic_kinehardening_init() result(myPlasticity) prm%h_sl_sl = crystal_interaction_SlipBySlip(N_sl,pl%get_as1dReal('h_sl-sl'), & phase_lattice(ph)) - xi_0 = math_expand(pl%get_as1dReal('xi_0', requiredSize=N_fam),N_sl) - prm%dot_gamma_0 = math_expand(pl%get_as1dReal('dot_gamma_0', requiredSize=N_fam),N_sl) - prm%n = math_expand(pl%get_as1dReal('n', requiredSize=N_fam),N_sl) - prm%xi_inf = math_expand(pl%get_as1dReal('xi_inf', requiredSize=N_fam),N_sl) - prm%chi_inf = math_expand(pl%get_as1dReal('chi_inf', requiredSize=N_fam),N_sl) - prm%h_0_xi = math_expand(pl%get_as1dReal('h_0_xi', requiredSize=N_fam),N_sl) - prm%h_0_chi = math_expand(pl%get_as1dReal('h_0_chi', requiredSize=N_fam),N_sl) - prm%h_inf_xi = math_expand(pl%get_as1dReal('h_inf_xi', requiredSize=N_fam),N_sl) - prm%h_inf_chi = math_expand(pl%get_as1dReal('h_inf_chi', requiredSize=N_fam),N_sl) + xi_0 = math_expand(pl%get_as1dReal('xi_0', requiredSize=size(N_sl)),N_sl) + prm%dot_gamma_0 = math_expand(pl%get_as1dReal('dot_gamma_0', requiredSize=size(N_sl)),N_sl) + prm%n = math_expand(pl%get_as1dReal('n', requiredSize=size(N_sl)),N_sl) + prm%xi_inf = math_expand(pl%get_as1dReal('xi_inf', requiredSize=size(N_sl)),N_sl) + prm%chi_inf = math_expand(pl%get_as1dReal('chi_inf', requiredSize=size(N_sl)),N_sl) + prm%h_0_xi = math_expand(pl%get_as1dReal('h_0_xi', requiredSize=size(N_sl)),N_sl) + prm%h_0_chi = math_expand(pl%get_as1dReal('h_0_chi', requiredSize=size(N_sl)),N_sl) + prm%h_inf_xi = math_expand(pl%get_as1dReal('h_inf_xi', requiredSize=size(N_sl)),N_sl) + prm%h_inf_chi = math_expand(pl%get_as1dReal('h_inf_chi', requiredSize=size(N_sl)),N_sl) !-------------------------------------------------------------------------------------------------- ! sanity checks diff --git a/src/phase_mechanical_plastic_phenopowerlaw.f90 b/src/phase_mechanical_plastic_phenopowerlaw.f90 index 94d00d342..b656b4a78 100644 --- a/src/phase_mechanical_plastic_phenopowerlaw.f90 +++ b/src/phase_mechanical_plastic_phenopowerlaw.f90 @@ -77,7 +77,6 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & ph, i, & - N_fam, & Nmembers, & sizeState, sizeDotState, & startIndex, endIndex @@ -142,17 +141,16 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! slip related parameters slipActive: if (prm%sum_N_sl > 0) then - N_fam = size(N_sl) - ones = [(1.0_pREAL,i=1,N_fam)] - prm%dot_gamma_0_sl = math_expand(pl%get_as1dReal('dot_gamma_0_sl',requiredSize=N_fam,defaultVal= ones), N_sl) - prm%n_sl = math_expand(pl%get_as1dReal('n_sl', requiredSize=N_fam,defaultVal= ones), N_sl) - prm%a_sl = math_expand(pl%get_as1dReal('a_sl', requiredSize=N_fam,defaultVal= ones), N_sl) - prm%h_0_sl_sl = math_expand(pl%get_as1dReal('h_0_sl-sl', requiredSize=N_fam,defaultVal=0.0_pREAL*ones), N_sl) - prm%c_1 = math_expand(pl%get_as1dReal('c_1', requiredSize=N_fam,defaultVal=0.0_pREAL*ones), N_sl) - prm%c_2 = math_expand(pl%get_as1dReal('c_2', requiredSize=N_fam,defaultVal= ones), N_sl) - xi_0_sl = math_expand(pl%get_as1dReal('xi_0_sl', requiredSize=N_fam), N_sl) - prm%xi_inf_sl = math_expand(pl%get_as1dReal('xi_inf_sl', requiredSize=N_fam), N_sl) - prm%f_sat_sl_tw = math_expand(pl%get_as1dReal('f_sat_sl-tw',requiredSize=N_fam,defaultVal=0.0_pREAL*ones), N_sl) + ones = [(1.0_pREAL,i=1,size(N_sl))] + prm%dot_gamma_0_sl = math_expand(pl%get_as1dReal('dot_gamma_0_sl',requiredSize=size(N_sl),defaultVal= ones), N_sl) + prm%n_sl = math_expand(pl%get_as1dReal('n_sl', requiredSize=size(N_sl),defaultVal= ones), N_sl) + prm%a_sl = math_expand(pl%get_as1dReal('a_sl', requiredSize=size(N_sl),defaultVal= ones), N_sl) + prm%h_0_sl_sl = math_expand(pl%get_as1dReal('h_0_sl-sl', requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) + prm%c_1 = math_expand(pl%get_as1dReal('c_1', requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) + prm%c_2 = math_expand(pl%get_as1dReal('c_2', requiredSize=size(N_sl),defaultVal= ones), N_sl) + xi_0_sl = math_expand(pl%get_as1dReal('xi_0_sl', requiredSize=size(N_sl)), N_sl) + prm%xi_inf_sl = math_expand(pl%get_as1dReal('xi_inf_sl', requiredSize=size(N_sl)), N_sl) + prm%f_sat_sl_tw = math_expand(pl%get_as1dReal('f_sat_sl-tw',requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) prm%h_sl_sl = crystal_interaction_SlipBySlip(N_sl,pl%get_as1dReal('h_sl-sl'),phase_lattice(ph)) @@ -194,14 +192,13 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! twin related parameters twinActive: if (prm%sum_N_tw > 0) then - N_fam = size(N_tw) - ones = [(1.0_pREAL,i=1,N_fam)] - prm%dot_gamma_0_tw = math_expand(pl%get_as1dReal('dot_gamma_0_tw', requiredSize=N_fam,defaultVal=ones), N_tw) - prm%n_tw = math_expand(pl%get_as1dReal('n_tw', requiredSize=N_fam,defaultVal= ones), N_tw) - prm%c_3 = math_expand(pl%get_as1dReal('c_3', requiredSize=N_fam,defaultVal=0.0_pREAL*ones), N_tw) - prm%c_4 = math_expand(pl%get_as1dReal('c_4', requiredSize=N_fam,defaultVal=0.0_pREAL*ones), N_tw) - prm%h_0_tw_tw = math_expand(pl%get_as1dReal('h_0_tw-tw', requiredSize=N_fam,defaultVal=0.0_pReal*ones), N_tw) - xi_0_tw = math_expand(pl%get_as1dReal('xi_0_tw', requiredSize=N_fam), N_tw) + ones = [(1.0_pREAL,i=1,size(N_tw))] + prm%dot_gamma_0_tw = math_expand(pl%get_as1dReal('dot_gamma_0_tw', requiredSize=size(N_tw),defaultVal=ones), N_tw) + prm%n_tw = math_expand(pl%get_as1dReal('n_tw', requiredSize=size(N_tw),defaultVal= ones), N_tw) + prm%c_3 = math_expand(pl%get_as1dReal('c_3', requiredSize=size(N_tw),defaultVal=0.0_pREAL*ones), N_tw) + prm%c_4 = math_expand(pl%get_as1dReal('c_4', requiredSize=size(N_tw),defaultVal=0.0_pREAL*ones), N_tw) + prm%h_0_tw_tw = math_expand(pl%get_as1dReal('h_0_tw-tw', requiredSize=size(N_tw),defaultVal=0.0_pReal*ones), N_tw) + xi_0_tw = math_expand(pl%get_as1dReal('xi_0_tw', requiredSize=size(N_tw)), N_tw) prm%gamma_char = crystal_characteristicShear_twin(N_tw,phase_lattice(ph),phase_cOverA(ph)) prm%h_tw_tw = crystal_interaction_TwinByTwin(N_tw,pl%get_as1dReal('h_tw-tw'),phase_lattice(ph)) @@ -230,7 +227,7 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! slip-twin related parameters slipAndTwinActive: if (prm%sum_N_sl > 0 .and. prm%sum_N_tw > 0) then - prm%h_0_tw_sl = math_expand(pl%get_as1dReal('h_0_tw-sl',requiredSize=N_fam,defaultVal=0.0_pReal*ones), N_tw) + prm%h_0_tw_sl = math_expand(pl%get_as1dReal('h_0_tw-sl',requiredSize=size(N_tw),defaultVal=0.0_pReal*ones), N_tw) prm%h_sl_tw = crystal_interaction_SlipByTwin(N_sl,N_tw,pl%get_as1dReal('h_sl-tw'),phase_lattice(ph)) prm%h_tw_sl = crystal_interaction_TwinBySlip(N_tw,N_sl,pl%get_as1dReal('h_tw-sl'),phase_lattice(ph)) else slipAndTwinActive From 1621473d99fe3a34a61151f1435e1ee96fb87c13 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 4 Dec 2023 16:35:25 +0100 Subject: [PATCH 2/5] should not have default values --- ...phase_mechanical_plastic_phenopowerlaw.f90 | 33 +++++++++---------- 1 file changed, 16 insertions(+), 17 deletions(-) diff --git a/src/phase_mechanical_plastic_phenopowerlaw.f90 b/src/phase_mechanical_plastic_phenopowerlaw.f90 index b656b4a78..77f74805c 100644 --- a/src/phase_mechanical_plastic_phenopowerlaw.f90 +++ b/src/phase_mechanical_plastic_phenopowerlaw.f90 @@ -141,16 +141,15 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! slip related parameters slipActive: if (prm%sum_N_sl > 0) then - ones = [(1.0_pREAL,i=1,size(N_sl))] - prm%dot_gamma_0_sl = math_expand(pl%get_as1dReal('dot_gamma_0_sl',requiredSize=size(N_sl),defaultVal= ones), N_sl) - prm%n_sl = math_expand(pl%get_as1dReal('n_sl', requiredSize=size(N_sl),defaultVal= ones), N_sl) - prm%a_sl = math_expand(pl%get_as1dReal('a_sl', requiredSize=size(N_sl),defaultVal= ones), N_sl) - prm%h_0_sl_sl = math_expand(pl%get_as1dReal('h_0_sl-sl', requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) - prm%c_1 = math_expand(pl%get_as1dReal('c_1', requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) - prm%c_2 = math_expand(pl%get_as1dReal('c_2', requiredSize=size(N_sl),defaultVal= ones), N_sl) - xi_0_sl = math_expand(pl%get_as1dReal('xi_0_sl', requiredSize=size(N_sl)), N_sl) - prm%xi_inf_sl = math_expand(pl%get_as1dReal('xi_inf_sl', requiredSize=size(N_sl)), N_sl) - prm%f_sat_sl_tw = math_expand(pl%get_as1dReal('f_sat_sl-tw',requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) + prm%dot_gamma_0_sl = math_expand(pl%get_as1dReal('dot_gamma_0_sl',requiredSize=size(N_sl) ), N_sl) + prm%n_sl = math_expand(pl%get_as1dReal('n_sl', requiredSize=size(N_sl) ), N_sl) + prm%a_sl = math_expand(pl%get_as1dReal('a_sl', requiredSize=size(N_sl) ), N_sl) + prm%h_0_sl_sl = math_expand(pl%get_as1dReal('h_0_sl-sl', requiredSize=size(N_sl) ), N_sl) + prm%c_1 = math_expand(pl%get_as1dReal('c_1', requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) + prm%c_2 = math_expand(pl%get_as1dReal('c_2', requiredSize=size(N_sl),defaultVal= ones), N_sl) + xi_0_sl = math_expand(pl%get_as1dReal('xi_0_sl', requiredSize=size(N_sl)), N_sl) + prm%xi_inf_sl = math_expand(pl%get_as1dReal('xi_inf_sl', requiredSize=size(N_sl)), N_sl) + prm%f_sat_sl_tw = math_expand(pl%get_as1dReal('f_sat_sl-tw', requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) prm%h_sl_sl = crystal_interaction_SlipBySlip(N_sl,pl%get_as1dReal('h_sl-sl'),phase_lattice(ph)) @@ -193,12 +192,12 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) ! twin related parameters twinActive: if (prm%sum_N_tw > 0) then ones = [(1.0_pREAL,i=1,size(N_tw))] - prm%dot_gamma_0_tw = math_expand(pl%get_as1dReal('dot_gamma_0_tw', requiredSize=size(N_tw),defaultVal=ones), N_tw) - prm%n_tw = math_expand(pl%get_as1dReal('n_tw', requiredSize=size(N_tw),defaultVal= ones), N_tw) - prm%c_3 = math_expand(pl%get_as1dReal('c_3', requiredSize=size(N_tw),defaultVal=0.0_pREAL*ones), N_tw) - prm%c_4 = math_expand(pl%get_as1dReal('c_4', requiredSize=size(N_tw),defaultVal=0.0_pREAL*ones), N_tw) - prm%h_0_tw_tw = math_expand(pl%get_as1dReal('h_0_tw-tw', requiredSize=size(N_tw),defaultVal=0.0_pReal*ones), N_tw) - xi_0_tw = math_expand(pl%get_as1dReal('xi_0_tw', requiredSize=size(N_tw)), N_tw) + prm%dot_gamma_0_tw = math_expand(pl%get_as1dReal('dot_gamma_0_tw', requiredSize=size(N_tw) ), N_tw) + prm%n_tw = math_expand(pl%get_as1dReal('n_tw', requiredSize=size(N_tw) ), N_tw) + prm%c_3 = math_expand(pl%get_as1dReal('c_3', requiredSize=size(N_tw),defaultVal=0.0_pREAL*ones), N_tw) + prm%c_4 = math_expand(pl%get_as1dReal('c_4', requiredSize=size(N_tw),defaultVal=0.0_pREAL*ones), N_tw) + prm%h_0_tw_tw = math_expand(pl%get_as1dReal('h_0_tw-tw', requiredSize=size(N_tw) ), N_tw) + xi_0_tw = math_expand(pl%get_as1dReal('xi_0_tw', requiredSize=size(N_tw)), N_tw) prm%gamma_char = crystal_characteristicShear_twin(N_tw,phase_lattice(ph),phase_cOverA(ph)) prm%h_tw_tw = crystal_interaction_TwinByTwin(N_tw,pl%get_as1dReal('h_tw-tw'),phase_lattice(ph)) @@ -227,7 +226,7 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! slip-twin related parameters slipAndTwinActive: if (prm%sum_N_sl > 0 .and. prm%sum_N_tw > 0) then - prm%h_0_tw_sl = math_expand(pl%get_as1dReal('h_0_tw-sl',requiredSize=size(N_tw),defaultVal=0.0_pReal*ones), N_tw) + prm%h_0_tw_sl = math_expand(pl%get_as1dReal('h_0_tw-sl',requiredSize=size(N_tw)), N_tw) prm%h_sl_tw = crystal_interaction_SlipByTwin(N_sl,N_tw,pl%get_as1dReal('h_sl-tw'),phase_lattice(ph)) prm%h_tw_sl = crystal_interaction_TwinBySlip(N_tw,N_sl,pl%get_as1dReal('h_tw-sl'),phase_lattice(ph)) else slipAndTwinActive From a5b5972e9ae0bad312dd76b37d91d4b2df764635 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 4 Dec 2023 16:48:09 +0100 Subject: [PATCH 3/5] already documented --- ...phase_mechanical_plastic_phenopowerlaw.f90 | 20 +++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/src/phase_mechanical_plastic_phenopowerlaw.f90 b/src/phase_mechanical_plastic_phenopowerlaw.f90 index 77f74805c..b074e2a67 100644 --- a/src/phase_mechanical_plastic_phenopowerlaw.f90 +++ b/src/phase_mechanical_plastic_phenopowerlaw.f90 @@ -176,14 +176,14 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) else slipActive xi_0_sl = emptyRealArray - allocate(prm%dot_gamma_0_sl, & !< reference shear strain rate for slip + allocate(prm%dot_gamma_0_sl, & prm%a_sl, & - prm%n_sl, & !< stress exponent for slip - prm%xi_inf_sl, & !< maximum critical shear stress for slip - prm%f_sat_sl_tw, & !< push-up factor for slip saturation due to twinning + prm%n_sl, & + prm%xi_inf_sl, & + prm%f_sat_sl_tw, & prm%c_1, & prm%c_2, & - prm%h_0_sl_sl, & !< reference hardening slip - slip + prm%h_0_sl_sl, & source=emptyRealArray) allocate(prm%h_sl_sl(0,0)) end if slipActive @@ -212,13 +212,13 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) else twinActive xi_0_tw = emptyRealArray - allocate(prm%dot_gamma_0_tw, & !< reference shear strain rate for twin - prm%n_tw, & !< stress exponent for twin + allocate(prm%dot_gamma_0_tw, & + prm%n_tw, & prm%c_3, & prm%c_4, & - prm%gamma_char, & !< characteristic shear for twins - prm%h_0_tw_sl, & !< reference hardening twin - slip - prm%h_0_tw_tw, & !< reference hardening twin - twin + prm%gamma_char, & + prm%h_0_tw_sl, & + prm%h_0_tw_tw, & source=emptyRealArray) allocate(prm%h_tw_tw(0,0)) end if twinActive From 660b4727470b1090b4064d63685719b6e85efdc9 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 4 Dec 2023 23:38:58 +0100 Subject: [PATCH 4/5] centralized function for ones/zeros --- src/misc.f90 | 40 ++++++++++++++++++- ...phase_mechanical_plastic_phenopowerlaw.f90 | 39 +++++++++--------- 2 files changed, 59 insertions(+), 20 deletions(-) diff --git a/src/misc.f90 b/src/misc.f90 index 902a5b731..1d546dd30 100644 --- a/src/misc.f90 +++ b/src/misc.f90 @@ -21,7 +21,9 @@ module misc misc_init, & misc_selfTest, & misc_optional, & - misc_prefixOptions + misc_prefixOptions, & + misc_ones, & + misc_zeros contains @@ -136,6 +138,34 @@ pure function misc_prefixOptions(string,prefix) result(prefixed) end function misc_prefixOptions +!-------------------------------------------------------------------------------------------------- +!> @brief 1D array of zeros. +!-------------------------------------------------------------------------------------------------- +pure function misc_zeros(N) + + integer, intent(in) :: N !< number of ones + real(pREAL), dimension(N) :: misc_zeros + + + misc_zeros = 0._pREAL + +end function misc_zeros + + +!-------------------------------------------------------------------------------------------------- +!> @brief 1D array of ones. +!-------------------------------------------------------------------------------------------------- +pure function misc_ones(N) + + integer, intent(in) :: N !< number of ones + real(pREAL), dimension(N) :: misc_ones + + + misc_ones = 1._pREAL + +end function misc_ones + + !-------------------------------------------------------------------------------------------------- !> @brief Check correctness of some misc functions. !-------------------------------------------------------------------------------------------------- @@ -143,7 +173,7 @@ subroutine misc_selfTest() real(pREAL) :: r character(len=:), allocatable :: str,out - + integer :: N call random_number(r) if (test_str('DAMASK') /= 'DAMASK') error stop 'optional_str, present' @@ -163,6 +193,12 @@ subroutine misc_selfTest() out=misc_prefixOptions(str,'p_') if (out /= '-p_a -1 -p_more 123 -p_flag -') error stop 'misc_prefixOptions' + N = int(r*99._pReal) + if (size(misc_zeros(N)) /= N) error stop 'shape zeros' + if (size(misc_ones(N)) /= N) error stop 'shape ones' + if (any(dNeq(misc_zeros(N),0.0_pReal))) error stop 'value zeros' + if (any(dNeq(misc_ones(N),1.0_pReal))) error stop 'value ones' + contains function test_str(str_in) result(str_out) diff --git a/src/phase_mechanical_plastic_phenopowerlaw.f90 b/src/phase_mechanical_plastic_phenopowerlaw.f90 index b074e2a67..2368e7e07 100644 --- a/src/phase_mechanical_plastic_phenopowerlaw.f90 +++ b/src/phase_mechanical_plastic_phenopowerlaw.f90 @@ -85,8 +85,7 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) N_tw !< number of twin-systems for a given twin family real(pREAL), dimension(:), allocatable :: & xi_0_sl, & !< initial critical shear stress for slip - xi_0_tw, & !< initial critical shear stress for twin - ones + xi_0_tw !< initial critical shear stress for twin real(pREAL), dimension(:,:), allocatable :: & a_nS !< non-Schmid coefficients character(len=:), allocatable :: & @@ -141,15 +140,18 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! slip related parameters slipActive: if (prm%sum_N_sl > 0) then - prm%dot_gamma_0_sl = math_expand(pl%get_as1dReal('dot_gamma_0_sl',requiredSize=size(N_sl) ), N_sl) - prm%n_sl = math_expand(pl%get_as1dReal('n_sl', requiredSize=size(N_sl) ), N_sl) - prm%a_sl = math_expand(pl%get_as1dReal('a_sl', requiredSize=size(N_sl) ), N_sl) - prm%h_0_sl_sl = math_expand(pl%get_as1dReal('h_0_sl-sl', requiredSize=size(N_sl) ), N_sl) - prm%c_1 = math_expand(pl%get_as1dReal('c_1', requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) - prm%c_2 = math_expand(pl%get_as1dReal('c_2', requiredSize=size(N_sl),defaultVal= ones), N_sl) - xi_0_sl = math_expand(pl%get_as1dReal('xi_0_sl', requiredSize=size(N_sl)), N_sl) - prm%xi_inf_sl = math_expand(pl%get_as1dReal('xi_inf_sl', requiredSize=size(N_sl)), N_sl) - prm%f_sat_sl_tw = math_expand(pl%get_as1dReal('f_sat_sl-tw', requiredSize=size(N_sl),defaultVal=0.0_pREAL*ones), N_sl) + prm%dot_gamma_0_sl = math_expand(pl%get_as1dReal('dot_gamma_0_sl',requiredSize=size(N_sl)), N_sl) + prm%n_sl = math_expand(pl%get_as1dReal('n_sl', requiredSize=size(N_sl)), N_sl) + prm%a_sl = math_expand(pl%get_as1dReal('a_sl', requiredSize=size(N_sl)), N_sl) + prm%h_0_sl_sl = math_expand(pl%get_as1dReal('h_0_sl-sl', requiredSize=size(N_sl)), N_sl) + xi_0_sl = math_expand(pl%get_as1dReal('xi_0_sl', requiredSize=size(N_sl)), N_sl) + prm%xi_inf_sl = math_expand(pl%get_as1dReal('xi_inf_sl', requiredSize=size(N_sl)), N_sl) + prm%c_1 = math_expand(pl%get_as1dReal('c_1', requiredSize=size(N_sl), & + defaultVal=misc_ones(size(N_sl))), N_sl) + prm%c_2 = math_expand(pl%get_as1dReal('c_2', requiredSize=size(N_sl), & + defaultVal=misc_ones(size(N_sl))), N_sl) + prm%f_sat_sl_tw = math_expand(pl%get_as1dReal('f_sat_sl-tw', requiredSize=size(N_sl), & + defaultVal=misc_zeros(size(N_sl))), N_sl) prm%h_sl_sl = crystal_interaction_SlipBySlip(N_sl,pl%get_as1dReal('h_sl-sl'),phase_lattice(ph)) @@ -191,13 +193,14 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! twin related parameters twinActive: if (prm%sum_N_tw > 0) then - ones = [(1.0_pREAL,i=1,size(N_tw))] - prm%dot_gamma_0_tw = math_expand(pl%get_as1dReal('dot_gamma_0_tw', requiredSize=size(N_tw) ), N_tw) - prm%n_tw = math_expand(pl%get_as1dReal('n_tw', requiredSize=size(N_tw) ), N_tw) - prm%c_3 = math_expand(pl%get_as1dReal('c_3', requiredSize=size(N_tw),defaultVal=0.0_pREAL*ones), N_tw) - prm%c_4 = math_expand(pl%get_as1dReal('c_4', requiredSize=size(N_tw),defaultVal=0.0_pREAL*ones), N_tw) - prm%h_0_tw_tw = math_expand(pl%get_as1dReal('h_0_tw-tw', requiredSize=size(N_tw) ), N_tw) - xi_0_tw = math_expand(pl%get_as1dReal('xi_0_tw', requiredSize=size(N_tw)), N_tw) + prm%dot_gamma_0_tw = math_expand(pl%get_as1dReal('dot_gamma_0_tw', requiredSize=size(N_tw)), N_tw) + prm%n_tw = math_expand(pl%get_as1dReal('n_tw', requiredSize=size(N_tw)), N_tw) + prm%h_0_tw_tw = math_expand(pl%get_as1dReal('h_0_tw-tw', requiredSize=size(N_tw)), N_tw) + xi_0_tw = math_expand(pl%get_as1dReal('xi_0_tw', requiredSize=size(N_tw)), N_tw) + prm%c_3 = math_expand(pl%get_as1dReal('c_3', requiredSize=size(N_tw), & + defaultVal=misc_ones(size(N_tw))), N_tw) + prm%c_4 = math_expand(pl%get_as1dReal('c_4', requiredSize=size(N_tw), & + defaultVal=misc_zeros(size(N_tw))), N_tw) prm%gamma_char = crystal_characteristicShear_twin(N_tw,phase_lattice(ph),phase_cOverA(ph)) prm%h_tw_tw = crystal_interaction_TwinByTwin(N_tw,pl%get_as1dReal('h_tw-tw'),phase_lattice(ph)) From 8ace3411b4c673daf9fa57589e57e45d6c0ca1cd Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 5 Dec 2023 15:58:06 +0100 Subject: [PATCH 5/5] restore correct behavior --- src/misc.f90 | 2 +- src/phase_mechanical_plastic_phenopowerlaw.f90 | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/src/misc.f90 b/src/misc.f90 index 1d546dd30..ac8f37839 100644 --- a/src/misc.f90 +++ b/src/misc.f90 @@ -143,7 +143,7 @@ end function misc_prefixOptions !-------------------------------------------------------------------------------------------------- pure function misc_zeros(N) - integer, intent(in) :: N !< number of ones + integer, intent(in) :: N !< number of zeros real(pREAL), dimension(N) :: misc_zeros diff --git a/src/phase_mechanical_plastic_phenopowerlaw.f90 b/src/phase_mechanical_plastic_phenopowerlaw.f90 index 2368e7e07..095f72b74 100644 --- a/src/phase_mechanical_plastic_phenopowerlaw.f90 +++ b/src/phase_mechanical_plastic_phenopowerlaw.f90 @@ -147,7 +147,7 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) xi_0_sl = math_expand(pl%get_as1dReal('xi_0_sl', requiredSize=size(N_sl)), N_sl) prm%xi_inf_sl = math_expand(pl%get_as1dReal('xi_inf_sl', requiredSize=size(N_sl)), N_sl) prm%c_1 = math_expand(pl%get_as1dReal('c_1', requiredSize=size(N_sl), & - defaultVal=misc_ones(size(N_sl))), N_sl) + defaultVal=misc_zeros(size(N_sl))), N_sl) prm%c_2 = math_expand(pl%get_as1dReal('c_2', requiredSize=size(N_sl), & defaultVal=misc_ones(size(N_sl))), N_sl) prm%f_sat_sl_tw = math_expand(pl%get_as1dReal('f_sat_sl-tw', requiredSize=size(N_sl), &