separation by instance does not add any value

src/phase_mechanical.f90

@@ -60,7 +60,7 @@ submodule(phase) mechanical
     module subroutine plastic_init
     end subroutine plastic_init
 
-    module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,me)
+    module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,ph,me)
       real(pReal), dimension(3,3),     intent(out) :: &
         Li                                                                                          !< inleastic velocity gradient
       real(pReal), dimension(3,3,3,3), intent(out)  :: &
@@ -68,7 +68,7 @@ submodule(phase) mechanical
       real(pReal), dimension(3,3),     intent(in) :: &
         Mi                                                                                          !< Mandel stress
       integer,                         intent(in) :: &
-        instance, &
+        ph, &
         me
     end subroutine plastic_isotropic_LiAndItsTangent
 
@@ -130,8 +130,8 @@ submodule(phase) mechanical
     end subroutine plastic_LpAndItsTangents
 
 
-    module subroutine plastic_isotropic_results(instance,group)
+    module subroutine plastic_isotropic_results(ph,group)
-      integer,          intent(in) :: instance
+      integer,          intent(in) :: ph
       character(len=*), intent(in) :: group
     end subroutine plastic_isotropic_results
 
@@ -145,13 +145,13 @@ submodule(phase) mechanical
       character(len=*), intent(in) :: group
     end subroutine plastic_kinehardening_results
 
-    module subroutine plastic_dislotwin_results(instance,group)
+    module subroutine plastic_dislotwin_results(ph,group)
-      integer,          intent(in) :: instance
+      integer,          intent(in) :: ph
       character(len=*), intent(in) :: group
     end subroutine plastic_dislotwin_results
 
-    module subroutine plastic_dislotungsten_results(instance,group)
+    module subroutine plastic_dislotungsten_results(ph,group)
-      integer,          intent(in) :: instance
+      integer,          intent(in) :: ph
       character(len=*), intent(in) :: group
     end subroutine plastic_dislotungsten_results
 
@@ -403,7 +403,7 @@ module subroutine mechanical_results(group,ph)
   select case(phase_plasticity(ph))
 
     case(PLASTICITY_ISOTROPIC_ID)
