diff --git a/code/IO.f90 b/code/IO.f90 index 0da6d9cb9..f18c42c29 100644 --- a/code/IO.f90 +++ b/code/IO.f90 @@ -134,7 +134,6 @@ function IO_hybridIA(Nast,ODFfileName) use prec, only: pReal, pInt - use math, only: inRad implicit none @@ -150,6 +149,9 @@ real(pReal), dimension(:,:,:), allocatable :: dV_V real(pReal), dimension(3,Nast) :: IO_hybridIA + real(pReal), parameter :: pi = 3.14159265358979323846264338327950288419716939937510_pReal + real(pReal), parameter :: inRad = pi/180.0_pReal + if (.not. IO_open_file(999,ODFfileName)) goto 100 !--- parse header of ODF file --- diff --git a/code/math.f90 b/code/math.f90 index a25f0de6d..0dfae69c7 100644 --- a/code/math.f90 +++ b/code/math.f90 @@ -213,13 +213,22 @@ real(pReal), dimension(132,3), parameter :: sym = & SUBROUTINE math_init () use prec, only: pReal,pInt + use numerics, only: fixedSeed implicit none + integer (pInt), dimension(1) :: randInit integer (pInt) seed - call random_seed() + if (fixedSeed > 0_pInt) then + randInit = fixedSeed + call random_seed(put=randInit) + else + call random_seed() + endif + call get_seed(seed) - seed = 1 + if (fixedSeed > 0_pInt) seed = int(dble(fixedSeed)/2.0) + 1_pInt + call halton_seed_set(seed) call halton_ndim_set(3) diff --git a/code/mpie_cpfem_marc.f90 b/code/mpie_cpfem_marc.f90 index c3efcff53..3623f201f 100644 --- a/code/mpie_cpfem_marc.f90 +++ b/code/mpie_cpfem_marc.f90 @@ -35,9 +35,9 @@ !******************************************************************** ! include "prec.f90" ! uses nothing else - include "math.f90" ! uses prec - include "IO.f90" ! uses prec, math + include "IO.f90" ! uses prec include "numerics.f90" ! uses prec, IO + include "math.f90" ! uses prec, numerics include "debug.f90" ! uses prec, numerics include "FEsolving.f90" ! uses prec, IO include "mesh.f90" ! uses prec, math, IO, FEsolving diff --git a/code/numerics.config b/code/numerics.config index ccf994419..6d94fa4b3 100644 --- a/code/numerics.config +++ b/code/numerics.config @@ -1,22 +1,24 @@ ### numerical parameters ### -relevantStrain 1.0e-7 # strain increment considered significant -iJacoStiffness 1 # frequency of stiffness update -iJacoLpresiduum 1 # frequency of Jacobian update of residuum in Lp -pert_Fg 1.0e-6 # strain perturbation for FEM Jacobi -nHomog 20 # homogenization loop limit -nMPstate 10 # material point state loop limit -nCryst 20 # crystallite loop limit (only for debugging info, real loop limit is "subStepMin") -nState 10 # state loop limit -nStress 40 # stress loop limit -subStepMin 1.0e-3 # minimum (relative) size of sub-step allowed during cutback in crystallite -rTol_crystalliteState 1.0e-6 # relative tolerance in crystallite state loop -rTol_crystalliteTemperature 1.0e-6 # relative tolerance in crystallite temperature loop -rTol_crystalliteStress 1.0e-6 # relative tolerance in crystallite stress loop -aTol_crystalliteStress 1.0e-8 # absolute tolerance in crystallite stress loop +relevantStrain 1.0e-9 # strain increment considered significant +iJacoStiffness 1 # frequency of stiffness update +iJacoLpresiduum 1 # frequency of Jacobian update of residuum in Lp +pert_Fg 1.0e-7 # strain perturbation for FEM Jacobi +nHomog 25 # homogenization loop limit +nCryst 20 # crystallite loop limit (only for debugging info, real loop limit is "subStepMin") +nState 40 # state loop limit +nStress 80 # stress loop limit +subStepMin 1.0e-3 # minimum (relative) size of sub-step allowed during cutback in crystallite +rTol_crystalliteState 1.0e-6 # relative tolerance in crystallite state loop +rTol_crystalliteStress 1.0e-6 # relative tolerance in crystallite stress loop +aTol_crystalliteStress 1.0e-8 # absolute tolerance in crystallite stress loop -resToler 1.0e-4 # relative tolerance of residual in GIA iteration -resAbsol 1.0e+2 # absolute tolerance of residual in GIA iteration (corresponds to ~1 Pa) -resBound 1.0e+1 # relative maximum value (upper bound) for GIA residual -NRiterMax 24 # maximum number of GIA iteration \ No newline at end of file +aTol_RGC 1.0e+4 # absolute tolerance of RGC residuum (in Pa) +rTol_RGC 1.0e-3 # relative ... +aMax_RGC 1.0e+12 # absolute upper-limit of RGC residuum (in Pa) +rMax_RGC 1.0e+3 # relative ... +perturbPenalty_RGC 1.0e-8 # perturbation for computing penalty tangent +relevantMismatch_RGC 1.0e-5 # minimum threshold of mismatch + +fixed_seed 1234 # put any number larger than zero, integer, if you want to have a pseudo random distribution \ No newline at end of file diff --git a/code/numerics.f90 b/code/numerics.f90 index 31e798d50..29480fec0 100644 --- a/code/numerics.f90 +++ b/code/numerics.f90 @@ -29,6 +29,9 @@ real(pReal) relevantStrain, & ! strain pPert_RGC, & ! perturbation for computing RGC penalty tangent xSmoo_RGC ! RGC penalty smoothing parameter (hyperbolic tangent) +!* Random seeding parameters: added <<>> +integer(pInt) fixedSeed ! fixed seeding for pseudo-random number generator + CONTAINS !******************************************* @@ -89,6 +92,9 @@ subroutine numerics_init() pPert_RGC = 1.0e-8 xSmoo_RGC = 1.0e-5 +!* Random seeding parameters: added <<>> + fixedSeed = 0_pInt + ! try to open the config file if(IO_open_file(fileunit,numerics_configFile)) then @@ -146,6 +152,9 @@ subroutine numerics_init() case ('relevantmismatch_rgc') xSmoo_RGC = IO_floatValue(line,positions,2) +!* Random seeding parameters: added <<>> + case ('fixed_seed') + fixedSeed = IO_floatValue(line,positions,2) endselect enddo 100 close(fileunit) @@ -181,6 +190,9 @@ subroutine numerics_init() write(6,'(a24,x,e8.1)') 'rMax_RGC: ',relMax_RGC write(6,'(a24,x,e8.1)') 'perturbPenalty_RGC: ',pPert_RGC write(6,'(a24,x,e8.1)') 'relevantMismatch_RGC: ',xSmoo_RGC + +!* Random seeding parameters: added <<>> + write(6,'(a24,x,i8)') 'fixed_seed: ',fixedSeed write(6,*) ! sanity check @@ -207,6 +219,7 @@ subroutine numerics_init() if (pPert_RGC <= 0.0_pReal) call IO_error(276) !! oops !! if (xSmoo_RGC <= 0.0_pReal) call IO_error(277) + if (fixedSeed <= 0_pInt) write(6,'(a)') 'Random is random!' endsubroutine END MODULE numerics