Merge branch '2024-10-24_polishing' into 'development'

polishing 2024-10-24

See merge request damask/DAMASK!993

Author: Javier Velo

examples/config/phase/mechanical/elastic/Hooke_X6CrNiMo17-12-2.yaml

@@ -6,13 +6,13 @@ references:
     https://doi.org/10.1007/s11661-006-0129-4,
     fit to Fig. 6 (solid lines, T_min=0°C, T_max=600°C)
 
-C_11: 118.6e+9
-C_11,T: -5.491e+7
+C_11: 229.0e+9
+C_11,T: -1.058e+8
 
 C_12: 154.2e+9
 C_12,T: -6.059e+7
 
-C_44: 229.0e+9
-C_44,T: -1.058e+8
+C_44: 118.6e+9
+C_44,T: -5.491e+7
 
 T_ref: 293.15

examples/config/phase/mechanical/plastic/dislotwin_IF-steel.yaml

@@ -14,7 +14,7 @@ N_sl: [12, 12]
 
 f_edge: [1.0, 1.0]
 b_sl: [2.49e-10, 2.49e-10]
-rho_mob_0: [2.81e+12, 2.8e+12]
+rho_mob_0: [2.81e+12, 2.81e+12]
 rho_dip_0: [1.0, 1.0] # not given
 v_0: [1.4e+3, 1.4e+3]
 Q_sl: [1.57e-19, 1.57e-19] # Delta_F

src/constants.f90

@@ -27,6 +27,9 @@ module constants
     STATUS_FAIL_PHASE_MECHANICAL_STATE, &
     STATUS_FAIL_PHASE_MECHANICAL_DELTASTATE, &
     STATUS_FAIL_PHASE_MECHANICAL_STRESS, &
+    STATUS_FAIL_PHASE_MECHANICAL_STRESS_LP_MAXIT, &
+    STATUS_FAIL_PHASE_MECHANICAL_STRESS_LP_NAN, &
+    STATUS_FAIL_PHASE_MECHANICAL_STRESS_LP_DGESV, &
     STATUS_FAIL_PHASE_DAMAGE, &
     STATUS_FAIL_PHASE_DAMAGE_STATE, &
     STATUS_FAIL_PHASE_DAMAGE_DELTASTATE, &

src/grid/grid_thermal_spectral.f90

@@ -237,7 +237,7 @@ function grid_thermal_spectral_solution(Delta_t) result(solution)
 
   if (solution%converged) &
     print'(/,1x,a)', '... thermal conduction converged ..................................'
-  print'(/,1x,a,f8.4,2x,f8.4,2x,f8.4)', 'Minimum|Maximum|Delta Temperature / K = ', T_min, T_max, stagNorm
+  print'(/,1x,a,f8.3,2x,f8.3,2x,f8.4)', 'Minimum|Maximum|Delta Temperature / K = ', T_min, T_max, stagNorm
   print'(/,1x,a)', '==========================================================================='
   flush(IO_STDOUT)
 

src/phase_mechanical.f90

@@ -449,7 +449,7 @@ function integrateStress(F,Fp0,Fi0,Delta_t,ph,en) result(status)
 
   jacoCounterLi  = 0
   steplengthLi   = 1.0_pREAL
-  residuumLi_old = 0.0_pREAL
+  residuumLi_old = huge(1.0_pREAL)
   Liguess_old    = Liguess
 
   NiterationStressLi = 0
@@ -462,32 +462,35 @@ function integrateStress(F,Fp0,Fi0,Delta_t,ph,en) result(status)
 
     jacoCounterLp  = 0
     steplengthLp   = 1.0_pREAL
-    residuumLp_old = 0.0_pREAL
+    residuumLp_old = huge(1.0_pREAL)
     Lpguess_old    = Lpguess
 
     NiterationStressLp = 0
     LpLoop: do
       NiterationStressLp = NiterationStressLp + 1
-      if (NiterationStressLp>num%nStress_Lp) return ! error
-
+      if (NiterationStressLp>num%nStress_Lp) then
+        status = STATUS_FAIL_PHASE_MECHANICAL_STRESS_LP_MAXIT
+        return
+      endif
       B  = math_I3 - Delta_t*Lpguess
       Fe = matmul(matmul(A,B), invFi_new)
       call phase_hooke_SandItsTangents(S, dS_dFe, dS_dFi, &
-                                        Fe, Fi_new, ph, en)
+                                       Fe, Fi_new, ph, en)
 
       call plastic_LpAndItsTangents(Lp_constitutive, dLp_dS, dLp_dFi, &
-                                         S, Fi_new, ph,en)
+                                    S, Fi_new, ph,en)
 
       !* update current residuum and check for convergence of loop
       atol_Lp = max(num%rtol_Lp * max(norm2(Lpguess),norm2(Lp_constitutive)), &                     ! absolute tolerance from largest acceptable relative error
                     num%atol_Lp)                                                                    ! minimum lower cutoff
       residuumLp = Lpguess - Lp_constitutive
 
       if (any(IEEE_is_NaN(residuumLp))) then
-        return ! error
+        status = STATUS_FAIL_PHASE_MECHANICAL_STRESS_LP_NAN
+        return
       elseif (norm2(residuumLp) < atol_Lp) then                                                     ! converged if below absolute tolerance
         exit LpLoop
-      elseif (NiterationStressLp == 1 .or. norm2(residuumLp) < norm2(residuumLp_old)) then          ! not converged, but improved norm of residuum (always proceed in first iteration)...
+      elseif (norm2(residuumLp) < norm2(residuumLp_old)) then                                       ! not converged, but improved norm of residuum...
         residuumLp_old = residuumLp                                                                 ! ...remember old values and...
         Lpguess_old    = Lpguess
         steplengthLp   = 1.0_pREAL                                                                  ! ...proceed with normal step length (calculate new search direction)
@@ -508,7 +511,10 @@ function integrateStress(F,Fp0,Fi0,Delta_t,ph,en) result(status)
                  - math_3333to99(math_mul3333xx3333(math_mul3333xx3333(dLp_dS,dS_dFe),dFe_dLp))
         temp_9 = math_33to9(residuumLp)
         call dgesv(9,1,dRLp_dLp,9,devNull_9,temp_9,9,ierr)                                          ! solve dRLp/dLp * delta Lp = -res for delta Lp
-        if (ierr /= 0) return ! error
+        if (ierr /= 0) then
+          status = STATUS_FAIL_PHASE_MECHANICAL_STRESS_LP_DGESV
+          return
+        end if
         deltaLp = - math_9to33(temp_9)
       end if calculateJacobiLp
 
@@ -527,7 +533,7 @@ function integrateStress(F,Fp0,Fi0,Delta_t,ph,en) result(status)
       return ! error
     elseif (norm2(residuumLi) < atol_Li) then                                                       ! converged if below absolute tolerance
       exit LiLoop
-    elseif (NiterationStressLi == 1 .or. norm2(residuumLi) < norm2(residuumLi_old)) then            ! not converged, but improved norm of residuum (always proceed in first iteration)...
+    elseif (norm2(residuumLi) < norm2(residuumLi_old)) then                                         ! not converged, but improved norm of residuum ...
       residuumLi_old = residuumLi                                                                   ! ...remember old values and...
       Liguess_old    = Liguess
       steplengthLi   = 1.0_pREAL                                                                    ! ...proceed with normal step length (calculate new search direction)