merged more recent development

This commit is contained in:
Zhuowen Zhao 2018-03-09 18:19:41 -05:00
commit 361cf680a5
89 changed files with 2831 additions and 1085 deletions

View File

@ -9,9 +9,23 @@ OUTFILE="system_report.txt"
echo =========================================== echo ===========================================
echo + Generating $OUTFILE echo + Generating $OUTFILE
echo + Send to damask@mpie.de for support echo + Send to damask@mpie.de for support
echo + view with \'cat $OUTFILE\'
echo =========================================== echo ===========================================
function getDetails {
if which $1 &> /dev/null; then
echo ----------------------------------------------------------------------------------------------
echo $1:
echo ----------------------------------------------------------------------------------------------
echo + location:
which $1
echo + $1 $2:
$1 $2
echo -e '\n'
else
echo $ does not exist
fi
}
# redirect STDOUT and STDERR to logfile # redirect STDOUT and STDERR to logfile
# https://stackoverflow.com/questions/11229385/redirect-all-output-in-a-bash-script-when-using-set-x^ # https://stackoverflow.com/questions/11229385/redirect-all-output-in-a-bash-script-when-using-set-x^
@ -28,13 +42,19 @@ echo
echo ============================================================================================== echo ==============================================================================================
echo DAMASK settings echo DAMASK settings
echo ============================================================================================== echo ==============================================================================================
echo ----------------------------------------------------------------------------------------------
echo DAMASK_ROOT: echo DAMASK_ROOT:
echo ----------------------------------------------------------------------------------------------
echo $DAMASK_ROOT echo $DAMASK_ROOT
echo echo
echo ----------------------------------------------------------------------------------------------
echo Version: echo Version:
echo ----------------------------------------------------------------------------------------------
cat VERSION cat VERSION
echo echo
echo ----------------------------------------------------------------------------------------------
echo Settings in CONFIG: echo Settings in CONFIG:
echo ----------------------------------------------------------------------------------------------
cat CONFIG cat CONFIG
echo echo
echo ============================================================================================== echo ==============================================================================================
@ -46,6 +66,8 @@ echo PATH: $PATH
echo LD_LIBRARY_PATH: $LD_LIBRARY_PATH echo LD_LIBRARY_PATH: $LD_LIBRARY_PATH
echo PYTHONPATH: $PYTHONPATH echo PYTHONPATH: $PYTHONPATH
echo SHELL: $SHELL echo SHELL: $SHELL
echo PETSC_ARCH: $PETSC_ARCH
echo PETSC_DIR: $PETSC_DIR
echo echo
echo ============================================================================================== echo ==============================================================================================
echo Python echo Python
@ -53,16 +75,14 @@ echo ===========================================================================
DEFAULT_PYTHON=python2.7 DEFAULT_PYTHON=python2.7
for executable in python python2 python3 python2.7; do for executable in python python2 python3 python2.7; do
if which $executable &> /dev/null; then getDetails $executable '--version'
echo $executable version: $($executable --version 2>&1)
else
echo $executable does not exist
fi
done done
echo Details on $DEFAULT_PYTHON: $(ls -la $(which $DEFAULT_PYTHON))
echo
for module in numpy scipy;do
echo ---------------------------------------------------------------------------------------------- echo ----------------------------------------------------------------------------------------------
echo Details on $DEFAULT_PYTHON:
echo ----------------------------------------------------------------------------------------------
echo $(ls -la $(which $DEFAULT_PYTHON))
for module in numpy scipy;do
echo -e '\n----------------------------------------------------------------------------------------------'
echo $module echo $module
echo ---------------------------------------------------------------------------------------------- echo ----------------------------------------------------------------------------------------------
$DEFAULT_PYTHON -c "import $module; \ $DEFAULT_PYTHON -c "import $module; \
@ -86,31 +106,36 @@ echo ===========================================================================
echo GCC echo GCC
echo ============================================================================================== echo ==============================================================================================
for executable in gcc g++ gfortran ;do for executable in gcc g++ gfortran ;do
if which $executable &> /dev/null; then getDetails $executable '--version'
echo $(which $executable) version: $($executable --version 2>&1)
else
echo $executable does not exist
fi
done done
echo echo
echo ============================================================================================== echo ==============================================================================================
echo Intel Compiler Suite echo Intel Compiler Suite
echo ============================================================================================== echo ==============================================================================================
for executable in icc icpc ifort ;do for executable in icc icpc ifort ;do
if which $executable &> /dev/null; then getDetails $executable '--version'
echo $(which $executable) version: $($executable --version 2>&1)
else
echo $executable does not exist
fi
done done
echo echo
echo ============================================================================================== echo ==============================================================================================
echo MPI Wrappers echo MPI Wrappers
echo ============================================================================================== echo ==============================================================================================
for executable in mpicc mpiCC mpicxx mpicxx mpifort mpif90 mpif77; do for executable in mpicc mpiCC mpic++ mpicpc mpicxx mpifort mpif90 mpif77; do
if which $executable &> /dev/null; then getDetails $executable '-show'
echo $(which $executable) version: $($executable --show 2>&1)
else
echo $executable does not exist
fi
done done
echo
echo ==============================================================================================
echo MPI Launchers
echo ==============================================================================================
for executable in mpirun mpiexec; do
getDetails $executable '--version'
done
echo
echo ==============================================================================================
echo Abaqus
echo ==============================================================================================
cd installation/mods_Abaqus # to have the right environment file
for executable in abaqus abq2016 abq2017; do
getDetails $executable 'information=all'
done
cd ../..

View File

@ -1,4 +1,4 @@
Copyright 2011-17 Max-Planck-Institut für Eisenforschung GmbH Copyright 2011-18 Max-Planck-Institut für Eisenforschung GmbH
DAMASK is free software: you can redistribute it and/or modify DAMASK is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by it under the terms of the GNU General Public License as published by

View File

@ -1 +1 @@
v2.0.1-1035-gd80a255 v2.0.1-1107-g2c3eae6

View File

@ -0,0 +1,7 @@
elasticity hooke
c11 160.0e9
c12 90.0e9
c13 66.0e9
c33 181.7e9
c44 46.5e9

View File

@ -0,0 +1,11 @@
elasticity hooke
c11 100.0e9
c22 100.0e9
c33 100.0e9
c12 0.0e9
c13 0.0e9
c23 0.0e9
c44 50.0e9
c55 50.0e9
c66 50.0e9

View File

@ -0,0 +1,4 @@
elasticity hooke
c11 100.0e9
c12 0.0e9

View File

@ -0,0 +1,3 @@
reference_temperature 300.0
specific_heat 1
mass_density 1

View File

@ -0,0 +1,8 @@
(kinematics) thermal_expansion
thermal_expansion11 9.5e-6
thermal_expansion22 9.5e-6
thermal_expansion33 5.6e-6
(source) thermal_externalheat
externalheat_time 0 500 500.001 1000 # 500 secs supplying 1 Watt, then removing 1 Watt
externalheat_rate 1 1 -1 -1

View File

@ -0,0 +1,9 @@
(kinematics) thermal_expansion
thermal_expansion11 5e-6
thermal_expansion22 10e-6
thermal_expansion33 20e-6
lattice_structure orthorhombic
(source) thermal_externalheat
externalheat_time 0 500 500.001 1000 # 500 secs supplying 1 Watt, then removing 1 Watt
externalheat_rate 1 1 -1 -1

View File

@ -0,0 +1,6 @@
(kinematics) thermal_expansion
thermal_expansion11 10e-6
(source) thermal_externalheat
externalheat_time 0 500 500.001 1000 # 500 secs supplying 1 Watt, then removing 1 Watt
externalheat_rate 1 1 -1 -1

View File

@ -0,0 +1,427 @@
9 header
geom_addPrimitive v2.0.1-1073-gc544fa1 -c 32 32 32 -d -16 -16 -16 inclusion.geom
geom_canvas v2.0.1-1073-gc544fa1 -o 16 16 16 -g 32 32 32
geom_translate v2.0.1-1073-gc544fa1 -s 1,2,2,6 Ti_Ti.geom
geom_pack v2.0.1-1073-gc544fa1 Ti_Ti.geom isotropic_anisotropic.geom isotropic_isotropic.geom isotropic_rotated.geom
grid a 32 b 32 c 32
size x 0.5 y 0.5 z 0.5
origin x 0.25 y 0.25 z 0.25
homogenization 1
microstructures 2
8623 of 2
2 of 6
29 of 2
4 of 6
27 of 2
6 of 6
26 of 2
6 of 6
27 of 2
4 of 6
29 of 2
2 of 6
798 of 2
2 of 6
28 of 2
6 of 6
25 of 2
8 of 6
24 of 2
8 of 6
23 of 2
10 of 6
22 of 2
10 of 6
23 of 2
8 of 6
24 of 2
8 of 6
25 of 2
6 of 6
28 of 2
2 of 6
701 of 2
4 of 6
26 of 2
8 of 6
23 of 2
10 of 6
22 of 2
10 of 6
21 of 2
12 of 6
20 of 2
12 of 6
20 of 2
12 of 6
20 of 2
12 of 6
21 of 2
10 of 6
22 of 2
10 of 6
23 of 2
8 of 6
26 of 2
4 of 6
637 of 2
2 of 6
27 of 2
8 of 6
23 of 2
10 of 6
21 of 2
12 of 6
20 of 2
12 of 6
20 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
20 of 2
12 of 6
20 of 2
12 of 6
21 of 2
10 of 6
23 of 2
8 of 6
27 of 2
2 of 6
604 of 2
6 of 6
24 of 2
10 of 6
21 of 2
12 of 6
20 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
20 of 2
12 of 6
21 of 2
10 of 6
24 of 2
6 of 6
572 of 2
2 of 6
27 of 2
8 of 6
23 of 2
10 of 6
21 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
17 of 2
16 of 6
16 of 2
16 of 6
17 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
21 of 2
10 of 6
23 of 2
8 of 6
27 of 2
2 of 6
541 of 2
4 of 6
26 of 2
8 of 6
22 of 2
12 of 6
20 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
17 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
17 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
20 of 2
12 of 6
22 of 2
8 of 6
26 of 2
4 of 6
539 of 2
6 of 6
24 of 2
10 of 6
21 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
17 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
17 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
21 of 2
10 of 6
24 of 2
6 of 6
538 of 2
6 of 6
24 of 2
10 of 6
21 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
17 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
17 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
21 of 2
10 of 6
24 of 2
6 of 6
539 of 2
4 of 6
26 of 2
8 of 6
22 of 2
12 of 6
20 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
17 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
16 of 2
16 of 6
17 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
20 of 2
12 of 6
22 of 2
8 of 6
26 of 2
4 of 6
541 of 2
2 of 6
27 of 2
8 of 6
23 of 2
10 of 6
21 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
17 of 2
16 of 6
16 of 2
16 of 6
17 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
21 of 2
10 of 6
23 of 2
8 of 6
27 of 2
2 of 6
572 of 2
6 of 6
24 of 2
10 of 6
21 of 2
12 of 6
20 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
20 of 2
12 of 6
21 of 2
10 of 6
24 of 2
6 of 6
604 of 2
2 of 6
27 of 2
8 of 6
23 of 2
10 of 6
21 of 2
12 of 6
20 of 2
12 of 6
20 of 2
12 of 6
19 of 2
14 of 6
18 of 2
14 of 6
19 of 2
12 of 6
20 of 2
12 of 6
20 of 2
12 of 6
21 of 2
10 of 6
23 of 2
8 of 6
27 of 2
2 of 6
637 of 2
4 of 6
26 of 2
8 of 6
23 of 2
10 of 6
22 of 2
10 of 6
21 of 2
12 of 6
20 of 2
12 of 6
20 of 2
12 of 6
20 of 2
12 of 6
21 of 2
10 of 6
22 of 2
10 of 6
23 of 2
8 of 6
26 of 2
4 of 6
701 of 2
2 of 6
28 of 2
6 of 6
25 of 2
8 of 6
24 of 2
8 of 6
23 of 2
10 of 6
22 of 2
10 of 6
23 of 2
8 of 6
24 of 2
8 of 6
25 of 2
6 of 6
28 of 2
2 of 6
798 of 2
2 of 6
29 of 2
4 of 6
27 of 2
6 of 6
26 of 2
6 of 6
27 of 2
4 of 6
29 of 2
2 of 6
8623 of 2

View File

@ -0,0 +1,427 @@
9 header
geom_addPrimitive v2.0.1-1073-gc544fa1 -c 32 32 32 -d -16 -16 -16 inclusion.geom
geom_canvas v2.0.1-1073-gc544fa1 -o 16 16 16 -g 32 32 32
geom_translate v2.0.1-1073-gc544fa1 -s 2,4 isotropic_anisotropic.geom
geom_pack v2.0.1-1073-gc544fa1 Ti_Ti.geom isotropic_anisotropic.geom isotropic_isotropic.geom isotropic_rotated.geom
grid a 32 b 32 c 32
size x 0.5 y 0.5 z 0.5
origin x 0.25 y 0.25 z 0.25
homogenization 1
microstructures 2
8623 of 1
2 of 4
29 of 1
4 of 4
27 of 1
6 of 4
26 of 1
6 of 4
27 of 1
4 of 4
29 of 1
2 of 4
798 of 1
2 of 4
28 of 1
6 of 4
25 of 1
8 of 4
24 of 1
8 of 4
23 of 1
10 of 4
22 of 1
10 of 4
23 of 1
8 of 4
24 of 1
8 of 4
25 of 1
6 of 4
28 of 1
2 of 4
701 of 1
4 of 4
26 of 1
8 of 4
23 of 1
10 of 4
22 of 1
10 of 4
21 of 1
12 of 4
20 of 1
12 of 4
20 of 1
12 of 4
20 of 1
12 of 4
21 of 1
10 of 4
22 of 1
10 of 4
23 of 1
8 of 4
26 of 1
4 of 4
637 of 1
2 of 4
27 of 1
8 of 4
23 of 1
10 of 4
21 of 1
12 of 4
20 of 1
12 of 4
20 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
20 of 1
12 of 4
20 of 1
12 of 4
21 of 1
10 of 4
23 of 1
8 of 4
27 of 1
2 of 4
604 of 1
6 of 4
24 of 1
10 of 4
21 of 1
12 of 4
20 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
20 of 1
12 of 4
21 of 1
10 of 4
24 of 1
6 of 4
572 of 1
2 of 4
27 of 1
8 of 4
23 of 1
10 of 4
21 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
17 of 1
16 of 4
16 of 1
16 of 4
17 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
21 of 1
10 of 4
23 of 1
8 of 4
27 of 1
2 of 4
541 of 1
4 of 4
26 of 1
8 of 4
22 of 1
12 of 4
20 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
17 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
17 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
20 of 1
12 of 4
22 of 1
8 of 4
26 of 1
4 of 4
539 of 1
6 of 4
24 of 1
10 of 4
21 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
17 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
17 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
21 of 1
10 of 4
24 of 1
6 of 4
538 of 1
6 of 4
24 of 1
10 of 4
21 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
17 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
17 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
21 of 1
10 of 4
24 of 1
6 of 4
539 of 1
4 of 4
26 of 1
8 of 4
22 of 1
12 of 4
20 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
17 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
16 of 1
16 of 4
17 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
20 of 1
12 of 4
22 of 1
8 of 4
26 of 1
4 of 4
541 of 1
2 of 4
27 of 1
8 of 4
23 of 1
10 of 4
21 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
17 of 1
16 of 4
16 of 1
16 of 4
17 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
21 of 1
10 of 4
23 of 1
8 of 4
27 of 1
2 of 4
572 of 1
6 of 4
24 of 1
10 of 4
21 of 1
12 of 4
20 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
20 of 1
12 of 4
21 of 1
10 of 4
24 of 1
6 of 4
604 of 1
2 of 4
27 of 1
8 of 4
23 of 1
10 of 4
21 of 1
12 of 4
20 of 1
12 of 4
20 of 1
12 of 4
19 of 1
14 of 4
18 of 1
14 of 4
19 of 1
12 of 4
20 of 1
12 of 4
20 of 1
12 of 4
21 of 1
10 of 4
23 of 1
8 of 4
27 of 1
2 of 4
637 of 1
4 of 4
26 of 1
8 of 4
23 of 1
10 of 4
22 of 1
10 of 4
21 of 1
12 of 4
20 of 1
12 of 4
20 of 1
12 of 4
20 of 1
12 of 4
21 of 1
10 of 4
22 of 1
10 of 4
23 of 1
8 of 4
26 of 1
4 of 4
701 of 1
2 of 4
28 of 1
6 of 4
25 of 1
8 of 4
24 of 1
8 of 4
23 of 1
10 of 4
22 of 1
10 of 4
23 of 1
8 of 4
24 of 1
8 of 4
25 of 1
6 of 4
28 of 1
2 of 4
798 of 1
2 of 4
29 of 1
4 of 4
27 of 1
6 of 4
26 of 1
6 of 4
27 of 1
4 of 4
29 of 1
2 of 4
8623 of 1

View File

@ -0,0 +1,427 @@
9 header
geom_addPrimitive v2.0.1-1073-gc544fa1 -c 32 32 32 -d -16 -16 -16 inclusion.geom
geom_canvas v2.0.1-1073-gc544fa1 -o 16 16 16 -g 32 32 32
geom_translate v2.0.1-1073-gc544fa1 -s 2,3 isotropic_isotropic.geom
geom_pack v2.0.1-1073-gc544fa1 Ti_Ti.geom isotropic_anisotropic.geom isotropic_isotropic.geom isotropic_rotated.geom
grid a 32 b 32 c 32
size x 0.5 y 0.5 z 0.5
origin x 0.25 y 0.25 z 0.25
homogenization 1
microstructures 2
8623 of 1
2 of 3
29 of 1
4 of 3
27 of 1
6 of 3
26 of 1
6 of 3
27 of 1
4 of 3
29 of 1
2 of 3
798 of 1
2 of 3
28 of 1
6 of 3
25 of 1
8 of 3
24 of 1
8 of 3
23 of 1
10 of 3
22 of 1
10 of 3
23 of 1
8 of 3
24 of 1
8 of 3
25 of 1
6 of 3
28 of 1
2 of 3
701 of 1
4 of 3
26 of 1
8 of 3
23 of 1
10 of 3
22 of 1
10 of 3
21 of 1
12 of 3
20 of 1
12 of 3
20 of 1
12 of 3
20 of 1
12 of 3
21 of 1
10 of 3
22 of 1
10 of 3
23 of 1
8 of 3
26 of 1
4 of 3
637 of 1
2 of 3
27 of 1
8 of 3
23 of 1
10 of 3
21 of 1
12 of 3
20 of 1
12 of 3
20 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
20 of 1
12 of 3
20 of 1
12 of 3
21 of 1
10 of 3
23 of 1
8 of 3
27 of 1
2 of 3
604 of 1
6 of 3
24 of 1
10 of 3
21 of 1
12 of 3
20 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
20 of 1
12 of 3
21 of 1
10 of 3
24 of 1
6 of 3
572 of 1
2 of 3
27 of 1
8 of 3
23 of 1
10 of 3
21 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
17 of 1
16 of 3
16 of 1
16 of 3
17 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
21 of 1
10 of 3
23 of 1
8 of 3
27 of 1
2 of 3
541 of 1
4 of 3
26 of 1
8 of 3
22 of 1
12 of 3
20 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
17 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
17 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
20 of 1
12 of 3
22 of 1
8 of 3
26 of 1
4 of 3
539 of 1
6 of 3
24 of 1
10 of 3
21 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
17 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
17 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
21 of 1
10 of 3
24 of 1
6 of 3
538 of 1
6 of 3
24 of 1
10 of 3
21 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
17 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
17 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
21 of 1
10 of 3
24 of 1
6 of 3
539 of 1
4 of 3
26 of 1
8 of 3
22 of 1
12 of 3
20 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
17 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
16 of 1
16 of 3
17 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
20 of 1
12 of 3
22 of 1
8 of 3
26 of 1
4 of 3
541 of 1
2 of 3
27 of 1
8 of 3
23 of 1
10 of 3
21 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
17 of 1
16 of 3
16 of 1
16 of 3
17 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
21 of 1
10 of 3
23 of 1
8 of 3
27 of 1
2 of 3
572 of 1
6 of 3
24 of 1
10 of 3
21 of 1
12 of 3
20 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
20 of 1
12 of 3
21 of 1
10 of 3
24 of 1
6 of 3
604 of 1
2 of 3
27 of 1
8 of 3
23 of 1
10 of 3
21 of 1
12 of 3
20 of 1
12 of 3
20 of 1
12 of 3
19 of 1
14 of 3
18 of 1
14 of 3
19 of 1
12 of 3
20 of 1
12 of 3
20 of 1
12 of 3
21 of 1
10 of 3
23 of 1
8 of 3
27 of 1
2 of 3
637 of 1
4 of 3
26 of 1
8 of 3
23 of 1
10 of 3
22 of 1
10 of 3
21 of 1
12 of 3
20 of 1
12 of 3
20 of 1
12 of 3
20 of 1
12 of 3
21 of 1
10 of 3
22 of 1
10 of 3
23 of 1
8 of 3
26 of 1
4 of 3
701 of 1
2 of 3
28 of 1
6 of 3
25 of 1
8 of 3
24 of 1
8 of 3
23 of 1
10 of 3
22 of 1
10 of 3
23 of 1
8 of 3
24 of 1
8 of 3
25 of 1
6 of 3
28 of 1
2 of 3
798 of 1
2 of 3
29 of 1
4 of 3
27 of 1
6 of 3
26 of 1
6 of 3
27 of 1
4 of 3
29 of 1
2 of 3
8623 of 1

View File

@ -0,0 +1,427 @@
9 header
geom_addPrimitive v2.0.1-1073-gc544fa1 -c 32 32 32 -d -16 -16 -16 inclusion.geom
geom_canvas v2.0.1-1073-gc544fa1 -o 16 16 16 -g 32 32 32
geom_translate v2.0.1-1073-gc544fa1 -s 2,5 isotropic_rotated.geom
geom_pack v2.0.1-1073-gc544fa1 Ti_Ti.geom isotropic_anisotropic.geom isotropic_isotropic.geom isotropic_rotated.geom
grid a 32 b 32 c 32
size x 0.5 y 0.5 z 0.5
origin x 0.25 y 0.25 z 0.25
homogenization 1
microstructures 2
8623 of 1
2 of 5
29 of 1
4 of 5
27 of 1
6 of 5
26 of 1
6 of 5
27 of 1
4 of 5
29 of 1
2 of 5
798 of 1
2 of 5
28 of 1
6 of 5
25 of 1
8 of 5
24 of 1
8 of 5
23 of 1
10 of 5
22 of 1
10 of 5
23 of 1
8 of 5
24 of 1
8 of 5
25 of 1
6 of 5
28 of 1
2 of 5
701 of 1
4 of 5
26 of 1
8 of 5
23 of 1
10 of 5
22 of 1
10 of 5
21 of 1
12 of 5
20 of 1
12 of 5
20 of 1
12 of 5
20 of 1
12 of 5
21 of 1
10 of 5
22 of 1
10 of 5
23 of 1
8 of 5
26 of 1
4 of 5
637 of 1
2 of 5
27 of 1
8 of 5
23 of 1
10 of 5
21 of 1
12 of 5
20 of 1
12 of 5
20 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
20 of 1
12 of 5
20 of 1
12 of 5
21 of 1
10 of 5
23 of 1
8 of 5
27 of 1
2 of 5
604 of 1
6 of 5
24 of 1
10 of 5
21 of 1
12 of 5
20 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
20 of 1
12 of 5
21 of 1
10 of 5
24 of 1
6 of 5
572 of 1
2 of 5
27 of 1
8 of 5
23 of 1
10 of 5
21 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
17 of 1
16 of 5
16 of 1
16 of 5
17 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
21 of 1
10 of 5
23 of 1
8 of 5
27 of 1
2 of 5
541 of 1
4 of 5
26 of 1
8 of 5
22 of 1
12 of 5
20 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
17 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
17 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
20 of 1
12 of 5
22 of 1
8 of 5
26 of 1
4 of 5
539 of 1
6 of 5
24 of 1
10 of 5
21 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
17 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
17 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
21 of 1
10 of 5
24 of 1
6 of 5
538 of 1
6 of 5
24 of 1
10 of 5
21 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
17 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
17 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
21 of 1
10 of 5
24 of 1
6 of 5
539 of 1
4 of 5
26 of 1
8 of 5
22 of 1
12 of 5
20 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
17 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
16 of 1
16 of 5
17 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
20 of 1
12 of 5
22 of 1
8 of 5
26 of 1
4 of 5
541 of 1
2 of 5
27 of 1
8 of 5
23 of 1
10 of 5
21 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
17 of 1
16 of 5
16 of 1
16 of 5
17 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
21 of 1
10 of 5
23 of 1
8 of 5
27 of 1
2 of 5
572 of 1
6 of 5
24 of 1
10 of 5
21 of 1
12 of 5
20 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
20 of 1
12 of 5
21 of 1
10 of 5
24 of 1
6 of 5
604 of 1
2 of 5
27 of 1
8 of 5
23 of 1
10 of 5
21 of 1
12 of 5
20 of 1
12 of 5
20 of 1
12 of 5
19 of 1
14 of 5
18 of 1
14 of 5
19 of 1
12 of 5
20 of 1
12 of 5
20 of 1
12 of 5
21 of 1
10 of 5
23 of 1
8 of 5
27 of 1
2 of 5
637 of 1
4 of 5
26 of 1
8 of 5
23 of 1
10 of 5
22 of 1
10 of 5
21 of 1
12 of 5
20 of 1
12 of 5
20 of 1
12 of 5
20 of 1
12 of 5
21 of 1
10 of 5
22 of 1
10 of 5
23 of 1
8 of 5
26 of 1
4 of 5
701 of 1
2 of 5
28 of 1
6 of 5
25 of 1
8 of 5
24 of 1
8 of 5
23 of 1
10 of 5
22 of 1
10 of 5
23 of 1
8 of 5
24 of 1
8 of 5
25 of 1
6 of 5
28 of 1
2 of 5
798 of 1
2 of 5
29 of 1
4 of 5
27 of 1
6 of 5
26 of 1
6 of 5
27 of 1
4 of 5
29 of 1
2 of 5
8623 of 1