-      call plastic_isotropic_results(phase_plasticInstance(ph),group//'plastic/')
+      call plastic_isotropic_results(ph,group//'plastic/')
 
     case(PLASTICITY_PHENOPOWERLAW_ID)
       call plastic_phenopowerlaw_results(ph,group//'plastic/')
@@ -412,10 +412,10 @@ module subroutine mechanical_results(group,ph)
       call plastic_kinehardening_results(ph,group//'plastic/')
 
     case(PLASTICITY_DISLOTWIN_ID)
-      call plastic_dislotwin_results(phase_plasticInstance(ph),group//'plastic/')
+      call plastic_dislotwin_results(ph,group//'plastic/')
 
     case(PLASTICITY_DISLOTUNGSTEN_ID)
-      call plastic_dislotungsten_results(phase_plasticInstance(ph),group//'plastic/')
+      call plastic_dislotungsten_results(ph,group//'plastic/')
 
     case(PLASTICITY_NONLOCAL_ID)
       call plastic_nonlocal_results(phase_plasticInstance(ph),group//'plastic/')

src/phase_mechanical_eigen.f90

@@ -180,7 +180,7 @@ module subroutine phase_LiAndItsTangents(Li, dLi_dS, dLi_dFi, &
 
   plasticType: select case (phase_plasticity(ph))
     case (PLASTICITY_isotropic_ID) plasticType
-      call plastic_isotropic_LiAndItsTangent(my_Li, my_dLi_dS, S ,phase_plasticInstance(ph),me)
+      call plastic_isotropic_LiAndItsTangent(my_Li, my_dLi_dS, S ,ph,me)
       Li = Li + my_Li
       dLi_dS = dLi_dS + my_dLi_dS
       active = .true.

src/phase_mechanical_plastic.f90

@@ -179,17 +179,17 @@ submodule(phase:mechanical) plastic
         el                                                                                          !< current element number
     end subroutine nonlocal_dotState
 
-    module subroutine dislotwin_dependentState(T,instance,me)
+    module subroutine dislotwin_dependentState(T,ph,me)
       integer,       intent(in) :: &
-        instance, &
+        ph, &
         me
       real(pReal),   intent(in) :: &
         T
     end subroutine dislotwin_dependentState
 
-    module subroutine dislotungsten_dependentState(instance,me)
+    module subroutine dislotungsten_dependentState(ph,me)
       integer,       intent(in) :: &
-        instance, &
+        ph, &
         me
     end subroutine dislotungsten_dependentState
 
@@ -374,10 +374,10 @@ module subroutine plastic_dependentState(co, ip, el)
   plasticType: select case (phase_plasticity(material_phaseAt(co,el)))
 
     case (PLASTICITY_DISLOTWIN_ID) plasticType
-      call dislotwin_dependentState(thermal_T(ph,me),instance,me)
+      call dislotwin_dependentState(thermal_T(ph,me),ph,me)
 
     case (PLASTICITY_DISLOTUNGSTEN_ID) plasticType
-      call dislotungsten_dependentState(instance,me)
+      call dislotungsten_dependentState(ph,me)
 
     case (PLASTICITY_NONLOCAL_ID) plasticType
       call nonlocal_dependentState(instance,me,ip,el)

src/phase_mechanical_plastic_dislotungsten.f90

@@ -78,7 +78,6 @@ module function plastic_dislotungsten_init() result(myPlasticity)
 
   logical, dimension(:), allocatable :: myPlasticity
   integer :: &
-    Ninstances, &
     p, i, &
     Nconstituents, &
     sizeState, sizeDotState, &
@@ -97,29 +96,29 @@ module function plastic_dislotungsten_init() result(myPlasticity)
     mech, &
     pl
 
+
   myPlasticity = plastic_active('dislotungsten')
-  Ninstances = count(myPlasticity)
+  if(count(myPlasticity) == 0) return
-  if(Ninstances == 0) return
 
   print'(/,a)', ' <<<+-  phase:mechanical:plastic:dislotungsten init  -+>>>'
-  print'(a,i0)', ' # phases: ',Ninstances; flush(IO_STDOUT)
+  print'(a,i0)', ' # phases: ',count(myPlasticity); flush(IO_STDOUT)
-
 
   print*, 'Cereceda et al., International Journal of Plasticity 78:242–256, 2016'
   print*, 'https://dx.doi.org/10.1016/j.ijplas.2015.09.002'
 
-  allocate(param(Ninstances))
-  allocate(state(Ninstances))
-  allocate(dotState(Ninstances))
-  allocate(dependentState(Ninstances))
 
   phases => config_material%get('phase')
-  i = 0
+  allocate(param(phases%length))
+  allocate(state(phases%length))
+  allocate(dotState(phases%length))
+  allocate(dependentState(phases%length))
+
+
   do p = 1, phases%length
+    if(.not. myPlasticity(p)) cycle
     phase => phases%get(p)
     mech  => phase%get('mechanics')
-    if(.not. myPlasticity(p)) cycle
+    i = p
-    i = i + 1
     associate(prm => param(i), &
               dot => dotState(i), &
               stt => state(i), &
@@ -290,16 +289,16 @@ pure module subroutine dislotungsten_LpAndItsTangent(Lp,dLp_dMp, &
 
   integer :: &
     i,k,l,m,n
-  real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl) :: &
     dot_gamma_pos,dot_gamma_neg, &
     ddot_gamma_dtau_pos,ddot_gamma_dtau_neg
 
   Lp = 0.0_pReal
   dLp_dMp = 0.0_pReal
 
-  associate(prm => param(phase_plasticInstance(ph)))
+  associate(prm => param(ph))
 
-  call kinetics(Mp,T,phase_plasticInstance(ph),me,dot_gamma_pos,dot_gamma_neg,ddot_gamma_dtau_pos,ddot_gamma_dtau_neg)
+  call kinetics(Mp,T,ph,me,dot_gamma_pos,dot_gamma_neg,ddot_gamma_dtau_pos,ddot_gamma_dtau_neg)
   do i = 1, prm%sum_N_sl
     Lp = Lp + (dot_gamma_pos(i)+dot_gamma_neg(i))*prm%P_sl(1:3,1:3,i)
     forall (k=1:3,l=1:3,m=1:3,n=1:3) &
@@ -328,7 +327,7 @@ module subroutine dislotungsten_dotState(Mp,T,ph,me)
 
   real(pReal) :: &
     VacancyDiffusion
-  real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl) :: &
     gdot_pos, gdot_neg,&
     tau_pos,&
     tau_neg, &
@@ -337,10 +336,10 @@ module subroutine dislotungsten_dotState(Mp,T,ph,me)
     dot_rho_dip_climb, &
     dip_distance
 
-  associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph)),&
+  associate(prm => param(ph), stt => state(ph),&
-            dot => dotState(phase_plasticInstance(ph)), dst => dependentState(phase_plasticInstance(ph)))
+            dot => dotState(ph), dst => dependentState(ph))
 
-  call kinetics(Mp,T,phase_plasticInstance(ph),me,&
+  call kinetics(Mp,T,ph,me,&
                 gdot_pos,gdot_neg, &
                 tau_pos_out = tau_pos,tau_neg_out = tau_neg)
 
@@ -377,16 +376,16 @@ end subroutine dislotungsten_dotState
 !--------------------------------------------------------------------------------------------------
 !> @brief Calculate derived quantities from state.
 !--------------------------------------------------------------------------------------------------
-module subroutine dislotungsten_dependentState(instance,me)
+module subroutine dislotungsten_dependentState(ph,me)
 
   integer,      intent(in) :: &
-    instance, &
+    ph, &
     me
 
-  real(pReal), dimension(param(instance)%sum_N_sl) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl) :: &
     dislocationSpacing
 
-  associate(prm => param(instance), stt => state(instance),dst => dependentState(instance))
+  associate(prm => param(ph), stt => state(ph),dst => dependentState(ph))
 
   dislocationSpacing = sqrt(matmul(prm%forestProjection,stt%rho_mob(:,me)+stt%rho_dip(:,me)))
   dst%threshold_stress(:,me) = prm%mu*prm%b_sl &
@@ -402,14 +401,14 @@ end subroutine dislotungsten_dependentState
 !--------------------------------------------------------------------------------------------------
 !> @brief Write results to HDF5 output file.
 !--------------------------------------------------------------------------------------------------
-module subroutine plastic_dislotungsten_results(instance,group)
+module subroutine plastic_dislotungsten_results(ph,group)
 
-  integer,          intent(in) :: instance
+  integer,          intent(in) :: ph
   character(len=*), intent(in) :: group
 
   integer :: o
 
-  associate(prm => param(instance), stt => state(instance), dst => dependentState(instance))
+  associate(prm => param(ph), stt => state(ph), dst => dependentState(ph))
   outputsLoop: do o = 1,size(prm%output)
     select case(trim(prm%output(o)))
       case('rho_mob')
@@ -441,7 +440,7 @@ end subroutine plastic_dislotungsten_results
 ! NOTE: Against the common convention, the result (i.e. intent(out)) variables are the last to
 ! have the optional arguments at the end
 !--------------------------------------------------------------------------------------------------
-pure subroutine kinetics(Mp,T,instance,me, &
+pure subroutine kinetics(Mp,T,ph,me, &
                  dot_gamma_pos,dot_gamma_neg,ddot_gamma_dtau_pos,ddot_gamma_dtau_neg,tau_pos_out,tau_neg_out)
 
   real(pReal), dimension(3,3),  intent(in) :: &
@@ -449,18 +448,18 @@ pure subroutine kinetics(Mp,T,instance,me, &
   real(pReal),                  intent(in) :: &
     T                                                                                               !< temperature
   integer,                      intent(in) :: &
-    instance, &
+    ph, &
     me
 
-  real(pReal),                  intent(out), dimension(param(instance)%sum_N_sl) :: &
+  real(pReal),                  intent(out), dimension(param(ph)%sum_N_sl) :: &
     dot_gamma_pos, &
     dot_gamma_neg
-  real(pReal),                  intent(out), optional, dimension(param(instance)%sum_N_sl) :: &
+  real(pReal),                  intent(out), optional, dimension(param(ph)%sum_N_sl) :: &
     ddot_gamma_dtau_pos, &
     ddot_gamma_dtau_neg, &
     tau_pos_out, &
     tau_neg_out
-  real(pReal), dimension(param(instance)%sum_N_sl) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl) :: &
     StressRatio, &
     StressRatio_p,StressRatio_pminus1, &
     dvel, vel, &
@@ -469,7 +468,7 @@ pure subroutine kinetics(Mp,T,instance,me, &
     needsGoodName                                                                                   ! ToDo: @Karo: any idea?
   integer :: j
 
-  associate(prm => param(instance), stt => state(instance), dst => dependentState(instance))
+  associate(prm => param(ph), stt => state(ph), dst => dependentState(ph))
 
   do j = 1, prm%sum_N_sl
     tau_pos(j) = math_tensordot(Mp,prm%nonSchmid_pos(1:3,1:3,j))

src/phase_mechanical_plastic_dislotwin.f90

@@ -48,7 +48,7 @@ submodule(phase:plastic) dislotwin
       dot_N_0_tr, &                                                                                 !< trans nucleation rate [1/m³s] for each trans system
       t_tw, &                                                                                       !< twin thickness [m] for each twin system
       i_sl, &                                                                                       !< Adj. parameter for distance between 2 forest dislocations for each slip system
-      t_tr, &                                                                                       !< martensite lamellar thickness [m] for each trans system and instance
+      t_tr, &                                                                                       !< martensite lamellar thickness [m] for each trans system
       p, &                                                                                          !< p-exponent in glide velocity
       q, &                                                                                          !< q-exponent in glide velocity
       r, &                                                                                          !< r-exponent in twin nucleation rate
@@ -126,7 +126,6 @@ module function plastic_dislotwin_init() result(myPlasticity)
 
   logical, dimension(:), allocatable :: myPlasticity
   integer :: &
-    Ninstances, &
     p, i, &
     Nconstituents, &
     sizeState, sizeDotState, &
@@ -144,12 +143,12 @@ module function plastic_dislotwin_init() result(myPlasticity)
     mech, &
     pl
 
+
   myPlasticity = plastic_active('dislotwin')
-  Ninstances = count(myPlasticity)
+  if(count(myPlasticity) == 0) return
-  if(Ninstances == 0) return
 
   print'(/,a)', ' <<<+-  phase:mechanical:plastic:dislotwin init  -+>>>'
-  print'(a,i0)', ' # phases: ',Ninstances; flush(IO_STDOUT)
+  print'(a,i0)', ' # phases: ',count(myPlasticity); flush(IO_STDOUT)
 
   print*, 'Ma and Roters, Acta Materialia 52(12):3603–3612, 2004'
   print*, 'https://doi.org/10.1016/j.actamat.2004.04.012'//IO_EOL
@@ -160,18 +159,19 @@ module function plastic_dislotwin_init() result(myPlasticity)
   print*, 'Wong et al., Acta Materialia 118:140–151, 2016'
   print*, 'https://doi.org/10.1016/j.actamat.2016.07.032'
 
-  allocate(param(Ninstances))
-  allocate(state(Ninstances))
-  allocate(dotState(Ninstances))
-  allocate(dependentState(Ninstances))
 
   phases => config_material%get('phase')
-  i = 0
+  allocate(param(phases%length))
+  allocate(state(phases%length))
+  allocate(dotState(phases%length))
+  allocate(dependentState(phases%length))
+
+
   do p = 1, phases%length
+    if(.not. myPlasticity(p)) cycle
     phase => phases%get(p)
     mech  => phase%get('mechanics')
-    if(.not. myPlasticity(p)) cycle
+    i = p
-    i = i + 1
     associate(prm => param(i), &
               dot => dotState(i), &
               stt => state(i), &
@@ -496,8 +496,8 @@ module function plastic_dislotwin_homogenizedC(ph,me) result(homogenizedC)
   real(pReal) :: f_unrotated
 
 
-  associate(prm => param(phase_plasticInstance(ph)),&
+  associate(prm => param(ph),&
-            stt => state(phase_plasticInstance(ph)))
+            stt => state(ph))
 
   f_unrotated = 1.0_pReal &
               - sum(stt%f_tw(1:prm%sum_N_tw,me)) &
@@ -535,11 +535,11 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me)
      BoltzmannRatio, &
      ddot_gamma_dtau, &
      tau
-  real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl) :: &
     dot_gamma_sl,ddot_gamma_dtau_slip
-  real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_tw) :: &
+  real(pReal), dimension(param(ph)%sum_N_tw) :: &
     dot_gamma_twin,ddot_gamma_dtau_twin
-  real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_tr) :: &
+  real(pReal), dimension(param(ph)%sum_N_tr) :: &
     dot_gamma_tr,ddot_gamma_dtau_trans
   real(pReal):: dot_gamma_sb
   real(pReal), dimension(3,3) :: eigVectors, P_sb
@@ -564,7 +564,7 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me)
          0, 1, 1  &
          ],pReal),[ 3,6])
 
-  associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph)))
+  associate(prm => param(ph), stt => state(ph))
 
   f_unrotated = 1.0_pReal &
               - sum(stt%f_tw(1:prm%sum_N_tw,me)) &
