diff --git a/code/lattice.f90 b/code/lattice.f90
index 37df1b9b3..26ff927b3 100644
--- a/code/lattice.f90
+++ b/code/lattice.f90
@@ -657,14 +657,14 @@ function lattice_initializeStructure(struct,CoverA)
      if (lattice_fcc_Nstructure == 1_pInt) then
        processMe = .true.
        do i = 1,myNslip
-         sn(:,i) = lattice_fcc_systemSlip(1:3,i)/dsqrt(math_mul3x3(lattice_fcc_systemSlip(1:3,i),lattice_fcc_systemSlip(1:3,i)))
-         sd(:,i) = lattice_fcc_systemSlip(4:6,i)/dsqrt(math_mul3x3(lattice_fcc_systemSlip(4:6,i),lattice_fcc_systemSlip(4:6,i)))
-         st(:,i) = math_vectorproduct(sn(:,i),sd(:,i))
+         sd(:,i) = lattice_fcc_systemSlip(1:3,i)/dsqrt(math_mul3x3(lattice_fcc_systemSlip(1:3,i),lattice_fcc_systemSlip(1:3,i)))
+         sn(:,i) = lattice_fcc_systemSlip(4:6,i)/dsqrt(math_mul3x3(lattice_fcc_systemSlip(4:6,i),lattice_fcc_systemSlip(4:6,i)))
+         st(:,i) = math_vectorproduct(sd(:,i),sn(:,i))
        enddo
        do i = 1,myNtwin
-         tn(:,i) = lattice_fcc_systemTwin(1:3,i)/dsqrt(math_mul3x3(lattice_fcc_systemTwin(1:3,i),lattice_fcc_systemTwin(1:3,i)))
-         td(:,i) = lattice_fcc_systemTwin(4:6,i)/dsqrt(math_mul3x3(lattice_fcc_systemTwin(4:6,i),lattice_fcc_systemTwin(4:6,i)))
-         tt(:,i) = math_vectorproduct(tn(:,i),td(:,i))
+         td(:,i) = lattice_fcc_systemTwin(1:3,i)/dsqrt(math_mul3x3(lattice_fcc_systemTwin(1:3,i),lattice_fcc_systemTwin(1:3,i)))
+         tn(:,i) = lattice_fcc_systemTwin(4:6,i)/dsqrt(math_mul3x3(lattice_fcc_systemTwin(4:6,i),lattice_fcc_systemTwin(4:6,i)))
+         tt(:,i) = math_vectorproduct(td(:,i),tn(:,i))
          ts(i)   = lattice_fcc_shearTwin(i)
        enddo
        interactionSlipSlip => lattice_fcc_interactionSlipSlip
@@ -681,14 +681,14 @@ function lattice_initializeStructure(struct,CoverA)
      if (lattice_bcc_Nstructure == 1_pInt) then
        processMe = .true.
        do i = 1,myNslip
-         sn(:,i) = lattice_bcc_systemSlip(1:3,i)/dsqrt(math_mul3x3(lattice_bcc_systemSlip(1:3,i),lattice_bcc_systemSlip(1:3,i)))
-         sd(:,i) = lattice_bcc_systemSlip(4:6,i)/dsqrt(math_mul3x3(lattice_bcc_systemSlip(4:6,i),lattice_bcc_systemSlip(4:6,i)))
-         st(:,i) = math_vectorproduct(sn(:,i),sd(:,i))
+         sd(:,i) = lattice_bcc_systemSlip(1:3,i)/dsqrt(math_mul3x3(lattice_bcc_systemSlip(1:3,i),lattice_bcc_systemSlip(1:3,i)))
+         sn(:,i) = lattice_bcc_systemSlip(4:6,i)/dsqrt(math_mul3x3(lattice_bcc_systemSlip(4:6,i),lattice_bcc_systemSlip(4:6,i)))
+         st(:,i) = math_vectorproduct(sd(:,i),sn(:,i))
        enddo
        do i = 1,myNtwin
-         tn(:,i) = lattice_bcc_systemTwin(1:3,i)/dsqrt(math_mul3x3(lattice_bcc_systemTwin(1:3,i),lattice_bcc_systemTwin(1:3,i)))
-         td(:,i) = lattice_bcc_systemTwin(4:6,i)/dsqrt(math_mul3x3(lattice_bcc_systemTwin(4:6,i),lattice_bcc_systemTwin(4:6,i)))
-         tt(:,i) = math_vectorproduct(tn(:,i),td(:,i))
+         td(:,i) = lattice_bcc_systemTwin(1:3,i)/dsqrt(math_mul3x3(lattice_bcc_systemTwin(1:3,i),lattice_bcc_systemTwin(1:3,i)))
+         tn(:,i) = lattice_bcc_systemTwin(4:6,i)/dsqrt(math_mul3x3(lattice_bcc_systemTwin(4:6,i),lattice_bcc_systemTwin(4:6,i)))
+         tt(:,i) = math_vectorproduct(td(:,i),tn(:,i))
          ts(i)   = lattice_bcc_shearTwin(i)
        enddo
        interactionSlipSlip => lattice_bcc_interactionSlipSlip
