Merge branch 'initial-eigenstrain' into 'development'

option to initialize F_i

See merge request damask/DAMASK!455

PRIVATE

@@ -1 +1 @@
-Subproject commit f730bcb8ddd7224011e70c57d0a5c03068532d2d
+Subproject commit ee4b1a03b443a2c497a13c45c9498313442d731e

python/damask/_configmaterial.py

@@ -440,9 +440,10 @@ class ConfigMaterial(Config):
         """
         N,n,shaped = 1,1,{}
 
+        map_dim = {'O':-1,'F_i':-2}
         for k,v in kwargs.items():
             shaped[k] = np.array(v)
-            s = shaped[k].shape[:-1] if k=='O' else shaped[k].shape
+            s = shaped[k].shape[:map_dim.get(k,None)]
             N = max(N,s[0]) if len(s)>0 else N
             n = max(n,s[1]) if len(s)>1 else n
 
@@ -451,11 +452,12 @@ class ConfigMaterial(Config):
         if 'v' not in kwargs:
             shaped['v'] = np.broadcast_to(1/n,(N,n))
 
+        map_shape = {'O':(N,n,4),'F_i':(N,n,3,3)}
         for k,v in shaped.items():
-            target = (N,n,4) if k=='O' else (N,n)
+            target = map_shape.get(k,(N,n))
             obj = np.broadcast_to(v.reshape(util.shapeshifter(v.shape,target,mode='right')),target)
             for i in range(N):
-                if k in ['phase','O','v']:
+                if k in ['phase','O','v','F_i']:
                     for j in range(n):
                         mat[i]['constituents'][j][k] = obj[i,j].item() if isinstance(obj[i,j],np.generic) else obj[i,j]
                 else:

python/tests/test_ConfigMaterial.py

@@ -99,12 +99,15 @@ class TestConfigMaterial:
     @pytest.mark.parametrize('N,n,kw',[
                                         (1,1,{'phase':'Gold',
                                               'O':[1,0,0,0],
+                                              'F_i':np.eye(3),
                                               'homogenization':'SX'}),
                                         (3,1,{'phase':'Gold',
                                               'O':Rotation.from_random(3),
+                                              'F_i':np.broadcast_to(np.eye(3),(3,3,3)),
                                               'homogenization':'SX'}),
                                         (2,3,{'phase':np.broadcast_to(['a','b','c'],(2,3)),
                                               'O':Rotation.from_random((2,3)),
+                                              'F_i':np.broadcast_to(np.eye(3),(2,3,3,3)),
                                               'homogenization':['SX','PX']}),
                                         ])
     def test_material_add(self,kw,N,n):

src/material.f90

@@ -26,6 +26,7 @@ module material
 
 
   type(tRotationContainer), dimension(:), allocatable :: material_O_0
+  type(tTensorContainer),   dimension(:), allocatable :: material_F_i_0
 
   integer, dimension(:), allocatable, public, protected :: &
     homogenization_Nconstituents                                                                    !< number of grains in each homogenization
@@ -48,6 +49,7 @@ module material
   public :: &
     tTensorContainer, &
     tRotationContainer, &
+    material_F_i_0, &
     material_O_0, &
     material_init
 
@@ -159,14 +161,17 @@ subroutine parse()
   enddo
 
   allocate(material_O_0(materials%length))
+  allocate(material_F_i_0(materials%length))
 
   do ma = 1, materials%length
     material     => materials%get(ma)
     constituents => material%get('constituents')
     allocate(material_O_0(ma)%data(constituents%length))
+    allocate(material_F_i_0(ma)%data(1:3,1:3,constituents%length))
     do co = 1, constituents%length
       constituent => constituents%get(co)
       call material_O_0(ma)%data(co)%fromQuaternion(constituent%get_as1dFloat('O',requiredSize=4))
+      material_F_i_0(ma)%data(1:3,1:3,co) = constituent%get_as2dFloat('F_i',defaultVal=math_I3,requiredShape=[3,3])
     enddo
  enddo
 

src/phase_mechanical.f90

@@ -205,6 +205,9 @@ module subroutine mechanical_init(phases)
     phases
 
   integer :: &
+    ce, &
+    co, &
+    ma, &
     ph, &
     en, &
     Nmembers
@@ -261,15 +264,21 @@ module subroutine mechanical_init(phases)
 #endif
   enddo
 
+  do ce = 1, size(material_phaseID,2)
+    ma = discretization_materialAt((ce-1)/discretization_nIPs+1)
+    do co = 1,homogenization_Nconstituents(material_homogenizationID(ce))
+      ph = material_phaseID(co,ce)
+      phase_mechanical_Fi0(ph)%data(1:3,1:3,material_phaseEntry(co,ce)) = material_F_i_0(ma)%data(1:3,1:3,co)
+    enddo
+  enddo
+
   do ph = 1, phases%length
     do en = 1, count(material_phaseID == ph)
 
       phase_mechanical_Fp0(ph)%data(1:3,1:3,en) = phase_O_0(ph)%data(en)%asMatrix()                 ! Fp reflects initial orientation (see 10.1016/j.actamat.2006.01.005)
       phase_mechanical_Fp0(ph)%data(1:3,1:3,en) = phase_mechanical_Fp0(ph)%data(1:3,1:3,en) &
                                                 / math_det33(phase_mechanical_Fp0(ph)%data(1:3,1:3,en))**(1.0_pReal/3.0_pReal)
-      phase_mechanical_Fi0(ph)%data(1:3,1:3,en) = math_I3
-      phase_mechanical_F0(ph)%data(1:3,1:3,en)  = math_I3
-
+      phase_mechanical_F0(ph)%data(1:3,1:3,en) = math_I3
       phase_mechanical_Fe(ph)%data(1:3,1:3,en) = math_inv33(matmul(phase_mechanical_Fi0(ph)%data(1:3,1:3,en), &
                                                                    phase_mechanical_Fp0(ph)%data(1:3,1:3,en)))  ! assuming that euler angles are given in internal strain free configuration
     enddo