@@ -573,7 +573,7 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me)
   Lp = 0.0_pReal
   dLp_dMp = 0.0_pReal
 
-  call kinetics_slip(Mp,T,phase_plasticInstance(ph),me,dot_gamma_sl,ddot_gamma_dtau_slip)
+  call kinetics_slip(Mp,T,ph,me,dot_gamma_sl,ddot_gamma_dtau_slip)
   slipContribution: do i = 1, prm%sum_N_sl
     Lp = Lp + dot_gamma_sl(i)*prm%P_sl(1:3,1:3,i)
     forall (k=1:3,l=1:3,m=1:3,n=1:3) &
@@ -581,7 +581,7 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me)
                        + ddot_gamma_dtau_slip(i) * prm%P_sl(k,l,i) * prm%P_sl(m,n,i)
   enddo slipContribution
 
-  call kinetics_twin(Mp,T,dot_gamma_sl,phase_plasticInstance(ph),me,dot_gamma_twin,ddot_gamma_dtau_twin)
+  call kinetics_twin(Mp,T,dot_gamma_sl,ph,me,dot_gamma_twin,ddot_gamma_dtau_twin)
   twinContibution: do i = 1, prm%sum_N_tw
     Lp = Lp + dot_gamma_twin(i)*prm%P_tw(1:3,1:3,i)
     forall (k=1:3,l=1:3,m=1:3,n=1:3) &
