Skip to content

Commit

Permalink
Merge branch 'restructure-numerics' into 'development'
Browse files Browse the repository at this point in the history 
Restructure numerical parameters RGC

See merge request damask/DAMASK!788
  • Loading branch information
@eisenlohr
eisenlohr committed Aug 1, 2023
2 parents 4d03d1a + 485839b commit d73f769
Show file tree
Hide file tree
Showing 3 changed files with 44 additions and 69 deletions.
36 changes: 15 additions & 21 deletions examples/config/numerics.yaml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,28 +2,22 @@
# Case sensitive keys

homogenization:
mech:
mechanical:
RGC:
atol: 1.0e+4 # absolute tolerance of RGC residuum (in Pa)
rtol: 1.0e-3 # relative ...
amax: 1.0e+10 # absolute upper-limit of RGC residuum (in Pa)
rmax: 1.0e+2 # relative ...
perturbpenalty: 1.0e-7 # perturbation for computing penalty tangent
relevantmismatch: 1.0e-5 # minimum threshold of mismatch
viscositypower: 1.0e+0 # power (sensitivity rate) of numerical viscosity in RGC scheme
viscositymodulus: 0.0e+0 # stress modulus of RGC numerical viscosity (zero = without numerical viscosity)
# suggestion: larger than the aTol_RGC but still far below the expected flow stress of material
refrelaxationrate: 1.0e-3 # reference rate of relaxation (about the same magnitude as straining rate, possibly a bit higher)
maxrelaxationrate: 1.0e+0 # threshold of maximum relaxation vector increment (if exceed this then cutback)
maxvoldiscrepancy: 1.0e-5 # maximum allowable relative volume discrepancy
voldiscrepancymod: 1.0e+12
discrepancypower: 5.0

generic:
subStepMin: 1.0e-3 # minimum (relative) size of sub-step allowed during cutback in homogenization
subStepSize: 0.25 # size of substep when cutback introduced in homogenization (value between 0 and 1)
stepIncrease: 1.5 # increase of next substep size when previous substep converged in homogenization (value higher than 1)
nMPstate: 10 # materialpoint state loop limit
eps_abs_P: 1.0e+4 # absolute tolerance of RGC residuum (in Pa)
eps_rel_P: 1.0e-3 # relative ...
eps_abs_max: 1.0e+10 # absolute upper-limit of RGC residuum (in Pa)
eps_rel_max: 1.0e+2 # relative ...
Delta_a: 1.0e-7 # perturbation for computing penalty tangent
relevant_mismatch: 1.0e-5 # minimum threshold of mismatch
viscosity_exponent: 1.0e+0 # power (sensitivity rate) of numerical viscosity in RGC scheme
viscosity_modulus: 0.0e+0 # stress modulus of RGC numerical viscosity (zero = without numerical viscosity)
# suggestion: larger than the aTol_RGC but still far below the expected flow stress of material
dot_a_ref: 1.0e-3 # reference rate of relaxation (about the same magnitude as straining rate, possibly a bit higher)
dot_a_max: 1.0e+0 # threshold of maximum relaxation vector increment (if exceed this then cutback)
Delta_V_max: 1.0e-5 # maximum allowable relative volume discrepancy
Delta_V_modulus: 1.0e+12
Delta_V_exponent: 5.0

solver:
grid:
Expand Down
23 changes: 2 additions & 21 deletions src/homogenization.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -59,15 +59,6 @@ module homogenization
real(pREAL), dimension(:,:,:,:,:), allocatable, public :: & !, protected :: &
homogenization_dPdF !< tangent of first P--K stress at IP


!--------------------------------------------------------------------------------------------------
type :: tNumerics
integer :: &
nMPstate !< materialpoint state loop limit
end type tNumerics

type(tNumerics) :: num

!--------------------------------------------------------------------------------------------------
interface

Expand Down Expand Up @@ -214,12 +205,6 @@ subroutine homogenization_init()
allocate(damageState_h (size(material_name_homogenization)))
call parseHomogenization()

