Merge branch '382-missing-damage-type-error' into 'development'

better consistency checks for thermal and damage

Closes #382

See merge request damask/DAMASK!1017

Author: hyharry

src/IO.f90

@@ -373,9 +373,7 @@ subroutine IO_error(error_ID,ext_msg,label1,ID1,label2,ID2)
 !--------------------------------------------------------------------------------------------------
 ! plasticity error messages
     case (200)
-      msg = 'unknown elasticity specified:'
-    case (201)
-      msg = 'unknown plasticity specified:'
+      msg = 'unknown type specified:'
 
     case (211)
       msg = 'material parameter out of bounds:'

src/phase_damage.f90

@@ -72,8 +72,10 @@ module subroutine damage_init()
   type(tDict), pointer :: &
     phases, &
     phase, &
-    source
-  character(len=:), allocatable :: refs
+    damage
+  character(len=:), allocatable :: &
+    refs, &
+    extmsg
   logical:: damage_active
 
 
@@ -85,40 +87,59 @@ module subroutine damage_init()
   allocate(damageState(phases%length))
   allocate(param(phases%length))
 
+  extmsg = ''
   damage_active = .false.
   do ph = 1,phases%length
 
     Nmembers = count(material_ID_phase == ph)
-
     allocate(current(ph)%phi(Nmembers),source=1.0_pREAL)
 
     phase => phases%get_dict(ph)
-    source => phase%get_dict('damage',defaultVal=emptyDict)
-    if (source%length > 0) then
+    if (damage_active) then
+      damage => phase%get_dict('damage')
+    else
+      damage => phase%get_dict('damage',defaultVal=emptyDict)
+    end if
+
+    if (damage%length > 0) then
+      damage_active = .true.
+
       print'(/,1x,a,i0,a)', 'phase ',ph,': '//phases%key(ph)
-      refs = config_listReferences(source,indent=3)
+      refs = config_listReferences(damage,indent=3)
       if (len(refs) > 0) print'(/,1x,a)', refs
-      damage_active = .true.
-      param(ph)%mu = source%get_asReal('mu')
-      param(ph)%l_c = source%get_asReal('l_c')
+
+      param(ph)%mu = damage%get_asReal('mu')
+      param(ph)%l_c = damage%get_asReal('l_c')
+
+      ! sanity checks
+      if (param(ph)%mu  < 0.0_pREAL)  extmsg = trim(extmsg)//' mu'
+      if (param(ph)%l_c < 0.0_pREAL)  extmsg = trim(extmsg)//' l_c'
+      if (extmsg /= '') call IO_error(211,ext_msg=trim(extmsg))
+
     end if
 
   end do
 
   allocate(damage_type(phases%length), source = UNDEFINED)
 
+  where(isobrittle_init()  ) damage_type = DAMAGE_ISOBRITTLE
+  where(anisobrittle_init()) damage_type = DAMAGE_ANISOBRITTLE
+  phase_damage_maxSizeDotState = maxval(damageState%sizeDotState)
+
   if (damage_active) then
-    where(isobrittle_init()  ) damage_type = DAMAGE_ISOBRITTLE
-    where(anisobrittle_init()) damage_type = DAMAGE_ANISOBRITTLE
+    do ph = 1,phases%length
+      phase => phases%get_dict(ph)
+      damage => phase%get_dict('damage')
+      if (any(damage%keys() == 'type') .and. damage_type(ph) == UNDEFINED) &
+        call IO_error(200,label1='damage',ext_msg=damage%get_asStr('type'))
+    end do
   end if
 
-  phase_damage_maxSizeDotState = maxval(damageState%sizeDotState)
-
 end subroutine damage_init
 
 
 !--------------------------------------------------------------------------------------------------
-!> @brief calculate stress (P)
+!> @brief t.b.d.
 !--------------------------------------------------------------------------------------------------
 module function phase_damage_constitutive(Delta_t,co,ce) result(status)
 
@@ -141,7 +162,7 @@ end function phase_damage_constitutive
 
 
 !--------------------------------------------------------------------------------------------------