@@ -589,7 +589,7 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,me)
                        + ddot_gamma_dtau_twin(i)* prm%P_tw(k,l,i)*prm%P_tw(m,n,i)
   enddo twinContibution
 
-  call kinetics_trans(Mp,T,dot_gamma_sl,phase_plasticInstance(ph),me,dot_gamma_tr,ddot_gamma_dtau_trans)
+  call kinetics_trans(Mp,T,dot_gamma_sl,ph,me,dot_gamma_tr,ddot_gamma_dtau_trans)
   transContibution: do i = 1, prm%sum_N_tr
     Lp = Lp + dot_gamma_tr(i)*prm%P_tr(1:3,1:3,i)
     forall (k=1:3,l=1:3,m=1:3,n=1:3) &
@@ -653,24 +653,24 @@ module subroutine dislotwin_dotState(Mp,T,ph,me)
     tau, &
     sigma_cl, &                                                                                     !< climb stress
     b_d                                                                                             !< ratio of Burgers vector to stacking fault width
-  real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_sl) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl) :: &
     dot_rho_dip_formation, &
     dot_rho_dip_climb, &
     rho_dip_distance_min, &
     dot_gamma_sl
-  real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_tw) :: &
+  real(pReal), dimension(param(ph)%sum_N_tw) :: &
     dot_gamma_twin
