does not make sense to store and use the 6-vector version of the Schmid

matrix

src/plastic_phenopowerlaw.f90

@@ -68,9 +68,6 @@ module plastic_phenopowerlaw
      interaction_SlipTwin, &                                                                        !< slip resistance from twin activity
      interaction_TwinSlip, &                                                                        !< twin resistance from slip activity
      interaction_TwinTwin                                                                           !< twin resistance from twin activity
-   real(pReal),   dimension(:,:),   allocatable :: &
-     Schmid_slip6, &
-     Schmid_twin6
    real(pReal),   dimension(:,:,:),   allocatable :: &
      Schmid_slip, &
      Schmid_twin
@@ -363,7 +360,6 @@ subroutine plastic_phenopowerlaw_init
    allocate(temp1(prm%totalNslip,prm%totalNslip),source = 0.0_pReal)
    allocate(temp2(prm%totalNslip,prm%totalNtwin),source = 0.0_pReal)
    allocate(prm%Schmid_slip(3,3,prm%totalNslip),source   = 0.0_pReal)
-   allocate(prm%Schmid_slip6(6,prm%totalNslip),source   = 0.0_pReal)
    allocate(prm%nonSchmid_pos(3,3,size(prm%nonSchmidCoeff)+1,prm%totalNslip),source = 0.0_pReal)
    allocate(prm%nonSchmid_neg(3,3,size(prm%nonSchmidCoeff)+1,prm%totalNslip),source = 0.0_pReal)
    i = 0_pInt
@@ -373,7 +369,6 @@ subroutine plastic_phenopowerlaw_init
      mySlipSystems: do j = 1_pInt,prm%Nslip(f)
        i = i + 1_pInt
        prm%Schmid_slip(1:3,1:3,i) = lattice_Sslip(1:3,1:3,1,sum(lattice_Nslipsystem(1:f-1,p))+j,p)
-       prm%Schmid_slip6(1:6,i)    = lattice_Sslip_v(1:6,1,sum(lattice_Nslipsystem(1:f-1,p))+j,p)
        !prm%nonSchmid_pos(1:3,1:3,1,i) = lattice_Sslip(1:3,1:3,1,sum(lattice_Nslipsystem(1:f-1,p))+j,p)
        !prm%nonSchmid_neg(1:3,1:3,1,i) = lattice_Sslip(1:3,1:3,1,sum(lattice_Nslipsystem(1:f-1,p))+j,p)
        !do k = 1,size(prm%nonSchmidCoeff)
@@ -410,7 +405,6 @@ subroutine plastic_phenopowerlaw_init
    allocate(temp1(prm%totalNtwin,prm%totalNslip),source = 0.0_pReal)
    allocate(temp2(prm%totalNtwin,prm%totalNtwin),source = 0.0_pReal)
    allocate(prm%Schmid_twin(3,3,prm%totalNtwin),source  = 0.0_pReal)
-   allocate(prm%Schmid_twin6(6,prm%totalNtwin),source  = 0.0_pReal)
    allocate(prm%shear_twin(prm%totalNtwin),source       = 0.0_pReal)
    i = 0_pInt
    myTwinFamilies: do f = 1_pInt,size(prm%Ntwin,1)                                    ! >>> interaction twin -- X
@@ -418,7 +412,6 @@ subroutine plastic_phenopowerlaw_init
      myTwinSystems: do j = 1_pInt,prm%Ntwin(f)
        i = i + 1_pInt
        prm%Schmid_twin(1:3,1:3,i) = lattice_Stwin(1:3,1:3,sum(lattice_NTwinsystem(1:f-1,p))+j,p)
-       prm%Schmid_twin6(1:6,i)    = lattice_Stwin_v(1:6,sum(lattice_Ntwinsystem(1:f-1,p))+j,p)
        prm%shear_twin(i)          = lattice_shearTwin(sum(lattice_Ntwinsystem(1:f-1,p))+j,p)
        slipFamilies: do o = 1_pInt,size(prm%Nslip,1)
          index_otherFamily = sum(prm%Nslip(1:o-1_pInt))
@@ -505,7 +498,7 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dMstar99,Mstar_v,ipc,ip,
    dNeq0
  use math, only: &
    math_mul33xx33,&
-   math_Mandel33to6, &
+   math_Mandel6to33, &
    math_Plain3333to99
  use material, only: &
    phasememberAt, &
