223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity
DAMASK_EICMD/src/crystallite.f90

2365 lines
127 KiB
Fortran
Raw Normal View History

doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH
!> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH
removed empty line remaining from old ID string
2016-02-29 18:56:06 +05:30
!> @author Franz Roters, Max-Planck-Institut für Eisenforschung GmbH
doxygen comments
2013-02-22 04:38:36 +05:30
!> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH
removed empty line remaining from old ID string
2016-02-29 18:56:06 +05:30
!> @author Christoph Kords, Max-Planck-Institut für Eisenforschung GmbH
!> @author Chen Zhang, Michigan State University
doxygen comments
2013-02-22 04:38:36 +05:30
!> @brief crystallite state integration functions and reporting of results
!--------------------------------------------------------------------------------------------------
* now able to choose method for state integration (integrator and integratorStiffness in numerics.config); standard is "Fixed Point Iteration", which is basically the same as before; others available are "Explicit Euler", "AdaptiveEuler", "Classical Runge-Kutta" and "Runge-Kutta Cash-Karp" * now remembering stiffness similar to how we do it for Lp etc.; avoids undefined stiffness values for nonconverged stiffness calculation * non-local stuff: * changed non-local kinetics (Gilman2002) * enforce zero shearrate for overall carrrier density below relevant density * enforce zero density for those states that become negative and were below relevant density before * dislocation velocity is not limited by V^(1/3) / dt anymore
2010-10-01 17:48:49 +05:30
some minor improvement on precision detection: checking only once (in prec and no longer in math and crystallite), added one more 4/8 switch for LAPACK, as there is no single precision FFTW, stopping compilation for spectral method if FLOAT=4 new function in IO to print integers without leading zeros, implemented it at some places in the new spectral solver (reporting still needs some serious polishing) updated preprocessing for documentation to handle precision correctly
2012-08-31 01:56:28 +05:30
module crystallite
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
use prec
use IO
use config
use debug
use numerics
use rotations
use math
use FEsolving
use material
use constitutive
mesh_elem and theMesh are deprecated
2019-06-07 02:19:17 +05:30
use discretization
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
use lattice
use plastic_nonlocal
separating
2019-06-06 12:04:01 +05:30
use geometry_plastic_nonlocal, only: &
further removal of mesh
2019-06-07 13:50:56 +05:30
nIPneighbors => geometry_plastic_nonlocal_nIPneighbors, &
separating
2019-06-06 12:04:01 +05:30
IPneighborhood => geometry_plastic_nonlocal_IPneighborhood
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
use HDF5_utilities
use results
fixed indentation
2019-04-03 16:24:07 +05:30
repetition not needed - implicit none at the beginning of the module is enough
2019-04-11 11:05:58 +05:30
implicit none
fixed indentation
2019-04-03 16:24:07 +05:30
private
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
character(len=64), dimension(:,:), allocatable :: &
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_output !< name of each post result output
integer, public, protected :: &
crystallite_maxSizePostResults !< description not available
integer, dimension(:), allocatable, public, protected :: &
crystallite_sizePostResults !< description not available
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
integer, dimension(:,:), allocatable :: &
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_sizePostResult !< description not available
real(pReal), dimension(:,:,:), allocatable, public :: &
crystallite_dt !< requested time increment of each grain
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
real(pReal), dimension(:,:,:), allocatable :: &
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_subdt, & !< substepped time increment of each grain
crystallite_subFrac, & !< already calculated fraction of increment
crystallite_subStep !< size of next integration step
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
type(rotation), dimension(:,:,:), allocatable :: &
rotation class has consistent set of conversions crytallite_oriention0 was essentially a copy of material_EulerAngles
2019-09-20 08:12:28 +05:30
crystallite_orientation !< current orientation
fixed indentation
2019-04-03 16:24:07 +05:30
real(pReal), dimension(:,:,:,:,:), allocatable, public, protected :: &
crystallite_Fe, & !< current "elastic" def grad (end of converged time step)
crystallite_P !< 1st Piola-Kirchhoff stress per grain
real(pReal), dimension(:,:,:,:,:), allocatable, public :: &
crystallite_S, & !< current 2nd Piola-Kirchhoff stress vector (end of converged time step)
crystallite_S0, & !< 2nd Piola-Kirchhoff stress vector at start of FE inc
crystallite_partionedS0, & !< 2nd Piola-Kirchhoff stress vector at start of homog inc
crystallite_Fp, & !< current plastic def grad (end of converged time step)
crystallite_Fp0, & !< plastic def grad at start of FE inc
crystallite_partionedFp0,& !< plastic def grad at start of homog inc
crystallite_Fi, & !< current intermediate def grad (end of converged time step)
crystallite_Fi0, & !< intermediate def grad at start of FE inc
crystallite_partionedFi0,& !< intermediate def grad at start of homog inc
crystallite_F0, & !< def grad at start of FE inc
crystallite_partionedF, & !< def grad to be reached at end of homog inc
crystallite_partionedF0, & !< def grad at start of homog inc
crystallite_Lp, & !< current plastic velocitiy grad (end of converged time step)
crystallite_Lp0, & !< plastic velocitiy grad at start of FE inc
crystallite_partionedLp0, & !< plastic velocity grad at start of homog inc
crystallite_Li, & !< current intermediate velocitiy grad (end of converged time step)
crystallite_Li0, & !< intermediate velocitiy grad at start of FE inc
crystallite_partionedLi0 !< intermediate velocity grad at start of homog inc
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
real(pReal), dimension(:,:,:,:,:), allocatable :: &
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_subS0, & !< 2nd Piola-Kirchhoff stress vector at start of crystallite inc
crystallite_invFp, & !< inverse of current plastic def grad (end of converged time step)
crystallite_subFp0,& !< plastic def grad at start of crystallite inc
crystallite_invFi, & !< inverse of current intermediate def grad (end of converged time step)
crystallite_subFi0,& !< intermediate def grad at start of crystallite inc
crystallite_subF, & !< def grad to be reached at end of crystallite inc
crystallite_subF0, & !< def grad at start of crystallite inc
crystallite_subLp0,& !< plastic velocity grad at start of crystallite inc
crystallite_subLi0 !< intermediate velocity grad at start of crystallite inc
real(pReal), dimension(:,:,:,:,:,:,:), allocatable, public :: &
crystallite_dPdF !< current individual dPdF per grain (end of converged time step)
logical, dimension(:,:,:), allocatable, public :: &
crystallite_requested !< used by upper level (homogenization) to request crystallite calculation
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
logical, dimension(:,:,:), allocatable :: &
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_converged, & !< convergence flag
crystallite_todo, & !< flag to indicate need for further computation
crystallite_localPlasticity !< indicates this grain to have purely local constitutive law
enum, bind(c)
enumerator :: undefined_ID, &
phase_ID, &
texture_ID, &
orientation_ID, &
grainrotation_ID, &
defgrad_ID, &
fe_ID, &
fp_ID, &
fi_ID, &
lp_ID, &
li_ID, &
p_ID, &
s_ID, &
elasmatrix_ID, &
neighboringip_ID, &
neighboringelement_ID
end enum
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
integer(kind(undefined_ID)),dimension(:,:), allocatable :: &
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_outputID !< ID of each post result output
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
type :: tOutput !< new requested output (per phase)
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
character(len=65536), allocatable, dimension(:) :: &
label
end type tOutput
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
type(tOutput), allocatable, dimension(:) :: output_constituent
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
type :: tNumerics
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
integer :: &
iJacoLpresiduum, & !< frequency of Jacobian update of residuum in Lp
nState, & !< state loop limit
nStress !< stress loop limit
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
real(pReal) :: &
subStepMinCryst, & !< minimum (relative) size of sub-step allowed during cutback
subStepSizeCryst, & !< size of first substep when cutback
subStepSizeLp, & !< size of first substep when cutback in Lp calculation
subStepSizeLi, & !< size of first substep when cutback in Li calculation
stepIncreaseCryst, & !< increase of next substep size when previous substep converged
rTol_crystalliteState, & !< relative tolerance in state loop
rTol_crystalliteStress, & !< relative tolerance in stress loop
aTol_crystalliteStress !< absolute tolerance in stress loop
end type tNumerics
type(tNumerics) :: num ! numerics parameters. Better name?
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
fixed indentation
2019-04-03 16:24:07 +05:30
procedure(), pointer :: integrateState
public :: &
crystallite_init, &
crystallite_stress, &
crystallite_stressTangent, &
crystallite_orientations, &
crystallite_push33ToRef, &
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
crystallite_postResults, &
crystallite_results
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
some minor improvement on precision detection: checking only once (in prec and no longer in math and crystallite), added one more 4/8 switch for LAPACK, as there is no single precision FFTW, stopping compilation for spectral method if FLOAT=4 new function in IO to print integers without leading zeros, implemented it at some places in the new spectral solver (reporting still needs some serious polishing) updated preprocessing for documentation to handle precision correctly
2012-08-31 01:56:28 +05:30
contains
major restructuring of code. homogenization as well as constitutive are now free to choose. the runtime got somewhat longer (25% on simple tests) compared to a hardcoded isostrain homogenization. this might be a point of further optimization at a later stage... please use homogenization_isostrain.f90 as starting point / example for future developments of homog-schemes. the homogenization scheme now can additionally output certain results. hence, the userdata structure at each integration point now looks like this: - sizeHomogPostResults - block of that size containing homogPostResults then for each grain: - sizeGrainPostResults - block of that size containing crystallitePostResults, which consist of: + phaseID + volFrac + Eulers (3) + any constitutive post results requested
2009-05-07 21:57:36 +05:30
* Introduced preprocessor directives in order to suppress compilation of most write statements when using openmp. This tremendously improves efficiency of parallelization. * Also added some more openmp directives to increase percentage of parallelized code. * "implicit none" was missing in two subroutines of homogenization and constitutive.
2011-03-29 12:57:19 +05:30
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief allocates and initialize per grain variables
!--------------------------------------------------------------------------------------------------
removed redundant crystallite_temperature/heat and materialpoint_heat and using temperature from the thermal module
2014-10-10 17:58:57 +05:30
subroutine crystallite_init
fixed indentation
2019-04-03 16:24:07 +05:30
integer, parameter :: FILEUNIT=434
logical, dimension(:,:), allocatable :: devNull
integer :: &
c, & !< counter in integration point component loop
i, & !< counter in integration point loop
e, & !< counter in element loop
o = 0, & !< counter in output loop
r, &
cMax, & !< maximum number of integration point components
iMax, & !< maximum number of integration points
eMax, & !< maximum number of elements
myNcomponents, & !< number of components at current IP
mySize
character(len=65536), dimension(:), allocatable :: str
write(6,'(/,a)') ' <<<+- crystallite init -+>>>'
cMax = homogenization_maxNgrains
mesh_elem and theMesh are deprecated
2019-06-07 02:19:17 +05:30
iMax = discretization_nIP
eMax = discretization_nElem
fixed indentation
2019-04-03 16:24:07 +05:30
allocate(crystallite_S0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_partionedS0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_S(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subS0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_P(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_F0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_partionedF0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_partionedF(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subF0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subF(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Fp0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_partionedFp0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subFp0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Fp(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_invFp(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Fi0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_partionedFi0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subFi0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Fi(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_invFi(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Fe(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Lp0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_partionedLp0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subLp0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Lp(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Li0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_partionedLi0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subLi0(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Li(3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_dPdF(3,3,3,3,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_dt(cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subdt(cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subFrac(cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_subStep(cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_orientation(cMax,iMax,eMax))
allocate(crystallite_localPlasticity(cMax,iMax,eMax), source=.true.)
allocate(crystallite_requested(cMax,iMax,eMax), source=.false.)
allocate(crystallite_todo(cMax,iMax,eMax), source=.false.)
allocate(crystallite_converged(cMax,iMax,eMax), source=.true.)
allocate(crystallite_output(maxval(crystallite_Noutput), &
size(config_crystallite))) ; crystallite_output = ''
allocate(crystallite_outputID(maxval(crystallite_Noutput), &
size(config_crystallite)), source=undefined_ID)
allocate(crystallite_sizePostResults(size(config_crystallite)),source=0)
allocate(crystallite_sizePostResult(maxval(crystallite_Noutput), &
size(config_crystallite)), source=0)
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
num%subStepMinCryst = config_numerics%getFloat('substepmincryst', defaultVal=1.0e-3_pReal)
num%subStepSizeCryst = config_numerics%getFloat('substepsizecryst', defaultVal=0.25_pReal)
num%stepIncreaseCryst = config_numerics%getFloat('stepincreasecryst', defaultVal=1.5_pReal)
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
num%subStepSizeLp = config_numerics%getFloat('substepsizelp', defaultVal=0.5_pReal)
num%subStepSizeLi = config_numerics%getFloat('substepsizeli', defaultVal=0.5_pReal)
num%rTol_crystalliteState = config_numerics%getFloat('rtol_crystallitestate', defaultVal=1.0e-6_pReal)
num%rTol_crystalliteStress = config_numerics%getFloat('rtol_crystallitestress',defaultVal=1.0e-6_pReal)
num%aTol_crystalliteStress = config_numerics%getFloat('atol_crystallitestress',defaultVal=1.0e-8_pReal)
num%iJacoLpresiduum = config_numerics%getInt ('ijacolpresiduum', defaultVal=1)
num%nState = config_numerics%getInt ('nstate', defaultVal=20)
num%nStress = config_numerics%getInt ('nstress', defaultVal=40)
if(num%subStepMinCryst <= 0.0_pReal) call IO_error(301,ext_msg='subStepMinCryst')
if(num%subStepSizeCryst <= 0.0_pReal) call IO_error(301,ext_msg='subStepSizeCryst')
if(num%stepIncreaseCryst <= 0.0_pReal) call IO_error(301,ext_msg='stepIncreaseCryst')
if(num%subStepSizeLp <= 0.0_pReal) call IO_error(301,ext_msg='subStepSizeLp')
if(num%subStepSizeLi <= 0.0_pReal) call IO_error(301,ext_msg='subStepSizeLi')
if(num%rTol_crystalliteState <= 0.0_pReal) call IO_error(301,ext_msg='rTol_crystalliteState')
if(num%rTol_crystalliteStress <= 0.0_pReal) call IO_error(301,ext_msg='rTol_crystalliteStress')
if(num%aTol_crystalliteStress <= 0.0_pReal) call IO_error(301,ext_msg='aTol_crystalliteStress')
if(num%iJacoLpresiduum < 1) call IO_error(301,ext_msg='iJacoLpresiduum')
if(num%nState < 1) call IO_error(301,ext_msg='nState')
if(num%nStress< 1) call IO_error(301,ext_msg='nStress')
fixed indentation
2019-04-03 16:24:07 +05:30
select case(numerics_integrator)
case(1)
integrateState => integrateStateFPI
case(2)
integrateState => integrateStateEuler
case(3)
integrateState => integrateStateAdaptiveEuler
case(4)
integrateState => integrateStateRK4
case(5)
integrateState => integrateStateRKCK45
end select
do c = 1, size(config_crystallite)
gfortran does not recognice an empty array as 'present'
2018-06-22 11:33:22 +05:30
#if defined(__GFORTRAN__)
fixed indentation
2019-04-03 16:24:07 +05:30
str = ['GfortranBug86277']
str = config_crystallite(c)%getStrings('(output)',defaultVal=str)
if (str(1) == 'GfortranBug86277') str = [character(len=65536)::]
gfortran does not recognice an empty array as 'present'
2018-06-22 11:33:22 +05:30
#else
fixed indentation
2019-04-03 16:24:07 +05:30
str = config_crystallite(c)%getStrings('(output)',defaultVal=[character(len=65536)::])
gfortran does not recognice an empty array as 'present'
2018-06-22 11:33:22 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
do o = 1, size(str)
crystallite_output(o,c) = str(o)
outputName: select case(str(o))
case ('phase') outputName
crystallite_outputID(o,c) = phase_ID
case ('texture') outputName
crystallite_outputID(o,c) = texture_ID
case ('orientation') outputName
crystallite_outputID(o,c) = orientation_ID
case ('grainrotation') outputName
crystallite_outputID(o,c) = grainrotation_ID
case ('defgrad','f') outputName ! ToDo: no alias (f only)
crystallite_outputID(o,c) = defgrad_ID
case ('fe') outputName
crystallite_outputID(o,c) = fe_ID
case ('fp') outputName
crystallite_outputID(o,c) = fp_ID
case ('fi') outputName
crystallite_outputID(o,c) = fi_ID
case ('lp') outputName
crystallite_outputID(o,c) = lp_ID
case ('li') outputName
crystallite_outputID(o,c) = li_ID
case ('p','firstpiola','1stpiola') outputName ! ToDo: no alias (p only)
crystallite_outputID(o,c) = p_ID
case ('s','tstar','secondpiola','2ndpiola') outputName ! ToDo: no alias (s only)
crystallite_outputID(o,c) = s_ID
case ('elasmatrix') outputName
crystallite_outputID(o,c) = elasmatrix_ID
case ('neighboringip') outputName ! ToDo: this is not a result, it is static. Should be written out by mesh
crystallite_outputID(o,c) = neighboringip_ID
case ('neighboringelement') outputName ! ToDo: this is not a result, it is static. Should be written out by mesh
crystallite_outputID(o,c) = neighboringelement_ID
case default outputName
call IO_error(105,ext_msg=trim(str(o))//' (Crystallite)')
end select outputName
enddo
enddo
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
allocate(output_constituent(size(config_phase)))
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
do c = 1, size(config_phase)
#if defined(__GFORTRAN__)
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
allocate(output_constituent(c)%label(1))
output_constituent(c)%label(1)= 'GfortranBug86277'
output_constituent(c)%label = config_phase(c)%getStrings('(output)',defaultVal=output_constituent(c)%label )
if (output_constituent(c)%label (1) == 'GfortranBug86277') output_constituent(c)%label = [character(len=pStringLen)::]
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
#else
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
output_constituent(c)%label = config_phase(c)%getStrings('(output)',defaultVal=[character(len=pStringLen)::])
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
#endif
enddo
fixed indentation
2019-04-03 16:24:07 +05:30
do r = 1,size(config_crystallite)
do o = 1,crystallite_Noutput(r)
select case(crystallite_outputID(o,r))
specific for nonlocal, can be calculated during post processing
2019-05-14 10:52:29 +05:30
case(phase_ID,texture_ID)
fixed indentation
2019-04-03 16:24:07 +05:30
mySize = 1
case(orientation_ID,grainrotation_ID)
mySize = 4
case(defgrad_ID,fe_ID,fp_ID,fi_ID,lp_ID,li_ID,p_ID,s_ID)
mySize = 9
case(elasmatrix_ID)
mySize = 36
case(neighboringip_ID,neighboringelement_ID)
further removal of mesh
2019-06-07 13:50:56 +05:30
mySize = nIPneighbors
fixed indentation
2019-04-03 16:24:07 +05:30
case default
mySize = 0
end select
crystallite_sizePostResult(o,r) = mySize
crystallite_sizePostResults(r) = crystallite_sizePostResults(r) + mySize
enddo
enddo
crystallite_maxSizePostResults = &
maxval(crystallite_sizePostResults(microstructure_crystallite),microstructure_active)
using array size not to rely on correctly set variable, cleaning crystallite data after use
2018-06-27 00:03:41 +05:30
Tstar_v is now restored at different stages of cutbacking; we need this because preguess of state relies on consistent Tstar_v
2009-06-16 14:33:30 +05:30
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
!--------------------------------------------------------------------------------------------------
reworked crystallite part to allow for flexible user output --> new "crystallite" part in config file --> new "crystallite" option for microstructures --> new output file "...job.outputCrystallite" to be used in conjunction with marc_addUserOutput for meaningful naming of User Defined Vars.
2010-02-25 23:09:11 +05:30
! write description file for crystallite output
fixed indentation
2019-04-03 16:24:07 +05:30
if (worldrank == 0) then
call IO_write_jobFile(FILEUNIT,'outputCrystallite')
do r = 1,size(config_crystallite)
microstructure is only a property of the discretization will be used only by material.f90 once crystallite is removed
2019-06-07 11:14:34 +05:30
if (any(microstructure_crystallite(discretization_microstructureAt) == r)) then
consistent names
2019-06-15 17:56:32 +05:30
write(FILEUNIT,'(/,a,/)') '['//trim(config_name_crystallite(r))//']'
fixed indentation
2019-04-03 16:24:07 +05:30
do o = 1,crystallite_Noutput(r)
write(FILEUNIT,'(a,i4)') trim(crystallite_output(o,r))//char(9),crystallite_sizePostResult(o,r)
enddo
endif
enddo
close(FILEUNIT)
endif
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
call config_deallocate('material.config/phase')
fixed indentation
2019-04-03 16:24:07 +05:30
call config_deallocate('material.config/crystallite')
using array size not to rely on correctly set variable, cleaning crystallite data after use
2018-06-27 00:03:41 +05:30
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
!--------------------------------------------------------------------------------------------------
! initialize
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
!$OMP PARALLEL DO PRIVATE(myNcomponents,i,c)
fixed indentation
2019-04-03 16:24:07 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
myNcomponents = homogenization_Ngrains(material_homogenizationAt(e))
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
do i = FEsolving_execIP(1,e), FEsolving_execIP(2,e); do c = 1, myNcomponents
this fix restores the previous behavior of dislotwin very strange behavior in plasticity_detect changes for dislotwin. This fix restores the old behavior with respect to number of iterations. Using orientation%fromMatrix() results in much more iterations and an increase in runtime for the test from approx 6 min to 40 min. Results still match in the end. Also, crystallite_Fp0 differs by approx 1e-15 only between the two methods. I assume that something is wrong with either dislotwin or the state/stress integration
2019-09-22 20:16:30 +05:30
crystallite_Fp0(1:3,1:3,c,i,e) = math_EulerToR(material_Eulers(1:3,c,i,e)) ! plastic def gradient reflects init orientation
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_Fi0(1:3,1:3,c,i,e) = constitutive_initialFi(c,i,e)
crystallite_F0(1:3,1:3,c,i,e) = math_I3
using new names
2019-06-14 12:19:05 +05:30
crystallite_localPlasticity(c,i,e) = phase_localPlasticity(material_phaseAt(c,e))
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_Fe(1:3,1:3,c,i,e) = math_inv33(matmul(crystallite_Fi0(1:3,1:3,c,i,e), &
preparing for PGI compiler
2019-04-03 21:54:15 +05:30
crystallite_Fp0(1:3,1:3,c,i,e))) ! assuming that euler angles are given in internal strain free configuration
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_Fp(1:3,1:3,c,i,e) = crystallite_Fp0(1:3,1:3,c,i,e)
crystallite_Fi(1:3,1:3,c,i,e) = crystallite_Fi0(1:3,1:3,c,i,e)
crystallite_requested(c,i,e) = .true.
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
enddo; enddo
fixed indentation
2019-04-03 16:24:07 +05:30
enddo
!$OMP END PARALLEL DO
if(any(.not. crystallite_localPlasticity) .and. .not. usePingPong) call IO_error(601) ! exit if nonlocal but no ping-pong ToDo: Why not check earlier? or in nonlocal?
crystallite_partionedFp0 = crystallite_Fp0
crystallite_partionedFi0 = crystallite_Fi0
crystallite_partionedF0 = crystallite_F0
crystallite_partionedF = crystallite_F0
call crystallite_orientations()
!$OMP PARALLEL DO
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do c = 1,homogenization_Ngrains(material_homogenizationAt(e))
fixed indentation
2019-04-03 16:24:07 +05:30
call constitutive_microstructure(crystallite_Fe(1:3,1:3,c,i,e), &
crystallite_Fp(1:3,1:3,c,i,e), &
c,i,e) ! update dependent state variables to be consistent with basic states
enddo
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
enddo
fixed indentation
2019-04-03 16:24:07 +05:30
enddo
!$OMP END PARALLEL DO
devNull = crystallite_stress()
call crystallite_stressTangent
further separation still using old combined function
2019-01-18 19:12:44 +05:30
define debug variables only if needed
2018-09-20 09:39:02 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) then
write(6,'(a42,1x,i10)') ' # of elements: ', eMax
write(6,'(a42,1x,i10)') 'max # of integration points/element: ', iMax
write(6,'(a42,1x,i10)') 'max # of constituents/integration point: ', cMax
write(6,'(a42,1x,i10)') ' # of nonlocal constituents: ',count(.not. crystallite_localPlasticity)
flush(6)
endif
call debug_info
call debug_reset
define debug variables only if needed
2018-09-20 09:39:02 +05:30
#endif
* checked and corrected parallelization of code. compiles now successfully, but simply aborts computation with first parallel directive without any comment :-((( * also put a call to constitutive_microstructure at the start of each crystallite_integration subroutine like it was before. need that for nonlocal model in case of crystallite cutback
2010-11-03 22:52:48 +05:30
added new, flexible debugging scheme. now all modules have their own debug specification. compiles and runs, I hope nothing is broken did a lot of polishing
2012-03-09 01:55:28 +05:30
end subroutine crystallite_init
major restructuring of code. homogenization as well as constitutive are now free to choose. the runtime got somewhat longer (25% on simple tests) compared to a hardcoded isostrain homogenization. this might be a point of further optimization at a later stage... please use homogenization_isostrain.f90 as starting point / example for future developments of homog-schemes. the homogenization scheme now can additionally output certain results. hence, the userdata structure at each integration point now looks like this: - sizeHomogPostResults - block of that size containing homogPostResults then for each grain: - sizeGrainPostResults - block of that size containing crystallitePostResults, which consist of: + phaseID + volFrac + Eulers (3) + any constitutive post results requested
2009-05-07 21:57:36 +05:30
further separation still using old combined function
2019-01-18 19:12:44 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief calculate stress (P)
!--------------------------------------------------------------------------------------------------
use parameter structure and avoid conversion 33<->6
2019-02-17 16:45:46 +05:30
function crystallite_stress(dummyArgumentToPreventInternalCompilerErrorWithGCC)
fixed indentation
2019-04-03 16:24:07 +05:30
mesh_elem and theMesh are deprecated
2019-06-07 02:19:17 +05:30
logical, dimension(discretization_nIP,discretization_nElem) :: crystallite_stress
fixed indentation
2019-04-03 16:24:07 +05:30
real(pReal), intent(in), optional :: &
dummyArgumentToPreventInternalCompilerErrorWithGCC
real(pReal) :: &
formerSubStep
integer :: &
NiterationCrystallite, & ! number of iterations in crystallite loop
c, & !< counter in integration point component loop
i, & !< counter in integration point loop
e, & !< counter in element loop
startIP, endIP, &
s
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite),debug_levelSelective) /= 0 &
.and. FEsolving_execElem(1) <= debug_e &
.and. debug_e <= FEsolving_execElem(2)) then
write(6,'(/,a,i8,1x,i2,1x,i3)') '<< CRYST stress >> boundary and initial values at el ip ipc ', &
debug_e,debug_i, debug_g
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST stress >> F ', &
transpose(crystallite_partionedF(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST stress >> F0 ', &
transpose(crystallite_partionedF0(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST stress >> Fp0', &
transpose(crystallite_partionedFp0(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST stress >> Fi0', &
transpose(crystallite_partionedFi0(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST stress >> Lp0', &
transpose(crystallite_partionedLp0(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST stress >> Li0', &
transpose(crystallite_partionedLi0(1:3,1:3,debug_g,debug_i,debug_e))
endif
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#endif
!--------------------------------------------------------------------------------------------------
! initialize to starting condition
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_subStep = 0.0_pReal
!$OMP PARALLEL DO
elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do c = 1,homogenization_Ngrains(material_homogenizationAt(e))
fixed indentation
2019-04-03 16:24:07 +05:30
homogenizationRequestsCalculation: if (crystallite_requested(c,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
plasticState (material_phaseAt(c,e))%subState0( :,material_phaseMemberAt(c,i,e)) = &
plasticState (material_phaseAt(c,e))%partionedState0(:,material_phaseMemberAt(c,i,e))
fixed indentation
2019-04-03 16:24:07 +05:30
using new names
2019-06-14 12:19:05 +05:30
do s = 1, phase_Nsources(material_phaseAt(c,e))
sourceState(material_phaseAt(c,e))%p(s)%subState0( :,material_phaseMemberAt(c,i,e)) = &
sourceState(material_phaseAt(c,e))%p(s)%partionedState0(:,material_phaseMemberAt(c,i,e))
fixed indentation
2019-04-03 16:24:07 +05:30
enddo
crystallite_subFp0(1:3,1:3,c,i,e) = crystallite_partionedFp0(1:3,1:3,c,i,e)
crystallite_subLp0(1:3,1:3,c,i,e) = crystallite_partionedLp0(1:3,1:3,c,i,e)
crystallite_subFi0(1:3,1:3,c,i,e) = crystallite_partionedFi0(1:3,1:3,c,i,e)
crystallite_subLi0(1:3,1:3,c,i,e) = crystallite_partionedLi0(1:3,1:3,c,i,e)
crystallite_subF0(1:3,1:3,c,i,e) = crystallite_partionedF0(1:3,1:3,c,i,e)
crystallite_subS0(1:3,1:3,c,i,e) = crystallite_partionedS0(1:3,1:3,c,i,e)
crystallite_subFrac(c,i,e) = 0.0_pReal
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
crystallite_subStep(c,i,e) = 1.0_pReal/num%subStepSizeCryst
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_todo(c,i,e) = .true.
crystallite_converged(c,i,e) = .false. ! pretend failed step of 1/subStepSizeCryst
endif homogenizationRequestsCalculation
enddo; enddo
enddo elementLooping1
!$OMP END PARALLEL DO
singleRun: if (FEsolving_execELem(1) == FEsolving_execElem(2) .and. &
FEsolving_execIP(1,FEsolving_execELem(1))==FEsolving_execIP(2,FEsolving_execELem(1))) then
startIP = FEsolving_execIP(1,FEsolving_execELem(1))
endIP = startIP
else singleRun
startIP = 1
mesh_elem and theMesh are deprecated
2019-06-07 02:19:17 +05:30
endIP = discretization_nIP
fixed indentation
2019-04-03 16:24:07 +05:30
endif singleRun
NiterationCrystallite = 0
cutbackLooping: do while (any(crystallite_todo(:,startIP:endIP,FEsolving_execELem(1):FEsolving_execElem(2))))
NiterationCrystallite = NiterationCrystallite + 1
more elaborate debug labelling
2019-02-27 00:17:46 +05:30
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite),debug_levelExtensive) /= 0) &
write(6,'(a,i6)') '<< CRYST stress >> crystallite iteration ',NiterationCrystallite
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
!$OMP PARALLEL DO PRIVATE(formerSubStep)
elementLooping3: do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do c = 1,homogenization_Ngrains(material_homogenizationAt(e))
further separation still using old combined function
2019-01-18 19:12:44 +05:30
!--------------------------------------------------------------------------------------------------
! wind forward
fixed indentation
2019-04-03 16:24:07 +05:30
if (crystallite_converged(c,i,e)) then
formerSubStep = crystallite_subStep(c,i,e)
crystallite_subFrac(c,i,e) = crystallite_subFrac(c,i,e) + crystallite_subStep(c,i,e)
crystallite_subStep(c,i,e) = min(1.0_pReal - crystallite_subFrac(c,i,e), &
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
num%stepIncreaseCryst * crystallite_subStep(c,i,e))
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_todo(c,i,e) = crystallite_subStep(c,i,e) > 0.0_pReal ! still time left to integrate on?
if (crystallite_todo(c,i,e)) then
crystallite_subF0 (1:3,1:3,c,i,e) = crystallite_subF(1:3,1:3,c,i,e)
crystallite_subLp0(1:3,1:3,c,i,e) = crystallite_Lp (1:3,1:3,c,i,e)
crystallite_subLi0(1:3,1:3,c,i,e) = crystallite_Li (1:3,1:3,c,i,e)
crystallite_subFp0(1:3,1:3,c,i,e) = crystallite_Fp (1:3,1:3,c,i,e)
crystallite_subFi0(1:3,1:3,c,i,e) = crystallite_Fi (1:3,1:3,c,i,e)
crystallite_subS0 (1:3,1:3,c,i,e) = crystallite_S (1:3,1:3,c,i,e)
!if abbrevation, make c and p private in omp
using new names
2019-06-14 12:19:05 +05:30
plasticState( material_phaseAt(c,e))%subState0(:,material_phaseMemberAt(c,i,e)) &
= plasticState(material_phaseAt(c,e))%state( :,material_phaseMemberAt(c,i,e))
do s = 1, phase_Nsources(material_phaseAt(c,e))
sourceState( material_phaseAt(c,e))%p(s)%subState0(:,material_phaseMemberAt(c,i,e)) &
= sourceState(material_phaseAt(c,e))%p(s)%state( :,material_phaseMemberAt(c,i,e))
fixed indentation
2019-04-03 16:24:07 +05:30
enddo
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite),debug_levelBasic) /= 0 &
.and. ((e == debug_e .and. i == debug_i .and. c == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) &
write(6,'(a,f12.8,a,f12.8,a,i8,1x,i2,1x,i3,/)') '<< CRYST stress >> winding forward from ', &
crystallite_subFrac(c,i,e)-formerSubStep,' to current crystallite_subfrac ', &
crystallite_subFrac(c,i,e),' in crystallite_stress at el ip ipc ',e,i,c
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
endif
further separation still using old combined function
2019-01-18 19:12:44 +05:30
!--------------------------------------------------------------------------------------------------
! cut back (reduced time and restore)
fixed indentation
2019-04-03 16:24:07 +05:30
else
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
crystallite_subStep(c,i,e) = num%subStepSizeCryst * crystallite_subStep(c,i,e)
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_Fp (1:3,1:3,c,i,e) = crystallite_subFp0(1:3,1:3,c,i,e)
crystallite_invFp(1:3,1:3,c,i,e) = math_inv33(crystallite_Fp (1:3,1:3,c,i,e))
crystallite_Fi (1:3,1:3,c,i,e) = crystallite_subFi0(1:3,1:3,c,i,e)
crystallite_invFi(1:3,1:3,c,i,e) = math_inv33(crystallite_Fi (1:3,1:3,c,i,e))
crystallite_S (1:3,1:3,c,i,e) = crystallite_S0 (1:3,1:3,c,i,e)
if (crystallite_subStep(c,i,e) < 1.0_pReal) then ! actual (not initial) cutback
crystallite_Lp (1:3,1:3,c,i,e) = crystallite_subLp0(1:3,1:3,c,i,e)
crystallite_Li (1:3,1:3,c,i,e) = crystallite_subLi0(1:3,1:3,c,i,e)
endif
using new names
2019-06-14 12:19:05 +05:30
plasticState (material_phaseAt(c,e))%state( :,material_phaseMemberAt(c,i,e)) &
= plasticState(material_phaseAt(c,e))%subState0(:,material_phaseMemberAt(c,i,e))
do s = 1, phase_Nsources(material_phaseAt(c,e))
sourceState( material_phaseAt(c,e))%p(s)%state( :,material_phaseMemberAt(c,i,e)) &
= sourceState(material_phaseAt(c,e))%p(s)%subState0(:,material_phaseMemberAt(c,i,e))
fixed indentation
2019-04-03 16:24:07 +05:30
enddo
! cant restore dotState here, since not yet calculated in first cutback after initialization
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
crystallite_todo(c,i,e) = crystallite_subStep(c,i,e) > num%subStepMinCryst ! still on track or already done (beyond repair)
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((e == debug_e .and. i == debug_i .and. c == debug_g) &
.or. .not. iand(debug_level(debug_crystallite),debug_levelSelective) /= 0)) then
if (crystallite_todo(c,i,e)) then
write(6,'(a,f12.8,a,i8,1x,i2,1x,i3,/)') '<< CRYST stress >> cutback with new crystallite_subStep: ', &
crystallite_subStep(c,i,e),' at el ip ipc ',e,i,c
else
write(6,'(a,i8,1x,i2,1x,i3,/)') '<< CRYST stress >> reached minimum step size at el ip ipc ',e,i,c
endif
endif
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
endif
further separation still using old combined function
2019-01-18 19:12:44 +05:30
!--------------------------------------------------------------------------------------------------
! prepare for integration
fixed indentation
2019-04-03 16:24:07 +05:30
if (crystallite_todo(c,i,e)) then
crystallite_subF(1:3,1:3,c,i,e) = crystallite_subF0(1:3,1:3,c,i,e) &
+ crystallite_subStep(c,i,e) * (crystallite_partionedF (1:3,1:3,c,i,e) &
- crystallite_partionedF0(1:3,1:3,c,i,e))
crystallite_Fe(1:3,1:3,c,i,e) = matmul(matmul(crystallite_subF (1:3,1:3,c,i,e), &
preparing for PGI compiler
2019-04-03 21:54:15 +05:30
crystallite_invFp(1:3,1:3,c,i,e)), &
crystallite_invFi(1:3,1:3,c,i,e))
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_subdt(c,i,e) = crystallite_subStep(c,i,e) * crystallite_dt(c,i,e)
crystallite_converged(c,i,e) = .false.
endif
enddo
enddo
enddo elementLooping3
!$OMP END PARALLEL DO
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite),debug_levelExtensive) /= 0) then
write(6,'(/,a,f8.5,a,f8.5,/)') '<< CRYST stress >> ',minval(crystallite_subStep), &
' ≤ subStep ≤ ',maxval(crystallite_subStep)
write(6,'(/,a,f8.5,a,f8.5,/)') '<< CRYST stress >> ',minval(crystallite_subFrac), &
' ≤ subFrac ≤ ',maxval(crystallite_subFrac)
flush(6)
if (iand(debug_level(debug_crystallite),debug_levelSelective) /= 0) then
write(6,'(/,a,f8.5,1x,a,1x,f8.5,1x,a)') '<< CRYST stress >> subFrac + subStep = ',&
crystallite_subFrac(debug_g,debug_i,debug_e),'+',crystallite_subStep(debug_g,debug_i,debug_e),'@selective'
flush(6)
endif
endif
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#endif
!--------------------------------------------------------------------------------------------------
! integrate --- requires fully defined state array (basic + dependent state)
fixed indentation
2019-04-03 16:24:07 +05:30
if (any(crystallite_todo)) call integrateState ! TODO: unroll into proper elementloop to avoid N^2 for single point evaluation
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
where(.not. crystallite_converged .and. crystallite_subStep > num%subStepMinCryst) & ! do not try non-converged but fully cutbacked any further
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_todo = .true. ! TODO: again unroll this into proper elementloop to avoid N^2 for single point evaluation
further separation still using old combined function
2019-01-18 19:12:44 +05:30
fixed indentation
2019-04-03 16:24:07 +05:30
enddo cutbackLooping
further separation still using old combined function
2019-01-18 19:12:44 +05:30
! return whether converged or not
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_stress = .false.
elementLooping5: do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
crystallite_stress(i,e) = all(crystallite_converged(:,i,e))
enddo
enddo elementLooping5
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
elementLooping6: do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do c = 1,homogenization_Ngrains(material_homogenizationAt(e))
fixed indentation
2019-04-03 16:24:07 +05:30
if (.not. crystallite_converged(c,i,e)) then
if(iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) &
write(6,'(a,i8,1x,i2,1x,i3,/)') '<< CRYST stress >> no convergence at el ip ipc ', &
e,i,c
endif
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((e == debug_e .and. i == debug_i .and. c == debug_g) &
.or. .not. iand(debug_level(debug_crystallite),debug_levelSelective) /= 0)) then
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST stress >> solution at el ip ipc ',e,i,c
write(6,'(/,a,/,3(12x,3(f12.4,1x)/))') '<< CRYST stress >> P / MPa', &
transpose(crystallite_P(1:3,1:3,c,i,e))*1.0e-6_pReal
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST stress >> Fp', &
transpose(crystallite_Fp(1:3,1:3,c,i,e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST stress >> Fi', &
transpose(crystallite_Fi(1:3,1:3,c,i,e))
write(6,'(a,/,3(12x,3(f14.9,1x)/),/)') '<< CRYST stress >> Lp', &
transpose(crystallite_Lp(1:3,1:3,c,i,e))
write(6,'(a,/,3(12x,3(f14.9,1x)/),/)') '<< CRYST stress >> Li', &
transpose(crystallite_Li(1:3,1:3,c,i,e))
flush(6)
endif
enddo
enddo
enddo elementLooping6
further separation still using old combined function
2019-01-18 19:12:44 +05:30
#endif
end function crystallite_stress
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief calculate tangent (dPdF)
!--------------------------------------------------------------------------------------------------
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
subroutine crystallite_stressTangent
fixed indentation
2019-04-03 16:24:07 +05:30
integer :: &
c, & !< counter in integration point component loop
i, & !< counter in integration point loop
e, & !< counter in element loop
o, &
p
real(pReal), dimension(3,3) :: temp_33_1, devNull,invSubFi0, temp_33_2, temp_33_3, temp_33_4
real(pReal), dimension(3,3,3,3) :: dSdFe, &
dSdF, &
dSdFi, &
dLidS, &
dLidFi, &
dLpdS, &
dLpdFi, &
dFidS, &
dFpinvdF, &
rhs_3333, &
lhs_3333, &
temp_3333
real(pReal), dimension(9,9):: temp_99
logical :: error
!$OMP PARALLEL DO PRIVATE(dSdF,dSdFe,dSdFi,dLpdS,dLpdFi,dFpinvdF,dLidS,dLidFi,dFidS,invSubFi0,o,p, &
!$OMP rhs_3333,lhs_3333,temp_99,temp_33_1,temp_33_2,temp_33_3,temp_33_4,temp_3333,error)
elementLooping: do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do c = 1,homogenization_Ngrains(material_homogenizationAt(e))
fixed indentation
2019-04-03 16:24:07 +05:30
call constitutive_SandItsTangents(devNull,dSdFe,dSdFi, &
crystallite_Fe(1:3,1:3,c,i,e), &
crystallite_Fi(1:3,1:3,c,i,e),c,i,e) ! call constitutive law to calculate elastic stress tangent
call constitutive_LiAndItsTangents(devNull,dLidS,dLidFi, &
crystallite_S (1:3,1:3,c,i,e), &
crystallite_Fi(1:3,1:3,c,i,e), &
c,i,e) ! call constitutive law to calculate Li tangent in lattice configuration
if (sum(abs(dLidS)) < tol_math_check) then
dFidS = 0.0_pReal
else
invSubFi0 = math_inv33(crystallite_subFi0(1:3,1:3,c,i,e))
lhs_3333 = 0.0_pReal; rhs_3333 = 0.0_pReal
do o=1,3; do p=1,3
lhs_3333(1:3,1:3,o,p) = lhs_3333(1:3,1:3,o,p) &
+ crystallite_subdt(c,i,e)*matmul(invSubFi0,dLidFi(1:3,1:3,o,p))
lhs_3333(1:3,o,1:3,p) = lhs_3333(1:3,o,1:3,p) &
+ crystallite_invFi(1:3,1:3,c,i,e)*crystallite_invFi(p,o,c,i,e)
rhs_3333(1:3,1:3,o,p) = rhs_3333(1:3,1:3,o,p) &
- crystallite_subdt(c,i,e)*matmul(invSubFi0,dLidS(1:3,1:3,o,p))
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
enddo; enddo
using newer interface
2019-09-21 06:50:33 +05:30
call math_invert(temp_99,error,math_3333to99(lhs_3333))
fixed indentation
2019-04-03 16:24:07 +05:30
if (error) then
call IO_warning(warning_ID=600,el=e,ip=i,g=c, &
ext_msg='inversion error in analytic tangent calculation')
dFidS = 0.0_pReal
else
dFidS = math_mul3333xx3333(math_99to3333(temp_99),rhs_3333)
endif
dLidS = math_mul3333xx3333(dLidFi,dFidS) + dLidS
endif
call constitutive_LpAndItsTangents(devNull,dLpdS,dLpdFi, &
crystallite_S (1:3,1:3,c,i,e), &
crystallite_Fi(1:3,1:3,c,i,e),c,i,e) ! call constitutive law to calculate Lp tangent in lattice configuration
dLpdS = math_mul3333xx3333(dLpdFi,dFidS) + dLpdS
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
!--------------------------------------------------------------------------------------------------
! calculate dSdF
fixed indentation
2019-04-03 16:24:07 +05:30
temp_33_1 = transpose(matmul(crystallite_invFp(1:3,1:3,c,i,e), &
crystallite_invFi(1:3,1:3,c,i,e)))
temp_33_2 = matmul( crystallite_subF (1:3,1:3,c,i,e), &
math_inv33(crystallite_subFp0(1:3,1:3,c,i,e)))
temp_33_3 = matmul(matmul(crystallite_subF (1:3,1:3,c,i,e), &
crystallite_invFp (1:3,1:3,c,i,e)), &
math_inv33(crystallite_subFi0(1:3,1:3,c,i,e)))
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
do o=1,3; do p=1,3
fixed indentation
2019-04-03 16:24:07 +05:30
rhs_3333(p,o,1:3,1:3) = matmul(dSdFe(p,o,1:3,1:3),temp_33_1)
temp_3333(1:3,1:3,p,o) = matmul(matmul(temp_33_2,dLpdS(1:3,1:3,p,o)), &
crystallite_invFi(1:3,1:3,c,i,e)) &
+ matmul(temp_33_3,dLidS(1:3,1:3,p,o))
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
enddo; enddo
fixed indentation
2019-04-03 16:24:07 +05:30
lhs_3333 = crystallite_subdt(c,i,e)*math_mul3333xx3333(dSdFe,temp_3333) &
+ math_mul3333xx3333(dSdFi,dFidS)
using newer interface
2019-09-21 06:50:33 +05:30
call math_invert(temp_99,error,math_identity2nd(9)+math_3333to99(lhs_3333))
fixed indentation
2019-04-03 16:24:07 +05:30
if (error) then
call IO_warning(warning_ID=600,el=e,ip=i,g=c, &
ext_msg='inversion error in analytic tangent calculation')
dSdF = rhs_3333
else
dSdF = math_mul3333xx3333(math_99to3333(temp_99),rhs_3333)
endif
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
!--------------------------------------------------------------------------------------------------
! calculate dFpinvdF
fixed indentation
2019-04-03 16:24:07 +05:30
temp_3333 = math_mul3333xx3333(dLpdS,dSdF)
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
do o=1,3; do p=1,3
fixed indentation
2019-04-03 16:24:07 +05:30
dFpinvdF(1:3,1:3,p,o) &
= -crystallite_subdt(c,i,e) &
* matmul(math_inv33(crystallite_subFp0(1:3,1:3,c,i,e)), &
matmul(temp_3333(1:3,1:3,p,o),crystallite_invFi(1:3,1:3,c,i,e)))
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
enddo; enddo
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
!--------------------------------------------------------------------------------------------------
! assemble dPdF
fixed indentation
2019-04-03 16:24:07 +05:30
temp_33_1 = matmul(crystallite_invFp(1:3,1:3,c,i,e), &
matmul(crystallite_S(1:3,1:3,c,i,e), &
transpose(crystallite_invFp(1:3,1:3,c,i,e))))
temp_33_2 = matmul(crystallite_S(1:3,1:3,c,i,e), &
transpose(crystallite_invFp(1:3,1:3,c,i,e)))
temp_33_3 = matmul(crystallite_subF(1:3,1:3,c,i,e), &
crystallite_invFp(1:3,1:3,c,i,e))
temp_33_4 = matmul(matmul(crystallite_subF(1:3,1:3,c,i,e), &
crystallite_invFp(1:3,1:3,c,i,e)), &
crystallite_S(1:3,1:3,c,i,e))
crystallite_dPdF(1:3,1:3,1:3,1:3,c,i,e) = 0.0_pReal
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
do p=1,3
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_dPdF(p,1:3,p,1:3,c,i,e) = transpose(temp_33_1)
enddo
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
do o=1,3; do p=1,3
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_dPdF(1:3,1:3,p,o,c,i,e) = crystallite_dPdF(1:3,1:3,p,o,c,i,e) + &
matmul(matmul(crystallite_subF(1:3,1:3,c,i,e),dFpinvdF(1:3,1:3,p,o)),temp_33_2) + &
matmul(matmul(temp_33_3,dSdF(1:3,1:3,p,o)),transpose(crystallite_invFp(1:3,1:3,c,i,e))) + &
matmul(temp_33_4,transpose(dFpinvdF(1:3,1:3,p,o)))
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
enddo; enddo
fixed indentation
2019-04-03 16:24:07 +05:30
enddo; enddo
enddo elementLooping
!$OMP END PARALLEL DO
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
end subroutine crystallite_stressTangent
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
!> @brief calculates orientations
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
subroutine crystallite_orientations
fixed indentation
2019-04-03 16:24:07 +05:30
integer &
c, & !< counter in integration point component loop
i, & !< counter in integration point loop
e !< counter in element loop
!$OMP PARALLEL DO
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do c = 1,homogenization_Ngrains(material_homogenizationAt(e))
same names as in python
2019-09-21 05:48:09 +05:30
call crystallite_orientation(c,i,e)%fromMatrix(transpose(math_rotationalPart33(crystallite_Fe(1:3,1:3,c,i,e))))
fixed indentation
2019-04-03 16:24:07 +05:30
enddo; enddo; enddo
!$OMP END PARALLEL DO
nonlocalPresent: if (any(plasticState%nonLocal)) then
!$OMP PARALLEL DO
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
using new names
2019-06-14 12:19:05 +05:30
if (plasticState(material_phaseAt(1,e))%nonLocal) & ! if nonlocal model
fixed indentation
2019-04-03 16:24:07 +05:30
call plastic_nonlocal_updateCompatibility(crystallite_orientation,i,e)
enddo; enddo
!$OMP END PARALLEL DO
endif nonlocalPresent
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
end subroutine crystallite_orientations
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief Map 2nd order tensor to reference config
!--------------------------------------------------------------------------------------------------
function crystallite_push33ToRef(ipc,ip,el, tensor33)
fixed indentation
2019-04-03 16:24:07 +05:30
real(pReal), dimension(3,3) :: crystallite_push33ToRef
real(pReal), dimension(3,3), intent(in) :: tensor33
real(pReal), dimension(3,3) :: T
integer, intent(in):: &
el, &
ip, &
ipc
same names as in python
2019-09-21 05:48:09 +05:30
T = matmul(material_orientation0(ipc,ip,el)%asMatrix(), & ! ToDo: initial orientation correct?
rotation class has consistent set of conversions crytallite_oriention0 was essentially a copy of material_EulerAngles
2019-09-20 08:12:28 +05:30
transpose(math_inv33(crystallite_subF(1:3,1:3,ipc,ip,el))))
fixed indentation
2019-04-03 16:24:07 +05:30
crystallite_push33ToRef = matmul(transpose(T),matmul(tensor33,T))
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
end function crystallite_push33ToRef
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief return results of particular grain
!--------------------------------------------------------------------------------------------------
function crystallite_postResults(ipc, ip, el)
white space adjustments
2019-06-15 19:14:15 +05:30
integer, intent(in):: &
el, & !< element index
ip, & !< integration point index
ipc !< grain index
real(pReal), dimension(1+crystallite_sizePostResults(microstructure_crystallite(discretization_microstructureAt(el))) + &
1+plasticState(material_phaseAt(ipc,el))%sizePostResults + &
sum(sourceState(material_phaseAt(ipc,el))%p(:)%sizePostResults)) :: &
crystallite_postResults
integer :: &
o, &
c, &
crystID, &
mySize, &
n
type(rotation) :: rot
crystID = microstructure_crystallite(discretization_microstructureAt(el))
crystallite_postResults = 0.0_pReal
crystallite_postResults(1) = real(crystallite_sizePostResults(crystID),pReal) ! header-like information (length)
c = 1
do o = 1,crystallite_Noutput(crystID)
mySize = 0
select case(crystallite_outputID(o,crystID))
case (phase_ID)
mySize = 1
crystallite_postResults(c+1) = real(material_phaseAt(ipc,el),pReal) ! phaseID of grain
case (texture_ID)
mySize = 1
crystallite_postResults(c+1) = real(material_texture(ipc,ip,el),pReal) ! textureID of grain
case (orientation_ID)
mySize = 4
crystallite_postResults(c+1:c+mySize) = crystallite_orientation(ipc,ip,el)%asQuaternion()
case (grainrotation_ID)
rotation class has consistent set of conversions crytallite_oriention0 was essentially a copy of material_EulerAngles
2019-09-20 08:12:28 +05:30
rot = material_orientation0(ipc,ip,el)%misorientation(crystallite_orientation(ipc,ip,el))
white space adjustments
2019-06-15 19:14:15 +05:30
mySize = 4
same names as in python
2019-09-21 05:48:09 +05:30
crystallite_postResults(c+1:c+mySize) = rot%asAxisAngle()
white space adjustments
2019-06-15 19:14:15 +05:30
crystallite_postResults(c+4) = inDeg * crystallite_postResults(c+4) ! angle in degree
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
! remark: tensor output is of the form 11,12,13, 21,22,23, 31,32,33
! thus row index i is slow, while column index j is fast. reminder: "row is slow"
white space adjustments
2019-06-15 19:14:15 +05:30
case (defgrad_ID)
mySize = 9
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_partionedF(1:3,1:3,ipc,ip,el)),[mySize])
case (fe_ID)
mySize = 9
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Fe(1:3,1:3,ipc,ip,el)),[mySize])
case (fp_ID)
mySize = 9
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Fp(1:3,1:3,ipc,ip,el)),[mySize])
case (fi_ID)
mySize = 9
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Fi(1:3,1:3,ipc,ip,el)),[mySize])
case (lp_ID)
mySize = 9
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Lp(1:3,1:3,ipc,ip,el)),[mySize])
case (li_ID)
mySize = 9
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Li(1:3,1:3,ipc,ip,el)),[mySize])
case (p_ID)
mySize = 9
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_P(1:3,1:3,ipc,ip,el)),[mySize])
case (s_ID)
mySize = 9
crystallite_postResults(c+1:c+mySize) = &
reshape(crystallite_S(1:3,1:3,ipc,ip,el),[mySize])
case (elasmatrix_ID)
mySize = 36
crystallite_postResults(c+1:c+mySize) = reshape(constitutive_homogenizedC(ipc,ip,el),[mySize])
case(neighboringelement_ID)
mySize = nIPneighbors
crystallite_postResults(c+1:c+mySize) = 0.0_pReal
forall (n = 1:mySize) &
crystallite_postResults(c+n) = real(IPneighborhood(1,n,ip,el),pReal)
case(neighboringip_ID)
mySize = nIPneighbors
crystallite_postResults(c+1:c+mySize) = 0.0_pReal
forall (n = 1:mySize) &
crystallite_postResults(c+n) = real(IPneighborhood(2,n,ip,el),pReal)
end select
c = c + mySize
enddo
crystallite_postResults(c+1) = real(plasticState(material_phaseAt(ipc,el))%sizePostResults,pReal) ! size of constitutive results
c = c + 1
if (size(crystallite_postResults)-c > 0) &
crystallite_postResults(c+1:size(crystallite_postResults)) = &
constitutive_postResults(crystallite_S(1:3,1:3,ipc,ip,el), crystallite_Fi(1:3,1:3,ipc,ip,el), &
ipc, ip, el)
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
end function crystallite_postResults
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
!--------------------------------------------------------------------------------------------------
leaner group structure, centrally handled
2019-04-13 13:17:56 +05:30
!> @brief writes crystallite results to HDF5 output file
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
!--------------------------------------------------------------------------------------------------
subroutine crystallite_results
#if defined(PETSc) || defined(DAMASK_HDF5)
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
integer :: p,o
polishing
2019-04-07 17:58:08 +05:30
real(pReal), allocatable, dimension(:,:,:) :: selected_tensors
type(rotation), allocatable, dimension(:) :: selected_rotations
small polishing
2019-04-07 18:17:21 +05:30
character(len=256) :: group,lattice_label
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
do p=1,size(config_name_phase)
adjustments for easier access to output data
2019-04-18 15:25:50 +05:30
group = trim('current/constituent')//'/'//trim(config_name_phase(p))//'/generic'
call HDF5_closeGroup(results_addGroup(group))
leaner group structure, centrally handled
2019-04-13 13:17:56 +05:30
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
do o = 1, size(output_constituent(p)%label)
select case (output_constituent(p)%label(o))
case('f')
polishing
2019-04-07 17:58:08 +05:30
selected_tensors = select_tensors(crystallite_partionedF,p)
call results_writeDataset(group,selected_tensors,'F',&
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
'deformation gradient','1')
case('fe')
polishing
2019-04-07 17:58:08 +05:30
selected_tensors = select_tensors(crystallite_Fe,p)
call results_writeDataset(group,selected_tensors,'Fe',&
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
'elastic deformation gradient','1')
case('fp')
polishing
2019-04-07 17:58:08 +05:30
selected_tensors = select_tensors(crystallite_Fp,p)
call results_writeDataset(group,selected_tensors,'Fp',&
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
'plastic deformation gradient','1')
case('fi')
polishing
2019-04-07 17:58:08 +05:30
selected_tensors = select_tensors(crystallite_Fi,p)
call results_writeDataset(group,selected_tensors,'Fi',&
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
'inelastic deformation gradient','1')
case('lp')
polishing
2019-04-07 17:58:08 +05:30
selected_tensors = select_tensors(crystallite_Lp,p)
call results_writeDataset(group,selected_tensors,'Lp',&
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
'plastic velocity gradient','1/s')
case('li')
polishing
2019-04-07 17:58:08 +05:30
selected_tensors = select_tensors(crystallite_Li,p)
call results_writeDataset(group,selected_tensors,'Li',&
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
'inelastic velocity gradient','1/s')
case('p')
polishing
2019-04-07 17:58:08 +05:30
selected_tensors = select_tensors(crystallite_P,p)
call results_writeDataset(group,selected_tensors,'P',&
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
'1st Piola-Kirchoff stress','Pa')
case('s')
polishing
2019-04-07 17:58:08 +05:30
selected_tensors = select_tensors(crystallite_S,p)
call results_writeDataset(group,selected_tensors,'S',&
'2nd Piola-Kirchoff stress','Pa')
case('orientation')
small polishing
2019-04-07 18:17:21 +05:30
select case(lattice_structure(p))
case(LATTICE_iso_ID)
lattice_label = 'iso'
case(LATTICE_fcc_ID)
lattice_label = 'fcc'
case(LATTICE_bcc_ID)
lattice_label = 'bcc'
case(LATTICE_bct_ID)
lattice_label = 'bct'
case(LATTICE_hex_ID)
lattice_label = 'hex'
case(LATTICE_ort_ID)
lattice_label = 'ort'
end select
polishing
2019-04-07 17:58:08 +05:30
selected_rotations = select_rotations(crystallite_orientation,p)
call results_writeDataset(group,selected_rotations,'orientation',&
small polishing
2019-04-07 18:17:21 +05:30
'crystal orientation as quaternion',lattice_label)
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
end select
enddo
white space adjustments
2019-06-15 19:14:15 +05:30
enddo
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
white space adjustments
2019-06-15 19:14:15 +05:30
contains
polishing
2019-04-07 17:58:08 +05:30
white space adjustments
2019-06-15 19:14:15 +05:30
!------------------------------------------------------------------------------------------------
!> @brief select tensors for output
!------------------------------------------------------------------------------------------------
better readable
2019-04-13 19:04:51 +05:30
function select_tensors(dataset,instance)
polishing
2019-04-07 17:58:08 +05:30
better readable
2019-04-13 19:04:51 +05:30
integer, intent(in) :: instance
real(pReal), dimension(:,:,:,:,:), intent(in) :: dataset
real(pReal), allocatable, dimension(:,:,:) :: select_tensors
integer :: e,i,c,j
correct output for multipe integration points
2019-10-19 23:25:00 +05:30
allocate(select_tensors(3,3,count(material_phaseAt==instance)*homogenization_maxNgrains*discretization_nIP))
better readable
2019-04-13 19:04:51 +05:30
j=0
do e = 1, size(material_phaseAt,2)
correct output for multipe integration points
2019-10-19 23:25:00 +05:30
do i = 1, discretization_nIP
do c = 1, size(material_phaseAt,1) !ToDo: this needs to be changed for varying Ngrains
if (material_phaseAt(c,e) == instance) then
j = j + 1
select_tensors(1:3,1:3,j) = dataset(1:3,1:3,c,i,e)
endif
better readable
2019-04-13 19:04:51 +05:30
enddo
enddo
polishing
2019-04-07 17:58:08 +05:30
enddo
better readable
2019-04-13 19:04:51 +05:30
end function select_tensors
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
polishing
2019-04-07 17:58:08 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief select rotations for output
!--------------------------------------------------------------------------------------------------
better readable
2019-04-13 19:04:51 +05:30
function select_rotations(dataset,instance)
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
better readable
2019-04-13 19:04:51 +05:30
integer, intent(in) :: instance
type(rotation), dimension(:,:,:), intent(in) :: dataset
type(rotation), allocatable, dimension(:) :: select_rotations
integer :: e,i,c,j
correct output for multipe integration points
2019-10-19 23:25:00 +05:30
allocate(select_rotations(count(material_phaseAt==instance)*homogenization_maxNgrains*discretization_nIP))
better readable
2019-04-13 19:04:51 +05:30
j=0
do e = 1, size(material_phaseAt,2)
correct output for multipe integration points
2019-10-19 23:25:00 +05:30
do i = 1, discretization_nIP
do c = 1, size(material_phaseAt,1) !ToDo: this needs to be changed for varying Ngrains
better readable
2019-04-13 19:04:51 +05:30
if (material_phaseAt(c,e) == instance) then
constituent (ex crystallite) results are stored in HDF5 currently, not the best code but new structure for crystallite data will fix that output of orientations still missing
2019-04-06 15:36:34 +05:30
j = j + 1
better readable
2019-04-13 19:04:51 +05:30
select_rotations(j) = dataset(c,i,e)
endif
enddo
enddo
enddo
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
polishing
2019-04-07 17:58:08 +05:30
end function select_rotations
WIP: crystallite HDF5 results will be stored according to the phase sections
2019-04-06 10:01:02 +05:30
#endif
end subroutine crystallite_results
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief calculation of stress (P) with time integration based on a residuum in Lp and
!> intermediate acceleration of the Newton-Raphson correction
!--------------------------------------------------------------------------------------------------
fixed indentation
2019-04-03 16:24:07 +05:30
logical function integrateStress(ipc,ip,el,timeFraction)
integer, intent(in):: el, & ! element index
ip, & ! integration point index
ipc ! grain index
real(pReal), optional, intent(in) :: timeFraction ! fraction of timestep
real(pReal), dimension(3,3):: Fg_new, & ! deformation gradient at end of timestep
Fp_new, & ! plastic deformation gradient at end of timestep
Fe_new, & ! elastic deformation gradient at end of timestep
invFp_new, & ! inverse of Fp_new
Fi_new, & ! gradient of intermediate deformation stages
invFi_new, &
invFp_current, & ! inverse of Fp_current
invFi_current, & ! inverse of Fp_current
Lpguess, & ! current guess for plastic velocity gradient
Lpguess_old, & ! known last good guess for plastic velocity gradient
Lp_constitutive, & ! plastic velocity gradient resulting from constitutive law
residuumLp, & ! current residuum of plastic velocity gradient
residuumLp_old, & ! last residuum of plastic velocity gradient
deltaLp, & ! direction of next guess
Liguess, & ! current guess for intermediate velocity gradient
Liguess_old, & ! known last good guess for intermediate velocity gradient
Li_constitutive, & ! intermediate velocity gradient resulting from constitutive law
residuumLi, & ! current residuum of intermediate velocity gradient
residuumLi_old, & ! last residuum of intermediate velocity gradient
deltaLi, & ! direction of next guess
S, & ! 2nd Piola-Kirchhoff Stress in plastic (lattice) configuration
A, &
B, &
Fe, & ! elastic deformation gradient
temp_33
real(pReal), dimension(9) :: work ! needed for matrix inversion by LAPACK
integer, dimension(9) :: devNull ! needed for matrix inversion by LAPACK
real(pReal), dimension(9,9) :: dRLp_dLp, & ! partial derivative of residuum (Jacobian for Newton-Raphson scheme)
dRLp_dLp2, & ! working copy of dRdLp
dRLi_dLi ! partial derivative of residuumI (Jacobian for Newton-Raphson scheme)
real(pReal), dimension(3,3,3,3):: dS_dFe, & ! partial derivative of 2nd Piola-Kirchhoff stress
dS_dFi, &
dFe_dLp, & ! partial derivative of elastic deformation gradient
dFe_dLi, &
dFi_dLi, &
dLp_dFi, &
dLi_dFi, &
dLp_dS, &
dLi_dS
real(pReal) detInvFi, & ! determinant of InvFi
steplengthLp, &
steplengthLi, &
dt, & ! time increment
aTolLp, &
aTolLi
integer NiterationStressLp, & ! number of stress integrations
NiterationStressLi, & ! number of inner stress integrations
ierr, & ! error indicator for LAPACK
o, &
p, &
jacoCounterLp, &
jacoCounterLi ! counters to check for Jacobian update
external :: &
dgesv
!* be pessimistic
integrateStress = .false.
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) &
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST integrateStress >> at el ip ipc ',el,ip,ipc
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
changed to new state, please report bugs to Luv or Martin
2014-07-02 17:57:39 +05:30
fixed indentation
2019-04-03 16:24:07 +05:30
if (present(timeFraction)) then
dt = crystallite_subdt(ipc,ip,el) * timeFraction
Fg_new = crystallite_subF0(1:3,1:3,ipc,ip,el) &
+ (crystallite_subF(1:3,1:3,ipc,ip,el) - crystallite_subF0(1:3,1:3,ipc,ip,el)) * timeFraction
else
dt = crystallite_subdt(ipc,ip,el)
Fg_new = crystallite_subF(1:3,1:3,ipc,ip,el)
endif
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
polishing
2019-07-01 10:39:51 +05:30
Lpguess = crystallite_Lp(1:3,1:3,ipc,ip,el) ! take as first guess
Liguess = crystallite_Li(1:3,1:3,ipc,ip,el) ! take as first guess
fixed indentation
2019-04-03 16:24:07 +05:30
Liguess_old = Liguess
invFp_current = math_inv33(crystallite_subFp0(1:3,1:3,ipc,ip,el))
failedInversionFp: if (all(dEq0(invFp_current))) then
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) &
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST integrateStress >> failed on inversion of current Fp at el ip ipc ',&
el,ip,ipc
if (iand(debug_level(debug_crystallite), debug_levelExtensive) > 0) &
write(6,'(/,a,/,3(12x,3(f12.7,1x)/))') '<< CRYST >> current Fp ',transpose(crystallite_subFp0(1:3,1:3,ipc,ip,el))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
return
endif failedInversionFp
A = matmul(Fg_new,invFp_current) ! intermediate tensor needed later to calculate dFe_dLp
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
fixed indentation
2019-04-03 16:24:07 +05:30
invFi_current = math_inv33(crystallite_subFi0(1:3,1:3,ipc,ip,el))
failedInversionFi: if (all(dEq0(invFi_current))) then
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) &
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST integrateStress >> failed on inversion of current Fi at el ip ipc ',&
el,ip,ipc
if (iand(debug_level(debug_crystallite), debug_levelExtensive) > 0) &
write(6,'(/,a,/,3(12x,3(f12.7,1x)/))') '<< CRYST integrateStress >> current Fi ', &
transpose(crystallite_subFi0(1:3,1:3,ipc,ip,el))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
return
endif failedInversionFi
!* start Li loop with normal step length
NiterationStressLi = 0
jacoCounterLi = 0
steplengthLi = 1.0_pReal
residuumLi_old = 0.0_pReal
LiLoop: do
NiterationStressLi = NiterationStressLi + 1
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
LiLoopLimit: if (NiterationStressLi > num%nStress) then
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) &
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
write(6,'(a,i3,a,i8,1x,i2,1x,i3,/)') '<< CRYST integrateStress >> reached Li loop limit',num%nStress, &
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
' at el ip ipc ', el,ip,ipc
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
return
endif LiLoopLimit
invFi_new = matmul(invFi_current,math_I3 - dt*Liguess)
Fi_new = math_inv33(invFi_new)
detInvFi = math_det33(invFi_new)
!* start Lp loop with normal step length
NiterationStressLp = 0
jacoCounterLp = 0
steplengthLp = 1.0_pReal
residuumLp_old = 0.0_pReal
Lpguess_old = Lpguess
LpLoop: do
NiterationStressLp = NiterationStressLp + 1
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
LpLoopLimit: if (NiterationStressLp > num%nStress) then
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) &
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
write(6,'(a,i3,a,i8,1x,i2,1x,i3,/)') '<< CRYST integrateStress >> reached Lp loop limit',num%nStress, &
fixed indentation
2019-04-03 16:24:07 +05:30
' at el ip ipc ', el,ip,ipc
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
return
endif LpLoopLimit
!* calculate (elastic) 2nd Piola--Kirchhoff stress tensor and its tangent from constitutive law
B = math_I3 - dt*Lpguess
Fe = matmul(matmul(A,B), invFi_new)
call constitutive_SandItsTangents(S, dS_dFe, dS_dFi, &
Fe, Fi_new, ipc, ip, el) ! call constitutive law to calculate 2nd Piola-Kirchhoff stress and its derivative in unloaded configuration
!* calculate plastic velocity gradient and its tangent from constitutive law
call constitutive_LpAndItsTangents(Lp_constitutive, dLp_dS, dLp_dFi, &
S, Fi_new, ipc, ip, el)
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
added some lines to help while debugging
2019-09-07 01:15:49 +05:30
write(6,'(a,i3,/)') '<< CRYST integrateStress >> Lp iteration ', NiterationStressLp
fixed indentation
2019-04-03 16:24:07 +05:30
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST integrateStress >> Lpguess', transpose(Lpguess)
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST integrateStress >> Lp_constitutive', transpose(Lp_constitutive)
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST integrateStress >> Fi', transpose(Fi_new)
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST integrateStress >> Fe', transpose(Fe)
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST integrateStress >> S', transpose(S)
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
fixed indentation
2019-04-03 16:24:07 +05:30
!* update current residuum and check for convergence of loop
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
aTolLp = max(num%rTol_crystalliteStress * max(norm2(Lpguess),norm2(Lp_constitutive)), & ! absolute tolerance from largest acceptable relative error
num%aTol_crystalliteStress) ! minimum lower cutoff
fixed indentation
2019-04-03 16:24:07 +05:30
residuumLp = Lpguess - Lp_constitutive
if (any(IEEE_is_NaN(residuumLp))) then
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) &
write(6,'(a,i8,1x,i2,1x,i3,a,i3,a)') '<< CRYST integrateStress >> encountered NaN for Lp-residuum at el ip ipc ', &
el,ip,ipc, &
' ; iteration ', NiterationStressLp,&
' >> returning..!'
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
return ! ...me = .false. to inform integrator about problem
elseif (norm2(residuumLp) < aTolLp) then ! converged if below absolute tolerance
exit LpLoop ! ...leave iteration loop
elseif ( NiterationStressLp == 1 &
.or. norm2(residuumLp) < norm2(residuumLp_old)) then ! not converged, but improved norm of residuum (always proceed in first iteration)...
residuumLp_old = residuumLp ! ...remember old values and...
Lpguess_old = Lpguess
steplengthLp = 1.0_pReal ! ...proceed with normal step length (calculate new search direction)
else ! not converged and residuum not improved...
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
steplengthLp = num%subStepSizeLp * steplengthLp ! ...try with smaller step length in same direction
fixed indentation
2019-04-03 16:24:07 +05:30
Lpguess = Lpguess_old + steplengthLp * deltaLp
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
write(6,'(a,1x,f7.4)') '<< CRYST integrateStress >> linear search for Lpguess with step', steplengthLp
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
cycle LpLoop
endif
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!* calculate Jacobian for correction term
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
if (mod(jacoCounterLp, num%iJacoLpresiduum) == 0) then
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
do o=1,3; do p=1,3
dFe_dLp(o,1:3,p,1:3) = A(o,p)*transpose(invFi_new) ! dFe_dLp(i,j,k,l) = -dt * A(i,k) invFi(l,j)
enddo; enddo
fixed indentation
2019-04-03 16:24:07 +05:30
dFe_dLp = - dt * dFe_dLp
dRLp_dLp = math_identity2nd(9) &
- math_3333to99(math_mul3333xx3333(math_mul3333xx3333(dLp_dS,dS_dFe),dFe_dLp))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
write(6,'(a,/,9(12x,9(e12.4,1x)/))') '<< CRYST integrateStress >> dLp_dS', math_3333to99(dLp_dS)
write(6,'(a,1x,e20.10)') '<< CRYST integrateStress >> dLp_dS norm', norm2(math_3333to99(dLp_dS))
write(6,'(a,/,9(12x,9(e12.4,1x)/))') '<< CRYST integrateStress >> dRLp_dLp', dRLp_dLp-math_identity2nd(9)
write(6,'(a,1x,e20.10)') '<< CRYST integrateStress >> dRLp_dLp norm', norm2(dRLp_dLp-math_identity2nd(9))
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
dRLp_dLp2 = dRLp_dLp ! will be overwritten in first call to LAPACK routine
work = math_33to9(residuumLp)
call dgesv(9,1,dRLp_dLp2,9,devNull,work,9,ierr) ! solve dRLp/dLp * delta Lp = -res for delta Lp
if (ierr /= 0) then
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) then
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST integrateStress >> failed on dR/dLp inversion at el ip ipc ', &
el,ip,ipc
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g)&
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
write(6,*)
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST integrateStress >> dR_dLp',transpose(dRLp_dLp)
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST integrateStress >> dFe_dLp',transpose(math_3333to99(dFe_dLp))
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST integrateStress >> dS_dFe (cnst)',transpose(math_3333to99(dS_dFe))
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST integrateStress >> dLp_dS (cnst)',transpose(math_3333to99(dLp_dS))
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST integrateStress >> A',transpose(A)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST integrateStress >> B',transpose(B)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST integrateStress >> Lp_constitutive',transpose(Lp_constitutive)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST integrateStress >> Lpguess',transpose(Lpguess)
endif
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
return
endif
deltaLp = - math_9to33(work)
endif
jacoCounterLp = jacoCounterLp + 1
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
fixed indentation
2019-04-03 16:24:07 +05:30
Lpguess = Lpguess + steplengthLp * deltaLp
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
fixed indentation
2019-04-03 16:24:07 +05:30
enddo LpLoop
general polishing; simplification of RKCK45 dotstate assembly for different stages.
2014-09-03 01:41:57 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!* calculate intermediate velocity gradient and its tangent from constitutive law
fixed indentation
2019-04-03 16:24:07 +05:30
call constitutive_LiAndItsTangents(Li_constitutive, dLi_dS, dLi_dFi, &
S, Fi_new, ipc, ip, el)
general polishing; simplification of RKCK45 dotstate assembly for different stages.
2014-09-03 01:41:57 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
added some lines to help while debugging
2019-09-07 01:15:49 +05:30
write(6,'(a,i3,/)') '<< CRYST integrateStress >> Li iteration ', NiterationStressLi
fixed indentation
2019-04-03 16:24:07 +05:30
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST integrateStress >> Li_constitutive', transpose(Li_constitutive)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST integrateStress >> Liguess', transpose(Liguess)
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
fixed indentation
2019-04-03 16:24:07 +05:30
!* update current residuum and check for convergence of loop
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
aTolLi = max(num%rTol_crystalliteStress * max(norm2(Liguess),norm2(Li_constitutive)), & ! absolute tolerance from largest acceptable relative error
num%aTol_crystalliteStress) ! minimum lower cutoff
fixed indentation
2019-04-03 16:24:07 +05:30
residuumLi = Liguess - Li_constitutive
if (any(IEEE_is_NaN(residuumLi))) then ! NaN in residuum...
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) &
write(6,'(a,i8,1x,i2,1x,i3,a,i3,a)') '<< CRYST integrateStress >> encountered NaN for Li-residuum at el ip ipc ', &
el,ip,ipc, &
' ; iteration ', NiterationStressLi,&
' >> returning..!'
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
return ! ...me = .false. to inform integrator about problem
elseif (norm2(residuumLi) < aTolLi) then ! converged if below absolute tolerance
exit LiLoop ! ...leave iteration loop
elseif ( NiterationStressLi == 1 &
.or. norm2(residuumLi) < norm2(residuumLi_old)) then ! not converged, but improved norm of residuum (always proceed in first iteration)...
residuumLi_old = residuumLi ! ...remember old values and...
Liguess_old = Liguess
steplengthLi = 1.0_pReal ! ...proceed with normal step length (calculate new search direction)
else ! not converged and residuum not improved...
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
steplengthLi = num%subStepSizeLi * steplengthLi ! ...try with smaller step length in same direction
fixed indentation
2019-04-03 16:24:07 +05:30
Liguess = Liguess_old + steplengthLi * deltaLi
added some lines to help while debugging
2019-09-07 01:15:49 +05:30
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
write(6,'(a,1x,f7.4)') '<< CRYST integrateStress >> linear search for Liguess with step', steplengthLi
endif
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
cycle LiLoop
endif
!* calculate Jacobian for correction term
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
if (mod(jacoCounterLi, num%iJacoLpresiduum) == 0) then
fixed indentation
2019-04-03 16:24:07 +05:30
temp_33 = matmul(matmul(A,B),invFi_current)
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
do o=1,3; do p=1,3
fixed indentation
2019-04-03 16:24:07 +05:30
dFe_dLi(1:3,o,1:3,p) = -dt*math_I3(o,p)*temp_33 ! dFe_dLp(i,j,k,l) = -dt * A(i,k) invFi(l,j)
dFi_dLi(1:3,o,1:3,p) = -dt*math_I3(o,p)*invFi_current
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
enddo; enddo
do o=1,3; do p=1,3
fixed indentation
2019-04-03 16:24:07 +05:30
dFi_dLi(1:3,1:3,o,p) = matmul(matmul(Fi_new,dFi_dLi(1:3,1:3,o,p)),Fi_new)
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
enddo; enddo
fixed indentation
2019-04-03 16:24:07 +05:30
dRLi_dLi = math_identity2nd(9) &
- math_3333to99(math_mul3333xx3333(dLi_dS, math_mul3333xx3333(dS_dFe, dFe_dLi) + &
math_mul3333xx3333(dS_dFi, dFi_dLi))) &
- math_3333to99(math_mul3333xx3333(dLi_dFi, dFi_dLi))
work = math_33to9(residuumLi)
call dgesv(9,1,dRLi_dLi,9,devNull,work,9,ierr) ! solve dRLi/dLp * delta Li = -res for delta Li
if (ierr /= 0) then
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) then
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST integrateStress >> failed on dR/dLi inversion at el ip ipc ', &
el,ip,ipc
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g)&
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
write(6,*)
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST integrateStress >> dR_dLi',transpose(dRLi_dLi)
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST integrateStress >> dFe_dLi',transpose(math_3333to99(dFe_dLi))
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST integrateStress >> dS_dFi (cnst)',transpose(math_3333to99(dS_dFi))
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST integrateStress >> dLi_dS (cnst)',transpose(math_3333to99(dLi_dS))
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST integrateStress >> Li_constitutive',transpose(Li_constitutive)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST integrateStress >> Liguess',transpose(Liguess)
endif
endif
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
return
endif
deltaLi = - math_9to33(work)
endif
jacoCounterLi = jacoCounterLi + 1
Liguess = Liguess + steplengthLi * deltaLi
added some lines to help while debugging
2019-09-07 01:15:49 +05:30
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST integrateStress >> corrected Liguess by', transpose(deltaLi)
endif
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
enddo LiLoop
!* calculate new plastic and elastic deformation gradient
invFp_new = matmul(invFp_current,B)
invFp_new = invFp_new / math_det33(invFp_new)**(1.0_pReal/3.0_pReal) ! regularize
Fp_new = math_inv33(invFp_new)
failedInversionInvFp: if (all(dEq0(Fp_new))) then
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0) then
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST integrateStress >> failed on invFp_new inversion at el ip ipc ', &
el,ip,ipc
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) &
write(6,'(/,a,/,3(12x,3(f12.7,1x)/))') '<< CRYST integrateStress >> invFp_new',transpose(invFp_new)
endif
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
#endif
fixed indentation
2019-04-03 16:24:07 +05:30
return
endif failedInversionInvFp
Fe_new = matmul(matmul(Fg_new,invFp_new),invFi_new)
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
!--------------------------------------------------------------------------------------------------
! stress integration was successful
fixed indentation
2019-04-03 16:24:07 +05:30
integrateStress = .true.
crystallite_P (1:3,1:3,ipc,ip,el) = matmul(matmul(Fg_new,invFp_new),matmul(S,transpose(invFp_new)))
crystallite_S (1:3,1:3,ipc,ip,el) = S
crystallite_Lp (1:3,1:3,ipc,ip,el) = Lpguess
crystallite_Li (1:3,1:3,ipc,ip,el) = Liguess
crystallite_Fp (1:3,1:3,ipc,ip,el) = Fp_new
crystallite_Fi (1:3,1:3,ipc,ip,el) = Fi_new
crystallite_Fe (1:3,1:3,ipc,ip,el) = Fe_new
crystallite_invFp(1:3,1:3,ipc,ip,el) = invFp_new
crystallite_invFi(1:3,1:3,ipc,ip,el) = invFi_new
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
#ifdef DEBUG
fixed indentation
2019-04-03 16:24:07 +05:30
if (iand(debug_level(debug_crystallite),debug_levelExtensive) /= 0 &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
write(6,'(a,/)') '<< CRYST integrateStress >> successful integration'
write(6,'(a,/,3(12x,3(f12.7,1x)/))') '<< CRYST integrateStress >> P / MPa', &
transpose(crystallite_P(1:3,1:3,ipc,ip,el))*1.0e-6_pReal
write(6,'(a,/,3(12x,3(f12.7,1x)/))') '<< CRYST integrateStress >> Cauchy / MPa', &
matmul(crystallite_P(1:3,1:3,ipc,ip,el), transpose(Fg_new)) * 1.0e-6_pReal / math_det33(Fg_new)
write(6,'(a,/,3(12x,3(f12.7,1x)/))') '<< CRYST integrateStress >> Fe Lp Fe^-1', &
transpose(matmul(Fe_new, matmul(crystallite_Lp(1:3,1:3,ipc,ip,el), math_inv33(Fe_new))))
write(6,'(a,/,3(12x,3(f12.7,1x)/))') '<< CRYST integrateStress >> Fp',transpose(crystallite_Fp(1:3,1:3,ipc,ip,el))
write(6,'(a,/,3(12x,3(f12.7,1x)/))') '<< CRYST integrateStress >> Fi',transpose(crystallite_Fi(1:3,1:3,ipc,ip,el))
endif
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
#endif
add pheno+ module in
2015-10-14 00:22:01 +05:30
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
end function integrateStress
* now able to choose method for state integration (integrator and integratorStiffness in numerics.config); standard is "Fixed Point Iteration", which is basically the same as before; others available are "Explicit Euler", "AdaptiveEuler", "Classical Runge-Kutta" and "Runge-Kutta Cash-Karp" * now remembering stiffness similar to how we do it for Lp etc.; avoids undefined stiffness values for nonconverged stiffness calculation * non-local stuff: * changed non-local kinetics (Gilman2002) * enforce zero shearrate for overall carrrier density below relevant density * enforce zero density for those states that become negative and were below relevant density before * dislocation velocity is not limited by V^(1/3) / dt anymore
2010-10-01 17:48:49 +05:30
instantaneous jumps in the state by constitutive_deltaState are now incorporated for all state integrators. still they (should) not influence the result, since all constitutive laws simply return zero for the deltaState
2012-05-17 20:55:21 +05:30
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!> @brief integrate stress, state with adaptive 1st order explicit Euler method
!> using Fixed Point Iteration to adapt the stepsize
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
repetition not needed - implicit none at the beginning of the module is enough
2019-04-11 11:05:58 +05:30
subroutine integrateStateFPI
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
NiterationState, & !< number of iterations in state loop
e, & !< element index in element loop
i, & !< integration point index in ip loop
g, & !< grain index in grain loop
modified state integrators to perform staggered integration of plasticity damage and thermal internal variables
2014-06-25 04:51:25 +05:30
p, &
c, &
more "building block" separation
2019-01-24 03:32:21 +05:30
s, &
no need for two variables only resulted in confusing code
2019-01-29 10:09:01 +05:30
sizeDotState
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal) :: &
clearer logic
2019-01-30 02:59:36 +05:30
zeta
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal), dimension(constitutive_plasticity_maxSizeDotState) :: &
taking over from old branch
2019-01-30 02:49:38 +05:30
residuum_plastic ! residuum for plastic state
clearer logic
2019-01-30 02:59:36 +05:30
real(pReal), dimension(constitutive_source_maxSizeDotState) :: &
taking over from old branch
2019-01-30 02:49:38 +05:30
residuum_source ! residuum for source state
in preguess of FPintegrator: resetting of previous dotstates can be done in same loop as collectdotstate polishing: indentation level and capitalization
2013-11-21 16:28:41 +05:30
logical :: &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
doneWithIntegration
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --+>> PREGUESS FOR STATE <<+--
further cleaning
2019-01-28 16:19:24 +05:30
call update_dotState(1.0_pReal)
call update_state(1.0_pReal)
add pheno+ module in
2015-10-14 00:22:01 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
NiterationState = 0
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
doneWithIntegration = .false.
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
crystalliteLooping: do while (.not. doneWithIntegration .and. NiterationState < num%nState)
pInt not needed
2019-04-03 16:02:30 +05:30
NiterationState = NiterationState + 1
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
shortened too long lines
2019-02-27 02:01:47 +05:30
#ifdef DEBUG
pInt not needed
2019-04-03 16:02:30 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0) &
more elaborate debug labelling
2019-02-27 00:17:46 +05:30
write(6,'(a,i6)') '<< CRYST stateFPI >> state iteration ',NiterationState
shortened too long lines
2019-02-27 02:01:47 +05:30
#endif
more elaborate debug labelling
2019-02-27 00:17:46 +05:30
more "building block" separation
2019-01-24 03:32:21 +05:30
! store previousDotState and previousDotState2
layout adjustments
2019-01-29 12:59:19 +05:30
more "building block" separation
2019-01-24 03:32:21 +05:30
!$OMP PARALLEL DO PRIVATE(p,c)
layout adjustments
2019-01-29 12:59:19 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
layout adjustments
2019-01-29 12:59:19 +05:30
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
layout adjustments
2019-01-29 12:59:19 +05:30
plasticState(p)%previousDotState2(:,c) = merge(plasticState(p)%previousDotState(:,c),&
0.0_pReal,&
pInt not needed
2019-04-03 16:02:30 +05:30
NiterationState > 1)
layout adjustments
2019-01-29 12:59:19 +05:30
plasticState(p)%previousDotState (:,c) = plasticState(p)%dotState(:,c)
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
layout adjustments
2019-01-29 12:59:19 +05:30
sourceState(p)%p(s)%previousDotState2(:,c) = merge(sourceState(p)%p(s)%previousDotState(:,c),&
0.0_pReal, &
pInt not needed
2019-04-03 16:02:30 +05:30
NiterationState > 1)
layout adjustments
2019-01-29 12:59:19 +05:30
sourceState(p)%p(s)%previousDotState (:,c) = sourceState(p)%p(s)%dotState(:,c)
enddo
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo
more "building block" separation
2019-01-24 03:32:21 +05:30
enddo
enddo
!$OMP END PARALLEL DO
call update_dependentState
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
call update_stress(1.0_pReal)
common function to update dot state
2019-01-23 10:44:19 +05:30
call update_dotState(1.0_pReal)
layout adjustments
2019-01-29 12:59:19 +05:30
!$OMP PARALLEL
clearer logic
2019-01-30 02:59:36 +05:30
!$OMP DO PRIVATE(sizeDotState,residuum_plastic,residuum_source,zeta,p,c)
more elaborate debug labelling
2019-02-27 00:17:46 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
further cleaning
2019-01-28 16:19:24 +05:30
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
layout adjustments
2019-01-29 12:59:19 +05:30
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
now readable (kind of)
2019-01-30 04:01:26 +05:30
sizeDotState = plasticState(p)%sizeDotState
more work on the new state
2014-05-27 20:16:03 +05:30
clearer logic
2019-01-30 02:59:36 +05:30
zeta = damper(plasticState(p)%dotState (:,c), &
plasticState(p)%previousDotState (:,c), &
plasticState(p)%previousDotState2(:,c))
now readable (kind of)
2019-01-30 04:01:26 +05:30
taking over from old branch
2019-01-30 02:49:38 +05:30
residuum_plastic(1:SizeDotState) = plasticState(p)%state (1:sizeDotState,c) &
layout adjustments
2019-01-29 12:59:19 +05:30
- plasticState(p)%subState0(1:sizeDotState,c) &
clearer logic
2019-01-30 02:59:36 +05:30
- ( plasticState(p)%dotState (:,c) * zeta &
+ plasticState(p)%previousDotState(:,c) * (1.0_pReal-zeta) &
layout adjustments
2019-01-29 12:59:19 +05:30
) * crystallite_subdt(g,i,e)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
layout adjustments
2019-01-29 12:59:19 +05:30
plasticState(p)%state(1:sizeDotState,c) = plasticState(p)%state(1:sizeDotState,c) &
taking over from old branch
2019-01-30 02:49:38 +05:30
- residuum_plastic(1:sizeDotState)
clearer logic
2019-01-30 02:59:36 +05:30
plasticState(p)%dotState(:,c) = plasticState(p)%dotState(:,c) * zeta &
+ plasticState(p)%previousDotState(:,c) * (1.0_pReal - zeta)
use the same formulation for convergence every where
2019-01-31 13:42:44 +05:30
crystallite_converged(g,i,e) = converged(residuum_plastic(1:sizeDotState), &
plasticState(p)%state(1:sizeDotState,c), &
plasticState(p)%aTolState(1:sizeDotState))
convergence of plastic state can be done earlier
2019-01-29 09:50:16 +05:30
more work on the new state
2014-05-27 20:16:03 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
now readable (kind of)
2019-01-30 04:01:26 +05:30
sizeDotState = sourceState(p)%p(s)%sizeDotState
clearer logic
2019-01-30 02:59:36 +05:30
zeta = damper(sourceState(p)%p(s)%dotState (:,c), &
sourceState(p)%p(s)%previousDotState (:,c), &
sourceState(p)%p(s)%previousDotState2(:,c))
now readable (kind of)
2019-01-30 04:01:26 +05:30
clearer logic
2019-01-30 02:59:36 +05:30
residuum_source(1:sizeDotState) = sourceState(p)%p(s)%state (1:sizeDotState,c) &
- sourceState(p)%p(s)%subState0(1:sizeDotState,c) &
- ( sourceState(p)%p(s)%dotState (:,c) * zeta &
+ sourceState(p)%p(s)%previousDotState(:,c) * (1.0_pReal - zeta) &
) * crystallite_subdt(g,i,e)
sourceState(p)%p(s)%state(1:sizeDotState,c) = sourceState(p)%p(s)%state(1:sizeDotState,c) &
- residuum_source(1:sizeDotState)
sourceState(p)%p(s)%dotState(:,c) = sourceState(p)%p(s)%dotState(:,c) * zeta &
+ sourceState(p)%p(s)%previousDotState(:,c)* (1.0_pReal - zeta)
merged and added correction to FreeSurface.config example
2019-02-27 00:26:49 +05:30
crystallite_converged(g,i,e) = &
crystallite_converged(g,i,e) .and. converged(residuum_source(1:sizeDotState), &
sourceState(p)%p(s)%state(1:sizeDotState,c), &
sourceState(p)%p(s)%aTolState(1:sizeDotState))
clearer logic
2019-01-30 02:59:36 +05:30
enddo
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
endif
clearer logic
2019-01-30 02:59:36 +05:30
enddo; enddo; enddo
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP ENDDO
!$OMP DO
further cleaning
2019-01-28 16:19:24 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e) .and. crystallite_converged(g,i,e)) then ! converged and still alive...
crystallite_todo(g,i,e) = stateJump(g,i,e)
!$OMP FLUSH(crystallite_todo)
if (.not. crystallite_todo(g,i,e)) then ! if state jump fails, then convergence is broken
crystallite_converged(g,i,e) = .false.
if (.not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local...
!$OMP CRITICAL (checkTodo)
crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped
!$OMP END CRITICAL (checkTodo)
endif
endif
doxygen comments
2013-02-22 04:38:36 +05:30
endif
enddo; enddo; enddo
!$OMP ENDDO
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP END PARALLEL
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
if (any(plasticState(:)%nonlocal)) call nonlocalConvergenceCheck
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
more work on the new state
2014-05-27 20:16:03 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- CHECK IF DONE WITH INTEGRATION ---
doneWithIntegration = .true.
further cleaning
2019-01-28 16:19:24 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
doneWithIntegration = .false.
further cleaning
2019-01-28 16:19:24 +05:30
exit
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
endif
enddo; enddo
further cleaning
2019-01-28 16:19:24 +05:30
enddo
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo crystalliteLooping
further cleaning
2019-01-28 16:19:24 +05:30
cleaning
2019-01-28 15:58:46 +05:30
contains
polishing
2019-01-29 04:57:58 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief calculate the damping for correction of state and dot state
!--------------------------------------------------------------------------------------------------
cleaning
2019-01-28 15:58:46 +05:30
real(pReal) pure function damper(current,previous,previous2)
real(pReal), dimension(:), intent(in) ::&
current, previous, previous2
real(pReal) :: dot_prod12, dot_prod22
wrong dot product in state damper
2019-01-29 15:22:00 +05:30
dot_prod12 = dot_product(current - previous, previous - previous2)
dot_prod22 = dot_product(previous - previous2, previous - previous2)
sorted according to importance
2019-01-30 11:17:36 +05:30
if ((dot_product(current,previous) < 0.0_pReal .or. dot_prod12 < 0.0_pReal) .and. dot_prod22 > 0.0_pReal) then
cleaning
2019-01-28 15:58:46 +05:30
damper = 0.75_pReal + 0.25_pReal * tanh(2.0_pReal + 4.0_pReal * dot_prod12 / dot_prod22)
else
damper = 1.0_pReal
endif
end function damper
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
end subroutine integrateStateFPI
* now able to choose method for state integration (integrator and integratorStiffness in numerics.config); standard is "Fixed Point Iteration", which is basically the same as before; others available are "Explicit Euler", "AdaptiveEuler", "Classical Runge-Kutta" and "Runge-Kutta Cash-Karp" * now remembering stiffness similar to how we do it for Lp etc.; avoids undefined stiffness values for nonconverged stiffness calculation * non-local stuff: * changed non-local kinetics (Gilman2002) * enforce zero shearrate for overall carrrier density below relevant density * enforce zero density for those states that become negative and were below relevant density before * dislocation velocity is not limited by V^(1/3) / dt anymore
2010-10-01 17:48:49 +05:30
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
!> @brief integrate state with 1st order explicit Euler method
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
subroutine integrateStateEuler
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
common function to update dot state
2019-01-23 10:44:19 +05:30
call update_dotState(1.0_pReal)
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
call update_state(1.0_pReal)
delta state update for all points replaced stateJump, which works only on one point
2019-01-24 16:03:04 +05:30
call update_deltaState
call update_dependentState
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
call update_stress(1.0_pReal)
fuction for convergence check avoid code repetition
2019-01-25 11:50:05 +05:30
call setConvergenceFlag
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
if (any(plasticState(:)%nonlocal)) call nonlocalConvergenceCheck
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
end subroutine integrateStateEuler
* now able to choose method for state integration (integrator and integratorStiffness in numerics.config); standard is "Fixed Point Iteration", which is basically the same as before; others available are "Explicit Euler", "AdaptiveEuler", "Classical Runge-Kutta" and "Runge-Kutta Cash-Karp" * now remembering stiffness similar to how we do it for Lp etc.; avoids undefined stiffness values for nonconverged stiffness calculation * non-local stuff: * changed non-local kinetics (Gilman2002) * enforce zero shearrate for overall carrrier density below relevant density * enforce zero density for those states that become negative and were below relevant density before * dislocation velocity is not limited by V^(1/3) / dt anymore
2010-10-01 17:48:49 +05:30
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!> @brief integrate stress, state with 1st order Euler method with adaptive step size
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
subroutine integrateStateAdaptiveEuler
changed to new state, please report bugs to Luv or Martin
2014-07-02 17:57:39 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
e, & ! element index in element loop
i, & ! integration point index in ip loop
g, & ! grain index in grain loop
modified state integrators to perform staggered integration of plasticity damage and thermal internal variables
2014-06-25 04:51:25 +05:30
p, &
c, &
consistent variable names
2019-01-30 03:34:50 +05:30
s, &
one variable is enough
2019-01-30 03:16:21 +05:30
sizeDotState
now readable (kind of)
2019-01-30 04:01:26 +05:30
no need to store relative residual for all points
2019-01-30 14:58:47 +05:30
! ToDo: MD: once all constitutives use allocate state, attach residuum arrays to the state in case of adaptive Euler
real(pReal), dimension(constitutive_plasticity_maxSizeDotState, &
mesh_elem and theMesh are deprecated
2019-06-07 02:19:17 +05:30
homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: &
no need to store relative residual for all points
2019-01-30 14:58:47 +05:30
residuum_plastic
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal), dimension(constitutive_source_maxSizeDotState,&
maxval(phase_Nsources), &
mesh_elem and theMesh are deprecated
2019-06-07 02:19:17 +05:30
homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: &
still improving readability
2019-01-30 04:15:41 +05:30
residuum_source
more work on the new state
2014-05-27 20:16:03 +05:30
common function to update dot state
2019-01-23 10:44:19 +05:30
!--------------------------------------------------------------------------------------------------
! contribution to state and relative residui and from Euler integration
call update_dotState(1.0_pReal)
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
one variable is enough
2019-01-30 03:16:21 +05:30
!$OMP PARALLEL DO PRIVATE(sizeDotState,p,c)
now readable (kind of)
2019-01-30 04:01:26 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (crystallite_todo(g,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
one variable is enough
2019-01-30 03:16:21 +05:30
sizeDotState = plasticState(p)%sizeDotState
now readable (kind of)
2019-01-30 04:01:26 +05:30
still improving readability
2019-01-30 04:15:41 +05:30
residuum_plastic(1:sizeDotState,g,i,e) = plasticState(p)%dotstate(1:sizeDotState,c) &
* (- 0.5_pReal * crystallite_subdt(g,i,e))
shortened too long lines
2019-02-27 02:01:47 +05:30
plasticState(p)%state(1:sizeDotState,c) = &
plasticState(p)%state(1:sizeDotState,c) + plasticState(p)%dotstate(1:sizeDotState,c) * crystallite_subdt(g,i,e) !ToDo: state, partitioned state?
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
consistent variable names
2019-01-30 03:34:50 +05:30
sizeDotState = sourceState(p)%p(s)%sizeDotState
now readable (kind of)
2019-01-30 04:01:26 +05:30
still improving readability
2019-01-30 04:15:41 +05:30
residuum_source(1:sizeDotState,s,g,i,e) = sourceState(p)%p(s)%dotstate(1:sizeDotState,c) &
* (- 0.5_pReal * crystallite_subdt(g,i,e))
shortened too long lines
2019-02-27 02:01:47 +05:30
sourceState(p)%p(s)%state(1:sizeDotState,c) = &
sourceState(p)%p(s)%state(1:sizeDotState,c) + sourceState(p)%p(s)%dotstate(1:sizeDotState,c) * crystallite_subdt(g,i,e) !ToDo: state, partitioned state?
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo
endif
enddo; enddo; enddo
further simplifications
2019-01-29 05:24:02 +05:30
!$OMP END PARALLEL DO
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
no need to store relative residual for all points
2019-01-30 14:58:47 +05:30
call update_deltaState
call update_dependentState
call update_stress(1.0_pReal)
call update_dotState(1.0_pReal)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
still improving readability
2019-01-30 04:15:41 +05:30
!$OMP PARALLEL DO PRIVATE(sizeDotState,p,c)
now readable (kind of)
2019-01-30 04:01:26 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (crystallite_todo(g,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
now readable (kind of)
2019-01-30 04:01:26 +05:30
sizeDotState = plasticState(p)%sizeDotState
still improving readability
2019-01-30 04:15:41 +05:30
residuum_plastic(1:sizeDotState,g,i,e) = residuum_plastic(1:sizeDotState,g,i,e) &
+ 0.5_pReal * plasticState(p)%dotState(:,c) * crystallite_subdt(g,i,e)
use the same formulation for convergence every where
2019-01-31 13:42:44 +05:30
crystallite_converged(g,i,e) = converged(residuum_plastic(1:sizeDotState,g,i,e), &
shortened too long lines
2019-02-27 02:01:47 +05:30
plasticState(p)%state(1:sizeDotState,c), &
plasticState(p)%aTolState(1:sizeDotState))
major restructuring of multi field handling in DAMASK and added some example config files for multi field simulations. please report bugs
2015-05-28 22:32:23 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
consistent variable names
2019-01-30 03:34:50 +05:30
sizeDotState = sourceState(p)%p(s)%sizeDotState
now readable (kind of)
2019-01-30 04:01:26 +05:30
shortened too long lines
2019-02-27 02:01:47 +05:30
residuum_source(1:sizeDotState,s,g,i,e) = &
residuum_source(1:sizeDotState,s,g,i,e) + 0.5_pReal * sourceState(p)%p(s)%dotState(:,c) * crystallite_subdt(g,i,e)
improved readability
2019-01-30 17:06:02 +05:30
shortened too long lines
2019-02-27 02:01:47 +05:30
crystallite_converged(g,i,e) = &
crystallite_converged(g,i,e) .and. converged(residuum_source(1:sizeDotState,s,g,i,e), &
sourceState(p)%p(s)%state(1:sizeDotState,c), &
sourceState(p)%p(s)%aTolState(1:sizeDotState))
improved readability
2019-01-30 17:06:02 +05:30
enddo
doxygen comments
2013-02-22 04:38:36 +05:30
endif
now readable (kind of)
2019-01-30 04:01:26 +05:30
enddo; enddo; enddo
further simplifications
2019-01-29 05:24:02 +05:30
!$OMP END PARALLEL DO
crystallite_integrateStateFPI: convergence check always ran about all e,i,g; now checking only those specified in FEsolving_execElem, FEsolving_execIP
2014-01-22 00:15:41 +05:30
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
if (any(plasticState(:)%nonlocal)) call nonlocalConvergenceCheck
improved readability
2019-01-30 17:06:02 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
end subroutine integrateStateAdaptiveEuler
* now able to choose method for state integration (integrator and integratorStiffness in numerics.config); standard is "Fixed Point Iteration", which is basically the same as before; others available are "Explicit Euler", "AdaptiveEuler", "Classical Runge-Kutta" and "Runge-Kutta Cash-Karp" * now remembering stiffness similar to how we do it for Lp etc.; avoids undefined stiffness values for nonconverged stiffness calculation * non-local stuff: * changed non-local kinetics (Gilman2002) * enforce zero shearrate for overall carrrier density below relevant density * enforce zero density for those states that become negative and were below relevant density before * dislocation velocity is not limited by V^(1/3) / dt anymore
2010-10-01 17:48:49 +05:30
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!> @brief integrate stress, state with 4th order explicit Runge Kutta method
just figured out that RK4 integrator is totally broken readable code helps ;)
2019-01-30 04:47:04 +05:30
! ToDo: This is totally BROKEN: RK4dotState is never used!!!
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
subroutine integrateStateRK4
added funtions to get averged properties at integration points.
2014-09-10 14:07:12 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal), dimension(4), parameter :: &
TIMESTEPFRACTION = [0.5_pReal, 0.5_pReal, 1.0_pReal, 1.0_pReal] ! factor giving the fraction of the original timestep used for Runge Kutta Integration
real(pReal), dimension(4), parameter :: &
WEIGHT = [1.0_pReal, 2.0_pReal, 2.0_pReal, 1.0_pReal/6.0_pReal] ! weight of slope used for Runge Kutta integration (final weight divided by 6)
added funtions to get averged properties at integration points.
2014-09-10 14:07:12 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: e, & ! element index in element loop
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
i, & ! integration point index in ip loop
g, & ! grain index in grain loop
p, & ! phase loop
c, &
n, &
common variable name
2019-01-30 04:31:40 +05:30
s
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
common function to update dot state
2019-01-23 10:44:19 +05:30
call update_dotState(1.0_pReal)
did only compile with BUILD_TYPE=DEBUG
2018-05-09 20:24:06 +05:30
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
do n = 1,4
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
on-the-fly initialization
2019-01-30 04:41:10 +05:30
!$OMP PARALLEL DO PRIVATE(p,c)
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (crystallite_todo(g,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
common variable name
2019-01-30 04:31:40 +05:30
on-the-fly initialization
2019-01-30 04:41:10 +05:30
plasticState(p)%RK4dotState(:,c) = WEIGHT(n)*plasticState(p)%dotState(:,c) &
pInt not needed
2019-04-03 16:02:30 +05:30
+ merge(plasticState(p)%RK4dotState(:,c),0.0_pReal,n>1)
do s = 1, phase_Nsources(p)
on-the-fly initialization
2019-01-30 04:41:10 +05:30
sourceState(p)%p(s)%RK4dotState(:,c) = WEIGHT(n)*sourceState(p)%p(s)%dotState(:,c) &
pInt not needed
2019-04-03 16:02:30 +05:30
+ merge(sourceState(p)%p(s)%RK4dotState(:,c),0.0_pReal,n>1)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo
endif
enddo; enddo; enddo
on-the-fly initialization
2019-01-30 04:41:10 +05:30
!$OMP END PARALLEL DO
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
call update_state(TIMESTEPFRACTION(n))
delta state update for all points replaced stateJump, which works only on one point
2019-01-24 16:03:04 +05:30
call update_deltaState
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
call update_dependentState
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
call update_stress(TIMESTEPFRACTION(n))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- dot state and RK dot state---
first3steps: if (n < 4) then
on-the-fly initialization
2019-01-30 04:41:10 +05:30
call update_dotState(TIMESTEPFRACTION(n))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
endif first3steps
enddo
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
call setConvergenceFlag
if (any(plasticState(:)%nonlocal)) call nonlocalConvergenceCheck
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
end subroutine integrateStateRK4
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief integrate stress, state with 5th order Runge-Kutta Cash-Karp method with
!> adaptive step size (use 5th order solution to advance = "local extrapolation")
!--------------------------------------------------------------------------------------------------
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
subroutine integrateStateRKCK45
switched order of Lp-Li nesting in stress integration loop for better convergence. temperature is integrated during stress integration so it need not be a dot state
2014-11-12 22:10:50 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal), dimension(5,5), parameter :: &
A = reshape([&
shortened too long lines
2019-02-27 02:01:47 +05:30
.2_pReal, .075_pReal, .3_pReal, -11.0_pReal/54.0_pReal, 1631.0_pReal/55296.0_pReal, &
.0_pReal, .225_pReal, -.9_pReal, 2.5_pReal, 175.0_pReal/512.0_pReal, &
.0_pReal, .0_pReal, 1.2_pReal, -70.0_pReal/27.0_pReal, 575.0_pReal/13824.0_pReal, &
.0_pReal, .0_pReal, .0_pReal, 35.0_pReal/27.0_pReal, 44275.0_pReal/110592.0_pReal, &
.0_pReal, .0_pReal, .0_pReal, .0_pReal, 253.0_pReal/4096.0_pReal], &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
[5,5], order=[2,1]) !< coefficients in Butcher tableau (used for preliminary integration in stages 2 to 6)
switched order of Lp-Li nesting in stress integration loop for better convergence. temperature is integrated during stress integration so it need not be a dot state
2014-11-12 22:10:50 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal), dimension(6), parameter :: &
B = &
[37.0_pReal/378.0_pReal, .0_pReal, 250.0_pReal/621.0_pReal, &
125.0_pReal/594.0_pReal, .0_pReal, 512.0_pReal/1771.0_pReal], & !< coefficients in Butcher tableau (used for final integration and error estimate)
DB = B - &
shortened too long lines
2019-02-27 02:01:47 +05:30
[2825.0_pReal/27648.0_pReal, .0_pReal, 18575.0_pReal/48384.0_pReal,&
13525.0_pReal/55296.0_pReal, 277.0_pReal/14336.0_pReal, 0.25_pReal] !< coefficients in Butcher tableau (used for final integration and error estimate)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal), dimension(5), parameter :: &
C = [0.2_pReal, 0.3_pReal, 0.6_pReal, 1.0_pReal, 0.875_pReal] !< coefficients in Butcher tableau (fractions of original time step in stages 2 to 6)
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
e, & ! element index in element loop
i, & ! integration point index in ip loop
g, & ! grain index in grain loop
stage, & ! stage index in integration stage loop
n, &
p, &
cc, &
same names everywhere if possible
2019-01-30 15:07:18 +05:30
s, &
only one variable needed
2019-01-30 15:18:59 +05:30
sizeDotState
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
no need to store relative residual pointwise
2019-01-30 15:34:49 +05:30
! ToDo: MD: once all constitutives use allocate state, attach residuum arrays to the state in case of RKCK45
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal), dimension(constitutive_plasticity_maxSizeDotState, &
mesh_elem and theMesh are deprecated
2019-06-07 02:19:17 +05:30
homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: &
no need to store relative residual pointwise
2019-01-30 15:34:49 +05:30
residuum_plastic ! relative residuum from evolution in microstructure
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal), dimension(constitutive_source_maxSizeDotState, &
maxval(phase_Nsources), &
mesh_elem and theMesh are deprecated
2019-06-07 02:19:17 +05:30
homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: &
no need to store relative residual pointwise
2019-01-30 15:34:49 +05:30
residuum_source ! relative residuum from evolution in microstructure
overhaul of kinematics from FeFp to FeFiFp decomposition
2014-11-01 00:33:08 +05:30
switched order of Lp-Li nesting in stress integration loop for better convergence. temperature is integrated during stress integration so it need not be a dot state
2014-11-12 22:10:50 +05:30
common function to update dot state
2019-01-23 10:44:19 +05:30
call update_dotState(1.0_pReal)
switched order of Lp-Li nesting in stress integration loop for better convergence. temperature is integrated during stress integration so it need not be a dot state
2014-11-12 22:10:50 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- SECOND TO SIXTH RUNGE KUTTA STEP ---
overhaul of kinematics from FeFp to FeFiFp decomposition
2014-11-01 00:33:08 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
do stage = 1,5
switched order of Lp-Li nesting in stress integration loop for better convergence. temperature is integrated during stress integration so it need not be a dot state
2014-11-12 22:10:50 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- state update ---
only one loop needed
2019-01-30 13:26:16 +05:30
!$OMP PARALLEL DO PRIVATE(p,cc)
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (crystallite_todo(g,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); cc = material_phaseMemberAt(g,i,e)
switched order of Lp-Li nesting in stress integration loop for better convergence. temperature is integrated during stress integration so it need not be a dot state
2014-11-12 22:10:50 +05:30
only one loop needed
2019-01-30 13:26:16 +05:30
plasticState(p)%RKCK45dotState(stage,:,cc) = plasticState(p)%dotState(:,cc)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
plasticState(p)%dotState(:,cc) = A(1,stage) * plasticState(p)%RKCK45dotState(1,:,cc)
only one loop needed
2019-01-30 13:26:16 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
same names everywhere if possible
2019-01-30 15:07:18 +05:30
sourceState(p)%p(s)%RKCK45dotState(stage,:,cc) = sourceState(p)%p(s)%dotState(:,cc)
sourceState(p)%p(s)%dotState(:,cc) = A(1,stage) * sourceState(p)%p(s)%RKCK45dotState(1,:,cc)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo
only one loop needed
2019-01-30 13:26:16 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
do n = 2, stage
only one loop needed
2019-01-30 13:26:16 +05:30
plasticState(p)%dotState(:,cc) = plasticState(p)%dotState(:,cc) &
+ A(n,stage) * plasticState(p)%RKCK45dotState(n,:,cc)
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
same names everywhere if possible
2019-01-30 15:07:18 +05:30
sourceState(p)%p(s)%dotState(:,cc) = sourceState(p)%p(s)%dotState(:,cc) &
avoid loops
2019-01-30 15:16:53 +05:30
+ A(n,stage) * sourceState(p)%p(s)%RKCK45dotState(n,:,cc)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo
enddo
only one loop needed
2019-01-30 13:26:16 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
endif
enddo; enddo; enddo
only one loop needed
2019-01-30 13:26:16 +05:30
!$OMP END PARALLEL DO
overhaul of kinematics from FeFp to FeFiFp decomposition
2014-11-01 00:33:08 +05:30
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
call update_state(1.0_pReal) !MD: 1.0 correct?
delta state update for all points replaced stateJump, which works only on one point
2019-01-24 16:03:04 +05:30
call update_deltaState
call update_dependentState
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
call update_stress(C(stage))
call update_dotState(C(stage))
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
changes related to intermediate configuration kinematics: - switched Fi and Li from state variables to crystallite variables - Lp and Li are now work conjugate with the corresponding mandel stresses defined in their respective configuration - T, Lp and Li need to return tangents wrt Fi arising from the convection of the material frame due to Fi - Updated analytic jacobian to take into account tangents wrt Fi - Updated Lp and Li residual jacobians to take into account tangents wrt Fi
2015-03-06 18:39:00 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!--------------------------------------------------------------------------------------------------
! --- STATE UPDATE WITH ERROR ESTIMATE FOR STATE ---
only one variable needed
2019-01-30 15:18:59 +05:30
!$OMP PARALLEL DO PRIVATE(sizeDotState,p,cc)
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (crystallite_todo(g,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); cc = material_phaseMemberAt(g,i,e)
no need for 2 loops
2019-01-30 13:51:33 +05:30
only one variable needed
2019-01-30 15:18:59 +05:30
sizeDotState = plasticState(p)%sizeDotState
no need for 2 loops
2019-01-30 13:51:33 +05:30
plasticState(p)%RKCK45dotState(6,:,cc) = plasticState (p)%dotState(:,cc)
unifying name
2019-01-30 15:21:24 +05:30
residuum_plastic(1:sizeDotState,g,i,e) = &
just polishing
2019-01-31 09:47:46 +05:30
matmul(transpose(plasticState(p)%RKCK45dotState(1:6,1:sizeDotState,cc)),DB) & ! why transpose? Better to transpose constant DB
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
* crystallite_subdt(g,i,e)
no need for 2 loops
2019-01-30 13:51:33 +05:30
plasticState(p)%dotState(:,cc) = &
just polishing
2019-01-31 09:47:46 +05:30
matmul(transpose(plasticState(p)%RKCK45dotState(1:6,1:sizeDotState,cc)), B) ! why transpose? Better to transpose constant B
no need for 2 loops
2019-01-30 13:51:33 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
only one variable needed
2019-01-30 15:18:59 +05:30
sizeDotState = sourceState(p)%p(s)%sizeDotState
no need for 2 loops
2019-01-30 13:51:33 +05:30
[skip ci] cleaning
2019-01-30 19:22:12 +05:30
sourceState(p)%p(s)%RKCK45dotState(6,:,cc) = sourceState(p)%p(s)%dotState(:,cc)
unifying name
2019-01-30 15:21:24 +05:30
residuum_source(1:sizeDotState,s,g,i,e) = &
only one variable needed
2019-01-30 15:18:59 +05:30
matmul(transpose(sourceState(p)%p(s)%RKCK45dotState(1:6,1:sizeDotState,cc)),DB) &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
* crystallite_subdt(g,i,e)
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
same names everywhere if possible
2019-01-30 15:07:18 +05:30
sourceState(p)%p(s)%dotState(:,cc) = &
only one variable needed
2019-01-30 15:18:59 +05:30
matmul(transpose(sourceState(p)%p(s)%RKCK45dotState(1:6,1:sizeDotState,cc)),B)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo
no need for 2 loops
2019-01-30 13:51:33 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
endif
enddo; enddo; enddo
no need for 2 loops
2019-01-30 13:51:33 +05:30
!$OMP END PARALLEL DO
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
call update_state(1.0_pReal)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- relative residui and state convergence ---
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
no need to store relative residual pointwise
2019-01-30 15:34:49 +05:30
!$OMP PARALLEL DO PRIVATE(sizeDotState,p,cc)
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
no need to store relative residual pointwise
2019-01-30 15:34:49 +05:30
if (crystallite_todo(g,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); cc = material_phaseMemberAt(g,i,e)
same names everywhere if possible
2019-01-30 15:07:18 +05:30
no need to store relative residual pointwise
2019-01-30 15:34:49 +05:30
sizeDotState = plasticState(p)%sizeDotState
improved readability
2019-01-30 17:06:02 +05:30
use the same formulation for convergence every where
2019-01-31 13:42:44 +05:30
crystallite_todo(g,i,e) = converged(residuum_plastic(1:sizeDotState,g,i,e), &
calculated convergence criteria wrongly
2019-01-30 21:34:58 +05:30
plasticState(p)%state(1:sizeDotState,cc), &
use the same formulation for convergence every where
2019-01-31 13:42:44 +05:30
plasticState(p)%aTolState(1:sizeDotState))
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
no need to store relative residual pointwise
2019-01-30 15:34:49 +05:30
sizeDotState = sourceState(p)%p(s)%sizeDotState
shortened too long lines
2019-02-27 02:01:47 +05:30
crystallite_todo(g,i,e) = &
crystallite_todo(g,i,e) .and. converged(residuum_source(1:sizeDotState,s,g,i,e), &
sourceState(p)%p(s)%state(1:sizeDotState,cc), &
sourceState(p)%p(s)%aTolState(1:sizeDotState))
[skip ci] cleaning
2019-01-30 19:22:12 +05:30
enddo
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
endif
enddo; enddo; enddo
no need to store relative residual pointwise
2019-01-30 15:34:49 +05:30
!$OMP END PARALLEL DO
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
delta state update for all points replaced stateJump, which works only on one point
2019-01-24 16:03:04 +05:30
call update_deltaState
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
call update_dependentState
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
call update_stress(1.0_pReal)
fuction for convergence check avoid code repetition
2019-01-25 11:50:05 +05:30
call setConvergenceFlag
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
if (any(plasticState(:)%nonlocal)) call nonlocalConvergenceCheck
improved readability
2019-01-30 17:06:02 +05:30
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
end subroutine integrateStateRKCK45
new substructure with 'include' and 'config' directories renaming of mpie_... to DAMASK for main usersubroutines extension of element outputs from 5 to 8 digits (FFT!!!)
2011-05-11 22:08:45 +05:30
change the debug level of cutback reporting to non-selective, of central solution from basic to extensive, of stiffness perturbation to selective. fixed few typos and deleted excess whitespace.
2014-08-26 20:14:32 +05:30
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief sets convergence flag for nonlocal calculations
!> @detail one non-converged nonlocal sets all other nonlocals to non-converged to trigger cut back
!--------------------------------------------------------------------------------------------------
repetition not needed - implicit none at the beginning of the module is enough
2019-04-11 11:05:58 +05:30
subroutine nonlocalConvergenceCheck
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
if (any(.not. crystallite_converged .and. .not. crystallite_localPlasticity)) & ! any non-local not yet converged (or broken)...
where( .not. crystallite_localPlasticity) crystallite_converged = .false.
modified state integrators to perform staggered integration of plasticity damage and thermal internal variables
2014-06-25 04:51:25 +05:30
nonlocal convergence check in function
2019-01-29 09:41:29 +05:30
end subroutine nonlocalConvergenceCheck
major restructuring of code. homogenization as well as constitutive are now free to choose. the runtime got somewhat longer (25% on simple tests) compared to a hardcoded isostrain homogenization. this might be a point of further optimization at a later stage... please use homogenization_isostrain.f90 as starting point / example for future developments of homog-schemes. the homogenization scheme now can additionally output certain results. hence, the userdata structure at each integration point now looks like this: - sizeHomogPostResults - block of that size containing homogPostResults then for each grain: - sizeGrainPostResults - block of that size containing crystallitePostResults, which consist of: + phaseID + volFrac + Eulers (3) + any constitutive post results requested
2009-05-07 21:57:36 +05:30
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
fuction for convergence check avoid code repetition
2019-01-25 11:50:05 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief Sets convergence flag based on "todo": every point that survived the integration (todo is
! still .true. is considered as converged
!> @details: For explicitEuler, RK4 and RKCK45, adaptive Euler and FPI have their on criteria
!--------------------------------------------------------------------------------------------------
repetition not needed - implicit none at the beginning of the module is enough
2019-04-11 11:05:58 +05:30
subroutine setConvergenceFlag
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: &
fuction for convergence check avoid code repetition
2019-01-25 11:50:05 +05:30
e, & !< element index in element loop
i, & !< integration point index in ip loop
g !< grain index in grain loop
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
!OMP DO PARALLEL PRIVATE
fuction for convergence check avoid code repetition
2019-01-25 11:50:05 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
forall is deprecated - performance-wise, there should be no difference for the small loops we have - still, the on-liner syntax was much nicer
2019-04-03 16:47:21 +05:30
crystallite_converged(g,i,e) = crystallite_todo(g,i,e) .or. crystallite_converged(g,i,e) ! if still "to do" then converged per definition
enddo; enddo; enddo
fuction for convergence check avoid code repetition
2019-01-25 11:50:05 +05:30
!OMP END DO PARALLEL
end subroutine setConvergenceFlag
no need to repeat the same code
2019-01-31 13:44:02 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief determines whether a point is converged
!--------------------------------------------------------------------------------------------------
logical pure function converged(residuum,state,aTol)
real(pReal), intent(in), dimension(:) ::&
residuum, state, aTol
keep connected data together - avoids dependencies - easier to read and modify
2019-04-11 10:54:04 +05:30
real(pReal) :: &
rTol
rTol = num%rTol_crystalliteState
no need to repeat the same code
2019-01-31 13:44:02 +05:30
converged = all(abs(residuum) <= max(aTol, rTol*abs(state)))
end function converged
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
!--------------------------------------------------------------------------------------------------
[skip ci] only corrections to comments made
2019-10-11 18:51:29 +05:30
!> @brief Standard forwarding of state as state = state0 + dotState * (delta t) comment seems wrong!
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
!--------------------------------------------------------------------------------------------------
subroutine update_stress(timeFraction)
repetition not needed - implicit none at the beginning of the module is enough
2019-04-11 11:05:58 +05:30
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
real(pReal), intent(in) :: &
timeFraction
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: &
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
e, & !< element index in element loop
i, & !< integration point index in ip loop
g
!$OMP PARALLEL DO
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
crystallite_todo(g,i,e) = integrateStress(g,i,e,timeFraction)
!$OMP FLUSH(crystallite_todo)
if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local...
!$OMP CRITICAL (checkTodo)
crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped
!$OMP END CRITICAL (checkTodo)
endif
endif
enddo; enddo; enddo
!$OMP END PARALLEL DO
end subroutine update_stress
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief tbd
!--------------------------------------------------------------------------------------------------
key words don't contain small letters moved more variables from numerics to crystallite
2019-04-11 14:57:03 +05:30
subroutine update_dependentState
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
use constitutive, only: &
constitutive_dependentState => constitutive_microstructure
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: e, & ! element index in element loop
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
i, & ! integration point index in ip loop
g ! grain index in grain loop
variables were not used
2019-01-24 03:48:14 +05:30
!$OMP PARALLEL DO
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) &
using new rotation class
2019-02-02 00:58:07 +05:30
call constitutive_dependentState(crystallite_Fe(1:3,1:3,g,i,e), &
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
crystallite_Fp(1:3,1:3,g,i,e), &
g, i, e)
enddo; enddo; enddo
variables were not used
2019-01-24 03:48:14 +05:30
!$OMP END PARALLEL DO
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
end subroutine update_dependentState
stress integration for all points in one function
2019-01-24 22:29:38 +05:30
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief Standard forwarding of state as state = state0 + dotState * (delta t)
!--------------------------------------------------------------------------------------------------
subroutine update_state(timeFraction)
real(pReal), intent(in) :: &
timeFraction
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: &
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
e, & !< element index in element loop
i, & !< integration point index in ip loop
g, & !< grain index in grain loop
p, &
c, &
s, &
mySize
!$OMP PARALLEL DO PRIVATE(mySize,p,c)
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
mySize = plasticState(p)%sizeDotState
plasticState(p)%state(1:mySize,c) = plasticState(p)%subState0(1:mySize,c) &
+ plasticState(p)%dotState (1:mySize,c) &
* crystallite_subdt(g,i,e) * timeFraction
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
mySize = sourceState(p)%p(s)%sizeDotState
sourceState(p)%p(s)%state(1:mySize,c) = sourceState(p)%p(s)%subState0(1:mySize,c) &
+ sourceState(p)%p(s)%dotState (1:mySize,c) &
* crystallite_subdt(g,i,e) * timeFraction
enddo
endif
enddo; enddo; enddo
!$OMP END PARALLEL DO
end subroutine update_state
common function to update dot state
2019-01-23 10:44:19 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief triggers calculation of all new rates
!> if NaN occurs, crystallite_todo is set to FALSE. Any NaN in a nonlocal propagates to all others
!--------------------------------------------------------------------------------------------------
subroutine update_dotState(timeFraction)
real(pReal), intent(in) :: &
timeFraction
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: &
common function to update dot state
2019-01-23 10:44:19 +05:30
e, & !< element index in element loop
i, & !< integration point index in ip loop
g, & !< grain index in grain loop
p, &
c, &
s
logical :: &
avoid flush of full array more clear logic
2019-01-24 18:45:26 +05:30
NaN, &
nonlocalStop
nonlocalStop = .false.
common function to update dot state
2019-01-23 10:44:19 +05:30
avoid flush of full array more clear logic
2019-01-24 18:45:26 +05:30
!$OMP PARALLEL DO PRIVATE (p,c,NaN)
common function to update dot state
2019-01-23 10:44:19 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
avoid flush of full array more clear logic
2019-01-24 18:45:26 +05:30
!$OMP FLUSH(nonlocalStop)
logic error for nonlocal
2019-01-30 20:36:14 +05:30
if ((crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) .and. .not. nonlocalStop) then
Tstar renamed to S (following the DAMASK paper)
2019-03-09 05:03:11 +05:30
call constitutive_collectDotState(crystallite_S(1:3,1:3,g,i,e), &
common function to update dot state
2019-01-23 10:44:19 +05:30
crystallite_Fe, &
crystallite_Fi(1:3,1:3,g,i,e), &
crystallite_Fp, &
Merge branch 'development' into New-Thermal
2019-02-26 12:24:03 +05:30
crystallite_subdt(g,i,e)*timeFraction, g,i,e)
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
common function to update dot state
2019-01-23 10:44:19 +05:30
NaN = any(IEEE_is_NaN(plasticState(p)%dotState(:,c)))
pInt not needed
2019-04-03 16:02:30 +05:30
do s = 1, phase_Nsources(p)
common function to update dot state
2019-01-23 10:44:19 +05:30
NaN = NaN .or. any(IEEE_is_NaN(sourceState(p)%p(s)%dotState(:,c)))
enddo
if (NaN) then
crystallite_todo(g,i,e) = .false. ! this one done (and broken)
avoid flush of full array more clear logic
2019-01-24 18:45:26 +05:30
if (.not. crystallite_localPlasticity(g,i,e)) nonlocalStop = .True.
common function to update dot state
2019-01-23 10:44:19 +05:30
endif
endif
enddo; enddo; enddo
avoid flush of full array more clear logic
2019-01-24 18:45:26 +05:30
!$OMP END PARALLEL DO
if (nonlocalStop) crystallite_todo = crystallite_todo .and. crystallite_localPlasticity
common function to update dot state
2019-01-23 10:44:19 +05:30
end subroutine update_DotState
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
subroutine update_deltaState
don't clutter the code with useless stuff we only need to be more strict about prefixing functions/subroutines/variables to see in which module they reside
2019-05-17 02:26:48 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: &
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
e, & !< element index in element loop
i, & !< integration point index in ip loop
g, & !< grain index in grain loop
p, &
on-the-fly initialization
2019-01-30 04:41:10 +05:30
mySize, &
WIP: update_deltaState
2019-01-24 12:04:30 +05:30
myOffset, &
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
c, &
s
avoid flush of full array more clear logic
2019-01-24 18:45:26 +05:30
logical :: &
NaN, &
nonlocalStop
nonlocalStop = .false.
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
on-the-fly initialization
2019-01-30 04:41:10 +05:30
!$OMP PARALLEL DO PRIVATE(p,c,myOffset,mySize,NaN)
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
mesh_element is deprecated (meaningless name)
2019-06-05 13:35:59 +05:30
do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
avoid flush of full array more clear logic
2019-01-24 18:45:26 +05:30
!$OMP FLUSH(nonlocalStop)
logic error for nonlocal
2019-01-30 20:36:14 +05:30
if ((crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) .and. .not. nonlocalStop) then
Tstar renamed to S (following the DAMASK paper)
2019-03-09 05:03:11 +05:30
call constitutive_collectDeltaState(crystallite_S(1:3,1:3,g,i,e), &
crystallite_Fe(1:3,1:3,g,i,e), &
crystallite_Fi(1:3,1:3,g,i,e), &
g,i,e)
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
WIP: update_deltaState
2019-01-24 12:04:30 +05:30
myOffset = plasticState(p)%offsetDeltaState
mySize = plasticState(p)%sizeDeltaState
NaN = any(IEEE_is_NaN(plasticState(p)%deltaState(1:mySize,c)))
if (.not. NaN) then
pInt not needed
2019-04-03 16:02:30 +05:30
plasticState(p)%state(myOffset + 1: myOffset + mySize,c) = &
plasticState(p)%state(myOffset + 1: myOffset + mySize,c) + plasticState(p)%deltaState(1:mySize,c)
do s = 1, phase_Nsources(p)
on-the-fly initialization
2019-01-30 04:41:10 +05:30
myOffset = sourceState(p)%p(s)%offsetDeltaState
mySize = sourceState(p)%p(s)%sizeDeltaState
NaN = NaN .or. any(IEEE_is_NaN(sourceState(p)%p(s)%deltaState(1:mySize,c)))
delta state update for all points replaced stateJump, which works only on one point
2019-01-24 16:03:04 +05:30
if (.not. NaN) then
pInt not needed
2019-04-03 16:02:30 +05:30
sourceState(p)%p(s)%state(myOffset + 1:myOffset + mySize,c) = &
sourceState(p)%p(s)%state(myOffset + 1:myOffset + mySize,c) + sourceState(p)%p(s)%deltaState(1:mySize,c)
delta state update for all points replaced stateJump, which works only on one point
2019-01-24 16:03:04 +05:30
endif
enddo
WIP: update_deltaState
2019-01-24 12:04:30 +05:30
endif
delta state update for all points replaced stateJump, which works only on one point
2019-01-24 16:03:04 +05:30
crystallite_todo(g,i,e) = .not. NaN
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
if (.not. crystallite_todo(g,i,e)) then ! if state jump fails, then convergence is broken
crystallite_converged(g,i,e) = .false.
avoid flush of full array more clear logic
2019-01-24 18:45:26 +05:30
if (.not. crystallite_localPlasticity(g,i,e)) nonlocalStop = .true.
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
endif
endif
enddo; enddo; enddo
!$OMP END PARALLEL DO
avoid flush of full array more clear logic
2019-01-24 18:45:26 +05:30
if (nonlocalStop) crystallite_todo = crystallite_todo .and. crystallite_localPlasticity
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
end subroutine update_deltaState
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief calculates a jump in the state according to the current state and the current stress
!> returns true, if state jump was successfull or not needed. false indicates NaN in delta state
!--------------------------------------------------------------------------------------------------
logical function stateJump(ipc,ip,el)
pInt not needed
2019-04-03 16:02:30 +05:30
integer, intent(in):: &
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
el, & ! element index
ip, & ! integration point index
ipc ! grain index
pInt not needed
2019-04-03 16:02:30 +05:30
integer :: &
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
c, &
p, &
mySource, &
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
myOffset, &
mySize
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
clearer name
2019-06-14 12:12:12 +05:30
c = material_phaseMemberAt(ipc,ip,el)
using new names
2019-06-14 12:19:05 +05:30
p = material_phaseAt(ipc,el)
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
Tstar renamed to S (following the DAMASK paper)
2019-03-09 05:03:11 +05:30
call constitutive_collectDeltaState(crystallite_S(1:3,1:3,ipc,ip,el), &
crystallite_Fe(1:3,1:3,ipc,ip,el), &
crystallite_Fi(1:3,1:3,ipc,ip,el), &
ipc,ip,el)
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
myOffset = plasticState(p)%offsetDeltaState
mySize = plasticState(p)%sizeDeltaState
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
if( any(IEEE_is_NaN(plasticState(p)%deltaState(1:mySize,c)))) then ! NaN occured in deltaState
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
stateJump = .false.
return
endif
pInt not needed
2019-04-03 16:02:30 +05:30
plasticState(p)%state(myOffset + 1:myOffset + mySize,c) = &
plasticState(p)%state(myOffset + 1:myOffset + mySize,c) + plasticState(p)%deltaState(1:mySize,c)
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
pInt not needed
2019-04-03 16:02:30 +05:30
do mySource = 1, phase_Nsources(p)
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
myOffset = sourceState(p)%p(mySource)%offsetDeltaState
fuction for convergence check avoid code repetition
2019-01-25 11:50:05 +05:30
mySize = sourceState(p)%p(mySource)%sizeDeltaState
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
if (any(IEEE_is_NaN(sourceState(p)%p(mySource)%deltaState(1:mySize,c)))) then ! NaN occured in deltaState
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
stateJump = .false.
return
endif
pInt not needed
2019-04-03 16:02:30 +05:30
sourceState(p)%p(mySource)%state(myOffset + 1: myOffset + mySize,c) = &
sourceState(p)%p(mySource)%state(myOffset + 1: myOffset + mySize,c) + sourceState(p)%p(mySource)%deltaState(1:mySize,c)
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
enddo
#ifdef DEBUG
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
if (any(dNeq0(plasticState(p)%deltaState(1:mySize,c))) &
pInt not needed
2019-04-03 16:02:30 +05:30
.and. iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0 &
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
pInt not needed
2019-04-03 16:02:30 +05:30
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0)) then
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
write(6,'(a,i8,1x,i2,1x,i3, /)') '<< CRYST >> update state at el ip ipc ',el,ip,ipc
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> deltaState', plasticState(p)%deltaState(1:mySize,c)
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> new state', &
pInt not needed
2019-04-03 16:02:30 +05:30
plasticState(p)%state(myOffset + 1 : &
preparing function for deltaState essential a "stateJump" over all points
2019-01-24 11:42:20 +05:30
myOffset + mySize,c)
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
endif
#endif
stateJump = .true.
end function stateJump
doxygen comments
2013-02-22 04:38:36 +05:30
end module crystallite