@@ -706,28 +706,28 @@ function lattice_initializeStructure(struct,CoverA)
        processMe = .true.
 ! converting from 4 axes coordinate system (a1=a2=a3=c) to ortho-hexgonal system (a, b, c)
        do i = 1,myNslip
-         hex_n(1) =  lattice_hex_systemSlip(1,i)           ! plane (hkil)->(h (h+2k)/sqrt(3) l/(c/a))
-         hex_n(2) = (lattice_hex_systemSlip(1,i)+2.0_pReal*lattice_hex_systemSlip(2,i))/dsqrt(3.0_pReal)
-         hex_n(3) =  lattice_hex_systemSlip(4,i)/CoverA
-         hex_d(1) =  lattice_hex_systemSlip(5,i)*1.5_pReal ! direction [uvtw]->[3u/2 (u+2v)*sqrt(3)/2 w*(c/a)]
-         hex_d(2) = (lattice_hex_systemSlip(5,i)+2.0_pReal*lattice_hex_systemSlip(6,i))*(0.5_pReal*dsqrt(3.0_pReal))
-         hex_d(3) =  lattice_hex_systemSlip(8,i)*CoverA
+         hex_d(1) =  lattice_hex_systemSlip(1,i)            ! direction [uvtw]->[3u/2 (u+2v)*sqrt(3)/2 w*(c/a)]
+         hex_d(2) = (lattice_hex_systemSlip(1,i)+2.0_pReal*lattice_hex_systemSlip(2,i))/dsqrt(3.0_pReal)
+         hex_d(3) =  lattice_hex_systemSlip(4,i)/CoverA
+         hex_n(1) =  lattice_hex_systemSlip(5,i)*1.5_pReal  ! plane (hkil)->(h (h+2k)/sqrt(3) l/(c/a))
+         hex_n(2) = (lattice_hex_systemSlip(5,i)+2.0_pReal*lattice_hex_systemSlip(6,i))*(0.5_pReal*dsqrt(3.0_pReal))
+         hex_n(3) =  lattice_hex_systemSlip(8,i)*CoverA
 
-         sn(:,i) = hex_n/dsqrt(math_mul3x3(hex_n,hex_n))
          sd(:,i) = hex_d/dsqrt(math_mul3x3(hex_d,hex_d))
-         st(:,i) = math_vectorproduct(sn(:,i),sd(:,i))
+         sn(:,i) = hex_n/dsqrt(math_mul3x3(hex_n,hex_n))
+         st(:,i) = math_vectorproduct(sd(:,i),sn(:,i))
        enddo
        do i = 1,myNtwin
-         hex_n(1) =  lattice_hex_systemTwin(1,i)
-         hex_n(2) = (lattice_hex_systemTwin(1,i)+2.0_pReal*lattice_hex_systemTwin(2,i))/dsqrt(3.0_pReal)
-         hex_n(3) =  lattice_hex_systemTwin(4,i)/CoverA
-         hex_d(1) =  lattice_hex_systemTwin(5,i)*1.5_pReal
-         hex_d(2) = (lattice_hex_systemTwin(5,i)+2.0_pReal*lattice_hex_systemTwin(6,i))*(0.5_pReal*dsqrt(3.0_pReal))
-         hex_d(3) =  lattice_hex_systemTwin(8,i)*CoverA
+         hex_d(1) =  lattice_hex_systemTwin(1,i)
+         hex_d(2) = (lattice_hex_systemTwin(1,i)+2.0_pReal*lattice_hex_systemTwin(2,i))/dsqrt(3.0_pReal)
+         hex_d(3) =  lattice_hex_systemTwin(4,i)/CoverA
+         hex_n(1) =  lattice_hex_systemTwin(5,i)*1.5_pReal
+         hex_n(2) = (lattice_hex_systemTwin(5,i)+2.0_pReal*lattice_hex_systemTwin(6,i))*(0.5_pReal*dsqrt(3.0_pReal))
+         hex_n(3) =  lattice_hex_systemTwin(8,i)*CoverA
 
-         tn(:,i) = hex_n/dsqrt(math_mul3x3(hex_n,hex_n))
          td(:,i) = hex_d/dsqrt(math_mul3x3(hex_d,hex_d))
-         tt(:,i) = math_vectorproduct(tn(:,i),td(:,i))
+         tn(:,i) = hex_n/dsqrt(math_mul3x3(hex_n,hex_n))
+         tt(:,i) = math_vectorproduct(td(:,i),tn(:,i))
          ts(i)   = lattice_hex_shearTwin(i)
        enddo
        interactionSlipSlip => lattice_hex_interactionSlipSlip