Merge branch 'development' of magit1.mpie.de:damask/DAMASK into miscImprovements
This commit is contained in:
Martin Diehl
commit
48a9aea9ff
11
CONFIG
11
CONFIG
|
|
@ -1,10 +1,11 @@
|
||||||
# "set"-syntax needed only for tcsh (but works with bash and zsh)
|
# "set"-syntax needed only for tcsh (but works with bash and zsh)
|
||||||
# DAMASK_ROOT will be expanded
|
# DAMASK_ROOT will be expanded
|
||||||
set DAMASK_BIN=${DAMASK_ROOT}/bin
|
|
||||||
|
|
||||||
set DAMASK_NUM_THREADS = 4
|
set DAMASK_BIN = ${DAMASK_ROOT}/bin
|
||||||
|
|
||||||
set MSC_ROOT=/opt/MSC
|
set DAMASK_NUM_THREADS = 4
|
||||||
set MARC_VERSION=2015
|
|
||||||
|
|
||||||
set ABAQUS_VERSION=6.14-5
|
set MSC_ROOT = /opt/MSC
|
||||||
|
set MARC_VERSION = 2015
|
||||||
|
|
||||||
|
set ABAQUS_VERSION = 6.14-5
|
||||||
|
|
|
||||||
|
|
@ -23,25 +23,26 @@ endif
|
||||||
|
|
||||||
# according to http://software.intel.com/en-us/forums/topic/501500
|
# according to http://software.intel.com/en-us/forums/topic/501500
|
||||||
# this seems to make sense for the stack size
|
# this seems to make sense for the stack size
|
||||||
set FREE=`which free`
|
if ( `which free` != "free: Command not found." ) then
|
||||||
set FREE=''
|
|
||||||
if ( "x$FREE" != "x" ) then
|
|
||||||
set freeMem=`free -k | grep -E '(Mem|Speicher):' | awk '{print $4;}'`
|
set freeMem=`free -k | grep -E '(Mem|Speicher):' | awk '{print $4;}'`
|
||||||
set heap=`expr $freeMem / 2`
|
|
||||||
set stack=`expr $freeMem / $DAMASK_NUM_THREADS / 2`
|
set stack=`expr $freeMem / $DAMASK_NUM_THREADS / 2`
|
||||||
|
set heap=` expr $freeMem / 2`
|
||||||
# http://superuser.com/questions/220059/what-parameters-has-ulimit
|
# http://superuser.com/questions/220059/what-parameters-has-ulimit
|
||||||
limit stacksize $stack # maximum stack size (kB)
|
limit stacksize $stack # maximum stack size (kB)
|
||||||
limit datasize $heap # maximum heap size (kB)
|
limit datasize $heap # maximum heap size (kB)
|
||||||
|
endif
|
||||||
|
if ( `limit | grep memoryuse` != "" ) then
|
||||||
|
limit memoryuse unlimited # maximum physical memory size
|
||||||
|
endif
|
||||||
|
if ( `limit | grep vmemoryuse` != "" ) then
|
||||||
|
limit vmemoryuse unlimited # maximum virtual memory size
|
||||||
endif
|
endif
|
||||||
limit coredumpsize 2000 # core file size (512-byte blocks)
|
|
||||||
limit vmemoryuse unlimited # maximum virtual memory size
|
|
||||||
limit memoryuse unlimited # maximum physical memory size
|
|
||||||
|
|
||||||
# disable output in case of scp
|
# disable output in case of scp
|
||||||
if($?prompt) then
|
if ( $?prompt ) then
|
||||||
echo ''
|
echo ''
|
||||||
echo Düsseldorf Advanced Materials Simulation Kit - DAMASK
|
echo Düsseldorf Advanced Materials Simulation Kit --- DAMASK
|
||||||
echo Max-Planck-Institut für Eisenforschung, Düsseldorf
|
echo Max-Planck-Institut für Eisenforschung GmbH, Düsseldorf
|
||||||
echo https://damask.mpie.de
|
echo https://damask.mpie.de
|
||||||
echo
|
echo
|
||||||
echo Using environment with ...
|
echo Using environment with ...
|
||||||
|
|
@ -59,8 +60,8 @@ if($?prompt) then
|
||||||
echo "MSC.Marc/Mentat $MSC_ROOT"
|
echo "MSC.Marc/Mentat $MSC_ROOT"
|
||||||
endif
|
endif
|
||||||
echo
|
echo
|
||||||
echo `limit stacksize`
|
|
||||||
echo `limit datasize`
|
echo `limit datasize`
|
||||||
|
echo `limit stacksize`
|
||||||
endif
|
endif
|
||||||
|
|
||||||
setenv DAMASK_NUM_THREADS $DAMASK_NUM_THREADS
|
setenv DAMASK_NUM_THREADS $DAMASK_NUM_THREADS
|
||||||
|
|
|
||||||
|
|
@ -5,8 +5,8 @@
|
||||||
if [ "$OSTYPE" == "linux-gnu" ] || [ "$OSTYPE" == 'linux' ]; then
|
if [ "$OSTYPE" == "linux-gnu" ] || [ "$OSTYPE" == 'linux' ]; then
|
||||||
DAMASK_ROOT=$(python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$(dirname $BASH_SOURCE)")
|
DAMASK_ROOT=$(python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$(dirname $BASH_SOURCE)")
|
||||||
else
|
else
|
||||||
[[ "${BASH_SOURCE::1}" == "/" ]] && BASE="" || BASE="`pwd`/"
|
[[ "${BASH_SOURCE::1}" == "/" ]] && BASE="" || BASE="$(pwd)/"
|
||||||
STAT=$(stat "`dirname $BASE$BASH_SOURCE`")
|
STAT=$(stat "$(dirname $BASE$BASH_SOURCE)")
|
||||||
DAMASK_ROOT=${STAT##* }
|
DAMASK_ROOT=${STAT##* }
|
||||||
fi
|
fi
|
||||||
|
|
||||||
|
|
@ -17,16 +17,16 @@ set() {
|
||||||
source $DAMASK_ROOT/CONFIG
|
source $DAMASK_ROOT/CONFIG
|
||||||
unset -f set
|
unset -f set
|
||||||
|
|
||||||
# if DAMASK_BIN is present and not in $PATH, add it
|
# add DAMASK_BIN if present but not yet in $PATH
|
||||||
if [[ "x$DAMASK_BIN" != "x" && ! `echo ":$PATH:" | grep $DAMASK_BIN:` ]]; then
|
if [[ "x$DAMASK_BIN" != "x" && ! $(echo ":$PATH:" | grep $DAMASK_BIN:) ]]; then
|
||||||
export PATH=$DAMASK_BIN:$PATH
|
export PATH=$DAMASK_BIN:$PATH
|
||||||
fi
|
fi
|
||||||
|
|
||||||
SOLVER=`which DAMASK_spectral 2>/dev/null`
|
SOLVER=$(which DAMASK_spectral 2>/dev/null)
|
||||||
if [ "x$SOLVER" == "x" ]; then
|
if [ "x$SOLVER" == "x" ]; then
|
||||||
SOLVER='Not found!'
|
SOLVER='Not found!'
|
||||||
fi
|
fi
|
||||||
PROCESSING=`which postResults 2>/dev/null`
|
PROCESSING=$(which postResults 2>/dev/null)
|
||||||
if [ "x$PROCESSING" == "x" ]; then
|
if [ "x$PROCESSING" == "x" ]; then
|
||||||
PROCESSING='Not found!'
|
PROCESSING='Not found!'
|
||||||
fi
|
fi
|
||||||
|
|
@ -36,25 +36,21 @@ fi
|
||||||
|
|
||||||
# according to http://software.intel.com/en-us/forums/topic/501500
|
# according to http://software.intel.com/en-us/forums/topic/501500
|
||||||
# this seems to make sense for the stack size
|
# this seems to make sense for the stack size
|
||||||
FREE=`which free 2>/dev/null`
|
FREE=$(which free 2>/dev/null)
|
||||||
if [ "x$FREE" != "x" ]; then
|
if [ "x$FREE" != "x" ]; then
|
||||||
freeMem=`free -k | grep -E '(Mem|Speicher):' | awk '{print $4;}'`
|
freeMem=$(free -k | grep -E '(Mem|Speicher):' | awk '{print $4;}')
|
||||||
heap=`expr $freeMem / 2`
|
|
||||||
stack=`expr $freeMem / $DAMASK_NUM_THREADS / 2`
|
|
||||||
|
|
||||||
# http://superuser.com/questions/220059/what-parameters-has-ulimit
|
# http://superuser.com/questions/220059/what-parameters-has-ulimit
|
||||||
ulimit -s $stack 2>/dev/null # maximum stack size (kB)
|
ulimit -s $(expr $freeMem / $DAMASK_NUM_THREADS / 2) 2>/dev/null # maximum stack size (kB)
|
||||||
ulimit -d $heap 2>/dev/null # maximum heap size (kB)
|
ulimit -d $(expr $freeMem / 2) 2>/dev/null # maximum heap size (kB)
|
||||||
fi
|
fi
|
||||||
ulimit -c 2000 2>/dev/null # core file size (512-byte blocks)
|
|
||||||
ulimit -v unlimited 2>/dev/null # maximum virtual memory size
|
ulimit -v unlimited 2>/dev/null # maximum virtual memory size
|
||||||
ulimit -m unlimited 2>/dev/null # maximum physical memory size
|
ulimit -m unlimited 2>/dev/null # maximum physical memory size
|
||||||
|
|
||||||
# disable output in case of scp
|
# disable output in case of scp
|
||||||
if [ ! -z "$PS1" ]; then
|
if [ ! -z "$PS1" ]; then
|
||||||
echo
|
echo
|
||||||
echo Düsseldorf Advanced Materials Simulation Kit - DAMASK
|
echo Düsseldorf Advanced Materials Simulation Kit --- DAMASK
|
||||||
echo Max-Planck-Institut für Eisenforschung, Düsseldorf
|
echo Max-Planck-Institut für Eisenforschung GmbH, Düsseldorf
|
||||||
echo https://damask.mpie.de
|
echo https://damask.mpie.de
|
||||||
echo
|
echo
|
||||||
echo Using environment with ...
|
echo Using environment with ...
|
||||||
|
|
@ -64,14 +60,29 @@ if [ ! -z "$PS1" ]; then
|
||||||
echo "Multithreading DAMASK_NUM_THREADS=$DAMASK_NUM_THREADS"
|
echo "Multithreading DAMASK_NUM_THREADS=$DAMASK_NUM_THREADS"
|
||||||
if [ "x$PETSC_DIR" != "x" ]; then
|
if [ "x$PETSC_DIR" != "x" ]; then
|
||||||
echo "PETSc location $PETSC_DIR"
|
echo "PETSc location $PETSC_DIR"
|
||||||
[[ `python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$PETSC_DIR"` == $PETSC_DIR ]] \
|
[[ $(python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$PETSC_DIR") == $PETSC_DIR ]] \
|
||||||
|| echo " ~~> "`python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$PETSC_DIR"`
|
|| echo " ~~> "$(python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$PETSC_DIR")
|
||||||
fi
|
fi
|
||||||
[[ "x$PETSC_ARCH" != "x" ]] && echo "PETSc architecture $PETSC_ARCH"
|
[[ "x$PETSC_ARCH" == "x" ]] \
|
||||||
|
|| echo "PETSc architecture $PETSC_ARCH"
|
||||||
echo "MSC.Marc/Mentat $MSC_ROOT"
|
echo "MSC.Marc/Mentat $MSC_ROOT"
|
||||||
echo
|
echo
|
||||||
echo -n "heap size/MiB "; echo "`ulimit -d`/1024" | bc
|
echo -n "heap size "
|
||||||
echo -n "stack size/MiB "; echo "`ulimit -s`/1024" | bc
|
[[ "$(ulimit -d)" == "unlimited" ]] \
|
||||||
|
&& echo "unlimited" \
|
||||||
|
|| echo $(python -c \
|
||||||
|
"import math; \
|
||||||
|
size=$(( 1024*$(ulimit -d) )); \
|
||||||
|
print('{:.4g} {}'.format(size / (1 << ((int(math.log(size,2) / 10) if size else 0) * 10)), \
|
||||||
|
['bytes','KiB','MiB','GiB','TiB','EiB','ZiB'][int(math.log(size,2) / 10) if size else 0]))")
|
||||||
|
echo -n "stack size "
|
||||||
|
[[ "$(ulimit -s)" == "unlimited" ]] \
|
||||||
|
&& echo "unlimited" \
|
||||||
|
|| echo $(python -c \
|
||||||
|
"import math; \
|
||||||
|
size=$(( 1024*$(ulimit -s) )); \
|
||||||
|
print('{:.4g} {}'.format(size / (1 << ((int(math.log(size,2) / 10) if size else 0) * 10)), \
|
||||||
|
['bytes','KiB','MiB','GiB','TiB','EiB','ZiB'][int(math.log(size,2) / 10) if size else 0]))")
|
||||||
fi
|
fi
|
||||||
|
|
||||||
export DAMASK_NUM_THREADS
|
export DAMASK_NUM_THREADS
|
||||||
|
|
|
||||||
|
|
@ -1,12 +1,7 @@
|
||||||
# sets up an environment for DAMASK on zsh
|
# sets up an environment for DAMASK on zsh
|
||||||
# usage: source DAMASK_env.zsh
|
# usage: source DAMASK_env.zsh
|
||||||
|
|
||||||
|
DAMASK_ROOT=${0:a:h}
|
||||||
if [ "$OSTYPE" = "linux-gnu" ] || [ "$OSTYPE" = 'linux' ]; then
|
|
||||||
DAMASK_ROOT=$(readlink -f "`dirname ${(%):-%N}`")
|
|
||||||
else
|
|
||||||
print 'Not done yet'
|
|
||||||
fi
|
|
||||||
|
|
||||||
# defining set() allows to source the same file for tcsh and zsh, with and without space around =
|
# defining set() allows to source the same file for tcsh and zsh, with and without space around =
|
||||||
set() {
|
set() {
|
||||||
|
|
@ -15,45 +10,35 @@ set() {
|
||||||
source $DAMASK_ROOT/CONFIG
|
source $DAMASK_ROOT/CONFIG
|
||||||
unset -f set
|
unset -f set
|
||||||
|
|
||||||
# if DAMASK_BIN is present and not in $PATH, add it
|
# add DAMASK_BIN if present but not yet in $PATH
|
||||||
MATCH=`echo ":$PATH:" | grep $DAMASK_BIN:`
|
MATCH=`echo ":$PATH:" | grep $DAMASK_BIN:`
|
||||||
if [[ ( "x$DAMASK_BIN" != "x" ) && ( "x$MATCH" = "x" ) ]]; then
|
if [[ ( "x$DAMASK_BIN" != "x" ) && ( "x$MATCH" = "x" ) ]]; then
|
||||||
export PATH=$DAMASK_BIN:$PATH
|
export PATH=$DAMASK_BIN:$PATH
|
||||||
fi
|
fi
|
||||||
|
|
||||||
SOLVER=`which DAMASK_spectral 2>/dev/null`
|
SOLVER=`which DAMASK_spectral 2>/dev/null`
|
||||||
if [ "x$SOLVER" = "x" ]; then
|
|
||||||
export SOLVER='Not found!'
|
|
||||||
fi
|
|
||||||
PROCESSING=`which postResults 2>/dev/null`
|
PROCESSING=`which postResults 2>/dev/null`
|
||||||
if [ "x$PROCESSING" = "x" ]; then
|
|
||||||
export PROCESSING='Not found!'
|
|
||||||
fi
|
|
||||||
if [ "x$DAMASK_NUM_THREADS" = "x" ]; then
|
if [ "x$DAMASK_NUM_THREADS" = "x" ]; then
|
||||||
DAMASK_NUM_THREADS=1
|
DAMASK_NUM_THREADS=1
|
||||||
fi
|
fi
|
||||||
|
|
||||||
# according to http://software.intel.com/en-us/forums/topic/501500
|
# according to http://software.intel.com/en-us/forums/topic/501500
|
||||||
# this seems to make sense for the stack size
|
# this seems to make sense for the stack size
|
||||||
FREE=`which free 2>/dev/null`
|
if [ "`which free 2>/dev/null`" != "free not found" ]; then
|
||||||
if [ "x$FREE" != "x" ]; then
|
|
||||||
freeMem=`free -k | grep -E '(Mem|Speicher):' | awk '{print $4;}'`
|
freeMem=`free -k | grep -E '(Mem|Speicher):' | awk '{print $4;}'`
|
||||||
heap=`expr $freeMem / 2`
|
|
||||||
stack=`expr $freeMem / 2`
|
|
||||||
|
|
||||||
# http://superuser.com/questions/220059/what-parameters-has-ulimit
|
# http://superuser.com/questions/220059/what-parameters-has-ulimit
|
||||||
ulimit -s $stack 2>/dev/null # maximum stack size (kB)
|
ulimit -s `expr $freeMem / $DAMASK_NUM_THREADS / 2` 2>/dev/null # maximum stack size (kB)
|
||||||
ulimit -d $heap 2>/dev/null # maximum heap size (kB)
|
ulimit -d `expr $freeMem / 2` 2>/dev/null # maximum heap size (kB)
|
||||||
fi
|
fi
|
||||||
ulimit -c 2000 2>/dev/null # core file size (512-byte blocks)
|
|
||||||
ulimit -v unlimited 2>/dev/null # maximum virtual memory size
|
ulimit -v unlimited 2>/dev/null # maximum virtual memory size
|
||||||
ulimit -m unlimited 2>/dev/null # maximum physical memory size
|
ulimit -m unlimited 2>/dev/null # maximum physical memory size
|
||||||
|
|
||||||
# disable output in case of scp
|
# disable output in case of scp
|
||||||
if [ ! -z "$PS1" ]; then
|
if [ ! -z "$PS1" ]; then
|
||||||
echo
|
echo
|
||||||
echo Düsseldorf Advanced Materials Simulation Kit - DAMASK
|
echo Düsseldorf Advanced Materials Simulation Kit --- DAMASK
|
||||||
echo Max-Planck-Institut für Eisenforschung, Düsseldorf
|
echo Max-Planck-Institut für Eisenforschung GmbH, Düsseldorf
|
||||||
echo https://damask.mpie.de
|
echo https://damask.mpie.de
|
||||||
echo
|
echo
|
||||||
echo Using environment with ...
|
echo Using environment with ...
|
||||||
|
|
@ -63,13 +48,29 @@ if [ ! -z "$PS1" ]; then
|
||||||
echo "Multithreading DAMASK_NUM_THREADS=$DAMASK_NUM_THREADS"
|
echo "Multithreading DAMASK_NUM_THREADS=$DAMASK_NUM_THREADS"
|
||||||
if [ "x$PETSC_DIR" != "x" ]; then
|
if [ "x$PETSC_DIR" != "x" ]; then
|
||||||
echo "PETSc location $PETSC_DIR"
|
echo "PETSc location $PETSC_DIR"
|
||||||
[[ `readlink -f $PETSC_DIR` == $PETSC_DIR ]] || echo " ~~> "`readlink -f $PETSC_DIR`
|
[[ $(python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$PETSC_DIR") == $PETSC_DIR ]] \
|
||||||
|
|| echo " ~~> "$(python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$PETSC_DIR")
|
||||||
fi
|
fi
|
||||||
[[ "x$PETSC_ARCH" != "x" ]] && echo "PETSc architecture $PETSC_ARCH"
|
[[ "x$PETSC_ARCH" == "x" ]] \
|
||||||
|
|| echo "PETSc architecture $PETSC_ARCH"
|
||||||
echo "MSC.Marc/Mentat $MSC_ROOT"
|
echo "MSC.Marc/Mentat $MSC_ROOT"
|
||||||
echo
|
echo
|
||||||
echo -n "heap size/MiB "; echo "`ulimit -d`/1024" | bc
|
echo -n "heap size "
|
||||||
echo -n "stack size/MiB "; echo "`ulimit -s`/1024" | bc
|
[[ "$(ulimit -d)" == "unlimited" ]] \
|
||||||
|
&& echo "unlimited" \
|
||||||
|
|| echo $(python -c \
|
||||||
|
"import math; \
|
||||||
|
size=$(( 1024*$(ulimit -d) )); \
|
||||||
|
print('{:.4g} {}'.format(size / (1 << ((int(math.log(size,2) / 10) if size else 0) * 10)), \
|
||||||
|
['bytes','KiB','MiB','GiB','TiB','EiB','ZiB'][int(math.log(size,2) / 10) if size else 0]))")
|
||||||
|
echo -n "stack size "
|
||||||
|
[[ "$(ulimit -s)" == "unlimited" ]] \
|
||||||
|
&& echo "unlimited" \
|
||||||
|
|| echo $(python -c \
|
||||||
|
"import math; \
|
||||||
|
size=$(( 1024*$(ulimit -s) )); \
|
||||||
|
print('{:.4g} {}'.format(size / (1 << ((int(math.log(size,2) / 10) if size else 0) * 10)), \
|
||||||
|
['bytes','KiB','MiB','GiB','TiB','EiB','ZiB'][int(math.log(size,2) / 10) if size else 0]))")
|
||||||
fi
|
fi
|
||||||
|
|
||||||
export DAMASK_NUM_THREADS
|
export DAMASK_NUM_THREADS
|
||||||
|
|
|
||||||
|
|
@ -156,11 +156,9 @@ program DAMASK_spectral
|
||||||
!--------------------------------------------------------------------------------------------------
|
!--------------------------------------------------------------------------------------------------
|
||||||
! init DAMASK (all modules)
|
! init DAMASK (all modules)
|
||||||
call CPFEM_initAll(el = 1_pInt, ip = 1_pInt)
|
call CPFEM_initAll(el = 1_pInt, ip = 1_pInt)
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- DAMASK_spectral init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- DAMASK_spectral init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
!--------------------------------------------------------------------------------------------------
|
!--------------------------------------------------------------------------------------------------
|
||||||
! initialize field solver information
|
! initialize field solver information
|
||||||
|
|
@ -199,14 +197,14 @@ program DAMASK_spectral
|
||||||
allocate(loadCases(i)%ID(nActiveFields))
|
allocate(loadCases(i)%ID(nActiveFields))
|
||||||
field = 1
|
field = 1
|
||||||
loadCases(i)%ID(field) = FIELD_MECH_ID ! mechanical active by default
|
loadCases(i)%ID(field) = FIELD_MECH_ID ! mechanical active by default
|
||||||
if (any(thermal_type == THERMAL_conduction_ID)) then ! thermal field active
|
thermalActive: if (any(thermal_type == THERMAL_conduction_ID)) then
|
||||||
field = field + 1
|
field = field + 1
|
||||||
loadCases(i)%ID(field) = FIELD_THERMAL_ID
|
loadCases(i)%ID(field) = FIELD_THERMAL_ID
|
||||||
endif
|
endif thermalActive
|
||||||
if (any(damage_type == DAMAGE_nonlocal_ID)) then ! damage field active
|
damageActive: if (any(damage_type == DAMAGE_nonlocal_ID)) then
|
||||||
field = field + 1
|
field = field + 1
|
||||||
loadCases(i)%ID(field) = FIELD_DAMAGE_ID
|
loadCases(i)%ID(field) = FIELD_DAMAGE_ID
|
||||||
endif
|
endif damageActive
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
!--------------------------------------------------------------------------------------------------
|
!--------------------------------------------------------------------------------------------------
|
||||||
|
|
|
||||||
|
|
@ -518,8 +518,7 @@ subroutine crystallite_stressAndItsTangent(updateJaco)
|
||||||
nCryst, &
|
nCryst, &
|
||||||
numerics_integrator, &
|
numerics_integrator, &
|
||||||
numerics_integrationMode, &
|
numerics_integrationMode, &
|
||||||
numerics_timeSyncing, &
|
numerics_timeSyncing
|
||||||
analyticJaco
|
|
||||||
use debug, only: &
|
use debug, only: &
|
||||||
debug_level, &
|
debug_level, &
|
||||||
debug_crystallite, &
|
debug_crystallite, &
|
||||||
|
|
@ -582,23 +581,6 @@ subroutine crystallite_stressAndItsTangent(updateJaco)
|
||||||
invFp, & ! inverse of the plastic deformation gradient
|
invFp, & ! inverse of the plastic deformation gradient
|
||||||
Fe_guess, & ! guess for elastic deformation gradient
|
Fe_guess, & ! guess for elastic deformation gradient
|
||||||
Tstar ! 2nd Piola-Kirchhoff stress tensor
|
Tstar ! 2nd Piola-Kirchhoff stress tensor
|
||||||
real(pReal), allocatable, dimension(:,:,:,:,:,:,:) :: &
|
|
||||||
dPdF_perturbation1, &
|
|
||||||
dPdF_perturbation2
|
|
||||||
real(pReal), allocatable, dimension(:,:,:,:,:) :: &
|
|
||||||
F_backup, &
|
|
||||||
Fp_backup, &
|
|
||||||
InvFp_backup, &
|
|
||||||
Fi_backup, &
|
|
||||||
InvFi_backup, &
|
|
||||||
Fe_backup, &
|
|
||||||
Lp_backup, &
|
|
||||||
Li_backup, &
|
|
||||||
P_backup
|
|
||||||
real(pReal), allocatable, dimension(:,:,:,:) :: &
|
|
||||||
Tstar_v_backup
|
|
||||||
logical, allocatable, dimension(:,:,:) :: &
|
|
||||||
convergenceFlag_backup
|
|
||||||
integer(pInt) :: &
|
integer(pInt) :: &
|
||||||
NiterationCrystallite, & ! number of iterations in crystallite loop
|
NiterationCrystallite, & ! number of iterations in crystallite loop
|
||||||
c, & !< counter in integration point component loop
|
c, & !< counter in integration point component loop
|
||||||
|
|
@ -1137,371 +1119,115 @@ subroutine crystallite_stressAndItsTangent(updateJaco)
|
||||||
! --+>> STIFFNESS CALCULATION <<+--
|
! --+>> STIFFNESS CALCULATION <<+--
|
||||||
|
|
||||||
computeJacobian: if(updateJaco) then
|
computeJacobian: if(updateJaco) then
|
||||||
jacobianMethod: if (analyticJaco) then
|
!$OMP PARALLEL DO PRIVATE(dSdF,dSdFe,dSdFi,dLpdS,dLpdFi,dFpinvdF,dLidS,dLidFi,dFidS,&
|
||||||
|
!$OMP rhs_3333,lhs_3333,temp_99,temp_33,temp_3333,myNcomponents,error)
|
||||||
|
elementLooping6: do e = FEsolving_execElem(1),FEsolving_execElem(2)
|
||||||
|
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
|
||||||
|
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) ! iterate over IPs of this element to be processed
|
||||||
|
do c = 1_pInt,myNcomponents
|
||||||
|
call constitutive_TandItsTangent(temp_33,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
|
||||||
|
|
||||||
! --- ANALYTIC JACOBIAN ---
|
call constitutive_LiAndItsTangent(temp_33,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
|
||||||
|
temp_33 = 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(temp_33,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(temp_33,dLidS(1:3,1:3,o,p))
|
||||||
|
enddo; enddo
|
||||||
|
call math_invert(9_pInt,math_Plain3333to99(lhs_3333),temp_99,error)
|
||||||
|
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
|
||||||
|
dFidS = math_mul3333xx3333(math_Plain99to3333(temp_99),rhs_3333)
|
||||||
|
endif
|
||||||
|
dLidS = math_mul3333xx3333(dLidFi,dFidS) + dLidS
|
||||||
|
endif
|
||||||
|
|
||||||
!$OMP PARALLEL DO PRIVATE(dSdF,dSdFe,dSdFi,dLpdS,dLpdFi,dFpinvdF,dLidS,dLidFi,dFidS,&
|
call constitutive_LpAndItsTangent(temp_33,dLpdS,dLpdFi,crystallite_Tstar_v(1:6,c,i,e), &
|
||||||
!$OMP rhs_3333,lhs_3333,temp_99,temp_33,temp_3333,myNcomponents,error)
|
crystallite_Fi(1:3,1:3,c,i,e),c,i,e) ! call constitutive law to calculate Lp tangent in lattice configuration
|
||||||
elementLooping6: do e = FEsolving_execElem(1),FEsolving_execElem(2)
|
dLpdS = math_mul3333xx3333(dLpdFi,dFidS) + dLpdS
|
||||||
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
|
|
||||||
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) ! iterate over IPs of this element to be processed
|
|
||||||
do c = 1_pInt,myNcomponents
|
|
||||||
call constitutive_TandItsTangent(temp_33,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_LiAndItsTangent(temp_33,dLidS,dLidFi,crystallite_Tstar_v(1:6,c,i,e), &
|
temp_33 = math_transpose33(math_mul33x33(crystallite_invFp(1:3,1:3,c,i,e), &
|
||||||
crystallite_Fi(1:3,1:3,c,i,e), &
|
crystallite_invFi(1:3,1:3,c,i,e)))
|
||||||
c,i,e) ! call constitutive law to calculate Li tangent in lattice configuration
|
rhs_3333 = 0.0_pReal
|
||||||
if (sum(abs(dLidS)) < tol_math_check) then
|
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
||||||
dFidS = 0.0_pReal
|
rhs_3333(p,o,1:3,1:3) = math_mul33x33(dSdFe(p,o,1:3,1:3),temp_33)
|
||||||
else
|
|
||||||
temp_33 = 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(temp_33,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(temp_33,dLidS(1:3,1:3,o,p))
|
|
||||||
enddo; enddo
|
|
||||||
call math_invert(9_pInt,math_Plain3333to99(lhs_3333),temp_99,error)
|
|
||||||
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
|
|
||||||
dFidS = math_mul3333xx3333(math_Plain99to3333(temp_99),rhs_3333)
|
|
||||||
endif
|
|
||||||
dLidS = math_mul3333xx3333(dLidFi,dFidS) + dLidS
|
|
||||||
endif
|
|
||||||
|
|
||||||
call constitutive_LpAndItsTangent(temp_33,dLpdS,dLpdFi,crystallite_Tstar_v(1:6,c,i,e), &
|
temp_3333 = 0.0_pReal
|
||||||
crystallite_Fi(1:3,1:3,c,i,e),c,i,e) ! call constitutive law to calculate Lp tangent in lattice configuration
|
temp_33 = math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
|
||||||
dLpdS = math_mul3333xx3333(dLpdFi,dFidS) + dLpdS
|
math_inv33(crystallite_subFp0(1:3,1:3,c,i,e)))
|
||||||
|
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
||||||
|
temp_3333(1:3,1:3,p,o) = math_mul33x33(math_mul33x33(temp_33,dLpdS(1:3,1:3,p,o)), &
|
||||||
|
crystallite_invFi(1:3,1:3,c,i,e))
|
||||||
|
|
||||||
temp_33 = math_transpose33(math_mul33x33(crystallite_invFp(1:3,1:3,c,i,e), &
|
temp_33 = math_mul33x33(math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
|
||||||
crystallite_invFi(1:3,1:3,c,i,e)))
|
crystallite_invFp(1:3,1:3,c,i,e)), &
|
||||||
rhs_3333 = 0.0_pReal
|
math_inv33(crystallite_subFi0(1:3,1:3,c,i,e)))
|
||||||
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
||||||
rhs_3333(p,o,1:3,1:3) = math_mul33x33(dSdFe(p,o,1:3,1:3),temp_33)
|
temp_3333(1:3,1:3,p,o) = temp_3333(1:3,1:3,p,o) + math_mul33x33(temp_33,dLidS(1:3,1:3,p,o))
|
||||||
|
|
||||||
temp_3333 = 0.0_pReal
|
lhs_3333 = crystallite_subdt(c,i,e)*math_mul3333xx3333(dSdFe,temp_3333) + &
|
||||||
temp_33 = math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
|
math_mul3333xx3333(dSdFi,dFidS)
|
||||||
math_inv33(crystallite_subFp0(1:3,1:3,c,i,e)))
|
|
||||||
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
|
||||||
temp_3333(1:3,1:3,p,o) = math_mul33x33(math_mul33x33(temp_33,dLpdS(1:3,1:3,p,o)), &
|
|
||||||
crystallite_invFi(1:3,1:3,c,i,e))
|
|
||||||
|
|
||||||
temp_33 = math_mul33x33(math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
|
call math_invert(9_pInt,math_identity2nd(9_pInt)+math_Plain3333to99(lhs_3333),temp_99,error)
|
||||||
crystallite_invFp(1:3,1:3,c,i,e)), &
|
if (error) then
|
||||||
math_inv33(crystallite_subFi0(1:3,1:3,c,i,e)))
|
call IO_warning(warning_ID=600_pInt,el=e,ip=i,g=c, &
|
||||||
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
ext_msg='inversion error in analytic tangent calculation')
|
||||||
temp_3333(1:3,1:3,p,o) = temp_3333(1:3,1:3,p,o) + math_mul33x33(temp_33,dLidS(1:3,1:3,p,o))
|
dSdF = rhs_3333
|
||||||
|
else
|
||||||
|
dSdF = math_mul3333xx3333(math_Plain99to3333(temp_99),rhs_3333)
|
||||||
|
endif
|
||||||
|
|
||||||
lhs_3333 = crystallite_subdt(c,i,e)*math_mul3333xx3333(dSdFe,temp_3333) + &
|
dFpinvdF = 0.0_pReal
|
||||||
math_mul3333xx3333(dSdFi,dFidS)
|
temp_3333 = math_mul3333xx3333(dLpdS,dSdF)
|
||||||
|
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
||||||
|
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)))
|
||||||
|
|
||||||
call math_invert(9_pInt,math_identity2nd(9_pInt)+math_Plain3333to99(lhs_3333),temp_99,error)
|
crystallite_dPdF(1:3,1:3,1:3,1:3,c,i,e) = 0.0_pReal
|
||||||
if (error) then
|
temp_33 = math_mul33x33(crystallite_invFp(1:3,1:3,c,i,e), &
|
||||||
call IO_warning(warning_ID=600_pInt,el=e,ip=i,g=c, &
|
math_mul33x33(math_Mandel6to33(crystallite_Tstar_v(1:6,c,i,e)), &
|
||||||
ext_msg='inversion error in analytic tangent calculation')
|
math_transpose33(crystallite_invFp(1:3,1:3,c,i,e))))
|
||||||
dSdF = rhs_3333
|
forall(p=1_pInt:3_pInt) &
|
||||||
else
|
crystallite_dPdF(p,1:3,p,1:3,c,i,e) = math_transpose33(temp_33)
|
||||||
dSdF = math_mul3333xx3333(math_Plain99to3333(temp_99),rhs_3333)
|
|
||||||
endif
|
|
||||||
|
|
||||||
dFpinvdF = 0.0_pReal
|
temp_33 = math_mul33x33(math_Mandel6to33(crystallite_Tstar_v(1:6,c,i,e)), &
|
||||||
temp_3333 = math_mul3333xx3333(dLpdS,dSdF)
|
math_transpose33(crystallite_invFp(1:3,1:3,c,i,e)))
|
||||||
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
||||||
dFpinvdF(1:3,1:3,p,o) = -crystallite_subdt(c,i,e)* &
|
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_inv33(crystallite_subFp0(1:3,1:3,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)
|
||||||
math_mul33x33(temp_3333(1:3,1:3,p,o), &
|
|
||||||
crystallite_invFi(1:3,1:3,c,i,e)))
|
|
||||||
|
|
||||||
crystallite_dPdF(1:3,1:3,1:3,1:3,c,i,e) = 0.0_pReal
|
temp_33 = math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
|
||||||
temp_33 = math_mul33x33(crystallite_invFp(1:3,1:3,c,i,e), &
|
crystallite_invFp(1:3,1:3,c,i,e))
|
||||||
math_mul33x33(math_Mandel6to33(crystallite_Tstar_v(1:6,c,i,e)), &
|
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
||||||
math_transpose33(crystallite_invFp(1:3,1:3,c,i,e))))
|
crystallite_dPdF(1:3,1:3,p,o,c,i,e) = crystallite_dPdF(1:3,1:3,p,o,c,i,e) + &
|
||||||
forall(p=1_pInt:3_pInt) &
|
math_mul33x33(math_mul33x33(temp_33,dSdF(1:3,1:3,p,o)), &
|
||||||
crystallite_dPdF(p,1:3,p,1:3,c,i,e) = math_transpose33(temp_33)
|
math_transpose33(crystallite_invFp(1:3,1:3,c,i,e)))
|
||||||
|
|
||||||
temp_33 = math_mul33x33(math_Mandel6to33(crystallite_Tstar_v(1:6,c,i,e)), &
|
temp_33 = math_mul33x33(math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
|
||||||
math_transpose33(crystallite_invFp(1:3,1:3,c,i,e)))
|
crystallite_invFp(1:3,1:3,c,i,e)), &
|
||||||
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
math_Mandel6to33(crystallite_Tstar_v(1:6,c,i,e)))
|
||||||
crystallite_dPdF(1:3,1:3,p,o,c,i,e) = crystallite_dPdF(1:3,1:3,p,o,c,i,e) + &
|
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
||||||
math_mul33x33(math_mul33x33(crystallite_subF(1:3,1:3,c,i,e),dFpinvdF(1:3,1:3,p,o)),temp_33)
|
crystallite_dPdF(1:3,1:3,p,o,c,i,e) = crystallite_dPdF(1:3,1:3,p,o,c,i,e) + &
|
||||||
|
math_mul33x33(temp_33,math_transpose33(dFpinvdF(1:3,1:3,p,o)))
|
||||||
|
|
||||||
temp_33 = math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
|
enddo; enddo
|
||||||
crystallite_invFp(1:3,1:3,c,i,e))
|
enddo elementLooping6
|
||||||
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
!$OMP END PARALLEL DO
|
||||||
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(temp_33,dSdF(1:3,1:3,p,o)), &
|
|
||||||
math_transpose33(crystallite_invFp(1:3,1:3,c,i,e)))
|
|
||||||
|
|
||||||
temp_33 = math_mul33x33(math_mul33x33(crystallite_subF(1:3,1:3,c,i,e), &
|
|
||||||
crystallite_invFp(1:3,1:3,c,i,e)), &
|
|
||||||
math_Mandel6to33(crystallite_Tstar_v(1:6,c,i,e)))
|
|
||||||
forall(p=1_pInt:3_pInt, o=1_pInt:3_pInt) &
|
|
||||||
crystallite_dPdF(1:3,1:3,p,o,c,i,e) = crystallite_dPdF(1:3,1:3,p,o,c,i,e) + &
|
|
||||||
math_mul33x33(temp_33,math_transpose33(dFpinvdF(1:3,1:3,p,o)))
|
|
||||||
|
|
||||||
enddo; enddo
|
|
||||||
enddo elementLooping6
|
|
||||||
!$OMP END PARALLEL DO
|
|
||||||
|
|
||||||
else jacobianMethod
|
|
||||||
|
|
||||||
! --- STANDARD (PERTURBATION METHOD) FOR JACOBIAN ---
|
|
||||||
|
|
||||||
numerics_integrationMode = 2_pInt
|
|
||||||
|
|
||||||
! --- BACKUP ---
|
|
||||||
allocate(dPdF_perturbation1(3,3,3,3,homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(dPdF_perturbation2(3,3,3,3,homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(F_backup (3,3, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(Fp_backup (3,3, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(InvFp_backup (3,3, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(Fi_backup (3,3, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(InvFi_backup (3,3, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(Fe_backup (3,3, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(Lp_backup (3,3, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(Li_backup (3,3, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(P_backup (3,3, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(Tstar_v_backup (6, homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = 0.0_pReal)
|
|
||||||
allocate(convergenceFlag_backup (homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source = .false.)
|
|
||||||
|
|
||||||
!$OMP PARALLEL DO PRIVATE(myNcomponents)
|
|
||||||
elementLooping7: do e = FEsolving_execElem(1),FEsolving_execElem(2)
|
|
||||||
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
|
|
||||||
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do c = 1,myNcomponents
|
|
||||||
|
|
||||||
plasticState (phaseAt(c,i,e))%state_backup(:,phasememberAt(c,i,e)) = &
|
|
||||||
plasticState (phaseAt(c,i,e))%state( :,phasememberAt(c,i,e))
|
|
||||||
do mySource = 1_pInt, phase_Nsources(phaseAt(c,i,e))
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%state_backup(:,phasememberAt(c,i,e)) = &
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%state( :,phasememberAt(c,i,e))
|
|
||||||
enddo
|
|
||||||
|
|
||||||
plasticState (phaseAt(c,i,e))%dotState_backup(:,phasememberAt(c,i,e)) = &
|
|
||||||
plasticState (phaseAt(c,i,e))%dotState( :,phasememberAt(c,i,e))
|
|
||||||
do mySource = 1_pInt, phase_Nsources(phaseAt(c,i,e))
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%dotState_backup(:,phasememberAt(c,i,e)) = &
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%dotState( :,phasememberAt(c,i,e))
|
|
||||||
enddo
|
|
||||||
|
|
||||||
F_backup(1:3,1:3,c,i,e) = crystallite_subF(1:3,1:3,c,i,e) ! ... and kinematics
|
|
||||||
Fp_backup(1:3,1:3,c,i,e) = crystallite_Fp(1:3,1:3,c,i,e)
|
|
||||||
InvFp_backup(1:3,1:3,c,i,e) = crystallite_invFp(1:3,1:3,c,i,e)
|
|
||||||
Fi_backup(1:3,1:3,c,i,e) = crystallite_Fi(1:3,1:3,c,i,e)
|
|
||||||
InvFi_backup(1:3,1:3,c,i,e) = crystallite_invFi(1:3,1:3,c,i,e)
|
|
||||||
Fe_backup(1:3,1:3,c,i,e) = crystallite_Fe(1:3,1:3,c,i,e)
|
|
||||||
Lp_backup(1:3,1:3,c,i,e) = crystallite_Lp(1:3,1:3,c,i,e)
|
|
||||||
Li_backup(1:3,1:3,c,i,e) = crystallite_Li(1:3,1:3,c,i,e)
|
|
||||||
Tstar_v_backup(1:6,c,i,e) = crystallite_Tstar_v(1:6,c,i,e)
|
|
||||||
P_backup(1:3,1:3,c,i,e) = crystallite_P(1:3,1:3,c,i,e)
|
|
||||||
convergenceFlag_backup(c,i,e) = crystallite_converged(c,i,e)
|
|
||||||
enddo; enddo
|
|
||||||
enddo elementLooping7
|
|
||||||
!$END PARALLEL DO
|
|
||||||
! --- CALCULATE STATE AND STRESS FOR PERTURBATION ---
|
|
||||||
|
|
||||||
dPdF_perturbation1 = crystallite_dPdF0 ! initialize stiffness with known good values from last increment
|
|
||||||
dPdF_perturbation2 = crystallite_dPdF0 ! initialize stiffness with known good values from last increment
|
|
||||||
pertubationLoop: do perturbation = 1,2 ! forward and backward perturbation
|
|
||||||
if (iand(pert_method,perturbation) > 0_pInt) then ! mask for desired direction
|
|
||||||
myPert = -pert_Fg * (-1.0_pReal)**perturbation ! set perturbation step
|
|
||||||
do k = 1,3; do l = 1,3 ! ...alter individual components
|
|
||||||
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)) &
|
|
||||||
write(6,'(a,2(1x,i1),1x,a,/)') '<< CRYST >> [[[[[[ Stiffness perturbation',k,l,']]]]]]'
|
|
||||||
! --- INITIALIZE UNPERTURBED STATE ---
|
|
||||||
|
|
||||||
select case(numerics_integrator(numerics_integrationMode))
|
|
||||||
case(1_pInt)
|
|
||||||
!why not OMP? ! Fix-point method: restore to last converged state at end of subinc, since this is probably closest to perturbed state
|
|
||||||
do e = FEsolving_execElem(1),FEsolving_execElem(2)
|
|
||||||
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
|
|
||||||
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do c = 1,myNcomponents
|
|
||||||
|
|
||||||
plasticState (phaseAt(c,i,e))%state( :,phasememberAt(c,i,e)) = &
|
|
||||||
plasticState (phaseAt(c,i,e))%state_backup(:,phasememberAt(c,i,e))
|
|
||||||
do mySource = 1_pInt, phase_Nsources(phaseAt(c,i,e))
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%state( :,phasememberAt(c,i,e)) = &
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%state_backup(:,phasememberAt(c,i,e))
|
|
||||||
enddo
|
|
||||||
|
|
||||||
plasticState (phaseAt(c,i,e))%dotState( :,phasememberAt(c,i,e)) = &
|
|
||||||
plasticState (phaseAt(c,i,e))%dotState_backup(:,phasememberAt(c,i,e))
|
|
||||||
do mySource = 1_pInt, phase_Nsources(phaseAt(c,i,e))
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%dotState( :,phasememberAt(c,i,e)) = &
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%dotState_backup(:,phasememberAt(c,i,e))
|
|
||||||
enddo
|
|
||||||
|
|
||||||
crystallite_Fp(1:3,1:3,c,i,e) = Fp_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_invFp(1:3,1:3,c,i,e) = InvFp_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Fi(1:3,1:3,c,i,e) = Fi_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_invFi(1:3,1:3,c,i,e) = InvFi_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Fe(1:3,1:3,c,i,e) = Fe_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Lp(1:3,1:3,c,i,e) = Lp_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Li(1:3,1:3,c,i,e) = Li_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Tstar_v(1:6,c,i,e) = Tstar_v_backup(1:6,c,i,e)
|
|
||||||
enddo; enddo
|
|
||||||
enddo
|
|
||||||
case(2_pInt,3_pInt) ! explicit Euler methods: nothing to restore (except for F), since we are only doing a stress integration step
|
|
||||||
case(4_pInt,5_pInt)
|
|
||||||
!why not OMP? ! explicit Runge-Kutta methods: restore to start of subinc, since we are doing a full integration of state and stress
|
|
||||||
do e = FEsolving_execElem(1),FEsolving_execElem(2)
|
|
||||||
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
|
|
||||||
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do c = 1,myNcomponents
|
|
||||||
|
|
||||||
plasticState (phaseAt(c,i,e))%state( :,phasememberAt(c,i,e)) = &
|
|
||||||
plasticState (phaseAt(c,i,e))%subState0(:,phasememberAt(c,i,e))
|
|
||||||
do mySource = 1_pInt, phase_Nsources(phaseAt(c,i,e))
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%state( :,phasememberAt(c,i,e)) = &
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%subState0(:,phasememberAt(c,i,e))
|
|
||||||
enddo
|
|
||||||
|
|
||||||
plasticState (phaseAt(c,i,e))%dotState( :,phasememberAt(c,i,e)) = &
|
|
||||||
plasticState (phaseAt(c,i,e))%dotState_backup(:,phasememberAt(c,i,e))
|
|
||||||
do mySource = 1_pInt, phase_Nsources(phaseAt(c,i,e))
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%dotState( :,phasememberAt(c,i,e)) = &
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%dotState_backup(:,phasememberAt(c,i,e))
|
|
||||||
enddo
|
|
||||||
|
|
||||||
crystallite_Fp(1:3,1:3,c,i,e) = crystallite_subFp0(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_Fe(1:3,1:3,c,i,e) = crystallite_subFe0(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)
|
|
||||||
crystallite_Tstar_v(1:6,c,i,e) = crystallite_subTstar0_v(1:6,c,i,e)
|
|
||||||
enddo; enddo
|
|
||||||
enddo
|
|
||||||
end select
|
|
||||||
|
|
||||||
! --- PERTURB EITHER FORWARD OR BACKWARD ---
|
|
||||||
!why not OMP?
|
|
||||||
do e = FEsolving_execElem(1),FEsolving_execElem(2)
|
|
||||||
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
|
|
||||||
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
|
|
||||||
do c = 1,myNcomponents
|
|
||||||
crystallite_subF(1:3,1:3,c,i,e) = F_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_subF(k,l,c,i,e) = crystallite_subF(k,l,c,i,e) + myPert
|
|
||||||
crystallite_todo(c,i,e) = crystallite_requested(c,i,e) &
|
|
||||||
.and. convergenceFlag_backup(c,i,e)
|
|
||||||
if (crystallite_todo(c,i,e)) crystallite_converged(c,i,e) = .false. ! start out non-converged
|
|
||||||
enddo; enddo; enddo
|
|
||||||
|
|
||||||
|
|
||||||
select case(numerics_integrator(numerics_integrationMode))
|
|
||||||
case(1_pInt)
|
|
||||||
call crystallite_integrateStateFPI()
|
|
||||||
case(2_pInt)
|
|
||||||
call crystallite_integrateStateEuler()
|
|
||||||
case(3_pInt)
|
|
||||||
call crystallite_integrateStateAdaptiveEuler()
|
|
||||||
case(4_pInt)
|
|
||||||
call crystallite_integrateStateRK4()
|
|
||||||
case(5_pInt)
|
|
||||||
call crystallite_integrateStateRKCK45()
|
|
||||||
end select
|
|
||||||
!why not OMP?
|
|
||||||
elementLooping8: do e = FEsolving_execElem(1),FEsolving_execElem(2)
|
|
||||||
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
|
|
||||||
select case(perturbation)
|
|
||||||
case(1_pInt)
|
|
||||||
forall (i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), c = 1:myNcomponents, &
|
|
||||||
crystallite_requested(c,i,e) .and. crystallite_converged(c,i,e)) & ! converged state warrants stiffness update
|
|
||||||
dPdF_perturbation1(1:3,1:3,k,l,c,i,e) = &
|
|
||||||
(crystallite_P(1:3,1:3,c,i,e) - P_backup(1:3,1:3,c,i,e)) / myPert ! tangent dP_ij/dFg_kl
|
|
||||||
case(2_pInt)
|
|
||||||
forall (i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), c = 1:myNcomponents, &
|
|
||||||
crystallite_requested(c,i,e) .and. crystallite_converged(c,i,e)) & ! converged state warrants stiffness update
|
|
||||||
dPdF_perturbation2(1:3,1:3,k,l,c,i,e) = &
|
|
||||||
(crystallite_P(1:3,1:3,c,i,e) - P_backup(1:3,1:3,c,i,e)) / myPert ! tangent dP_ij/dFg_kl
|
|
||||||
end select
|
|
||||||
enddo elementLooping8
|
|
||||||
|
|
||||||
enddo; enddo ! k,l component perturbation loop
|
|
||||||
|
|
||||||
endif
|
|
||||||
enddo pertubationLoop
|
|
||||||
|
|
||||||
! --- STIFFNESS ACCORDING TO PERTURBATION METHOD AND CONVERGENCE ---
|
|
||||||
|
|
||||||
elementLooping9: do e = FEsolving_execElem(1),FEsolving_execElem(2)
|
|
||||||
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
|
|
||||||
select case(pert_method)
|
|
||||||
case(1_pInt)
|
|
||||||
forall (i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), c = 1:myNcomponents, &
|
|
||||||
crystallite_requested(c,i,e) .and. convergenceFlag_backup(c,i,e)) & ! perturbation mode 1: central solution converged
|
|
||||||
crystallite_dPdF(1:3,1:3,1:3,1:3,c,i,e) = dPdF_perturbation1(1:3,1:3,1:3,1:3,c,i,e)
|
|
||||||
case(2_pInt)
|
|
||||||
forall (i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), c = 1:myNcomponents, &
|
|
||||||
crystallite_requested(c,i,e) .and. convergenceFlag_backup(c,i,e)) & ! perturbation mode 2: central solution converged
|
|
||||||
crystallite_dPdF(1:3,1:3,1:3,1:3,c,i,e) = dPdF_perturbation2(1:3,1:3,1:3,1:3,c,i,e)
|
|
||||||
case(3_pInt)
|
|
||||||
forall (i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), c = 1:myNcomponents, &
|
|
||||||
crystallite_requested(c,i,e) .and. convergenceFlag_backup(c,i,e)) & ! perturbation mode 3: central solution converged
|
|
||||||
crystallite_dPdF(1:3,1:3,1:3,1:3,c,i,e) = 0.5_pReal* ( dPdF_perturbation1(1:3,1:3,1:3,1:3,c,i,e) &
|
|
||||||
+ dPdF_perturbation2(1:3,1:3,1:3,1:3,c,i,e))
|
|
||||||
end select
|
|
||||||
forall (i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), c = 1:myNcomponents, &
|
|
||||||
crystallite_requested(c,i,e) .and. .not. convergenceFlag_backup(c,i,e)) & ! for any pertubation mode: if central solution did not converge...
|
|
||||||
crystallite_dPdF(1:3,1:3,1:3,1:3,c,i,e) = crystallite_fallbackdPdF(1:3,1:3,1:3,1:3,c,i,e) ! ...use (elastic) fallback
|
|
||||||
enddo elementLooping9
|
|
||||||
|
|
||||||
! --- RESTORE ---
|
|
||||||
!why not OMP?
|
|
||||||
elementLooping10: do e = FEsolving_execElem(1),FEsolving_execElem(2)
|
|
||||||
myNcomponents = homogenization_Ngrains(mesh_element(3,e))
|
|
||||||
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do c = 1,myNcomponents
|
|
||||||
|
|
||||||
plasticState (phaseAt(c,i,e))%state( :,phasememberAt(c,i,e)) = &
|
|
||||||
plasticState (phaseAt(c,i,e))%state_backup(:,phasememberAt(c,i,e))
|
|
||||||
do mySource = 1_pInt, phase_Nsources(phaseAt(c,i,e))
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%state( :,phasememberAt(c,i,e)) = &
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%state_backup(:,phasememberAt(c,i,e))
|
|
||||||
enddo
|
|
||||||
|
|
||||||
plasticState (phaseAt(c,i,e))%dotState( :,phasememberAt(c,i,e)) = &
|
|
||||||
plasticState (phaseAt(c,i,e))%dotState_backup(:,phasememberAt(c,i,e))
|
|
||||||
do mySource = 1_pInt, phase_Nsources(phaseAt(c,i,e))
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%dotState( :,phasememberAt(c,i,e)) = &
|
|
||||||
sourceState(phaseAt(c,i,e))%p(mySource)%dotState_backup(:,phasememberAt(c,i,e))
|
|
||||||
enddo
|
|
||||||
|
|
||||||
crystallite_subF(1:3,1:3,c,i,e) = F_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Fp(1:3,1:3,c,i,e) = Fp_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_invFp(1:3,1:3,c,i,e) = InvFp_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Fi(1:3,1:3,c,i,e) = Fi_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_invFi(1:3,1:3,c,i,e) = InvFi_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Fe(1:3,1:3,c,i,e) = Fe_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Lp(1:3,1:3,c,i,e) = Lp_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Li(1:3,1:3,c,i,e) = Li_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_Tstar_v(1:6,c,i,e) = Tstar_v_backup(1:6,c,i,e)
|
|
||||||
crystallite_P(1:3,1:3,c,i,e) = P_backup(1:3,1:3,c,i,e)
|
|
||||||
crystallite_converged(c,i,e) = convergenceFlag_backup(c,i,e)
|
|
||||||
enddo; enddo
|
|
||||||
enddo elementLooping10
|
|
||||||
|
|
||||||
deallocate(dPdF_perturbation1)
|
|
||||||
deallocate(dPdF_perturbation2)
|
|
||||||
deallocate(F_backup )
|
|
||||||
deallocate(Fp_backup )
|
|
||||||
deallocate(InvFp_backup )
|
|
||||||
deallocate(Fi_backup )
|
|
||||||
deallocate(InvFi_backup )
|
|
||||||
deallocate(Fe_backup )
|
|
||||||
deallocate(Lp_backup )
|
|
||||||
deallocate(Li_backup )
|
|
||||||
deallocate(P_backup )
|
|
||||||
deallocate(Tstar_v_backup )
|
|
||||||
deallocate(convergenceFlag_backup)
|
|
||||||
|
|
||||||
endif jacobianMethod
|
|
||||||
endif computeJacobian
|
endif computeJacobian
|
||||||
!why not OMP?
|
!why not OMP?
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -64,7 +64,6 @@ module numerics
|
||||||
charLength = 1.0_pReal, & !< characteristic length scale for gradient problems
|
charLength = 1.0_pReal, & !< characteristic length scale for gradient problems
|
||||||
residualStiffness = 1.0e-6_pReal !< non-zero residual damage
|
residualStiffness = 1.0e-6_pReal !< non-zero residual damage
|
||||||
logical, protected, public :: &
|
logical, protected, public :: &
|
||||||
analyticJaco = .true., & !< use analytic Jacobian or perturbation, Default for Spectral solver .true.:
|
|
||||||
usePingPong = .true., &
|
usePingPong = .true., &
|
||||||
numerics_timeSyncing = .false. !< flag indicating if time synchronization in crystallite is used for nonlocal plasticity
|
numerics_timeSyncing = .false. !< flag indicating if time synchronization in crystallite is used for nonlocal plasticity
|
||||||
|
|
||||||
|
|
@ -315,8 +314,6 @@ subroutine numerics_init
|
||||||
numerics_integrator(1) = IO_intValue(line,chunkPos,2_pInt)
|
numerics_integrator(1) = IO_intValue(line,chunkPos,2_pInt)
|
||||||
case ('integratorstiffness')
|
case ('integratorstiffness')
|
||||||
numerics_integrator(2) = IO_intValue(line,chunkPos,2_pInt)
|
numerics_integrator(2) = IO_intValue(line,chunkPos,2_pInt)
|
||||||
case ('analyticjaco')
|
|
||||||
analyticJaco = IO_intValue(line,chunkPos,2_pInt) > 0_pInt
|
|
||||||
case ('usepingpong')
|
case ('usepingpong')
|
||||||
usepingpong = IO_intValue(line,chunkPos,2_pInt) > 0_pInt
|
usepingpong = IO_intValue(line,chunkPos,2_pInt) > 0_pInt
|
||||||
case ('timesyncing')
|
case ('timesyncing')
|
||||||
|
|
@ -528,7 +525,6 @@ subroutine numerics_init
|
||||||
write(6,'(a24,1x,es8.1)') ' aTol_crystalliteStress: ',aTol_crystalliteStress
|
write(6,'(a24,1x,es8.1)') ' aTol_crystalliteStress: ',aTol_crystalliteStress
|
||||||
write(6,'(a24,2(1x,i8))') ' integrator: ',numerics_integrator
|
write(6,'(a24,2(1x,i8))') ' integrator: ',numerics_integrator
|
||||||
write(6,'(a24,1x,L8)') ' timeSyncing: ',numerics_timeSyncing
|
write(6,'(a24,1x,L8)') ' timeSyncing: ',numerics_timeSyncing
|
||||||
write(6,'(a24,1x,L8)') ' analytic Jacobian: ',analyticJaco
|
|
||||||
write(6,'(a24,1x,L8)') ' use ping pong scheme: ',usepingpong
|
write(6,'(a24,1x,L8)') ' use ping pong scheme: ',usepingpong
|
||||||
write(6,'(a24,1x,es8.1,/)')' unitlength: ',numerics_unitlength
|
write(6,'(a24,1x,es8.1,/)')' unitlength: ',numerics_unitlength
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -152,8 +152,6 @@ subroutine plastic_disloUCLA_init(fileUnit)
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use lattice
|
use lattice
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -173,11 +171,9 @@ subroutine plastic_disloUCLA_init(fileUnit)
|
||||||
line = ''
|
line = ''
|
||||||
real(pReal), dimension(:), allocatable :: tempPerSlip
|
real(pReal), dimension(:), allocatable :: tempPerSlip
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_DISLOUCLA_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_DISLOUCLA_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_plasticity == PLASTICITY_DISLOUCLA_ID),pInt)
|
maxNinstance = int(count(phase_plasticity == PLASTICITY_DISLOUCLA_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -498,10 +494,6 @@ subroutine plastic_disloUCLA_init(fileUnit)
|
||||||
|
|
||||||
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(plasticState(phase)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(plasticState(phase)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -238,8 +238,6 @@ subroutine plastic_dislotwin_init(fileUnit)
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use lattice
|
use lattice
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -261,11 +259,9 @@ subroutine plastic_dislotwin_init(fileUnit)
|
||||||
line = ''
|
line = ''
|
||||||
real(pReal), dimension(:), allocatable :: tempPerSlip, tempPerTwin, tempPerTrans
|
real(pReal), dimension(:), allocatable :: tempPerSlip, tempPerTwin, tempPerTrans
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_DISLOTWIN_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_DISLOTWIN_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_plasticity == PLASTICITY_DISLOTWIN_ID),pInt)
|
maxNinstance = int(count(phase_plasticity == PLASTICITY_DISLOTWIN_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -930,10 +926,6 @@ subroutine plastic_dislotwin_init(fileUnit)
|
||||||
|
|
||||||
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(plasticState(phase)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(plasticState(phase)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -94,8 +94,6 @@ subroutine plastic_isotropic_init(fileUnit)
|
||||||
debug_constitutive, &
|
debug_constitutive, &
|
||||||
debug_levelBasic
|
debug_levelBasic
|
||||||
use numerics, only: &
|
use numerics, only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
use math, only: &
|
use math, only: &
|
||||||
math_Mandel3333to66, &
|
math_Mandel3333to66, &
|
||||||
|
|
@ -145,11 +143,9 @@ subroutine plastic_isotropic_init(fileUnit)
|
||||||
outputtag = ''
|
outputtag = ''
|
||||||
integer(pInt) :: NipcMyPhase
|
integer(pInt) :: NipcMyPhase
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_ISOTROPIC_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_ISOTROPIC_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_plasticity == PLASTICITY_ISOTROPIC_ID),pInt)
|
maxNinstance = int(count(phase_plasticity == PLASTICITY_ISOTROPIC_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -316,10 +312,6 @@ subroutine plastic_isotropic_init(fileUnit)
|
||||||
allocate(plasticState(phase)%state ( sizeState,NipcMyPhase),source=0.0_pReal)
|
allocate(plasticState(phase)%state ( sizeState,NipcMyPhase),source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%dotState (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
allocate(plasticState(phase)%dotState (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%deltaState (sizeDeltaState,NipcMyPhase),source=0.0_pReal)
|
allocate(plasticState(phase)%deltaState (sizeDeltaState,NipcMyPhase),source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(plasticState(phase)%state_backup ( sizeState,NipcMyPhase),source=0.0_pReal)
|
|
||||||
allocate(plasticState(phase)%dotState_backup (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(plasticState(phase)%previousDotState (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%previousDotState2(sizeDotState,NipcMyPhase),source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState2(sizeDotState,NipcMyPhase),source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -34,7 +34,6 @@ subroutine plastic_none_init
|
||||||
use IO, only: &
|
use IO, only: &
|
||||||
IO_timeStamp
|
IO_timeStamp
|
||||||
use numerics, only: &
|
use numerics, only: &
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
use material, only: &
|
use material, only: &
|
||||||
phase_plasticity, &
|
phase_plasticity, &
|
||||||
|
|
@ -53,11 +52,9 @@ subroutine plastic_none_init
|
||||||
sizeDotState, &
|
sizeDotState, &
|
||||||
sizeDeltaState
|
sizeDeltaState
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_NONE_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_NONE_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_plasticity == PLASTICITY_none_ID),pInt)
|
maxNinstance = int(count(phase_plasticity == PLASTICITY_none_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -84,11 +81,9 @@ subroutine plastic_none_init
|
||||||
allocate(plasticState(phase)%partionedState0 (sizeState,NofMyPhase))
|
allocate(plasticState(phase)%partionedState0 (sizeState,NofMyPhase))
|
||||||
allocate(plasticState(phase)%subState0 (sizeState,NofMyPhase))
|
allocate(plasticState(phase)%subState0 (sizeState,NofMyPhase))
|
||||||
allocate(plasticState(phase)%state (sizeState,NofMyPhase))
|
allocate(plasticState(phase)%state (sizeState,NofMyPhase))
|
||||||
allocate(plasticState(phase)%state_backup (sizeState,NofMyPhase))
|
|
||||||
|
|
||||||
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase))
|
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase))
|
||||||
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase))
|
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase))
|
||||||
allocate(plasticState(phase)%dotState_backup (sizeDotState,NofMyPhase))
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase))
|
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase))
|
||||||
allocate(plasticState(phase)%previousDotState2(sizeDotState,NofMyPhase))
|
allocate(plasticState(phase)%previousDotState2(sizeDotState,NofMyPhase))
|
||||||
|
|
|
||||||
|
|
@ -295,8 +295,6 @@ use material, only: phase_plasticity, &
|
||||||
material_phase
|
material_phase
|
||||||
use lattice
|
use lattice
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
|
|
||||||
|
|
@ -332,11 +330,9 @@ integer(pInt) :: phase, &
|
||||||
|
|
||||||
integer(pInt) :: NofMyPhase
|
integer(pInt) :: NofMyPhase
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_NONLOCAL_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_NONLOCAL_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstances = int(count(phase_plasticity == PLASTICITY_NONLOCAL_ID),pInt)
|
maxNinstances = int(count(phase_plasticity == PLASTICITY_NONLOCAL_ID),pInt)
|
||||||
if (maxNinstances == 0) return ! we don't have to do anything if there's no instance for this constitutive law
|
if (maxNinstances == 0) return ! we don't have to do anything if there's no instance for this constitutive law
|
||||||
|
|
@ -1310,10 +1306,6 @@ allocate(nonSchmidProjection(3,3,4,maxTotalNslip,maxNinstances),
|
||||||
|
|
||||||
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(plasticState(phase)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(plasticState(phase)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -145,8 +145,6 @@ subroutine plastic_phenoplus_init(fileUnit)
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use lattice
|
use lattice
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -168,11 +166,9 @@ subroutine plastic_phenoplus_init(fileUnit)
|
||||||
line = ''
|
line = ''
|
||||||
real(pReal), dimension(:), allocatable :: tempPerSlip
|
real(pReal), dimension(:), allocatable :: tempPerSlip
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_PHENOPLUS_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_PHENOPLUS_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_plasticity == PLASTICITY_PHENOPLUS_ID),pInt)
|
maxNinstance = int(count(phase_plasticity == PLASTICITY_PHENOPLUS_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -589,10 +585,6 @@ subroutine plastic_phenoplus_init(fileUnit)
|
||||||
allocate(plasticState(phase)%state ( sizeState,NipcMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%state ( sizeState,NipcMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%dotState (sizeDotState,NipcMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%dotState (sizeDotState,NipcMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%deltaState (sizeDeltaState,NipcMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%deltaState (sizeDeltaState,NipcMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(plasticState(phase)%state_backup ( sizeState,NipcMyPhase),source=0.0_pReal)
|
|
||||||
allocate(plasticState(phase)%dotState_backup (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(plasticState(phase)%previousDotState (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%previousDotState2(sizeDotState,NipcMyPhase),source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState2(sizeDotState,NipcMyPhase),source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -157,8 +157,6 @@ subroutine plastic_phenopowerlaw_init(fileUnit)
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use lattice
|
use lattice
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -181,11 +179,9 @@ subroutine plastic_phenopowerlaw_init(fileUnit)
|
||||||
line = ''
|
line = ''
|
||||||
real(pReal), dimension(:), allocatable :: tempPerSlip
|
real(pReal), dimension(:), allocatable :: tempPerSlip
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_PHENOPOWERLAW_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_PHENOPOWERLAW_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_plasticity == PLASTICITY_PHENOPOWERLAW_ID),pInt)
|
maxNinstance = int(count(phase_plasticity == PLASTICITY_PHENOPOWERLAW_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -587,10 +583,6 @@ subroutine plastic_phenopowerlaw_init(fileUnit)
|
||||||
allocate(plasticState(phase)%state ( sizeState,NipcMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%state ( sizeState,NipcMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%dotState (sizeDotState,NipcMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%dotState (sizeDotState,NipcMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%deltaState (sizeDeltaState,NipcMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%deltaState (sizeDeltaState,NipcMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(plasticState(phase)%state_backup ( sizeState,NipcMyPhase),source=0.0_pReal)
|
|
||||||
allocate(plasticState(phase)%dotState_backup (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(plasticState(phase)%previousDotState (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState (sizeDotState,NipcMyPhase),source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%previousDotState2(sizeDotState,NipcMyPhase),source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState2(sizeDotState,NipcMyPhase),source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -216,8 +216,6 @@ subroutine plastic_titanmod_init(fileUnit)
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use lattice
|
use lattice
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -241,11 +239,9 @@ subroutine plastic_titanmod_init(fileUnit)
|
||||||
tag = '', &
|
tag = '', &
|
||||||
line = ''
|
line = ''
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_TITANMOD_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_TITANMOD_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_plasticity == PLASTICITY_TITANMOD_ID),pInt)
|
maxNinstance = int(count(phase_plasticity == PLASTICITY_TITANMOD_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -859,10 +855,6 @@ subroutine plastic_titanmod_init(fileUnit)
|
||||||
|
|
||||||
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(plasticState(phase)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(plasticState(phase)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(plasticState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(plasticState(phase)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -67,11 +67,9 @@ module prec
|
||||||
real(pReal), allocatable, dimension(:,:) :: &
|
real(pReal), allocatable, dimension(:,:) :: &
|
||||||
partionedState0, &
|
partionedState0, &
|
||||||
subState0, &
|
subState0, &
|
||||||
state_backup, &
|
|
||||||
deltaState, &
|
deltaState, &
|
||||||
previousDotState, & !< state rate of previous xxxx
|
previousDotState, & !< state rate of previous xxxx
|
||||||
previousDotState2, & !< state rate two xxxx ago
|
previousDotState2, & !< state rate two xxxx ago
|
||||||
dotState_backup, & !< backup of state rate
|
|
||||||
RK4dotState
|
RK4dotState
|
||||||
real(pReal), allocatable, dimension(:,:,:) :: &
|
real(pReal), allocatable, dimension(:,:,:) :: &
|
||||||
RKCK45dotState
|
RKCK45dotState
|
||||||
|
|
|
||||||
|
|
@ -92,8 +92,6 @@ subroutine source_damage_anisoBrittle_init(fileUnit)
|
||||||
sourceState, &
|
sourceState, &
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
use lattice, only: &
|
use lattice, only: &
|
||||||
lattice_maxNcleavageFamily, &
|
lattice_maxNcleavageFamily, &
|
||||||
|
|
@ -111,11 +109,9 @@ subroutine source_damage_anisoBrittle_init(fileUnit)
|
||||||
tag = '', &
|
tag = '', &
|
||||||
line = ''
|
line = ''
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- source_'//SOURCE_damage_anisoBrittle_LABEL//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- source_'//SOURCE_damage_anisoBrittle_LABEL//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_source == SOURCE_damage_anisoBrittle_ID),pInt)
|
maxNinstance = int(count(phase_source == SOURCE_damage_anisoBrittle_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -268,10 +264,6 @@ subroutine source_damage_anisoBrittle_init(fileUnit)
|
||||||
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -96,8 +96,6 @@ subroutine source_damage_anisoDuctile_init(fileUnit)
|
||||||
sourceState, &
|
sourceState, &
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
use lattice, only: &
|
use lattice, only: &
|
||||||
lattice_maxNslipFamily, &
|
lattice_maxNslipFamily, &
|
||||||
|
|
@ -115,11 +113,9 @@ subroutine source_damage_anisoDuctile_init(fileUnit)
|
||||||
tag = '', &
|
tag = '', &
|
||||||
line = ''
|
line = ''
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- source_'//SOURCE_damage_anisoDuctile_LABEL//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- source_'//SOURCE_damage_anisoDuctile_LABEL//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_source == SOURCE_damage_anisoDuctile_ID),pInt)
|
maxNinstance = int(count(phase_source == SOURCE_damage_anisoDuctile_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -270,10 +266,6 @@ subroutine source_damage_anisoDuctile_init(fileUnit)
|
||||||
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -82,8 +82,6 @@ subroutine source_damage_isoBrittle_init(fileUnit)
|
||||||
sourceState, &
|
sourceState, &
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -97,11 +95,9 @@ subroutine source_damage_isoBrittle_init(fileUnit)
|
||||||
tag = '', &
|
tag = '', &
|
||||||
line = ''
|
line = ''
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- source_'//SOURCE_damage_isoBrittle_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- source_'//SOURCE_damage_isoBrittle_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_source == SOURCE_damage_isoBrittle_ID),pInt)
|
maxNinstance = int(count(phase_source == SOURCE_damage_isoBrittle_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -222,10 +218,6 @@ subroutine source_damage_isoBrittle_init(fileUnit)
|
||||||
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -82,8 +82,6 @@ subroutine source_damage_isoDuctile_init(fileUnit)
|
||||||
sourceState, &
|
sourceState, &
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -97,11 +95,9 @@ subroutine source_damage_isoDuctile_init(fileUnit)
|
||||||
tag = '', &
|
tag = '', &
|
||||||
line = ''
|
line = ''
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- source_'//SOURCE_damage_isoDuctile_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- source_'//SOURCE_damage_isoDuctile_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_source == SOURCE_damage_isoDuctile_ID),pInt)
|
maxNinstance = int(count(phase_source == SOURCE_damage_isoDuctile_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -222,10 +218,6 @@ subroutine source_damage_isoDuctile_init(fileUnit)
|
||||||
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -68,8 +68,6 @@ subroutine source_thermal_dissipation_init(fileUnit)
|
||||||
sourceState, &
|
sourceState, &
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -83,11 +81,9 @@ subroutine source_thermal_dissipation_init(fileUnit)
|
||||||
tag = '', &
|
tag = '', &
|
||||||
line = ''
|
line = ''
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- source_'//SOURCE_thermal_dissipation_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- source_'//SOURCE_thermal_dissipation_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_source == SOURCE_thermal_dissipation_ID),pInt)
|
maxNinstance = int(count(phase_source == SOURCE_thermal_dissipation_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -163,10 +159,6 @@ subroutine source_thermal_dissipation_init(fileUnit)
|
||||||
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -73,8 +73,6 @@ subroutine source_thermal_externalheat_init(fileUnit)
|
||||||
sourceState, &
|
sourceState, &
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -89,11 +87,9 @@ subroutine source_thermal_externalheat_init(fileUnit)
|
||||||
line = ''
|
line = ''
|
||||||
real(pReal), allocatable, dimension(:,:) :: temp_time, temp_rate
|
real(pReal), allocatable, dimension(:,:) :: temp_time, temp_rate
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- source_'//SOURCE_thermal_externalheat_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- source_'//SOURCE_thermal_externalheat_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_source == SOURCE_thermal_externalheat_ID),pInt)
|
maxNinstance = int(count(phase_source == SOURCE_thermal_externalheat_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -189,10 +185,6 @@ subroutine source_thermal_externalheat_init(fileUnit)
|
||||||
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -70,8 +70,6 @@ subroutine source_vacancy_irradiation_init(fileUnit)
|
||||||
sourceState, &
|
sourceState, &
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -85,11 +83,9 @@ subroutine source_vacancy_irradiation_init(fileUnit)
|
||||||
tag = '', &
|
tag = '', &
|
||||||
line = ''
|
line = ''
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- source_'//SOURCE_vacancy_irradiation_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- source_'//SOURCE_vacancy_irradiation_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_source == SOURCE_vacancy_irradiation_ID),pInt)
|
maxNinstance = int(count(phase_source == SOURCE_vacancy_irradiation_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -169,10 +165,6 @@ subroutine source_vacancy_irradiation_init(fileUnit)
|
||||||
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -68,8 +68,6 @@ subroutine source_vacancy_phenoplasticity_init(fileUnit)
|
||||||
sourceState, &
|
sourceState, &
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -83,11 +81,9 @@ subroutine source_vacancy_phenoplasticity_init(fileUnit)
|
||||||
tag = '', &
|
tag = '', &
|
||||||
line = ''
|
line = ''
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- source_'//SOURCE_vacancy_phenoplasticity_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- source_'//SOURCE_vacancy_phenoplasticity_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_source == SOURCE_vacancy_phenoplasticity_ID),pInt)
|
maxNinstance = int(count(phase_source == SOURCE_vacancy_phenoplasticity_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -163,10 +159,6 @@ subroutine source_vacancy_phenoplasticity_init(fileUnit)
|
||||||
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -72,8 +72,6 @@ subroutine source_vacancy_thermalfluc_init(fileUnit)
|
||||||
sourceState, &
|
sourceState, &
|
||||||
MATERIAL_partPhase
|
MATERIAL_partPhase
|
||||||
use numerics,only: &
|
use numerics,only: &
|
||||||
analyticJaco, &
|
|
||||||
worldrank, &
|
|
||||||
numerics_integrator
|
numerics_integrator
|
||||||
|
|
||||||
implicit none
|
implicit none
|
||||||
|
|
@ -87,11 +85,9 @@ subroutine source_vacancy_thermalfluc_init(fileUnit)
|
||||||
tag = '', &
|
tag = '', &
|
||||||
line = ''
|
line = ''
|
||||||
|
|
||||||
mainProcess: if (worldrank == 0) then
|
write(6,'(/,a)') ' <<<+- source_'//SOURCE_vacancy_thermalfluc_label//' init -+>>>'
|
||||||
write(6,'(/,a)') ' <<<+- source_'//SOURCE_vacancy_thermalfluc_label//' init -+>>>'
|
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
||||||
write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
|
|
||||||
#include "compilation_info.f90"
|
#include "compilation_info.f90"
|
||||||
endif mainProcess
|
|
||||||
|
|
||||||
maxNinstance = int(count(phase_source == SOURCE_vacancy_thermalfluc_ID),pInt)
|
maxNinstance = int(count(phase_source == SOURCE_vacancy_thermalfluc_ID),pInt)
|
||||||
if (maxNinstance == 0_pInt) return
|
if (maxNinstance == 0_pInt) return
|
||||||
|
|
@ -170,10 +166,6 @@ subroutine source_vacancy_thermalfluc_init(fileUnit)
|
||||||
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
|
||||||
if (.not. analyticJaco) then
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%state_backup (sizeState,NofMyPhase), source=0.0_pReal)
|
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%dotState_backup (sizeDotState,NofMyPhase), source=0.0_pReal)
|
|
||||||
endif
|
|
||||||
if (any(numerics_integrator == 1_pInt)) then
|
if (any(numerics_integrator == 1_pInt)) then
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
|
||||||
|
|
|
||||||
|
|
@ -1,2 +1 @@
|
||||||
fixed_seed 1697667030
|
fixed_seed 1697667030
|
||||||
analyticJaco 1
|
|
||||||
|
|
|
||||||
|
|
@ -8,7 +8,6 @@ pert_Fg 1.0e-7 # deformation gradient perturbation for g
|
||||||
pert_method 1 # perturbation method (1 = forward, 2 = backward or 3 = central)
|
pert_method 1 # perturbation method (1 = forward, 2 = backward or 3 = central)
|
||||||
integrator 1 # integration method (1 = Fixed Point Iteration, 2 = Euler, 3 = Adaptive Euler, 4 = classical 4th order Runge-Kutta, 5 = 5th order Runge-Kutta Cash-Karp)
|
integrator 1 # integration method (1 = Fixed Point Iteration, 2 = Euler, 3 = Adaptive Euler, 4 = classical 4th order Runge-Kutta, 5 = 5th order Runge-Kutta Cash-Karp)
|
||||||
integratorStiffness 1 # integration method used for stiffness (1 = Fixed Point Iteration, 2 = Euler, 3 = Adaptive Euler, 4 = classical 4th order Runge-Kutta, 5 = 5th order Runge-Kutta Cash-Karp)
|
integratorStiffness 1 # integration method used for stiffness (1 = Fixed Point Iteration, 2 = Euler, 3 = Adaptive Euler, 4 = classical 4th order Runge-Kutta, 5 = 5th order Runge-Kutta Cash-Karp)
|
||||||
analyticJaco 1 # use analytic Jacobian or perturbation (0 = perturbations, 1 = analytic)
|
|
||||||
unitlength 1 # physical length of one computational length unit
|
unitlength 1 # physical length of one computational length unit
|
||||||
usepingpong 1 # use the ping pong (collect <-> calc) scheme (always off for Abaqus exp, must be on for Spectral Solver)
|
usepingpong 1 # use the ping pong (collect <-> calc) scheme (always off for Abaqus exp, must be on for Spectral Solver)
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -1,7 +1,6 @@
|
||||||
# -*- coding: UTF-8 no BOM -*-
|
# -*- coding: UTF-8 no BOM -*-
|
||||||
|
|
||||||
|
import os,subprocess,shlex,re,string
|
||||||
import os,subprocess,shlex,re
|
|
||||||
|
|
||||||
class Environment():
|
class Environment():
|
||||||
__slots__ = [ \
|
__slots__ = [ \
|
||||||
|
|
@ -23,7 +22,7 @@ class Environment():
|
||||||
for line in configFile:
|
for line in configFile:
|
||||||
l = re.sub('^set ', '', line).strip() # remove "set" (tcsh) when setting variables
|
l = re.sub('^set ', '', line).strip() # remove "set" (tcsh) when setting variables
|
||||||
if l and not l.startswith('#'):
|
if l and not l.startswith('#'):
|
||||||
items = l.split('=')
|
items = map(string.strip,l.split('='))
|
||||||
if len(items) == 2:
|
if len(items) == 2:
|
||||||
self.options[items[0].upper()] = \
|
self.options[items[0].upper()] = \
|
||||||
re.sub('\$\{*DAMASK_ROOT\}*',self.rootDir(),os.path.expandvars(items[1])) # expand all shell variables and DAMASK_ROOT
|
re.sub('\$\{*DAMASK_ROOT\}*',self.rootDir(),os.path.expandvars(items[1])) # expand all shell variables and DAMASK_ROOT
|
||||||
|
|
|
||||||
|
|
@ -413,11 +413,15 @@ class Quaternion:
|
||||||
|
|
||||||
|
|
||||||
@classmethod
|
@classmethod
|
||||||
def fromAngleAxis(cls, angle, axis):
|
def fromAngleAxis(cls,
|
||||||
|
angle,
|
||||||
|
axis,
|
||||||
|
degrees = False):
|
||||||
if not isinstance(axis, np.ndarray): axis = np.array(axis,dtype='d')
|
if not isinstance(axis, np.ndarray): axis = np.array(axis,dtype='d')
|
||||||
axis = axis.astype(float)/np.linalg.norm(axis)
|
axis = axis.astype(float)/np.linalg.norm(axis)
|
||||||
s = math.sin(angle / 2.0)
|
angle = np.radians(angle) if degrees else angle
|
||||||
w = math.cos(angle / 2.0)
|
s = math.sin(0.5 * angle)
|
||||||
|
w = math.cos(0.5 * angle)
|
||||||
x = axis[0] * s
|
x = axis[0] * s
|
||||||
y = axis[1] * s
|
y = axis[1] * s
|
||||||
z = axis[2] * s
|
z = axis[2] * s
|
||||||
|
|
@ -425,27 +429,26 @@ class Quaternion:
|
||||||
|
|
||||||
|
|
||||||
@classmethod
|
@classmethod
|
||||||
def fromEulers(cls, eulers, type = 'Bunge'):
|
def fromEulers(cls,
|
||||||
|
eulers,
|
||||||
|
type = 'Bunge',
|
||||||
|
degrees = False):
|
||||||
|
if not isinstance(eulers, np.ndarray): eulers = np.array(eulers,dtype='d')
|
||||||
|
eulers = np.radians(eulers) if degrees else eulers
|
||||||
|
|
||||||
eulers *= 0.5 # reduce to half angles
|
c = np.cos(0.5 * eulers)
|
||||||
|
s = np.sin(0.5 * eulers)
|
||||||
c1 = math.cos(eulers[0])
|
|
||||||
s1 = math.sin(eulers[0])
|
|
||||||
c2 = math.cos(eulers[1])
|
|
||||||
s2 = math.sin(eulers[1])
|
|
||||||
c3 = math.cos(eulers[2])
|
|
||||||
s3 = math.sin(eulers[2])
|
|
||||||
|
|
||||||
if type.lower() == 'bunge' or type.lower() == 'zxz':
|
if type.lower() == 'bunge' or type.lower() == 'zxz':
|
||||||
w = c1 * c2 * c3 - s1 * c2 * s3
|
w = c[0] * c[1] * c[2] - s[0] * c[1] * s[2]
|
||||||
x = c1 * s2 * c3 + s1 * s2 * s3
|
x = c[0] * s[1] * c[2] + s[0] * s[1] * s[2]
|
||||||
y = - c1 * s2 * s3 + s1 * s2 * c3
|
y = - c[0] * s[1] * s[2] + s[0] * s[1] * c[2]
|
||||||
z = c1 * c2 * s3 + s1 * c2 * c3
|
z = c[0] * c[1] * s[2] + s[0] * c[1] * c[2]
|
||||||
else:
|
else:
|
||||||
w = c1 * c2 * c3 - s1 * s2 * s3
|
w = c[0] * c[1] * c[2] - s[0] * s[1] * s[2]
|
||||||
x = s1 * s2 * c3 + c1 * c2 * s3
|
x = s[0] * s[1] * c[2] + c[0] * c[1] * s[2]
|
||||||
y = s1 * c2 * c3 + c1 * s2 * s3
|
y = s[0] * c[1] * c[2] + c[0] * s[1] * s[2]
|
||||||
z = c1 * s2 * c3 - s1 * c2 * s3
|
z = c[0] * s[1] * c[2] - s[0] * c[1] * s[2]
|
||||||
return cls([w,x,y,z])
|
return cls([w,x,y,z])
|
||||||
|
|
||||||
|
|
||||||
|
|
@ -819,6 +822,7 @@ class Orientation:
|
||||||
Eulers = None,
|
Eulers = None,
|
||||||
random = False, # integer to have a fixed seed or True for real random
|
random = False, # integer to have a fixed seed or True for real random
|
||||||
symmetry = None,
|
symmetry = None,
|
||||||
|
degrees = False,
|
||||||
):
|
):
|
||||||
if random: # produce random orientation
|
if random: # produce random orientation
|
||||||
if isinstance(random, bool ):
|
if isinstance(random, bool ):
|
||||||
|
|
@ -826,16 +830,16 @@ class Orientation:
|
||||||
else:
|
else:
|
||||||
self.quaternion = Quaternion.fromRandom(randomSeed=random)
|
self.quaternion = Quaternion.fromRandom(randomSeed=random)
|
||||||
elif isinstance(Eulers, np.ndarray) and Eulers.shape == (3,): # based on given Euler angles
|
elif isinstance(Eulers, np.ndarray) and Eulers.shape == (3,): # based on given Euler angles
|
||||||
self.quaternion = Quaternion.fromEulers(Eulers,'bunge')
|
self.quaternion = Quaternion.fromEulers(Eulers,type='bunge',degrees=degrees)
|
||||||
elif isinstance(matrix, np.ndarray) : # based on given rotation matrix
|
elif isinstance(matrix, np.ndarray) : # based on given rotation matrix
|
||||||
self.quaternion = Quaternion.fromMatrix(matrix)
|
self.quaternion = Quaternion.fromMatrix(matrix)
|
||||||
elif isinstance(angleAxis, np.ndarray) and angleAxis.shape == (4,): # based on given angle and rotation axis
|
elif isinstance(angleAxis, np.ndarray) and angleAxis.shape == (4,): # based on given angle and rotation axis
|
||||||
self.quaternion = Quaternion.fromAngleAxis(angleAxis[0],angleAxis[1:4])
|
self.quaternion = Quaternion.fromAngleAxis(angleAxis[0],angleAxis[1:4],degrees=degrees)
|
||||||
elif isinstance(Rodrigues, np.ndarray) and Rodrigues.shape == (3,): # based on given Rodrigues vector
|
elif isinstance(Rodrigues, np.ndarray) and Rodrigues.shape == (3,): # based on given Rodrigues vector
|
||||||
self.quaternion = Quaternion.fromRodrigues(Rodrigues)
|
self.quaternion = Quaternion.fromRodrigues(Rodrigues)
|
||||||
elif isinstance(quaternion, Quaternion): # based on given quaternion
|
elif isinstance(quaternion, Quaternion): # based on given quaternion
|
||||||
self.quaternion = quaternion.homomorphed()
|
self.quaternion = quaternion.homomorphed()
|
||||||
elif isinstance(quaternion, np.ndarray) and quaternion.shape == (4,): # based on given quaternion
|
elif isinstance(quaternion, np.ndarray) and quaternion.shape == (4,): # based on given quaternion-like array
|
||||||
self.quaternion = Quaternion(quaternion).homomorphed()
|
self.quaternion = Quaternion(quaternion).homomorphed()
|
||||||
|
|
||||||
self.symmetry = Symmetry(symmetry)
|
self.symmetry = Symmetry(symmetry)
|
||||||
|
|
|
||||||
|
|
@ -1,7 +1,7 @@
|
||||||
#!/usr/bin/env python2.7
|
#!/usr/bin/env python2.7
|
||||||
# -*- coding: UTF-8 no BOM -*-
|
# -*- coding: UTF-8 no BOM -*-
|
||||||
|
|
||||||
import os,re,sys
|
import os,re,sys,collections
|
||||||
import math # noqa
|
import math # noqa
|
||||||
import numpy as np
|
import numpy as np
|
||||||
from optparse import OptionParser
|
from optparse import OptionParser
|
||||||
|
|
@ -10,6 +10,10 @@ import damask
|
||||||
scriptName = os.path.splitext(os.path.basename(__file__))[0]
|
scriptName = os.path.splitext(os.path.basename(__file__))[0]
|
||||||
scriptID = ' '.join([scriptName,damask.version])
|
scriptID = ' '.join([scriptName,damask.version])
|
||||||
|
|
||||||
|
def listify(x):
|
||||||
|
return x if isinstance(x, collections.Iterable) else [x]
|
||||||
|
|
||||||
|
|
||||||
# --------------------------------------------------------------------
|
# --------------------------------------------------------------------
|
||||||
# MAIN
|
# MAIN
|
||||||
# --------------------------------------------------------------------
|
# --------------------------------------------------------------------
|
||||||
|
|
@ -17,10 +21,12 @@ scriptID = ' '.join([scriptName,damask.version])
|
||||||
parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [file[s]]', description = """
|
parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [file[s]]', description = """
|
||||||
Add or alter column(s) with derived values according to user-defined arithmetic operation between column(s).
|
Add or alter column(s) with derived values according to user-defined arithmetic operation between column(s).
|
||||||
Column labels are tagged by '#label#' in formulas. Use ';' for ',' in functions.
|
Column labels are tagged by '#label#' in formulas. Use ';' for ',' in functions.
|
||||||
Numpy is available as np.
|
Numpy is available as 'np'.
|
||||||
|
|
||||||
Special variables: #_row_# -- row index
|
Special variables: #_row_# -- row index
|
||||||
Examples: (1) magnitude of vector -- "np.linalg.norm(#vec#)" (2) rounded root of row number -- "round(math.sqrt(#_row_#);3)"
|
Examples:
|
||||||
|
(1) magnitude of vector -- "np.linalg.norm(#vec#)"
|
||||||
|
(2) rounded root of row number -- "round(math.sqrt(#_row_#);3)"
|
||||||
|
|
||||||
""", version = scriptID)
|
""", version = scriptID)
|
||||||
|
|
||||||
|
|
@ -50,30 +56,26 @@ if len(options.labels) != len(options.formulas):
|
||||||
for i in xrange(len(options.formulas)):
|
for i in xrange(len(options.formulas)):
|
||||||
options.formulas[i] = options.formulas[i].replace(';',',')
|
options.formulas[i] = options.formulas[i].replace(';',',')
|
||||||
|
|
||||||
# --- loop over input files -------------------------------------------------------------------------
|
# ------------------------------------- loop over input files --------------------------------------
|
||||||
|
|
||||||
if filenames == []: filenames = [None]
|
if filenames == []: filenames = [None]
|
||||||
|
|
||||||
for name in filenames:
|
for name in filenames:
|
||||||
try:
|
try: table = damask.ASCIItable(name = name,
|
||||||
table = damask.ASCIItable(name = name,
|
buffered = False)
|
||||||
buffered = False)
|
except: continue
|
||||||
output = damask.ASCIItable(name = name,
|
|
||||||
buffered = False)
|
|
||||||
except:
|
|
||||||
continue
|
|
||||||
damask.util.report(scriptName,name)
|
damask.util.report(scriptName,name)
|
||||||
|
|
||||||
# ------------------------------------------ read header -------------------------------------------
|
# ------------------------------------------ read header -------------------------------------------
|
||||||
|
|
||||||
table.head_read()
|
table.head_read()
|
||||||
|
|
||||||
# -----------------------------------------------------------------------------------------------------
|
# --------------------------------------------------------------------------------------------------
|
||||||
specials = { \
|
specials = { \
|
||||||
'_row_': 0,
|
'_row_': 0,
|
||||||
}
|
}
|
||||||
|
|
||||||
# ------------------------------------------ Evaluate condition ---------------------------------------
|
# --------------------------------------- evaluate condition ---------------------------------------
|
||||||
if options.condition is not None:
|
if options.condition is not None:
|
||||||
interpolator = []
|
interpolator = []
|
||||||
condition = options.condition # copy per file, since might be altered inline
|
condition = options.condition # copy per file, since might be altered inline
|
||||||
|
|
@ -98,7 +100,7 @@ for name in filenames:
|
||||||
|
|
||||||
evaluator_condition = "'" + condition + "'.format(" + ','.join(interpolator) + ")"
|
evaluator_condition = "'" + condition + "'.format(" + ','.join(interpolator) + ")"
|
||||||
|
|
||||||
# ------------------------------------------ build formulae ----------------------------------------
|
# ------------------------------------------ build formulas ----------------------------------------
|
||||||
|
|
||||||
evaluator = {}
|
evaluator = {}
|
||||||
|
|
||||||
|
|
@ -121,6 +123,10 @@ for name in filenames:
|
||||||
|
|
||||||
evaluator[label] = formula
|
evaluator[label] = formula
|
||||||
|
|
||||||
|
# ---------------------------- separate requested labels into old and new --------------------------
|
||||||
|
|
||||||
|
veterans = list(set(options.labels)&set(table.labels(raw=False)+table.labels(raw=True)) ) # intersection of requested and existing
|
||||||
|
newbies = list(set(options.labels)-set(table.labels(raw=False)+table.labels(raw=True)) ) # requested but not existing
|
||||||
|
|
||||||
# ------------------------------------------ process data ------------------------------------------
|
# ------------------------------------------ process data ------------------------------------------
|
||||||
|
|
||||||
|
|
@ -129,53 +135,45 @@ for name in filenames:
|
||||||
|
|
||||||
while outputAlive and table.data_read(): # read next data line of ASCII table
|
while outputAlive and table.data_read(): # read next data line of ASCII table
|
||||||
specials['_row_'] += 1 # count row
|
specials['_row_'] += 1 # count row
|
||||||
output.data_clear()
|
|
||||||
|
|
||||||
# ------------------------------------------ calculate one result to get length of labels ---------
|
|
||||||
|
|
||||||
if firstLine:
|
if firstLine:
|
||||||
firstLine = False
|
firstLine = False
|
||||||
labelDim = {}
|
|
||||||
for label in [x for x in options.labels]:
|
|
||||||
labelDim[label] = np.size(eval(evaluator[label]))
|
|
||||||
if labelDim[label] == 0: options.labels.remove(label)
|
|
||||||
|
|
||||||
# ------------------------------------------ assemble header ---------------------------------------
|
# ---------------------------- line 1: determine dimension of formulas -----------------------------
|
||||||
|
|
||||||
output.labels_clear()
|
resultDim = {}
|
||||||
tabLabels = table.labels()
|
for label in list(options.labels): # iterate over stable copy
|
||||||
for label in tabLabels:
|
resultDim[label] = np.size(eval(evaluator[label])) # get dimension of formula[label]
|
||||||
dim = labelDim[label] if label in options.labels \
|
if resultDim[label] == 0: options.labels.remove(label) # remove label if invalid result
|
||||||
else table.label_dimension(label)
|
|
||||||
output.labels_append(['{}_{}'.format(i+1,label) for i in xrange(dim)] if dim > 1 else label)
|
for veteran in list(veterans):
|
||||||
|
if resultDim[veteran] != table.label_dimension(veteran):
|
||||||
|
damask.util.croak('skipping {} due to inconsistent dimension...'.format(veteran))
|
||||||
|
veterans.remove(veteran) # discard culprit
|
||||||
|
|
||||||
for label in options.labels:
|
# ----------------------------------- line 1: assemble header --------------------------------------
|
||||||
if label in tabLabels: continue
|
|
||||||
output.labels_append(['{}_{}'.format(i+1,label) for i in xrange(labelDim[label])]
|
|
||||||
if labelDim[label] > 1
|
|
||||||
else label)
|
|
||||||
|
|
||||||
output.info = table.info
|
for newby in newbies:
|
||||||
output.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
|
table.labels_append(['{}_{}'.format(i+1,newby) for i in xrange(resultDim[newby])]
|
||||||
output.head_write()
|
if resultDim[newby] > 1 else newby)
|
||||||
|
|
||||||
# ------------------------------------------ process data ------------------------------------------
|
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
|
||||||
|
table.head_write()
|
||||||
|
|
||||||
for label in output.labels():
|
# -------------------------------------- evaluate formulas -----------------------------------------
|
||||||
oldIndices = table.label_indexrange(label)
|
|
||||||
Nold = max(1,len(oldIndices)) # Nold could be zero for new columns
|
|
||||||
Nnew = len(output.label_indexrange(label))
|
|
||||||
output.data_append(eval(evaluator[label]) if label in options.labels and
|
|
||||||
(options.condition is None or eval(eval(evaluator_condition)))
|
|
||||||
else np.tile([table.data[i] for i in oldIndices]
|
|
||||||
if label in tabLabels
|
|
||||||
else np.nan,
|
|
||||||
np.ceil(float(Nnew)/Nold))[:Nnew]) # spread formula result into given number of columns
|
|
||||||
|
|
||||||
outputAlive = output.data_write() # output processed line
|
if options.condition is None or eval(eval(evaluator_condition)): # condition for veteran replacement fulfilled
|
||||||
|
for veteran in veterans: # evaluate formulas that overwrite
|
||||||
|
table.data[table.label_index(veteran):
|
||||||
|
table.label_index(veteran)+table.label_dimension(veteran)] = \
|
||||||
|
listify(eval(evaluator[veteran]))
|
||||||
|
|
||||||
|
for newby in newbies: # evaluate formulas that append
|
||||||
|
table.data_append(listify(eval(evaluator[newby])))
|
||||||
|
|
||||||
# ------------------------------------------ output finalization -----------------------------------
|
outputAlive = table.data_write() # output processed line
|
||||||
|
|
||||||
|
# ------------------------------------- output finalization ----------------------------------------
|
||||||
|
|
||||||
|
table.close() # close ASCII table
|
||||||
|
|
||||||
table.input_close() # close ASCII tables
|
|
||||||
output.close() # close ASCII tables
|
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -68,7 +68,7 @@ for name in filenames:
|
||||||
# ------------------------------------------ assemble header --------------------------------------
|
# ------------------------------------------ assemble header --------------------------------------
|
||||||
|
|
||||||
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
|
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
|
||||||
table.labels_append(['%i_Cauchy'%(i+1) for i in xrange(9)]) # extend ASCII header with new labels
|
table.labels_append(['{}_Cauchy'.format(i+1) for i in xrange(9)]) # extend ASCII header with new labels
|
||||||
table.head_write()
|
table.head_write()
|
||||||
|
|
||||||
# ------------------------------------------ process data ------------------------------------------
|
# ------------------------------------------ process data ------------------------------------------
|
||||||
|
|
|
||||||
|
|
@ -17,6 +17,7 @@ parser = OptionParser(option_class=damask.extendableOption, usage='%prog options
|
||||||
Add quaternion and/or Bunge Euler angle representation of crystal lattice orientation.
|
Add quaternion and/or Bunge Euler angle representation of crystal lattice orientation.
|
||||||
Orientation is given by quaternion, Euler angles, rotation matrix, or crystal frame coordinates
|
Orientation is given by quaternion, Euler angles, rotation matrix, or crystal frame coordinates
|
||||||
(i.e. component vectors of rotation matrix).
|
(i.e. component vectors of rotation matrix).
|
||||||
|
Additional (globally fixed) rotations of the lab frame and/or crystal frame can be applied.
|
||||||
|
|
||||||
""", version = scriptID)
|
""", version = scriptID)
|
||||||
|
|
||||||
|
|
@ -24,24 +25,27 @@ outputChoices = ['quaternion','rodrigues','eulers']
|
||||||
parser.add_option('-o', '--output',
|
parser.add_option('-o', '--output',
|
||||||
dest = 'output',
|
dest = 'output',
|
||||||
action = 'extend', metavar = '<string LIST>',
|
action = 'extend', metavar = '<string LIST>',
|
||||||
help = 'output orientation formats {%s}'%(','.join(outputChoices)))
|
help = 'output orientation formats {{{}}}'.format(', '.join(outputChoices)))
|
||||||
parser.add_option('-r', '--rotation',
|
|
||||||
dest='rotation',
|
|
||||||
type = 'float', nargs = 4, metavar = ' '.join(['float']*4),
|
|
||||||
help = 'angle and axis to (pre)rotate orientation')
|
|
||||||
|
|
||||||
parser.add_option('-s', '--symmetry',
|
parser.add_option('-s', '--symmetry',
|
||||||
dest = 'symmetry',
|
dest = 'symmetry',
|
||||||
type = 'choice', choices = damask.Symmetry.lattices[1:], metavar='string',
|
type = 'choice', choices = damask.Symmetry.lattices[1:], metavar='string',
|
||||||
help = 'crystal symmetry [%default] {{{}}} '.format(', '.join(damask.Symmetry.lattices[1:])))
|
help = 'crystal symmetry [%default] {{{}}} '.format(', '.join(damask.Symmetry.lattices[1:])))
|
||||||
|
parser.add_option('-d', '--degrees',
|
||||||
|
dest = 'degrees',
|
||||||
|
action = 'store_true',
|
||||||
|
help = 'angles are given in degrees [%default]')
|
||||||
|
parser.add_option('-R', '--labrotation',
|
||||||
|
dest='labrotation',
|
||||||
|
type = 'float', nargs = 4, metavar = ' '.join(['float']*4),
|
||||||
|
help = 'angle and axis of additional lab frame rotation')
|
||||||
|
parser.add_option('-r', '--crystalrotation',
|
||||||
|
dest='crystalrotation',
|
||||||
|
type = 'float', nargs = 4, metavar = ' '.join(['float']*4),
|
||||||
|
help = 'angle and axis of additional crystal frame rotation')
|
||||||
parser.add_option('-e', '--eulers',
|
parser.add_option('-e', '--eulers',
|
||||||
dest = 'eulers',
|
dest = 'eulers',
|
||||||
type = 'string', metavar = 'string',
|
type = 'string', metavar = 'string',
|
||||||
help = 'Euler angles label')
|
help = 'Euler angles label')
|
||||||
parser.add_option('-d', '--degrees',
|
|
||||||
dest = 'degrees',
|
|
||||||
action = 'store_true',
|
|
||||||
help = 'Euler angles are given in degrees [%default]')
|
|
||||||
parser.add_option('-m', '--matrix',
|
parser.add_option('-m', '--matrix',
|
||||||
dest = 'matrix',
|
dest = 'matrix',
|
||||||
type = 'string', metavar = 'string',
|
type = 'string', metavar = 'string',
|
||||||
|
|
@ -65,7 +69,8 @@ parser.add_option('-q', '--quaternion',
|
||||||
|
|
||||||
parser.set_defaults(output = [],
|
parser.set_defaults(output = [],
|
||||||
symmetry = damask.Symmetry.lattices[-1],
|
symmetry = damask.Symmetry.lattices[-1],
|
||||||
rotation = (0.,1.,1.,1.), # no rotation about 1,1,1
|
labrotation = (0.,1.,1.,1.), # no rotation about 1,1,1
|
||||||
|
crystalrotation = (0.,1.,1.,1.), # no rotation about 1,1,1
|
||||||
degrees = False,
|
degrees = False,
|
||||||
)
|
)
|
||||||
|
|
||||||
|
|
@ -91,16 +96,16 @@ if np.sum(input) != 1: parser.error('needs exactly one input format.')
|
||||||
(options.quaternion,4,'quaternion'),
|
(options.quaternion,4,'quaternion'),
|
||||||
][np.where(input)[0][0]] # select input label that was requested
|
][np.where(input)[0][0]] # select input label that was requested
|
||||||
toRadians = math.pi/180.0 if options.degrees else 1.0 # rescale degrees to radians
|
toRadians = math.pi/180.0 if options.degrees else 1.0 # rescale degrees to radians
|
||||||
r = damask.Quaternion().fromAngleAxis(toRadians*options.rotation[0],options.rotation[1:]) # pre-rotation
|
r = damask.Quaternion().fromAngleAxis(toRadians*options.crystalrotation[0],options.crystalrotation[1:]) # crystal frame rotation
|
||||||
|
R = damask.Quaternion().fromAngleAxis(toRadians*options. labrotation[0],options. labrotation[1:]) # lab frame rotation
|
||||||
|
|
||||||
# --- loop over input files ------------------------------------------------------------------------
|
# --- loop over input files ------------------------------------------------------------------------
|
||||||
|
|
||||||
if filenames == []: filenames = [None]
|
if filenames == []: filenames = [None]
|
||||||
|
|
||||||
for name in filenames:
|
for name in filenames:
|
||||||
try:
|
try: table = damask.ASCIItable(name = name,
|
||||||
table = damask.ASCIItable(name = name,
|
buffered = False)
|
||||||
buffered = False)
|
|
||||||
except: continue
|
except: continue
|
||||||
damask.util.report(scriptName,name)
|
damask.util.report(scriptName,name)
|
||||||
|
|
||||||
|
|
@ -126,9 +131,9 @@ for name in filenames:
|
||||||
|
|
||||||
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
|
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
|
||||||
for output in options.output:
|
for output in options.output:
|
||||||
if output == 'quaternion': table.labels_append(['{}_quat({})'.format( i+1,options.symmetry) for i in xrange(4)])
|
if output == 'quaternion': table.labels_append(['{}_{}_{}({})'.format(i+1,'quat',options.symmetry,label) for i in xrange(4)])
|
||||||
if output == 'rodrigues': table.labels_append(['{}_rodrigues({})'.format(i+1,options.symmetry) for i in xrange(3)])
|
elif output == 'rodrigues': table.labels_append(['{}_{}_{}({})'.format(i+1,'rodr',options.symmetry,label) for i in xrange(3)])
|
||||||
if output == 'eulers': table.labels_append(['{}_eulers({})'.format( i+1,options.symmetry) for i in xrange(3)])
|
elif output == 'eulers': table.labels_append(['{}_{}_{}({})'.format(i+1,'eulr',options.symmetry,label) for i in xrange(3)])
|
||||||
table.head_write()
|
table.head_write()
|
||||||
|
|
||||||
# ------------------------------------------ process data ------------------------------------------
|
# ------------------------------------------ process data ------------------------------------------
|
||||||
|
|
@ -150,12 +155,12 @@ for name in filenames:
|
||||||
o = damask.Orientation(quaternion = np.array(map(float,table.data[column:column+4])),
|
o = damask.Orientation(quaternion = np.array(map(float,table.data[column:column+4])),
|
||||||
symmetry = options.symmetry).reduced()
|
symmetry = options.symmetry).reduced()
|
||||||
|
|
||||||
o.quaternion = r*o.quaternion
|
o.quaternion = r*o.quaternion*R # apply additional lab and crystal frame rotations
|
||||||
|
|
||||||
for output in options.output:
|
for output in options.output:
|
||||||
if output == 'quaternion': table.data_append(o.asQuaternion())
|
if output == 'quaternion': table.data_append(o.asQuaternion())
|
||||||
if output == 'rodrigues': table.data_append(o.asRodrigues())
|
elif output == 'rodrigues': table.data_append(o.asRodrigues())
|
||||||
if output == 'eulers': table.data_append(o.asEulers('Bunge', degrees=options.degrees))
|
elif output == 'eulers': table.data_append(o.asEulers('Bunge', degrees=options.degrees))
|
||||||
outputAlive = table.data_write() # output processed line
|
outputAlive = table.data_write() # output processed line
|
||||||
|
|
||||||
# ------------------------------------------ output finalization -----------------------------------
|
# ------------------------------------------ output finalization -----------------------------------
|
||||||
|
|
|
||||||
|
|
@ -49,7 +49,7 @@ parser.set_defaults(right = (0.0,0.0,0.0))
|
||||||
(options,filename) = parser.parse_args()
|
(options,filename) = parser.parse_args()
|
||||||
|
|
||||||
if options.format not in outtypes:
|
if options.format not in outtypes:
|
||||||
parser.error('invalid format: "%s" (can be %s).'%(options.format,', '.join(outtypes)))
|
parser.error('invalid format: "{}" (choices: {}).'.format(options.format,', '.join(outtypes)))
|
||||||
|
|
||||||
if options.N < 2:
|
if options.N < 2:
|
||||||
parser.error('too few steps (need at least 2).')
|
parser.error('too few steps (need at least 2).')
|
||||||
|
|
@ -59,10 +59,9 @@ if options.trim[0] < -1.0 or \
|
||||||
options.trim[0] >= options.trim[1]:
|
options.trim[0] >= options.trim[1]:
|
||||||
parser.error('invalid trim range (-1 +1).')
|
parser.error('invalid trim range (-1 +1).')
|
||||||
|
|
||||||
|
name = options.format if filename == [] \
|
||||||
name = options.format if filename[0] is None\
|
|
||||||
else filename[0]
|
else filename[0]
|
||||||
output = sys.stdout if filename[0] is None\
|
output = sys.stdout if filename == [] \
|
||||||
else open(os.path.basename(filename[0])+extensions[outtypes.index(options.format)],'w')
|
else open(os.path.basename(filename[0])+extensions[outtypes.index(options.format)],'w')
|
||||||
|
|
||||||
colorLeft = damask.Color(options.colormodel.upper(), list(options.left))
|
colorLeft = damask.Color(options.colormodel.upper(), list(options.left))
|
||||||
|
|
|
||||||
Loading…
Reference in New Issue