-  real(pReal), dimension(param(phase_plasticInstance(ph))%sum_N_tr) :: &
+  real(pReal), dimension(param(ph)%sum_N_tr) :: &
     dot_gamma_tr
 
-  associate(prm => param(phase_plasticInstance(ph)),    stt => state(phase_plasticInstance(ph)), &
+  associate(prm => param(ph),    stt => state(ph), &
-            dot => dotState(phase_plasticInstance(ph)), dst => dependentState(phase_plasticInstance(ph)))
+            dot => dotState(ph), dst => dependentState(ph))
 
   f_unrotated = 1.0_pReal &
               - sum(stt%f_tw(1:prm%sum_N_tw,me)) &
               - sum(stt%f_tr(1:prm%sum_N_tr,me))
 
-  call kinetics_slip(Mp,T,phase_plasticInstance(ph),me,dot_gamma_sl)
+  call kinetics_slip(Mp,T,ph,me,dot_gamma_sl)
   dot%gamma_sl(:,me) = abs(dot_gamma_sl)
 
   rho_dip_distance_min = prm%D_a*prm%b_sl
@@ -721,10 +721,10 @@ module subroutine dislotwin_dotState(Mp,T,ph,me)
                     - 2.0_pReal*rho_dip_distance_min/prm%b_sl * stt%rho_dip(:,me)*abs(dot_gamma_sl) &
                     - dot_rho_dip_climb
 
-  call kinetics_twin(Mp,T,dot_gamma_sl,phase_plasticInstance(ph),me,dot_gamma_twin)
+  call kinetics_twin(Mp,T,dot_gamma_sl,ph,me,dot_gamma_twin)
   dot%f_tw(:,me) = f_unrotated*dot_gamma_twin/prm%gamma_char
 
-  call kinetics_trans(Mp,T,dot_gamma_sl,phase_plasticInstance(ph),me,dot_gamma_tr)
+  call kinetics_trans(Mp,T,dot_gamma_sl,ph,me,dot_gamma_tr)
   dot%f_tr(:,me) = f_unrotated*dot_gamma_tr
 
   end associate
@@ -735,33 +735,33 @@ end subroutine dislotwin_dotState
 !--------------------------------------------------------------------------------------------------
 !> @brief Calculate derived quantities from state.
 !--------------------------------------------------------------------------------------------------
