From 06e76c1a81a3e46b295e774a0f89b562e4da3c84 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 25 Jun 2021 06:57:08 +0200 Subject: [PATCH] consistent default absolute tolerances better use conservative values, users can easily relax if needed --- src/grid/grid_damage_spectral.f90 | 1 + src/phase_damage_anisobrittle.f90 | 2 +- src/phase_damage_isobrittle.f90 | 2 +- src/phase_damage_isoductile.f90 | 2 +- src/phase_mechanical_plastic_dislotungsten.f90 | 5 +++-- src/phase_mechanical_plastic_dislotwin.f90 | 7 ++++--- src/phase_mechanical_plastic_nonlocal.f90 | 9 +++++---- src/phase_mechanical_plastic_phenopowerlaw.f90 | 2 -- 8 files changed, 16 insertions(+), 14 deletions(-) diff --git a/src/grid/grid_damage_spectral.f90 b/src/grid/grid_damage_spectral.f90 index 967ddb73a..26e0a909e 100644 --- a/src/grid/grid_damage_spectral.f90 +++ b/src/grid/grid_damage_spectral.f90 @@ -261,6 +261,7 @@ subroutine formResidual(in,x_scal,f_scal,dummy,ierr) PetscErrorCode :: ierr integer :: i, j, k, ce + phi_current = x_scal !-------------------------------------------------------------------------------------------------- ! evaluate polarization field diff --git a/src/phase_damage_anisobrittle.f90 b/src/phase_damage_anisobrittle.f90 index 721c27ce2..4f2dbfdf6 100644 --- a/src/phase_damage_anisobrittle.f90 +++ b/src/phase_damage_anisobrittle.f90 @@ -94,7 +94,7 @@ module function anisobrittle_init() result(mySources) Nmembers = count(material_phaseID==p) call phase_allocateState(damageState(p),Nmembers,1,1,0) - damageState(p)%atol = src%get_asFloat('anisobrittle_atol',defaultVal=1.0e-3_pReal) + damageState(p)%atol = src%get_asFloat('anisobrittle_atol',defaultVal=1.0e-9_pReal) if(any(damageState(p)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisobrittle_atol' end associate diff --git a/src/phase_damage_isobrittle.f90 b/src/phase_damage_isobrittle.f90 index 296f7d80e..7229cb22e 100644 --- a/src/phase_damage_isobrittle.f90 +++ b/src/phase_damage_isobrittle.f90 @@ -66,7 +66,7 @@ module function isobrittle_init() result(mySources) Nmembers = count(material_phaseID==ph) call phase_allocateState(damageState(ph),Nmembers,1,1,1) - damageState(ph)%atol = src%get_asFloat('isobrittle_atol',defaultVal=1.0e-3_pReal) + damageState(ph)%atol = src%get_asFloat('isobrittle_atol',defaultVal=1.0e-9_pReal) if(any(damageState(ph)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isobrittle_atol' end associate diff --git a/src/phase_damage_isoductile.f90 b/src/phase_damage_isoductile.f90 index ea6d8acdc..0e5489b8c 100644 --- a/src/phase_damage_isoductile.f90 +++ b/src/phase_damage_isoductile.f90 @@ -70,7 +70,7 @@ module function isoductile_init() result(mySources) Nmembers=count(material_phaseID==ph) call phase_allocateState(damageState(ph),Nmembers,1,1,0) - damageState(ph)%atol = src%get_asFloat('isoductile_atol',defaultVal=1.0e-3_pReal) + damageState(ph)%atol = src%get_asFloat('isoductile_atol',defaultVal=1.0e-9_pReal) if(any(damageState(ph)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isoductile_atol' end associate diff --git a/src/phase_mechanical_plastic_dislotungsten.f90 b/src/phase_mechanical_plastic_dislotungsten.f90 index 7c9accc8d..8d0dbe948 100644 --- a/src/phase_mechanical_plastic_dislotungsten.f90 +++ b/src/phase_mechanical_plastic_dislotungsten.f90 @@ -247,9 +247,10 @@ module function plastic_dislotungsten_init() result(myPlasticity) endIndex = endIndex + prm%sum_N_sl stt%gamma_sl => plasticState(ph)%state(startIndex:endIndex,:) dot%gamma_sl => plasticState(ph)%dotState(startIndex:endIndex,:) - plasticState(ph)%atol(startIndex:endIndex) = 1.0e-2_pReal + plasticState(ph)%atol(startIndex:endIndex) = pl%get_asFloat('atol_gamma',defaultVal=1.0e-6_pReal) + if(any(plasticState(ph)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_gamma' ! global alias - plasticState(ph)%slipRate => plasticState(ph)%dotState(startIndex:endIndex,:) + plasticState(ph)%slipRate => plasticState(ph)%dotState(startIndex:endIndex,:) allocate(dst%Lambda_sl(prm%sum_N_sl,Nmembers), source=0.0_pReal) allocate(dst%threshold_stress(prm%sum_N_sl,Nmembers), source=0.0_pReal) diff --git a/src/phase_mechanical_plastic_dislotwin.f90 b/src/phase_mechanical_plastic_dislotwin.f90 index f1f39c71f..93b1b84a0 100644 --- a/src/phase_mechanical_plastic_dislotwin.f90 +++ b/src/phase_mechanical_plastic_dislotwin.f90 @@ -436,15 +436,16 @@ module function plastic_dislotwin_init() result(myPlasticity) endIndex = endIndex + prm%sum_N_sl stt%gamma_sl=>plasticState(ph)%state(startIndex:endIndex,:) dot%gamma_sl=>plasticState(ph)%dotState(startIndex:endIndex,:) - plasticState(ph)%atol(startIndex:endIndex) = 1.0e-2_pReal + plasticState(ph)%atol(startIndex:endIndex) = pl%get_asFloat('atol_gamma',defaultVal=1.0e-6_pReal) + if(any(plasticState(ph)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_gamma' ! global alias - plasticState(ph)%slipRate => plasticState(ph)%dotState(startIndex:endIndex,:) + plasticState(ph)%slipRate => plasticState(ph)%dotState(startIndex:endIndex,:) startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_tw stt%f_tw=>plasticState(ph)%state(startIndex:endIndex,:) dot%f_tw=>plasticState(ph)%dotState(startIndex:endIndex,:) - plasticState(ph)%atol(startIndex:endIndex) = pl%get_asFloat('atol_f_tw',defaultVal=1.0e-7_pReal) + plasticState(ph)%atol(startIndex:endIndex) = pl%get_asFloat('atol_f_tw',defaultVal=1.0e-6_pReal) if (any(plasticState(ph)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_f_tw' startIndex = endIndex + 1 diff --git a/src/phase_mechanical_plastic_nonlocal.f90 b/src/phase_mechanical_plastic_nonlocal.f90 index 8f38a9301..7c16b1bc7 100644 --- a/src/phase_mechanical_plastic_nonlocal.f90 +++ b/src/phase_mechanical_plastic_nonlocal.f90 @@ -231,14 +231,14 @@ module function plastic_nonlocal_init() result(myPlasticity) mech => phase%get('mechanical') pl => mech%get('plastic') - plasticState(ph)%nonlocal = pl%get_asBool('nonlocal') + plasticState(ph)%nonlocal = pl%get_asBool('flux') #if defined (__GFORTRAN__) prm%output = output_as1dString(pl) #else prm%output = pl%get_as1dString('output',defaultVal=emptyStringArray) #endif - prm%atol_rho = pl%get_asFloat('atol_rho',defaultVal=1.0e4_pReal) + prm%atol_rho = pl%get_asFloat('atol_rho',defaultVal=1.0_pReal) ! This data is read in already in lattice prm%mu = lattice_mu(ph) @@ -488,10 +488,11 @@ module function plastic_nonlocal_init() result(myPlasticity) stt%gamma => plasticState(ph)%state (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nmembers) dot%gamma => plasticState(ph)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nmembers) del%gamma => plasticState(ph)%deltaState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nmembers) - plasticState(ph)%atol(10*prm%sum_N_sl+1:11*prm%sum_N_sl ) = pl%get_asFloat('atol_gamma', defaultVal = 1.0e-2_pReal) + plasticState(ph)%atol(10*prm%sum_N_sl+1:11*prm%sum_N_sl ) = pl%get_asFloat('atol_gamma', defaultVal = 1.0e-6_pReal) if(any(plasticState(ph)%atol(10*prm%sum_N_sl+1:11*prm%sum_N_sl) < 0.0_pReal)) & extmsg = trim(extmsg)//' atol_gamma' - plasticState(ph)%slipRate => plasticState(ph)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nmembers) + ! global alias + plasticState(ph)%slipRate => plasticState(ph)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nmembers) stt%rho_forest => plasticState(ph)%state (11*prm%sum_N_sl + 1:12*prm%sum_N_sl,1:Nmembers) stt%v => plasticState(ph)%state (12*prm%sum_N_sl + 1:16*prm%sum_N_sl,1:Nmembers) diff --git a/src/phase_mechanical_plastic_phenopowerlaw.f90 b/src/phase_mechanical_plastic_phenopowerlaw.f90 index af6734add..5f5279237 100644 --- a/src/phase_mechanical_plastic_phenopowerlaw.f90 +++ b/src/phase_mechanical_plastic_phenopowerlaw.f90 @@ -247,7 +247,6 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) stt%xi_twin = spread(xi_0_tw, 2, Nmembers) dot%xi_twin => plasticState(ph)%dotState(startIndex:endIndex,:) plasticState(ph)%atol(startIndex:endIndex) = pl%get_asFloat('atol_xi',defaultVal=1.0_pReal) - if(any(plasticState(ph)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_xi' startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_sl @@ -263,7 +262,6 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) stt%gamma_twin => plasticState(ph)%state (startIndex:endIndex,:) dot%gamma_twin => plasticState(ph)%dotState(startIndex:endIndex,:) plasticState(ph)%atol(startIndex:endIndex) = pl%get_asFloat('atol_gamma',defaultVal=1.0e-6_pReal) - if(any(plasticState(ph)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_gamma' end associate