@@ -534,6 +527,8 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dMstar99,Mstar_v,ipc,ip,
    gdot_slip_pos,gdot_slip_neg, &
    dgdot_dtauslip_pos,dgdot_dtauslip_neg, &
    gdot_twin,dgdot_dtautwin,tau_twin
+ real(pReal), dimension(3,3) :: &
+   S                                                                                                !< Second-Piola Kirchhoff stress
  real(pReal), dimension(3,3,3,3) :: &
    dLp_dMstar                                                                                       !< derivative of Lp with respect to Mstar as 4th order tensor
  type(tParameters)          :: prm
@@ -547,11 +542,12 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dMstar99,Mstar_v,ipc,ip,
  Lp = 0.0_pReal
  dLp_dMstar = 0.0_pReal
 
+ S = math_Mandel6to33(Mstar_v)
 !--------------------------------------------------------------------------------------------------
 ! Slip part
  do j = 1_pInt, prm%totalNslip
 
-   tau_slip_pos  = dot_product(Mstar_v,prm%Schmid_slip6(1:6,j))
+   tau_slip_pos  = math_mul33xx33(S,prm%Schmid_slip(1:3,1:3,j))
    tau_slip_neg  = tau_slip_pos
    !do k = 1,size(prm%nonSchmidCoeff)
    !  tau_slip_pos = tau_slip_pos &
@@ -585,7 +581,7 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dMstar99,Mstar_v,ipc,ip,
 ! Twinning part
  do j = 1_pInt, prm%totalNtwin
 
-   tau_twin  = dot_product(Mstar_v,prm%Schmid_twin6(1:6,j))
+   tau_twin  = math_mul33xx33(S,prm%Schmid_twin(1:3,1:3,j))
    gdot_twin = (1.0_pReal-stt%sumF(of))*prm%gdot0_twin*(abs(tau_twin)/stt%s_twin(j,of))**prm%n_twin&
              * max(0.0_pReal,sign(1.0_pReal,tau_twin))
    Lp = Lp + gdot_twin*prm%Schmid_twin(1:3,1:3,j)
@@ -634,8 +630,8 @@ subroutine plastic_phenopowerlaw_dotState(Mstar6,ipc,ip,el)
    ssat_offset, &
    tau_slip_pos,tau_slip_neg,tau_twin
 
- !real(pReal), dimension(3,3) :: &
- !  Mstar
+ real(pReal), dimension(3,3) :: &
+   S                                                                                                !< Second-Piola Kirchhoff stress
  real(pReal), dimension(param(phase_plasticityInstance(material_phase(ipc,ip,el)))%totalNslip) :: &
    gdot_slip,left_SlipSlip,right_SlipSlip
  real(pReal), dimension(param(phase_plasticityInstance(material_phase(ipc,ip,el)))%totalNtwin) :: &
@@ -650,7 +646,7 @@ subroutine plastic_phenopowerlaw_dotState(Mstar6,ipc,ip,el)
            dst => dotState(phase_plasticityInstance(material_phase(ipc,ip,el))))
 
  dst%whole(:,of) = 0.0_pReal
- !Mstar = math_Mandel6to33(Mstar6)
+ S = math_Mandel6to33(Mstar6)
 
 !--------------------------------------------------------------------------------------------------
 ! system-independent (nonlinear) prefactors to M_Xx (X influenced by x) matrices
@@ -666,7 +662,7 @@ subroutine plastic_phenopowerlaw_dotState(Mstar6,ipc,ip,el)
    right_SlipSlip(j) = abs(1.0_pReal-stt%s_slip(j,of) / (prm%tausat_slip(j)+ssat_offset)) **prm%a_slip &
                      * sign(1.0_pReal,1.0_pReal-stt%s_slip(j,of) / (prm%tausat_slip(j)+ssat_offset))
 
-   tau_slip_pos  = dot_product(Mstar6,prm%Schmid_slip6(1:6,j))
+   tau_slip_pos  = math_mul33xx33(S,prm%Schmid_slip(1:3,1:3,j))
    tau_slip_neg  = tau_slip_pos
    !nonSchmidSystems: do k = 1,size(prm%nonSchmidCoeff)
    !  tau_slip_pos = tau_slip_pos + math_mul33xx33(Mstar,prm%nonSchmid_pos(1:3,1:3,k,j))