View File

@ -0,0 +1,111 @@
#-------------------#
<homogenization>
#-------------------#
[direct]
type none # isostrain 1 grain
thermal adiabatic # thermal strain (stress) induced mass transport
initialT 300.0
(output) temperature
#-------------------#
<crystallite>
#-------------------#
[aLittleSomething]
(output) texture
(output) f
(output) p
(output) fe
(output) fi
(output) fp
#-------------------#
<phase>
#-------------------#
#.................
[isotropic matrix]
lattice_structure isotropic
plasticity none
{config/elastic_isotropic.config}
{config/thermal.config}
#.................
[Ti matrix]
lattice_structure hex
covera_ratio 1.587
plasticity none
{config/elastic_Ti.config}
{config/thermal.config}
#.................
[isotropic inclusion]
lattice_structure isotropic
plasticity none
{config/elastic_isotropic.config}
{config/thermal.config}
{config/thermalExpansion_isotropic.config}
#.................
[anisotropic inclusion]
lattice_structure orthorhombic
plasticity none
{config/elastic_fullyAnisotropic.config}
{config/thermal.config}
{config/thermalExpansion_fullyAnisotropic.config}
#.................
[Ti inclusion]
lattice_structure hex
covera_ratio 1.587
plasticity none
{config/elastic_Ti.config}
{config/thermal.config}
{config/thermalExpansion_Ti.config}
#--------------------------#
<microstructure>
#--------------------------#
[isotropic matrix]
crystallite 1
(constituent) phase 1 texture 1 fraction 1.0
[Ti matrix]
crystallite 1
(constituent) phase 2 texture 1 fraction 1.0
[isotropic inclusion]
crystallite 1
(constituent) phase 3 texture 1 fraction 1.0
[anisotropic inclusion]
crystallite 1
(constituent) phase 4 texture 1 fraction 1.0
[rotated inclusion]
crystallite 1
(constituent) phase 4 texture 2 fraction 1.0
[Ti inclusion]
crystallite 1
(constituent) phase 5 texture 1 fraction 1.0
#--------------------------#
<texture>
#--------------------------#
[cube]
(gauss) phi1 0.0 Phi 0.0 phi2 0.0 scatter 0.0 fraction 1.0
[rotated]
(gauss) phi1 0.0 Phi 45.0 phi2 0.0 scatter 0.0 fraction 1.0

View File

@ -0,0 +1,5 @@
#spectralsolver polarisation
spectralderivative fwbw_difference
err_div_tolrel 1e-3
itmin 2
petsc_options -mech_snes_type anderson -mech_snes_anderson_beta 1.0 -mech_snes_anderson_restart 10 -thermal_snes_type anderson -thermal_snes_anderson_beta 1.0

View File