-module subroutine dislotwin_dependentState(T,instance,me)
+module subroutine dislotwin_dependentState(T,ph,me)
 
   integer,       intent(in) :: &
-    instance, &
+    ph, &
     me
   real(pReal),   intent(in) :: &
     T
 
   real(pReal) :: &
     sumf_twin,Gamma,sumf_trans
-  real(pReal), dimension(param(instance)%sum_N_sl) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl) :: &
     inv_lambda_sl_sl, &                                                                             !< 1/mean free distance between 2 forest dislocations seen by a moving dislocation
     inv_lambda_sl_tw, &                                                                             !< 1/mean free distance between 2 twin stacks from different systems seen by a moving dislocation
     inv_lambda_sl_tr                                                                                !< 1/mean free distance between 2 martensite lamellar from different systems seen by a moving dislocation
-  real(pReal), dimension(param(instance)%sum_N_tw) :: &
+  real(pReal), dimension(param(ph)%sum_N_tw) :: &
     inv_lambda_tw_tw, &                                                                             !< 1/mean free distance between 2 twin stacks from different systems seen by a growing twin
     f_over_t_tw
-   real(pReal), dimension(param(instance)%sum_N_tr) :: &
+   real(pReal), dimension(param(ph)%sum_N_tr) :: &
     inv_lambda_tr_tr, &                                                                             !< 1/mean free distance between 2 martensite stacks from different systems seen by a growing martensite
     f_over_t_tr
   real(pReal), dimension(:), allocatable :: &
     x0
 
 
-  associate(prm => param(instance),&
+  associate(prm => param(ph),&
-            stt => state(instance),&
+            stt => state(ph),&
-            dst => dependentState(instance))
+            dst => dependentState(ph))
 
   sumf_twin  = sum(stt%f_tw(1:prm%sum_N_tw,me))
   sumf_trans = sum(stt%f_tr(1:prm%sum_N_tr,me))
@@ -827,14 +827,14 @@ end subroutine dislotwin_dependentState
 !--------------------------------------------------------------------------------------------------
 !> @brief Write results to HDF5 output file.
 !--------------------------------------------------------------------------------------------------
-module subroutine plastic_dislotwin_results(instance,group)
+module subroutine plastic_dislotwin_results(ph,group)
 
-  integer,          intent(in) :: instance
+  integer,          intent(in) :: ph
   character(len=*), intent(in) :: group
 
   integer :: o
 
- associate(prm => param(instance), stt => state(instance), dst => dependentState(instance))
+ associate(prm => param(ph), stt => state(ph), dst => dependentState(ph))
   outputsLoop: do o = 1,size(prm%output)
     select case(trim(prm%output(o)))
 
@@ -882,7 +882,7 @@ end subroutine plastic_dislotwin_results
 ! NOTE: Against the common convention, the result (i.e. intent(out)) variables are the last to
 ! have the optional arguments at the end
 !--------------------------------------------------------------------------------------------------
-pure subroutine kinetics_slip(Mp,T,instance,me, &
+pure subroutine kinetics_slip(Mp,T,ph,me, &
                               dot_gamma_sl,ddot_gamma_dtau_slip,tau_slip)
 
   real(pReal), dimension(3,3),  intent(in) :: &
@@ -890,18 +890,18 @@ pure subroutine kinetics_slip(Mp,T,instance,me, &
   real(pReal),                  intent(in) :: &
     T                                                                                               !< temperature
   integer,                      intent(in) :: &
-    instance, &
+    ph, &
     me
 
-  real(pReal), dimension(param(instance)%sum_N_sl), intent(out) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl), intent(out) :: &
     dot_gamma_sl
-  real(pReal), dimension(param(instance)%sum_N_sl), optional, intent(out) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl), optional, intent(out) :: &
     ddot_gamma_dtau_slip, &
     tau_slip
-  real(pReal), dimension(param(instance)%sum_N_sl) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl) :: &
     ddot_gamma_dtau
 
-  real(pReal), dimension(param(instance)%sum_N_sl) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl) :: &
     tau, &
     stressRatio, &
     StressRatio_p, &
