Separating functionality

lattice should become a module with static data and functions

PRIVATE

@@ -1 +1 @@
-Subproject commit 185cb53be76eded17565c5fa91bd9b4499cda4b8
+Subproject commit 5ecf7d74a79da155aed31e776963c7bc10eb9890

src/lattice.f90

@@ -365,25 +365,11 @@ module lattice
        1, 1, 1,      1,-2, 1  &
        ],pReal),shape(BCT_SYSTEMSLIP))                                                              !< bct slip systems for c/a = 0.5456 (Sn), sorted by Bieler 2009 (https://doi.org/10.1007/s11664-009-0909-x)
 
-
-  enum, bind(c); enumerator :: &
-    lattice_UNDEFINED_ID, &
-    lattice_FCC_ID, &
-    lattice_BCC_ID, &
-    lattice_HEX_ID, &
-    lattice_BCT_ID
-  end enum
-
 ! SHOULD NOT BE PART OF LATTICE BEGIN
   real(pReal),                        dimension(:),     allocatable, public, protected :: &
-    lattice_mu, lattice_nu, &
+    lattice_mu, lattice_nu
-    lattice_mu_phi, &
-    lattice_rho, &
-    lattice_c_p
    real(pReal),                       dimension(:,:,:), allocatable, public, protected :: &
     lattice_C66
- integer(kind(lattice_UNDEFINED_ID)), dimension(:),     allocatable, public, protected :: &
-    lattice_structure
 ! SHOULD NOT BE PART OF LATTICE END
 
   interface lattice_forestProjection_edge
@@ -396,10 +382,6 @@ module lattice
 
   public :: &
     lattice_init, &
-    lattice_FCC_ID, &
-    lattice_BCC_ID, &
-    lattice_HEX_ID, &
-    lattice_BCT_ID, &
     lattice_equivalent_nu, &
     lattice_equivalent_mu, &
     lattice_applyLatticeSymmetry33, &
@@ -446,11 +428,9 @@ subroutine lattice_init
   phases => config_material%get('phase')
   Nphases = phases%length
 
-  allocate(lattice_structure(Nphases),source = lattice_UNDEFINED_ID)
   allocate(lattice_C66(6,6,Nphases), source=0.0_pReal)
 
-  allocate(lattice_rho, &
+  allocate(lattice_mu, lattice_nu,&
-           lattice_mu, lattice_nu,&
            source=[(0.0_pReal,i=1,Nphases)])
 
   do ph = 1, phases%length
@@ -466,19 +446,6 @@ subroutine lattice_init
     lattice_C66(3,3,ph) = elasticity%get_asFloat('C_33',defaultVal=0.0_pReal)
     lattice_C66(6,6,ph) = elasticity%get_asFloat('C_66',defaultVal=0.0_pReal)
 
-    select case(phase%get_asString('lattice'))
-      case('cF')
-        lattice_structure(ph) = lattice_FCC_ID
-      case('cI')
-        lattice_structure(ph) = lattice_BCC_ID
-      case('hP')
-        lattice_structure(ph) = lattice_HEX_ID
-      case('tI')
-        lattice_structure(ph) = lattice_BCT_ID
-      case default
-        call IO_error(130,ext_msg='lattice_init: '//phase%get_asString('lattice'))
-    end select
-
     lattice_C66(1:6,1:6,ph) = applyLatticeSymmetryC66(lattice_C66(1:6,1:6,ph),phase%get_asString('lattice'))
 
     lattice_nu(ph) = lattice_equivalent_nu(lattice_C66(1:6,1:6,ph),'voigt')
@@ -489,8 +456,6 @@ subroutine lattice_init
       if (abs(lattice_C66(i,i,ph))<tol_math_check) &
         call IO_error(135,el=i,ip=ph,ext_msg='matrix diagonal "el"ement of phase "ip"')
     enddo
-
-    lattice_rho(ph) = phase%get_asFloat('rho', defaultVal=0.0_pReal)
 ! SHOULD NOT BE PART OF LATTICE END
 
     call selfTest

src/phase.f90

@@ -19,6 +19,11 @@ module phase
   implicit none
   private
 
+
+  character(len=2), allocatable, dimension(:) :: phase_lattice
+  real(pReal),      allocatable, dimension(:) :: phase_cOverA
+  real(pReal),      allocatable, dimension(:) :: phase_rho
+
   type(tRotationContainer), dimension(:), allocatable :: &
     phase_orientation0, &
     phase_orientation
@@ -348,7 +353,8 @@ subroutine phase_init
   class (tNode), pointer :: &
     debug_constitutive, &
     materials, &
-    phases
+    phases, &
+    phase
 
 
   print'(/,a)', ' <<<+-  phase init  -+>>>'; flush(IO_STDOUT)
@@ -365,9 +371,19 @@ subroutine phase_init
   materials => config_material%get('material')
   phases    => config_material%get('phase')
 
-
+  allocate(phase_lattice(phases%length))
+  allocate(phase_cOverA(phases%length),source=-1.0_pReal)
+  allocate(phase_rho(phases%length))
   allocate(phase_orientation0(phases%length))
+
   do ph = 1,phases%length
+    phase => phases%get(ph)
+    phase_lattice(ph) = phase%get_asString('lattice')
+    if (all(phase_lattice(ph) /= ['cF','cI','hP','tI'])) &
+      call IO_error(130,ext_msg='phase_init: '//phase%get_asString('lattice'))
+    if (any(phase_lattice(ph) == ['hP','tI'])) &
+      phase_cOverA(ph) = phase%get_asFloat('c/a')
+    phase_rho(ph) = phase%get_asFloat('rho',defaultVal=0.0_pReal)
     allocate(phase_orientation0(ph)%data(count(material_phaseID==ph)))
   enddo
 

src/phase_mechanical.f90

@@ -896,7 +896,6 @@ subroutine crystallite_results(group,ph)
 
   integer :: ou
   real(pReal), allocatable, dimension(:,:)   :: selected_rotations
-  character(len=:), allocatable              :: structureLabel
 
 
   call results_closeGroup(results_addGroup(group//'/mechanical'))
@@ -929,20 +928,12 @@ subroutine crystallite_results(group,ph)
         call results_writeDataset(group//'/mechanical/',phase_mechanical_S(ph)%data,'S', &
                                  'second Piola-Kirchhoff stress','Pa')
       case('O')
-          select case(lattice_structure(ph))
-            case(lattice_FCC_ID)
-              structureLabel = 'cF'
-            case(lattice_BCC_ID)
-              structureLabel = 'cI'
-            case(lattice_HEX_ID)
-              structureLabel = 'hP'
-            case(lattice_BCT_ID)
-              structureLabel = 'tI'
-          end select
         selected_rotations = select_rotations(phase_orientation(ph)%data)
         call results_writeDataset(group//'/mechanical',selected_rotations,output_constituent(ph)%label(ou),&
                                  'crystal orientation as quaternion','q_0 (q_1 q_2 q_3)')
-          call results_addAttribute('lattice',structureLabel,group//'/mechanical/'//output_constituent(ph)%label(ou))
+        call results_addAttribute('lattice',phase_lattice(ph),group//'/mechanical/'//output_constituent(ph)%label(ou))
+        if (any(phase_lattice(ph) == ['hP', 'tI'])) &
+          call results_addAttribute('c/a',phase_cOverA(ph),group//'/mechanical/'//output_constituent(ph)%label(ou))
     end select
   enddo
 

src/phase_mechanical_plastic_dislotwin.f90

@@ -202,7 +202,7 @@ module function plastic_dislotwin_init() result(myPlasticity)
 
       prm%n0_sl            = lattice_slip_normal(N_sl,phase%get_asString('lattice'),&
                                                  phase%get_asFloat('c/a',defaultVal=0.0_pReal))
-      prm%fccTwinTransNucleation = lattice_structure(ph) == lattice_FCC_ID .and. (N_sl(1) == 12)
+      prm%fccTwinTransNucleation = phase_lattice(ph) == 'cF' .and. (N_sl(1) == 12)
       if(prm%fccTwinTransNucleation) prm%fcc_twinNucleationSlipPair = lattice_FCC_TWINNUCLEATIONSLIPPAIR
 
       rho_mob_0                = pl%get_as1dFloat('rho_mob_0',   requiredSize=size(N_sl))
@@ -231,7 +231,7 @@ module function plastic_dislotwin_init() result(myPlasticity)
       ! multiplication factor according to crystal structure (nearest neighbors bcc vs fcc/hex)
       ! details: Argon & Moffat, Acta Metallurgica, Vol. 29, pg 293 to 299, 1981
       prm%omega = pl%get_asFloat('omega',  defaultVal = 1000.0_pReal) &
-                * merge(12.0_pReal,8.0_pReal,any(lattice_structure(ph) == [lattice_FCC_ID,lattice_HEX_ID]))
+                * merge(12.0_pReal,8.0_pReal,any(phase_lattice(ph) == ['cF','hP']))
 
       ! expand: family => system
       rho_mob_0        = math_expand(rho_mob_0,       N_sl)
@@ -340,7 +340,7 @@ module function plastic_dislotwin_init() result(myPlasticity)
                                                pl%get_asFloat('a_cI', defaultVal=0.0_pReal), &
                                                pl%get_asFloat('a_cF', defaultVal=0.0_pReal))
 
-      if (lattice_structure(ph) /= lattice_FCC_ID) then
+      if (phase_lattice(ph) /= 'cF') then
         prm%dot_N_0_tr = pl%get_as1dFloat('dot_N_0_tr')
         prm%dot_N_0_tr = math_expand(prm%dot_N_0_tr,N_tr)
       endif
@@ -355,7 +355,7 @@ module function plastic_dislotwin_init() result(myPlasticity)
       if (    prm%i_tr          < 0.0_pReal)  extmsg = trim(extmsg)//' i_tr'
       if (any(prm%t_tr          < 0.0_pReal)) extmsg = trim(extmsg)//' t_tr'
       if (any(prm%s             < 0.0_pReal)) extmsg = trim(extmsg)//' p_tr'
-      if (lattice_structure(ph) /= lattice_FCC_ID) then
+      if (phase_lattice(ph) /= 'cF') then
         if (any(prm%dot_N_0_tr  < 0.0_pReal)) extmsg = trim(extmsg)//' dot_N_0_tr'
       endif
     else transActive

src/phase_mechanical_plastic_nonlocal.f90

@@ -622,7 +622,7 @@ module subroutine nonlocal_dependentState(ph, en, ip, el)
 
   ! coefficients are corrected for the line tension effect
   ! (see Kubin,Devincre,Hoc; 2008; Modeling dislocation storage rates and mean free paths in face-centered cubic crystals)
-  if (any(lattice_structure(ph) == [LATTICE_bcc_ID,LATTICE_fcc_ID])) then
+  if (any(phase_lattice(ph) == ['cI','cF'])) then
     myInteractionMatrix = prm%h_sl_sl &
                         * spread((  1.0_pReal - prm%f_F &
                                    + prm%f_F &
@@ -1055,7 +1055,7 @@ module subroutine nonlocal_dotState(Mp, Temperature,timestep, &
   !****************************************************************************
   !*** dislocation multiplication
   rhoDotMultiplication = 0.0_pReal
-  isBCC: if (lattice_structure(ph) == LATTICE_bcc_ID) then
+  isBCC: if (phase_lattice(ph) == 'cI') then
     forall (s = 1:ns, sum(abs(v(s,1:4))) > 0.0_pReal)
       rhoDotMultiplication(s,1:2) = sum(abs(gdot(s,3:4))) / prm%b_sl(s) &                           ! assuming double-cross-slip of screws to be decisive for multiplication
                                   * sqrt(stt%rho_forest(s,en)) / prm%i_sl(s) ! &                    ! mean free path
@@ -1111,7 +1111,7 @@ module subroutine nonlocal_dotState(Mp, Temperature,timestep, &
           + rhoDip(:,c) * (abs(gdot(:,2*c-1)) + abs(gdot(:,2*c))))                                  ! single knocks dipole constituent
 
   ! annihilated screw dipoles leave edge jogs behind on the colinear system
-  if (lattice_structure(ph) == LATTICE_fcc_ID) &
+  if (phase_lattice(ph) == 'cF') &
     forall (s = 1:ns, prm%colinearSystem(s) > 0) &
       rhoDotAthermalAnnihilation(prm%colinearSystem(s),1:2) = - rhoDotAthermalAnnihilation(s,10) &
         * 0.25_pReal * sqrt(stt%rho_forest(s,en)) * (dUpper(s,2) + dLower(s,2)) * prm%f_ed

src/phase_thermal.f90

@@ -196,7 +196,7 @@ module function phase_mu_T(co,ce) result(mu)
   real(pReal) :: mu
 
 
-  mu = lattice_rho(material_phaseID(co,ce)) &
+  mu = phase_rho(material_phaseID(co,ce)) &
      * param(material_phaseID(co,ce))%C_p
 
 end function phase_mu_T