Elasticity matrices for transformed phase

code/IO.f90

@@ -1552,6 +1552,8 @@ subroutine IO_error(error_ID,el,ip,g,ext_msg)
    msg = 'transformed hex lattice structure with invalid c/a ratio'
  case (135_pInt)
    msg = 'zero entry on stiffness diagonal'
+ case (136_pInt)
+   msg = 'zero entry on stiffness diagonal for transformed phase'
 
 !--------------------------------------------------------------------------------------------------
 ! material error messages and related messages in mesh

code/lattice.f90

@@ -1642,6 +1642,37 @@ subroutine lattice_initializeStructure(myPhase,CoverA,CoverA_trans,a_fcc,a_bcc)
    if (abs(lattice_C66(i,i,myPhase))<tol_math_check) &
      call IO_error(135_pInt,el=i,ip=myPhase,ext_msg='matrix diagonal "el"ement of phase "ip"')
  enddo
+
+ ! Elasticity matrices for transformed phase
+ select case(lattice_structure(myPhase))
+   case (LATTICE_fcc_ID)
+     select case(trans_lattice_structure(myPhase))
+       case (LATTICE_bcc_ID)
+         lattice_trans_C66(1:6,1:6,myPhase) = lattice_C66(1:6,1:6,myPhase)
+         lattice_trans_mu(myPhase) = lattice_mu(myPhase)
+         lattice_trans_nu(myPhase) = lattice_nu(myPhase)
+         lattice_trans_C3333(1:3,1:3,1:3,1:3,myPhase) = lattice_C3333(1:3,1:3,1:3,1:3,myPhase)
+       case (LATTICE_hex_ID)
+         lattice_trans_C66(1:6,1:6,myPhase) = lattice_symmetrizeC66(trans_lattice_structure(myPhase),&
+                                                                   lattice_trans_C66(1:6,1:6,myPhase))
+         lattice_trans_mu(myPhase) = 0.2_pReal *(  lattice_trans_C66(1,1,myPhase) &
+                                                 - lattice_trans_C66(1,2,myPhase) &
+                                                 + 3.0_pReal*lattice_trans_C66(4,4,myPhase))
+         lattice_trans_nu(myPhase) = (  lattice_trans_C66(1,1,myPhase) &
+                                      + 4.0_pReal*lattice_trans_C66(1,2,myPhase) &
+                                      - 2.0_pReal*lattice_trans_C66(4,4,myPhase)) &
+                                                     /(  4.0_pReal*lattice_trans_C66(1,1,myPhase) &
+                                                       + 6.0_pReal*lattice_trans_C66(1,2,myPhase) &
+                                                       + 2.0_pReal*lattice_trans_C66(4,4,myPhase))
+         lattice_trans_C3333(1:3,1:3,1:3,1:3,myPhase) = math_Voigt66to3333(lattice_trans_C66(1:6,1:6,myPhase))
+         lattice_trans_C66(1:6,1:6,myPhase) = math_Mandel3333to66(lattice_trans_C3333(1:3,1:3,1:3,1:3,myPhase))
+         do i = 1_pInt, 6_pInt
+           if (abs(lattice_trans_C66(i,i,myPhase))<tol_math_check) &
+           call IO_error(136_pInt,el=i,ip=myPhase,ext_msg='matrix diagonal "el"ement of phase "ip"')
+         enddo
+     end select
+ end select
+
  lattice_thermalConductivity33(1:3,1:3,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),&
                                                                        lattice_thermalConductivity33(1:3,1:3,myPhase))
  lattice_thermalExpansion33(1:3,1:3,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),&

code/plastic_dislotwin.f90

@@ -982,7 +982,7 @@ subroutine plastic_dislotwin_init(fileUnit)
            do p = 1_pInt,3_pInt; do q = 1_pInt,3_pInt; do r = 1_pInt,3_pInt; do s = 1_pInt,3_pInt
              plastic_dislotwin_Ctrans3333(l,m,n,o,index_myFamily+j,instance) = &
              plastic_dislotwin_Ctrans3333(l,m,n,o,index_myFamily+j,instance) + &
-               lattice_C3333(p,q,r,s,instance) * &
+               lattice_trans_C3333(p,q,r,s,instance) * &
                lattice_Qtrans(l,p,index_otherFamily+j,phase) * &
                lattice_Qtrans(m,q,index_otherFamily+j,phase) * &
                lattice_Qtrans(n,r,index_otherFamily+j,phase) * &
@@ -993,6 +993,8 @@ subroutine plastic_dislotwin_init(fileUnit)
            math_Mandel3333to66(plastic_dislotwin_Ctrans3333(1:3,1:3,1:3,1:3,index_myFamily+j,instance))
 
          !* Projection matrices for shear from slip systems to fault-band (twin) systems for strain-induced martensite nucleation
+         select case(trans_lattice_structure(phase))
+           case (LATTICE_bcc_ID)
              do o = 1_pInt,lattice_maxNtransFamily
                index_otherFamily = sum(plastic_dislotwin_Nslip(1:o-1_pInt,instance))
                do k = 1_pInt,plastic_dislotwin_Nslip(o,instance)                                   ! loop over (active) systems in other family (trans)
@@ -1000,6 +1002,7 @@ subroutine plastic_dislotwin_init(fileUnit)
                        lattice_projectionTrans( sum(lattice_NtransSystem(1:f-1,phase))+j, &
                                                 sum(lattice_NslipSystem(1:o-1,phase))+k, phase)
              enddo; enddo
+         end select 
 
        enddo transSystemsLoop
      enddo transFamiliesLoop