@@ -914,7 +914,7 @@ pure subroutine kinetics_slip(Mp,T,instance,me, &
     tau_eff                                                                                         !< effective resolved stress
   integer :: i
 
-  associate(prm => param(instance), stt => state(instance), dst => dependentState(instance))
+  associate(prm => param(ph), stt => state(ph), dst => dependentState(ph))
 
   do i = 1, prm%sum_N_sl
     tau(i) = math_tensordot(Mp,prm%P_sl(1:3,1:3,i))
@@ -959,7 +959,7 @@ end subroutine kinetics_slip
 ! NOTE: Against the common convention, the result (i.e. intent(out)) variables are the last to
 ! have the optional arguments at the end.
 !--------------------------------------------------------------------------------------------------
-pure subroutine kinetics_twin(Mp,T,dot_gamma_sl,instance,me,&
+pure subroutine kinetics_twin(Mp,T,dot_gamma_sl,ph,me,&
                               dot_gamma_twin,ddot_gamma_dtau_twin)
 
   real(pReal), dimension(3,3),  intent(in) :: &
@@ -967,17 +967,17 @@ pure subroutine kinetics_twin(Mp,T,dot_gamma_sl,instance,me,&
   real(pReal),                  intent(in) :: &
     T                                                                                               !< temperature
   integer,                      intent(in) :: &
-    instance, &
+    ph, &
     me
-  real(pReal), dimension(param(instance)%sum_N_sl), intent(in) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl), intent(in) :: &
     dot_gamma_sl
 
-  real(pReal), dimension(param(instance)%sum_N_tw), intent(out) :: &
+  real(pReal), dimension(param(ph)%sum_N_tw), intent(out) :: &
     dot_gamma_twin
-  real(pReal), dimension(param(instance)%sum_N_tw), optional, intent(out) :: &
+  real(pReal), dimension(param(ph)%sum_N_tw), optional, intent(out) :: &
     ddot_gamma_dtau_twin
 
-  real, dimension(param(instance)%sum_N_tw) :: &
+  real, dimension(param(ph)%sum_N_tw) :: &
     tau, &
     Ndot0, &
     stressRatio_r, &
@@ -985,7 +985,7 @@ pure subroutine kinetics_twin(Mp,T,dot_gamma_sl,instance,me,&
 
   integer :: i,s1,s2
 
-  associate(prm => param(instance), stt => state(instance), dst => dependentState(instance))
+  associate(prm => param(ph), stt => state(ph), dst => dependentState(ph))
 
   do i = 1, prm%sum_N_tw
     tau(i) = math_tensordot(Mp,prm%P_tw(1:3,1:3,i))
@@ -1028,7 +1028,7 @@ end subroutine kinetics_twin
 ! NOTE: Against the common convention, the result (i.e. intent(out)) variables are the last to
 ! have the optional arguments at the end.
 !--------------------------------------------------------------------------------------------------
-pure subroutine kinetics_trans(Mp,T,dot_gamma_sl,instance,me,&
+pure subroutine kinetics_trans(Mp,T,dot_gamma_sl,ph,me,&
                               dot_gamma_tr,ddot_gamma_dtau_trans)
 
   real(pReal), dimension(3,3),  intent(in) :: &
@@ -1036,24 +1036,24 @@ pure subroutine kinetics_trans(Mp,T,dot_gamma_sl,instance,me,&
   real(pReal),                  intent(in) :: &
     T                                                                                               !< temperature
   integer,                      intent(in) :: &
-    instance, &
+    ph, &
     me
-  real(pReal), dimension(param(instance)%sum_N_sl), intent(in) :: &
+  real(pReal), dimension(param(ph)%sum_N_sl), intent(in) :: &
     dot_gamma_sl
 
-  real(pReal), dimension(param(instance)%sum_N_tr), intent(out) :: &
+  real(pReal), dimension(param(ph)%sum_N_tr), intent(out) :: &
     dot_gamma_tr
-  real(pReal), dimension(param(instance)%sum_N_tr), optional, intent(out) :: &
+  real(pReal), dimension(param(ph)%sum_N_tr), optional, intent(out) :: &
     ddot_gamma_dtau_trans
 
-  real, dimension(param(instance)%sum_N_tr) :: &
+  real, dimension(param(ph)%sum_N_tr) :: &
     tau, &
     Ndot0, &
     stressRatio_s, &
     ddot_gamma_dtau
 
   integer :: i,s1,s2
-  associate(prm => param(instance), stt => state(instance), dst => dependentState(instance))
+  associate(prm => param(ph), stt => state(ph), dst => dependentState(ph))
 
   do i = 1, prm%sum_N_tr
     tau(i) = math_tensordot(Mp,prm%P_tr(1:3,1:3,i))

src/phase_mechanical_plastic_isotropic.f90

@@ -53,7 +53,6 @@ module function plastic_isotropic_init() result(myPlasticity)
 
   logical, dimension(:), allocatable :: myPlasticity
   integer :: &
-    Ninstances, &
     p, &
     i, &
     Nconstituents, &
@@ -68,28 +67,26 @@ module function plastic_isotropic_init() result(myPlasticity)
     mech, &
     pl
 
+
   myPlasticity = plastic_active('isotropic')
-  Ninstances = count(myPlasticity)
+  if(count(myPlasticity) == 0) return
-  if(Ninstances == 0) return
 
   print'(/,a)', ' <<<+-  phase:mechanical:plastic:isotropic init  -+>>>'
-  print'(a,i0)', ' # phses: ',Ninstances; flush(IO_STDOUT)
+  print'(a,i0)', ' # phases: ',count(myPlasticity); flush(IO_STDOUT)
-
 
   print*, 'Maiti and Eisenlohr, Scripta Materialia 145:37–40, 2018'
   print*, 'https://doi.org/10.1016/j.scriptamat.2017.09.047'
 
-  allocate(param(Ninstances))
-  allocate(state(Ninstances))
-  allocate(dotState(Ninstances))
-
   phases => config_material%get('phase')
-  i = 0
+  allocate(param(phases%length))
+  allocate(state(phases%length))
+  allocate(dotState(phases%length))
+
   do p = 1, phases%length
+    if(.not. myPlasticity(p)) cycle
     phase => phases%get(p)
     mech  => phase%get('mechanics')
-    if(.not. myPlasticity(p)) cycle
+    i = p
-    i = i + 1
     associate(prm => param(i), &
               dot => dotState(i), &
               stt => state(i))
@@ -191,7 +188,7 @@ module subroutine isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,ph,me)
   integer :: &
     k, l, m, n
 
-  associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph)))
+  associate(prm => param(ph), stt => state(ph))
 
   Mp_dev = math_deviatoric33(Mp)
   squarenorm_Mp_dev = math_tensordot(Mp_dev,Mp_dev)
@@ -221,7 +218,7 @@ end subroutine isotropic_LpAndItsTangent
 !--------------------------------------------------------------------------------------------------
 !> @brief Calculate inelastic velocity gradient and its tangent.
 !--------------------------------------------------------------------------------------------------
-module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,me)
+module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,ph,me)
 
   real(pReal), dimension(3,3),     intent(out) :: &
     Li                                                                                              !< inleastic velocity gradient