num_homog => config_numerics%get_dict('homogenization',defaultVal=emptyDict)
num_homogGeneric => num_homog%get_dict('generic',defaultVal=emptyDict)

num%nMPstate = num_homogGeneric%get_asInt('nMPstate',defaultVal=10)
if (num%nMPstate < 1) call IO_error(301,ext_msg='nMPstate')

call mechanical_init()
call thermal_init()
call damage_init()
Expand All @@ -236,15 +221,14 @@ subroutine homogenization_mechanical_response(Delta_t,cell_start,cell_end)
integer, intent(in) :: &
cell_start, cell_end
integer :: &
NiterationMPstate, &
co, ce, ho, en
logical :: &
converged
logical, dimension(2) :: &
doneAndHappy


!$OMP PARALLEL DO PRIVATE(en,ho,co,NiterationMPstate,converged,doneAndHappy)
!$OMP PARALLEL DO PRIVATE(en,ho,co,converged,doneAndHappy)
do ce = cell_start, cell_end

en = material_entry_homogenization(ce)
Expand All @@ -258,10 +242,7 @@ subroutine homogenization_mechanical_response(Delta_t,cell_start,cell_end)

doneAndHappy = [.false.,.true.]

NiterationMPstate = 0
convergenceLooping: do while (.not. (terminallyIll .or. doneAndHappy(1)) &
.and. NiterationMPstate < num%nMPstate)
NiterationMPstate = NiterationMPstate + 1
convergenceLooping: do while (.not. (terminallyIll .or. doneAndHappy(1)))

call mechanical_partition(homogenization_F(1:3,1:3,ce),ce)
converged = all([(phase_mechanical_constitutive(Delta_t,co,ce),co=1,homogenization_Nconstituents(ho))])
Expand Down
54 changes: 27 additions & 27 deletions src/homogenization_mechanical_RGC.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -108,33 +108,33 @@ module subroutine RGC_init()
num_mechanical => num_homogenization%get_dict('mechanical',defaultVal=emptyDict)
num_RGC => num_mechanical%get_dict('RGC',defaultVal=emptyDict)

num%atol = num_RGC%get_asReal('atol', defaultVal=1.0e+4_pREAL)
num%rtol = num_RGC%get_asReal('rtol', defaultVal=1.0e-3_pREAL)
num%absMax = num_RGC%get_asReal('amax', defaultVal=1.0e+10_pREAL)
num%relMax = num_RGC%get_asReal('rmax', defaultVal=1.0e+2_pREAL)
num%pPert = num_RGC%get_asReal('perturbpenalty', defaultVal=1.0e-7_pREAL)
num%xSmoo = num_RGC%get_asReal('relvantmismatch', defaultVal=1.0e-5_pREAL)
num%viscPower = num_RGC%get_asReal('viscositypower', defaultVal=1.0e+0_pREAL)
num%viscModus = num_RGC%get_asReal('viscositymodulus', defaultVal=0.0e+0_pREAL)
num%refRelaxRate = num_RGC%get_asReal('refrelaxationrate', defaultVal=1.0e-3_pREAL)
num%maxdRelax = num_RGC%get_asReal('maxrelaxationrate', defaultVal=1.0e+0_pREAL)
num%maxVolDiscr = num_RGC%get_asReal('maxvoldiscrepancy', defaultVal=1.0e-5_pREAL)
num%volDiscrMod = num_RGC%get_asReal('voldiscrepancymod', defaultVal=1.0e+12_pREAL)
num%volDiscrPow = num_RGC%get_asReal('dicrepancypower', defaultVal=5.0_pREAL)

