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

3048 lines
153 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
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
use FEsolving, only: &
FEsolving_execElem, &
FEsolving_execIP
use mesh, only: &
mesh_element
use material, only: &
homogenization_Ngrains
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
use prec, only: &
pReal, &
pInt
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
implicit none
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
doxygen comments
2013-02-22 04:38:36 +05:30
private
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
character(len=64), dimension(:,:), allocatable, private :: &
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_output !< name of each post result output
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
integer(pInt), public, protected :: &
introduced plain mode (no ping pong) again and added test for MSC.MArc 2012 as prove that its working
2013-05-17 23:22:46 +05:30
crystallite_maxSizePostResults !< description not available
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
integer(pInt), dimension(:), allocatable, public, protected :: &
introduced plain mode (no ping pong) again and added test for MSC.MArc 2012 as prove that its working
2013-05-17 23:22:46 +05:30
crystallite_sizePostResults !< description not available
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
integer(pInt), dimension(:,:), allocatable, private :: &
introduced plain mode (no ping pong) again and added test for MSC.MArc 2012 as prove that its working
2013-05-17 23:22:46 +05:30
crystallite_sizePostResult !< description not available
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
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
real(pReal), dimension(:,:,:), allocatable, public :: &
removed temperature integration and corresponding data structures and debugging options temperature is stored in crystallite, but homogeneous on one IP (not an component (grain) quantity and an input value parsed in by the BVP solver. introduced heat, a component (grain) quantity which is homogenized before returned to the heat transfer solver. went ahead with removal of dummy functions in homogenization and constitutive, this time mainly reduced function signatures to reflect actually needed quantities.
2013-10-16 18:34:59 +05:30
crystallite_dt !< requested time increment of each grain
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
real(pReal), dimension(:,:,:), allocatable, private :: &
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_subdt, & !< substepped time increment of each grain
crystallite_subFrac, & !< already calculated fraction of increment
removed temperature integration and corresponding data structures and debugging options temperature is stored in crystallite, but homogeneous on one IP (not an component (grain) quantity and an input value parsed in by the BVP solver. introduced heat, a component (grain) quantity which is homogenized before returned to the heat transfer solver. went ahead with removal of dummy functions in homogenization and constitutive, this time mainly reduced function signatures to reflect actually needed quantities.
2013-10-16 18:34:59 +05:30
crystallite_subStep !< size of next integration step
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
real(pReal), dimension(:,:,:,:), allocatable, public :: &
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
crystallite_Tstar_v, & !< current 2nd Piola-Kirchhoff stress vector (end of converged time step) ToDo: Should be called S, 3x3
crystallite_Tstar0_v, & !< 2nd Piola-Kirchhoff stress vector at start of FE inc ToDo: Should be called S, 3x3
crystallite_partionedTstar0_v !< 2nd Piola-Kirchhoff stress vector at start of homog inc ToDo: Should be called S, 3x3
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
real(pReal), dimension(:,:,:,:), allocatable, private :: &
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_orientation, & !< orientation as quaternion
crystallite_orientation0, & !< initial orientation as quaternion
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
crystallite_rotation !< grain rotation away from initial orientation as axis-angle (in degrees) in crystal reference frame
only read access
2018-09-19 20:34:12 +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
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
real(pReal), dimension(:,:,:,:,:), allocatable, public :: &
doxygen comments
2013-02-22 04:38:36 +05:30
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
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
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
doxygen comments
2013-02-22 04:38:36 +05:30
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
only read access
2018-09-19 20:34:12 +05:30
crystallite_partionedLp0, & !< plastic velocity grad at start of homog inc
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
crystallite_Li, & !< current intermediate velocitiy grad (end of converged time step)
crystallite_Li0, & !< intermediate velocitiy grad at start of FE inc
only read access
2018-09-19 20:34:12 +05:30
crystallite_partionedLi0 !< intermediate velocity grad at start of homog inc
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
real(pReal), dimension(:,:,:,:,:), allocatable, private :: &
further separation still using old combined function
2019-01-18 19:12:44 +05:30
crystallite_subS0, & !< 2nd Piola-Kirchhoff stress vector at start of crystallite inc
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_invFp, & !< inverse of current plastic def grad (end of converged time step)
crystallite_subFp0,& !< plastic def grad at start of crystallite inc
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
crystallite_invFi, & !< inverse of current intermediate def grad (end of converged time step)
crystallite_subFi0,& !< intermediate def grad at start of crystallite inc
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
crystallite_subF, & !< def grad to be reached at end of crystallite inc
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_subF0, & !< def grad at start of crystallite inc
crystallite_subLp0,& !< plastic velocity grad at start of crystallite inc
disorientation was never use not even for nonlocal, but it slows down calculation a lot
2019-01-14 12:14:36 +05:30
crystallite_subLi0 !< intermediate velocity grad at start of crystallite inc
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
real(pReal), dimension(:,:,:,:,:,:,:), allocatable, public :: &
no need to store converged tangent
2019-01-15 12:34:50 +05:30
crystallite_dPdF !< current individual dPdF per grain (end of converged time step)
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
logical, dimension(:,:,:), allocatable, public :: &
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_requested !< flag to request crystallite calculation
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
logical, dimension(:,:,:), allocatable, public, protected :: &
also not needed was a remainder from time syncinc
2019-01-15 12:22:01 +05:30
crystallite_converged !< convergence flag
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
logical, dimension(:,:,:), allocatable, private :: &
also not needed was a remainder from time syncinc
2019-01-15 12:22:01 +05:30
crystallite_localPlasticity, & !< indicates this grain to have purely local constitutive law
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_todo !< flag to indicate need for further computation
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
enum, bind(c)
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
enumerator :: undefined_ID, &
phase_ID, &
texture_ID, &
volume_ID, &
orientation_ID, &
grainrotation_ID, &
eulerangles_ID, &
defgrad_ID, &
fe_ID, &
fp_ID, &
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
fi_ID, &
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
lp_ID, &
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
li_ID, &
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
p_ID, &
s_ID, &
elasmatrix_ID, &
neighboringip_ID, &
neighboringelement_ID
end enum
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
integer(kind(undefined_ID)),dimension(:,:), allocatable, private :: &
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
crystallite_outputID !< ID of each post result output
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
procedure(), pointer :: integrateState
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
doxygen comments
2013-02-22 04:38:36 +05:30
public :: &
crystallite_init, &
using separate functions for stress and tangent extensively tested in 46-simplification-of-crystallite-f90-NEW3 already
2019-01-18 20:00:50 +05:30
crystallite_stress, &
crystallite_stressTangent, &
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_orientations, &
added funtions to get averged properties at integration points.
2014-09-10 14:07:12 +05:30
crystallite_push33ToRef, &
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_postResults
private :: &
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
integrateStress, &
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
integrateState, &
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
integrateStateFPI, &
integrateStateEuler, &
integrateStateAdaptiveEuler, &
integrateStateRK4, &
integrateStateRKCK45, &
stateJump
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
compiler_options and compiler_version supported by new Intel compiler
2018-02-02 17:06:09 +05:30
#if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
suppress warning when compiling with gfortran
2017-10-05 20:05:34 +05:30
use, intrinsic :: iso_fortran_env, only: &
compiler_version, &
compiler_options
#endif
define debug variables only if needed
2018-09-20 09:39:02 +05:30
#ifdef DEBUG
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
use debug, only: &
debug_info, &
debug_reset, &
debug_level, &
debug_crystallite, &
debug_levelBasic
define debug variables only if needed
2018-09-20 09:39:02 +05:30
#endif
further improvements on implementing the plain mode
2013-04-16 22:37:27 +05:30
use numerics, only: &
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
numerics_integrator, &
only output from the root processor for parallel runs
2014-10-10 01:53:06 +05:30
worldrank, &
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
usePingPong
introduced plain mode (no ping pong) again and added test for MSC.MArc 2012 as prove that its working
2013-05-17 23:22:46 +05:30
use math, only: &
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
math_I3, &
math_EulerToR, &
math_inv33, &
math_mul33x33
use mesh, only: &
mesh_element, &
mesh_NcpElems, &
mesh_maxNips, &
mesh_maxNipNeighbors
use IO, only: &
IO_timeStamp, &
IO_stringValue, &
IO_write_jobFile, &
homogenization and crystallite use new structure math_transpose33 has no advantage over transpose intrinsic
2018-06-02 22:58:08 +05:30
IO_error
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
use material
using array size not to rely on correctly set variable, cleaning crystallite data after use
2018-06-27 00:03:41 +05:30
use config, only: &
seems to work now anyway, nicer code
2018-08-22 18:00:51 +05:30
config_deallocate, &
consistent naming
2018-06-27 00:24:54 +05:30
config_crystallite, &
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
crystallite_name
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
use constitutive, only: &
initialise Fi correctly for initial field values away from equilibrium
2015-07-24 20:17:18 +05:30
constitutive_initialFi, &
fixed missing propagation of dependent state variables
2018-09-07 02:13:29 +05:30
constitutive_microstructure ! derived (shortcut) quantities of given state
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
implicit none
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
homogenization and crystallite use new structure math_transpose33 has no advantage over transpose intrinsic
2018-06-02 22:58:08 +05:30
integer(pInt), parameter :: FILEUNIT=434_pInt
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
logical, dimension(:,:), allocatable :: devNull
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
integer(pInt) :: &
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
c, & !< counter in integration point component loop
i, & !< counter in integration point loop
e, & !< counter in element loop
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
o = 0_pInt, & !< counter in output loop
fixed missing propagation of dependent state variables
2018-09-07 02:13:29 +05:30
r, &
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
cMax, & !< maximum number of integration point components
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
iMax, & !< maximum number of integration points
eMax, & !< maximum number of elements
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
myNcomponents, & !< number of components at current IP
moved reading in of lattice type and elastic constants to lattice module removed structure type for hex, fcc, bcc, now defining slip/twin systems for each phase found in material.config constitutive modules will only be initialized if needed homogenizedC function is only needed for models incorporating twinning in a physical way (titanmod and dislotwin)
2014-03-09 02:20:31 +05:30
mySize
gfortran does not recognice an empty array as 'present'
2018-06-22 11:33:22 +05:30
character(len=65536), dimension(:), allocatable :: str
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 simplifications
2016-06-30 02:57:22 +05:30
write(6,'(/,a)') ' <<<+- crystallite init -+>>>'
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
made code standard conform to Fortran 2008 (ignoring warning concerning comments beyond character 132). Basically, changing "x" format specifier to "Nx" ("1x") plus removing $ format specifier added compiler switches for gfortran and ifort to check for standard conformity old gnu compilers <4.4 are not longer supported because they don't provide the c binding for fftw
2012-02-01 00:48:55 +05:30
#include "compilation_info.f90"
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
renamed g -> c/ipc
2016-01-17 20:20:33 +05:30
cMax = homogenization_maxNgrains
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
iMax = mesh_maxNips
eMax = mesh_NcpElems
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
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
! ---------------------------------------------------------------------------
! ToDo (when working on homogenization): should be 3x3 tensor called S
renamed g -> c/ipc
2016-01-17 20:20:33 +05:30
allocate(crystallite_Tstar0_v(6,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_partionedTstar0_v(6,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_Tstar_v(6,cMax,iMax,eMax), source=0.0_pReal)
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
! ---------------------------------------------------------------------------
further separation still using old combined function
2019-01-18 19:12:44 +05:30
allocate(crystallite_subS0(3,3,cMax,iMax,eMax), source=0.0_pReal)
renamed g -> c/ipc
2016-01-17 20:20:33 +05:30
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(4,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_orientation0(4,cMax,iMax,eMax), source=0.0_pReal)
allocate(crystallite_rotation(4,cMax,iMax,eMax), source=0.0_pReal)
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.)
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
allocate(crystallite_output(maxval(crystallite_Noutput), &
consistent naming
2018-06-27 00:24:54 +05:30
size(config_crystallite))) ; crystallite_output = ''
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
allocate(crystallite_outputID(maxval(crystallite_Noutput), &
consistent naming
2018-06-27 00:24:54 +05:30
size(config_crystallite)), source=undefined_ID)
allocate(crystallite_sizePostResults(size(config_crystallite)),source=0_pInt)
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
allocate(crystallite_sizePostResult(maxval(crystallite_Noutput), &
consistent naming
2018-06-27 00:24:54 +05:30
size(config_crystallite)), source=0_pInt)
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
select case(numerics_integrator(1))
case(1_pInt)
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
integrateState => integrateStateFPI
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
case(2_pInt)
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
integrateState => integrateStateEuler
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
case(3_pInt)
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
integrateState => integrateStateAdaptiveEuler
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
case(4_pInt)
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
integrateState => integrateStateRK4
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
case(5_pInt)
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
integrateState => integrateStateRKCK45
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
end select
using separate functions for stress and tangent extensively tested in 46-simplification-of-crystallite-f90-NEW3 already
2019-01-18 20:00:50 +05:30
consistent naming
2018-06-27 00:24:54 +05:30
do c = 1_pInt, size(config_crystallite)
gfortran does not recognice an empty array as 'present'
2018-06-22 11:33:22 +05:30
#if defined(__GFORTRAN__)
str = ['GfortranBug86277']
consistent naming
2018-06-27 00:24:54 +05:30
str = config_crystallite(c)%getStrings('(output)',defaultVal=str)
gfortran does not recognice an empty array as 'present'
2018-06-22 11:33:22 +05:30
if (str(1) == 'GfortranBug86277') str = [character(len=65536)::]
#else
consistent naming
2018-06-27 00:24:54 +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
homogenization and crystallite use new structure math_transpose33 has no advantage over transpose intrinsic
2018-06-02 22:58:08 +05:30
do o = 1_pInt, size(str)
output tag was empty:
2018-06-03 00:30:47 +05:30
crystallite_output(o,c) = str(o)
homogenization and crystallite use new structure math_transpose33 has no advantage over transpose intrinsic
2018-06-02 22:58:08 +05:30
outputName: select case(str(o))
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
case ('phase') outputName
crystallite_outputID(o,c) = phase_ID
case ('texture') outputName
crystallite_outputID(o,c) = texture_ID
case ('volume') outputName
crystallite_outputID(o,c) = volume_ID
case ('orientation') outputName
crystallite_outputID(o,c) = orientation_ID
case ('grainrotation') outputName
crystallite_outputID(o,c) = grainrotation_ID
case ('eulerangles') outputName
crystallite_outputID(o,c) = eulerangles_ID
case ('defgrad','f') outputName
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
crystallite_outputID(o,c) = p_ID
case ('s','tstar','secondpiola','2ndpiola') outputName
crystallite_outputID(o,c) = s_ID
case ('elasmatrix') outputName
crystallite_outputID(o,c) = elasmatrix_ID
case ('neighboringip') outputName
crystallite_outputID(o,c) = neighboringip_ID
case ('neighboringelement') outputName
crystallite_outputID(o,c) = neighboringelement_ID
case default outputName
call IO_error(105_pInt,ext_msg=trim(str(o))//' (Crystallite)')
end select outputName
having no output in crystallite caused trouble
2018-06-20 02:28:46 +05:30
enddo
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
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
consistent naming
2018-06-27 00:24:54 +05:30
do r = 1_pInt,size(config_crystallite)
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
do o = 1_pInt,crystallite_Noutput(r)
select case(crystallite_outputID(o,r))
can be easily computed during post processing
2018-09-20 10:05:30 +05:30
case(phase_ID,texture_ID,volume_ID)
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
mySize = 1_pInt
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
case(orientation_ID,grainrotation_ID)
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
mySize = 4_pInt
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
case(eulerangles_ID)
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
mySize = 3_pInt
can be easily computed during post processing
2018-09-20 10:05:30 +05:30
case(defgrad_ID,fe_ID,fp_ID,fi_ID,lp_ID,li_ID,p_ID,s_ID)
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
mySize = 9_pInt
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
case(elasmatrix_ID)
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
mySize = 36_pInt
indroduced sourced allocation, enums where applicable (some parts still missing). fixed bug when having recursive file input
2013-12-12 22:39:59 +05:30
case(neighboringip_ID,neighboringelement_ID)
added crystallite output "neighboringelement" and "neighboringip" some polishing in crystallite_postResults
2013-05-08 17:32:30 +05:30
mySize = mesh_maxNipNeighbors
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
case default
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
mySize = 0_pInt
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
end select
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
crystallite_sizePostResult(o,r) = mySize
crystallite_sizePostResults(r) = crystallite_sizePostResults(r) + mySize
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
enddo
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
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
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
multi physics output now working for all solvers
2015-03-25 21:32:30 +05:30
if (worldrank == 0_pInt) then
call IO_write_jobFile(FILEUNIT,'outputCrystallite')
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
consistent naming
2018-06-27 00:24:54 +05:30
do r = 1_pInt,size(config_crystallite)
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
if (any(microstructure_crystallite(mesh_element(4,:)) == r)) then
write(FILEUNIT,'(/,a,/)') '['//trim(crystallite_name(r))//']'
do o = 1_pInt,crystallite_Noutput(r)
write(FILEUNIT,'(a,i4)') trim(crystallite_output(o,r))//char(9),crystallite_sizePostResult(o,r)
write output information only for active crystallites. removed unused variables. fixed some real == 0 comparisons
2015-04-21 17:53:00 +05:30
enddo
endif
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
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
multi physics output now working for all solvers
2015-03-25 21:32:30 +05:30
close(FILEUNIT)
add pheno+ module in
2015-10-14 00:22:01 +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
using array size not to rely on correctly set variable, cleaning crystallite data after use
2018-06-27 00:03:41 +05:30
call config_deallocate('material.config/crystallite')
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)
do e = FEsolving_execElem(1),FEsolving_execElem(2)
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
forall (i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), c = 1_pInt:myNcomponents)
crystallite_Fp0(1:3,1:3,c,i,e) = math_EulerToR(material_EulerAngles(1:3,c,i,e)) ! plastic def gradient reflects init orientation
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
crystallite_localPlasticity(c,i,e) = phase_localPlasticity(material_phase(c,i,e))
crystallite_Fe(1:3,1:3,c,i,e) = math_inv33(math_mul33x33(crystallite_Fi0(1:3,1:3,c,i,e), &
crystallite_Fp0(1:3,1:3,c,i,e))) ! assuming that euler angles are given in internal strain free configuration
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.
endforall
enddo
more work on the new state
2014-05-27 20:16:03 +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
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
if(any(.not. crystallite_localPlasticity) .and. .not. usePingPong) call IO_error(601_pInt) ! exit if nonlocal but no ping-pong ToDo: Why not check earlier? or in nonlocal?
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
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
crystallite_partionedFp0 = crystallite_Fp0
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
crystallite_partionedFi0 = crystallite_Fi0
more rename g -> c/ipc, wrong ipc iteration for crystallite_integrateStateRK4 and single run fixed
2016-01-17 23:26:24 +05:30
crystallite_partionedF0 = crystallite_F0
passing the entire defomration gradient information through microstructure
2016-04-13 23:36:04 +05:30
crystallite_partionedF = crystallite_F0
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
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
call crystallite_orientations()
crystallite_orientation0 = crystallite_orientation ! store initial orientations for calculation of grain rotations
homogenization_RGC.f90 >> just switching-off one of the debugging statement. crystallite.f90 >> optimizing the perturbation algorithm for local tangent.
2009-12-22 17:58:02 +05:30
flushes not needed + further cleaning
2018-09-19 23:15:57 +05:30
!$OMP PARALLEL DO
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do c = 1_pInt,homogenization_Ngrains(mesh_element(3,e))
call constitutive_microstructure(crystallite_orientation, &
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
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
enddo
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
enddo
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
!$OMP END PARALLEL DO
fixed two severe bugs in new state, now a little bit faster (399 vs. 378 sec) and more memory efficicient (251312 vs. 236832 kbytes) for Phenopowerlaw 16x16x16 example (short load only)
2014-06-17 00:49:38 +05:30
further separation still using old combined function
2019-01-18 19:12:44 +05:30
devNull = crystallite_stress()
call crystallite_stressTangent
define debug variables only if needed
2018-09-20 09:39:02 +05:30
#ifdef DEBUG
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) then
only print out the essential information
2018-09-20 10:28:31 +05:30
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
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
write(6,'(a42,1x,i10)') 'max # of neigbours/integration point: ', mesh_maxNipNeighbors
only print out the essential information
2018-09-20 10:28:31 +05:30
write(6,'(a42,1x,i10)') ' # of nonlocal constituents: ',count(.not. crystallite_localPlasticity)
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
flush(6)
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
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
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)
!--------------------------------------------------------------------------------------------------
function crystallite_stress()
use prec, only: &
tol_math_check, &
dNeq0
use numerics, only: &
subStepMinCryst, &
subStepSizeCryst, &
stepIncreaseCryst
#ifdef DEBUG
use debug, only: &
debug_level, &
debug_crystallite, &
debug_levelBasic, &
debug_levelExtensive, &
debug_levelSelective, &
debug_e, &
debug_i, &
debug_g
#endif
use IO, only: &
IO_warning, &
IO_error
use math, only: &
math_inv33, &
math_mul33x33, &
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
math_6toSym33, &
math_sym33to6
further separation still using old combined function
2019-01-18 19:12:44 +05:30
use mesh, only: &
mesh_NcpElems, &
mesh_element, &
mesh_maxNips, &
FE_geomtype
use material, only: &
homogenization_Ngrains, &
plasticState, &
sourceState, &
phase_Nsources, &
phaseAt, phasememberAt
use constitutive, only: &
constitutive_SandItsTangents, &
constitutive_LpAndItsTangents, &
constitutive_LiAndItsTangents
implicit none
logical, dimension(mesh_maxNips,mesh_NcpElems) :: crystallite_stress
real(pReal) :: &
formerSubStep
integer(pInt) :: &
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
#ifdef DEBUG
if (iand(debug_level(debug_crystallite),debug_levelSelective) /= 0_pInt &
.and. FEsolving_execElem(1) <= debug_e &
.and. debug_e <= FEsolving_execElem(2)) then
write(6,'(/,a,i8,1x,i2,1x,i3)') '<< CRYST >> boundary values at el ip ipc ', &
debug_e,debug_i, debug_g
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST >> F ', &
transpose(crystallite_partionedF(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST >> F0 ', &
transpose(crystallite_partionedF0(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST >> Fp0', &
transpose(crystallite_partionedFp0(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST >> Fi0', &
transpose(crystallite_partionedFi0(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST >> Lp0', &
transpose(crystallite_partionedLp0(1:3,1:3,debug_g,debug_i,debug_e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST >> Li0', &
transpose(crystallite_partionedLi0(1:3,1:3,debug_g,debug_i,debug_e))
endif
#endif
!--------------------------------------------------------------------------------------------------
! initialize to starting condition
crystallite_subStep = 0.0_pReal
!$OMP PARALLEL DO
elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do c = 1_pInt,homogenization_Ngrains(mesh_element(3,e))
homogenizationRequestsCalculation: if (crystallite_requested(c,i,e)) then
plasticState (phaseAt(c,i,e))%subState0( :,phasememberAt(c,i,e)) = &
plasticState (phaseAt(c,i,e))%partionedState0(:,phasememberAt(c,i,e))
do s = 1_pInt, phase_Nsources(phaseAt(c,i,e))
sourceState(phaseAt(c,i,e))%p(s)%subState0( :,phasememberAt(c,i,e)) = &
sourceState(phaseAt(c,i,e))%p(s)%partionedState0(:,phasememberAt(c,i,e))
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)
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
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) = math_6toSym33(crystallite_partionedTstar0_v(1:6,c,i,e))
further separation still using old combined function
2019-01-18 19:12:44 +05:30
crystallite_subFrac(c,i,e) = 0.0_pReal
crystallite_subStep(c,i,e) = 1.0_pReal/subStepSizeCryst
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_pInt
endIP = mesh_maxNips
endif singleRun
NiterationCrystallite = 0_pInt
cutbackLooping: do while (any(crystallite_todo(:,startIP:endIP,FEsolving_execELem(1):FEsolving_execElem(2))))
#ifdef DEBUG
if (iand(debug_level(debug_crystallite),debug_levelExtensive) /= 0_pInt) &
write(6,'(a,i6)') '<< CRYST >> crystallite iteration ',NiterationCrystallite
#endif
!$OMP PARALLEL DO PRIVATE(formerSubStep)
elementLooping3: do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do c = 1,homogenization_Ngrains(mesh_element(3,e))
!--------------------------------------------------------------------------------------------------
! wind forward
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), &
stepIncreaseCryst * crystallite_subStep(c,i,e))
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)
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
crystallite_subS0 (1:3,1:3,c,i,e) = math_6toSym33(crystallite_Tstar_v(1:6,c,i,e))
further separation still using old combined function
2019-01-18 19:12:44 +05:30
!if abbrevation, make c and p private in omp
plasticState( phaseAt(c,i,e))%subState0(:,phasememberAt(c,i,e)) &
= plasticState(phaseAt(c,i,e))%state( :,phasememberAt(c,i,e))
do s = 1_pInt, phase_Nsources(phaseAt(c,i,e))
sourceState( phaseAt(c,i,e))%p(s)%subState0(:,phasememberAt(c,i,e)) &
= sourceState(phaseAt(c,i,e))%p(s)%state( :,phasememberAt(c,i,e))
enddo
#ifdef DEBUG
if (iand(debug_level(debug_crystallite),debug_levelBasic) /= 0_pInt &
.and. ((e == debug_e .and. i == debug_i .and. c == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) &
write(6,'(a,f12.8,a,f12.8,a,i8,1x,i2,1x,i3,/)') '<< CRYST >> 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
#endif
endif
!--------------------------------------------------------------------------------------------------
! cut back (reduced time and restore)
else
crystallite_subStep(c,i,e) = subStepSizeCryst * crystallite_subStep(c,i,e)
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_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)
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
crystallite_Tstar_v(1:6,c,i,e) = math_sym33to6(crystallite_subS0(1:3,1:3,c,i,e))
further separation still using old combined function
2019-01-18 19:12:44 +05:30
plasticState (phaseAt(c,i,e))%state( :,phasememberAt(c,i,e)) &
= plasticState(phaseAt(c,i,e))%subState0(:,phasememberAt(c,i,e))
do s = 1_pInt, phase_Nsources(phaseAt(c,i,e))
sourceState( phaseAt(c,i,e))%p(s)%state( :,phasememberAt(c,i,e)) &
= sourceState(phaseAt(c,i,e))%p(s)%subState0(:,phasememberAt(c,i,e))
enddo
! cant restore dotState here, since not yet calculated in first cutback after initialization
crystallite_todo(c,i,e) = crystallite_subStep(c,i,e) > subStepMinCryst ! still on track or already done (beyond repair)
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((e == debug_e .and. i == debug_i .and. c == debug_g) &
.or. .not. iand(debug_level(debug_crystallite),debug_levelSelective) /= 0_pInt)) then
if (crystallite_todo(c,i,e)) then
write(6,'(a,f12.8,a,i8,1x,i2,1x,i3,/)') '<< CRYST >> cutback step in crystallite_stress &
&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 >> reached minimum step size &
&in crystallite_stress at el ip ipc ',e,i,c
endif
endif
#endif
endif
!--------------------------------------------------------------------------------------------------
! prepare for integration
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) = math_mul33x33(math_mul33x33(crystallite_subF (1:3,1:3,c,i,e), &
crystallite_invFp(1:3,1:3,c,i,e)), &
crystallite_invFi(1:3,1:3,c,i,e))
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
#ifdef DEBUG
if (iand(debug_level(debug_crystallite),debug_levelExtensive) /= 0_pInt) then
write(6,'(/,a,f8.5,a,f8.5,/)') '<< CRYST >> ',minval(crystallite_subStep),' ≤ subStep ≤ ',maxval(crystallite_subStep)
write(6,'(/,a,f8.5,a,f8.5,/)') '<< CRYST >> ',minval(crystallite_subFrac),' ≤ subFrac ≤ ',maxval(crystallite_subFrac)
flush(6)
if (iand(debug_level(debug_crystallite),debug_levelSelective) /= 0_pInt) then
write(6,'(/,a,f8.5,1x,a,1x,f8.5,1x,a)') '<< CRYST >> subFrac + subStep = ',&
crystallite_subFrac(debug_g,debug_i,debug_e),'+',crystallite_subStep(debug_g,debug_i,debug_e),'@selective'
flush(6)
endif
endif
#endif
!--------------------------------------------------------------------------------------------------
! integrate --- requires fully defined state array (basic + dependent state)
more "building block" separation
2019-01-24 03:32:21 +05:30
if (any(crystallite_todo)) call integrateState() ! TODO: unroll into proper elementloop to avoid N^2 for single point evaluation
where(.not. crystallite_converged .and. crystallite_subStep > subStepMinCryst) & ! do not try non-converged but fully cutbacked any further
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
NiterationCrystallite = NiterationCrystallite + 1_pInt
enddo cutbackLooping
! return whether converged or not
crystallite_stress = .false.
elementLooping5: do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) ! iterate over IPs of this element to be processed
crystallite_stress(i,e) = all(crystallite_converged(:,i,e))
enddo
enddo elementLooping5
#ifdef DEBUG
elementLooping6: do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) ! iterate over IPs of this element to be processed
do c = 1,homogenization_Ngrains(mesh_element(3,e))
if (.not. crystallite_converged(c,i,e)) then
if(iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) &
write(6,'(a,i8,1x,i2,1x,i3,/)') '<< CRYST >> no convergence at el ip ipc ', &
e,i,c
endif
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((e == debug_e .and. i == debug_i .and. c == debug_g) &
.or. .not. iand(debug_level(debug_crystallite),debug_levelSelective) /= 0_pInt)) then
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST >> central solution of cryst_StressAndTangent at el ip ipc ',e,i,c
write(6,'(/,a,/,3(12x,3(f12.4,1x)/))') '<< CRYST >> 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 >> Fp', &
transpose(crystallite_Fp(1:3,1:3,c,i,e))
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< CRYST >> Fi', &
transpose(crystallite_Fi(1:3,1:3,c,i,e))
write(6,'(a,/,3(12x,3(f14.9,1x)/),/)') '<< CRYST >> Lp', &
transpose(crystallite_Lp(1:3,1:3,c,i,e))
write(6,'(a,/,3(12x,3(f14.9,1x)/),/)') '<< CRYST >> Li', &
transpose(crystallite_Li(1:3,1:3,c,i,e))
flush(6)
endif
enddo
enddo
enddo elementLooping6
#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)
!--------------------------------------------------------------------------------------------------
subroutine crystallite_stressTangent()
use prec, only: &
tol_math_check, &
dNeq0
use IO, only: &
IO_warning, &
IO_error
use math, only: &
math_inv33, &
math_identity2nd, &
math_mul33x33, &
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
math_6toSym33, &
math_3333to99, &
math_99to3333, &
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
math_I3, &
math_mul3333xx3333, &
math_mul33xx33, &
using separate functions for stress and tangent extensively tested in 46-simplification-of-crystallite-f90-NEW3 already
2019-01-18 20:00:50 +05:30
math_invert2, &
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
math_det33
use mesh, only: &
mesh_element, &
FE_geomtype
use material, only: &
homogenization_Ngrains
use constitutive, only: &
constitutive_SandItsTangents, &
constitutive_LpAndItsTangents, &
constitutive_LiAndItsTangents
implicit none
integer(pInt) :: &
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)
do c = 1_pInt,homogenization_Ngrains(mesh_element(3,e))
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_Tstar_v(1:6,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_pInt,3_pInt; do p=1_pInt,3_pInt
lhs_3333(1:3,1:3,o,p) = lhs_3333(1:3,1:3,o,p) &
+ crystallite_subdt(c,i,e)*math_mul33x33(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)*math_mul33x33(invSubFi0,dLidS(1:3,1:3,o,p))
enddo;enddo
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
call math_invert2(temp_99,error,math_3333to99(lhs_3333))
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
if (error) then
call IO_warning(warning_ID=600_pInt,el=e,ip=i,g=c, &
ext_msg='inversion error in analytic tangent calculation')
dFidS = 0.0_pReal
else
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
dFidS = math_mul3333xx3333(math_99to3333(temp_99),rhs_3333)
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
endif
dLidS = math_mul3333xx3333(dLidFi,dFidS) + dLidS
endif
call constitutive_LpAndItsTangents(devNull,dLpdS,dLpdFi, &
crystallite_Tstar_v(1:6,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
!--------------------------------------------------------------------------------------------------
! calculate dSdF
temp_33_1 = transpose(math_mul33x33(crystallite_invFp(1:3,1:3,c,i,e), &
crystallite_invFi(1:3,1:3,c,i,e)))
temp_33_2 = math_mul33x33( crystallite_subF (1:3,1:3,c,i,e), &
math_inv33(crystallite_subFp0(1:3,1:3,c,i,e)))
temp_33_3 = math_mul33x33(math_mul33x33(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)))
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt)
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
rhs_3333(p,o,1:3,1:3) = math_mul33x33(dSdFe(p,o,1:3,1:3),temp_33_1)
temp_3333(1:3,1:3,p,o) = math_mul33x33(math_mul33x33(temp_33_2,dLpdS(1:3,1:3,p,o)), &
crystallite_invFi(1:3,1:3,c,i,e)) &
+ math_mul33x33(temp_33_3,dLidS(1:3,1:3,p,o))
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
end forall
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
lhs_3333 = crystallite_subdt(c,i,e)*math_mul3333xx3333(dSdFe,temp_3333) + &
math_mul3333xx3333(dSdFi,dFidS)
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
call math_invert2(temp_99,error,math_identity2nd(9_pInt)+math_3333to99(lhs_3333))
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
if (error) then
call IO_warning(warning_ID=600_pInt,el=e,ip=i,g=c, &
ext_msg='inversion error in analytic tangent calculation')
dSdF = rhs_3333
else
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
dSdF = math_mul3333xx3333(math_99to3333(temp_99),rhs_3333)
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
endif
!--------------------------------------------------------------------------------------------------
! calculate dFpinvdF
temp_3333 = math_mul3333xx3333(dLpdS,dSdF)
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt)
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
dFpinvdF(1:3,1:3,p,o) &
= -crystallite_subdt(c,i,e) &
* math_mul33x33(math_inv33(crystallite_subFp0(1:3,1:3,c,i,e)), &
math_mul33x33(temp_3333(1:3,1:3,p,o),crystallite_invFi(1:3,1:3,c,i,e)))
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
end forall
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
!--------------------------------------------------------------------------------------------------
! assemble dPdF
temp_33_1 = math_mul33x33(crystallite_invFp(1:3,1:3,c,i,e), &
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
math_mul33x33(math_6toSym33(crystallite_Tstar_v(1:6,c,i,e)), &
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
transpose(crystallite_invFp(1:3,1:3,c,i,e))))
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
temp_33_2 = math_mul33x33(math_6toSym33(crystallite_Tstar_v(1:6,c,i,e)), &
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
transpose(crystallite_invFp(1:3,1:3,c,i,e)))
temp_33_3 = math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
crystallite_invFp(1:3,1:3,c,i,e))
temp_33_4 = math_mul33x33(math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
crystallite_invFp(1:3,1:3,c,i,e)), &
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
math_6toSym33(crystallite_Tstar_v(1:6,c,i,e)))
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
crystallite_dPdF(1:3,1:3,1:3,1:3,c,i,e) = 0.0_pReal
do p=1_pInt, 3_pInt
crystallite_dPdF(p,1:3,p,1:3,c,i,e) = transpose(temp_33_1)
enddo
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt)
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
crystallite_dPdF(1:3,1:3,p,o,c,i,e) = crystallite_dPdF(1:3,1:3,p,o,c,i,e) + &
math_mul33x33(math_mul33x33(crystallite_subF(1:3,1:3,c,i,e),dFpinvdF(1:3,1:3,p,o)),temp_33_2) + &
math_mul33x33(math_mul33x33(temp_33_3,dSdF(1:3,1:3,p,o)),transpose(crystallite_invFp(1:3,1:3,c,i,e))) + &
math_mul33x33(temp_33_4,transpose(dFpinvdF(1:3,1:3,p,o)))
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
end forall
preparing for separation of stress calculation and tangent calculatin
2019-01-18 16:46:26 +05:30
enddo; enddo
enddo elementLooping
!$OMP END PARALLEL DO
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
use math, only: &
math_rotationalPart33, &
math_RtoQ
introduced dummy heat calculation, available as an crystallite output temperature, avgP and avgF are available as homogenization output. move crystallite output of ipcoordinates to homogenization
2013-10-19 00:27:28 +05:30
use material, only: &
fixed small allocation flaws for new state
2014-05-12 18:30:37 +05:30
plasticState, &
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
material_phase, &
homogenization_Ngrains
use mesh, only: &
mesh_element
use lattice, only: &
lattice_qDisorientation
use plastic_nonlocal, only: &
plastic_nonlocal_updateCompatibility
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
doxygen comments
2013-02-22 04:38:36 +05:30
implicit none
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
integer(pInt) &
c, & !< counter in integration point component loop
i, & !< counter in integration point loop
cleaning no reason to assume that the math functions are not thread safe
2019-01-19 14:32:04 +05:30
e !< counter in element loop
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
!$OMP PARALLEL DO
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do c = 1_pInt,homogenization_Ngrains(mesh_element(3,e))
crystallite_orientation(1:4,c,i,e) = math_RtoQ(transpose(math_rotationalPart33(crystallite_Fe(1:3,1:3,c,i,e))))
crystallite_rotation(1:4,c,i,e) = lattice_qDisorientation(crystallite_orientation0(1:4,c,i,e), &! active rotation from initial
cleaning no reason to assume that the math functions are not thread safe
2019-01-19 14:32:04 +05:30
crystallite_orientation(1:4,c,i,e)) ! to current orientation (with no symmetry)
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
enddo; enddo; enddo
!$OMP END PARALLEL DO
! --- we use crystallite_orientation from above, so need a separate loop
nonlocalPresent: if (any(plasticState%nonLocal)) then
cleaning no reason to assume that the math functions are not thread safe
2019-01-19 14:32:04 +05:30
!$OMP PARALLEL DO
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
cleaning no reason to assume that the math functions are not thread safe
2019-01-19 14:32:04 +05:30
if (plasticState(material_phase(1,i,e))%nonLocal) & ! if nonlocal model
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +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)
use math, only: &
math_mul33x33, &
math_inv33, &
math_EulerToR
use material, only: &
material_EulerAngles
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
implicit none
real(pReal), dimension(3,3) :: crystallite_push33ToRef
real(pReal), dimension(3,3), intent(in) :: tensor33
real(pReal), dimension(3,3) :: T
integer(pInt), intent(in):: &
el, & ! element index
ip, & ! integration point index
ipc ! grain index
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
T = math_mul33x33(math_EulerToR(material_EulerAngles(1:3,ipc,ip,el)), &
transpose(math_inv33(crystallite_subF(1:3,1:3,ipc,ip,el))))
crystallite_push33ToRef = math_mul33x33(transpose(T),math_mul33x33(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)
use math, only: &
math_qToEuler, &
math_qToEulerAxisAngle, &
math_mul33x33, &
math_det33, &
math_I3, &
inDeg, &
math_6toSym33
use mesh, only: &
mesh_element, &
mesh_ipVolume, &
mesh_maxNipNeighbors, &
mesh_ipNeighborhood, &
FE_NipNeighbors, &
FE_geomtype, &
FE_celltype
use material, only: &
plasticState, &
sourceState, &
microstructure_crystallite, &
crystallite_Noutput, &
material_phase, &
material_texture, &
homogenization_Ngrains
use constitutive, only: &
constitutive_homogenizedC, &
constitutive_postResults
implicit none
integer(pInt), intent(in):: &
el, & !< element index
ip, & !< integration point index
ipc !< grain index
real(pReal), dimension(1+crystallite_sizePostResults(microstructure_crystallite(mesh_element(4,el))) + &
1+plasticState(material_phase(ipc,ip,el))%sizePostResults + &
sum(sourceState(material_phase(ipc,ip,el))%p(:)%sizePostResults)) :: &
crystallite_postResults
real(pReal) :: &
detF
integer(pInt) :: &
o, &
c, &
crystID, &
mySize, &
n
crystID = microstructure_crystallite(mesh_element(4,el))
crystallite_postResults = 0.0_pReal
c = 0_pInt
crystallite_postResults(c+1) = real(crystallite_sizePostResults(crystID),pReal) ! size of results from cryst
c = c + 1_pInt
do o = 1_pInt,crystallite_Noutput(crystID)
mySize = 0_pInt
select case(crystallite_outputID(o,crystID))
case (phase_ID)
mySize = 1_pInt
crystallite_postResults(c+1) = real(material_phase(ipc,ip,el),pReal) ! phaseID of grain
case (texture_ID)
mySize = 1_pInt
crystallite_postResults(c+1) = real(material_texture(ipc,ip,el),pReal) ! textureID of grain
case (volume_ID)
mySize = 1_pInt
detF = math_det33(crystallite_partionedF(1:3,1:3,ipc,ip,el)) ! V_current = det(F) * V_reference
crystallite_postResults(c+1) = detF * mesh_ipVolume(ip,el) &
/ real(homogenization_Ngrains(mesh_element(3,el)),pReal) ! grain volume (not fraction but absolute)
case (orientation_ID)
mySize = 4_pInt
crystallite_postResults(c+1:c+mySize) = crystallite_orientation(1:4,ipc,ip,el) ! grain orientation as quaternion
case (eulerangles_ID)
mySize = 3_pInt
crystallite_postResults(c+1:c+mySize) = inDeg &
* math_qToEuler(crystallite_orientation(1:4,ipc,ip,el)) ! grain orientation as Euler angles in degree
case (grainrotation_ID)
mySize = 4_pInt
crystallite_postResults(c+1:c+mySize) = &
math_qToEulerAxisAngle(crystallite_rotation(1:4,ipc,ip,el)) ! grain rotation away from initial orientation as axis-angle in sample reference coordinates
crystallite_postResults(c+4) = inDeg * crystallite_postResults(c+4) ! angle in degree
! 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"
case (defgrad_ID)
mySize = 9_pInt
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_partionedF(1:3,1:3,ipc,ip,el)),[mySize])
case (fe_ID)
mySize = 9_pInt
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Fe(1:3,1:3,ipc,ip,el)),[mySize])
case (fp_ID)
mySize = 9_pInt
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Fp(1:3,1:3,ipc,ip,el)),[mySize])
case (fi_ID)
mySize = 9_pInt
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Fi(1:3,1:3,ipc,ip,el)),[mySize])
case (lp_ID)
mySize = 9_pInt
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Lp(1:3,1:3,ipc,ip,el)),[mySize])
case (li_ID)
mySize = 9_pInt
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_Li(1:3,1:3,ipc,ip,el)),[mySize])
case (p_ID)
mySize = 9_pInt
crystallite_postResults(c+1:c+mySize) = &
reshape(transpose(crystallite_P(1:3,1:3,ipc,ip,el)),[mySize])
case (s_ID)
mySize = 9_pInt
crystallite_postResults(c+1:c+mySize) = &
reshape(math_6toSym33(crystallite_Tstar_v(1:6,ipc,ip,el)),[mySize])
case (elasmatrix_ID)
mySize = 36_pInt
crystallite_postResults(c+1:c+mySize) = reshape(constitutive_homogenizedC(ipc,ip,el),[mySize])
case(neighboringelement_ID)
mySize = mesh_maxNipNeighbors
crystallite_postResults(c+1:c+mySize) = 0.0_pReal
forall (n = 1_pInt:FE_NipNeighbors(FE_celltype(FE_geomtype(mesh_element(2,el))))) &
crystallite_postResults(c+n) = real(mesh_ipNeighborhood(1,n,ip,el),pReal)
case(neighboringip_ID)
mySize = mesh_maxNipNeighbors
crystallite_postResults(c+1:c+mySize) = 0.0_pReal
forall (n = 1_pInt:FE_NipNeighbors(FE_celltype(FE_geomtype(mesh_element(2,el))))) &
crystallite_postResults(c+n) = real(mesh_ipNeighborhood(2,n,ip,el),pReal)
end select
c = c + mySize
enddo
crystallite_postResults(c+1) = real(plasticState(material_phase(ipc,ip,el))%sizePostResults,pReal) ! size of constitutive results
c = c + 1_pInt
if (size(crystallite_postResults)-c > 0_pInt) &
crystallite_postResults(c+1:size(crystallite_postResults)) = &
constitutive_postResults(crystallite_Tstar_v(1:6,ipc,ip,el), crystallite_Fi(1:3,1:3,ipc,ip,el), &
crystallite_Fe, ipc, ip, el)
end function crystallite_postResults
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
!--------------------------------------------------------------------------------------------------
logical function integrateStress(&
ipc,& ! grain number
ip,& ! integration point number
el,& ! element number
timeFraction &
)
use, intrinsic :: &
IEEE_arithmetic
use prec, only: pLongInt, &
tol_math_check, &
dEq0
use numerics, only: nStress, &
aTol_crystalliteStress, &
rTol_crystalliteStress, &
iJacoLpresiduum, &
subStepSizeLp, &
subStepSizeLi
#ifdef DEBUG
use debug, only: debug_level, &
debug_e, &
debug_i, &
debug_g, &
debug_crystallite, &
debug_levelBasic, &
debug_levelExtensive, &
debug_levelSelective
#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
use constitutive, only: constitutive_LpAndItsTangents, &
constitutive_LiAndItsTangents, &
constitutive_SandItsTangents
use math, only: math_mul33x33, &
math_mul33xx33, &
math_mul3333xx3333, &
math_inv33, &
math_det33, &
math_I3, &
math_identity2nd, &
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
math_sym33to6, &
math_3333to99, &
math_33to9, &
math_9to33
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
implicit none
integer(pInt), intent(in):: el, & ! element index
ip, & ! integration point index
ipc ! grain index
real(pReal), optional, intent(in) :: timeFraction ! fraction of timestep
add pheno+ module in
2015-10-14 00:22:01 +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(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
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
S, & ! 2nd Piola-Kirchhoff Stress in plastic (lattice) configuration
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
A, &
B, &
Fe, & ! elastic deformation gradient
temp_33
real(pReal), dimension(9):: work ! needed for matrix inversion by LAPACK
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
integer(pInt), 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)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
dRLp_dLp2, & ! working copy of dRdLp
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
dRLi_dLi ! partial derivative of residuumI (Jacobian for Newton-Raphson scheme)
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
real(pReal), dimension(3,3,3,3):: dS_dFe, & ! partial derivative of 2nd Piola-Kirchhoff stress
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
dS_dFi, &
do concurrent causes problems on some Intel compilers use forall instead. Mandel/Plain prefixes not needed any more
2019-01-19 03:50:44 +05:30
dFe_dLp, & ! partial derivative of elastic deformation gradient
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
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(pInt) 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
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
!* be pessimistic
integrateStress = .false.
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) &
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST >> integrateStress at el ip ipc ',el,ip,ipc
#endif
changed to new state, please report bugs to Luv or Martin
2014-07-02 17:57:39 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!* feed local variables
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
Lpguess = crystallite_Lp(1:3,1:3,ipc,ip,el) ! ... and take it as first guess
Liguess = crystallite_Li(1:3,1:3,ipc,ip,el) ! ... and take it as first guess
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
Liguess_old = Liguess
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
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
invFp_current = math_inv33(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
failedInversionFp: if (all(dEq0(invFp_current))) then
#ifdef DEBUG
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) &
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_pInt) &
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
return
endif failedInversionFp
A = math_mul33x33(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
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
invFi_current = math_inv33(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
failedInversionFi: if (all(dEq0(invFi_current))) then
#ifdef DEBUG
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) &
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_pInt) &
write(6,'(/,a,/,3(12x,3(f12.7,1x)/))') '<< CRYST >> 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
return
endif failedInversionFi
* 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
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
!* start Li loop with normal step length
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
NiterationStressLi = 0_pInt
jacoCounterLi = 0_pInt
steplengthLi = 1.0_pReal
residuumLi_old = 0.0_pReal
LiLoop: do
NiterationStressLi = NiterationStressLi + 1_pInt
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
LiLoopLimit: if (NiterationStressLi > nStress) then
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
#ifdef DEBUG
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) &
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
write(6,'(a,i3,a,i8,1x,i2,1x,i3,/)') '<< CRYST >> integrateStress reached Li loop limit',nStress, &
' at el ip ipc ', el,ip,ipc
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
return
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
endif LiLoopLimit
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
invFi_new = math_mul33x33(invFi_current,math_I3 - dt*Liguess)
Fi_new = math_inv33(invFi_new)
detInvFi = math_det33(invFi_new)
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
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
!* start Lp loop with normal step length
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
NiterationStressLp = 0_pInt
jacoCounterLp = 0_pInt
steplengthLp = 1.0_pReal
residuumLp_old = 0.0_pReal
Lpguess_old = Lpguess
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
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
LpLoop: do
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
NiterationStressLp = NiterationStressLp + 1_pInt
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
LpLoopLimit: if (NiterationStressLp > 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
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) &
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
write(6,'(a,i3,a,i8,1x,i2,1x,i3,/)') '<< CRYST >> integrateStress reached Lp loop limit',nStress, &
' 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
return
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
endif LpLoopLimit
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!* calculate (elastic) 2nd Piola--Kirchhoff stress tensor and its tangent from constitutive law
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
B = math_I3 - dt*Lpguess
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
Fe = math_mul33x33(math_mul33x33(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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!* calculate plastic velocity gradient and its tangent from constitutive law
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
call constitutive_LpAndItsTangents(Lp_constitutive, dLp_dS, dLp_dFi, &
math_sym33to6(S), Fi_new, ipc, ip, el)
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
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#ifdef DEBUG
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
write(6,'(a,i3,/)') '<< CRYST >> stress iteration ', NiterationStressLp
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST >> Lpguess', transpose(Lpguess)
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST >> Fi', transpose(Fi_new)
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST >> Fe', transpose(Fe)
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST >> S', transpose(S)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
write(6,'(a,/,3(12x,3(e20.10,1x)/))') '<< CRYST >> Lp_constitutive', transpose(Lp_constitutive)
endif
#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
!* update current residuum and check for convergence of loop
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
aTolLp = max(rTol_crystalliteStress * max(norm2(Lpguess),norm2(Lp_constitutive)), & ! absolute tolerance from largest acceptable relative error
aTol_crystalliteStress) ! minimum lower cutoff
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
residuumLp = Lpguess - Lp_constitutive
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
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
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
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) &
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
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
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_pInt &
.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...
steplengthLp = subStepSizeLp * steplengthLp ! ...try with smaller step length in same direction
takeover from old 46-XXX branch
2019-01-19 03:39:46 +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
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
write(6,'(a,1x,f7.4)') '<< CRYST >> linear search for Lpguess with step', steplengthLp
more work on the new state
2014-05-27 20:16:03 +05:30
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
#endif
cycle LpLoop
doxygen comments
2013-02-22 04:38:36 +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
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
if (mod(jacoCounterLp, iJacoLpresiduum) == 0_pInt) then
forall(o=1_pInt:3_pInt,p=1_pInt:3_pInt) &
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
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)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
dFe_dLp = - dt * dFe_dLp
dRLp_dLp = math_identity2nd(9_pInt) &
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
- 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
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
write(6,'(a,/,9(12x,9(e12.4,1x)/))') '<< CRYST >> dLp_dS', math_3333to99(dLp_dS)
write(6,'(a,1x,e20.10)') '<< CRYST >> dLp_dS norm', norm2(math_3333to99(dLp_dS))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
write(6,'(a,/,9(12x,9(e12.4,1x)/))') '<< CRYST >> dRLp_dLp', dRLp_dLp - math_identity2nd(9_pInt)
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
write(6,'(a,1x,e20.10)') '<< CRYST >> dRLp_dLp norm', norm2(dRLp_dLp - math_identity2nd(9_pInt))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
endif
#endif
takeover from old 46-XXX branch
2019-01-19 03:39:46 +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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (ierr /= 0_pInt) then
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) then
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST >> integrateStress failed on dR/dLp inversion 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
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g)&
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
write(6,*)
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST >> dR_dLp',transpose(dRLp_dLp)
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST >> dFe_dLp',transpose(math_3333to99(dFe_dLp))
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST >> dS_dFe_constitutive',transpose(math_3333to99(dS_dFe))
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST >> dLp_dS_constitutive',transpose(math_3333to99(dLp_dS))
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST >> A',transpose(A)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST >> B',transpose(B)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST >> Lp_constitutive',transpose(Lp_constitutive)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST >> Lpguess',transpose(Lpguess)
endif
endif
#endif
return
endif
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
deltaLp = - math_9to33(work)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
endif
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
jacoCounterLp = jacoCounterLp + 1_pInt
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
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +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
call constitutive_LiAndItsTangents(Li_constitutive, dLi_dS, dLi_dFi, &
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
math_sym33to6(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
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST >> Li_constitutive', transpose(Li_constitutive)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST >> Liguess', transpose(Liguess)
endif
#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
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
!* update current residuum and check for convergence of loop
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
aTolLi = max(rTol_crystalliteStress * max(norm2(Liguess),norm2(Li_constitutive)), & ! absolute tolerance from largest acceptable relative error
aTol_crystalliteStress) ! minimum lower cutoff
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
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) &
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..!'
#endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +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_pInt &
.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...
steplengthLi = subStepSizeLi * steplengthLi ! ...try with smaller step length in same direction
Liguess = Liguess_old + steplengthLi * deltaLi
cycle LiLoop
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
if (mod(jacoCounterLi, iJacoLpresiduum) == 0_pInt) then
temp_33 = math_mul33x33(math_mul33x33(A,B),invFi_current)
forall(o=1_pInt:3_pInt,p=1_pInt:3_pInt)
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
end forall
forall(o=1_pInt:3_pInt,p=1_pInt:3_pInt) &
dFi_dLi(1:3,1:3,o,p) = math_mul33x33(math_mul33x33(Fi_new,dFi_dLi(1:3,1:3,o,p)),Fi_new)
dRLi_dLi = math_identity2nd(9_pInt) &
takeover from old 46-XXX branch
2019-01-19 03:39:46 +05:30
- math_3333to99(math_mul3333xx3333(dLi_dS, math_mul3333xx3333(dS_dFe, dFe_dLi) + &
math_mul3333xx3333(dS_dFi, dFi_dLi))) &
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
- 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_pInt) then
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) 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_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g)&
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
write(6,*)
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST >> dR_dLi',transpose(dRLi_dLi)
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST >> dFe_dLi',transpose(math_3333to99(dFe_dLi))
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST >> dS_dFi_constitutive',transpose(math_3333to99(dS_dFi))
write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST >> dLi_dS_constitutive',transpose(math_3333to99(dLi_dS))
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST >> Li_constitutive',transpose(Li_constitutive)
write(6,'(a,/,3(12x,3(e20.7,1x)/))') '<< CRYST >> Liguess',transpose(Liguess)
endif
endif
#endif
return
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
[skip ci] re-ordered functions internal/private functions at the end
2019-01-19 14:05:45 +05:30
deltaLi = - math_9to33(work)
endif
jacoCounterLi = jacoCounterLi + 1_pInt
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
Liguess = Liguess + steplengthLi * deltaLi
enddo LiLoop
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
!* calculate new plastic and elastic deformation gradient
invFp_new = math_mul33x33(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
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) 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_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) &
write(6,'(/,a,/,3(12x,3(f12.7,1x)/))') '<< CRYST >> invFp_new',transpose(invFp_new)
endif
#endif
return
endif failedInversionInvFp
Fe_new = math_mul33x33(math_mul33x33(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
integrateStress = .true.
crystallite_P (1:3,1:3,ipc,ip,el) = math_mul33x33(math_mul33x33(Fg_new,invFp_new), &
math_mul33x33(S,transpose(invFp_new)))
crystallite_Tstar_v (1:6,ipc,ip,el) = math_sym33to6(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
if (iand(debug_level(debug_crystallite),debug_levelExtensive) /= 0_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
write(6,'(a,/,3(12x,3(f12.7,1x)/))') '<< CRYST >> 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 >> Cauchy / MPa', &
math_mul33x33(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 >> Fe Lp Fe^-1', &
transpose(math_mul33x33(Fe_new, math_mul33x33(crystallite_Lp(1:3,1:3,ipc,ip,el), math_inv33(Fe_new))))
write(6,'(a,/,3(12x,3(f12.7,1x)/))') '<< CRYST >> Fp',transpose(crystallite_Fp(1:3,1:3,ipc,ip,el))
write(6,'(a,/,3(12x,3(f12.7,1x)/))') '<< CRYST >> Fi',transpose(crystallite_Fi(1:3,1:3,ipc,ip,el))
endif
#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
!--------------------------------------------------------------------------------------------------
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
subroutine integrateStateFPI()
using IEEE_is_NAN and IEEE_quiet_NaN instead of hand-written solution, will not work for gfortran < 5
2017-05-04 04:02:44 +05:30
use, intrinsic :: &
IEEE_arithmetic
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
#ifdef DEBUG
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
use debug, only: &
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
debug_e, &
debug_i, &
debug_g, &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
debug_level,&
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
debug_crystallite, &
debug_levelBasic, &
debug_levelExtensive, &
Deleted performance profiling code
2018-05-17 20:03:35 +05:30
debug_levelSelective
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#endif
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
use numerics, only: &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
nState, &
integrationMode was always 1, very obscure code was most likely the leftover of old functionality
2018-08-05 20:36:03 +05:30
rTol_crystalliteState
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
use mesh, only: &
mesh_element, &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
mesh_NcpElems
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
use material, only: &
more work on the new state
2014-05-27 20:16:03 +05:30
plasticState, &
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
sourceState, &
exchanged mappingConstitutive array for more descriptive "phaseAt and phasememberAt"
2016-01-15 05:49:44 +05:30
phaseAt, phasememberAt, &
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
phase_Nsources, &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
homogenization_Ngrains
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
use constitutive, only: &
more work on the new state
2014-05-27 20:16:03 +05:30
constitutive_collectDotState, &
constitutive_microstructure, &
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
constitutive_plasticity_maxSizeDotState, &
constitutive_source_maxSizeDotState
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
doxygen comments
2013-02-22 04:38:36 +05:30
implicit none
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
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
integer(pInt) :: &
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, &
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
mySource, &
Heavy bug fixing: Butcher tableau multiplication in RK45 integration scheme, state updating of adaptive Euler, wrong assignment to evolution rate instead of state in fixed-point iterator. Exchanged possibly non thread-safe matmul in RKCK45, simplified RK45 integration step 1--3 vs 4 handling. (Thanks to Pratheek, Luv, and Chen for their help!)
2014-09-03 01:16:52 +05:30
mySizePlasticDotState, & ! size of dot states
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
mySizeSourceDotState
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
integer(pInt), dimension(2) :: &
eIter ! bounds for element iteration
integer(pInt), dimension(2,mesh_NcpElems) :: &
iIter, & ! bounds for ip iteration
gIter ! bounds for grain iteration
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
real(pReal) :: &
dot_prod12, &
dot_prod22, &
plasticStateDamper, & ! damper for integration of state
sourceStateDamper
real(pReal), dimension(constitutive_plasticity_maxSizeDotState) :: &
plasticStateResiduum, &
tempPlasticState
real(pReal), dimension(constitutive_source_maxSizeDotState, maxval(phase_Nsources)) :: &
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
sourceStateResiduum, & ! residuum from evolution in micrstructure
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
tempSourceState
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 :: &
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
converged, &
name conflict
2015-09-10 03:22:00 +05:30
NaN, &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
singleRun, & ! flag indicating computation for single (g,i,e) triple
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
simplified logic for single run, now not checking for extra arguments but simply checking if start and end elem/IP are the same
2013-04-26 18:53:36 +05:30
eIter = FEsolving_execElem(1:2)
do e = eIter(1),eIter(2)
iIter(1:2,e) = FEsolving_execIP(1:2,e)
gIter(1:2,e) = [ 1_pInt,homogenization_Ngrains(mesh_element(3,e))]
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
fixed out of range bug introduced with singleRun changes in rev. 2356
2013-04-30 17:44:07 +05:30
singleRun = (eIter(1) == eIter(2) .and. iIter(1,eIter(1)) == iIter(2,eIter(2)))
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
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt) &
write(6,'(a,i8,a)') '<< CRYST >> ', count(crystallite_todo(:,:,:)),' grains todo at start of state integration'
#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
! --+>> PREGUESS FOR STATE <<+--
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
call update_dotState(1.0_pReal)
call update_state(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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --+>> STATE LOOP <<+--
add pheno+ module in
2015-10-14 00:22:01 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
NiterationState = 0_pInt
doneWithIntegration = .false.
crystalliteLooping: do while (.not. doneWithIntegration .and. NiterationState < nState)
NiterationState = NiterationState + 1_pInt
more "building block" separation
2019-01-24 03:32:21 +05:30
! store previousDotState and previousDotState2
!$OMP PARALLEL DO PRIVATE(p,c)
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(mesh_element(3,e))
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
p = phaseAt(g,i,e); c = phasememberAt(g,i,e)
plasticState(p)%previousDotState2(:,c) = merge(plasticState(p)%previousDotState(:,c),&
0.0_pReal,&
NiterationState > 1_pInt)
plasticState(p)%previousDotState (:,c) = plasticState(p)%dotState(:,c)
do s = 1_pInt, phase_Nsources(p)
sourceState(p)%p(s)%previousDotState2(:,c) = merge(sourceState(p)%p(s)%previousDotState(:,c),&
0.0_pReal, &
NiterationState > 1_pInt)
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
!$OMP 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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- STRESS INTEGRATION ---
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
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt) &
write(6,'(a,i8,a)') '<< CRYST >> ', count(crystallite_todo(:,:,:)),' grains todo before stress integration'
#endif
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
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
crystallite_todo(g,i,e) = integrateStress(g,i,e)
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP FLUSH(crystallite_todo)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! broken non-local...
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP CRITICAL (checkTodo)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ... then all non-locals skipped
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP END CRITICAL (checkTodo)
endif
endif
enddo; enddo; enddo
!$OMP ENDDO
common function to update dot state
2019-01-23 10:44:19 +05:30
!$OMP END PARALLEL
call update_dotState(1.0_pReal)
!$OMP PARALLEL
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- UPDATE STATE ---
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
!$OMP DO PRIVATE(dot_prod12,dot_prod22, &
!$OMP& mySizePlasticDotState,mySizeSourceDotState, &
!$OMP& plasticStateResiduum,sourceStateResiduum, &
!$OMP& plasticStatedamper,sourceStateDamper, &
!$OMP& tempPlasticState,tempSourceState,converged,p,c)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
more work on the new state
2014-05-27 20:16:03 +05:30
exchanged mappingConstitutive array for more descriptive "phaseAt and phasememberAt"
2016-01-15 05:49:44 +05:30
p = phaseAt(g,i,e)
c = phasememberAt(g,i,e)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
dot_prod12 = dot_product( plasticState(p)%dotState (:,c) &
- plasticState(p)%previousDotState (:,c), &
plasticState(p)%previousDotState (:,c) &
- plasticState(p)%previousDotState2(:,c))
dot_prod22 = dot_product( plasticState(p)%previousDotState (:,c) &
- plasticState(p)%previousDotState2(:,c), &
plasticState(p)%previousDotState (:,c) &
- plasticState(p)%previousDotState2(:,c))
if ( dot_prod22 > 0.0_pReal &
.and. ( dot_prod12 < 0.0_pReal &
.or. dot_product(plasticState(p)%dotState(:,c), &
plasticState(p)%previousDotState(:,c)) < 0.0_pReal) ) then
plasticStateDamper = 0.75_pReal + 0.25_pReal * tanh(2.0_pReal + 4.0_pReal * dot_prod12 / dot_prod22)
else
plasticStateDamper = 1.0_pReal
endif
! --- get residui ---
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
mySizePlasticDotState = plasticState(p)%sizeDotState
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
plasticStateResiduum(1:mySizePlasticDotState) = &
plasticState(p)%state(1:mySizePlasticDotState,c) &
- plasticState(p)%subState0(1:mySizePlasticDotState,c) &
- ( plasticState(p)%dotState(1:mySizePlasticDotState,c) * plasticStateDamper &
+ plasticState(p)%previousDotState(1:mySizePlasticDotState,c) &
* (1.0_pReal - plasticStateDamper)) * crystallite_subdt(g,i,e)
! --- correct state with residuum ---
tempPlasticState(1:mySizePlasticDotState) = &
plasticState(p)%state(1:mySizePlasticDotState,c) &
- plasticStateResiduum(1:mySizePlasticDotState) ! need to copy to local variable, since we cant flush a pointer in openmp
! --- store corrected dotState --- (cannot do this before state update, because not sure how to flush pointers in openmp)
plasticState(p)%dotState(:,c) = plasticState(p)%dotState(:,c) * plasticStateDamper &
+ plasticState(p)%previousDotState(:,c) &
* (1.0_pReal - plasticStateDamper)
add pheno+ module in
2015-10-14 00:22:01 +05:30
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
do mySource = 1_pInt, phase_Nsources(p)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
dot_prod12 = dot_product( sourceState(p)%p(mySource)%dotState (:,c) &
- sourceState(p)%p(mySource)%previousDotState (:,c), &
sourceState(p)%p(mySource)%previousDotState (:,c) &
- sourceState(p)%p(mySource)%previousDotState2(:,c))
dot_prod22 = dot_product( sourceState(p)%p(mySource)%previousDotState (:,c) &
- sourceState(p)%p(mySource)%previousDotState2(:,c), &
sourceState(p)%p(mySource)%previousDotState (:,c) &
- sourceState(p)%p(mySource)%previousDotState2(:,c))
add pheno+ module in
2015-10-14 00:22:01 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if ( dot_prod22 > 0.0_pReal &
.and. ( dot_prod12 < 0.0_pReal &
.or. dot_product(sourceState(p)%p(mySource)%dotState(:,c), &
sourceState(p)%p(mySource)%previousDotState(:,c)) < 0.0_pReal) ) then
sourceStateDamper = 0.75_pReal + 0.25_pReal * tanh(2.0_pReal + 4.0_pReal * dot_prod12 / dot_prod22)
else
sourceStateDamper = 1.0_pReal
endif
! --- get residui ---
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
sourceStateResiduum(1:mySizeSourceDotState,mySource) = &
sourceState(p)%p(mySource)%state(1:mySizeSourceDotState,c) &
- sourceState(p)%p(mySource)%subState0(1:mySizeSourceDotState,c) &
- ( sourceState(p)%p(mySource)%dotState(1:mySizeSourceDotState,c) * sourceStateDamper &
+ sourceState(p)%p(mySource)%previousDotState(1:mySizeSourceDotState,c) &
* (1.0_pReal - sourceStateDamper)) * crystallite_subdt(g,i,e)
! --- correct state with residuum ---
tempSourceState(1:mySizeSourceDotState,mySource) = &
sourceState(p)%p(mySource)%state(1:mySizeSourceDotState,c) &
- sourceStateResiduum(1:mySizeSourceDotState,mySource) ! need to copy to local variable, since we cant flush a pointer in openmp
! --- store corrected dotState --- (cannot do this before state update, because not sure how to flush pointers in openmp)
sourceState(p)%p(mySource)%dotState(:,c) = &
sourceState(p)%p(mySource)%dotState(:,c) * sourceStateDamper &
+ sourceState(p)%p(mySource)%previousDotState(:,c) &
* (1.0_pReal - sourceStateDamper)
enddo
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((e == debug_e .and. i == debug_i .and. g == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
write(6,'(a,i8,1x,i2,1x,i3,/)') '<< CRYST >> update state at el ip g ',e,i,g
write(6,'(a,f6.1,/)') '<< CRYST >> plasticstatedamper ',plasticStatedamper
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> plastic state residuum',&
abs(plasticStateResiduum(1:mySizePlasticDotState))
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> abstol dotstate',plasticState(p)%aTolState(1:mySizePlasticDotState)
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> reltol dotstate',rTol_crystalliteState* &
abs(tempPlasticState(1:mySizePlasticDotState))
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> new state',tempPlasticState(1:mySizePlasticDotState)
more work on the new state
2014-05-27 20:16:03 +05:30
endif
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
#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
doxygen comments
2013-02-22 04:38:36 +05:30
! --- converged ? ---
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
converged = all( abs(plasticStateResiduum(1:mySizePlasticDotState)) < &
plasticState(p)%aTolState(1:mySizePlasticDotState) &
.or. abs(plasticStateResiduum(1:mySizePlasticDotState)) < &
rTol_crystalliteState * abs(tempPlasticState(1:mySizePlasticDotState)))
add pheno+ module in
2015-10-14 00:22:01 +05:30
do mySource = 1_pInt, phase_Nsources(p)
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
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
converged = converged .and. &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
all( abs(sourceStateResiduum(1:mySizeSourceDotState,mySource)) < &
sourceState(p)%p(mySource)%aTolState(1:mySizeSourceDotState) &
.or. abs(sourceStateResiduum(1:mySizeSourceDotState,mySource)) < &
rTol_crystalliteState * abs(tempSourceState(1:mySizeSourceDotState,mySource)))
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
enddo
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (converged) crystallite_converged(g,i,e) = .true. ! ... converged per definition
plasticState(p)%state(1:mySizePlasticDotState,c) = &
tempPlasticState(1:mySizePlasticDotState)
do mySource = 1_pInt, phase_Nsources(p)
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
sourceState(p)%p(mySource)%state(1:mySizeSourceDotState,c) = &
tempSourceState(1:mySizeSourceDotState,mySource)
enddo
endif
enddo; enddo; enddo
!$OMP ENDDO
! --- STATE JUMP ---
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$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
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt) &
write(6,'(a,i8,a,i2)') '<< CRYST >> ', count(crystallite_converged(:,:,:)), &
' grains converged after state integration #', NiterationState
#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
! --- NON-LOCAL CONVERGENCE CHECK ---
if (.not. singleRun) then ! if not requesting Integration of just a single IP
if (any(.not. crystallite_converged .and. .not. crystallite_localPlasticity)) & ! any non-local not yet converged (or broken)...
crystallite_converged = crystallite_converged .and. crystallite_localPlasticity ! ...restart all non-local as not converged
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
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#ifdef DEBUG
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt) then
write(6,'(a,i8,a)') '<< CRYST >> ', count(crystallite_converged(:,:,:)), &
' grains converged after non-local check'
write(6,'(a,i8,a,i2,/)') '<< CRYST >> ', count(crystallite_todo(:,:,:)), &
' grains todo after state integration #', NiterationState
endif
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#endif
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 ---
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
doneWithIntegration = .true.
elemLoop: do e = eIter(1),eIter(2)
do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
doneWithIntegration = .false.
exit elemLoop
endif
enddo; enddo
enddo elemLoop
enddo crystalliteLooping
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
!--------------------------------------------------------------------------------------------------
removed temperature integration and corresponding data structures and debugging options temperature is stored in crystallite, but homogeneous on one IP (not an component (grain) quantity and an input value parsed in by the BVP solver. introduced heat, a component (grain) quantity which is homogenized before returned to the heat transfer solver. went ahead with removal of dummy functions in homogenization and constitutive, this time mainly reduced function signatures to reflect actually needed quantities.
2013-10-16 18:34:59 +05:30
!> @brief integrate stress, and state with 1st order explicit Euler method
doxygen comments
2013-02-22 04:38:36 +05:30
!--------------------------------------------------------------------------------------------------
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
subroutine integrateStateEuler()
using IEEE_is_NAN and IEEE_quiet_NaN instead of hand-written solution, will not work for gfortran < 5
2017-05-04 04:02:44 +05:30
use, intrinsic :: &
IEEE_arithmetic
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
#ifdef DEBUG
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
use debug, only: &
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
debug_e, &
debug_i, &
debug_g, &
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
debug_level, &
debug_crystallite, &
debug_levelBasic, &
debug_levelExtensive, &
Deleted performance profiling code
2018-05-17 20:03:35 +05:30
debug_levelSelective
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#endif
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
use mesh, only: &
mesh_element, &
mesh_NcpElems
use material, only: &
more work on the new state
2014-05-27 20:16:03 +05:30
plasticState, &
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
sourceState, &
exchanged mappingConstitutive array for more descriptive "phaseAt and phasememberAt"
2016-01-15 05:49:44 +05:30
phaseAt, phasememberAt, &
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
phase_Nsources, &
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
homogenization_Ngrains
use constitutive, only: &
more work on the new state
2014-05-27 20:16:03 +05:30
constitutive_collectDotState, &
moved mappings for new state layout to material
2014-06-23 00:28:29 +05:30
constitutive_microstructure
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
doxygen comments
2013-02-22 04:38:36 +05:30
implicit none
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
integer(pInt) :: &
e, & ! element index in element loop
i, & ! integration point index in ip loop
g, & ! grain index in grain loop
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
p, & ! phase loop
modified state integrators to perform staggered integration of plasticity damage and thermal internal variables
2014-06-25 04:51:25 +05:30
c, &
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
mySource, &
modified state integrators to perform staggered integration of plasticity damage and thermal internal variables
2014-06-25 04:51:25 +05:30
mySizePlasticDotState, &
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
mySizeSourceDotState
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
integer(pInt), dimension(2) :: &
eIter ! bounds for element iteration
integer(pInt), dimension(2,mesh_NcpElems) :: &
iIter, & ! bounds for ip iteration
gIter ! bounds for grain iteration
logical :: &
name conflict
2015-09-10 03:22:00 +05:30
NaN, &
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
singleRun ! flag indicating computation for single (g,i,e) triple
simplified logic for single run, now not checking for extra arguments but simply checking if start and end elem/IP are the same
2013-04-26 18:53:36 +05:30
eIter = FEsolving_execElem(1:2)
do e = eIter(1),eIter(2)
iIter(1:2,e) = FEsolving_execIP(1:2,e)
gIter(1:2,e) = [ 1_pInt,homogenization_Ngrains(mesh_element(3,e))]
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
fixed out of range bug introduced with singleRun changes in rev. 2356
2013-04-30 17:44:07 +05:30
singleRun = (eIter(1) == eIter(2) .and. iIter(1,eIter(1)) == iIter(2,eIter(2)))
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)
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)
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
!$OMP 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
add pheno+ module in
2015-10-14 00:22:01 +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
doxygen comments
2013-02-22 04:38:36 +05:30
! --- STATE JUMP ---
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
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP DO
name conflict
2015-09-10 03:22:00 +05:30
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
crystallite_todo(g,i,e) = stateJump(g,i,e)
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP FLUSH(crystallite_todo)
removed time syncing
2019-01-14 17:26:46 +05:30
if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local...
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP CRITICAL (checkTodo)
name conflict
2015-09-10 03:22:00 +05:30
crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP END CRITICAL (checkTodo)
endif
endif
enddo; enddo; enddo
!$OMP 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
more work on the new state
2014-05-27 20:16:03 +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
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
call update_dependentState
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
!$OMP PARALLEL
doxygen comments
2013-02-22 04:38:36 +05:30
! --- STRESS INTEGRATION ---
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
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
don't use unnecessarily long names
2018-09-20 09:57:53 +05:30
crystallite_todo(g,i,e) = integrateStress(g,i,e)
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP FLUSH(crystallite_todo)
removed time syncing
2019-01-14 17:26:46 +05:30
if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local...
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP CRITICAL (checkTodo)
crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped
!$OMP END CRITICAL (checkTodo)
endif
endif
Depending on which state integrator one uses for the stiffness calculation, the initial state has to be chosen accordingly: e.g. for FPI choose last converged state, but for explicit RK choose converged state from start of increment (in case of explicit euler no state integration at all, but only stress integration). For this purpose we also need to remember "Fe" which now follows the cutbacking procedure as it is used for "Fp".
2011-11-04 18:14:50 +05:30
enddo; enddo; enddo
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP 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
doxygen comments
2013-02-22 04:38:36 +05:30
! --- SET CONVERGENCE FLAG ---
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
doxygen comments
2013-02-22 04:38:36 +05:30
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
cut out more now useless pieces
2018-07-04 01:42:25 +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 definitionem
doxygen comments
2013-02-22 04:38:36 +05:30
enddo; enddo; enddo
!$OMP 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
doxygen comments
2013-02-22 04:38:36 +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
doxygen comments
2013-02-22 04:38:36 +05:30
! --- CHECK NON-LOCAL 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
if (.not. singleRun) then ! if not requesting Integration of just a single IP
removed time syncing
2019-01-14 17:26:46 +05:30
if (any(.not. crystallite_converged .and. .not. crystallite_localPlasticity) ) & ! any non-local not yet converged (or broken)...
doxygen comments
2013-02-22 04:38:36 +05:30
crystallite_converged = crystallite_converged .and. crystallite_localPlasticity ! ...restart all non-local as not converged
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
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
!--------------------------------------------------------------------------------------------------
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
subroutine integrateStateAdaptiveEuler()
using IEEE_is_NAN and IEEE_quiet_NaN instead of hand-written solution, will not work for gfortran < 5
2017-05-04 04:02:44 +05:30
use, intrinsic :: &
IEEE_arithmetic
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
#ifdef DEBUG
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
use debug, only: &
corrected buggy previous commit
2015-04-21 20:03:38 +05:30
debug_e, &
debug_i, &
debug_g, &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
debug_level, &
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
debug_crystallite, &
debug_levelBasic, &
debug_levelExtensive, &
Deleted performance profiling code
2018-05-17 20:03:35 +05:30
debug_levelSelective
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#endif
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
use numerics, only: &
rTol_crystalliteState
use mesh, only: &
mesh_element, &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
mesh_NcpElems, &
mesh_maxNips
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
use material, only: &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
homogenization_Ngrains, &
more work on the new state
2014-05-27 20:16:03 +05:30
plasticState, &
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
sourceState, &
exchanged mappingConstitutive array for more descriptive "phaseAt and phasememberAt"
2016-01-15 05:49:44 +05:30
phaseAt, phasememberAt, &
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
phase_Nsources, &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
homogenization_maxNgrains
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
use constitutive, only: &
more work on the new state
2014-05-27 20:16:03 +05:30
constitutive_collectDotState, &
constitutive_microstructure, &
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
constitutive_plasticity_maxSizeDotState, &
constitutive_source_maxSizeDotState
changed to new state, please report bugs to Luv or Martin
2014-07-02 17:57:39 +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
implicit none
introduced plain mode (no ping pong) again and added test for MSC.MArc 2012 as prove that its working
2013-05-17 23:22:46 +05:30
integer(pInt) :: &
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
s, & ! state index
modified state integrators to perform staggered integration of plasticity damage and thermal internal variables
2014-06-25 04:51:25 +05:30
p, &
c, &
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
mySource, &
Heavy bug fixing: Butcher tableau multiplication in RK45 integration scheme, state updating of adaptive Euler, wrong assignment to evolution rate instead of state in fixed-point iterator. Exchanged possibly non thread-safe matmul in RKCK45, simplified RK45 integration step 1--3 vs 4 handling. (Thanks to Pratheek, Luv, and Chen for their help!)
2014-09-03 01:16:52 +05:30
mySizePlasticDotState, & ! size of dot states
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
mySizeSourceDotState
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
integer(pInt), dimension(2) :: &
eIter ! bounds for element iteration
integer(pInt), dimension(2,mesh_NcpElems) :: &
iIter, & ! bounds for ip iteration
gIter ! bounds for grain iteration
real(pReal), dimension(constitutive_plasticity_maxSizeDotState, &
homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems) :: &
plasticStateResiduum, & ! residuum from evolution in micrstructure
relPlasticStateResiduum ! relative residuum from evolution in microstructure
real(pReal), dimension(constitutive_source_maxSizeDotState,&
maxval(phase_Nsources), &
homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems) :: &
sourceStateResiduum, & ! residuum from evolution in micrstructure
relSourceStateResiduum ! relative residuum from evolution in microstructure
disabled multiphysics output because homogenization output was not accessible anymore
2014-07-09 19:04:33 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
logical :: &
converged, &
NaN, &
singleRun ! flag indicating computation for single (g,i,e) triple
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
! --- LOOP ITERATOR FOR ELEMENT, GRAIN, IP ---
eIter = FEsolving_execElem(1:2)
do e = eIter(1),eIter(2)
iIter(1:2,e) = FEsolving_execIP(1:2,e)
gIter(1:2,e) = [ 1_pInt,homogenization_Ngrains(mesh_element(3,e))]
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
singleRun = (eIter(1) == eIter(2) .and. iIter(1,eIter(1)) == iIter(2,eIter(2)))
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
disabled multiphysics output because homogenization output was not accessible anymore
2014-07-09 19:04:33 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
plasticStateResiduum = 0.0_pReal
relPlasticStateResiduum = 0.0_pReal
sourceStateResiduum = 0.0_pReal
relSourceStateResiduum = 0.0_pReal
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP PARALLEL
more work on the new state
2014-05-27 20:16:03 +05:30
disabled multiphysics output because homogenization output was not accessible anymore
2014-07-09 19:04:33 +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 (EULER INTEGRATION) ---
removed temperature integration and corresponding data structures and debugging options temperature is stored in crystallite, but homogeneous on one IP (not an component (grain) quantity and an input value parsed in by the BVP solver. introduced heat, a component (grain) quantity which is homogenized before returned to the heat transfer solver. went ahead with removal of dummy functions in homogenization and constitutive, this time mainly reduced function signatures to reflect actually needed quantities.
2013-10-16 18:34:59 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP DO PRIVATE(mySizePlasticDotState,mySizeSourceDotState,p,c)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e)) then
exchanged mappingConstitutive array for more descriptive "phaseAt and phasememberAt"
2016-01-15 05:49:44 +05:30
p = phaseAt(g,i,e)
c = phasememberAt(g,i,e)
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
mySizePlasticDotState = plasticState(p)%sizeDotState
plasticStateResiduum(1:mySizePlasticDotState,g,i,e) = &
- 0.5_pReal &
* plasticState(p)%dotstate(1:mySizePlasticDotState,c) &
* crystallite_subdt(g,i,e) ! contribution to absolute residuum in state
plasticState(p)%state (1:mySizePlasticDotState,c) = &
plasticState(p)%state (1:mySizePlasticDotState,c) &
+ plasticState(p)%dotstate(1:mySizePlasticDotState,c) &
* crystallite_subdt(g,i,e)
do mySource = 1_pInt, phase_Nsources(p)
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
sourceStateResiduum(1:mySizeSourceDotState,mySource,g,i,e) = &
- 0.5_pReal &
* sourceState(p)%p(mySource)%dotstate(1:mySizeSourceDotState,c) &
* crystallite_subdt(g,i,e) ! contribution to absolute residuum in state
sourceState(p)%p(mySource)%state (1:mySizeSourceDotState,c) = &
sourceState(p)%p(mySource)%state (1:mySizeSourceDotState,c) &
+ sourceState(p)%p(mySource)%dotstate(1:mySizeSourceDotState,c) &
* crystallite_subdt(g,i,e)
enddo
endif
enddo; enddo; enddo
!$OMP ENDDO
! --- STATE JUMP ---
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e)) then
crystallite_todo(g,i,e) = stateJump(g,i,e)
!$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)
remved bug in phenopower law state initialization and a minor change in Makefile
2014-07-07 19:51:58 +05:30
endif
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
!$OMP ENDDO
!$OMP END PARALLEL
added vacancy state integrators
2014-12-08 16:15:12 +05:30
more "building block" separation
2019-01-24 03:32:21 +05:30
call update_dependentState
add pheno+ module in
2015-10-14 00:22:01 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- STRESS INTEGRATION (EULER INTEGRATION) ---
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP PARALLEL DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e)) then
crystallite_todo(g,i,e) = integrateStress(g,i,e)
!$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
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
call update_dotState(1.0_pReal)
more "building block" separation
2019-01-24 03:32:21 +05:30
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
!$OMP PARALLEL
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP DO PRIVATE(p,c,NaN)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e)
c = phasememberAt(g,i,e)
NaN = any(IEEE_is_NaN(plasticState(p)%dotState(:,c)))
do mySource = 1_pInt, phase_Nsources(p)
NaN = NaN .or. any(IEEE_is_NaN(sourceState(p)%p(mySource)%dotState(:,c)))
enddo
if (NaN) then ! NaN occured in any dotState
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)
else ! if broken local...
crystallite_todo(g,i,e) = .false. ! ... skip this one next time
endif
endif
endif
enddo; enddo; enddo
!$OMP 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
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- ERROR ESTIMATE FOR STATE (HEUN METHOD) ---
!$OMP SINGLE
relPlasticStateResiduum = 0.0_pReal
relSourceStateResiduum = 0.0_pReal
!$OMP END SINGLE
!$OMP DO PRIVATE(mySizePlasticDotState,mySizeSourceDotState,converged,p,c,s)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e)
c = phasememberAt(g,i,e)
! --- contribution of heun step to absolute residui ---
mySizePlasticDotState = plasticState(p)%sizeDotState
plasticStateResiduum(1:mySizePlasticDotState,g,i,e) = &
plasticStateResiduum(1:mySizePlasticDotState,g,i,e) &
+ 0.5_pReal * plasticState(p)%dotState(:,c) &
* crystallite_subdt(g,i,e) ! contribution to absolute residuum in state
do mySource = 1_pInt, phase_Nsources(p)
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
sourceStateResiduum(1:mySizeSourceDotState,mySource,g,i,e) = &
sourceStateResiduum(1:mySizeSourceDotState,mySource,g,i,e) &
+ 0.5_pReal * sourceState(p)%p(mySource)%dotState(:,c) &
* crystallite_subdt(g,i,e) ! contribution to absolute residuum in state
enddo
! --- relative residui ---
forall (s = 1_pInt:mySizePlasticDotState, abs(plasticState(p)%dotState(s,c)) > 0.0_pReal) &
relPlasticStateResiduum(s,g,i,e) = &
plasticStateResiduum(s,g,i,e) / plasticState(p)%dotState(s,c)
do mySource = 1_pInt, phase_Nsources(p)
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
forall (s = 1_pInt:mySizeSourceDotState,abs(sourceState(p)%p(mySource)%dotState(s,c)) > 0.0_pReal) &
relSourceStateResiduum(s,mySource,g,i,e) = &
sourceStateResiduum(s,mySource,g,i,e) / sourceState(p)%p(mySource)%dotState(s,c)
add pheno+ module in
2015-10-14 00:22:01 +05:30
enddo
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
added more debugging messages and switched to DEBUG as flag
2017-09-30 04:02:07 +05:30
#ifdef DEBUG
increased debug threshold for some outputs (decreased verbosity level at "basic")
2014-05-15 18:31:54 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
.and. ((e == debug_e .and. i == debug_i .and. g == debug_g)&
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
write(6,'(a,i8,1x,i2,1x,i3,/)') '<< CRYST >> updateState at el ip g ',e,i,g
write(6,'(a,/,(12x,12(f12.1,1x)),/)') '<< CRYST >> absolute residuum tolerance', &
plasticStateResiduum(1:mySizePlasticDotState,g,i,e) / plasticState(p)%aTolState(1:mySizePlasticDotState)
write(6,'(a,/,(12x,12(f12.1,1x)),/)') '<< CRYST >> relative residuum tolerance', &
relPlasticStateResiduum(1:mySizePlasticDotState,g,i,e) / rTol_crystalliteState
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> dotState', plasticState(p)%dotState(1:mySizePlasticDotState,c) &
- 2.0_pReal * plasticStateResiduum(1:mySizePlasticDotState,g,i,e) / crystallite_subdt(g,i,e) ! calculate former dotstate from higher order solution and state residuum
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> new state', plasticState(p)%state(1:mySizePlasticDotState,c)
doxygen comments
2013-02-22 04:38:36 +05:30
endif
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
#endif
more work on the new state
2014-05-27 20:16:03 +05:30
! --- converged ? ---
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
converged = all(abs(relPlasticStateResiduum(1:mySizePlasticDotState,g,i,e)) < &
rTol_crystalliteState .or. &
abs(plasticStateResiduum(1:mySizePlasticDotState,g,i,e)) < &
plasticState(p)%aTolState(1:mySizePlasticDotState))
add pheno+ module in
2015-10-14 00:22:01 +05:30
do mySource = 1_pInt, phase_Nsources(p)
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
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
converged = converged .and. &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
all(abs(relSourceStateResiduum(1:mySizeSourceDotState,mySource,g,i,e)) < &
rTol_crystalliteState .or. &
abs(sourceStateResiduum(1:mySizeSourceDotState,mySource,g,i,e)) < &
sourceState(p)%p(mySource)%aTolState(1:mySizeSourceDotState))
add pheno+ module in
2015-10-14 00:22:01 +05:30
enddo
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (converged) crystallite_converged(g,i,e) = .true. ! ... converged per definitionem
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- NONLOCAL CONVERGENCE CHECK ---
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#ifdef DEBUG
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt) &
write(6,'(a,i8,a,i2,/)') '<< CRYST >> ', count(crystallite_converged(:,:,:)), ' grains converged'
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
#endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if ((.not. singleRun) .and. any(.not. crystallite_converged .and. .not. crystallite_localPlasticity)) & ! any non-local not yet converged (or broken)...
crystallite_converged = crystallite_converged .and. crystallite_localPlasticity ! ...restart all non-local as not converged
use (import) debug variables only when needed
2018-09-20 09:54:03 +05:30
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
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
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
subroutine integrateStateRK4()
using IEEE_is_NAN and IEEE_quiet_NaN instead of hand-written solution, will not work for gfortran < 5
2017-05-04 04:02:44 +05:30
use, intrinsic :: &
IEEE_arithmetic
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
#ifdef DEBUG
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
use debug, only: &
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
debug_e, &
debug_i, &
debug_g, &
corrected buggy previous commit
2015-04-21 20:03:38 +05:30
debug_level, &
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
debug_crystallite, &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
debug_levelBasic, &
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
debug_levelExtensive, &
included changes to correct cutback issue of spectral solver
2018-02-16 20:06:18 +05:30
debug_levelSelective
#endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
use mesh, only: &
mesh_element, &
mesh_NcpElems
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
use material, only: &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
homogenization_Ngrains, &
more work on the new state
2014-05-27 20:16:03 +05:30
plasticState, &
changed up handling of delta states, and some bug fixes
2015-06-01 21:32:27 +05:30
sourceState, &
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
phase_Nsources, &
phaseAt, phasememberAt
use config, only: &
material_Nphase
use constitutive, only: &
constitutive_collectDotState, &
constitutive_microstructure
added funtions to get averged properties at integration points.
2014-09-10 14:07:12 +05:30
implicit none
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
integer(pInt) :: e, & ! element index in element loop
i, & ! integration point index in ip loop
g, & ! grain index in grain loop
p, & ! phase loop
c, &
n, &
mySource, &
mySizePlasticDotState, &
mySizeSourceDotState
integer(pInt), dimension(2) :: eIter ! bounds for element iteration
integer(pInt), dimension(2,mesh_NcpElems) :: iIter, & ! bounds for ip iteration
gIter ! bounds for grain iteration
logical :: NaN, &
singleRun ! flag indicating computation for single (g,i,e) triple
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
eIter = FEsolving_execElem(1:2)
do e = eIter(1),eIter(2)
iIter(1:2,e) = FEsolving_execIP(1:2,e)
gIter(1:2,e) = [ 1_pInt,homogenization_Ngrains(mesh_element(3,e))]
enddo
named 99 representation according to other constitutive laws
2014-11-06 17:19:37 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
singleRun = (eIter(1) == eIter(2) .and. iIter(1,eIter(1)) == iIter(2,eIter(2)))
named 99 representation according to other constitutive laws
2014-11-06 17:19:37 +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
! initialize dotState
if (.not. singleRun) then
do p = 1_pInt, material_Nphase
plasticState(p)%RK4dotState = 0.0_pReal
do mySource = 1_pInt, phase_Nsources(p)
sourceState(p)%p(mySource)%RK4dotState = 0.0_pReal
enddo
enddo
else
e = eIter(1)
i = iIter(1,e)
do g = gIter(1,e), gIter(2,e)
plasticState(phaseAt(g,i,e))%RK4dotState(:,phasememberAt(g,i,e)) = 0.0_pReal
do mySource = 1_pInt, phase_Nsources(phaseAt(g,i,e))
sourceState(phaseAt(g,i,e))%p(mySource)%RK4dotState(:,phasememberAt(g,i,e)) = 0.0_pReal
enddo
enddo
endif
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!--------------------------------------------------------------------------------------------------
! --- SECOND TO FOURTH RUNGE KUTTA STEP PLUS FINAL INTEGRATION ---
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
do n = 1_pInt,4_pInt
! --- state update ---
named some loops and made BUTCHER TABLEAUs parameters in crystallite_integrateStateRKCK45
2013-04-29 16:47:30 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP PARALLEL
!$OMP DO PRIVATE(p,c)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e)
c = phasememberAt(g,i,e)
plasticState(p)%RK4dotState(:,c) = plasticState(p)%RK4dotState(:,c) &
+ weight(n)*plasticState(p)%dotState(:,c)
do mySource = 1_pInt, phase_Nsources(p)
sourceState(p)%p(mySource)%RK4dotState(:,c) = sourceState(p)%p(mySource)%RK4dotState(:,c) &
+ weight(n)*sourceState(p)%p(mySource)%dotState(:,c)
enddo
endif
enddo; enddo; enddo
!$OMP ENDDO
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
!$OMP END PARALLEL
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))
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
!$OMP PARALLEL
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- state jump ---
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
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e)) then
crystallite_todo(g,i,e) = stateJump(g,i,e)
!$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 ENDDO
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
!$OMP END PARALLEL
call update_dependentState
!$OMP 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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- stress integration ---
cleaned up stress integrations and some documentation
2014-11-13 18:23:20 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e)) then
crystallite_todo(g,i,e) = integrateStress(g,i,e,timeStepFraction(n)) ! fraction of original times step
!$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 ENDDO
common function to update dot state
2019-01-23 10:44:19 +05:30
!$OMP END PARALLEL
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
common function to update dot state
2019-01-23 10:44:19 +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
common function to update dot state
2019-01-23 10:44:19 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo
! --- SET CONVERGENCE FLAG ---
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
crystallite_converged(g,i,e) = crystallite_todo(g,i,e) .or. crystallite_converged(g,i,e) ! if still "to do" then converged per definitionem
enddo; enddo; enddo
! --- CHECK NONLOCAL CONVERGENCE ---
if (.not. singleRun) then ! if not requesting Integration of just a single IP
if (any(.not. crystallite_converged .and. .not. crystallite_localPlasticity)) then ! any non-local not yet converged (or broken)...
crystallite_converged = crystallite_converged .and. crystallite_localPlasticity ! ...restart all non-local as not converged
overhaul of kinematics from FeFp to FeFiFp decomposition
2014-11-01 00:33:08 +05:30
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
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")
!--------------------------------------------------------------------------------------------------
subroutine integrateStateRKCK45()
use, intrinsic :: &
IEEE_arithmetic
added more debugging messages and switched to DEBUG as flag
2017-09-30 04:02:07 +05:30
#ifdef DEBUG
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
use debug, only: &
debug_e, &
debug_i, &
debug_g, &
debug_level, &
debug_crystallite, &
debug_levelBasic, &
debug_levelExtensive, &
debug_levelSelective
* 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
#endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
use numerics, only: &
rTol_crystalliteState
use mesh, only: &
mesh_element, &
mesh_NcpElems, &
mesh_maxNips
use material, only: &
homogenization_Ngrains, &
plasticState, &
sourceState, &
phase_Nsources, &
phaseAt, phasememberAt, &
homogenization_maxNgrains
use constitutive, only: &
constitutive_collectDotState, &
constitutive_plasticity_maxSizeDotState, &
constitutive_source_maxSizeDotState, &
constitutive_microstructure
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
implicit none
real(pReal), dimension(5,5), parameter :: &
A = reshape([&
.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], &
[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 - &
[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)
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)
integer(pInt) :: &
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
s, & ! state index
n, &
p, &
cc, &
mySource, &
mySizePlasticDotState, & ! size of dot States
mySizeSourceDotState
integer(pInt), dimension(2) :: &
eIter ! bounds for element iteration
integer(pInt), dimension(2,mesh_NcpElems) :: &
iIter, & ! bounds for ip iteration
gIter ! bounds for grain iteration
real(pReal), dimension(constitutive_plasticity_maxSizeDotState, &
homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems) :: &
plasticStateResiduum, & ! residuum from evolution in microstructure
relPlasticStateResiduum ! relative residuum from evolution in microstructure
real(pReal), dimension(constitutive_source_maxSizeDotState, &
maxval(phase_Nsources), &
homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems) :: &
sourceStateResiduum, & ! residuum from evolution in microstructure
relSourceStateResiduum ! relative residuum from evolution in microstructure
logical :: &
NaN, &
singleRun ! flag indicating computation for single (g,i,e) triple
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
eIter = FEsolving_execElem(1:2)
added more debugging messages and switched to DEBUG as flag
2017-09-30 04:02:07 +05:30
#ifdef DEBUG
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt) &
write(6,'(a,1x,i1)') '<< CRYST >> Runge--Kutta step',1
overhaul of kinematics from FeFp to FeFiFp decomposition
2014-11-01 00:33:08 +05:30
#endif
add pheno+ module in
2015-10-14 00:22:01 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- LOOP ITERATOR FOR ELEMENT, GRAIN, IP ---
do e = eIter(1),eIter(2)
iIter(1:2,e) = FEsolving_execIP(1:2,e)
gIter(1:2,e) = [ 1_pInt,homogenization_Ngrains(mesh_element(3,e))]
enddo
overhaul of kinematics from FeFp to FeFiFp decomposition
2014-11-01 00:33:08 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
singleRun = (eIter(1) == eIter(2) .and. iIter(1,eIter(1)) == iIter(2,eIter(2)))
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
do stage = 1_pInt,5_pInt
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 ---
!$OMP PARALLEL
!$OMP DO PRIVATE(p,cc)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e)
cc = phasememberAt(g,i,e)
plasticState(p)%RKCK45dotState(stage,:,cc) = plasticState(p)%dotState(:,cc) ! store Runge-Kutta dotState
do mySource = 1_pInt, phase_Nsources(p)
sourceState(p)%p(mySource)%RKCK45dotState(stage,:,cc) = sourceState(p)%p(mySource)%dotState(:,cc)
enddo
question marks on those files
2017-12-14 05:48:45 +05:30
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo; enddo; enddo
!$OMP 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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP DO PRIVATE(p,cc,n)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e)
cc = 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
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)
do mySource = 1_pInt, phase_Nsources(p)
sourceState(p)%p(mySource)%dotState(:,cc) = A(1,stage) * sourceState(p)%p(mySource)%RKCK45dotState(1,:,cc)
enddo
do n = 2_pInt, stage
plasticState(p)%dotState(:,cc) = &
plasticState(p)%dotState(:,cc) + A(n,stage) * plasticState(p)%RKCK45dotState(n,:,cc)
do mySource = 1_pInt, phase_Nsources(p)
sourceState(p)%p(mySource)%dotState(:,cc) = &
sourceState(p)%p(mySource)%dotState(:,cc) + A(n,stage) * sourceState(p)%p(mySource)%RKCK45dotState(n,:,cc)
enddo
enddo
endif
enddo; enddo; enddo
!$OMP ENDDO
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +05:30
!$OMP END PARALLEL
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?
!$OMP PARALLEL
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- state jump ---
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e)) then
crystallite_todo(g,i,e) = stateJump(g,i,e)
!$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)
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
endif
overhaul of kinematics from FeFp to FeFiFp decomposition
2014-11-01 00:33:08 +05:30
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo; enddo; enddo
!$OMP ENDDO
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
call update_dependentState
!$OMP PARALLEL
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- stress integration ---
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e)) then
crystallite_todo(g,i,e) = integrateStress(g,i,e,C(stage)) ! fraction of original time step
!$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 ENDDO
!$OMP END PARALLEL
common function to update dot state
2019-01-23 10:44:19 +05:30
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 ---
relPlasticStateResiduum = 0.0_pReal
relSourceStateResiduum = 0.0_pReal
!$OMP PARALLEL
!$OMP DO PRIVATE(p,cc)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e)
cc = phasememberAt(g,i,e)
plasticState(p)%RKCK45dotState(6,:,cc) = plasticState (p)%dotState(:,cc) ! store Runge-Kutta dotState
do mySource = 1_pInt, phase_Nsources(p)
sourceState(p)%p(mySource)%RKCK45dotState(6,:,cc) = sourceState(p)%p(mySource)%dotState(:,cc) ! store Runge-Kutta dotState
enddo
switched dNeq and dEq to relative tolerance, removed single precision (makes things complicated and was never used anyway)
2016-03-21 03:50:58 +05:30
endif
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
enddo; enddo; enddo
!$OMP ENDDO
improved debug statements
2015-08-05 02:56:22 +05:30
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP DO PRIVATE(mySizePlasticDotState,mySizeSourceDotState,p,cc)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e)
cc = phasememberAt(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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- absolute residuum in state ---
mySizePlasticDotState = plasticState(p)%sizeDotState
plasticStateResiduum(1:mySizePlasticDotState,g,i,e) = &
matmul(transpose(plasticState(p)%RKCK45dotState(1:6,1:mySizePlasticDotState,cc)),DB) &
* crystallite_subdt(g,i,e)
do mySource = 1_pInt, phase_Nsources(p)
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
sourceStateResiduum(1:mySizeSourceDotState,mySource,g,i,e) = &
matmul(transpose(sourceState(p)%p(mySource)%RKCK45dotState(1:6,1:mySizeSourceDotState,cc)),DB) &
* crystallite_subdt(g,i,e)
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- dot state ---
plasticState(p)%dotState(:,cc) = &
matmul(transpose(plasticState(p)%RKCK45dotState(1:6,1:mySizePlasticDotState,cc)), B)
do mySource = 1_pInt, phase_Nsources(p)
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
sourceState(p)%p(mySource)%dotState(:,cc) = &
matmul(transpose(sourceState(p)%p(mySource)%RKCK45dotState(1:6,1:mySizeSourceDotState,cc)),B)
enddo
endif
enddo; enddo; enddo
!$OMP ENDDO
using a function for state update avoids a lot of code repetition
2019-01-23 16:21:43 +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
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)
!$OMP PARALLEL
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!$OMP DO PRIVATE(mySizePlasticDotState,mySizeSourceDotState,p,cc,s)
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e)
cc = phasememberAt(g,i,e)
mySizePlasticDotState = plasticState(p)%sizeDotState
forall (s = 1_pInt:mySizePlasticDotState, abs(plasticState(p)%state(s,cc)) > 0.0_pReal) &
relPlasticStateResiduum(s,g,i,e) = &
plasticStateResiduum(s,g,i,e) / plasticState(p)%state(s,cc)
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
do mySource = 1_pInt, phase_Nsources(p)
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
forall (s = 1_pInt:mySizeSourceDotState,abs(sourceState(p)%p(mySource)%state(s,cc)) > 0.0_pReal) &
relSourceStateResiduum(s,mySource,g,i,e) = &
sourceStateResiduum(s,mySource,g,i,e) / sourceState(p)%p(mySource)%state(s,cc)
enddo
crystallite_todo(g,i,e) = all(abs(relPlasticStateResiduum(1:mySizePlasticDotState,g,i,e)) < &
rTol_crystalliteState .or. &
abs(plasticStateResiduum(1:mySizePlasticDotState,g,i,e)) < &
plasticState(p)%aTolState(1:mySizePlasticDotState))
do mySource = 1_pInt, phase_Nsources(p)
mySizeSourceDotState = sourceState(p)%p(mySource)%sizeDotState
crystallite_todo(g,i,e) = crystallite_todo(g,i,e) .and. &
all(abs(relSourceStateResiduum(1:mySizeSourceDotState,mySource,g,i,e)) < &
rTol_crystalliteState .or. &
abs(sourceStateResiduum(1:mySizeSourceDotState,mySource,g,i,e)) < &
sourceState(p)%p(mySource)%aTolState(1:mySizeSourceDotState))
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
added more debugging messages and switched to DEBUG as flag
2017-09-30 04:02:07 +05:30
#ifdef DEBUG
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt&
.and. ((e == debug_e .and. i == debug_i .and. g == debug_g)&
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
write(6,'(a,i8,1x,i3,1x,i3,/)') '<< CRYST >> updateState at el ip ipc ',e,i,g
write(6,'(a,/,(12x,12(f12.1,1x)),/)') '<< CRYST >> absolute residuum tolerance', &
plasticStateResiduum(1:mySizePlasticDotState,g,i,e) / plasticState(p)%aTolState(1:mySizePlasticDotState)
write(6,'(a,/,(12x,12(f12.1,1x)),/)') '<< CRYST >> relative residuum tolerance', &
relPlasticStateResiduum(1:mySizePlasticDotState,g,i,e) / rTol_crystalliteState
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> dotState', &
plasticState(p)%dotState(1:mySizePlasticDotState,cc)
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> new state', &
plasticState(p)%state(1:mySizePlasticDotState,cc)
endif
* 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
#endif
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
!$OMP 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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- STATE JUMP ---
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
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e)) then
crystallite_todo(g,i,e) = stateJump(g,i,e)
!$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)
doxygen comments
2013-02-22 04:38:36 +05:30
endif
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
!$OMP ENDDO
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +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
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
call update_dependentState
new subroutine crystallite_orientations: calculates the orientation matrix and the euler angles and - in case of single grain ips - the misorientation for each neighboring ip. The calculation of the misorientation is now done once in crystallite init and at the end of every FE increment. This saves a lot of time compared to doing it in dotState for every crystallite subinc.
2009-12-18 21:16:33 +05:30
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
!$OMP PARALLEL
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
!--------------------------------------------------------------------------------------------------
! --- FINAL STRESS INTEGRATION STEP IF RESIDUUM BELOW TOLERANCE ---
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e)) then
crystallite_todo(g,i,e) = integrateStress(g,i,e)
!$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 ENDDO
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
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
!--------------------------------------------------------------------------------------------------
! --- SET CONVERGENCE FLAG ---
!$OMP DO
do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains
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
!$OMP ENDDO
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
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
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
internal (private) functions at the end ordered state integrators according to their id
2019-01-15 08:57:57 +05:30
! --- nonlocal convergence check ---
#ifdef DEBUG
if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt) &
write(6,'(a,i8,a,i2,/)') '<< CRYST >> ', count(crystallite_converged(:,:,:)), ' grains converged' ! if not requesting Integration of just a single IP
#endif
if ((.not. singleRun) .and. any(.not. crystallite_converged .and. .not. crystallite_localPlasticity)) & ! any non-local not yet converged (or broken)...
crystallite_converged = crystallite_converged .and. crystallite_localPlasticity ! ...restart all non-local as not converged
modified state integrators to perform staggered integration of plasticity damage and thermal internal variables
2014-06-25 04:51:25 +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 integrateStateRKCK45
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
function for update of dependent state not introduced everywhere
2019-01-23 22:34:19 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief tbd
!--------------------------------------------------------------------------------------------------
subroutine update_dependentState()
use material, only: &
plasticState
use constitutive, only: &
constitutive_dependentState => constitutive_microstructure
implicit none
integer(pInt) :: e, & ! element index in element loop
i, & ! integration point index in ip loop
g ! grain index in grain loop
!$OMP PARALLEL
!$OMP DO
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(mesh_element(3,e))
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) &
call constitutive_dependentState(crystallite_orientation, &
crystallite_Fe(1:3,1:3,g,i,e), &
crystallite_Fp(1:3,1:3,g,i,e), &
g, i, e)
enddo; enddo; enddo
!$OMP ENDDO
!$OMP END PARALLEL
end subroutine update_dependentState
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)
use material, only: &
plasticState, &
sourceState, &
phase_Nsources, &
phaseAt, phasememberAt
implicit none
real(pReal), intent(in) :: &
timeFraction
integer(pInt) :: &
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)
do g = 1,homogenization_Ngrains(mesh_element(3,e))
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
p = phaseAt(g,i,e); c = phasememberAt(g,i,e)
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
do s = 1_pInt, phase_Nsources(p)
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)
use, intrinsic :: &
IEEE_arithmetic
use material, only: &
plasticState, &
sourceState, &
phaseAt, phasememberAt, &
phase_Nsources
use constitutive, only: &
constitutive_collectDotState
implicit none
real(pReal), intent(in) :: &
timeFraction
integer(pInt) :: &
e, & !< element index in element loop
i, & !< integration point index in ip loop
g, & !< grain index in grain loop
p, &
c, &
s
logical :: &
NaN
!$OMP PARALLEL
!$OMP DO PRIVATE (p,c,NaN)
do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(mesh_element(3,e))
!$OMP FLUSH(crystallite_todo)
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
call constitutive_collectDotState(crystallite_Tstar_v(1:6,g,i,e), &
crystallite_Fe, &
crystallite_Fi(1:3,1:3,g,i,e), &
crystallite_Fp, &
crystallite_subdt(g,i,e)*timeFraction, crystallite_subFrac, g,i,e)
p = phaseAt(g,i,e); c = phasememberAt(g,i,e)
NaN = any(IEEE_is_NaN(plasticState(p)%dotState(:,c)))
do s = 1_pInt, phase_Nsources(p)
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)
if (.not. crystallite_localPlasticity(g,i,e)) then ! if broken is a local...
!$OMP CRITICAL (checkTodo)
crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals done (and broken)
!$OMP END CRITICAL (checkTodo)
endif
endif
endif
enddo; enddo; enddo
!$OMP ENDDO
!$OMP END PARALLEL
end subroutine update_DotState
[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)
use, intrinsic :: &
IEEE_arithmetic
use prec, only: &
dNeq0
#ifdef DEBUG
use debug, only: &
debug_e, &
debug_i, &
debug_g, &
debug_level, &
debug_crystallite, &
debug_levelExtensive, &
debug_levelSelective
#endif
use material, only: &
plasticState, &
sourceState, &
phase_Nsources, &
phaseAt, phasememberAt
use constitutive, only: &
constitutive_collectDeltaState
implicit none
integer(pInt), intent(in):: &
el, & ! element index
ip, & ! integration point index
ipc ! grain index
integer(pInt) :: &
c, &
p, &
mySource, &
myOffsetPlasticDeltaState, &
myOffsetSourceDeltaState, &
mySizePlasticDeltaState, &
mySizeSourceDeltaState
c = phasememberAt(ipc,ip,el)
p = phaseAt(ipc,ip,el)
call constitutive_collectDeltaState(crystallite_Tstar_v(1:6,ipc,ip,el), &
crystallite_Fe(1:3,1:3,ipc,ip,el), &
crystallite_Fi(1:3,1:3,ipc,ip,el), &
ipc,ip,el)
myOffsetPlasticDeltaState = plasticState(p)%offsetDeltaState
mySizePlasticDeltaState = plasticState(p)%sizeDeltaState
if( any(IEEE_is_NaN(plasticState(p)%deltaState(1:mySizePlasticDeltaState,c)))) then ! NaN occured in deltaState
stateJump = .false.
return
endif
plasticState(p)%state(myOffsetPlasticDeltaState + 1_pInt : &
myOffsetPlasticDeltaState + mySizePlasticDeltaState,c) = &
plasticState(p)%state(myOffsetPlasticDeltaState + 1_pInt : &
myOffsetPlasticDeltaState + mySizePlasticDeltaState,c) + &
plasticState(p)%deltaState(1:mySizePlasticDeltaState,c)
do mySource = 1_pInt, phase_Nsources(p)
myOffsetSourceDeltaState = sourceState(p)%p(mySource)%offsetDeltaState
mySizeSourceDeltaState = sourceState(p)%p(mySource)%sizeDeltaState
if (any(IEEE_is_NaN(sourceState(p)%p(mySource)%deltaState(1:mySizeSourceDeltaState,c)))) then ! NaN occured in deltaState
stateJump = .false.
return
endif
sourceState(p)%p(mySource)%state(myOffsetSourceDeltaState + 1_pInt : &
myOffsetSourceDeltaState + mySizeSourceDeltaState,c) = &
sourceState(p)%p(mySource)%state(myOffsetSourceDeltaState + 1_pInt : &
myOffsetSourceDeltaState + mySizeSourceDeltaState,c) + &
sourceState(p)%p(mySource)%deltaState(1:mySizeSourceDeltaState,c)
enddo
#ifdef DEBUG
if (any(dNeq0(plasticState(p)%deltaState(1:mySizePlasticDeltaState,c))) &
.and. iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt &
.and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) &
.or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then
write(6,'(a,i8,1x,i2,1x,i3, /)') '<< CRYST >> update state at el ip ipc ',el,ip,ipc
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> deltaState', plasticState(p)%deltaState(1:mySizePlasticDeltaState,c)
write(6,'(a,/,(12x,12(e12.5,1x)),/)') '<< CRYST >> new state', &
plasticState(p)%state(myOffsetPlasticDeltaState + 1_pInt : &
myOffsetPlasticDeltaState + mySizePlasticDeltaState,c)
endif
#endif
stateJump = .true.
end function stateJump
doxygen comments
2013-02-22 04:38:36 +05:30
end module crystallite