-!> @brief returns the degraded/modified elasticity matrix
+!> @brief Return the degraded/modified elasticity matrix.
 !--------------------------------------------------------------------------------------------------
 module function phase_damage_C66(C66,ph,en) result(C66_degraded)
 

src/phase_mechanical_elastic.f90

@@ -46,7 +46,8 @@ module subroutine elastic_init(phases)
     print'(/,1x,a,1x,i0,a)', 'phase',ph,': '//phases%key(ph)
     refs = config_listReferences(elastic,indent=3)
     if (len(refs) > 0) print'(/,1x,a)', refs
-    if (elastic%get_asStr('type') /= 'Hooke') call IO_error(200,ext_msg=elastic%get_asStr('type'))
+    if (elastic%get_asStr('type') /= 'Hooke') &
+      call IO_error(200,label1='elasticity',ext_msg=elastic%get_asStr('type'))
 
     associate(prm => param(ph))
 

src/phase_mechanical_plastic.f90

@@ -226,7 +226,7 @@ module subroutine plastic_init
   where(plastic_dislotungsten_init())     mechanical_plasticity_type = MECHANICAL_PLASTICITY_DISLOTUNGSTEN
   where(plastic_nonlocal_init())          mechanical_plasticity_type = MECHANICAL_PLASTICITY_NONLOCAL
 
-  if (any(mechanical_plasticity_type == UNDEFINED)) call IO_error(201)
+  if (any(mechanical_plasticity_type == UNDEFINED)) call IO_error(200, label1='plasticity')
 
 end subroutine plastic_init
 

src/phase_thermal.f90

@@ -82,28 +82,37 @@ module subroutine thermal_init(phases)
   integer :: &
     ph, so, &
     Nmembers
-
+  logical :: thermal_active
 
   print'(/,1x,a)', '<<<+-  phase:thermal init  -+>>>'
 
   allocate(current(phases%length))
   allocate(thermalState(phases%length))
   allocate(thermal_Nsources(phases%length),source = 0)
   allocate(param(phases%length))
-  extmsg = ''
 
+  extmsg = ''
+  thermal_active = .false.
   do ph = 1, phases%length
     Nmembers = count(material_ID_phase == ph)
     allocate(current(ph)%T(Nmembers),source=T_ROOM)
     allocate(current(ph)%dot_T(Nmembers),source=0.0_pREAL)
+
     phase => phases%get_dict(ph)
-    thermal => phase%get_dict('thermal',defaultVal=emptyDict)
+    if (thermal_active) then
+      thermal => phase%get_dict('thermal')
+    else
+      thermal => phase%get_dict('thermal',defaultVal=emptyDict)
+    end if
 
     ! ToDo: temperature dependency of K and C_p
     if (thermal%length > 0) then
+      thermal_active = .true.
+
       print'(/,1x,a,i0,a)', 'phase ',ph,': '//phases%key(ph)
       refs = config_listReferences(thermal,indent=3)
       if (len(refs) > 0) print'(/,1x,a)', refs
+
       param(ph)%C_p = thermal%get_asReal('C_p')
       param(ph)%K(1,1) = thermal%get_asReal('K_11')
       if (any(phase_lattice(ph) == ['hP','tI'])) param(ph)%K(3,3) = thermal%get_asReal('K_33')
@@ -113,9 +122,6 @@ module subroutine thermal_init(phases)
       if (    param(ph)%C_p <= 0.0_pREAL )  extmsg = trim(extmsg)//' C_p'
       if (any(param(ph)%K   <  0.0_pREAL))  extmsg = trim(extmsg)//' K'
       if (    phase_rho(ph) <= 0.0_pREAL )  extmsg = trim(extmsg)//' rho'
-
-!--------------------------------------------------------------------------------------------------
-!  exit if any parameter is out of range
       if (extmsg /= '') call IO_error(211,ext_msg=trim(extmsg))
 
 #if defined(__GFORTRAN__)