if (num%atol <= 0.0_pREAL) call IO_error(301,ext_msg='absTol_RGC')
if (num%rtol <= 0.0_pREAL) call IO_error(301,ext_msg='relTol_RGC')
if (num%absMax <= 0.0_pREAL) call IO_error(301,ext_msg='absMax_RGC')
if (num%relMax <= 0.0_pREAL) call IO_error(301,ext_msg='relMax_RGC')
if (num%pPert <= 0.0_pREAL) call IO_error(301,ext_msg='pPert_RGC')
if (num%xSmoo <= 0.0_pREAL) call IO_error(301,ext_msg='xSmoo_RGC')
if (num%viscPower < 0.0_pREAL) call IO_error(301,ext_msg='viscPower_RGC')
if (num%viscModus < 0.0_pREAL) call IO_error(301,ext_msg='viscModus_RGC')
if (num%refRelaxRate <= 0.0_pREAL) call IO_error(301,ext_msg='refRelaxRate_RGC')
if (num%maxdRelax <= 0.0_pREAL) call IO_error(301,ext_msg='maxdRelax_RGC')
if (num%maxVolDiscr <= 0.0_pREAL) call IO_error(301,ext_msg='maxVolDiscr_RGC')
if (num%volDiscrMod < 0.0_pREAL) call IO_error(301,ext_msg='volDiscrMod_RGC')
if (num%volDiscrPow <= 0.0_pREAL) call IO_error(301,ext_msg='volDiscrPw_RGC')
num%atol = num_RGC%get_asReal('eps_abs_P', defaultVal=1.0e+4_pREAL)
num%rtol = num_RGC%get_asReal('eps_rel_P', defaultVal=1.0e-3_pREAL)
num%absMax = num_RGC%get_asReal('eps_abs_max', defaultVal=1.0e+10_pREAL)
num%relMax = num_RGC%get_asReal('eps_rel_max', defaultVal=1.0e+2_pREAL)
num%pPert = num_RGC%get_asReal('Delta_a', defaultVal=1.0e-7_pREAL)
num%xSmoo = num_RGC%get_asReal('relevant_mismatch', defaultVal=1.0e-5_pREAL)
num%viscPower = num_RGC%get_asReal('viscosity_exponent', defaultVal=1.0e+0_pREAL)
num%viscModus = num_RGC%get_asReal('viscosity_modulus', defaultVal=0.0e+0_pREAL)
num%refRelaxRate = num_RGC%get_asReal('dot_a_ref', defaultVal=1.0e-3_pREAL)
num%maxdRelax = num_RGC%get_asReal('dot_a_max', defaultVal=1.0e+0_pREAL)
num%maxVolDiscr = num_RGC%get_asReal('Delta_V_max', defaultVal=1.0e-5_pREAL)
num%volDiscrMod = num_RGC%get_asReal('Delta_V_modulus', defaultVal=1.0e+12_pREAL)
num%volDiscrPow = num_RGC%get_asReal('Delta_V_exponent', defaultVal=5.0_pREAL)

if (num%atol <= 0.0_pREAL) call IO_error(301,ext_msg='eps_abs_P')
if (num%rtol <= 0.0_pREAL) call IO_error(301,ext_msg='eps_rel_P')
if (num%absMax <= 0.0_pREAL) call IO_error(301,ext_msg='eps_abs_max')
if (num%relMax <= 0.0_pREAL) call IO_error(301,ext_msg='eps_rel_max')
if (num%pPert <= 0.0_pREAL) call IO_error(301,ext_msg='Delta_a')
if (num%xSmoo <= 0.0_pREAL) call IO_error(301,ext_msg='relevant_mismatch')
if (num%viscPower < 0.0_pREAL) call IO_error(301,ext_msg='viscosity_exponent')
if (num%viscModus < 0.0_pREAL) call IO_error(301,ext_msg='viscosity_modulus')
if (num%refRelaxRate <= 0.0_pREAL) call IO_error(301,ext_msg='dot_a_ref')
if (num%maxdRelax <= 0.0_pREAL) call IO_error(301,ext_msg='dot_a_max')
if (num%maxVolDiscr <= 0.0_pREAL) call IO_error(301,ext_msg='Delta_V_max')
if (num%volDiscrMod < 0.0_pREAL) call IO_error(301,ext_msg='Delta_V_modulus')
if (num%volDiscrPow <= 0.0_pREAL) call IO_error(301,ext_msg='Delta_V_exponent')


do ho = 1, size(mechanical_type)
Expand Down

0 comments on commit d73f769

Please sign in to comment.