@@ -231,7 +228,7 @@ module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,me)
   real(pReal), dimension(3,3), intent(in) :: &
     Mi                                                                                              !< Mandel stress
   integer,                     intent(in) :: &
-    instance, &
+    ph, &
     me
 
   real(pReal) :: &
@@ -239,7 +236,7 @@ module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,me)
   integer :: &
     k, l, m, n
 
-  associate(prm => param(instance), stt => state(instance))
+  associate(prm => param(ph), stt => state(ph))
 
   tr=math_trace33(math_spherical33(Mi))
 
@@ -276,8 +273,8 @@ module subroutine isotropic_dotState(Mp,ph,me)
     xi_inf_star, &                                                                                  !< saturation xi
     norm_Mp                                                                                         !< norm of the (deviatoric) Mandel stress
 
-  associate(prm => param(phase_plasticInstance(ph)), stt => state(phase_plasticInstance(ph)), &
+  associate(prm => param(ph), stt => state(ph), &
-   dot => dotState(phase_plasticInstance(ph)))
+   dot => dotState(ph))
 
   if (prm%dilatation) then
     norm_Mp = sqrt(math_tensordot(Mp,Mp))
@@ -313,14 +310,14 @@ end subroutine isotropic_dotState
 !--------------------------------------------------------------------------------------------------
 !> @brief Write results to HDF5 output file.
 !--------------------------------------------------------------------------------------------------
-module subroutine plastic_isotropic_results(instance,group)
+module subroutine plastic_isotropic_results(ph,group)
 
-  integer,          intent(in) :: instance
+  integer,          intent(in) :: ph
   character(len=*), intent(in) :: group
 
   integer :: o
 
-  associate(prm => param(instance), stt => state(instance))
+  associate(prm => param(ph), stt => state(ph))
   outputsLoop: do o = 1,size(prm%output)
     select case(trim(prm%output(o)))
       case ('xi')

src/phase_mechanical_plastic_nonlocal.f90

@@ -521,7 +521,7 @@ module function plastic_nonlocal_init() result(myPlasticity)
     if(.not. myPlasticity(p)) cycle
     i = i + 1
 
-    Nconstituents  = count(material_phaseAt==p) * discretization_nIPs
+    Nconstituents = count(material_phaseAt2 == p)
     l = 0
     do t = 1,4
       do s = 1,param(i)%sum_N_sl