@@ -678,7 +674,7 @@ subroutine plastic_phenopowerlaw_dotState(Mstar6,ipc,ip,el)
  enddo
 
  do j = 1_pInt, prm%totalNtwin
-   tau_twin  = dot_product(Mstar6,prm%Schmid_twin6(1:6,j))
+   tau_twin  = math_mul33xx33(S,prm%Schmid_twin(1:3,1:3,j))
    gdot_twin(j) = (1.0_pReal-stt%sumF(of))*prm%gdot0_twin* abs(tau_twin/stt%s_twin(j,of))**prm%n_twin & !ToDo: save to dotState
                 * max(0.0_pReal,sign(1.0_pReal,tau_twin))
  enddo
@@ -733,8 +729,8 @@ function plastic_phenopowerlaw_postResults(Mstar6,ipc,ip,el)
    ip, &                                                                                            !< integration point
    el                                                                                               !< element                                                                                        !< microstructure state
 
- !real(pReal), dimension(3,3) :: &
- !  Mstar
+ real(pReal), dimension(3,3) :: &
+   S                                                                                                !< Second-Piola Kirchhoff stress
  real(pReal), dimension(plasticState(material_phase(ipc,ip,el))%sizePostResults) :: &
    plastic_phenopowerlaw_postResults
 
@@ -753,6 +749,7 @@ function plastic_phenopowerlaw_postResults(Mstar6,ipc,ip,el)
 
  plastic_phenopowerlaw_postResults = 0.0_pReal
  c = 0_pInt
+ S = math_Mandel6to33(Mstar6)
 
  outputsLoop: do o = 1_pInt,size(prm%outputID)
    select case(prm%outputID(o))
@@ -766,7 +763,7 @@ function plastic_phenopowerlaw_postResults(Mstar6,ipc,ip,el)
 
      case (shearrate_slip_ID)
        do j = 1_pInt, prm%totalNslip
-         tau_slip_pos  = dot_product(Mstar6,prm%Schmid_slip6(1:6,j))
+         tau_slip_pos  = math_mul33xx33(S,prm%Schmid_slip(1:3,1:3,j))
          tau_slip_neg  = tau_slip_pos
          !nonSchmidSystems: do k = 1,size(prm%nonSchmidCoeff)
          !  tau_slip_pos = tau_slip_pos + math_mul33xx33(Mstar,prm%nonSchmid_pos(1:3,1:3,k,j))
@@ -780,7 +777,7 @@ function plastic_phenopowerlaw_postResults(Mstar6,ipc,ip,el)
 
      case (resolvedstress_slip_ID)
        do j = 1_pInt, prm%totalNslip
-         plastic_phenopowerlaw_postResults(c+j) = dot_product(Mstar6,prm%Schmid_slip6(1:6,j))
+         plastic_phenopowerlaw_postResults(c+j) = math_mul33xx33(S,prm%Schmid_slip(1:3,1:3,j))
        enddo
        c = c + prm%totalNslip
 
@@ -800,7 +797,7 @@ function plastic_phenopowerlaw_postResults(Mstar6,ipc,ip,el)
 
      case (shearrate_twin_ID)
        do j = 1_pInt, prm%totalNtwin
-         tau_twin  = dot_product(Mstar6,prm%Schmid_twin6(1:6,j))
+         tau_twin  = math_mul33xx33(S,prm%Schmid_twin(1:3,1:3,j))
          plastic_phenopowerlaw_postResults(c+j) = (1.0_pReal-stt%sumF(of))*&  ! 1-F
                                                  prm%gdot0_twin*(abs(tau_twin)/stt%s_twin(j,of))**&
                                            prm%n_twin*max(0.0_pReal,sign(1.0_pReal,tau_twin))
@@ -809,7 +806,7 @@ function plastic_phenopowerlaw_postResults(Mstar6,ipc,ip,el)
 
      case (resolvedstress_twin_ID)
        do j = 1_pInt, prm%totalNtwin
-         plastic_phenopowerlaw_postResults(c+j) = dot_product(Mstar6,prm%Schmid_twin6(1:6,j))
+         plastic_phenopowerlaw_postResults(c+j) = math_mul33xx33(S,prm%Schmid_twin(1:3,1:3,j))
        enddo
        c = c + prm%totalNtwin