@ -0,0 +1,56 @@
#!/usr/local/bin/bash
for geom in $(ls geom/*.geom)
do
base=${geom%.geom}
base=${base#geom/}
name=${base}_thermal
vtr=${base}.vtr
[[ -f ${name}.spectralOut ]] || \
DAMASK_spectral \
--workingdir ./ \
--load thermal.load \
--geom $geom \
> ${name}.out
if [ ! -f postProc/${name}_inc10.txt ]
then
postResults ${name}.spectralOut \
--ho temperature \
--cr f,fe,fi,fp,p \
--split \
--separation x,y,z \
addCauchy postProc/${name}_inc*.txt \
addDeviator postProc/${name}_inc*.txt \
--spherical \
--tensor p,Cauchy \
addDisplacement postProc/${name}_inc*.txt \
--nodal \
fi
geom_check ${geom}
for inc in {00..10}
do
echo "generating postProc/${name}_inc${inc}.vtr"
cp geom/${vtr} postProc/${name}_inc${inc}.vtr
vtk_addRectilinearGridData \
postProc/${name}_inc${inc}.txt \
--inplace \
--vtk postProc/${name}_inc${inc}.vtr \
--data 'sph(p)','sph(Cauchy)',temperature \
--tensor f,fe,fi,fp,p,Cauchy \
vtk_addRectilinearGridData \
postProc/${name}_inc${inc}_nodal.txt \
--inplace \
--vtk postProc/${name}_inc${inc}.vtr \
--data 'avg(f).pos','fluct(f).pos' \
done
done

View File

@ -0,0 +1 @@
Fdot 0 0 0 0 0 0 0 0 0 stress * * * * * * * * * time 1000 incs 10

View File

@ -2,7 +2,7 @@
# -*- coding: UTF-8 no BOM -*- # -*- coding: UTF-8 no BOM -*-
import os,re,sys,collections import os,re,sys,collections
import math # noqa import math,scipy,scipy.linalg # noqa
import numpy as np import numpy as np
from optparse import OptionParser from optparse import OptionParser
import damask import damask

View File

@ -29,9 +29,14 @@ parser.add_option('-m', '--microstructureoffset',
type = 'int', metavar = 'int', type = 'int', metavar = 'int',
help = 'offset (positive or negative) for tagged microstructure indices. '+ help = 'offset (positive or negative) for tagged microstructure indices. '+
'"0" selects maximum microstructure index [%default]') '"0" selects maximum microstructure index [%default]')
parser.add_option('-n', '--nonperiodic',
dest = 'mode',
action = 'store_const', const = 'nearest',
help = 'assume geometry to be non-periodic')
parser.set_defaults(vicinity = 1, parser.set_defaults(vicinity = 1,
offset = 0, offset = 0,
mode = 'wrap',
) )
(options, filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
@ -79,8 +84,8 @@ for name in filenames:
if options.offset == 0: options.offset = microstructure.max() if options.offset == 0: options.offset = microstructure.max()
microstructure = np.where(ndimage.filters.maximum_filter(microstructure,size=1+2*options.vicinity,mode='wrap') == microstructure = np.where(ndimage.filters.maximum_filter(microstructure,size=1+2*options.vicinity,mode=options.mode) ==
ndimage.filters.minimum_filter(microstructure,size=1+2*options.vicinity,mode='wrap'), ndimage.filters.minimum_filter(microstructure,size=1+2*options.vicinity,mode=options.mode),
microstructure, microstructure + options.offset) microstructure, microstructure + options.offset)
newInfo['microstructures'] = microstructure.max() newInfo['microstructures'] = microstructure.max()

View File

@ -251,7 +251,7 @@ for name in filenames:
table.head_read() table.head_read()
errors = [] errors = []
labels = ['phi1','Phi','phi2','intensity'] labels = ['1_euler','2_euler','3_euler','intensity']
for i,index in enumerate(table.label_index(labels)): for i,index in enumerate(table.label_index(labels)):
if index < 0: errors.append('label {} not present.'.format(labels[i])) if index < 0: errors.append('label {} not present.'.format(labels[i]))

View File

@ -166,28 +166,15 @@ for name in filenames:
grainEuler[2,:] *= 360.0 # phi_2 is uniformly distributed grainEuler[2,:] *= 360.0 # phi_2 is uniformly distributed
if not options.selective: if not options.selective:
seeds = np.array([]) n = np.maximum(np.ones(3),np.array(options.grid*options.fraction),
dtype=int,casting='unsafe') # find max grid indices within fraction
while len(seeds) < options.N: meshgrid = np.meshgrid(*map(np.arange,n),indexing='ij') # create a meshgrid within fraction
coords = np.vstack((meshgrid[0],meshgrid[1],meshgrid[2])).reshape(3,n.prod()).T # assemble list of 3D coordinates
theSeeds = np.zeros((options.N,3),dtype=float) # seed positions array seeds = ((random.sample(coords,options.N)+np.random.random(options.N*3).reshape(options.N,3))\
gridpoints = random.sample(range(gridSize),options.N) # choose first N from random permutation of grid positions / \
(n/options.fraction)).T # pick options.N of those, rattle position,
theSeeds[:,0] = (np.mod(gridpoints ,options.grid[0])\ # and rescale to fall within fraction
+np.random.random(options.N)) /options.grid[0]
theSeeds[:,1] = (np.mod(gridpoints// options.grid[0] ,options.grid[1])\
+np.random.random(options.N)) /options.grid[1]
theSeeds[:,2] = (np.mod(gridpoints//(options.grid[1]*options.grid[0]),options.grid[2])\
+np.random.random(options.N)) /options.grid[2]
goodSeeds = theSeeds[np.all(theSeeds<=options.fraction,axis=1)] # pick seeds within threshold fraction
seeds = goodSeeds if len(seeds) == 0 else np.vstack((seeds,goodSeeds))
if len(seeds) > options.N: seeds = seeds[:min(options.N,len(seeds))]
seeds = seeds.T # switch layout to point index as last index
else: else:
seeds = np.zeros((options.N,3),dtype=float) # seed positions array seeds = np.zeros((options.N,3),dtype=float) # seed positions array
seeds[0] = np.random.random(3)*options.grid/max(options.grid) seeds[0] = np.random.random(3)*options.grid/max(options.grid)
i = 1 # start out with one given point i = 1 # start out with one given point

View File

@ -113,7 +113,7 @@ end subroutine CPFEM_initAll
!> @brief allocate the arrays defined in module CPFEM and initialize them !> @brief allocate the arrays defined in module CPFEM and initialize them
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine CPFEM_init subroutine CPFEM_init
#ifdef __GFORTRAN__ #if __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -82,7 +82,7 @@ end subroutine CPFEM_initAll
!> @brief allocate the arrays defined in module CPFEM and initialize them !> @brief allocate the arrays defined in module CPFEM and initialize them
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine CPFEM_init subroutine CPFEM_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -140,36 +140,28 @@ subroutine CPFEM_init
write(rankStr,'(a1,i0)')'_',worldrank write(rankStr,'(a1,i0)')'_',worldrank
call IO_read_intFile(777,'recordedPhase'//trim(rankStr),modelName,size(material_phase)) call IO_read_intFile(777,'recordedPhase'//trim(rankStr),modelName,size(material_phase))
read (777,rec=1) material_phase read (777,rec=1) material_phase; close (777)
close (777)
call IO_read_realFile(777,'convergedF'//trim(rankStr),modelName,size(crystallite_F0)) call IO_read_realFile(777,'convergedF'//trim(rankStr),modelName,size(crystallite_F0))
read (777,rec=1) crystallite_F0 read (777,rec=1) crystallite_F0; close (777)
close (777)
call IO_read_realFile(777,'convergedFp'//trim(rankStr),modelName,size(crystallite_Fp0)) call IO_read_realFile(777,'convergedFp'//trim(rankStr),modelName,size(crystallite_Fp0))
read (777,rec=1) crystallite_Fp0 read (777,rec=1) crystallite_Fp0; close (777)
close (777)
call IO_read_realFile(777,'convergedFi'//trim(rankStr),modelName,size(crystallite_Fi0)) call IO_read_realFile(777,'convergedFi'//trim(rankStr),modelName,size(crystallite_Fi0))
read (777,rec=1) crystallite_Fi0 read (777,rec=1) crystallite_Fi0; close (777)
close (777)
call IO_read_realFile(777,'convergedLp'//trim(rankStr),modelName,size(crystallite_Lp0)) call IO_read_realFile(777,'convergedLp'//trim(rankStr),modelName,size(crystallite_Lp0))
read (777,rec=1) crystallite_Lp0 read (777,rec=1) crystallite_Lp0; close (777)
close (777)
call IO_read_realFile(777,'convergedLi'//trim(rankStr),modelName,size(crystallite_Li0)) call IO_read_realFile(777,'convergedLi'//trim(rankStr),modelName,size(crystallite_Li0))
read (777,rec=1) crystallite_Li0 read (777,rec=1) crystallite_Li0; close (777)
close (777)
call IO_read_realFile(777,'convergeddPdF'//trim(rankStr),modelName,size(crystallite_dPdF0)) call IO_read_realFile(777,'convergeddPdF'//trim(rankStr),modelName,size(crystallite_dPdF0))
read (777,rec=1) crystallite_dPdF0 read (777,rec=1) crystallite_dPdF0; close (777)
close (777)
call IO_read_realFile(777,'convergedTstar'//trim(rankStr),modelName,size(crystallite_Tstar0_v)) call IO_read_realFile(777,'convergedTstar'//trim(rankStr),modelName,size(crystallite_Tstar0_v))
read (777,rec=1) crystallite_Tstar0_v read (777,rec=1) crystallite_Tstar0_v; close (777)
close (777)
call IO_read_realFile(777,'convergedStateConst'//trim(rankStr),modelName) call IO_read_realFile(777,'convergedStateConst'//trim(rankStr),modelName)
m = 0_pInt m = 0_pInt

View File

@ -7,7 +7,7 @@
!> results !> results
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
program DAMASK_spectral program DAMASK_spectral
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -78,8 +78,7 @@ program DAMASK_spectral
FIELD_UNDEFINED_ID, & FIELD_UNDEFINED_ID, &
FIELD_MECH_ID, & FIELD_MECH_ID, &
FIELD_THERMAL_ID, & FIELD_THERMAL_ID, &
FIELD_DAMAGE_ID, & FIELD_DAMAGE_ID
utilities_calcPlasticity
use spectral_mech_Basic use spectral_mech_Basic
use spectral_mech_AL use spectral_mech_AL
use spectral_mech_Polarisation use spectral_mech_Polarisation
@ -157,19 +156,6 @@ program DAMASK_spectral
MPI_finalize, & MPI_finalize, &
MPI_allreduce, & MPI_allreduce, &
PETScFinalize PETScFinalize
!--------------------------------------------------------------------------------------------------
! variables related to stop criterion for yielding
real(pReal) :: plasticWorkOld, plasticWorkNew, & ! plastic work
eqTotalStrainOld, eqTotalStrainNew, & ! total equivalent strain
eqPlasticStrainOld, eqPlasticStrainNew, & ! total equivalent plastic strain
eqStressOld, eqStressNew , & ! equivalent stress
yieldStopValue
real(pReal), dimension(3,3) :: yieldStress,yieldStressOld,yieldStressNew, &
plasticStrainOld, plasticStrainNew, plasticStrainRate
integer(pInt) :: yieldResUnit = 0_pInt
integer(pInt) :: stressstrainUnit = 0_pInt
character(len=13) :: stopFlag
logical :: yieldStop, yieldStopSatisfied
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! init DAMASK (all modules) ! init DAMASK (all modules)
@ -227,8 +213,6 @@ program DAMASK_spectral
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! reading the load case and assign values to the allocated data structure ! reading the load case and assign values to the allocated data structure
yieldStop = .False.
yieldStopSatisfied = .False.
rewind(FILEUNIT) rewind(FILEUNIT)
do do
line = IO_read(FILEUNIT) line = IO_read(FILEUNIT)
@ -303,30 +287,10 @@ program DAMASK_spectral
temp_valueVector(j) = IO_floatValue(line,chunkPos,i+j) temp_valueVector(j) = IO_floatValue(line,chunkPos,i+j)
enddo enddo
loadCases(currentLoadCase)%rotation = math_plain9to33(temp_valueVector) loadCases(currentLoadCase)%rotation = math_plain9to33(temp_valueVector)
case('totalstrain')
yieldStop = .True.
stopFlag = 'totalStrain'
yieldStopValue = IO_floatValue(line,chunkPos,i+1_pInt)
case('plasticstrain')
yieldStop = .True.
stopFlag = 'plasticStrain'
yieldStopValue = IO_floatValue(line,chunkPos,i+1_pInt)
case('plasticwork')
yieldStop = .True.
stopFlag = 'plasticWork'
yieldStopValue = IO_floatValue(line,chunkPos,i+1_pInt)
end select end select
enddo; enddo enddo; enddo
close(FILEUNIT) close(FILEUNIT)
if(yieldStop) then ! initialize variables related to yield stop
yieldStressNew = 0.0_pReal
plasticStrainNew = 0.0_pReal
eqStressNew = 0.0_pReal
eqTotalStrainNew = 0.0_pReal
eqPlasticStrainNew = 0.0_pReal
plasticWorkNew = 0.0_pReal
endif
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! consistency checks and output of load case ! consistency checks and output of load case
loadCases(1)%followFormerTrajectory = .false. ! cannot guess along trajectory for first inc of first currentLoadCase loadCases(1)%followFormerTrajectory = .false. ! cannot guess along trajectory for first inc of first currentLoadCase
@ -443,7 +407,7 @@ program DAMASK_spectral
write(resUnit) 'times:', loadCases%time ! one entry per LoadCase write(resUnit) 'times:', loadCases%time ! one entry per LoadCase
write(resUnit) 'logscales:', loadCases%logscale write(resUnit) 'logscales:', loadCases%logscale
write(resUnit) 'increments:', loadCases%incs ! one entry per LoadCase write(resUnit) 'increments:', loadCases%incs ! one entry per LoadCase
write(resUnit) 'startingIncrement:', restartInc - 1_pInt ! start with writing out the previous inc write(resUnit) 'startingIncrement:', restartInc ! start with writing out the previous inc
write(resUnit) 'eoh' write(resUnit) 'eoh'
close(resUnit) ! end of header close(resUnit) ! end of header
open(newunit=statUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//& open(newunit=statUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//&
@ -518,14 +482,13 @@ program DAMASK_spectral
timeinc = time0 * & timeinc = time0 * &
( (1.0_pReal + loadCases(currentLoadCase)%time/time0 )**(real( inc ,pReal)/& ( (1.0_pReal + loadCases(currentLoadCase)%time/time0 )**(real( inc ,pReal)/&
real(loadCases(currentLoadCase)%incs ,pReal))& real(loadCases(currentLoadCase)%incs ,pReal))&
-(1.0_pReal + loadCases(currentLoadCase)%time/time0 )**(real( (inc-1_pInt),pReal)/& -(1.0_pReal + loadCases(currentLoadCase)%time/time0 )**(real( inc-1_pInt ,pReal)/&
real(loadCases(currentLoadCase)%incs ,pReal))) real(loadCases(currentLoadCase)%incs ,pReal)))
endif endif
endif endif
timeinc = timeinc * real(subStepFactor,pReal)**real(-cutBackLevel,pReal) ! depending on cut back level, decrease time step
timeinc = timeinc / real(subStepFactor,pReal)**real(cutBackLevel,pReal) ! depending on cut back level, decrease time step skipping: if (totalIncsCounter <= restartInc) then ! not yet at restart inc?
skipping: if (totalIncsCounter < restartInc) then ! not yet at restart inc?
time = time + timeinc ! just advance time, skip already performed calculation time = time + timeinc ! just advance time, skip already performed calculation
guess = .true. ! QUESTION:why forced guessing instead of inheriting loadcase preference guess = .true. ! QUESTION:why forced guessing instead of inheriting loadcase preference
else skipping else skipping
@ -549,8 +512,11 @@ program DAMASK_spectral
's: Increment ', inc,'/',loadCases(currentLoadCase)%incs,& 's: Increment ', inc,'/',loadCases(currentLoadCase)%incs,&
'-', stepFraction,'/',subStepFactor**cutBackLevel,& '-', stepFraction,'/',subStepFactor**cutBackLevel,&
' of load case ', currentLoadCase,'/',size(loadCases) ' of load case ', currentLoadCase,'/',size(loadCases)
write(incInfo,'(a,'//IO_intOut(totalIncsCounter)//',a,'//IO_intOut(sum(loadCases%incs))//& write(incInfo,&
',a,'//IO_intOut(stepFraction)//',a,'//IO_intOut(subStepFactor**cutBackLevel)//')') & '(a,'//IO_intOut(totalIncsCounter)//&
',a,'//IO_intOut(sum(loadCases%incs))//&
',a,'//IO_intOut(stepFraction)//&
',a,'//IO_intOut(subStepFactor**cutBackLevel)//')') &
'Increment ',totalIncsCounter,'/',sum(loadCases%incs),& 'Increment ',totalIncsCounter,'/',sum(loadCases%incs),&
'-', stepFraction,'/',subStepFactor**cutBackLevel '-', stepFraction,'/',subStepFactor**cutBackLevel
flush(6) flush(6)
@ -629,7 +595,6 @@ program DAMASK_spectral
stagIter = stagIter + 1_pInt stagIter = stagIter + 1_pInt
stagIterate = stagIter < stagItMax & stagIterate = stagIter < stagItMax &
.and. all(solres(:)%converged) & .and. all(solres(:)%converged) &
.and. .not. all(solres(:)%stagConverged) ! stationary with respect to staggered iteration .and. .not. all(solres(:)%stagConverged) ! stationary with respect to staggered iteration
enddo enddo
@ -675,12 +640,10 @@ program DAMASK_spectral
endif; flush(6) endif; flush(6)
if (mod(inc,loadCases(currentLoadCase)%outputFrequency) == 0_pInt) then ! at output frequency if (mod(inc,loadCases(currentLoadCase)%outputFrequency) == 0_pInt) then ! at output frequency
if (worldrank == 0) &
write(6,'(1/,a)') ' ... writing results to file ......................................' write(6,'(1/,a)') ' ... writing results to file ......................................'
flush(6) flush(6)
call materialpoint_postResults() call materialpoint_postResults()
call MPI_file_seek (resUnit,fileOffset,MPI_SEEK_SET,ierr) call MPI_file_seek (resUnit,fileOffset,MPI_SEEK_SET,ierr)
if (ierr /= 0_pInt) call IO_error(894_pInt, ext_msg='MPI_file_seek') if (ierr /= 0_pInt) call IO_error(894_pInt, ext_msg='MPI_file_seek')
do i=1, size(materialpoint_results,3)/(maxByteOut/(materialpoint_sizeResults*pReal))+1 ! slice the output of my process in chunks not exceeding the limit for one output do i=1, size(materialpoint_results,3)/(maxByteOut/(materialpoint_sizeResults*pReal))+1 ! slice the output of my process in chunks not exceeding the limit for one output
outputIndex=int([(i-1_pInt)*((maxRealOut)/materialpoint_sizeResults)+1_pInt, & outputIndex=int([(i-1_pInt)*((maxRealOut)/materialpoint_sizeResults)+1_pInt, &
@ -699,84 +662,18 @@ program DAMASK_spectral
lastRestartWritten = inc ! QUESTION: first call to CPFEM_general will write? lastRestartWritten = inc ! QUESTION: first call to CPFEM_general will write?
endif endif
else forwarding endif skipping
time = time + timeinc
guess = .true.
endif forwarding
yieldCheck: if(yieldStop) then ! check if it yields or satisfies the certain stop condition
yieldStressOld = yieldStressNew
plasticStrainOld = plasticStrainNew
eqStressOld = eqStressNew
eqTotalStrainOld = eqTotalStrainNew
eqPlasticStrainOld = eqPlasticStrainNew
plasticWorkOld = plasticWorkNew
call utilities_calcPlasticity(yieldStressNew, plasticStrainNew, eqStressNew, eqTotalStrainNew, &
eqPlasticStrainNew, plasticWorkNew, loadCases(currentLoadCase)%rotation)
if (worldrank == 0) then ! output the stress-strain curve to file if yield stop criterion is used
if ((currentLoadCase == 1_pInt) .and. (inc == 1_pInt)) then
open(newunit=stressstrainUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//&
'.stressstrain',form='FORMATTED',status='REPLACE')
write(stressstrainUnit,*) 0.0_pReal, 0.0_pReal
write(stressstrainUnit,*) eqTotalStrainNew, eqStressNew
close(stressstrainUnit)
else
open(newunit=stressstrainUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//&
'.stressstrain',form='FORMATTED', position='APPEND', status='OLD')
write(stressstrainUnit,*) eqTotalStrainNew, eqStressNew
close(stressstrainUnit)
endif
endif
if(stopFlag == 'totalStrain') then
if(eqTotalStrainNew > yieldStopValue) then
yieldStress = yieldStressOld * (eqTotalStrainNew - yieldStopValue)/(eqTotalStrainNew - eqTotalStrainOld) & ! linear interpolation of stress values
+ yieldStressNew * (yieldStopValue - eqTotalStrainOld)/(eqTotalStrainNew - eqTotalStrainOld)
plasticStrainRate = (plasticStrainNew - plasticStrainOld)/(time - time0) ! calculate plastic strain rate
yieldStopSatisfied = .True.
endif
elseif(stopFlag == 'plasticStrain') then
if(eqPlasticStrainNew > yieldStopValue) then
yieldStress = yieldStressOld * (eqPlasticStrainNew - yieldStopValue)/(eqPlasticStrainNew - eqPlasticStrainOld) &
+ yieldStressNew * (yieldStopValue - eqPlasticStrainOld)/(eqPlasticStrainNew - eqPlasticStrainOld)
plasticStrainRate = (plasticStrainNew - plasticStrainOld)/(time - time0)
yieldStopSatisfied = .True.
endif
elseif(stopFlag == 'plasticWork') then
if(plasticWorkNew > yieldStopValue) then
yieldStress = yieldStressOld * (plasticWorkNew - yieldStopValue)/(plasticWorkNew - plasticWorkOld) &
+ yieldStressNew * (yieldStopValue - plasticWorkOld)/(plasticWorkNew - plasticWorkOld)
plasticStrainRate = (plasticStrainNew - plasticStrainOld)/(time - time0)
yieldStopSatisfied = .True.
endif
endif
endif yieldCheck
if (yieldStopSatisfied) then ! when yield, write the yield stress and strain rate to file and quit the job
if (worldrank == 0) then
open(newunit=yieldResUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//&
'.yield',form='FORMATTED',status='REPLACE')
do i = 1_pInt,3_pInt
write(yieldResUnit,*) (yieldStress(i,j), j=1,3)
enddo
do i = 1_pInt,3_pInt
write(yieldResUnit,*) (plasticStrainRate(i,j), j=1,3)
enddo
close(yieldResUnit)
call quit(0_pInt)
endif
endif
enddo incLooping enddo incLooping
enddo loadCaseLooping enddo loadCaseLooping
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! report summary of whole calculation ! report summary of whole calculation
write(6,'(/,a)') ' ###########################################################################' write(6,'(/,a)') ' ###########################################################################'
write(6,'(1x,i6.6,a,i6.6,a,f5.1,a)') convergedCounter, ' out of ', & write(6,'(1x,'//IO_intOut(convergedCounter)//',a,'//IO_intOut(notConvergedCounter + convergedCounter)//',a,f5.1,a)') &
convergedCounter, ' out of ', &
notConvergedCounter + convergedCounter, ' (', & notConvergedCounter + convergedCounter, ' (', &
real(convergedCounter, pReal)/& real(convergedCounter, pReal)/&
real(notConvergedCounter + convergedCounter,pReal)*100.0_pReal, & real(notConvergedCounter + convergedCounter,pReal)*100.0_pReal, &

View File

@ -43,7 +43,7 @@ contains
!> solver the information is provided by the interface module !> solver the information is provided by the interface module
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine FE_init subroutine FE_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -79,20 +79,22 @@ subroutine FE_init
#include "compilation_info.f90" #include "compilation_info.f90"
modelName = getSolverJobName() modelName = getSolverJobName()
#if defined(Spectral) || defined(FEM)
#ifdef Spectral #ifdef Spectral
restartInc = spectralRestartInc restartInc = spectralRestartInc
if(restartInc <= 0_pInt) then #endif
call IO_warning(warning_ID=34_pInt) #ifdef FEM
restartInc = 1_pInt
endif
restartRead = restartInc > 1_pInt ! only read in if "true" restart requested
#elif defined FEM
restartInc = FEMRestartInc restartInc = FEMRestartInc
if(restartInc <= 0_pInt) then #endif
if(restartInc < 0_pInt) then
call IO_warning(warning_ID=34_pInt) call IO_warning(warning_ID=34_pInt)
restartInc = 1_pInt restartInc = 0_pInt
endif endif
restartRead = restartInc > 1_pInt restartRead = restartInc > 0_pInt ! only read in if "true" restart requested
#else #else
call IO_open_inputFile(FILEUNIT,modelName) call IO_open_inputFile(FILEUNIT,modelName)
rewind(FILEUNIT) rewind(FILEUNIT)
@ -131,13 +133,13 @@ subroutine FE_init
do do
read (FILEUNIT,'(a1024)',END=200) line read (FILEUNIT,'(a1024)',END=200) line
chunkPos = IO_stringPos(line) chunkPos = IO_stringPos(line)
if ( IO_lc(IO_stringValue(line,chunkPos,1_pInt)) == 'restart' .and. & if ( IO_lc(IO_stringValue(line,chunkPos,1_pInt)) == 'restart' &
IO_lc(IO_stringValue(line,chunkPos,2_pInt)) == 'file' .and. & .and. IO_lc(IO_stringValue(line,chunkPos,2_pInt)) == 'file' &
IO_lc(IO_stringValue(line,chunkPos,3_pInt)) == 'job' .and. & .and. IO_lc(IO_stringValue(line,chunkPos,3_pInt)) == 'job' &
IO_lc(IO_stringValue(line,chunkPos,4_pInt)) == 'id' ) & .and. IO_lc(IO_stringValue(line,chunkPos,4_pInt)) == 'id' ) &
modelName = IO_StringValue(line,chunkPos,6_pInt) modelName = IO_StringValue(line,chunkPos,6_pInt)
enddo enddo
#else #else ! QUESTION: is this meaningful for the spectral/FEM case?
call IO_open_inputFile(FILEUNIT,modelName) call IO_open_inputFile(FILEUNIT,modelName)
rewind(FILEUNIT) rewind(FILEUNIT)
do do

View File

@ -82,7 +82,7 @@ contains
!> @brief only outputs revision number !> @brief only outputs revision number
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine IO_init subroutine IO_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -1268,11 +1268,11 @@ integer(pInt) function IO_countNumericalDataLines(fileUnit)
line = IO_read(fileUnit) line = IO_read(fileUnit)
chunkPos = IO_stringPos(line) chunkPos = IO_stringPos(line)
tmp = IO_lc(IO_stringValue(line,chunkPos,1_pInt)) tmp = IO_lc(IO_stringValue(line,chunkPos,1_pInt))
if (verify(trim(tmp) ,"0123456789")/=0) then ! found keyword if (verify(trim(tmp),'0123456789') == 0) then ! numerical values
IO_countNumericalDataLines = IO_countNumericalDataLines + 1_pInt
else
line = IO_read(fileUnit, .true.) ! reset IO_read line = IO_read(fileUnit, .true.) ! reset IO_read
exit exit
else
IO_countNumericalDataLines = IO_countNumericalDataLines + 1_pInt
endif endif
enddo enddo
backspace(fileUnit) backspace(fileUnit)

View File

@ -1,8 +1,7 @@
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
write(6,*) 'Compiled with ', compiler_version() !not supported by and ifort <= 15 (and old gfortran) write(6,*) 'Compiled with ', compiler_version()
write(6,*) 'With options ', compiler_options() write(6,*) 'With options ', compiler_options()
#endif #else
#ifdef __INTEL_COMPILER
write(6,'(a,i4.4,a,i8.8)') ' Compiled with Intel fortran version ', __INTEL_COMPILER,& write(6,'(a,i4.4,a,i8.8)') ' Compiled with Intel fortran version ', __INTEL_COMPILER,&
', build date ', __INTEL_COMPILER_BUILD_DATE ', build date ', __INTEL_COMPILER_BUILD_DATE
#endif #endif

View File

@ -37,7 +37,7 @@ contains
!> @brief allocates arrays pointing to array of the various constitutive modules !> @brief allocates arrays pointing to array of the various constitutive modules
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine constitutive_init() subroutine constitutive_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -902,7 +902,8 @@ subroutine constitutive_collectDotState(Tstar_v, FeArray, FpArray, subdt, subfra
real(pReal), intent(in), dimension(6) :: & real(pReal), intent(in), dimension(6) :: &
Tstar_v !< 2nd Piola Kirchhoff stress tensor (Mandel) Tstar_v !< 2nd Piola Kirchhoff stress tensor (Mandel)
integer(pLongInt) :: & integer(pLongInt) :: &
tick, tock, & tick = 0_pLongInt, &
tock = 0_pLongInt, &
tickrate, & tickrate, &
maxticks maxticks
integer(pInt) :: & integer(pInt) :: &

View File

@ -55,14 +55,14 @@ module crystallite
crystallite_Li0, & !< intermediate velocitiy grad at start of FE inc crystallite_Li0, & !< intermediate velocitiy grad at start of FE inc
crystallite_partionedLi0,& !< intermediate velocity grad at start of homog inc crystallite_partionedLi0,& !< intermediate velocity grad at start of homog inc
crystallite_Fe, & !< current "elastic" def grad (end of converged time step) crystallite_Fe, & !< current "elastic" def grad (end of converged time step)
crystallite_P, & !< 1st Piola-Kirchhoff stress per grain crystallite_P !< 1st Piola-Kirchhoff stress per grain
crystallite_subF !< def grad to be reached at end of crystallite inc
real(pReal), dimension(:,:,:,:,:), allocatable, private :: & real(pReal), dimension(:,:,:,:,:), allocatable, private :: &
crystallite_subFe0,& !< "elastic" def grad at start of crystallite inc crystallite_subFe0,& !< "elastic" def grad at start of crystallite inc
crystallite_invFp, & !< inverse of current plastic def grad (end of converged time step) crystallite_invFp, & !< inverse of current plastic def grad (end of converged time step)
crystallite_subFp0,& !< plastic def grad at start of crystallite inc crystallite_subFp0,& !< plastic def grad at start of crystallite inc
crystallite_invFi, & !< inverse of current intermediate def grad (end of converged time step) crystallite_invFi, & !< inverse of current intermediate def grad (end of converged time step)
crystallite_subFi0,& !< intermediate def grad at start of crystallite inc crystallite_subFi0,& !< intermediate def grad at start of crystallite inc
crystallite_subF, & !< def grad to be reached at end of crystallite inc
crystallite_subF0, & !< def grad at start of crystallite inc crystallite_subF0, & !< def grad at start of crystallite inc
crystallite_subLp0,& !< plastic velocity grad at start of crystallite inc crystallite_subLp0,& !< plastic velocity grad at start of crystallite inc
crystallite_subLi0,& !< intermediate velocity grad at start of crystallite inc crystallite_subLi0,& !< intermediate velocity grad at start of crystallite inc
@ -137,7 +137,7 @@ contains
!> @brief allocates and initialize per grain variables !> @brief allocates and initialize per grain variables
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine crystallite_init subroutine crystallite_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -193,7 +193,7 @@ subroutine crystallite_init
c, & !< counter in integration point component loop c, & !< counter in integration point component loop
i, & !< counter in integration point loop i, & !< counter in integration point loop
e, & !< counter in element loop e, & !< counter in element loop
o, & !< counter in output loop o = 0_pInt, & !< counter in output loop
r, & !< counter in crystallite loop r, & !< counter in crystallite loop
cMax, & !< maximum number of integration point components cMax, & !< maximum number of integration point components
iMax, & !< maximum number of integration points iMax, & !< maximum number of integration points
@ -1250,14 +1250,16 @@ subroutine crystallite_integrateStateRK4()
use numerics, only: & use numerics, only: &
numerics_integrationMode numerics_integrationMode
use debug, only: & use debug, only: &
#ifdef DEBUG
debug_e, &
debug_i, &
debug_g, &
#endif
debug_level, & debug_level, &
debug_crystallite, & debug_crystallite, &
debug_levelBasic, & debug_levelBasic, &
debug_levelExtensive, & debug_levelExtensive, &
debug_levelSelective, & debug_levelSelective, &
debug_e, &
debug_i, &
debug_g, &
debug_StateLoopDistribution debug_StateLoopDistribution
use FEsolving, only: & use FEsolving, only: &
FEsolving_execElem, & FEsolving_execElem, &
@ -1544,14 +1546,16 @@ subroutine crystallite_integrateStateRKCK45()
use, intrinsic :: & use, intrinsic :: &
IEEE_arithmetic IEEE_arithmetic
use debug, only: & use debug, only: &
#ifdef DEBUG
debug_e, &
debug_i, &
debug_g, &
#endif
debug_level, & debug_level, &
debug_crystallite, & debug_crystallite, &
debug_levelBasic, & debug_levelBasic, &
debug_levelExtensive, & debug_levelExtensive, &
debug_levelSelective, & debug_levelSelective, &
debug_e, &
debug_i, &
debug_g, &
debug_StateLoopDistribution debug_StateLoopDistribution
use numerics, only: & use numerics, only: &
rTol_crystalliteState, & rTol_crystalliteState, &
@ -2047,14 +2051,16 @@ subroutine crystallite_integrateStateAdaptiveEuler()
use, intrinsic :: & use, intrinsic :: &
IEEE_arithmetic IEEE_arithmetic
use debug, only: & use debug, only: &
#ifdef DEBUG
debug_e, &
debug_i, &
debug_g, &
#endif
debug_level, & debug_level, &
debug_crystallite, & debug_crystallite, &
debug_levelBasic, & debug_levelBasic, &
debug_levelExtensive, & debug_levelExtensive, &
debug_levelSelective, & debug_levelSelective, &
debug_e, &
debug_i, &
debug_g, &
debug_StateLoopDistribution debug_StateLoopDistribution
use numerics, only: & use numerics, only: &
rTol_crystalliteState, & rTol_crystalliteState, &
@ -2407,14 +2413,16 @@ subroutine crystallite_integrateStateEuler()
use, intrinsic :: & use, intrinsic :: &
IEEE_arithmetic IEEE_arithmetic
use debug, only: & use debug, only: &
#ifdef DEBUG
debug_e, &
debug_i, &
debug_g, &
#endif
debug_level, & debug_level, &
debug_crystallite, & debug_crystallite, &
debug_levelBasic, & debug_levelBasic, &
debug_levelExtensive, & debug_levelExtensive, &
debug_levelSelective, & debug_levelSelective, &
debug_e, &
debug_i, &
debug_g, &
debug_StateLoopDistribution debug_StateLoopDistribution
use numerics, only: & use numerics, only: &
numerics_integrationMode, & numerics_integrationMode, &
@ -2630,9 +2638,11 @@ subroutine crystallite_integrateStateFPI()
use, intrinsic :: & use, intrinsic :: &
IEEE_arithmetic IEEE_arithmetic
use debug, only: & use debug, only: &
#ifdef DEBUG
debug_e, & debug_e, &
debug_i, & debug_i, &
debug_g, & debug_g, &
#endif
debug_level,& debug_level,&
debug_crystallite, & debug_crystallite, &
debug_levelBasic, & debug_levelBasic, &
@ -3094,14 +3104,16 @@ logical function crystallite_stateJump(ipc,ip,el)
IEEE_arithmetic IEEE_arithmetic
use prec, only: & use prec, only: &
dNeq0 dNeq0
#ifdef DEBUG
use debug, only: & use debug, only: &
debug_e, &
debug_i, &
debug_g, &
debug_level, & debug_level, &
debug_crystallite, & debug_crystallite, &
debug_levelExtensive, & debug_levelExtensive, &
debug_levelSelective, & debug_levelSelective
debug_e, & #endif
debug_i, &
debug_g
use material, only: & use material, only: &
plasticState, & plasticState, &
sourceState, & sourceState, &
@ -3231,9 +3243,11 @@ logical function crystallite_integrateStress(&
debug_levelBasic, & debug_levelBasic, &
debug_levelExtensive, & debug_levelExtensive, &
debug_levelSelective, & debug_levelSelective, &
#ifdef DEBUG
debug_e, & debug_e, &
debug_i, & debug_i, &
debug_g, & debug_g, &
#endif
debug_cumLpCalls, & debug_cumLpCalls, &
debug_cumLpTicks, & debug_cumLpTicks, &
debug_StressLoopLpDistribution, & debug_StressLoopLpDistribution, &
@ -3259,7 +3273,9 @@ logical function crystallite_integrateStress(&
math_Plain33to9, & math_Plain33to9, &
math_Plain9to33, & math_Plain9to33, &
math_Plain99to3333 math_Plain99to3333
#ifdef DEBUG
use mesh, only: mesh_element use mesh, only: mesh_element
#endif
implicit none implicit none
integer(pInt), intent(in):: el, & ! element index integer(pInt), intent(in):: el, & ! element index
@ -3323,8 +3339,8 @@ logical function crystallite_integrateStress(&
p, & p, &
jacoCounterLp, & jacoCounterLp, &
jacoCounterLi ! counters to check for Jacobian update jacoCounterLi ! counters to check for Jacobian update
integer(pLongInt) tick, & integer(pLongInt) :: tick = 0_pLongInt, &
tock, & tock = 0_pLongInt, &
tickrate, & tickrate, &
maxticks maxticks

View File

@ -41,7 +41,7 @@ contains
!> @brief allocates all neccessary fields, reads information from material configuration file !> @brief allocates all neccessary fields, reads information from material configuration file
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine damage_local_init(fileUnit) subroutine damage_local_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -16,7 +16,7 @@ contains
!> @brief allocates all neccessary fields, reads information from material configuration file !> @brief allocates all neccessary fields, reads information from material configuration file
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine damage_none_init() subroutine damage_none_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -46,7 +46,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine damage_nonlocal_init(fileUnit) subroutine damage_nonlocal_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -102,7 +102,7 @@ contains
!> @brief reads in parameters from debug.config and allocates arrays !> @brief reads in parameters from debug.config and allocates arrays
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine debug_init subroutine debug_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -71,7 +71,7 @@ contains
!> @brief module initialization !> @brief module initialization
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine homogenization_init subroutine homogenization_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -72,7 +72,7 @@ contains
!> @brief allocates all neccessary fields, reads information from material configuration file !> @brief allocates all neccessary fields, reads information from material configuration file
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine homogenization_RGC_init(fileUnit) subroutine homogenization_RGC_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -49,7 +49,7 @@ contains
!> @brief allocates all neccessary fields, reads information from material configuration file !> @brief allocates all neccessary fields, reads information from material configuration file
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine homogenization_isostrain_init(fileUnit) subroutine homogenization_isostrain_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -18,7 +18,7 @@ contains
!> @brief allocates all neccessary fields, reads information from material configuration file !> @brief allocates all neccessary fields, reads information from material configuration file
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine homogenization_none_init() subroutine homogenization_none_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -51,7 +51,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine hydrogenflux_cahnhilliard_init(fileUnit) subroutine hydrogenflux_cahnhilliard_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -16,7 +16,7 @@ contains
!> @brief allocates all neccessary fields, reads information from material configuration file !> @brief allocates all neccessary fields, reads information from material configuration file
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine hydrogenflux_isoconc_init() subroutine hydrogenflux_isoconc_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -51,7 +51,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine kinematics_cleavage_opening_init(fileUnit) subroutine kinematics_cleavage_opening_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -41,7 +41,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine kinematics_hydrogen_strain_init(fileUnit) subroutine kinematics_hydrogen_strain_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -51,7 +51,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine kinematics_slipplane_opening_init(fileUnit) subroutine kinematics_slipplane_opening_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -41,7 +41,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine kinematics_thermal_expansion_init(fileUnit) subroutine kinematics_thermal_expansion_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -174,8 +174,12 @@ pure function kinematics_thermal_expansion_initialStrain(ipc, ip, el)
offset = thermalMapping(homog)%p(ip,el) offset = thermalMapping(homog)%p(ip,el)
kinematics_thermal_expansion_initialStrain = & kinematics_thermal_expansion_initialStrain = &
(temperature(homog)%p(offset) - lattice_referenceTemperature(phase)) * & (temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**1 / 1. * &
lattice_thermalExpansion33(1:3,1:3,phase) lattice_thermalExpansion33(1:3,1:3,1,phase) + & ! constant coefficient
(temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**2 / 2. * &
lattice_thermalExpansion33(1:3,1:3,2,phase) + & ! linear coefficient
(temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**3 / 3. * &
lattice_thermalExpansion33(1:3,1:3,3,phase) ! quadratic coefficient
end function kinematics_thermal_expansion_initialStrain end function kinematics_thermal_expansion_initialStrain
@ -215,9 +219,16 @@ subroutine kinematics_thermal_expansion_LiAndItsTangent(Li, dLi_dTstar3333, ipc,
TDot = temperatureRate(homog)%p(offset) TDot = temperatureRate(homog)%p(offset)
TRef = lattice_referenceTemperature(phase) TRef = lattice_referenceTemperature(phase)
Li = TDot* & Li = TDot * ( &
lattice_thermalExpansion33(1:3,1:3,phase)/ & lattice_thermalExpansion33(1:3,1:3,1,phase)*(T - TRef)**0 & ! constant coefficient
(1.0_pReal + lattice_thermalExpansion33(1:3,1:3,phase)*(T - TRef)) + lattice_thermalExpansion33(1:3,1:3,2,phase)*(T - TRef)**1 & ! linear coefficient
+ lattice_thermalExpansion33(1:3,1:3,3,phase)*(T - TRef)**2 & ! quadratic coefficient
) / &
(1.0_pReal \
+ lattice_thermalExpansion33(1:3,1:3,1,phase)*(T - TRef)**1 / 1. &
+ lattice_thermalExpansion33(1:3,1:3,2,phase)*(T - TRef)**2 / 2. &
+ lattice_thermalExpansion33(1:3,1:3,3,phase)*(T - TRef)**3 / 3. &
)
dLi_dTstar3333 = 0.0_pReal dLi_dTstar3333 = 0.0_pReal
end subroutine kinematics_thermal_expansion_LiAndItsTangent end subroutine kinematics_thermal_expansion_LiAndItsTangent

View File

@ -41,7 +41,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine kinematics_vacancy_strain_init(fileUnit) subroutine kinematics_vacancy_strain_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -1082,9 +1082,10 @@ module lattice
lattice_nu, & lattice_nu, &
lattice_trans_mu, & lattice_trans_mu, &
lattice_trans_nu lattice_trans_nu
real(pReal), dimension(:,:,:,:), allocatable, public, protected :: & ! with higher-order parameters (e.g. temperature-dependent)
lattice_thermalExpansion33
real(pReal), dimension(:,:,:), allocatable, public, protected :: & real(pReal), dimension(:,:,:), allocatable, public, protected :: &
lattice_thermalConductivity33, & lattice_thermalConductivity33, &
lattice_thermalExpansion33, &
lattice_damageDiffusion33, & lattice_damageDiffusion33, &
lattice_vacancyfluxDiffusion33, & lattice_vacancyfluxDiffusion33, &
lattice_vacancyfluxMobility33, & lattice_vacancyfluxMobility33, &
@ -1243,7 +1244,7 @@ contains
!> @brief Module initialization !> @brief Module initialization
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine lattice_init subroutine lattice_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -1380,8 +1381,8 @@ subroutine lattice_init
allocate(lattice_C3333(3,3,3,3,Nphases), source=0.0_pReal) allocate(lattice_C3333(3,3,3,3,Nphases), source=0.0_pReal)
allocate(lattice_trans_C66(6,6,Nphases), source=0.0_pReal) allocate(lattice_trans_C66(6,6,Nphases), source=0.0_pReal)
allocate(lattice_trans_C3333(3,3,3,3,Nphases), source=0.0_pReal) allocate(lattice_trans_C3333(3,3,3,3,Nphases), source=0.0_pReal)
allocate(lattice_thermalExpansion33 (3,3,3,Nphases), source=0.0_pReal) ! constant, linear, quadratic coefficients
allocate(lattice_thermalConductivity33 (3,3,Nphases), source=0.0_pReal) allocate(lattice_thermalConductivity33 (3,3,Nphases), source=0.0_pReal)
allocate(lattice_thermalExpansion33 (3,3,Nphases), source=0.0_pReal)
allocate(lattice_damageDiffusion33 (3,3,Nphases), source=0.0_pReal) allocate(lattice_damageDiffusion33 (3,3,Nphases), source=0.0_pReal)
allocate(lattice_vacancyfluxDiffusion33 (3,3,Nphases), source=0.0_pReal) allocate(lattice_vacancyfluxDiffusion33 (3,3,Nphases), source=0.0_pReal)
allocate(lattice_vacancyfluxMobility33 (3,3,Nphases), source=0.0_pReal) allocate(lattice_vacancyfluxMobility33 (3,3,Nphases), source=0.0_pReal)
@ -1545,11 +1546,17 @@ subroutine lattice_init
case ('thermal_conductivity33') case ('thermal_conductivity33')
lattice_thermalConductivity33(3,3,section) = IO_floatValue(line,chunkPos,2_pInt) lattice_thermalConductivity33(3,3,section) = IO_floatValue(line,chunkPos,2_pInt)
case ('thermal_expansion11') case ('thermal_expansion11')
lattice_thermalExpansion33(1,1,section) = IO_floatValue(line,chunkPos,2_pInt) do i = 2_pInt, min(4,chunkPos(1)) ! read up to three parameters (constant, linear, quadratic with T)
lattice_thermalExpansion33(1,1,i-1_pInt,section) = IO_floatValue(line,chunkPos,i)
enddo
case ('thermal_expansion22') case ('thermal_expansion22')
lattice_thermalExpansion33(2,2,section) = IO_floatValue(line,chunkPos,2_pInt) do i = 2_pInt, min(4,chunkPos(1)) ! read up to three parameters (constant, linear, quadratic with T)
lattice_thermalExpansion33(2,2,i-1_pInt,section) = IO_floatValue(line,chunkPos,i)
enddo
case ('thermal_expansion33') case ('thermal_expansion33')
lattice_thermalExpansion33(3,3,section) = IO_floatValue(line,chunkPos,2_pInt) do i = 2_pInt, min(4,chunkPos(1)) ! read up to three parameters (constant, linear, quadratic with T)
lattice_thermalExpansion33(3,3,i-1_pInt,section) = IO_floatValue(line,chunkPos,i)
enddo
case ('specific_heat') case ('specific_heat')
lattice_specificHeat(section) = IO_floatValue(line,chunkPos,2_pInt) lattice_specificHeat(section) = IO_floatValue(line,chunkPos,2_pInt)
case ('vacancyformationenergy') case ('vacancyformationenergy')
@ -1756,10 +1763,12 @@ subroutine lattice_initializeStructure(myPhase,CoverA,CoverA_trans,a_fcc,a_bcc)
end select end select
end select end select
forall (i = 1_pInt:3_pInt) &
lattice_thermalExpansion33 (1:3,1:3,i,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),&
lattice_thermalExpansion33 (1:3,1:3,i,myPhase))
lattice_thermalConductivity33 (1:3,1:3,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),& lattice_thermalConductivity33 (1:3,1:3,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),&
lattice_thermalConductivity33 (1:3,1:3,myPhase)) lattice_thermalConductivity33 (1:3,1:3,myPhase))
lattice_thermalExpansion33(1:3,1:3,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),&
lattice_thermalExpansion33(1:3,1:3,myPhase))
lattice_DamageDiffusion33 (1:3,1:3,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),& lattice_DamageDiffusion33 (1:3,1:3,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),&
lattice_DamageDiffusion33 (1:3,1:3,myPhase)) lattice_DamageDiffusion33 (1:3,1:3,myPhase))
lattice_vacancyfluxDiffusion33 (1:3,1:3,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),& lattice_vacancyfluxDiffusion33 (1:3,1:3,myPhase) = lattice_symmetrize33(lattice_structure(myPhase),&

View File

@ -364,7 +364,7 @@ contains
!> material.config !> material.config
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine material_init() subroutine material_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -174,7 +174,7 @@ contains
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine math_init subroutine math_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -1743,13 +1743,12 @@ real(pReal) pure function math_EulerMisorientation(EulerA,EulerB)
implicit none implicit none
real(pReal), dimension(3), intent(in) :: EulerA,EulerB real(pReal), dimension(3), intent(in) :: EulerA,EulerB
real(pReal), dimension(3,3) :: r real(pReal) :: cosTheta
real(pReal) :: tr
r = math_mul33x33(math_EulerToR(EulerB),transpose(math_EulerToR(EulerA))) cosTheta = (math_trace33(math_mul33x33(math_EulerToR(EulerB), &
transpose(math_EulerToR(EulerA)))) - 1.0_pReal) * 0.5_pReal
tr = (math_trace33(r)-1.0_pReal)*0.4999999_pReal math_EulerMisorientation = acos(math_limit(cosTheta,-1.0_pReal,1.0_pReal))
math_EulerMisorientation = abs(0.5_pReal*PI-asin(tr))
end function math_EulerMisorientation end function math_EulerMisorientation
@ -1782,15 +1781,12 @@ function math_sampleGaussOri(center,noise)
real(pReal), dimension(3), intent(in) :: center real(pReal), dimension(3), intent(in) :: center
real(pReal) :: cosScatter,scatter real(pReal) :: cosScatter,scatter
real(pReal), dimension(3) :: math_sampleGaussOri, disturb real(pReal), dimension(3) :: math_sampleGaussOri, disturb
real(pReal), dimension(3), parameter :: ORIGIN = [0.0_pReal,0.0_pReal,0.0_pReal] real(pReal), dimension(3), parameter :: ORIGIN = 0.0_pReal
real(pReal), dimension(5) :: rnd real(pReal), dimension(5) :: rnd
integer(pInt) :: i
if (abs(noise) < tol_math_check) then noScatter: if (abs(noise) < tol_math_check) then
math_sampleGaussOri = center math_sampleGaussOri = center
return else noScatter
endif
! Helming uses different distribution with Bessel functions ! Helming uses different distribution with Bessel functions
! therefore the gauss scatter width has to be scaled differently ! therefore the gauss scatter width has to be scaled differently
scatter = 0.95_pReal * noise scatter = 0.95_pReal * noise
@ -1798,14 +1794,15 @@ function math_sampleGaussOri(center,noise)
do do
call halton(5_pInt,rnd) call halton(5_pInt,rnd)
forall (i=1_pInt:3_pInt) rnd(i) = 2.0_pReal*rnd(i)-1.0_pReal ! expand 1:3 to range [-1,+1] rnd(1:3) = 2.0_pReal*rnd(1:3)-1.0_pReal ! expand 1:3 to range [-1,+1]
disturb = [ scatter * rnd(1), & ! phi1 disturb = [ scatter * rnd(1), & ! phi1
sign(1.0_pReal,rnd(2))*acos(cosScatter+(1.0_pReal-cosScatter)*rnd(4)), & ! Phi sign(1.0_pReal,rnd(2))*acos(cosScatter+(1.0_pReal-cosScatter)*rnd(4)), & ! Phi
scatter * rnd(2)] ! phi2 scatter * rnd(3)] ! phi2
if (rnd(5) <= exp(-1.0_pReal*(math_EulerMisorientation(ORIGIN,disturb)/scatter)**2_pReal)) exit if (rnd(5) <= exp(-1.0_pReal*(math_EulerMisorientation(ORIGIN,disturb)/scatter)**2_pReal)) exit
enddo enddo
math_sampleGaussOri = math_RtoEuler(math_mul33x33(math_EulerToR(disturb),math_EulerToR(center))) math_sampleGaussOri = math_RtoEuler(math_mul33x33(math_EulerToR(disturb),math_EulerToR(center)))
endif noScatter
end function math_sampleGaussOri end function math_sampleGaussOri
@ -2759,8 +2756,7 @@ pure function math_rotate_forward33(tensor,rot_tensor)
real(pReal), dimension(3,3) :: math_rotate_forward33 real(pReal), dimension(3,3) :: math_rotate_forward33
real(pReal), dimension(3,3), intent(in) :: tensor, rot_tensor real(pReal), dimension(3,3), intent(in) :: tensor, rot_tensor
math_rotate_forward33 = math_mul33x33(rot_tensor,& math_rotate_forward33 = math_mul33x33(rot_tensor,math_mul33x33(tensor,transpose(rot_tensor)))
math_mul33x33(tensor,math_transpose33(rot_tensor)))
end function math_rotate_forward33 end function math_rotate_forward33
@ -2774,8 +2770,7 @@ pure function math_rotate_backward33(tensor,rot_tensor)
real(pReal), dimension(3,3) :: math_rotate_backward33 real(pReal), dimension(3,3) :: math_rotate_backward33
real(pReal), dimension(3,3), intent(in) :: tensor, rot_tensor real(pReal), dimension(3,3), intent(in) :: tensor, rot_tensor
math_rotate_backward33 = math_mul33x33(math_transpose33(rot_tensor),& math_rotate_backward33 = math_mul33x33(transpose(rot_tensor),math_mul33x33(tensor,rot_tensor))
math_mul33x33(tensor,rot_tensor))
end function math_rotate_backward33 end function math_rotate_backward33
@ -2796,8 +2791,8 @@ pure function math_rotate_forward3333(tensor,rot_tensor)
do i = 1_pInt,3_pInt; do j = 1_pInt,3_pInt; do k = 1_pInt,3_pInt; do l = 1_pInt,3_pInt do i = 1_pInt,3_pInt; do j = 1_pInt,3_pInt; do k = 1_pInt,3_pInt; do l = 1_pInt,3_pInt
do m = 1_pInt,3_pInt; do n = 1_pInt,3_pInt; do o = 1_pInt,3_pInt; do p = 1_pInt,3_pInt do m = 1_pInt,3_pInt; do n = 1_pInt,3_pInt; do o = 1_pInt,3_pInt; do p = 1_pInt,3_pInt
math_rotate_forward3333(i,j,k,l) = math_rotate_forward3333(i,j,k,l) & math_rotate_forward3333(i,j,k,l) = math_rotate_forward3333(i,j,k,l) &
+ rot_tensor(m,i) * rot_tensor(n,j) & + rot_tensor(i,m) * rot_tensor(j,n) &
* rot_tensor(o,k) * rot_tensor(p,l) * tensor(m,n,o,p) * rot_tensor(k,o) * rot_tensor(l,p) * tensor(m,n,o,p)
enddo; enddo; enddo; enddo; enddo; enddo; enddo; enddo enddo; enddo; enddo; enddo; enddo; enddo; enddo; enddo
end function math_rotate_forward3333 end function math_rotate_forward3333

View File

@ -379,30 +379,30 @@ module mesh
],pInt) ],pInt)
integer(pInt), dimension(FE_Nelemtypes), parameter, private :: MESH_VTKELEMTYPE = & ! integer(pInt), dimension(FE_Nelemtypes), parameter, private :: MESH_VTKELEMTYPE = &
int([ & ! int([ &
5, & ! element 6 (2D 3node 1ip) ! 5, & ! element 6 (2D 3node 1ip)
22, & ! element 125 (2D 6node 3ip) ! 22, & ! element 125 (2D 6node 3ip)
9, & ! element 11 (2D 4node 4ip) ! 9, & ! element 11 (2D 4node 4ip)
23, & ! element 27 (2D 8node 9ip) ! 23, & ! element 27 (2D 8node 9ip)
23, & ! element 54 (2D 8node 4ip) ! 23, & ! element 54 (2D 8node 4ip)
10, & ! element 134 (3D 4node 1ip) ! 10, & ! element 134 (3D 4node 1ip)
10, & ! element 157 (3D 5node 4ip) ! 10, & ! element 157 (3D 5node 4ip)
24, & ! element 127 (3D 10node 4ip) ! 24, & ! element 127 (3D 10node 4ip)
13, & ! element 136 (3D 6node 6ip) ! 13, & ! element 136 (3D 6node 6ip)
12, & ! element 117 (3D 8node 1ip) ! 12, & ! element 117 (3D 8node 1ip)
12, & ! element 7 (3D 8node 8ip) ! 12, & ! element 7 (3D 8node 8ip)
25, & ! element 57 (3D 20node 8ip) ! 25, & ! element 57 (3D 20node 8ip)
25 & ! element 21 (3D 20node 27ip) ! 25 & ! element 21 (3D 20node 27ip)
],pInt) ! ],pInt)
!
integer(pInt), dimension(FE_Ncelltypes), parameter, private :: MESH_VTKCELLTYPE = & ! integer(pInt), dimension(FE_Ncelltypes), parameter, private :: MESH_VTKCELLTYPE = &
int([ & ! int([ &
5, & ! (2D 3node) ! 5, & ! (2D 3node)
9, & ! (2D 4node) ! 9, & ! (2D 4node)
10, & ! (3D 4node) ! 10, & ! (3D 4node)
12 & ! (3D 8node) ! 12 & ! (3D 8node)
],pInt) ! ],pInt)
public :: & public :: &
@ -477,7 +477,7 @@ contains
!! Order and routines strongly depend on type of solver !! Order and routines strongly depend on type of solver
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine mesh_init(ip,el) subroutine mesh_init(ip,el)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -2848,16 +2848,18 @@ use IO, only: &
implicit none implicit none
integer(pInt), intent(in) :: fileUnit integer(pInt), intent(in) :: fileUnit
#ifndef Spectral #ifdef Spectral
mesh_periodicSurface = .true.
end subroutine mesh_get_damaskOptions
#else
integer(pInt), allocatable, dimension(:) :: chunkPos integer(pInt), allocatable, dimension(:) :: chunkPos
integer(pInt) chunk, Nchunks integer(pInt) chunk, Nchunks
character(len=300) :: line, damaskOption, v character(len=300) :: line, damaskOption, v
character(len=300) :: keyword character(len=300) :: keyword
#endif
#ifdef Spectral
mesh_periodicSurface = .true.
#else
mesh_periodicSurface = .false. mesh_periodicSurface = .false.
#ifdef Marc4DAMASK #ifdef Marc4DAMASK
keyword = '$damask' keyword = '$damask'
@ -2886,9 +2888,9 @@ use IO, only: &
enddo enddo
610 FORMAT(A300) 610 FORMAT(A300)
#endif
620 end subroutine mesh_get_damaskOptions 620 end subroutine mesh_get_damaskOptions
#endif
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------

View File

@ -198,7 +198,7 @@ contains
! a sanity check ! a sanity check
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine numerics_init subroutine numerics_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -119,7 +119,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine plastic_disloUCLA_init(fileUnit) subroutine plastic_disloUCLA_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -198,7 +198,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine plastic_dislotwin_init(fileUnit) subroutine plastic_dislotwin_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -90,7 +90,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine plastic_isotropic_init(fileUnit) subroutine plastic_isotropic_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -26,7 +26,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine plastic_none_init subroutine plastic_none_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -123,7 +123,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine plastic_phenopowerlaw_init(fileUnit) subroutine plastic_phenopowerlaw_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -16,7 +16,7 @@ contains
!> @brief allocates all neccessary fields, reads information from material configuration file !> @brief allocates all neccessary fields, reads information from material configuration file
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine porosity_none_init() subroutine porosity_none_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -48,7 +48,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine porosity_phasefield_init(fileUnit) subroutine porosity_phasefield_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -111,7 +111,7 @@ contains
!> @brief reporting precision !> @brief reporting precision
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine prec_init subroutine prec_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -63,7 +63,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_damage_anisoBrittle_init(fileUnit) subroutine source_damage_anisoBrittle_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -67,7 +67,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_damage_anisoDuctile_init(fileUnit) subroutine source_damage_anisoDuctile_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -53,7 +53,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_damage_isoBrittle_init(fileUnit) subroutine source_damage_isoBrittle_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -53,7 +53,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_damage_isoDuctile_init(fileUnit) subroutine source_damage_isoDuctile_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -39,7 +39,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_thermal_dissipation_init(fileUnit) subroutine source_thermal_dissipation_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -45,7 +45,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_thermal_externalheat_init(fileUnit) subroutine source_thermal_externalheat_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -41,7 +41,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_vacancy_irradiation_init(fileUnit) subroutine source_vacancy_irradiation_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -39,7 +39,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_vacancy_phenoplasticity_init(fileUnit) subroutine source_vacancy_phenoplasticity_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -41,7 +41,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_vacancy_thermalfluc_init(fileUnit) subroutine source_vacancy_thermalfluc_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -60,7 +60,7 @@ contains
!> @brief allocates all neccessary fields and fills them with data, potentially from restart info !> @brief allocates all neccessary fields and fills them with data, potentially from restart info
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine spectral_damage_init() subroutine spectral_damage_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -15,7 +15,7 @@ module DAMASK_interface
private private
#include <petsc/finclude/petscsys.h> #include <petsc/finclude/petscsys.h>
logical, public, protected :: appendToOutFile = .false. !< Append to existing spectralOut file (in case of restart, not in case of regridding) logical, public, protected :: appendToOutFile = .false. !< Append to existing spectralOut file (in case of restart, not in case of regridding)
integer(pInt), public, protected :: spectralRestartInc = 1_pInt !< Increment at which calculation starts integer(pInt), public, protected :: spectralRestartInc = 0_pInt !< Increment at which calculation starts
character(len=1024), public, protected :: & character(len=1024), public, protected :: &
geometryFile = '', & !< parameter given for geometry file geometryFile = '', & !< parameter given for geometry file
loadCaseFile = '' !< parameter given for load case file loadCaseFile = '' !< parameter given for load case file
@ -83,7 +83,9 @@ subroutine DAMASK_interface_init()
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! PETSc Init ! PETSc Init
#ifdef _OPENMP #ifdef _OPENMP
call MPI_Init_Thread(MPI_THREAD_FUNNELED,threadLevel,ierr);CHKERRQ(ierr) ! in case of OpenMP, don't rely on PETScInitialize doing MPI init ! If openMP is enabled, check if the MPI libary supports it and initialize accordingly.
! Otherwise, the first call to PETSc will do the initialization.
call MPI_Init_Thread(MPI_THREAD_FUNNELED,threadLevel,ierr);CHKERRQ(ierr)
if (threadLevel<MPI_THREAD_FUNNELED) then if (threadLevel<MPI_THREAD_FUNNELED) then
write(6,'(a)') ' MPI library does not support OpenMP' write(6,'(a)') ' MPI library does not support OpenMP'
call quit(1_pInt) call quit(1_pInt)
@ -165,12 +167,12 @@ subroutine DAMASK_interface_init()
write(6,'(a)') ' For further configuration place "numerics.config"' write(6,'(a)') ' For further configuration place "numerics.config"'
write(6,'(a)')' and "debug.config" in that directory.' write(6,'(a)')' and "debug.config" in that directory.'
write(6,'(/,a)')' --restart XX' write(6,'(/,a)')' --restart XX'
write(6,'(a)') ' Reads in total increment No. XX-1 and continues to' write(6,'(a)') ' Reads in increment XX and continues with calculating'
write(6,'(a)') ' calculate total increment No. XX.' write(6,'(a)') ' increment XX+1 based on this.'
write(6,'(a)') ' Appends to existing results file' write(6,'(a)') ' Appends to existing results file'
write(6,'(a)') ' "NameOfGeom_NameOfLoadFile.spectralOut".' write(6,'(a)') ' "NameOfGeom_NameOfLoadFile.spectralOut".'
write(6,'(a)') ' Works only if the restart information for total increment' write(6,'(a)') ' Works only if the restart information for increment XX'
write(6,'(a)') ' No. XX-1 is available in the working directory.' write(6,'(a)') ' is available in the working directory.'
write(6,'(/,a)')' -----------------------------------------------------------------------' write(6,'(/,a)')' -----------------------------------------------------------------------'
write(6,'(a)') ' Help:' write(6,'(a)') ' Help:'
write(6,'(/,a)')' --help' write(6,'(/,a)')' --help'
@ -210,8 +212,8 @@ subroutine DAMASK_interface_init()
write(6,'(a,a)') ' Geometry file: ', trim(geometryFile) write(6,'(a,a)') ' Geometry file: ', trim(geometryFile)
write(6,'(a,a)') ' Loadcase file: ', trim(loadCaseFile) write(6,'(a,a)') ' Loadcase file: ', trim(loadCaseFile)
write(6,'(a,a)') ' Solver job name: ', trim(getSolverJobName()) write(6,'(a,a)') ' Solver job name: ', trim(getSolverJobName())
if (SpectralRestartInc > 1_pInt) & if (SpectralRestartInc > 0_pInt) &
write(6,'(a,i6.6)') ' Restart at increment: ', spectralRestartInc write(6,'(a,i6.6)') ' Restart from increment: ', spectralRestartInc
write(6,'(a,l1,/)') ' Append to result file: ', appendToOutFile write(6,'(a,l1,/)') ' Append to result file: ', appendToOutFile
end subroutine DAMASK_interface_init end subroutine DAMASK_interface_init

View File

@ -43,17 +43,20 @@ module spectral_mech_AL
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! stress, stiffness and compliance average etc. ! stress, stiffness and compliance average etc.
real(pReal), private, dimension(3,3) :: & real(pReal), private, dimension(3,3) :: &
F_aimDot, & !< assumed rate of average deformation gradient F_aimDot = 0.0_pReal, & !< assumed rate of average deformation gradient
F_aim = math_I3, & !< current prescribed deformation gradient F_aim = math_I3, & !< current prescribed deformation gradient
F_aim_lastInc = math_I3, & !< previous average deformation gradient F_aim_lastInc = math_I3, & !< previous average deformation gradient
F_av = 0.0_pReal, & !< average incompatible def grad field F_av = 0.0_pReal, & !< average incompatible def grad field
P_av = 0.0_pReal, & !< average 1st Piola--Kirchhoff stress P_av = 0.0_pReal, & !< average 1st Piola--Kirchhoff stress
P_avLastEval = 0.0_pReal !< average 1st Piola--Kirchhoff stress last call of CPFEM_general P_avLastEval = 0.0_pReal !< average 1st Piola--Kirchhoff stress last call of CPFEM_general
character(len=1024), private :: incInfo !< time and increment information character(len=1024), private :: incInfo !< time and increment information
real(pReal), private, dimension(3,3,3,3) :: & real(pReal), private, dimension(3,3,3,3) :: &
C_volAvg = 0.0_pReal, & !< current volume average stiffness C_volAvg = 0.0_pReal, & !< current volume average stiffness
C_volAvgLastInc = 0.0_pReal, & !< previous volume average stiffness C_volAvgLastInc = 0.0_pReal, & !< previous volume average stiffness
C_minMaxAvg = 0.0_pReal, & !< current (min+max)/2 stiffness C_minMaxAvg = 0.0_pReal, & !< current (min+max)/2 stiffness
C_minMaxAvgLastInc = 0.0_pReal, & !< previous (min+max)/2 stiffness
S = 0.0_pReal, & !< current compliance (filled up with zeros) S = 0.0_pReal, & !< current compliance (filled up with zeros)
C_scale = 0.0_pReal, & C_scale = 0.0_pReal, &
S_scale = 0.0_pReal S_scale = 0.0_pReal
@ -62,7 +65,7 @@ module spectral_mech_AL
err_BC, & !< deviation from stress BC err_BC, & !< deviation from stress BC
err_curl, & !< RMS of curl of F err_curl, & !< RMS of curl of F
err_div !< RMS of div of P err_div !< RMS of div of P
logical, private :: ForwardData
integer(pInt), private :: & integer(pInt), private :: &
totalIter = 0_pInt !< total iteration in current increment totalIter = 0_pInt !< total iteration in current increment
@ -80,11 +83,11 @@ module spectral_mech_AL
contains contains
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief allocates all neccessary fields and fills them with data, potentially from restart info !> @brief allocates all necessary fields and fills them with data, potentially from restart info
!> @todo use sourced allocation, e.g. allocate(Fdot,source = F_lastInc) !> @todo use sourced allocation, e.g. allocate(Fdot,source = F_lastInc)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine AL_init subroutine AL_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -102,12 +105,15 @@ subroutine AL_init
use numerics, only: & use numerics, only: &
worldrank, & worldrank, &
worldsize worldsize
use homogenization, only: &
materialpoint_F0
use DAMASK_interface, only: & use DAMASK_interface, only: &
getSolverJobName getSolverJobName
use spectral_utilities, only: & use spectral_utilities, only: &
Utilities_constitutiveResponse, & Utilities_constitutiveResponse, &
Utilities_updateGamma, & Utilities_updateGamma, &
Utilities_updateIPcoords Utilities_updateIPcoords, &
wgt
use mesh, only: & use mesh, only: &
grid, & grid, &
grid3 grid3
@ -120,7 +126,11 @@ subroutine AL_init
temp33_Real = 0.0_pReal temp33_Real = 0.0_pReal
PetscErrorCode :: ierr PetscErrorCode :: ierr
PetscScalar, pointer, dimension(:,:,:,:) :: xx_psc, F, F_lambda PetscScalar, pointer, dimension(:,:,:,:) :: &
FandF_lambda, & ! overall pointer to solution data
F, & ! specific (sub)pointer
F_lambda ! specific (sub)pointer
integer(pInt), dimension(:), allocatable :: localK integer(pInt), dimension(:), allocatable :: localK
integer(pInt) :: proc integer(pInt) :: proc
character(len=1024) :: rankStr character(len=1024) :: rankStr
@ -164,79 +174,70 @@ subroutine AL_init
grid(1),grid(2),localK, & ! local grid grid(1),grid(2),localK, & ! local grid
da,ierr) ! handle, error da,ierr) ! handle, error
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetDM(snes,da,ierr); CHKERRQ(ierr) call SNESSetDM(snes,da,ierr); CHKERRQ(ierr) ! connect snes to da
call DMCreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr) call DMCreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr) ! global solution vector (grid x 9, i.e. every def grad tensor)
call DMDASNESSetFunctionLocal(da,INSERT_VALUES,AL_formResidual,PETSC_NULL_OBJECT,ierr) call DMDASNESSetFunctionLocal(da,INSERT_VALUES,AL_formResidual,PETSC_NULL_OBJECT,ierr) ! residual vector of same shape as solution vector
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetConvergenceTest(snes,AL_converged,PETSC_NULL_OBJECT,PETSC_NULL_FUNCTION,ierr) call SNESSetConvergenceTest(snes,AL_converged,PETSC_NULL_OBJECT,PETSC_NULL_FUNCTION,ierr) ! specify custom convergence check function "_converged"
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetFromOptions(snes,ierr); CHKERRQ(ierr) call SNESSetFromOptions(snes,ierr); CHKERRQ(ierr) ! pull it all together with additional cli arguments
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! init fields ! init fields
call DMDAVecGetArrayF90(da,solution_vec,xx_psc,ierr); CHKERRQ(ierr) ! places pointer xx_psc on PETSc data call DMDAVecGetArrayF90(da,solution_vec,FandF_lambda,ierr); CHKERRQ(ierr) ! places pointer on PETSc data
F => xx_psc(0:8,:,:,:) F => FandF_lambda( 0: 8,:,:,:)
F_lambda => xx_psc(9:17,:,:,:) F_lambda => FandF_lambda( 9:17,:,:,:)
restart: if (restartInc > 1_pInt) then
if (iand(debug_level(debug_spectral),debug_spectralRestart)/= 0 .and. worldrank == 0_pInt) & restart: if (restartInc > 0_pInt) then
write(6,'(/,a,'//IO_intOut(restartInc-1_pInt)//',a)') & if (iand(debug_level(debug_spectral),debug_spectralRestart) /= 0) then
'reading values of increment ', restartInc - 1_pInt, ' from file' write(6,'(/,a,'//IO_intOut(restartInc)//',a)') &
'reading values of increment ', restartInc, ' from file'
flush(6) flush(6)
endif
write(rankStr,'(a1,i0)')'_',worldrank write(rankStr,'(a1,i0)')'_',worldrank
call IO_read_realFile(777,'F'//trim(rankStr),trim(getSolverJobName()),size(F)) call IO_read_realFile(777,'F'//trim(rankStr),trim(getSolverJobName()),size(F))
read (777,rec=1) F read (777,rec=1) F; close (777)
close (777)
call IO_read_realFile(777,'F_lastInc'//trim(rankStr),trim(getSolverJobName()),size(F_lastInc)) call IO_read_realFile(777,'F_lastInc'//trim(rankStr),trim(getSolverJobName()),size(F_lastInc))
read (777,rec=1) F_lastInc read (777,rec=1) F_lastInc; close (777)
close (777)
call IO_read_realFile(777,'F_lambda'//trim(rankStr),trim(getSolverJobName()),size(F_lambda)) call IO_read_realFile(777,'F_lambda'//trim(rankStr),trim(getSolverJobName()),size(F_lambda))
read (777,rec=1) F_lambda read (777,rec=1) F_lambda; close (777)
close (777) call IO_read_realFile(777,'F_lambda_lastInc'//trim(rankStr),trim(getSolverJobName()),size(F_lambda_lastInc))
call IO_read_realFile(777,'F_lambda_lastInc'//trim(rankStr),& read (777,rec=1) F_lambda_lastInc; close (777)
trim(getSolverJobName()),size(F_lambda_lastInc)) call IO_read_realFile(777,'F_aimDot',trim(getSolverJobName()),size(F_aimDot))
read (777,rec=1) F_lambda_lastInc read (777,rec=1) F_aimDot; close (777)
close (777) F_aim = reshape(sum(sum(sum(F,dim=4),dim=3),dim=2) * wgt, [3,3]) ! average of F
call IO_read_realFile(777,'F_aim', trim(getSolverJobName()),size(F_aim)) F_aim_lastInc = sum(sum(sum(F_lastInc,dim=5),dim=4),dim=3) * wgt ! average of F_lastInc
read (777,rec=1) F_aim elseif (restartInc == 0_pInt) then restart
close (777)
call IO_read_realFile(777,'F_aim_lastInc', trim(getSolverJobName()),size(F_aim_lastInc))
read (777,rec=1) F_aim_lastInc
close (777)
call IO_read_realFile(777,'F_aimDot',trim(getSolverJobName()),size(f_aimDot))
read (777,rec=1) f_aimDot
close (777)
elseif (restartInc == 1_pInt) then restart
F_lastInc = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid3) ! initialize to identity F_lastInc = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid3) ! initialize to identity
F = reshape(F_lastInc,[9,grid(1),grid(2),grid3]) F = reshape(F_lastInc,[9,grid(1),grid(2),grid3])
F_lambda = F F_lambda = F
F_lambda_lastInc = F_lastInc F_lambda_lastInc = F_lastInc
endif restart endif restart
materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent
call Utilities_updateIPcoords(reshape(F,shape(F_lastInc))) call Utilities_updateIPcoords(reshape(F,shape(F_lastInc)))
call Utilities_constitutiveResponse(P,temp33_Real,C_volAvg,C_minMaxAvg, & call Utilities_constitutiveResponse(P,temp33_Real,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2
reshape(F,shape(F_lastInc)), 0.0_pReal, math_I3) reshape(F,shape(F_lastInc)), & ! target F
0.0_pReal, & ! time increment
math_I3) ! no rotation of boundary condition
nullify(F) nullify(F)
nullify(F_lambda) nullify(F_lambda)
call DMDAVecRestoreArrayF90(da,solution_vec,xx_psc,ierr); CHKERRQ(ierr) ! write data back to PETSc call DMDAVecRestoreArrayF90(da,solution_vec,FandF_lambda,ierr); CHKERRQ(ierr) ! write data back to PETSc
restartRead: if (restartInc > 1_pInt) then restartRead: if (restartInc > 0_pInt) then
if (iand(debug_level(debug_spectral),debug_spectralRestart)/= 0 .and. worldrank == 0_pInt) & if (iand(debug_level(debug_spectral),debug_spectralRestart)/= 0 .and. worldrank == 0_pInt) &
write(6,'(/,a,'//IO_intOut(restartInc-1_pInt)//',a)') & write(6,'(/,a,'//IO_intOut(restartInc)//',a)') &
'reading more values of increment', restartInc - 1_pInt, 'from file' 'reading more values of increment ', restartInc, ' from file'
flush(6) flush(6)
call IO_read_realFile(777,'C_volAvg',trim(getSolverJobName()),size(C_volAvg)) call IO_read_realFile(777,'C_volAvg',trim(getSolverJobName()),size(C_volAvg))
read (777,rec=1) C_volAvg read (777,rec=1) C_volAvg; close (777)
close (777)
call IO_read_realFile(777,'C_volAvgLastInc',trim(getSolverJobName()),size(C_volAvgLastInc)) call IO_read_realFile(777,'C_volAvgLastInc',trim(getSolverJobName()),size(C_volAvgLastInc))
read (777,rec=1) C_volAvgLastInc read (777,rec=1) C_volAvgLastInc; close (777)
close (777)
call IO_read_realFile(777,'C_ref',trim(getSolverJobName()),size(C_minMaxAvg)) call IO_read_realFile(777,'C_ref',trim(getSolverJobName()),size(C_minMaxAvg))
read (777,rec=1) C_minMaxAvg read (777,rec=1) C_minMaxAvg; close (777)
close (777)
endif restartRead endif restartRead
call Utilities_updateGamma(C_minMaxAvg,.True.) call Utilities_updateGamma(C_minMaxAvg,.true.)
C_scale = C_minMaxAvg C_scale = C_minMaxAvg
S_scale = math_invSym3333(C_minMaxAvg) S_scale = math_invSym3333(C_minMaxAvg)
@ -246,8 +247,7 @@ end subroutine AL_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief solution for the AL scheme with internal iterations !> @brief solution for the AL scheme with internal iterations
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
type(tSolutionState) function & type(tSolutionState) function AL_solution(incInfoIn,timeinc,timeinc_old,stress_BC,rotation_BC)
AL_solution(incInfoIn,timeinc,timeinc_old,stress_BC,rotation_BC)
use IO, only: & use IO, only: &
IO_error IO_error
use numerics, only: & use numerics, only: &
@ -266,13 +266,13 @@ type(tSolutionState) function &
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! input data for solution ! input data for solution
real(pReal), intent(in) :: &
timeinc, & !< increment in time for current solution
timeinc_old !< increment in time of last increment
type(tBoundaryCondition), intent(in) :: &
stress_BC
character(len=*), intent(in) :: & character(len=*), intent(in) :: &
incInfoIn incInfoIn
real(pReal), intent(in) :: &
timeinc, & !< increment time for current solution
timeinc_old !< increment time of last successful increment
type(tBoundaryCondition), intent(in) :: &
stress_BC
real(pReal), dimension(3,3), intent(in) :: rotation_BC real(pReal), dimension(3,3), intent(in) :: rotation_BC
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -305,18 +305,17 @@ type(tSolutionState) function &
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! solve BVP ! solve BVP
call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr) call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr); CHKERRQ(ierr)
CHKERRQ(ierr)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! check convergence ! check convergence
call SNESGetConvergedReason(snes,reason,ierr) call SNESGetConvergedReason(snes,reason,ierr); CHKERRQ(ierr)
CHKERRQ(ierr)
AL_solution%converged = reason > 0
AL_solution%iterationsNeeded = totalIter
AL_solution%termIll = terminallyIll AL_solution%termIll = terminallyIll
terminallyIll = .false. terminallyIll = .false.
if (reason == -4) call IO_error(893_pInt) if (reason == -4) call IO_error(893_pInt) ! MPI error
if (reason < 1) AL_solution%converged = .false.
AL_solution%iterationsNeeded = totalIter
end function AL_solution end function AL_solution
@ -331,8 +330,8 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
polarAlpha, & polarAlpha, &
polarBeta polarBeta
use mesh, only: & use mesh, only: &
grid3, & grid, &
grid grid3
use IO, only: & use IO, only: &
IO_intOut IO_intOut
use math, only: & use math, only: &
@ -341,6 +340,10 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
math_mul3333xx33, & math_mul3333xx33, &
math_invSym3333, & math_invSym3333, &
math_mul33x33 math_mul33x33
use debug, only: &
debug_level, &
debug_spectral, &
debug_spectralRotation
use spectral_utilities, only: & use spectral_utilities, only: &
wgt, & wgt, &
tensorField_real, & tensorField_real, &
@ -350,25 +353,17 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
Utilities_constitutiveResponse, & Utilities_constitutiveResponse, &
Utilities_divergenceRMS, & Utilities_divergenceRMS, &
Utilities_curlRMS Utilities_curlRMS
use debug, only: &
debug_level, &
debug_spectral, &
debug_spectralRotation
use homogenization, only: & use homogenization, only: &
materialpoint_dPdF materialpoint_dPdF
use FEsolving, only: & use FEsolving, only: &
terminallyIll terminallyIll
implicit none implicit none
!-------------------------------------------------------------------------------------------------- DMDALocalInfo, dimension(DMDA_LOCAL_INFO_SIZE) :: in
! strange syntax in the next line because otherwise macros expand beyond 132 character limit
DMDALocalInfo, dimension(&
DMDA_LOCAL_INFO_SIZE) :: &
in
PetscScalar, & PetscScalar, &
target, dimension(3,3,2, XG_RANGE,YG_RANGE,ZG_RANGE), intent(in) :: x_scal target, dimension(3,3,2, XG_RANGE,YG_RANGE,ZG_RANGE), intent(in) :: x_scal !< what is this?
PetscScalar, & PetscScalar, &
target, dimension(3,3,2, X_RANGE, Y_RANGE, Z_RANGE), intent(out) :: f_scal target, dimension(3,3,2, X_RANGE,Y_RANGE,Z_RANGE), intent(out) :: f_scal !< what is this?
PetscScalar, pointer, dimension(:,:,:,:,:) :: & PetscScalar, pointer, dimension(:,:,:,:,:) :: &
F, & F, &
F_lambda, & F_lambda, &
@ -395,24 +390,24 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
residual_F_lambda => f_scal(1:3,1:3,2,& residual_F_lambda => f_scal(1:3,1:3,2,&
X_RANGE, Y_RANGE, Z_RANGE) X_RANGE, Y_RANGE, Z_RANGE)
call SNESGetNumberFunctionEvals(snes,nfuncs,ierr); CHKERRQ(ierr)
call SNESGetIterationNumber(snes,PETScIter,ierr); CHKERRQ(ierr)
F_av = sum(sum(sum(F,dim=5),dim=4),dim=3) * wgt F_av = sum(sum(sum(F,dim=5),dim=4),dim=3) * wgt
call MPI_Allreduce(MPI_IN_PLACE,F_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,F_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
call SNESGetNumberFunctionEvals(snes,nfuncs,ierr); CHKERRQ(ierr)
call SNESGetIterationNumber(snes,PETScIter,ierr); CHKERRQ(ierr)
if (nfuncs == 0 .and. PETScIter == 0) totalIter = -1_pInt ! new increment if (nfuncs == 0 .and. PETScIter == 0) totalIter = -1_pInt ! new increment
newIteration: if(totalIter <= PETScIter) then
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! report begin of new iteration ! begin of new iteration
newIteration: if (totalIter <= PETScIter) then
totalIter = totalIter + 1_pInt totalIter = totalIter + 1_pInt
write(6,'(1x,a,3(a,'//IO_intOut(itmax)//'))') trim(incInfo), & write(6,'(1x,a,3(a,'//IO_intOut(itmax)//'))') &
' @ Iteration ', itmin, '≤',totalIter, '≤', itmax trim(incInfo), ' @ Iteration ', itmin, '≤',totalIter, '≤', itmax
if (iand(debug_level(debug_spectral),debug_spectralRotation) /= 0) & if (iand(debug_level(debug_spectral),debug_spectralRotation) /= 0) &
write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') ' deformation gradient aim (lab) =', & write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') &
math_transpose33(math_rotate_backward33(F_aim,params%rotation_BC)) ' deformation gradient aim (lab) =', math_transpose33(math_rotate_backward33(F_aim,params%rotation_BC))
write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') ' deformation gradient aim =', & write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') &
math_transpose33(F_aim) ' deformation gradient aim =', math_transpose33(F_aim)
flush(6) flush(6)
endif newIteration endif newIteration
@ -424,7 +419,6 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
polarBeta*math_mul3333xx33(C_scale,F(1:3,1:3,i,j,k) - math_I3) -& polarBeta*math_mul3333xx33(C_scale,F(1:3,1:3,i,j,k) - math_I3) -&
polarAlpha*math_mul33x33(F(1:3,1:3,i,j,k), & polarAlpha*math_mul33x33(F(1:3,1:3,i,j,k), &
math_mul3333xx33(C_scale,F_lambda(1:3,1:3,i,j,k) - math_I3)) math_mul3333xx33(C_scale,F_lambda(1:3,1:3,i,j,k) - math_I3))
enddo; enddo; enddo enddo; enddo; enddo
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -434,8 +428,8 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
call utilities_FFTtensorBackward() call utilities_FFTtensorBackward()
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! constructing residual ! constructing F_lambda residual
residual_F_lambda = polarBeta*F - tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) residual_F_lambda = polarBeta*F - tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) !< eq (16) in doi: 10.1016/j.ijplas.2014.02.006
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! evaluate constitutive response ! evaluate constitutive response
@ -444,15 +438,13 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
call Utilities_constitutiveResponse(residual_F,P_av,C_volAvg,C_minMaxAvg, & call Utilities_constitutiveResponse(residual_F,P_av,C_volAvg,C_minMaxAvg, &
F - residual_F_lambda/polarBeta,params%timeinc, params%rotation_BC) F - residual_F_lambda/polarBeta,params%timeinc, params%rotation_BC)
call MPI_Allreduce(MPI_IN_PLACE,terminallyIll,1,MPI_LOGICAL,MPI_LOR,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,terminallyIll,1,MPI_LOGICAL,MPI_LOR,PETSC_COMM_WORLD,ierr)
ForwardData = .False.
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! calculate divergence ! calculate divergence
tensorField_real = 0.0_pReal tensorField_real = 0.0_pReal
tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) = residual_F tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) = residual_F !< stress field in disguise
call utilities_FFTtensorForward() call utilities_FFTtensorForward()
err_div = Utilities_divergenceRMS() err_div = Utilities_divergenceRMS() !< root mean squared error in divergence of stress
call utilities_FFTtensorBackward()
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! constructing residual ! constructing residual
@ -463,7 +455,7 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
residual_F(1:3,1:3,i,j,k) - & residual_F(1:3,1:3,i,j,k) - &
math_mul33x33(F(1:3,1:3,i,j,k), & math_mul33x33(F(1:3,1:3,i,j,k), &
math_mul3333xx33(C_scale,F_lambda(1:3,1:3,i,j,k) - math_I3))) & math_mul3333xx33(C_scale,F_lambda(1:3,1:3,i,j,k) - math_I3))) &
+ residual_F_lambda(1:3,1:3,i,j,k) + residual_F_lambda(1:3,1:3,i,j,k) !< eq (16) in doi: 10.1016/j.ijplas.2014.02.006
enddo; enddo; enddo enddo; enddo; enddo
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -472,8 +464,11 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) = F tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) = F
call utilities_FFTtensorForward() call utilities_FFTtensorForward()
err_curl = Utilities_curlRMS() err_curl = Utilities_curlRMS()
call utilities_FFTtensorBackward()
nullify(F)
nullify(F_lambda)
nullify(residual_F)
nullify(residual_F_lambda)
end subroutine AL_formResidual end subroutine AL_formResidual
@ -488,8 +483,8 @@ subroutine AL_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,dummy,ierr
err_div_tolAbs, & err_div_tolAbs, &
err_curl_tolRel, & err_curl_tolRel, &
err_curl_tolAbs, & err_curl_tolAbs, &
err_stress_tolAbs, & err_stress_tolRel, &
err_stress_tolRel err_stress_tolAbs
use math, only: & use math, only: &
math_mul3333xx33 math_mul3333xx33
use FEsolving, only: & use FEsolving, only: &
@ -508,24 +503,24 @@ subroutine AL_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,dummy,ierr
real(pReal) :: & real(pReal) :: &
curlTol, & curlTol, &
divTol, & divTol, &
BC_tol BCTol
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! stress BC handling ! stress BC handling
F_aim = F_aim - math_mul3333xx33(S, ((P_av - params%stress_BC))) ! S = 0.0 for no bc F_aim = F_aim - math_mul3333xx33(S, ((P_av - params%stress_BC))) ! S = 0.0 for no bc
err_BC = maxval(abs((-mask_stress+1.0_pReal)*math_mul3333xx33(C_scale,F_aim-F_av) + & err_BC = maxval(abs((1.0_pReal-mask_stress) * math_mul3333xx33(C_scale,F_aim-F_av) + &
mask_stress * (P_av-params%stress_BC))) ! mask = 0.0 for no bc mask_stress * (P_av-params%stress_BC))) ! mask = 0.0 for no bc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! error calculation ! error calculation
curlTol = max(maxval(abs(F_aim-math_I3))*err_curl_tolRel ,err_curl_tolAbs) curlTol = max(maxval(abs(F_aim-math_I3))*err_curl_tolRel ,err_curl_tolAbs)
divTol = max(maxval(abs(P_av)) *err_div_tolRel ,err_div_tolAbs) divTol = max(maxval(abs(P_av)) *err_div_tolRel ,err_div_tolAbs)
BC_tol = max(maxval(abs(P_av)) *err_stress_tolrel,err_stress_tolabs) BCTol = max(maxval(abs(P_av)) *err_stress_tolRel,err_stress_tolAbs)
converged: if ((totalIter >= itmin .and. & converged: if ((totalIter >= itmin .and. &
all([ err_div /divTol, & all([ err_div /divTol, &
err_curl/curlTol, & err_curl/curlTol, &
err_BC/BC_tol ] < 1.0_pReal)) & err_BC /BCTol ] < 1.0_pReal)) &
.or. terminallyIll) then .or. terminallyIll) then
reason = 1 reason = 1
elseif (totalIter >= itmax) then converged elseif (totalIter >= itmax) then converged
@ -537,12 +532,12 @@ subroutine AL_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,dummy,ierr
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! report ! report
write(6,'(1/,a)') ' ... reporting .............................................................' write(6,'(1/,a)') ' ... reporting .............................................................'
write(6,'(/,a,f12.2,a,es8.2,a,es9.2,a)') ' error curl = ', & write(6,'(/,a,f12.2,a,es8.2,a,es9.2,a)') ' error divergence = ', &
err_curl/curlTol,' (',err_curl,' -, tol =',curlTol,')'
write(6,' (a,f12.2,a,es8.2,a,es9.2,a)') ' error divergence = ', &
err_div/divTol, ' (',err_div, ' / m, tol = ',divTol,')' err_div/divTol, ' (',err_div, ' / m, tol = ',divTol,')'
write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error curl = ', &
err_curl/curlTol,' (',err_curl,' -, tol = ',curlTol,')'
write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error BC = ', & write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error BC = ', &
err_BC/BC_tol, ' (',err_BC, ' Pa, tol =',BC_tol,')' err_BC/BCTol, ' (',err_BC, ' Pa, tol = ',BCTol,')'
write(6,'(/,a)') ' ===========================================================================' write(6,'(/,a)') ' ==========================================================================='
flush(6) flush(6)
@ -550,6 +545,8 @@ end subroutine AL_converged
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief forwarding routine !> @brief forwarding routine
!> @details find new boundary conditions and best F estimate for end of current timestep
!> possibly writing restart information, triggering of state increment in DAMASK, and updating of IPcoordinates
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine AL_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation_BC,stress_BC,rotation_BC) subroutine AL_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation_BC,stress_BC,rotation_BC)
use math, only: & use math, only: &
@ -559,9 +556,13 @@ subroutine AL_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation_BC,stre
math_rotate_backward33 math_rotate_backward33
use numerics, only: & use numerics, only: &
worldrank worldrank
use homogenization, only: &
materialpoint_F0
use mesh, only: & use mesh, only: &
grid3, & grid, &
grid grid3
use CPFEM2, only: &
CPFEM_age
use spectral_utilities, only: & use spectral_utilities, only: &
Utilities_calculateRate, & Utilities_calculateRate, &
Utilities_forwardField, & Utilities_forwardField, &
@ -574,6 +575,8 @@ subroutine AL_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation_BC,stre
restartWrite restartWrite
implicit none implicit none
logical, intent(in) :: &
guess
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
timeinc_old, & timeinc_old, &
timeinc, & timeinc, &
@ -581,100 +584,102 @@ subroutine AL_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation_BC,stre
type(tBoundaryCondition), intent(in) :: & type(tBoundaryCondition), intent(in) :: &
stress_BC, & stress_BC, &
deformation_BC deformation_BC
real(pReal), dimension(3,3), intent(in) :: rotation_BC real(pReal), dimension(3,3), intent(in) ::&
logical, intent(in) :: & rotation_BC
guess
PetscErrorCode :: ierr PetscErrorCode :: ierr
PetscScalar, dimension(:,:,:,:), pointer :: xx_psc, F, F_lambda PetscScalar, dimension(:,:,:,:), pointer :: FandF_lambda, F, F_lambda
integer(pInt) :: i, j, k integer(pInt) :: i, j, k
real(pReal), dimension(3,3) :: F_lambda33 real(pReal), dimension(3,3) :: F_lambda33
character(len=1024) :: rankStr character(len=32) :: rankStr
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! update coordinates and rate and forward last inc ! update coordinates and rate and forward last inc
call DMDAVecGetArrayF90(da,solution_vec,xx_psc,ierr) call DMDAVecGetArrayF90(da,solution_vec,FandF_lambda,ierr); CHKERRQ(ierr)
F => xx_psc(0:8,:,:,:) F => FandF_lambda( 0: 8,:,:,:)
F_lambda => xx_psc(9:17,:,:,:) F_lambda => FandF_lambda( 9:17,:,:,:)
if (restartWrite) then
write(6,'(/,a)') ' writing converged results for restart'
flush(6)
write(rankStr,'(a1,i0)')'_',worldrank
call IO_write_jobRealFile(777,'F'//trim(rankStr),size(F)) ! writing deformation gradient field to file
write (777,rec=1) F
close (777)
call IO_write_jobRealFile(777,'F_lastInc'//trim(rankStr),size(F_lastInc)) ! writing F_lastInc field to file
write (777,rec=1) F_lastInc
close (777)
call IO_write_jobRealFile(777,'F_lambda'//trim(rankStr),size(F_lambda)) ! writing deformation gradient field to file
write (777,rec=1) F_lambda
close (777)
call IO_write_jobRealFile(777,'F_lambda_lastInc'//trim(rankStr),size(F_lambda_lastInc)) ! writing F_lastInc field to file
write (777,rec=1) F_lambda_lastInc
close (777)
if (worldrank == 0_pInt) then
call IO_write_jobRealFile(777,'F_aim',size(F_aim))
write (777,rec=1) F_aim
close(777)
call IO_write_jobRealFile(777,'F_aim_lastInc',size(F_aim_lastInc))
write (777,rec=1) F_aim_lastInc
close(777)
call IO_write_jobRealFile(777,'F_aimDot',size(F_aimDot))
write (777,rec=1) F_aimDot
close(777)
call IO_write_jobRealFile(777,'C_volAvg',size(C_volAvg))
write (777,rec=1) C_volAvg
close(777)
call IO_write_jobRealFile(777,'C_volAvgLastInc',size(C_volAvgLastInc))
write (777,rec=1) C_volAvgLastInc
close(777)
endif
endif
call utilities_updateIPcoords(F)
if (cutBack) then if (cutBack) then
F_aim = F_aim_lastInc C_volAvg = C_volAvgLastInc ! QUESTION: where is this required?
F_lambda = reshape(F_lambda_lastInc,[9,grid(1),grid(2),grid3]) C_minMaxAvg = C_minMaxAvgLastInc ! QUESTION: where is this required?
F = reshape(F_lastInc, [9,grid(1),grid(2),grid3])
C_volAvg = C_volAvgLastInc
else else
ForwardData = .True. !--------------------------------------------------------------------------------------------------
! restart information for spectral solver
if (restartWrite) then ! QUESTION: where is this logical properly set?
write(6,'(/,a)') ' writing converged results for restart'
flush(6)
if (worldrank == 0_pInt) then
call IO_write_jobRealFile(777,'C_volAvg',size(C_volAvg))
write (777,rec=1) C_volAvg; close(777)
call IO_write_jobRealFile(777,'C_volAvgLastInc',size(C_volAvgLastInc))
write (777,rec=1) C_volAvgLastInc; close(777)
! call IO_write_jobRealFile(777,'C_minMaxAvg',size(C_volAvg))
! write (777,rec=1) C_minMaxAvg; close(777)
! call IO_write_jobRealFile(777,'C_minMaxAvgLastInc',size(C_volAvgLastInc))
! write (777,rec=1) C_minMaxAvgLastInc; close(777)
call IO_write_jobRealFile(777,'F_aimDot',size(F_aimDot))
write (777,rec=1) F_aimDot; close(777)
endif
write(rankStr,'(a1,i0)')'_',worldrank
call IO_write_jobRealFile(777,'F'//trim(rankStr),size(F)) ! writing deformation gradient field to file
write (777,rec=1) F; close (777)
call IO_write_jobRealFile(777,'F_lastInc'//trim(rankStr),size(F_lastInc)) ! writing F_lastInc field to file
write (777,rec=1) F_lastInc; close (777)
call IO_write_jobRealFile(777,'F_lambda'//trim(rankStr),size(F_lambda)) ! writing deformation gradient field to file
write (777,rec=1) F_lambda; close (777)
call IO_write_jobRealFile(777,'F_lambda_lastInc'//trim(rankStr),size(F_lambda_lastInc)) ! writing F_lastInc field to file
write (777,rec=1) F_lambda_lastInc; close (777)
endif
call CPFEM_age() ! age state and kinematics
call utilities_updateIPcoords(F)
C_volAvgLastInc = C_volAvg C_volAvgLastInc = C_volAvg
C_minMaxAvgLastInc = C_minMaxAvg
if (guess) then ! QUESTION: better with a = L ? x:y
F_aimDot = stress_BC%maskFloat * (F_aim - F_aim_lastInc)/timeinc_old ! initialize with correction based on last inc
else
F_aimDot = 0.0_pReal
endif
F_aim_lastInc = F_aim
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! calculate rate for aim ! calculate rate for aim
if (deformation_BC%myType=='l') then ! calculate f_aimDot from given L and current F if (deformation_BC%myType=='l') then ! calculate F_aimDot from given L and current F
f_aimDot = deformation_BC%maskFloat * math_mul33x33(deformation_BC%values, F_aim) F_aimDot = &
elseif(deformation_BC%myType=='fdot') then ! f_aimDot is prescribed F_aimDot + deformation_BC%maskFloat * math_mul33x33(deformation_BC%values, F_aim_lastInc)
f_aimDot = deformation_BC%maskFloat * deformation_BC%values elseif(deformation_BC%myType=='fdot') then ! F_aimDot is prescribed
F_aimDot = &
F_aimDot + deformation_BC%maskFloat * deformation_BC%values
elseif (deformation_BC%myType=='f') then ! aim at end of load case is prescribed elseif (deformation_BC%myType=='f') then ! aim at end of load case is prescribed
f_aimDot = deformation_BC%maskFloat * (deformation_BC%values -F_aim)/loadCaseTime F_aimDot = &
F_aimDot + deformation_BC%maskFloat * (deformation_BC%values - F_aim_lastInc)/loadCaseTime
endif endif
if (guess) f_aimDot = f_aimDot + stress_BC%maskFloat * (F_aim - F_aim_lastInc)/timeinc_old
F_aim_lastInc = F_aim
!--------------------------------------------------------------------------------------------------
! update coordinates and rate and forward last inc
call utilities_updateIPcoords(F)
Fdot = Utilities_calculateRate(guess, & Fdot = Utilities_calculateRate(guess, &
F_lastInc,reshape(F,[3,3,grid(1),grid(2),grid3]),timeinc_old, & F_lastInc,reshape(F,[3,3,grid(1),grid(2),grid3]),timeinc_old, &
math_rotate_backward33(f_aimDot,rotation_BC)) math_rotate_backward33(F_aimDot,rotation_BC))
F_lambdaDot = Utilities_calculateRate(guess, & F_lambdaDot = Utilities_calculateRate(guess, &
F_lambda_lastInc,reshape(F_lambda,[3,3,grid(1),grid(2),grid3]), timeinc_old, & F_lambda_lastInc,reshape(F_lambda,[3,3,grid(1),grid(2),grid3]), timeinc_old, &
math_rotate_backward33(f_aimDot,rotation_BC)) math_rotate_backward33(F_aimDot,rotation_BC))
F_lastInc = reshape(F, [3,3,grid(1),grid(2),grid3]) F_lastInc = reshape(F, [3,3,grid(1),grid(2),grid3]) ! winding F forward
F_lambda_lastInc = reshape(F_lambda,[3,3,grid(1),grid(2),grid3]) F_lambda_lastInc = reshape(F_lambda, [3,3,grid(1),grid(2),grid3]) ! winding F_lambda forward
materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent
endif endif
F_aim = F_aim + f_aimDot * timeinc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! update local deformation gradient ! update average and local deformation gradients
F = reshape(Utilities_forwardField(timeinc,F_lastInc,Fdot, & ! ensure that it matches rotated F_aim F_aim = F_aim_lastInc + F_aimDot * timeinc
F = reshape(Utilities_forwardField(timeinc,F_lastInc,Fdot, & ! estimate of F at end of time+timeinc that matches rotated F_aim on average
math_rotate_backward33(F_aim,rotation_BC)),& math_rotate_backward33(F_aim,rotation_BC)),&
[9,grid(1),grid(2),grid3]) [9,grid(1),grid(2),grid3])
if (guess) then
F_lambda = reshape(Utilities_forwardField(timeinc,F_lambda_lastInc,F_lambdadot), & F_lambda = reshape(Utilities_forwardField(timeinc,F_lambda_lastInc,F_lambdadot), &
[9,grid(1),grid(2),grid3]) ! does not have any average value as boundary condition [9,grid(1),grid(2),grid3]) ! does not have any average value as boundary condition
if (.not. guess) then ! large strain forwarding else
do k = 1_pInt, grid3; do j = 1_pInt, grid(2); do i = 1_pInt, grid(1) do k = 1_pInt, grid3; do j = 1_pInt, grid(2); do i = 1_pInt, grid(1)
F_lambda33 = reshape(F_lambda(1:9,i,j,k),[3,3]) F_lambda33 = reshape(F_lambda(1:9,i,j,k),[3,3])
F_lambda33 = math_mul3333xx33(S_scale,math_mul33x33(F_lambda33, & F_lambda33 = math_mul3333xx33(S_scale,math_mul33x33(F_lambda33, &
@ -685,7 +690,10 @@ subroutine AL_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation_BC,stre
F_lambda(1:9,i,j,k) = reshape(F_lambda33,[9]) F_lambda(1:9,i,j,k) = reshape(F_lambda33,[9])
enddo; enddo; enddo enddo; enddo; enddo
endif endif
call DMDAVecRestoreArrayF90(da,solution_vec,xx_psc,ierr); CHKERRQ(ierr)
nullify(F)
nullify(F_lambda)
call DMDAVecRestoreArrayF90(da,solution_vec,FandF_lambda,ierr); CHKERRQ(ierr)
end subroutine AL_forward end subroutine AL_forward

View File

@ -38,21 +38,27 @@ module spectral_mech_basic
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! stress, stiffness and compliance average etc. ! stress, stiffness and compliance average etc.
real(pReal), private, dimension(3,3) :: & real(pReal), private, dimension(3,3) :: &
F_aim = math_I3, & F_aimDot = 0.0_pReal, & !< assumed rate of average deformation gradient
F_aim_lastInc = math_I3, & F_aim = math_I3, & !< current prescribed deformation gradient
P_av = 0.0_pReal, & F_aim_lastInc = math_I3, & !< previous average deformation gradient
F_aimDot = 0.0_pReal P_av = 0.0_pReal !< average 1st Piola--Kirchhoff stress
character(len=1024), private :: incInfo
character(len=1024), private :: incInfo !< time and increment information
real(pReal), private, dimension(3,3,3,3) :: & real(pReal), private, dimension(3,3,3,3) :: &
C_volAvg = 0.0_pReal, & !< current volume average stiffness C_volAvg = 0.0_pReal, & !< current volume average stiffness
C_volAvgLastInc = 0.0_pReal, & !< previous volume average stiffness C_volAvgLastInc = 0.0_pReal, & !< previous volume average stiffness
C_minMaxAvg = 0.0_pReal, & !< current (min+max)/2 stiffness C_minMaxAvg = 0.0_pReal, & !< current (min+max)/2 stiffness
C_minMaxAvgLastInc = 0.0_pReal, & !< previous (min+max)/2 stiffness C_minMaxAvgLastInc = 0.0_pReal, & !< previous (min+max)/2 stiffness
S = 0.0_pReal !< current compliance (filled up with zeros) S = 0.0_pReal !< current compliance (filled up with zeros)
real(pReal), private :: err_stress, err_div
logical, private :: ForwardData real(pReal), private :: &
err_BC, & !< deviation from stress BC
err_div !< RMS of div of P
integer(pInt), private :: & integer(pInt), private :: &
totalIter = 0_pInt !< total iteration in current increment totalIter = 0_pInt !< total iteration in current increment
real(pReal), private, dimension(3,3) :: mask_stress = 0.0_pReal real(pReal), private, dimension(3,3) :: mask_stress = 0.0_pReal
public :: & public :: &
@ -72,7 +78,7 @@ contains
!> @brief allocates all necessary fields and fills them with data, potentially from restart info !> @brief allocates all necessary fields and fills them with data, potentially from restart info
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine basicPETSc_init subroutine basicPETSc_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -154,11 +160,10 @@ subroutine basicPETSc_init
grid(1),grid(2),localK, & ! local grid grid(1),grid(2),localK, & ! local grid
da,ierr) ! handle, error da,ierr) ! handle, error
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetDM(snes,da,ierr); CHKERRQ(ierr) call SNESSetDM(snes,da,ierr); CHKERRQ(ierr) ! connect snes to da
call DMCreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr) ! global solution vector (grid x 9, i.e. every def grad tensor) call DMCreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr) ! global solution vector (grid x 9, i.e. every def grad tensor)
call DMDASNESSetFunctionLocal(da,INSERT_VALUES,BasicPETSC_formResidual,PETSC_NULL_OBJECT,ierr) ! residual vector of same shape as solution vector call DMDASNESSetFunctionLocal(da,INSERT_VALUES,BasicPETSC_formResidual,PETSC_NULL_OBJECT,ierr) ! residual vector of same shape as solution vector
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetDM(snes,da,ierr); CHKERRQ(ierr) ! connect snes to da
call SNESSetConvergenceTest(snes,BasicPETSC_converged,PETSC_NULL_OBJECT,PETSC_NULL_FUNCTION,ierr) ! specify custom convergence check function "_converged" call SNESSetConvergenceTest(snes,BasicPETSC_converged,PETSC_NULL_OBJECT,PETSC_NULL_FUNCTION,ierr) ! specify custom convergence check function "_converged"
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetFromOptions(snes,ierr); CHKERRQ(ierr) ! pull it all together with additional cli arguments call SNESSetFromOptions(snes,ierr); CHKERRQ(ierr) ! pull it all together with additional cli arguments
@ -167,10 +172,10 @@ subroutine basicPETSc_init
! init fields ! init fields
call DMDAVecGetArrayF90(da,solution_vec,F,ierr); CHKERRQ(ierr) ! get the data out of PETSc to work with call DMDAVecGetArrayF90(da,solution_vec,F,ierr); CHKERRQ(ierr) ! get the data out of PETSc to work with
restart: if (restartInc > 1_pInt) then restart: if (restartInc > 0_pInt) then
if (iand(debug_level(debug_spectral),debug_spectralRestart) /= 0) then if (iand(debug_level(debug_spectral),debug_spectralRestart) /= 0) then
write(6,'(/,a,'//IO_intOut(restartInc-1_pInt)//',a)') & write(6,'(/,a,'//IO_intOut(restartInc)//',a)') &
'reading values of increment ', restartInc - 1_pInt, ' from file' 'reading values of increment ', restartInc, ' from file'
flush(6) flush(6)
endif endif
write(rankStr,'(a1,i0)')'_',worldrank write(rankStr,'(a1,i0)')'_',worldrank
@ -182,7 +187,7 @@ subroutine basicPETSc_init
read (777,rec=1) F_aimDot; close (777) read (777,rec=1) F_aimDot; close (777)
F_aim = reshape(sum(sum(sum(F,dim=4),dim=3),dim=2) * wgt, [3,3]) ! average of F F_aim = reshape(sum(sum(sum(F,dim=4),dim=3),dim=2) * wgt, [3,3]) ! average of F
F_aim_lastInc = sum(sum(sum(F_lastInc,dim=5),dim=4),dim=3) * wgt ! average of F_lastInc F_aim_lastInc = sum(sum(sum(F_lastInc,dim=5),dim=4),dim=3) * wgt ! average of F_lastInc
elseif (restartInc == 1_pInt) then restart elseif (restartInc == 0_pInt) then restart
F_lastInc = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid3) ! initialize to identity F_lastInc = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid3) ! initialize to identity
F = reshape(F_lastInc,[9,grid(1),grid(2),grid3]) F = reshape(F_lastInc,[9,grid(1),grid(2),grid3])
endif restart endif restart
@ -196,10 +201,10 @@ subroutine basicPETSc_init
call DMDAVecRestoreArrayF90(da,solution_vec,F,ierr); CHKERRQ(ierr) ! write data back to PETSc call DMDAVecRestoreArrayF90(da,solution_vec,F,ierr); CHKERRQ(ierr) ! write data back to PETSc
! QUESTION: why not writing back right after reading (l.189)? ! QUESTION: why not writing back right after reading (l.189)?
restartRead: if (restartInc > 1_pInt) then ! QUESTION: are those values not calc'ed by constitutiveResponse? why reading from file? restartRead: if (restartInc > 0_pInt) then ! QUESTION: are those values not calc'ed by constitutiveResponse? why reading from file?
if (iand(debug_level(debug_spectral),debug_spectralRestart) /= 0 .and. worldrank == 0_pInt) & if (iand(debug_level(debug_spectral),debug_spectralRestart) /= 0 .and. worldrank == 0_pInt) &
write(6,'(/,a,'//IO_intOut(restartInc-1_pInt)//',a)') & write(6,'(/,a,'//IO_intOut(restartInc)//',a)') &
'reading more values of increment', restartInc-1_pInt, 'from file' 'reading more values of increment ', restartInc, ' from file'
flush(6) flush(6)
call IO_read_realFile(777,'C_volAvg',trim(getSolverJobName()),size(C_volAvg)) call IO_read_realFile(777,'C_volAvg',trim(getSolverJobName()),size(C_volAvg))
read (777,rec=1) C_volAvg; close (777) read (777,rec=1) C_volAvg; close (777)
@ -209,7 +214,7 @@ subroutine basicPETSc_init
read (777,rec=1) C_minMaxAvg; close (777) read (777,rec=1) C_minMaxAvg; close (777)
endif restartRead endif restartRead
call Utilities_updateGamma(C_minmaxAvg,.true.) call Utilities_updateGamma(C_minMaxAvg,.true.)
end subroutine basicPETSc_init end subroutine basicPETSc_init
@ -256,7 +261,7 @@ type(tSolutionState) function basicPETSc_solution(incInfoIn,timeinc,timeinc_old,
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! update stiffness (and gamma operator) ! update stiffness (and gamma operator)
S = Utilities_maskedCompliance(rotation_BC,stress_BC%maskLogical,C_volAvg) S = Utilities_maskedCompliance(rotation_BC,stress_BC%maskLogical,C_volAvg)
if (update_gamma) call Utilities_updateGamma(C_minmaxAvg,restartWrite) if (update_gamma) call Utilities_updateGamma(C_minMaxAvg,restartWrite)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -269,8 +274,7 @@ type(tSolutionState) function basicPETSc_solution(incInfoIn,timeinc,timeinc_old,
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! solve BVP ! solve BVP
call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr) call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr); CHKERRQ(ierr)
CHKERRQ(ierr)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! check convergence ! check convergence
@ -282,7 +286,6 @@ type(tSolutionState) function basicPETSc_solution(incInfoIn,timeinc,timeinc_old,
terminallyIll = .false. terminallyIll = .false.
if (reason == -4) call IO_error(893_pInt) ! MPI error if (reason == -4) call IO_error(893_pInt) ! MPI error
end function BasicPETSc_solution end function BasicPETSc_solution
@ -298,7 +301,6 @@ subroutine BasicPETSC_formResidual(in,x_scal,f_scal,dummy,ierr)
grid3 grid3
use math, only: & use math, only: &
math_rotate_backward33, & math_rotate_backward33, &
math_transpose33, &
math_mul3333xx33 math_mul3333xx33
use debug, only: & use debug, only: &
debug_level, & debug_level, &
@ -346,15 +348,15 @@ subroutine BasicPETSC_formResidual(in,x_scal,f_scal,dummy,ierr)
trim(incInfo), ' @ Iteration ', itmin, '≤',totalIter, '≤', itmax trim(incInfo), ' @ Iteration ', itmin, '≤',totalIter, '≤', itmax
if (iand(debug_level(debug_spectral),debug_spectralRotation) /= 0) & if (iand(debug_level(debug_spectral),debug_spectralRotation) /= 0) &
write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') & write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') &
' deformation gradient aim (lab) =', math_transpose33(math_rotate_backward33(F_aim,params%rotation_BC)) ' deformation gradient aim (lab) =', transpose(math_rotate_backward33(F_aim,params%rotation_BC))
write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') & write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') &
' deformation gradient aim =', math_transpose33(F_aim) ' deformation gradient aim =', transpose(F_aim)
flush(6) flush(6)
endif newIteration endif newIteration
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! evaluate constitutive response ! evaluate constitutive response
call Utilities_constitutiveResponse(f_scal,P_av,C_volAvg,C_minmaxAvg, & call Utilities_constitutiveResponse(f_scal,P_av,C_volAvg,C_minMaxAvg, &
x_scal,params%timeinc, params%rotation_BC) x_scal,params%timeinc, params%rotation_BC)
call MPI_Allreduce(MPI_IN_PLACE,terminallyIll,1,MPI_LOGICAL,MPI_LOR,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,terminallyIll,1,MPI_LOGICAL,MPI_LOR,PETSC_COMM_WORLD,ierr)
@ -362,7 +364,7 @@ subroutine BasicPETSC_formResidual(in,x_scal,f_scal,dummy,ierr)
! stress BC handling ! stress BC handling
deltaF_aim = math_mul3333xx33(S, P_av - params%stress_BC) deltaF_aim = math_mul3333xx33(S, P_av - params%stress_BC)
F_aim = F_aim - deltaF_aim F_aim = F_aim - deltaF_aim
err_stress = maxval(abs(mask_stress * (P_av - params%stress_BC))) ! mask = 0.0 when no stress bc err_BC = maxval(abs(mask_stress * (P_av - params%stress_BC))) ! mask = 0.0 when no stress bc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! updated deformation gradient using fix point algorithm of basic scheme ! updated deformation gradient using fix point algorithm of basic scheme
@ -406,14 +408,14 @@ subroutine BasicPETSc_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,du
PetscErrorCode :: ierr PetscErrorCode :: ierr
real(pReal) :: & real(pReal) :: &
divTol, & divTol, &
stressTol BCTol
divTol = max(maxval(abs(P_av))*err_div_tolRel ,err_div_tolAbs) divTol = max(maxval(abs(P_av))*err_div_tolRel ,err_div_tolAbs)
stressTol = max(maxval(abs(P_av))*err_stress_tolrel,err_stress_tolabs) BCTol = max(maxval(abs(P_av))*err_stress_tolRel,err_stress_tolAbs)
converged: if ((totalIter >= itmin .and. & converged: if ((totalIter >= itmin .and. &
all([ err_div/divTol, & all([ err_div/divTol, &
err_stress/stressTol ] < 1.0_pReal)) & err_BC /BCTol ] < 1.0_pReal)) &
.or. terminallyIll) then .or. terminallyIll) then
reason = 1 reason = 1
elseif (totalIter >= itmax) then converged elseif (totalIter >= itmax) then converged
@ -428,7 +430,7 @@ subroutine BasicPETSc_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,du
write(6,'(1/,a,f12.2,a,es8.2,a,es9.2,a)') ' error divergence = ', & write(6,'(1/,a,f12.2,a,es8.2,a,es9.2,a)') ' error divergence = ', &
err_div/divTol, ' (',err_div,' / m, tol = ',divTol,')' err_div/divTol, ' (',err_div,' / m, tol = ',divTol,')'
write(6,'(a,f12.2,a,es8.2,a,es9.2,a)') ' error stress BC = ', & write(6,'(a,f12.2,a,es8.2,a,es9.2,a)') ' error stress BC = ', &
err_stress/stressTol, ' (',err_stress, ' Pa, tol =',stressTol,')' err_BC/BCTol, ' (',err_BC, ' Pa, tol = ',BCTol,')'
write(6,'(/,a)') ' ===========================================================================' write(6,'(/,a)') ' ==========================================================================='
flush(6) flush(6)
@ -476,7 +478,7 @@ subroutine BasicPETSc_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation
real(pReal), dimension(3,3), intent(in) ::& real(pReal), dimension(3,3), intent(in) ::&
rotation_BC rotation_BC
PetscErrorCode :: ierr PetscErrorCode :: ierr
PetscScalar, pointer :: F(:,:,:,:) PetscScalar, dimension(:,:,:,:), pointer :: F
character(len=32) :: rankStr character(len=32) :: rankStr
@ -497,10 +499,6 @@ subroutine BasicPETSc_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation
write (777,rec=1) C_volAvg; close(777) write (777,rec=1) C_volAvg; close(777)
call IO_write_jobRealFile(777,'C_volAvgLastInc',size(C_volAvgLastInc)) call IO_write_jobRealFile(777,'C_volAvgLastInc',size(C_volAvgLastInc))
write (777,rec=1) C_volAvgLastInc; close(777) write (777,rec=1) C_volAvgLastInc; close(777)
call IO_write_jobRealFile(777,'C_minMaxAvg',size(C_volAvg))
write (777,rec=1) C_minMaxAvg; close(777)
call IO_write_jobRealFile(777,'C_minMaxAvgLastInc',size(C_volAvgLastInc))
write (777,rec=1) C_minMaxAvgLastInc; close(777)
call IO_write_jobRealFile(777,'F_aimDot',size(F_aimDot)) call IO_write_jobRealFile(777,'F_aimDot',size(F_aimDot))
write (777,rec=1) F_aimDot; close(777) write (777,rec=1) F_aimDot; close(777)
endif endif
@ -518,12 +516,9 @@ subroutine BasicPETSc_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation
C_volAvgLastInc = C_volAvg C_volAvgLastInc = C_volAvg
C_minMaxAvgLastInc = C_minMaxAvg C_minMaxAvgLastInc = C_minMaxAvg
if (guess) then ! QUESTION: better with a = L ? x:y F_aimDot = merge(stress_BC%maskFloat*(F_aim-F_aim_lastInc)/timeinc_old, 0.0_pReal, guess)
F_aimDot = stress_BC%maskFloat * (F_aim - F_aim_lastInc)/timeinc_old ! initialize with correction based on last inc
else
F_aimDot = 0.0_pReal
endif
F_aim_lastInc = F_aim F_aim_lastInc = F_aim
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! calculate rate for aim ! calculate rate for aim
if (deformation_BC%myType=='l') then ! calculate F_aimDot from given L and current F if (deformation_BC%myType=='l') then ! calculate F_aimDot from given L and current F

View File

@ -43,7 +43,7 @@ module spectral_mech_Polarisation
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! stress, stiffness and compliance average etc. ! stress, stiffness and compliance average etc.
real(pReal), private, dimension(3,3) :: & real(pReal), private, dimension(3,3) :: &
F_aimDot, & !< assumed rate of average deformation gradient F_aimDot = 0.0_pReal, & !< assumed rate of average deformation gradient
F_aim = math_I3, & !< current prescribed deformation gradient F_aim = math_I3, & !< current prescribed deformation gradient
F_aim_lastInc = math_I3, & !< previous average deformation gradient F_aim_lastInc = math_I3, & !< previous average deformation gradient
F_av = 0.0_pReal, & !< average incompatible def grad field F_av = 0.0_pReal, & !< average incompatible def grad field
@ -54,6 +54,7 @@ module spectral_mech_Polarisation
C_volAvg = 0.0_pReal, & !< current volume average stiffness C_volAvg = 0.0_pReal, & !< current volume average stiffness
C_volAvgLastInc = 0.0_pReal, & !< previous volume average stiffness C_volAvgLastInc = 0.0_pReal, & !< previous volume average stiffness
C_minMaxAvg = 0.0_pReal, & !< current (min+max)/2 stiffness C_minMaxAvg = 0.0_pReal, & !< current (min+max)/2 stiffness
C_minMaxAvgLastInc = 0.0_pReal, & !< previous (min+max)/2 stiffness
S = 0.0_pReal, & !< current compliance (filled up with zeros) S = 0.0_pReal, & !< current compliance (filled up with zeros)
C_scale = 0.0_pReal, & C_scale = 0.0_pReal, &
S_scale = 0.0_pReal S_scale = 0.0_pReal
@ -62,7 +63,7 @@ module spectral_mech_Polarisation
err_BC, & !< deviation from stress BC err_BC, & !< deviation from stress BC
err_curl, & !< RMS of curl of F err_curl, & !< RMS of curl of F
err_div !< RMS of div of P err_div !< RMS of div of P
logical, private :: ForwardData
integer(pInt), private :: & integer(pInt), private :: &
totalIter = 0_pInt !< total iteration in current increment totalIter = 0_pInt !< total iteration in current increment
@ -80,11 +81,11 @@ module spectral_mech_Polarisation
contains contains
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief allocates all neccessary fields and fills them with data, potentially from restart info !> @brief allocates all necessary fields and fills them with data, potentially from restart info
!> @todo use sourced allocation, e.g. allocate(Fdot,source = F_lastInc) !> @todo use sourced allocation, e.g. allocate(Fdot,source = F_lastInc)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine Polarisation_init subroutine Polarisation_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -102,12 +103,15 @@ subroutine Polarisation_init
use numerics, only: & use numerics, only: &
worldrank, & worldrank, &
worldsize worldsize
use homogenization, only: &
materialpoint_F0
use DAMASK_interface, only: & use DAMASK_interface, only: &
getSolverJobName getSolverJobName
use spectral_utilities, only: & use spectral_utilities, only: &
Utilities_constitutiveResponse, & Utilities_constitutiveResponse, &
Utilities_updateGamma, & Utilities_updateGamma, &
Utilities_updateIPcoords Utilities_updateIPcoords, &
wgt
use mesh, only: & use mesh, only: &
grid, & grid, &
grid3 grid3
@ -120,7 +124,11 @@ subroutine Polarisation_init
temp33_Real = 0.0_pReal temp33_Real = 0.0_pReal
PetscErrorCode :: ierr PetscErrorCode :: ierr
PetscScalar, pointer, dimension(:,:,:,:) :: xx_psc, F, F_tau PetscScalar, pointer, dimension(:,:,:,:) :: &
FandF_tau, & ! overall pointer to solution data
F, & ! specific (sub)pointer
F_tau ! specific (sub)pointer
integer(pInt), dimension(:), allocatable :: localK integer(pInt), dimension(:), allocatable :: localK
integer(pInt) :: proc integer(pInt) :: proc
character(len=1024) :: rankStr character(len=1024) :: rankStr
@ -164,78 +172,69 @@ subroutine Polarisation_init
grid(1),grid(2),localK, & ! local grid grid(1),grid(2),localK, & ! local grid
da,ierr) ! handle, error da,ierr) ! handle, error
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetDM(snes,da,ierr); CHKERRQ(ierr) call SNESSetDM(snes,da,ierr); CHKERRQ(ierr) ! connect snes to da
call DMCreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr) call DMCreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr) ! global solution vector (grid x 9, i.e. every def grad tensor)
call DMDASNESSetFunctionLocal(da,INSERT_VALUES,Polarisation_formResidual,PETSC_NULL_OBJECT,ierr) call DMDASNESSetFunctionLocal(da,INSERT_VALUES,Polarisation_formResidual,PETSC_NULL_OBJECT,ierr) ! residual vector of same shape as solution vector
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetConvergenceTest(snes,Polarisation_converged,PETSC_NULL_OBJECT,PETSC_NULL_FUNCTION,ierr) call SNESSetConvergenceTest(snes,Polarisation_converged,PETSC_NULL_OBJECT,PETSC_NULL_FUNCTION,ierr) ! specify custom convergence check function "_converged"
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetFromOptions(snes,ierr); CHKERRQ(ierr) call SNESSetFromOptions(snes,ierr); CHKERRQ(ierr) ! pull it all together with additional cli arguments
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! init fields ! init fields
call DMDAVecGetArrayF90(da,solution_vec,xx_psc,ierr); CHKERRQ(ierr) ! places pointer xx_psc on PETSc data call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr) ! places pointer on PETSc data
F => xx_psc(0:8,:,:,:) F => FandF_tau( 0: 8,:,:,:)
F_tau => xx_psc(9:17,:,:,:) F_tau => FandF_tau( 9:17,:,:,:)
restart: if (restartInc > 1_pInt) then restart: if (restartInc > 0_pInt) then
if (iand(debug_level(debug_spectral),debug_spectralRestart)/= 0) & if (iand(debug_level(debug_spectral),debug_spectralRestart) /= 0) then
write(6,'(/,a,'//IO_intOut(restartInc-1_pInt)//',a)') & write(6,'(/,a,'//IO_intOut(restartInc)//',a)') &
'reading values of increment ', restartInc - 1_pInt, ' from file' 'reading values of increment ', restartInc, ' from file'
flush(6) flush(6)
endif
write(rankStr,'(a1,i0)')'_',worldrank write(rankStr,'(a1,i0)')'_',worldrank
call IO_read_realFile(777,'F'//trim(rankStr),trim(getSolverJobName()),size(F)) call IO_read_realFile(777,'F'//trim(rankStr),trim(getSolverJobName()),size(F))
read (777,rec=1) F read (777,rec=1) F; close (777)
close (777)
call IO_read_realFile(777,'F_lastInc'//trim(rankStr),trim(getSolverJobName()),size(F_lastInc)) call IO_read_realFile(777,'F_lastInc'//trim(rankStr),trim(getSolverJobName()),size(F_lastInc))
read (777,rec=1) F_lastInc read (777,rec=1) F_lastInc; close (777)
close (777)
call IO_read_realFile(777,'F_tau'//trim(rankStr),trim(getSolverJobName()),size(F_tau)) call IO_read_realFile(777,'F_tau'//trim(rankStr),trim(getSolverJobName()),size(F_tau))
read (777,rec=1) F_tau read (777,rec=1) F_tau; close (777)
close (777) call IO_read_realFile(777,'F_tau_lastInc'//trim(rankStr),trim(getSolverJobName()),size(F_tau_lastInc))
call IO_read_realFile(777,'F_tau_lastInc'//trim(rankStr),& read (777,rec=1) F_tau_lastInc; close (777)
trim(getSolverJobName()),size(F_tau_lastInc)) call IO_read_realFile(777,'F_aimDot',trim(getSolverJobName()),size(F_aimDot))
read (777,rec=1) F_tau_lastInc read (777,rec=1) F_aimDot; close (777)
close (777) F_aim = reshape(sum(sum(sum(F,dim=4),dim=3),dim=2) * wgt, [3,3]) ! average of F
call IO_read_realFile(777,'F_aim', trim(getSolverJobName()),size(F_aim)) F_aim_lastInc = sum(sum(sum(F_lastInc,dim=5),dim=4),dim=3) * wgt ! average of F_lastInc
read (777,rec=1) F_aim elseif (restartInc == 0_pInt) then restart
close (777)
call IO_read_realFile(777,'F_aim_lastInc', trim(getSolverJobName()),size(F_aim_lastInc))
read (777,rec=1) F_aim_lastInc
close (777)
call IO_read_realFile(777,'F_aimDot',trim(getSolverJobName()),size(f_aimDot))
read (777,rec=1) f_aimDot
close (777)
elseif (restartInc == 1_pInt) then restart
F_lastInc = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid3) ! initialize to identity F_lastInc = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid3) ! initialize to identity
F = reshape(F_lastInc,[9,grid(1),grid(2),grid3]) F = reshape(F_lastInc,[9,grid(1),grid(2),grid3])
F_tau = 2.0_pReal*F F_tau = 2.0_pReal*F
F_tau_lastInc = 2.0_pReal*F_lastInc F_tau_lastInc = 2.0_pReal*F_lastInc
endif restart endif restart
materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent
call Utilities_updateIPcoords(reshape(F,shape(F_lastInc))) call Utilities_updateIPcoords(reshape(F,shape(F_lastInc)))
call Utilities_constitutiveResponse(P,temp33_Real,C_volAvg,C_minMaxAvg, & call Utilities_constitutiveResponse(P,temp33_Real,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2
reshape(F,shape(F_lastInc)),0.0_pReal,math_I3) reshape(F,shape(F_lastInc)), & ! target F
0.0_pReal, & ! time increment
math_I3) ! no rotation of boundary condition
nullify(F) nullify(F)
nullify(F_tau) nullify(F_tau)
call DMDAVecRestoreArrayF90(da,solution_vec,xx_psc,ierr); CHKERRQ(ierr) ! write data back to PETSc call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr) ! write data back to PETSc
restartRead: if (restartInc > 1_pInt) then restartRead: if (restartInc > 0_pInt) then
if (iand(debug_level(debug_spectral),debug_spectralRestart)/= 0 .and. worldrank == 0_pInt) & if (iand(debug_level(debug_spectral),debug_spectralRestart)/= 0 .and. worldrank == 0_pInt) &
write(6,'(/,a,'//IO_intOut(restartInc-1_pInt)//',a)') & write(6,'(/,a,'//IO_intOut(restartInc)//',a)') &
'reading more values of increment', restartInc - 1_pInt, 'from file' 'reading more values of increment ', restartInc, ' from file'
flush(6) flush(6)
call IO_read_realFile(777,'C_volAvg',trim(getSolverJobName()),size(C_volAvg)) call IO_read_realFile(777,'C_volAvg',trim(getSolverJobName()),size(C_volAvg))
read (777,rec=1) C_volAvg read (777,rec=1) C_volAvg; close (777)
close (777)
call IO_read_realFile(777,'C_volAvgLastInc',trim(getSolverJobName()),size(C_volAvgLastInc)) call IO_read_realFile(777,'C_volAvgLastInc',trim(getSolverJobName()),size(C_volAvgLastInc))
read (777,rec=1) C_volAvgLastInc read (777,rec=1) C_volAvgLastInc; close (777)
close (777)
call IO_read_realFile(777,'C_ref',trim(getSolverJobName()),size(C_minMaxAvg)) call IO_read_realFile(777,'C_ref',trim(getSolverJobName()),size(C_minMaxAvg))
read (777,rec=1) C_minMaxAvg read (777,rec=1) C_minMaxAvg; close (777)
close (777)
endif restartRead endif restartRead
call Utilities_updateGamma(C_minMaxAvg,.True.) call Utilities_updateGamma(C_minMaxAvg,.true.)
C_scale = C_minMaxAvg C_scale = C_minMaxAvg
S_scale = math_invSym3333(C_minMaxAvg) S_scale = math_invSym3333(C_minMaxAvg)
@ -245,8 +244,7 @@ end subroutine Polarisation_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief solution for the Polarisation scheme with internal iterations !> @brief solution for the Polarisation scheme with internal iterations
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
type(tSolutionState) function & type(tSolutionState) function Polarisation_solution(incInfoIn,timeinc,timeinc_old,stress_BC,rotation_BC)
Polarisation_solution(incInfoIn,timeinc,timeinc_old,stress_BC,rotation_BC)
use IO, only: & use IO, only: &
IO_error IO_error
use numerics, only: & use numerics, only: &
@ -265,13 +263,13 @@ type(tSolutionState) function &
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! input data for solution ! input data for solution
real(pReal), intent(in) :: &
timeinc, & !< increment in time for current solution
timeinc_old !< increment in time of last increment
type(tBoundaryCondition), intent(in) :: &
stress_BC
character(len=*), intent(in) :: & character(len=*), intent(in) :: &
incInfoIn incInfoIn
real(pReal), intent(in) :: &
timeinc, & !< increment time for current solution
timeinc_old !< increment time of last successful increment
type(tBoundaryCondition), intent(in) :: &
stress_BC
real(pReal), dimension(3,3), intent(in) :: rotation_BC real(pReal), dimension(3,3), intent(in) :: rotation_BC
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -304,18 +302,17 @@ type(tSolutionState) function &
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! solve BVP ! solve BVP
call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr) call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr); CHKERRQ(ierr)
CHKERRQ(ierr)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! check convergence ! check convergence
call SNESGetConvergedReason(snes,reason,ierr) call SNESGetConvergedReason(snes,reason,ierr); CHKERRQ(ierr)
CHKERRQ(ierr)
Polarisation_solution%converged = reason > 0
Polarisation_solution%iterationsNeeded = totalIter
Polarisation_solution%termIll = terminallyIll Polarisation_solution%termIll = terminallyIll
terminallyIll = .false. terminallyIll = .false.
if (reason == -4) call IO_error(893_pInt) if (reason == -4) call IO_error(893_pInt) ! MPI error
if (reason < 1) Polarisation_solution%converged = .false.
Polarisation_solution%iterationsNeeded = totalIter
end function Polarisation_solution end function Polarisation_solution
@ -330,16 +327,19 @@ subroutine Polarisation_formResidual(in,x_scal,f_scal,dummy,ierr)
polarAlpha, & polarAlpha, &
polarBeta polarBeta
use mesh, only: & use mesh, only: &
grid3, & grid, &
grid grid3
use IO, only: & use IO, only: &
IO_intOut IO_intOut
use math, only: & use math, only: &
math_rotate_backward33, & math_rotate_backward33, &
math_transpose33, &
math_mul3333xx33, & math_mul3333xx33, &
math_invSym3333, & math_invSym3333, &
math_mul33x33 math_mul33x33
use debug, only: &
debug_level, &
debug_spectral, &
debug_spectralRotation
use spectral_utilities, only: & use spectral_utilities, only: &
wgt, & wgt, &
tensorField_real, & tensorField_real, &
@ -349,21 +349,13 @@ subroutine Polarisation_formResidual(in,x_scal,f_scal,dummy,ierr)
Utilities_constitutiveResponse, & Utilities_constitutiveResponse, &
Utilities_divergenceRMS, & Utilities_divergenceRMS, &
Utilities_curlRMS Utilities_curlRMS
use debug, only: &
debug_level, &
debug_spectral, &
debug_spectralRotation
use homogenization, only: & use homogenization, only: &
materialpoint_dPdF materialpoint_dPdF
use FEsolving, only: & use FEsolving, only: &
terminallyIll terminallyIll
implicit none implicit none
!-------------------------------------------------------------------------------------------------- DMDALocalInfo, dimension(DMDA_LOCAL_INFO_SIZE) :: in
! strange syntax in the next line because otherwise macros expand beyond 132 character limit
DMDALocalInfo, dimension(&
DMDA_LOCAL_INFO_SIZE) :: &
in
PetscScalar, & PetscScalar, &
target, dimension(3,3,2, XG_RANGE,YG_RANGE,ZG_RANGE), intent(in) :: x_scal target, dimension(3,3,2, XG_RANGE,YG_RANGE,ZG_RANGE), intent(in) :: x_scal
PetscScalar, & PetscScalar, &
@ -394,24 +386,24 @@ subroutine Polarisation_formResidual(in,x_scal,f_scal,dummy,ierr)
residual_F_tau => f_scal(1:3,1:3,2,& residual_F_tau => f_scal(1:3,1:3,2,&
X_RANGE, Y_RANGE, Z_RANGE) X_RANGE, Y_RANGE, Z_RANGE)
call SNESGetNumberFunctionEvals(snes,nfuncs,ierr); CHKERRQ(ierr)
call SNESGetIterationNumber(snes,PETScIter,ierr); CHKERRQ(ierr)
F_av = sum(sum(sum(F,dim=5),dim=4),dim=3) * wgt F_av = sum(sum(sum(F,dim=5),dim=4),dim=3) * wgt
call MPI_Allreduce(MPI_IN_PLACE,F_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,F_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
call SNESGetNumberFunctionEvals(snes,nfuncs,ierr); CHKERRQ(ierr)
call SNESGetIterationNumber(snes,PETScIter,ierr); CHKERRQ(ierr)
if (nfuncs == 0 .and. PETScIter == 0) totalIter = -1_pInt ! new increment if (nfuncs == 0 .and. PETScIter == 0) totalIter = -1_pInt ! new increment
newIteration: if(totalIter <= PETScIter) then
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! report begin of new iteration ! begin of new iteration
newIteration: if (totalIter <= PETScIter) then
totalIter = totalIter + 1_pInt totalIter = totalIter + 1_pInt
write(6,'(1x,a,3(a,'//IO_intOut(itmax)//'))') trim(incInfo), & write(6,'(1x,a,3(a,'//IO_intOut(itmax)//'))') &
' @ Iteration ', itmin, '≤',totalIter, '≤', itmax trim(incInfo), ' @ Iteration ', itmin, '≤',totalIter, '≤', itmax
if (iand(debug_level(debug_spectral),debug_spectralRotation) /= 0) & if (iand(debug_level(debug_spectral),debug_spectralRotation) /= 0) &
write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') ' deformation gradient aim (lab) =', & write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') &
math_transpose33(math_rotate_backward33(F_aim,params%rotation_BC)) ' deformation gradient aim (lab) =', transpose(math_rotate_backward33(F_aim,params%rotation_BC))
write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') ' deformation gradient aim =', & write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') &
math_transpose33(F_aim) ' deformation gradient aim =', transpose(F_aim)
flush(6) flush(6)
endif newIteration endif newIteration
@ -441,15 +433,13 @@ subroutine Polarisation_formResidual(in,x_scal,f_scal,dummy,ierr)
call Utilities_constitutiveResponse(residual_F,P_av,C_volAvg,C_minMaxAvg, & call Utilities_constitutiveResponse(residual_F,P_av,C_volAvg,C_minMaxAvg, &
F - residual_F_tau/polarBeta,params%timeinc,params%rotation_BC) F - residual_F_tau/polarBeta,params%timeinc,params%rotation_BC)
call MPI_Allreduce(MPI_IN_PLACE,terminallyIll,1,MPI_LOGICAL,MPI_LOR,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,terminallyIll,1,MPI_LOGICAL,MPI_LOR,PETSC_COMM_WORLD,ierr)
ForwardData = .False.
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! calculate divergence ! calculate divergence
tensorField_real = 0.0_pReal tensorField_real = 0.0_pReal
tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) = residual_F tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) = residual_F !< stress field in disguise
call utilities_FFTtensorForward() call utilities_FFTtensorForward()
err_div = Utilities_divergenceRMS() err_div = Utilities_divergenceRMS() !< root mean squared error in divergence of stress
call utilities_FFTtensorBackward()
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! constructing residual ! constructing residual
@ -469,8 +459,11 @@ subroutine Polarisation_formResidual(in,x_scal,f_scal,dummy,ierr)
tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) = F tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) = F
call utilities_FFTtensorForward() call utilities_FFTtensorForward()
err_curl = Utilities_curlRMS() err_curl = Utilities_curlRMS()
call utilities_FFTtensorBackward()
nullify(F)
nullify(F_tau)
nullify(residual_F)
nullify(residual_F_tau)
end subroutine Polarisation_formResidual end subroutine Polarisation_formResidual
@ -485,8 +478,8 @@ subroutine Polarisation_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,
err_div_tolAbs, & err_div_tolAbs, &
err_curl_tolRel, & err_curl_tolRel, &
err_curl_tolAbs, & err_curl_tolAbs, &
err_stress_tolAbs, & err_stress_tolRel, &
err_stress_tolRel err_stress_tolAbs
use math, only: & use math, only: &
math_mul3333xx33 math_mul3333xx33
use FEsolving, only: & use FEsolving, only: &
@ -505,24 +498,24 @@ subroutine Polarisation_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,
real(pReal) :: & real(pReal) :: &
curlTol, & curlTol, &
divTol, & divTol, &
BC_tol BCTol
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! stress BC handling ! stress BC handling
F_aim = F_aim - math_mul3333xx33(S, ((P_av - params%stress_BC))) ! S = 0.0 for no bc F_aim = F_aim - math_mul3333xx33(S, ((P_av - params%stress_BC))) ! S = 0.0 for no bc
err_BC = maxval(abs((-mask_stress+1.0_pReal)*math_mul3333xx33(C_scale,F_aim-F_av) + & err_BC = maxval(abs((1.0_pReal-mask_stress) * math_mul3333xx33(C_scale,F_aim-F_av) + &
mask_stress * (P_av-params%stress_BC))) ! mask = 0.0 for no bc mask_stress * (P_av-params%stress_BC))) ! mask = 0.0 for no bc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! error calculation ! error calculation
curlTol = max(maxval(abs(F_aim-math_I3))*err_curl_tolRel ,err_curl_tolAbs) curlTol = max(maxval(abs(F_aim-math_I3))*err_curl_tolRel ,err_curl_tolAbs)
divTol = max(maxval(abs(P_av)) *err_div_tolRel ,err_div_tolAbs) divTol = max(maxval(abs(P_av)) *err_div_tolRel ,err_div_tolAbs)
BC_tol = max(maxval(abs(P_av)) *err_stress_tolrel,err_stress_tolabs) BCTol = max(maxval(abs(P_av)) *err_stress_tolRel,err_stress_tolAbs)
converged: if ((totalIter >= itmin .and. & converged: if ((totalIter >= itmin .and. &
all([ err_div /divTol, & all([ err_div /divTol, &
err_curl/curlTol, & err_curl/curlTol, &
err_BC/BC_tol ] < 1.0_pReal)) & err_BC /BCTol ] < 1.0_pReal)) &
.or. terminallyIll) then .or. terminallyIll) then
reason = 1 reason = 1
elseif (totalIter >= itmax) then converged elseif (totalIter >= itmax) then converged
@ -534,12 +527,12 @@ subroutine Polarisation_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! report ! report
write(6,'(1/,a)') ' ... reporting .............................................................' write(6,'(1/,a)') ' ... reporting .............................................................'
write(6,'(/,a,f12.2,a,es8.2,a,es9.2,a)') ' error curl = ', & write(6,'(/,a,f12.2,a,es8.2,a,es9.2,a)') ' error divergence = ', &
err_curl/curlTol,' (',err_curl,' -, tol =',curlTol,')'
write(6,' (a,f12.2,a,es8.2,a,es9.2,a)') ' error divergence = ', &
err_div/divTol, ' (',err_div, ' / m, tol = ',divTol,')' err_div/divTol, ' (',err_div, ' / m, tol = ',divTol,')'
write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error curl = ', &
err_curl/curlTol,' (',err_curl,' -, tol = ',curlTol,')'
write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error BC = ', & write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error BC = ', &
err_BC/BC_tol, ' (',err_BC, ' Pa, tol =',BC_tol,')' err_BC/BCTol, ' (',err_BC, ' Pa, tol = ',BCTol,')'
write(6,'(/,a)') ' ===========================================================================' write(6,'(/,a)') ' ==========================================================================='
flush(6) flush(6)
@ -547,18 +540,23 @@ end subroutine Polarisation_converged
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief forwarding routine !> @brief forwarding routine
!> @details find new boundary conditions and best F estimate for end of current timestep
!> possibly writing restart information, triggering of state increment in DAMASK, and updating of IPcoordinates
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine Polarisation_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation_BC,stress_BC,rotation_BC) subroutine Polarisation_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation_BC,stress_BC,rotation_BC)
use math, only: & use math, only: &
math_mul33x33, & math_mul33x33, &
math_mul3333xx33, & math_mul3333xx33, &
math_transpose33, &
math_rotate_backward33 math_rotate_backward33
use numerics, only: & use numerics, only: &
worldrank worldrank
use homogenization, only: &
materialpoint_F0
use mesh, only: & use mesh, only: &
grid3, & grid, &
grid grid3
use CPFEM2, only: &
CPFEM_age
use spectral_utilities, only: & use spectral_utilities, only: &
Utilities_calculateRate, & Utilities_calculateRate, &
Utilities_forwardField, & Utilities_forwardField, &
@ -571,6 +569,8 @@ subroutine Polarisation_forward(guess,timeinc,timeinc_old,loadCaseTime,deformati
restartWrite restartWrite
implicit none implicit none
logical, intent(in) :: &
guess
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
timeinc_old, & timeinc_old, &
timeinc, & timeinc, &
@ -578,111 +578,108 @@ subroutine Polarisation_forward(guess,timeinc,timeinc_old,loadCaseTime,deformati
type(tBoundaryCondition), intent(in) :: & type(tBoundaryCondition), intent(in) :: &
stress_BC, & stress_BC, &
deformation_BC deformation_BC
real(pReal), dimension(3,3), intent(in) :: rotation_BC real(pReal), dimension(3,3), intent(in) ::&
logical, intent(in) :: & rotation_BC
guess
PetscErrorCode :: ierr PetscErrorCode :: ierr
PetscScalar, dimension(:,:,:,:), pointer :: xx_psc, F, F_tau PetscScalar, dimension(:,:,:,:), pointer :: FandF_tau, F, F_tau
integer(pInt) :: i, j, k integer(pInt) :: i, j, k
real(pReal), dimension(3,3) :: F_lambda33 real(pReal), dimension(3,3) :: F_lambda33
character(len=1024) :: rankStr character(len=32) :: rankStr
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! update coordinates and rate and forward last inc ! update coordinates and rate and forward last inc
call DMDAVecGetArrayF90(da,solution_vec,xx_psc,ierr) call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr)
F => xx_psc(0:8,:,:,:) F => FandF_tau( 0: 8,:,:,:)
F_tau => xx_psc(9:17,:,:,:) F_tau => FandF_tau( 9:17,:,:,:)
if (restartWrite) then
write(6,'(/,a)') ' writing converged results for restart'
flush(6)
write(rankStr,'(a1,i0)')'_',worldrank
call IO_write_jobRealFile(777,'F'//trim(rankStr),size(F)) ! writing deformation gradient field to file
write (777,rec=1) F
close (777)
call IO_write_jobRealFile(777,'F_lastInc'//trim(rankStr),size(F_lastInc)) ! writing F_lastInc field to file
write (777,rec=1) F_lastInc
close (777)
call IO_write_jobRealFile(777,'F_tau'//trim(rankStr),size(F_tau)) ! writing deformation gradient field to file
write (777,rec=1) F_tau
close (777)
call IO_write_jobRealFile(777,'F_tau_lastInc'//trim(rankStr),size(F_tau_lastInc)) ! writing F_lastInc field to file
write (777,rec=1) F_tau_lastInc
close (777)
if (worldrank == 0_pInt) then
call IO_write_jobRealFile(777,'F_aim',size(F_aim))
write (777,rec=1) F_aim
close(777)
call IO_write_jobRealFile(777,'F_aim_lastInc',size(F_aim_lastInc))
write (777,rec=1) F_aim_lastInc
close (777)
call IO_write_jobRealFile(777,'F_aimDot',size(F_aimDot))
write (777,rec=1) F_aimDot
close(777)
call IO_write_jobRealFile(777,'C_volAvg',size(C_volAvg))
write (777,rec=1) C_volAvg
close(777)
call IO_write_jobRealFile(777,'C_volAvgLastInc',size(C_volAvgLastInc))
write (777,rec=1) C_volAvgLastInc
close(777)
endif
endif
call utilities_updateIPcoords(F)
if (cutBack) then if (cutBack) then
F_aim = F_aim_lastInc C_volAvg = C_volAvgLastInc ! QUESTION: where is this required?
F_tau= reshape(F_tau_lastInc,[9,grid(1),grid(2),grid3]) C_minMaxAvg = C_minMaxAvgLastInc ! QUESTION: where is this required?
F = reshape(F_lastInc, [9,grid(1),grid(2),grid3])
C_volAvg = C_volAvgLastInc
else else
ForwardData = .True. !--------------------------------------------------------------------------------------------------
! restart information for spectral solver
if (restartWrite) then ! QUESTION: where is this logical properly set?
write(6,'(/,a)') ' writing converged results for restart'
flush(6)
if (worldrank == 0_pInt) then
call IO_write_jobRealFile(777,'C_volAvg',size(C_volAvg))
write (777,rec=1) C_volAvg; close(777)
call IO_write_jobRealFile(777,'C_volAvgLastInc',size(C_volAvgLastInc))
write (777,rec=1) C_volAvgLastInc; close(777)
call IO_write_jobRealFile(777,'F_aimDot',size(F_aimDot))
write (777,rec=1) F_aimDot; close(777)
endif
write(rankStr,'(a1,i0)')'_',worldrank
call IO_write_jobRealFile(777,'F'//trim(rankStr),size(F)) ! writing deformation gradient field to file
write (777,rec=1) F; close (777)
call IO_write_jobRealFile(777,'F_lastInc'//trim(rankStr),size(F_lastInc)) ! writing F_lastInc field to file
write (777,rec=1) F_lastInc; close (777)
call IO_write_jobRealFile(777,'F_tau'//trim(rankStr),size(F_tau)) ! writing deformation gradient field to file
write (777,rec=1) F_tau; close (777)
call IO_write_jobRealFile(777,'F_tau_lastInc'//trim(rankStr),size(F_tau_lastInc)) ! writing F_tau_lastInc field to file
write (777,rec=1) F_tau_lastInc; close (777)
endif
call CPFEM_age() ! age state and kinematics
call utilities_updateIPcoords(F)
C_volAvgLastInc = C_volAvg C_volAvgLastInc = C_volAvg
C_minMaxAvgLastInc = C_minMaxAvg
F_aimDot = merge(stress_BC%maskFloat*(F_aim-F_aim_lastInc)/timeinc_old, 0.0_pReal, guess)
F_aim_lastInc = F_aim
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! calculate rate for aim ! calculate rate for aim
if (deformation_BC%myType=='l') then ! calculate f_aimDot from given L and current F if (deformation_BC%myType=='l') then ! calculate F_aimDot from given L and current F
f_aimDot = deformation_BC%maskFloat * math_mul33x33(deformation_BC%values, F_aim) F_aimDot = &
elseif(deformation_BC%myType=='fdot') then ! f_aimDot is prescribed F_aimDot + deformation_BC%maskFloat * math_mul33x33(deformation_BC%values, F_aim_lastInc)
f_aimDot = deformation_BC%maskFloat * deformation_BC%values elseif(deformation_BC%myType=='fdot') then ! F_aimDot is prescribed
F_aimDot = &
F_aimDot + deformation_BC%maskFloat * deformation_BC%values
elseif (deformation_BC%myType=='f') then ! aim at end of load case is prescribed elseif (deformation_BC%myType=='f') then ! aim at end of load case is prescribed
f_aimDot = deformation_BC%maskFloat * (deformation_BC%values -F_aim)/loadCaseTime F_aimDot = &
F_aimDot + deformation_BC%maskFloat * (deformation_BC%values - F_aim_lastInc)/loadCaseTime
endif endif
if (guess) f_aimDot = f_aimDot + stress_BC%maskFloat * (F_aim - F_aim_lastInc)/timeinc_old
F_aim_lastInc = F_aim
!--------------------------------------------------------------------------------------------------
! update coordinates and rate and forward last inc
call utilities_updateIPcoords(F)
Fdot = Utilities_calculateRate(guess, & Fdot = Utilities_calculateRate(guess, &
F_lastInc,reshape(F,[3,3,grid(1),grid(2),grid3]),timeinc_old, & F_lastInc,reshape(F,[3,3,grid(1),grid(2),grid3]),timeinc_old, &
math_rotate_backward33( f_aimDot,rotation_BC)) math_rotate_backward33(F_aimDot,rotation_BC))
F_tauDot = Utilities_calculateRate(guess, & F_tauDot = Utilities_calculateRate(guess, &
F_tau_lastInc,reshape(F_tau,[3,3,grid(1),grid(2),grid3]), timeinc_old, & F_tau_lastInc,reshape(F_tau,[3,3,grid(1),grid(2),grid3]), timeinc_old, &
math_rotate_backward33(2.0_pReal*f_aimDot,rotation_BC)) math_rotate_backward33(F_aimDot,rotation_BC))
F_lastInc = reshape(F, [3,3,grid(1),grid(2),grid3]) F_lastInc = reshape(F, [3,3,grid(1),grid(2),grid3]) ! winding F forward
F_tau_lastInc = reshape(F_tau,[3,3,grid(1),grid(2),grid3]) F_tau_lastInc = reshape(F_tau, [3,3,grid(1),grid(2),grid3]) ! winding F_tau forward
materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent
endif endif
F_aim = F_aim + f_aimDot * timeinc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! update local deformation gradient ! update average and local deformation gradients
F = reshape(Utilities_forwardField(timeinc,F_lastInc,Fdot, & ! ensure that it matches rotated F_aim F_aim = F_aim_lastInc + F_aimDot * timeinc
F = reshape(Utilities_forwardField(timeinc,F_lastInc,Fdot, & ! estimate of F at end of time+timeinc that matches rotated F_aim on average
math_rotate_backward33(F_aim,rotation_BC)),& math_rotate_backward33(F_aim,rotation_BC)),&
[9,grid(1),grid(2),grid3]) [9,grid(1),grid(2),grid3])
F_tau = reshape(Utilities_forwardField(timeinc,F_tau_lastInc,F_taudot), & ! does not have any average value as boundary condition if (guess) then
[9,grid(1),grid(2),grid3]) F_tau = reshape(Utilities_forwardField(timeinc,F_tau_lastInc,F_taudot), &
if (.not. guess) then ! large strain forwarding [9,grid(1),grid(2),grid3]) ! does not have any average value as boundary condition
else
do k = 1_pInt, grid3; do j = 1_pInt, grid(2); do i = 1_pInt, grid(1) do k = 1_pInt, grid3; do j = 1_pInt, grid(2); do i = 1_pInt, grid(1)
F_lambda33 = reshape(F_tau(1:9,i,j,k)-F(1:9,i,j,k),[3,3]) F_lambda33 = reshape(F_tau(1:9,i,j,k)-F(1:9,i,j,k),[3,3])
F_lambda33 = math_mul3333xx33(S_scale,math_mul33x33(F_lambda33, & F_lambda33 = math_mul3333xx33(S_scale,math_mul33x33(F_lambda33, &
math_mul3333xx33(C_scale,& math_mul3333xx33(C_scale,&
math_mul33x33(math_transpose33(F_lambda33),& math_mul33x33(transpose(F_lambda33),&
F_lambda33)-math_I3))*0.5_pReal)& F_lambda33)-math_I3))*0.5_pReal)&
+ math_I3 + math_I3
F_tau(1:9,i,j,k) = reshape(F_lambda33,[9])+F(1:9,i,j,k) F_tau(1:9,i,j,k) = reshape(F_lambda33,[9])+F(1:9,i,j,k)
enddo; enddo; enddo enddo; enddo; enddo
endif endif
call DMDAVecRestoreArrayF90(da,solution_vec,xx_psc,ierr); CHKERRQ(ierr)
nullify(F)
nullify(F_tau)
call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr)
end subroutine Polarisation_forward end subroutine Polarisation_forward

View File

@ -60,7 +60,7 @@ contains
!> @brief allocates all neccessary fields and fills them with data, potentially from restart info !> @brief allocates all neccessary fields and fills them with data, potentially from restart info
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine spectral_thermal_init subroutine spectral_thermal_init
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -70,7 +70,6 @@ module spectral_utilities
! derived types ! derived types
type, public :: tSolutionState !< return type of solution from spectral solver variants type, public :: tSolutionState !< return type of solution from spectral solver variants
logical :: converged = .true. logical :: converged = .true.
logical :: regrid = .false.
logical :: stagConverged = .true. logical :: stagConverged = .true.
logical :: termIll = .false. logical :: termIll = .false.
integer(pInt) :: iterationsNeeded = 0_pInt integer(pInt) :: iterationsNeeded = 0_pInt
@ -145,8 +144,7 @@ module spectral_utilities
FIELD_UNDEFINED_ID, & FIELD_UNDEFINED_ID, &
FIELD_MECH_ID, & FIELD_MECH_ID, &
FIELD_THERMAL_ID, & FIELD_THERMAL_ID, &
FIELD_DAMAGE_ID, & FIELD_DAMAGE_ID
utilities_calcPlasticity
private :: & private :: &
utilities_getFreqDerivative utilities_getFreqDerivative
@ -154,14 +152,14 @@ contains
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief allocates all neccessary fields, sets debug flags, create plans for FFTW !> @brief allocates all neccessary fields, sets debug flags, create plans for FFTW
!> @details Sets the debug levels for general, divergence, restart and FFTW from the biwise coding !> @details Sets the debug levels for general, divergence, restart, and FFTW from the bitwise coding
!> provided by the debug module to logicals. !> provided by the debug module to logicals.
!> Allocates all fields used by FFTW and create the corresponding plans depending on the debug !> Allocate all fields used by FFTW and create the corresponding plans depending on the debug
!> level chosen. !> level chosen.
!> Initializes FFTW. !> Initializes FFTW.
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine utilities_init() subroutine utilities_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options
@ -378,10 +376,10 @@ end subroutine utilities_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief updates references stiffness and potentially precalculated gamma operator !> @brief updates reference stiffness and potentially precalculated gamma operator
!> @details Sets the current reference stiffness to the stiffness given as an argument. !> @details Sets the current reference stiffness to the stiffness given as an argument.
!> If the gamma operator is precalculated, it is calculated with this stiffness. !> If the gamma operator is precalculated, it is calculated with this stiffness.
!> In case of a on-the-fly calculation, only the reference stiffness is updated. !> In case of an on-the-fly calculation, only the reference stiffness is updated.
!> Also writes out the current reference stiffness for restart. !> Also writes out the current reference stiffness for restart.
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine utilities_updateGamma(C,saveReference) subroutine utilities_updateGamma(C,saveReference)
@ -414,8 +412,7 @@ subroutine utilities_updateGamma(C,saveReference)
write(6,'(/,a)') ' writing reference stiffness to file' write(6,'(/,a)') ' writing reference stiffness to file'
flush(6) flush(6)
call IO_write_jobRealFile(777,'C_ref',size(C_ref)) call IO_write_jobRealFile(777,'C_ref',size(C_ref))
write (777,rec=1) C_ref write (777,rec=1) C_ref; close(777)
close(777)
endif endif
endif endif
@ -784,7 +781,7 @@ function utilities_maskedCompliance(rot_BC,mask_stress,C)
if(debugGeneral) then if(debugGeneral) then
write(6,'(/,a)') ' ... updating masked compliance ............................................' write(6,'(/,a)') ' ... updating masked compliance ............................................'
write(6,'(/,a,/,9(9(2x,f12.7,1x)/))',advance='no') ' Stiffness C (load) / GPa =',& write(6,'(/,a,/,9(9(2x,f12.7,1x)/))',advance='no') ' Stiffness C (load) / GPa =',&
transpose(temp99_Real)/1.e9_pReal transpose(temp99_Real)*1.0e-9_pReal
flush(6) flush(6)
endif endif
k = 0_pInt ! calculate reduced stiffness k = 0_pInt ! calculate reduced stiffness
@ -823,7 +820,6 @@ function utilities_maskedCompliance(rot_BC,mask_stress,C)
.or. (m/=n .and. abs(sTimesC(m,n)) > 1.0e-12_pReal) ! off-diagonal elements of S*C should be 0 .or. (m/=n .and. abs(sTimesC(m,n)) > 1.0e-12_pReal) ! off-diagonal elements of S*C should be 0
enddo enddo
enddo enddo
if (debugGeneral .or. errmatinv) then if (debugGeneral .or. errmatinv) then
write(formatString, '(i2)') size_reduced write(formatString, '(i2)') size_reduced
formatString = '(/,a,/,'//trim(formatString)//'('//trim(formatString)//'(2x,es9.2,1x)/))' formatString = '(/,a,/,'//trim(formatString)//'('//trim(formatString)//'(2x,es9.2,1x)/))'
@ -838,13 +834,11 @@ function utilities_maskedCompliance(rot_BC,mask_stress,C)
else else
temp99_real = 0.0_pReal temp99_real = 0.0_pReal
endif endif
if(debugGeneral) then if(debugGeneral) then
write(6,'(/,a,/,9(9(2x,f10.5,1x)/),/)',advance='no') & write(6,'(/,a,/,9(9(2x,f10.5,1x)/),/)',advance='no') &
' Masked Compliance (load) / GPa =', transpose(temp99_Real*1.e-9_pReal) ' Masked Compliance (load) * GPa =', transpose(temp99_Real)*1.0e9_pReal
flush(6) flush(6)
endif endif
utilities_maskedCompliance = math_Plain99to3333(temp99_Real) utilities_maskedCompliance = math_Plain99to3333(temp99_Real)
end function utilities_maskedCompliance end function utilities_maskedCompliance
@ -940,7 +934,6 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,&
debug_reset, & debug_reset, &
debug_info debug_info
use math, only: & use math, only: &
math_transpose33, &
math_rotate_forward33, & math_rotate_forward33, &
math_det33 math_det33
use mesh, only: & use mesh, only: &
@ -999,10 +992,10 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,&
call MPI_Allreduce(MPI_IN_PLACE,P_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,P_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
if (debugRotation) & if (debugRotation) &
write(6,'(/,a,/,3(3(2x,f12.4,1x)/))',advance='no') ' Piola--Kirchhoff stress (lab) / MPa =',& write(6,'(/,a,/,3(3(2x,f12.4,1x)/))',advance='no') ' Piola--Kirchhoff stress (lab) / MPa =',&
math_transpose33(P_av)*1.e-6_pReal transpose(P_av)*1.e-6_pReal
P_av = math_rotate_forward33(P_av,rotation_BC) P_av = math_rotate_forward33(P_av,rotation_BC)
write(6,'(/,a,/,3(3(2x,f12.4,1x)/))',advance='no') ' Piola--Kirchhoff stress / MPa =',& write(6,'(/,a,/,3(3(2x,f12.4,1x)/))',advance='no') ' Piola--Kirchhoff stress / MPa =',&
math_transpose33(P_av)*1.e-6_pReal transpose(P_av)*1.e-6_pReal
flush(6) flush(6)
max_dPdF = 0.0_pReal max_dPdF = 0.0_pReal
@ -1034,125 +1027,6 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,&
end subroutine utilities_constitutiveResponse end subroutine utilities_constitutiveResponse
!--------------------------------------------------------------------------------------------------
!> @brief calculates yield stress, plastic strain, total strain and their equivalent values
!--------------------------------------------------------------------------------------------------
subroutine utilities_calcPlasticity(yieldStress, plasticStrain, eqStress, eqTotalStrain, &
eqPlasticStrain, plasticWork, rotation_BC)
use crystallite, only: &
crystallite_Fe, &
crystallite_P, &
crystallite_subF
use material, only: &
homogenization_maxNgrains
use mesh, only: &
mesh_maxNips,&
mesh_NcpElems
use math, only: &
math_det33, &
math_inv33, &
math_mul33x33, &
math_trace33, &
math_transpose33, &
math_equivStrain33, &
math_equivStress33, &
math_rotate_forward33, &
math_identity2nd, &
math_crossproduct, &
math_eigenvectorBasisSym, &
math_eigenvectorBasisSym33, &
math_eigenvectorBasisSym33_log, &
math_eigenValuesVectorsSym33
implicit none
real(pReal), intent(inout) :: eqStress, eqPlasticStrain, plasticWork
real(pReal), intent(out) :: eqTotalStrain
real(pReal), dimension(3,3),intent(out) :: yieldStress, plasticStrain
real(pReal), intent(in), dimension(3,3) :: rotation_BC !< rotation of load frame
real(pReal), dimension(3,3) :: cauchy, P_av, F_av, Ve_av !< average
real(pReal), dimension(3) :: Values, S
real(pReal), dimension(3,3) :: Vectors, diag
real(pReal), dimension(3,3) :: &
Vp, F_temp, U, VT, R, V, V_total
real(pReal), dimension(3,3,homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems) :: &
Be, Ve, Fe
real(pReal), dimension(15) :: WORK !< previous deformation gradient
integer(pInt) :: INFO, i, j, k, l, ierr
real(pReal) :: wgtm
real(pReal) :: eqStressOld, eqPlasticStrainOld, plasticWorkOld
external :: dgesvd
eqStressOld = eqStress
eqPlasticStrainOld = eqPlasticStrain
plasticWorkOld = plasticWork
wgtm = 1.0_pReal/real(mesh_NcpElems*mesh_maxNips*homogenization_maxNgrains,pReal)
diag = 0.0_pReal
P_av = sum(sum(sum(crystallite_P,dim=5),dim=4),dim=3) * wgtm
call MPI_Allreduce(MPI_IN_PLACE,P_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
P_av = math_rotate_forward33(P_av,rotation_BC)
F_av = sum(sum(sum(crystallite_subF,dim=5),dim=4),dim=3) * wgtm
call MPI_Allreduce(MPI_IN_PLACE,F_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
F_av = math_rotate_forward33(F_av,rotation_BC)
cauchy = 1.0_pReal/math_det33(F_av)*math_mul33x33(P_av,transpose(F_av))
yieldStress = cauchy
eqStress = math_equivStress33(cauchy)
F_temp = F_av
call dgesvd ('A', 'A', 3, 3, F_temp, 3, S, U, 3, VT, 3, WORK, 15, INFO) ! singular value decomposition
R = math_mul33x33(U, VT) ! rotation of polar decomposition
V = math_mul33x33(F_av,math_inv33(R))
call math_eigenValuesVectorsSym33(V,Values,Vectors)
do l = 1_pInt, 3_pInt
if (Values(l) < 0.0_pReal) then
Values(l) = -Values(l)
Vectors(1:3, l) = -Vectors(1:3, l)
endif
Values(l) = log(Values(l))
diag(l,l) = Values(l)
enddo
if (dot_product(Vectors(1:3,1),Vectors(1:3,2)) /= 0) then
Vectors(1:3,2) = math_crossproduct(Vectors(1:3,3), Vectors(1:3,1))
Vectors(1:3,2) = Vectors(1:3,2)/sqrt(dot_product(Vectors(1:3,2),Vectors(1:3,2)))
endif
if (dot_product(Vectors(1:3,2),Vectors(1:3,3)) /= 0) then
Vectors(1:3,3) = math_crossproduct(Vectors(1:3,1), Vectors(1:3,2))
Vectors(1:3,3) = Vectors(1:3,3)/sqrt(dot_product(Vectors(1:3,3),Vectors(1:3,3)))
endif
if (dot_product(Vectors(1:3,3),Vectors(1:3,1)) /= 0) then
Vectors(1:3,1) = math_crossproduct(Vectors(1:3,2), Vectors(1:3,3))
Vectors(1:3,1) = Vectors(1:3,1)/sqrt(dot_product(Vectors(1:3,1),Vectors(1:3,1)))
endif
V_total = REAL(math_mul33x33(Vectors, math_mul33x33(diag, transpose(Vectors))))
eqTotalStrain = math_equivStrain33(V_total)
do k = 1_pInt, mesh_NcpElems; do j = 1_pInt, mesh_maxNips; do i = 1_pInt,homogenization_maxNgrains
Fe(1:3,1:3,i,j,k) = crystallite_Fe(1:3,1:3,i,j,k)
Fe(1:3,1:3,i,j,k) = math_rotate_forward33(Fe(1:3,1:3,i,j,k),rotation_BC)
Be(1:3,1:3,i,j,k) = math_mul33x33(Fe(1:3,1:3,i,j,k),math_transpose33(Fe(1:3,1:3,i,j,k))) ! elastic part of left CauchyGreen deformation tensor
Ve(1:3,1:3,i,j,k) = math_eigenvectorBasisSym33_log(Be(1:3,1:3,i,j,k))
enddo; enddo; enddo
Ve_av = sum(sum(sum(Ve,dim=5),dim=4),dim=3) * wgtm
call MPI_Allreduce(MPI_IN_PLACE,Ve_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
Vp = V_total - Ve_av
eqPlasticStrain = math_equivStrain33(Vp)
plasticStrain = Vp
plasticWork = plasticWorkOld + 0.5*(eqStressOld + eqStress) * (eqPlasticStrain - eqPlasticStrainOld)
end subroutine utilities_calcPlasticity
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief calculates forward rate, either guessing or just add delta/timeinc !> @brief calculates forward rate, either guessing or just add delta/timeinc

View File

@ -46,7 +46,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine thermal_adiabatic_init(fileUnit) subroutine thermal_adiabatic_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -47,7 +47,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine thermal_conduction_init(fileUnit) subroutine thermal_conduction_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -16,7 +16,7 @@ contains
!> @brief allocates all neccessary fields, reads information from material configuration file !> @brief allocates all neccessary fields, reads information from material configuration file
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine thermal_isothermal_init() subroutine thermal_isothermal_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -61,7 +61,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine vacancyflux_cahnhilliard_init(fileUnit) subroutine vacancyflux_cahnhilliard_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -44,7 +44,7 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks !> @details reads in material parameters, allocates arrays, and does sanity checks
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine vacancyflux_isochempot_init(fileUnit) subroutine vacancyflux_isochempot_init(fileUnit)
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options

View File

@ -16,7 +16,7 @@ contains
!> @brief allocates all neccessary fields, reads information from material configuration file !> @brief allocates all neccessary fields, reads information from material configuration file
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine vacancyflux_isoconc_init() subroutine vacancyflux_isoconc_init()
#ifdef __GFORTRAN__ #if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800
use, intrinsic :: iso_fortran_env, only: & use, intrinsic :: iso_fortran_env, only: &
compiler_version, & compiler_version, &
compiler_options compiler_options