diff --git a/CMakeLists.txt b/CMakeLists.txt index 533b14e85..a6a133fe4 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,6 +1,12 @@ cmake_minimum_required (VERSION 3.10.0) include (FindPkgConfig REQUIRED) +if (DEFINED ENV{PETSC_DIR}) + message ("PETSC_DIR:\n$ENV{PETSC_DIR}\n") +else () + message (FATAL_ERROR "PETSc location (PETSC_DIR) is not defined") +endif () + # Dummy project to determine compiler names and version project (Prerequisites LANGUAGES) set(ENV{PKG_CONFIG_PATH} "$ENV{PETSC_DIR}/$ENV{PETSC_ARCH}/lib/pkgconfig") @@ -8,81 +14,8 @@ pkg_check_modules (PETSC REQUIRED PETSc>=3.12.0 PETSc<3.16.0) pkg_get_variable (CMAKE_Fortran_COMPILER PETSc fcompiler) pkg_get_variable (CMAKE_C_COMPILER PETSc ccompiler) -find_program (CAT_EXECUTABLE NAMES cat) -execute_process (COMMAND ${CAT_EXECUTABLE} ${PROJECT_SOURCE_DIR}/VERSION - RESULT_VARIABLE DAMASK_VERSION_RETURN - OUTPUT_VARIABLE DAMASK_VERSION - OUTPUT_STRIP_TRAILING_WHITESPACE) - -#--------------------------------------------------------------------------------------- -# Find PETSc from system environment -set(PETSC_DIR $ENV{PETSC_DIR}) -if (PETSC_DIR STREQUAL "") - message (FATAL_ERROR "PETSc location (PETSC_DIR) is not defined") -endif () - -set (petsc_conf_variables "${PETSC_DIR}/lib/petsc/conf/variables") -set (petsc_conf_rules "${PETSC_DIR}/lib/petsc/conf/rules" ) - -# Use existing variables from PETSc -# https://github.com/jedbrown/cmake-modules/blob/master/FindPETSc.cmake - -# Generate a temporary makefile to probe the PETSc configuration -# This file will be deleted once the settings from PETSc are parsed into CMake -exec_program (mktemp ARGS -d OUTPUT_VARIABLE TEMPDIR) -set (petsc_config_makefile "${TEMPDIR}/Makefile.petsc") -file (WRITE -"${petsc_config_makefile}" -"## This file was auto generated by CMake -# PETSC_DIR = ${PETSC_DIR} -SHELL = /bin/sh -include ${petsc_conf_rules} -include ${petsc_conf_variables} -INCLUDE_DIRS := \${PETSC_FC_INCLUDES} -LIBRARIES := \${PETSC_WITH_EXTERNAL_LIB} -includes: -\t@echo \${INCLUDE_DIRS} -extlibs: -\t@echo \${LIBRARIES} -") - -# CMake will execute each target in the ${petsc_config_makefile} -# to acquire corresponding PETSc Variables. -find_program (MAKE_EXECUTABLE NAMES gmake make) -# Find the PETSc includes directory settings -execute_process (COMMAND ${MAKE_EXECUTABLE} --no-print-directory -f ${petsc_config_makefile} "includes" - RESULT_VARIABLE PETSC_INCLUDES_RETURN - OUTPUT_VARIABLE petsc_includes - OUTPUT_STRIP_TRAILING_WHITESPACE) -# Find the PETSc external linking directory settings -execute_process (COMMAND ${MAKE_EXECUTABLE} --no-print-directory -f ${petsc_config_makefile} "extlibs" - RESULT_VARIABLE PETSC_EXTERNAL_LIB_RETURN - OUTPUT_VARIABLE petsc_external_lib - OUTPUT_STRIP_TRAILING_WHITESPACE) -# Remove temporary makefile, no need to keep it anymore. -file (REMOVE_RECURSE ${TEMPDIR}) - -# Remove duplicate compiler and linker flags -string (REGEX MATCHALL "-I([^\" ]+)" TMP_LIST "${petsc_includes}") -list (REMOVE_DUPLICATES TMP_LIST) -foreach (dir ${TMP_LIST}) - set (PETSC_INCLUDES "${PETSC_INCLUDES} ${dir}") -endforeach (dir) -string (REGEX MATCHALL "-[lLW]([^\" ]+)" TMP_LIST "${petsc_external_lib}") -list (REMOVE_DUPLICATES TMP_LIST) -foreach (exlib ${TMP_LIST}) - set (PETSC_EXTERNAL_LIB "${PETSC_EXTERNAL_LIB} ${exlib}") -endforeach (exlib) - -message ("Found PETSC_DIR:\n${PETSC_DIR}\n" ) -message ("Found PETSC_INCLUDES:\n${PETSC_INCLUDES}\n" ) -message ("Found PETSC_EXTERNAL_LIB:\n${PETSC_EXTERNAL_LIB}\n") - -#--------------------------------------------------------------------------------------- -# Now start to care about DAMASK - -# DAMASK solver defines project to build -string(TOLOWER ${DAMASK_SOLVER} DAMASK_SOLVER) +# Solver determines name of project +string(TOLOWER "${DAMASK_SOLVER}" DAMASK_SOLVER) if (DAMASK_SOLVER STREQUAL "grid") project (damask-grid HOMEPAGE_URL https://damask.mpie.de LANGUAGES Fortran C) add_definitions (-DGrid) @@ -92,10 +25,13 @@ elseif (DAMASK_SOLVER STREQUAL "mesh") else () message (FATAL_ERROR "Build target (DAMASK_SOLVER) is not defined") endif () -add_definitions (-DDAMASKVERSION="${DAMASK_VERSION}") -add_definitions (-DPETSc) -message ("\nBuilding ${CMAKE_PROJECT_NAME}\n") +file (STRINGS ${PROJECT_SOURCE_DIR}/VERSION DAMASK_VERSION) + +message ("\nBuilding ${CMAKE_PROJECT_NAME} ${DAMASK_VERSION}\n") + +add_definitions (-DPETSc) +add_definitions (-DDAMASKVERSION="${DAMASK_VERSION}") if (CMAKE_BUILD_TYPE STREQUAL "") set (CMAKE_BUILD_TYPE "RELEASE") @@ -139,12 +75,22 @@ if (CMAKE_Fortran_COMPILER_ID STREQUAL "Intel") include (Compiler-Intel) elseif(CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") include (Compiler-GNU) -elseif(CMAKE_Fortran_COMPILER_ID STREQUAL "PGI") - include (Compiler-PGI) else () message (FATAL_ERROR "Compiler type (CMAKE_Fortran_COMPILER_ID) not recognized") endif () +file (STRINGS "$ENV{PETSC_DIR}/$ENV{PETSC_ARCH}/lib/petsc/conf/petscvariables" PETSC_EXTERNAL_LIB REGEX "PETSC_WITH_EXTERNAL_LIB = .*$?") +string (REGEX MATCHALL "-[lLW]([^\" ]+)" PETSC_EXTERNAL_LIB "${PETSC_EXTERNAL_LIB}") +list (REMOVE_DUPLICATES PETSC_EXTERNAL_LIB) +string (REPLACE ";" " " PETSC_EXTERNAL_LIB "${PETSC_EXTERNAL_LIB}") +message ("PETSC_EXTERNAL_LIB:\n${PETSC_EXTERNAL_LIB}\n") + +file (STRINGS "$ENV{PETSC_DIR}/$ENV{PETSC_ARCH}/lib/petsc/conf/petscvariables" PETSC_INCLUDES REGEX "PETSC_FC_INCLUDES = .*$?") +string (REGEX MATCHALL "-I([^\" ]+)" PETSC_INCLUDES "${PETSC_INCLUDES}") +list (REMOVE_DUPLICATES PETSC_INCLUDES) +string (REPLACE ";" " " PETSC_INCLUDES "${PETSC_INCLUDES}") +message ("PETSC_INCLUDES:\n${PETSC_INCLUDES}\n") + set (CMAKE_Fortran_FLAGS_${CMAKE_BUILD_TYPE} "${BUILDCMD_PRE} ${OPENMP_FLAGS} ${STANDARD_CHECK} ${OPTIMIZATION_FLAGS} ${COMPILE_FLAGS} ${PRECISION_FLAGS}") set (CMAKE_Fortran_LINK_EXECUTABLE "${BUILDCMD_PRE} ${CMAKE_Fortran_COMPILER} ${OPENMP_FLAGS} ${OPTIMIZATION_FLAGS} ${LINKER_FLAGS}") diff --git a/PRIVATE b/PRIVATE index 95f7faea9..0744cf7f9 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 95f7faea920dd6956884e4a55f72e5d5b1ffcdc8 +Subproject commit 0744cf7f93dbd06423baaae97a67959f11c3e6a5 diff --git a/VERSION b/VERSION index eb33473c3..7567c1527 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v3.0.0-alpha3-51-g38ca85b58 +v3.0.0-alpha3-137-g0b80252d9 diff --git a/cmake/Compiler-PGI.cmake b/cmake/Compiler-PGI.cmake deleted file mode 100644 index 8ca8cfc5c..000000000 --- a/cmake/Compiler-PGI.cmake +++ /dev/null @@ -1,52 +0,0 @@ -################################################################################################### -# PGI Compiler -################################################################################################### - -if (OPENMP) - set (OPENMP_FLAGS "-mp") -else () - set (OPENMP_FLAGS "-nomp") -endif () - - -if (OPTIMIZATION STREQUAL "OFF") - set (OPTIMIZATION_FLAGS "-O0" ) -elseif (OPTIMIZATION STREQUAL "DEFENSIVE") - set (OPTIMIZATION_FLAGS "-O2 -fast") -elseif (OPTIMIZATION STREQUAL "AGGRESSIVE") - set (OPTIMIZATION_FLAGS "-O4 -fast -Mvect=sse") -endif () - -set (STANDARD_CHECK "-Mallocatable=03 -Mstandard") - -#------------------------------------------------------------------------------------------------ -# Fine tuning compilation options -set (COMPILE_FLAGS "${COMPILE_FLAGS} -Mpreprocess") -# preprocessor - -set (COMPILE_FLAGS "${COMPILE_FLAGS} -Minfo=all") -# instructs the compiler to produce information on standard error - -set (COMPILE_FLAGS "${COMPILE_FLAGS} -Minform=warn") -# instructs the compiler to display error messages at the specified and higher levels - -set (COMPILE_FLAGS "${COMPILE_FLAGS} -Mdclchk") -# instructs the compiler to require that all program variables be declared - -#------------------------------------------------------------------------------------------------O -# Runtime debugging -set (DEBUG_FLAGS "${DEBUG_FLAGS} -g") -# Includes debugging information in the object module; sets the optimization level to zero unless a -⁠O option is present on the command line -set (DEBUG_FLAGS "${DEBUG_FLAGS} -C") -# Generates code to check array bounds -set (DEBUG_FLAGS "${DEBUG_FLAGS} -Mchkptr") -# Check for NULL pointers (pgf95, pgfortran only) -set (DEBUG_FLAGS "${DEBUG_FLAGS} -Mchkstk") -# Check the stack for available space upon entry to and before the start of a parallel region. Useful when many private variables are declared -set (DEBUG_FLAGS "${DEBUG_FLAGS} -Mbounds") -# Specifies whether array bounds checking is enabled or disabled - -#------------------------------------------------------------------------------------------------ -# precision settings -set (PRECISION_FLAGS "${PRECISION_FLAGS} -r8") -# Determines whether the compiler promotes REAL variables and constants to DOUBLE PRECISION diff --git a/examples/Marc/material.yaml b/examples/Marc/material.yaml index 5c893ac52..b82fcc75c 100644 --- a/examples/Marc/material.yaml +++ b/examples/Marc/material.yaml @@ -17,7 +17,7 @@ phase: atol_xi: 1.0 dot_gamma_0_sl: 0.001 h_0_sl_sl: 75e6 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4, 1.4] n_sl: 20 output: [xi_sl] xi_0_sl: [31e6] diff --git a/examples/config/Phase_DisloUCLA_Tungsten.config b/examples/config/Phase_DisloUCLA_Tungsten.config index 2cdd5c58f..e680213d5 100644 --- a/examples/config/Phase_DisloUCLA_Tungsten.config +++ b/examples/config/Phase_DisloUCLA_Tungsten.config @@ -51,5 +51,5 @@ D0 4.0e-5 # Vacancy diffusion prefactor [m**2/ Qsd 4.5e-19 # Activation energy for climb [J] Catomicvolume 1.0 # Adj. parameter controlling the atomic volume [in b] Cedgedipmindistance 1.0 # Adj. parameter controlling the minimum dipole distance [in b] -interaction_slipslip 0.009 0.009 0.72 0.05 0.09 0.06 +interaction_slipslip 0.009 0.72 0.009 0.05 0.05 0.06 0.09 nonschmid_coefficients 0.938 0.71 4.43 0.0 0.0 0.0 diff --git a/examples/config/Phase_Dislotwin_TWIP-Steel-FeMnC.yaml b/examples/config/Phase_Dislotwin_TWIP-Steel-FeMnC.yaml index 2660f6924..610421886 100644 --- a/examples/config/Phase_Dislotwin_TWIP-Steel-FeMnC.yaml +++ b/examples/config/Phase_Dislotwin_TWIP-Steel-FeMnC.yaml @@ -19,7 +19,7 @@ TWIP_Steel_FeMnC: D_0: 4.0e-5 # Vacancy diffusion prefactor / m^2/s D_a: 1.0 # minimum dipole distance / b Q_cl: 4.5e-19 # Activation energy for climb / J - h_sl_sl: [0.122, 0.122, 0.625, 0.07, 0.137, 0.122] # Interaction coefficients (Kubin et al. 2008) + h_sl_sl: [0.122, 0.122, 0.625, 0.07, 0.137, 0.137, 0.122] # Interaction coefficients (Kubin et al. 2008) # shear band parameters xi_sb: 180.0e6 Q_sb: 3.7e-19 diff --git a/examples/config/Phase_Dislotwin_Tungsten.yaml b/examples/config/Phase_Dislotwin_Tungsten.yaml index c3f3e5864..d49149017 100644 --- a/examples/config/Phase_Dislotwin_Tungsten.yaml +++ b/examples/config/Phase_Dislotwin_Tungsten.yaml @@ -18,4 +18,4 @@ Tungsten: D_0: 4.0e-5 # Vacancy diffusion prefactor / m^2/s D_a: 1.0 # minimum dipole distance / b Q_cl: 4.5e-19 # Activation energy for climb / J - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1.4, 1, 1.4, 1.4, 1.4, 1.4] diff --git a/examples/config/Phase_Nonlocal_Aluminum.config b/examples/config/Phase_Nonlocal_Aluminum.config index c69d2fd67..0d835d5aa 100644 --- a/examples/config/Phase_Nonlocal_Aluminum.config +++ b/examples/config/Phase_Nonlocal_Aluminum.config @@ -52,7 +52,7 @@ q 1 # exponent for thermal barrier attackFrequency 50e9 # attack frequency in Hz surfaceTransmissivity 1.0 # transmissivity of free surfaces for dislocation flux grainboundaryTransmissivity 0.0 # transmissivity of grain boundaries for dislocation flux (grain bundaries are identified as interfaces with different textures on both sides); if not set or set to negative number, the subroutine automatically determines the transmissivity at the grain boundary -interaction_SlipSlip 0 0 0.625 0.07 0.137 0.122 # Dislocation interaction coefficient +interaction_SlipSlip 0 0 0.625 0.07 0.137 0.137 0.122 # Dislocation interaction coefficient linetension 0.8 # constant indicating the effect of the line tension on the hardening coefficients (0 to 1) edgejog 1.0 # fraction of annihilated screw dipoles that forms edge jogs (0 to 1) shortRangeStressCorrection 0 # switch for use of short range correction stress diff --git a/examples/config/Phase_Nonlocal_Nickel.config b/examples/config/Phase_Nonlocal_Nickel.config index c3814adeb..7bc2e1581 100644 --- a/examples/config/Phase_Nonlocal_Nickel.config +++ b/examples/config/Phase_Nonlocal_Nickel.config @@ -57,6 +57,6 @@ significantN 1 shortRangeStressCorrection 0 CFLfactor 1.1 # safety factor for CFL flux check (numerical parameter) r 1 -interaction_SlipSlip 0 0 0.625 0.07 0.137 0.122 # Dislocation interaction coefficient +interaction_SlipSlip 0 0 0.625 0.07 0.137 0.137 0.122 # Dislocation interaction coefficient linetension 0.8 edgejog 0.01 # 0.2 diff --git a/examples/config/Phase_Phenopowerlaw_Aluminum.yaml b/examples/config/Phase_Phenopowerlaw_Aluminum.yaml index a48550296..24da127fc 100644 --- a/examples/config/Phase_Phenopowerlaw_Aluminum.yaml +++ b/examples/config/Phase_Phenopowerlaw_Aluminum.yaml @@ -8,7 +8,7 @@ Aluminum: a_sl: 2.25 dot_gamma_0_sl: 0.001 h_0_sl_sl: 75e6 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4, 1.4] n_sl: 20 output: [xi_sl, gamma_sl] type: phenopowerlaw diff --git a/examples/config/Phase_Phenopowerlaw_BCC-Ferrite.yaml b/examples/config/Phase_Phenopowerlaw_BCC-Ferrite.yaml index 31c828b9f..03b7e2013 100644 --- a/examples/config/Phase_Phenopowerlaw_BCC-Ferrite.yaml +++ b/examples/config/Phase_Phenopowerlaw_BCC-Ferrite.yaml @@ -10,7 +10,7 @@ Ferrite: a_sl: 2.0 dot_gamma_0_sl: 0.001 h_0_sl_sl: 1000.0e6 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1.4, 1, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4] n_sl: 20 type: phenopowerlaw xi_0_sl: [95.e6, 96.e6] diff --git a/examples/config/Phase_Phenopowerlaw_BCC-Martensite.yaml b/examples/config/Phase_Phenopowerlaw_BCC-Martensite.yaml index 3c56f80aa..d2667f111 100644 --- a/examples/config/Phase_Phenopowerlaw_BCC-Martensite.yaml +++ b/examples/config/Phase_Phenopowerlaw_BCC-Martensite.yaml @@ -10,7 +10,7 @@ Martensite: a_sl: 2.0 dot_gamma_0_sl: 0.001 h_0_sl_sl: 563.0e9 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1.4, 1, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4, 1.4] n_sl: 20 type: phenopowerlaw xi_0_sl: [405.8e6, 456.7e6] diff --git a/examples/config/phase/damage/anisobrittle_cubic.yaml b/examples/config/phase/damage/anisobrittle_cubic.yaml index 410bcf2c3..3539a24d2 100644 --- a/examples/config/phase/damage/anisobrittle_cubic.yaml +++ b/examples/config/phase/damage/anisobrittle_cubic.yaml @@ -1,8 +1,11 @@ -N_cl: [3] -dot_o: 1e-3 -g_crit: [0.50e7] -q: 20 -s_crit: [0.006666] type: anisobrittle +N_cl: [3] +g_crit: [0.50e7] +s_crit: [0.006666] +dot_o: 1e-3 +q: 20 + +output: [f_phi] + D_11: 1.0 M: 0.001 diff --git a/examples/config/phase/damage/isobrittle_generic.yaml b/examples/config/phase/damage/isobrittle_generic.yaml new file mode 100644 index 000000000..b4416f995 --- /dev/null +++ b/examples/config/phase/damage/isobrittle_generic.yaml @@ -0,0 +1,9 @@ +type: isobrittle +W_crit: 1400000.0 +m: 1.0 +isoBrittle_atol: 0.01 + +output: [f_phi] + +D_11: 1.0 +M: 0.001 diff --git a/examples/config/phase/mechanical/plastic/dislotwin_IF-steel.yaml b/examples/config/phase/mechanical/plastic/dislotwin_IF-steel.yaml index bebe0f341..b5defb6db 100644 --- a/examples/config/phase/mechanical/plastic/dislotwin_IF-steel.yaml +++ b/examples/config/phase/mechanical/plastic/dislotwin_IF-steel.yaml @@ -18,7 +18,7 @@ q_sl: [1.55, 1.55] i_sl: [23.3, 23.3] D_a: 7.4 # C_anni B: [0.001, 0.001] -h_sl_sl: [0.1, 0.1, 0.72, 0.053, 0.137, 0.073] +h_sl_sl: [0.1, 0.72, 0.1, 0.053, 0.053, 0.073, 0.137, 0.72, 0.72, 0.053, 0.053, 0.053, 0.053, 0.073, 0.073, 0.073, 0.073, 0.073, 0.073, 0.137, 0.073, 0.073, 0.137, 0.073] D_0: 4.000E-05 Q_cl: 5.400E-19 # no recovery! D: 40e-6 # estimated diff --git a/examples/config/phase/mechanical/plastic/phenopowerlaw_Au.yaml b/examples/config/phase/mechanical/plastic/phenopowerlaw_Au.yaml index daaab878d..840f3ff30 100644 --- a/examples/config/phase/mechanical/plastic/phenopowerlaw_Au.yaml +++ b/examples/config/phase/mechanical/plastic/phenopowerlaw_Au.yaml @@ -8,7 +8,7 @@ N_sl: [12] n_sl: 83 dot_gamma_0_sl: 0.001 h_0_sl_sl: 75e6 -h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] +h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4, 1.4] a_sl: 1.0 xi_0_sl: [26e6] xi_inf_sl: [53e6] diff --git a/examples/grid/material.yaml b/examples/grid/material.yaml index b499e7413..47ab66d8b 100644 --- a/examples/grid/material.yaml +++ b/examples/grid/material.yaml @@ -17,7 +17,7 @@ phase: atol_xi: 1.0 dot_gamma_0_sl: 0.001 h_0_sl_sl: 75e6 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4, 1.4] n_sl: 20 output: [xi_sl] xi_0_sl: [31e6] diff --git a/examples/mesh/material.yaml b/examples/mesh/material.yaml index 5a135bff1..4c3026c21 100644 --- a/examples/mesh/material.yaml +++ b/examples/mesh/material.yaml @@ -1,22 +1,22 @@ homogenization: SX: N_constituents: 1 - mechanics: {type: pass} + mechanical: {type: pass} phase: Aluminum: lattice: cF - mechanics: + mechanical: output: [F, P, F_e, F_p, L_p] - elasticity: {type: Hooke, C_11: 106.75e9, C_12: 60.41e9, C_44: 28.34e9} - plasticity: + elastic: {type: Hooke, C_11: 106.75e9, C_12: 60.41e9, C_44: 28.34e9} + plastic: type: phenopowerlaw N_sl: [12] a_sl: 2.25 atol_xi: 1.0 dot_gamma_0_sl: 0.001 h_0_sl_sl: 75e6 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4, 1.4] n_sl: 20 output: [xi_sl] xi_0_sl: [31e6] diff --git a/installation/mods_MarcMentat/2020/Marc_tools/include_linux64 b/installation/mods_MarcMentat/2020/Marc_tools/include_linux64 index 23e7b589b..081042e84 100644 --- a/installation/mods_MarcMentat/2020/Marc_tools/include_linux64 +++ b/installation/mods_MarcMentat/2020/Marc_tools/include_linux64 @@ -503,7 +503,7 @@ then PROFILE=" $PROFILE -pg" fi -FORT_OPT="-c -implicitnone -stand f08 -standard-semantics -assume nostd_mod_proc_name -safe_cray_ptr -mp1 -WB -fp-model source" +FORT_OPT="-c -implicitnone -stand f18 -standard-semantics -assume nostd_mod_proc_name -safe_cray_ptr -mp1 -WB -fp-model source" if test "$MTHREAD" = "OPENMP" then FORT_OPT=" $FORT_OPT -qopenmp" @@ -541,15 +541,15 @@ fi # DAMASK compiler calls: additional flags are in line 2 OpenMP flags in line 3 -DFORTLOWMP="$FCOMP -c -O0 -qno-offload -implicitnone -stand f08 -standard-semantics -assume nostd_mod_proc_name -safe_cray_ptr $PROFILE -zero -mp1 -WB $I8FFLAGS -I$MARC_SOURCE/common \ +DFORTLOWMP="$FCOMP -c -O0 -qno-offload -implicitnone -stand f18 -standard-semantics -assume nostd_mod_proc_name -safe_cray_ptr $PROFILE -zero -mp1 -WB $I8FFLAGS -I$MARC_SOURCE/common \ -fpp -ftz -diag-disable 5268 -warn declarations -warn general -warn usage -warn interfaces -warn ignore_loc -warn alignments -DMarc4DAMASK=2020 -DDAMASKVERSION=$DAMASKVERSION \ -qopenmp -qopenmp-threadprivate=compat\ $MUMPS_INCLUDE $I8DEFINES -DLinux -DLINUX -DLinux_intel $FDEFINES $DDM $SOLVERFLAGS -I$KDTREE2_MOD" -DFORTRANMP="$FCOMP -c -O1 -qno-offload -implicitnone -stand f08 -standard-semantics -assume nostd_mod_proc_name -safe_cray_ptr $PROFILE -zero -mp1 -WB $I8FFLAGS -I$MARC_SOURCE/common \ +DFORTRANMP="$FCOMP -c -O1 -qno-offload -implicitnone -stand f18 -standard-semantics -assume nostd_mod_proc_name -safe_cray_ptr $PROFILE -zero -mp1 -WB $I8FFLAGS -I$MARC_SOURCE/common \ -fpp -ftz -diag-disable 5268 -warn declarations -warn general -warn usage -warn interfaces -warn ignore_loc -warn alignments -DMarc4DAMASK=2020 -DDAMASKVERSION=$DAMASKVERSION \ -qopenmp -qopenmp-threadprivate=compat\ $MUMPS_INCLUDE $I8DEFINES -DLinux -DLINUX -DLinux_intel $FDEFINES $DDM $SOLVERFLAGS -I$KDTREE2_MOD" -DFORTHIGHMP="$FCOMP -c -O3 -qno-offload -implicitnone -stand f08 -standard-semantics -assume nostd_mod_proc_name -safe_cray_ptr $PROFILE -zero -mp1 -WB $I8FFLAGS -I$MARC_SOURCE/common \ +DFORTHIGHMP="$FCOMP -c -O3 -qno-offload -implicitnone -stand f18 -standard-semantics -assume nostd_mod_proc_name -safe_cray_ptr $PROFILE -zero -mp1 -WB $I8FFLAGS -I$MARC_SOURCE/common \ -fpp -ftz -diag-disable 5268 -warn declarations -warn general -warn usage -warn interfaces -warn ignore_loc -warn alignments -DMarc4DAMASK=2020 -DDAMASKVERSION=$DAMASKVERSION \ -qopenmp -qopenmp-threadprivate=compat\ $MUMPS_INCLUDE $I8DEFINES -DLinux -DLINUX -DLinux_intel $FDEFINES $DDM $SOLVERFLAGS -I$KDTREE2_MOD" diff --git a/python/damask/__init__.py b/python/damask/__init__.py index ad46d454f..060c5e295 100644 --- a/python/damask/__init__.py +++ b/python/damask/__init__.py @@ -14,7 +14,6 @@ from . import tensor # noqa from . import mechanics # noqa from . import solver # noqa from . import grid_filters # noqa -from . import lattice # noqa #Modules that contain only one class (of the same name), are prefixed by a '_'. #For example, '_colormap' containsa class called 'Colormap' which is imported as 'damask.Colormap'. from ._rotation import Rotation # noqa diff --git a/python/damask/_colormap.py b/python/damask/_colormap.py index 02e2b922c..7e82853b4 100644 --- a/python/damask/_colormap.py +++ b/python/damask/_colormap.py @@ -110,13 +110,13 @@ class Colormap(mpl.colors.ListedColormap): low_,high_ = map(Colormap._rgb2msh,low_high) elif model.lower() == 'hsv': - if np.any(low_high<0) or np.any(low_high[:,1:3]>1) or np.any(low_high[:,0]>360): + if np.any(low_high<0) or np.any(low_high>[360,1,1]): raise ValueError(f'HSV color {low} | {high} are out of range.') low_,high_ = map(Colormap._hsv2msh,low_high) elif model.lower() == 'hsl': - if np.any(low_high<0) or np.any(low_high[:,1:3]>1) or np.any(low_high[:,0]>360): + if np.any(low_high<0) or np.any(low_high>[360,1,1]): raise ValueError(f'HSL color {low} | {high} are out of range.') low_,high_ = map(Colormap._hsl2msh,low_high) diff --git a/python/damask/lattice.py b/python/damask/_lattice.py similarity index 94% rename from python/damask/lattice.py rename to python/damask/_lattice.py index 768dc5ace..5b4621f3d 100644 --- a/python/damask/lattice.py +++ b/python/damask/_lattice.py @@ -63,6 +63,30 @@ kinematics = { [+1,-1,+1 , -1,+1,+2], [-1,+1,+1 , +1,-1,+2], [+1,+1,+1 , +1,+1,-2], + [+1,+1,-1 , +1,+2,+3], + [+1,-1,+1 , -1,+2,+3], + [-1,+1,+1 , +1,-2,+3], + [+1,+1,+1 , +1,+2,-3], + [+1,-1,+1 , +1,+3,+2], + [+1,+1,-1 , -1,+3,+2], + [+1,+1,+1 , +1,-3,+2], + [-1,+1,+1 , +1,+3,-2], + [+1,+1,-1 , +2,+1,+3], + [+1,-1,+1 , -2,+1,+3], + [-1,+1,+1 , +2,-1,+3], + [+1,+1,+1 , +2,+1,-3], + [+1,-1,+1 , +2,+3,+1], + [+1,+1,-1 , -2,+3,+1], + [+1,+1,+1 , +2,-3,+1], + [-1,+1,+1 , +2,+3,-1], + [-1,+1,+1 , +3,+1,+2], + [+1,+1,+1 , -3,+1,+2], + [+1,+1,-1 , +3,-1,+2], + [+1,-1,+1 , +3,+1,-2], + [-1,+1,+1 , +3,+2,+1], + [+1,+1,+1 , -3,+2,+1], + [+1,+1,-1 , +3,-2,+1], + [+1,-1,+1 , +3,+2,-1], ],'d'), 'twin' : _np.array([ [-1, 1, 1, 2, 1, 1], diff --git a/python/damask/_orientation.py b/python/damask/_orientation.py index 00e74d36c..919d9e518 100644 --- a/python/damask/_orientation.py +++ b/python/damask/_orientation.py @@ -5,6 +5,35 @@ import numpy as np from . import Rotation from . import util from . import tensor +from . import _lattice + +_crystal_families = ['triclinic', + 'monoclinic', + 'orthorhombic', + 'tetragonal', + 'hexagonal', + 'cubic'] + +_lattice_symmetries = { + 'aP': 'triclinic', + + 'mP': 'monoclinic', + 'mS': 'monoclinic', + + 'oP': 'orthorhombic', + 'oS': 'orthorhombic', + 'oI': 'orthorhombic', + 'oF': 'orthorhombic', + + 'tP': 'tetragonal', + 'tI': 'tetragonal', + + 'hP': 'hexagonal', + + 'cP': 'cubic', + 'cI': 'cubic', + 'cF': 'cubic', + } _parameter_doc = \ """lattice : str @@ -33,7 +62,7 @@ class Orientation(Rotation): """ Representation of crystallographic orientation as combination of rotation and either crystal family or Bravais lattice. - The crystal family is one of Orientation.crystal_families: + The crystal family is one of: - triclinic - monoclinic @@ -45,7 +74,7 @@ class Orientation(Rotation): and enables symmetry-related operations such as "equivalent", "reduced", "disorientation", "IPF_color", or "to_SST". - The Bravais lattice is one of Orientation.lattice_symmetries: + The Bravais lattice is given in the Pearson notation: - triclinic - aP : primitive @@ -85,35 +114,6 @@ class Orientation(Rotation): """ - crystal_families = ['triclinic', - 'monoclinic', - 'orthorhombic', - 'tetragonal', - 'hexagonal', - 'cubic'] - - lattice_symmetries = { - 'aP': 'triclinic', - - 'mP': 'monoclinic', - 'mS': 'monoclinic', - - 'oP': 'orthorhombic', - 'oS': 'orthorhombic', - 'oI': 'orthorhombic', - 'oF': 'orthorhombic', - - 'tP': 'tetragonal', - 'tI': 'tetragonal', - - 'hP': 'hexagonal', - - 'cP': 'cubic', - 'cI': 'cubic', - 'cF': 'cubic', - } - - @util.extend_docstring(_parameter_doc) def __init__(self, rotation = None, @@ -132,34 +132,17 @@ class Orientation(Rotation): Defaults to no rotation. """ - from damask.lattice import kinematics - Rotation.__init__(self) if rotation is None else Rotation.__init__(self,rotation=rotation) - if ( lattice not in self.lattice_symmetries - and lattice not in self.crystal_families): - raise KeyError(f'Lattice "{lattice}" is unknown') - - self.family = None - self.lattice = None - self.a = None - self.b = None - self.c = None - self.alpha = None - self.beta = None - self.gamma = None self.kinematics = None - if lattice in self.lattice_symmetries: - self.family = self.lattice_symmetries[lattice] + if lattice in _lattice_symmetries: + self.family = _lattice_symmetries[lattice] self.lattice = lattice + self.a = 1 if a is None else a self.b = b self.c = c - self.alpha = (np.radians(alpha) if degrees else alpha) if alpha is not None else None - self.beta = (np.radians(beta) if degrees else beta) if beta is not None else None - self.gamma = (np.radians(gamma) if degrees else gamma) if gamma is not None else None - self.a = float(self.a) if self.a is not None else \ (self.b / self.ratio['b'] if self.b is not None and self.ratio['b'] is not None else self.c / self.ratio['c'] if self.c is not None and self.ratio['c'] is not None else None) @@ -171,9 +154,13 @@ class Orientation(Rotation): (self.a * self.ratio['c'] if self.a is not None and self.ratio['c'] is not None else self.b / self.ratio['b'] * self.ratio['c'] if self.c is not None and self.ratio['b'] is not None and self.ratio['c'] is not None else None) - self.alpha = self.alpha if self.alpha is not None else self.immutable['alpha'] if 'alpha' in self.immutable else None - self.beta = self.beta if self.beta is not None else self.immutable['beta'] if 'beta' in self.immutable else None - self.gamma = self.gamma if self.gamma is not None else self.immutable['gamma'] if 'gamma' in self.immutable else None + + self.alpha = np.radians(alpha) if degrees and alpha is not None else alpha + self.beta = np.radians(beta) if degrees and beta is not None else beta + self.gamma = np.radians(gamma) if degrees and gamma is not None else gamma + if self.alpha is None and 'alpha' in self.immutable: self.alpha = self.immutable['alpha'] + if self.beta is None and 'beta' in self.immutable: self.beta = self.immutable['beta'] + if self.gamma is None and 'gamma' in self.immutable: self.gamma = self.immutable['gamma'] if \ (self.a is None) \ @@ -190,16 +177,22 @@ class Orientation(Rotation): > np.sum(np.roll([self.alpha,self.beta,self.gamma],r)[1:]) for r in range(3)]): raise ValueError ('Each lattice angle must be less than sum of others') - if self.lattice in kinematics: - master = kinematics[self.lattice] + if self.lattice in _lattice.kinematics: + master = _lattice.kinematics[self.lattice] self.kinematics = {} for m in master: self.kinematics[m] = {'direction':master[m][:,0:3],'plane':master[m][:,3:6]} \ if master[m].shape[-1] == 6 else \ {'direction':self.Bravais_to_Miller(uvtw=master[m][:,0:4]), 'plane': self.Bravais_to_Miller(hkil=master[m][:,4:8])} - elif lattice in self.crystal_families: - self.family = lattice + elif lattice in _crystal_families: + self.family = lattice + self.lattice = None + + self.a = self.b = self.c = None + self.alpha = self.beta = self.gamma = None + else: + raise KeyError(f'Lattice "{lattice}" is unknown') def __repr__(self): @@ -676,11 +669,9 @@ class Orientation(Rotation): https://doi.org/10.1016/j.actamat.2004.11.021 """ - from damask.lattice import relations - - if model not in relations: + if model not in _lattice.relations: raise KeyError(f'Orientation relationship "{model}" is unknown') - r = relations[model] + r = _lattice.relations[model] if self.lattice not in r: raise KeyError(f'Relationship "{model}" not supported for lattice "{self.lattice}"') diff --git a/python/damask/_rotation.py b/python/damask/_rotation.py index c329464ed..c76f95c8e 100644 --- a/python/damask/_rotation.py +++ b/python/damask/_rotation.py @@ -629,7 +629,7 @@ class Rotation: else: if np.any(qu[...,0] < 0.0): raise ValueError('Quaternion with negative first (real) component.') - if not np.all(np.isclose(np.linalg.norm(qu,axis=-1), 1.0)): + if not np.all(np.isclose(np.linalg.norm(qu,axis=-1), 1.0,rtol=0.0)): raise ValueError('Quaternion is not of unit length.') return Rotation(qu) diff --git a/python/damask/_table.py b/python/damask/_table.py index dd6a981b7..a65971f39 100644 --- a/python/damask/_table.py +++ b/python/damask/_table.py @@ -384,7 +384,7 @@ class Table: Returns ------- - udated : damask.Table + updated : damask.Table Updated table. """ @@ -412,7 +412,7 @@ class Table: Returns ------- - udated : damask.Table + updated : damask.Table Updated table. """ @@ -435,7 +435,7 @@ class Table: Returns ------- - udated : damask.Table + updated : damask.Table Updated table. """ @@ -461,7 +461,7 @@ class Table: Returns ------- - udated : damask.Table + updated : damask.Table Updated table. """ @@ -494,7 +494,7 @@ class Table: Returns ------- - udated : damask.Table + updated : damask.Table Updated table. """ @@ -519,7 +519,7 @@ class Table: Returns ------- - udated : damask.Table + updated : damask.Table Updated table. """ diff --git a/python/damask/_vtk.py b/python/damask/_vtk.py index 029abab19..a40a5990f 100644 --- a/python/damask/_vtk.py +++ b/python/damask/_vtk.py @@ -150,9 +150,8 @@ class VTK: ---------- fname : str or pathlib.Path Filename for reading. Valid extensions are .vtr, .vtu, .vtp, and .vtk. - dataset_type : str, optional + dataset_type : {'vtkRectilinearGrid', 'vtkUnstructuredGrid', 'vtkPolyData'}, optional Name of the vtk.vtkDataSet subclass when opening a .vtk file. - Valid types are vtkRectilinearGrid, vtkUnstructuredGrid, and vtkPolyData. Returns ------- diff --git a/python/damask/mechanics.py b/python/damask/mechanics.py index 6ba785df1..c8519abad 100644 --- a/python/damask/mechanics.py +++ b/python/damask/mechanics.py @@ -124,9 +124,6 @@ def strain(F,t,m): """ Calculate strain tensor (Seth–Hill family). - For details refer to https://en.wikipedia.org/wiki/Finite_strain_theory and - https://de.wikipedia.org/wiki/Verzerrungstensor - Parameters ---------- F : numpy.ndarray of shape (...,3,3) @@ -142,6 +139,11 @@ def strain(F,t,m): epsilon : numpy.ndarray of shape (...,3,3) Strain of F. + References + ---------- + https://en.wikipedia.org/wiki/Finite_strain_theory + https://de.wikipedia.org/wiki/Verzerrungstensor + """ if t == 'V': w,n = _np.linalg.eigh(deformation_Cauchy_Green_left(F)) @@ -150,7 +152,6 @@ def strain(F,t,m): if m > 0.0: eps = 1.0/(2.0*abs(m)) * (+ _np.einsum('...j,...kj,...lj',w**m,n,n) - _np.eye(3)) - elif m < 0.0: eps = 1.0/(2.0*abs(m)) * (- _np.einsum('...j,...kj,...lj',w**m,n,n) + _np.eye(3)) else: diff --git a/python/damask/seeds.py b/python/damask/seeds.py index 126f138aa..b2e91a352 100644 --- a/python/damask/seeds.py +++ b/python/damask/seeds.py @@ -3,8 +3,8 @@ from scipy import spatial as _spatial import numpy as _np -from . import util -from . import grid_filters +from . import util as _util +from . import grid_filters as _grid_filters def from_random(size,N_seeds,cells=None,rng_seed=None): @@ -34,7 +34,7 @@ def from_random(size,N_seeds,cells=None,rng_seed=None): if cells is None: coords = rng.random((N_seeds,3)) * size else: - grid_coords = grid_filters.coordinates0_point(cells,size).reshape(-1,3,order='F') + grid_coords = _grid_filters.coordinates0_point(cells,size).reshape(-1,3,order='F') coords = grid_coords[rng.choice(_np.prod(cells),N_seeds, replace=False)] \ + _np.broadcast_to(size/cells,(N_seeds,3))*(rng.random((N_seeds,3))*.5-.25) # wobble without leaving cells @@ -73,7 +73,7 @@ def from_Poisson_disc(size,N_seeds,N_candidates,distance,periodic=True,rng_seed= s = 1 i = 0 - progress = util._ProgressBar(N_seeds+1,'',50) + progress = _util._ProgressBar(N_seeds+1,'',50) while s < N_seeds: candidates = rng.random((N_candidates,3))*_np.broadcast_to(size,(N_candidates,3)) tree = _spatial.cKDTree(coords[:s],boxsize=size) if periodic else \ @@ -120,7 +120,7 @@ def from_grid(grid,selection=None,invert=False,average=False,periodic=True): material = grid.material.reshape((-1,1),order='F') mask = _np.full(grid.cells.prod(),True,dtype=bool) if selection is None else \ _np.isin(material,selection,invert=invert).flatten() - coords = grid_filters.coordinates0_point(grid.cells,grid.size).reshape(-1,3,order='F') + coords = _grid_filters.coordinates0_point(grid.cells,grid.size).reshape(-1,3,order='F') if not average: return (coords[mask],material[mask]) diff --git a/python/damask/util.py b/python/damask/util.py index fcf41c251..529b63d20 100644 --- a/python/damask/util.py +++ b/python/damask/util.py @@ -1,3 +1,5 @@ +"""Miscellaneous helper functionality.""" + import sys import datetime import os @@ -177,26 +179,36 @@ def execute(cmd,wd='./',env=None): def natural_sort(key): + """ + Natural sort. + + For use in python's 'sorted'. + + References + ---------- + https://en.wikipedia.org/wiki/Natural_sort_order + + """ convert = lambda text: int(text) if text.isdigit() else text return [ convert(c) for c in re.split('([0-9]+)', key) ] def show_progress(iterable,N_iter=None,prefix='',bar_length=50): """ - Decorate a loop with a status bar. + Decorate a loop with a progress bar. Use similar like enumerate. Parameters ---------- - iterable : iterable/function with yield statement + iterable : iterable or function with yield statement Iterable (or function with yield statement) to be decorated. - N_iter : int - Total # of iterations. Needed if number of iterations can not be obtained as len(iterable). - prefix : str, optional. + N_iter : int, optional + Total number of iterations. Required unless obtainable as len(iterable). + prefix : str, optional Prefix string. bar_length : int, optional - Character length of bar. Defaults to 50. + Length of progress bar in characters. Defaults to 50. """ if N_iter in [0,1] or (hasattr(iterable,'__len__') and len(iterable) <= 1): @@ -509,6 +521,7 @@ def dict_prune(d): v = dict_prune(v) if not isinstance(v,dict) or v != {}: new[k] = v + return new diff --git a/python/setup.py b/python/setup.py index e590a1197..5d66d3c88 100644 --- a/python/setup.py +++ b/python/setup.py @@ -12,7 +12,7 @@ setuptools.setup( author='The DAMASK team', author_email='damask@mpie.de', description='DAMASK library', - long_description='Python library for pre and post processing of DAMASK simulations', + long_description='Python library for managing DAMASK simulations', url='https://damask.mpie.de', packages=setuptools.find_packages(), include_package_data=True, diff --git a/python/tests/reference/Orientation/cI_slip.txt b/python/tests/reference/Orientation/cI_slip.txt index 5401392c8..0f8db6c21 100644 --- a/python/tests/reference/Orientation/cI_slip.txt +++ b/python/tests/reference/Orientation/cI_slip.txt @@ -23,3 +23,27 @@ -0.2357022603955159 0.2357022603955159 0.4714045207910318 0.23570226039551587 -0.23570226039551587 -0.47140452079103173 -0.2357022603955159 0.2357022603955159 0.4714045207910318 -0.2357022603955158 0.2357022603955158 -0.4714045207910316 0.2357022603955159 -0.2357022603955159 0.4714045207910318 0.2357022603955159 -0.2357022603955159 0.4714045207910318 0.2357022603955159 0.23570226039551587 -0.4714045207910318 0.2357022603955159 0.23570226039551587 -0.4714045207910318 0.2357022603955159 0.23570226039551587 -0.4714045207910318 +0.1543033499620919 0.3086066999241838 0.4629100498862757 0.1543033499620919 0.3086066999241838 0.4629100498862757 -0.15430334996209194 -0.3086066999241839 -0.4629100498862758 +-0.15430334996209194 0.3086066999241839 0.4629100498862758 0.1543033499620919 -0.3086066999241838 -0.4629100498862757 -0.15430334996209194 0.3086066999241839 0.4629100498862758 +-0.15430334996209188 0.30860669992418377 -0.46291004988627565 0.15430334996209194 -0.3086066999241839 0.4629100498862758 0.15430334996209194 -0.3086066999241839 0.4629100498862758 +0.15430334996209194 0.3086066999241839 -0.4629100498862758 0.15430334996209194 0.3086066999241839 -0.4629100498862758 0.15430334996209194 0.3086066999241839 -0.4629100498862758 +0.15430334996209194 0.4629100498862758 0.3086066999241838 -0.1543033499620919 -0.4629100498862757 -0.30860669992418377 0.15430334996209194 0.4629100498862758 0.3086066999241838 +-0.1543033499620919 0.4629100498862757 0.30860669992418377 -0.1543033499620919 0.4629100498862757 0.30860669992418377 0.15430334996209194 -0.4629100498862758 -0.3086066999241838 +0.15430334996209194 -0.4629100498862758 0.3086066999241839 0.15430334996209194 -0.4629100498862758 0.3086066999241839 0.15430334996209194 -0.4629100498862758 0.3086066999241839 +-0.15430334996209188 -0.46291004988627565 0.30860669992418377 0.15430334996209194 0.4629100498862758 -0.3086066999241839 0.15430334996209194 0.4629100498862758 -0.3086066999241839 +0.3086066999241838 0.15430334996209188 0.4629100498862757 0.3086066999241838 0.15430334996209188 0.4629100498862757 -0.3086066999241839 -0.1543033499620919 -0.4629100498862758 +-0.3086066999241839 0.15430334996209197 0.4629100498862758 0.3086066999241838 -0.15430334996209194 -0.4629100498862757 -0.3086066999241839 0.15430334996209197 0.4629100498862758 +-0.30860669992418377 0.1543033499620919 -0.46291004988627565 0.3086066999241839 -0.15430334996209197 0.4629100498862758 0.3086066999241839 -0.15430334996209197 0.4629100498862758 +0.3086066999241839 0.1543033499620919 -0.4629100498862758 0.3086066999241839 0.1543033499620919 -0.4629100498862758 0.3086066999241839 0.1543033499620919 -0.4629100498862758 +0.3086066999241839 0.4629100498862758 0.15430334996209188 -0.3086066999241838 -0.4629100498862757 -0.15430334996209186 0.3086066999241839 0.4629100498862758 0.15430334996209188 +-0.3086066999241838 0.4629100498862757 0.15430334996209188 -0.3086066999241838 0.4629100498862757 0.15430334996209188 0.3086066999241839 -0.4629100498862758 -0.1543033499620919 +0.3086066999241839 -0.4629100498862758 0.15430334996209194 0.3086066999241839 -0.4629100498862758 0.15430334996209194 0.3086066999241839 -0.4629100498862758 0.15430334996209194 +-0.30860669992418377 -0.46291004988627565 0.15430334996209194 0.3086066999241839 0.4629100498862758 -0.154303349962092 0.3086066999241839 0.4629100498862758 -0.154303349962092 +-0.46291004988627565 -0.15430334996209186 -0.3086066999241837 0.4629100498862758 0.1543033499620919 0.3086066999241838 0.4629100498862758 0.1543033499620919 0.3086066999241838 +-0.4629100498862758 0.15430334996209197 0.3086066999241839 -0.4629100498862758 0.15430334996209197 0.3086066999241839 -0.4629100498862758 0.15430334996209197 0.3086066999241839 +0.4629100498862757 -0.15430334996209194 0.30860669992418377 0.4629100498862757 -0.15430334996209194 0.30860669992418377 -0.4629100498862758 0.15430334996209197 -0.3086066999241838 +0.4629100498862758 0.1543033499620919 -0.3086066999241839 -0.4629100498862757 -0.15430334996209188 0.3086066999241838 0.4629100498862758 0.1543033499620919 -0.3086066999241839 +-0.46291004988627565 -0.3086066999241837 -0.1543033499620918 0.4629100498862758 0.3086066999241838 0.15430334996209188 0.4629100498862758 0.3086066999241838 0.15430334996209188 +-0.4629100498862758 0.3086066999241839 0.15430334996209194 -0.4629100498862758 0.3086066999241839 0.15430334996209194 -0.4629100498862758 0.3086066999241839 0.15430334996209194 +0.4629100498862757 -0.3086066999241838 0.1543033499620919 0.4629100498862757 -0.3086066999241838 0.1543033499620919 -0.4629100498862758 0.3086066999241839 -0.15430334996209194 +0.4629100498862758 0.3086066999241838 -0.154303349962092 -0.4629100498862757 -0.30860669992418377 0.15430334996209197 0.4629100498862758 0.3086066999241838 -0.154303349962092 diff --git a/python/tests/reference/Result/12grains6x7x8.material.yaml b/python/tests/reference/Result/12grains6x7x8.material.yaml index 436b0094e..8d3f53cfb 100644 --- a/python/tests/reference/Result/12grains6x7x8.material.yaml +++ b/python/tests/reference/Result/12grains6x7x8.material.yaml @@ -16,7 +16,7 @@ phase: atol_xi: 1.0 dot_gamma_0_sl: 0.001 h_0_sl_sl: 75e6 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4, 1.4] n_sl: 20 output: [xi_sl] type: phenopowerlaw @@ -33,7 +33,7 @@ phase: atol_xi: 1.0 dot_gamma_0_sl: 0.001 h_0_sl_sl: 75e6 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1.4, 1, 1.4, 1.4, 1.4, 1.4] n_sl: 20 output: [xi_sl] type: phenopowerlaw diff --git a/python/tests/reference/Result/4grains2x4x3.material.yaml b/python/tests/reference/Result/4grains2x4x3.material.yaml index 4cff681c0..ff2c90164 100644 --- a/python/tests/reference/Result/4grains2x4x3.material.yaml +++ b/python/tests/reference/Result/4grains2x4x3.material.yaml @@ -665,7 +665,7 @@ phase: atol_xi: 1.0 dot_gamma_0_sl: 0.001 h_0_sl_sl: 75e6 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1.4, 1, 1.4, 1.4, 1.4, 1.4] n_sl: 20 output: [xi_sl] type: phenopowerlaw diff --git a/python/tests/reference/Result/6grains6x7x8_single_phase.material.yaml b/python/tests/reference/Result/6grains6x7x8_single_phase.material.yaml index 993aa6cee..a06d679db 100644 --- a/python/tests/reference/Result/6grains6x7x8_single_phase.material.yaml +++ b/python/tests/reference/Result/6grains6x7x8_single_phase.material.yaml @@ -16,7 +16,7 @@ phase: atol_xi: 1.0 dot_gamma_0_sl: 0.001 h_0_sl_sl: 75e6 - h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4] + h_sl_sl: [1, 1, 1.4, 1.4, 1.4, 1.4, 1.4] n_sl: 20 output: [xi_sl] type: phenopowerlaw diff --git a/python/tests/test_Orientation.py b/python/tests/test_Orientation.py index dce0cf40f..3c4e5bbff 100644 --- a/python/tests/test_Orientation.py +++ b/python/tests/test_Orientation.py @@ -5,9 +5,10 @@ from itertools import permutations from damask import Rotation from damask import Orientation from damask import Table -from damask import lattice from damask import util from damask import grid_filters +from damask import _lattice as lattice +from damask._orientation import _crystal_families as crystal_families @pytest.fixture @@ -22,7 +23,7 @@ def set_of_rodrigues(set_of_quaternions): class TestOrientation: - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('shape',[None,5,(4,6)]) def test_equal(self,lattice,shape): R = Rotation.from_random(shape) @@ -30,14 +31,14 @@ class TestOrientation: (Orientation(R,lattice) == Orientation(R,lattice)).all() - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('shape',[None,5,(4,6)]) def test_unequal(self,lattice,shape): R = Rotation.from_random(shape) assert not ( Orientation(R,lattice) != Orientation(R,lattice) if shape is None else \ (Orientation(R,lattice) != Orientation(R,lattice)).any()) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('shape',[None,5,(4,6)]) def test_close(self,lattice,shape): R = Orientation.from_random(lattice=lattice,shape=shape) @@ -182,14 +183,14 @@ class TestOrientation: with pytest.raises(ValueError): Orientation(lattice='aP',a=1,b=2,c=3,alpha=45,beta=45,gamma=90.0001,degrees=True) # noqa - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('angle',[10,20,30,40]) def test_average(self,angle,lattice): o = Orientation.from_axis_angle(lattice=lattice,axis_angle=[[0,0,1,10],[0,0,1,angle]],degrees=True) avg_angle = o.average().as_axis_angle(degrees=True,pair=True)[1] assert np.isclose(avg_angle,10+(angle-10)/2.) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) def test_reduced_equivalent(self,lattice): i = Orientation(lattice=lattice) o = Orientation.from_random(lattice=lattice) @@ -197,7 +198,7 @@ class TestOrientation: FZ = np.argmin(abs(eq.misorientation(i.broadcast_to(len(eq))).as_axis_angle(pair=True)[1])) assert o.reduced == eq[FZ] - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('N',[1,8,32]) def test_disorientation(self,lattice,N): o = Orientation.from_random(lattice=lattice,shape=N) @@ -215,7 +216,7 @@ class TestOrientation: .misorientation(p[n].equivalent[ops[n][1]]) .as_quaternion()) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('a,b',[ ((2,3,2),(2,3,2)), ((2,2),(4,4)), @@ -230,20 +231,20 @@ class TestOrientation: assert o[tuple(loc[:len(o.shape)])].disorientation(p[tuple(loc[-len(p.shape):])]) \ .isclose(o.disorientation(p)[tuple(loc)]) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) def test_disorientation360(self,lattice): o_1 = Orientation(Rotation(),lattice) o_2 = Orientation.from_Euler_angles(lattice=lattice,phi=[360,0,0],degrees=True) assert np.allclose((o_1.disorientation(o_2)).as_matrix(),np.eye(3)) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('shape',[(1),(2,3),(4,3,2)]) def test_reduced_vectorization(self,lattice,shape): o = Orientation.from_random(lattice=lattice,shape=shape) for r, theO in zip(o.reduced.flatten(),o.flatten()): assert r == theO.reduced - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) def test_reduced_corner_cases(self,lattice): # test whether there is always a sym-eq rotation that falls into the FZ N = np.random.randint(10,40) @@ -253,7 +254,7 @@ class TestOrientation: assert evenly_distributed.shape == evenly_distributed.reduced.shape - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('shape',[(1),(2,3),(4,3,2)]) @pytest.mark.parametrize('vector',np.array([[1,0,0],[1,2,3],[-1,1,-1]])) @pytest.mark.parametrize('proper',[True,False]) @@ -262,7 +263,7 @@ class TestOrientation: for r, theO in zip(o.to_SST(vector=vector,proper=proper).reshape((-1,3)),o.flatten()): assert np.allclose(r,theO.to_SST(vector=vector,proper=proper)) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('shape',[(1),(2,3),(4,3,2)]) @pytest.mark.parametrize('vector',np.array([[1,0,0],[1,2,3],[-1,1,-1]])) @pytest.mark.parametrize('proper',[True,False]) @@ -272,7 +273,7 @@ class TestOrientation: for r, theO in zip(o.IPF_color(vector,in_SST=in_SST,proper=proper).reshape((-1,3)),o.flatten()): assert np.allclose(r,theO.IPF_color(vector,in_SST=in_SST,proper=proper)) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('a,b',[ ((2,3,2),(2,3,2)), ((2,2),(4,4)), @@ -300,7 +301,7 @@ class TestOrientation: assert np.allclose(np.array(color['RGB']), cube.IPF_color(vector=np.array(direction),proper=proper)) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('proper',[True,False]) def test_IPF_equivalent(self,set_of_quaternions,lattice,proper): direction = np.random.random(3)*2.0-1.0 @@ -308,13 +309,13 @@ class TestOrientation: color = o.IPF_color(vector=direction,proper=proper) assert np.allclose(np.broadcast_to(color[0,...],color.shape),color) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) def test_in_FZ_vectorization(self,set_of_rodrigues,lattice): result = Orientation.from_Rodrigues_vector(rho=set_of_rodrigues.reshape((-1,4,4)),lattice=lattice).in_FZ.reshape(-1) for r,rho in zip(result,set_of_rodrigues[:len(result)]): assert r == Orientation.from_Rodrigues_vector(rho=rho,lattice=lattice).in_FZ - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) def test_in_disorientation_FZ_vectorization(self,set_of_rodrigues,lattice): result = Orientation.from_Rodrigues_vector(rho=set_of_rodrigues.reshape((-1,4,4)), lattice=lattice).in_disorientation_FZ.reshape(-1) @@ -322,7 +323,7 @@ class TestOrientation: assert r == Orientation.from_Rodrigues_vector(rho=rho,lattice=lattice).in_disorientation_FZ @pytest.mark.parametrize('proper',[True,False]) - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) def test_in_SST_vectorization(self,lattice,proper): vecs = np.random.rand(20,4,3) result = Orientation(lattice=lattice).in_SST(vecs,proper).flatten() @@ -393,7 +394,7 @@ class TestOrientation: a=a,b=b,c=c, alpha=alpha,beta=beta,gamma=gamma).related(relation) # noqa - @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('lattice',crystal_families) @pytest.mark.parametrize('proper',[True,False]) def test_in_SST(self,lattice,proper): assert Orientation(lattice=lattice).in_SST(np.zeros(3),proper) diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index c20754a67..af0c75223 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -14,23 +14,23 @@ list(FILTER damask-sources EXCLUDE REGEX ".*commercialFEM_fileList.*.f90") if (PROJECT_NAME STREQUAL "damask-grid") - file(GLOB grid-sources grid/*.f90) + file(GLOB grid-sources grid/*.f90) - if(NOT CMAKE_BUILD_TYPE STREQUAL "SYNTAXONLY") + if (NOT CMAKE_BUILD_TYPE STREQUAL "SYNTAXONLY") add_executable(DAMASK_grid ${damask-sources} ${grid-sources}) install (TARGETS DAMASK_grid RUNTIME DESTINATION bin) - else() + else () add_library(DAMASK_grid OBJECT ${damask-sources} ${grid-sources}) exec_program (mktemp OUTPUT_VARIABLE nothing) exec_program (mktemp ARGS -d OUTPUT_VARIABLE black_hole) install (PROGRAMS ${nothing} DESTINATION ${black_hole}) - endif() + endif () elseif (PROJECT_NAME STREQUAL "damask-mesh") - file(GLOB mesh-sources mesh/*.f90) + file(GLOB mesh-sources mesh/*.f90) - add_executable(DAMASK_mesh ${damask-sources} ${mesh-sources}) - install (TARGETS DAMASK_mesh RUNTIME DESTINATION bin) + add_executable(DAMASK_mesh ${damask-sources} ${mesh-sources}) + install (TARGETS DAMASK_mesh RUNTIME DESTINATION bin) -endif() +endif () diff --git a/src/DAMASK_interface.f90 b/src/DAMASK_interface.f90 index e01c050fa..43398d0ad 100644 --- a/src/DAMASK_interface.f90 +++ b/src/DAMASK_interface.f90 @@ -96,15 +96,10 @@ subroutine DAMASK_interface_init print'(/,a)', ' Version: '//DAMASKVERSION - ! https://github.com/jeffhammond/HPCInfo/blob/master/docs/Preprocessor-Macros.md -#if defined(__PGI) - print'(/,a,i4.4,a,i8.8)', ' Compiled with PGI fortran version :', __PGIC__,& - '.', __PGIC_MINOR__ -#else print'(/,a)', ' Compiled with: '//compiler_version() print'(a)', ' Compiler options: '//compiler_options() -#endif + ! https://github.com/jeffhammond/HPCInfo/blob/master/docs/Preprocessor-Macros.md print'(/,a)', ' Compiled on: '//__DATE__//' at '//__TIME__ print'(/,a,i0,a,i0,a,i0)', & diff --git a/src/IO.f90 b/src/IO.f90 index 9ea35503b..81745441a 100644 --- a/src/IO.f90 +++ b/src/IO.f90 @@ -113,8 +113,7 @@ end function IO_readlines !-------------------------------------------------------------------------------------------------- !> @brief Read whole file. -!> @details ensures that the string ends with a new line (expected UNIX behavior) and rejects -! windows (CRLF) line endings +!> @details ensures that the string ends with a new line (expected UNIX behavior) !-------------------------------------------------------------------------------------------------- function IO_read(fileName) result(fileContent) @@ -124,8 +123,7 @@ function IO_read(fileName) result(fileContent) integer :: & fileLength, & fileUnit, & - myStat, & - firstEOL + myStat character, parameter :: CR = achar(13) @@ -143,12 +141,17 @@ function IO_read(fileName) result(fileContent) if(myStat /= 0) call IO_error(102,ext_msg=trim(fileName)) close(fileUnit) + foundCRLF: if (scan(fileContent(:index(fileContent,IO_EOL)),CR) /= 0) then + CRLF2LF: block + integer :: c + do c=1, len(fileContent) + if (fileContent(c:c) == CR) fileContent(c:c) = ' ' + enddo + end block CRLF2LF + endif foundCRLF if(fileContent(fileLength:fileLength) /= IO_EOL) fileContent = fileContent//IO_EOL ! ensure EOL@EOF - firstEOL = index(fileContent,IO_EOL) - if(scan(fileContent(firstEOL:firstEOL),CR) /= 0) call IO_error(115) - end function IO_read @@ -400,9 +403,6 @@ subroutine IO_error(error_ID,el,ip,g,instance,ext_msg) msg = 'invalid character for logical:' case (114) msg = 'cannot decode base64 string:' - case (115) - msg = 'found CR. Windows file endings (CRLF) are not supported.' - !-------------------------------------------------------------------------------------------------- ! lattice error messages @@ -644,11 +644,16 @@ subroutine selfTest if(dNeq(1.0_pReal, IO_stringAsFloat('1.0'))) error stop 'IO_stringAsFloat' if(dNeq(1.0_pReal, IO_stringAsFloat('1e0'))) error stop 'IO_stringAsFloat' if(dNeq(0.1_pReal, IO_stringAsFloat('1e-1'))) error stop 'IO_stringAsFloat' + if(dNeq(0.1_pReal, IO_stringAsFloat('1.0e-1'))) error stop 'IO_stringAsFloat' + if(dNeq(0.1_pReal, IO_stringAsFloat('1.00e-1'))) error stop 'IO_stringAsFloat' + if(dNeq(10._pReal, IO_stringAsFloat(' 1.0e+1 '))) error stop 'IO_stringAsFloat' if(3112019 /= IO_stringAsInt( '3112019')) error stop 'IO_stringAsInt' if(3112019 /= IO_stringAsInt(' 3112019')) error stop 'IO_stringAsInt' if(-3112019 /= IO_stringAsInt('-3112019')) error stop 'IO_stringAsInt' if(3112019 /= IO_stringAsInt('+3112019 ')) error stop 'IO_stringAsInt' + if(3112019 /= IO_stringAsInt('03112019 ')) error stop 'IO_stringAsInt' + if(3112019 /= IO_stringAsInt('+03112019')) error stop 'IO_stringAsInt' if(.not. IO_stringAsBool(' true')) error stop 'IO_stringAsBool' if(.not. IO_stringAsBool(' True ')) error stop 'IO_stringAsBool' diff --git a/src/commercialFEM_fileList.f90 b/src/commercialFEM_fileList.f90 index 5ddc7219e..b0bb96402 100644 --- a/src/commercialFEM_fileList.f90 +++ b/src/commercialFEM_fileList.f90 @@ -21,6 +21,7 @@ #include "lattice.f90" #include "phase.f90" #include "phase_mechanical.f90" +#include "phase_mechanical_elastic.f90" #include "phase_mechanical_plastic.f90" #include "phase_mechanical_plastic_none.f90" #include "phase_mechanical_plastic_isotropic.f90" diff --git a/src/lattice.f90 b/src/lattice.f90 index 06986b42b..055f2d2f1 100644 --- a/src/lattice.f90 +++ b/src/lattice.f90 @@ -31,44 +31,38 @@ module lattice FCC_NCLEAVAGESYSTEM = [3] !< # of cleavage systems per family for fcc integer, parameter :: & -#ifndef __PGI FCC_NSLIP = sum(FCC_NSLIPSYSTEM), & !< total # of slip systems for fcc FCC_NTWIN = sum(FCC_NTWINSYSTEM), & !< total # of twin systems for fcc FCC_NTRANS = sum(FCC_NTRANSSYSTEM), & !< total # of transformation systems for fcc FCC_NCLEAVAGE = sum(FCC_NCLEAVAGESYSTEM) !< total # of cleavage systems for fcc -#else - FCC_NSLIP = 18, & - FCC_NTWIN = 12, & - FCC_NTRANS = 12, & - FCC_NCLEAVAGE = 3 -#endif real(pReal), dimension(3+3,FCC_NSLIP), parameter :: & FCC_SYSTEMSLIP = reshape(real([& - ! Slip direction Plane normal ! SCHMID-BOAS notation - 0, 1,-1, 1, 1, 1, & ! B2 - -1, 0, 1, 1, 1, 1, & ! B4 - 1,-1, 0, 1, 1, 1, & ! B5 - 0,-1,-1, -1,-1, 1, & ! C1 - 1, 0, 1, -1,-1, 1, & ! C3 - -1, 1, 0, -1,-1, 1, & ! C5 - 0,-1, 1, 1,-1,-1, & ! A2 - -1, 0,-1, 1,-1,-1, & ! A3 - 1, 1, 0, 1,-1,-1, & ! A6 - 0, 1, 1, -1, 1,-1, & ! D1 - 1, 0,-1, -1, 1,-1, & ! D4 - -1,-1, 0, -1, 1,-1, & ! D6 - ! Slip system <110>{110} + ! <110>{111} systems + 0, 1,-1, 1, 1, 1, & ! B2 + -1, 0, 1, 1, 1, 1, & ! B4 + 1,-1, 0, 1, 1, 1, & ! B5 + 0,-1,-1, -1,-1, 1, & ! C1 + 1, 0, 1, -1,-1, 1, & ! C3 + -1, 1, 0, -1,-1, 1, & ! C5 + 0,-1, 1, 1,-1,-1, & ! A2 + -1, 0,-1, 1,-1,-1, & ! A3 + 1, 1, 0, 1,-1,-1, & ! A6 + 0, 1, 1, -1, 1,-1, & ! D1 + 1, 0,-1, -1, 1,-1, & ! D4 + -1,-1, 0, -1, 1,-1, & ! D6 + ! <110>{110}/non-octahedral systems 1, 1, 0, 1,-1, 0, & 1,-1, 0, 1, 1, 0, & 1, 0, 1, 1, 0,-1, & 1, 0,-1, 1, 0, 1, & 0, 1, 1, 0, 1,-1, & 0, 1,-1, 0, 1, 1 & - ],pReal),shape(FCC_SYSTEMSLIP)) !< Slip system <110>{111} directions. Sorted according to Eisenlohr & Hantcherli + ],pReal),shape(FCC_SYSTEMSLIP)) !< fcc slip systems real(pReal), dimension(3+3,FCC_NTWIN), parameter :: & FCC_SYSTEMTWIN = reshape(real( [& + ! <112>{111} systems -2, 1, 1, 1, 1, 1, & 1,-2, 1, 1, 1, 1, & 1, 1,-2, 1, 1, 1, & @@ -81,7 +75,7 @@ module lattice 2, 1,-1, -1, 1,-1, & -1,-2,-1, -1, 1,-1, & -1, 1, 2, -1, 1,-1 & - ],pReal),shape(FCC_SYSTEMTWIN)) !< Twin system <112>{111} directions. Sorted according to Eisenlohr & Hantcherli + ],pReal),shape(FCC_SYSTEMTWIN)) !< fcc twin systems integer, dimension(2,FCC_NTWIN), parameter, public :: & lattice_FCC_TWINNUCLEATIONSLIPPAIR = reshape( [& @@ -101,16 +95,16 @@ module lattice real(pReal), dimension(3+3,FCC_NCLEAVAGE), parameter :: & FCC_SYSTEMCLEAVAGE = reshape(real([& - ! Cleavage direction Plane normal + ! <001>{001} systems 0, 1, 0, 1, 0, 0, & 0, 0, 1, 0, 1, 0, & 1, 0, 0, 0, 0, 1 & - ],pReal),shape(FCC_SYSTEMCLEAVAGE)) + ],pReal),shape(FCC_SYSTEMCLEAVAGE)) !< fcc cleavage systems !-------------------------------------------------------------------------------------------------- ! body centered cubic (cI) integer, dimension(*), parameter :: & - BCC_NSLIPSYSTEM = [12, 12] !< # of slip systems per family for bcc + BCC_NSLIPSYSTEM = [12, 12, 24] !< # of slip systems per family for bcc integer, dimension(*), parameter :: & BCC_NTWINSYSTEM = [12] !< # of twin systems per family for bcc @@ -119,50 +113,68 @@ module lattice BCC_NCLEAVAGESYSTEM = [3] !< # of cleavage systems per family for bcc integer, parameter :: & -#ifndef __PGI BCC_NSLIP = sum(BCC_NSLIPSYSTEM), & !< total # of slip systems for bcc BCC_NTWIN = sum(BCC_NTWINSYSTEM), & !< total # of twin systems for bcc BCC_NCLEAVAGE = sum(BCC_NCLEAVAGESYSTEM) !< total # of cleavage systems for bcc -#else - BCC_NSLIP = 24, & - BCC_NTWIN = 12, & - BCC_NCLEAVAGE = 3 -#endif real(pReal), dimension(3+3,BCC_NSLIP), parameter :: & BCC_SYSTEMSLIP = reshape(real([& - ! Slip direction Plane normal - ! Slip system <111>{110} - 1,-1, 1, 0, 1, 1, & - -1,-1, 1, 0, 1, 1, & - 1, 1, 1, 0,-1, 1, & - -1, 1, 1, 0,-1, 1, & - -1, 1, 1, 1, 0, 1, & - -1,-1, 1, 1, 0, 1, & - 1, 1, 1, -1, 0, 1, & - 1,-1, 1, -1, 0, 1, & - -1, 1, 1, 1, 1, 0, & - -1, 1,-1, 1, 1, 0, & - 1, 1, 1, -1, 1, 0, & - 1, 1,-1, -1, 1, 0, & - ! Slip system <111>{112} - -1, 1, 1, 2, 1, 1, & - 1, 1, 1, -2, 1, 1, & - 1, 1,-1, 2,-1, 1, & - 1,-1, 1, 2, 1,-1, & - 1,-1, 1, 1, 2, 1, & - 1, 1,-1, -1, 2, 1, & - 1, 1, 1, 1,-2, 1, & - -1, 1, 1, 1, 2,-1, & - 1, 1,-1, 1, 1, 2, & - 1,-1, 1, -1, 1, 2, & - -1, 1, 1, 1,-1, 2, & - 1, 1, 1, 1, 1,-2 & - ],pReal),shape(BCC_SYSTEMSLIP)) + ! <111>{110} systems + 1,-1, 1, 0, 1, 1, & ! D1 + -1,-1, 1, 0, 1, 1, & ! C1 + 1, 1, 1, 0,-1, 1, & ! B2 + -1, 1, 1, 0,-1, 1, & ! A2 + -1, 1, 1, 1, 0, 1, & ! A3 + -1,-1, 1, 1, 0, 1, & ! C3 + 1, 1, 1, -1, 0, 1, & ! B4 + 1,-1, 1, -1, 0, 1, & ! D4 + -1, 1, 1, 1, 1, 0, & ! A6 + -1, 1,-1, 1, 1, 0, & ! D6 + 1, 1, 1, -1, 1, 0, & ! B5 + 1, 1,-1, -1, 1, 0, & ! C5 + ! <111>{112} systems + -1, 1, 1, 2, 1, 1, & ! A-4 + 1, 1, 1, -2, 1, 1, & ! B-3 + 1, 1,-1, 2,-1, 1, & ! C-10 + 1,-1, 1, 2, 1,-1, & ! D-9 + 1,-1, 1, 1, 2, 1, & ! D-6 + 1, 1,-1, -1, 2, 1, & ! C-5 + 1, 1, 1, 1,-2, 1, & ! B-12 + -1, 1, 1, 1, 2,-1, & ! A-11 + 1, 1,-1, 1, 1, 2, & ! C-2 + 1,-1, 1, -1, 1, 2, & ! D-1 + -1, 1, 1, 1,-1, 2, & ! A-8 + 1, 1, 1, 1, 1,-2, & ! B-7 + ! Slip system <111>{123} + 1, 1,-1, 1, 2, 3, & + 1,-1, 1, -1, 2, 3, & + -1, 1, 1, 1,-2, 3, & + 1, 1, 1, 1, 2,-3, & + 1,-1, 1, 1, 3, 2, & + 1, 1,-1, -1, 3, 2, & + 1, 1, 1, 1,-3, 2, & + -1, 1, 1, 1, 3,-2, & + 1, 1,-1, 2, 1, 3, & + 1,-1, 1, -2, 1, 3, & + -1, 1, 1, 2,-1, 3, & + 1, 1, 1, 2, 1,-3, & + 1,-1, 1, 2, 3, 1, & + 1, 1,-1, -2, 3, 1, & + 1, 1, 1, 2,-3, 1, & + -1, 1, 1, 2, 3,-1, & + -1, 1, 1, 3, 1, 2, & + 1, 1, 1, -3, 1, 2, & + 1, 1,-1, 3,-1, 2, & + 1,-1, 1, 3, 1,-2, & + -1, 1, 1, 3, 2, 1, & + 1, 1, 1, -3, 2, 1, & + 1, 1,-1, 3,-2, 1, & + 1,-1, 1, 3, 2,-1 & + ],pReal),shape(BCC_SYSTEMSLIP)) !< bcc slip systems real(pReal), dimension(3+3,BCC_NTWIN), parameter :: & BCC_SYSTEMTWIN = reshape(real([& - ! Twin system <111>{112} + ! <111>{112} systems -1, 1, 1, 2, 1, 1, & 1, 1, 1, -2, 1, 1, & 1, 1,-1, 2,-1, 1, & @@ -175,15 +187,15 @@ module lattice 1,-1, 1, -1, 1, 2, & -1, 1, 1, 1,-1, 2, & 1, 1, 1, 1, 1,-2 & - ],pReal),shape(BCC_SYSTEMTWIN)) + ],pReal),shape(BCC_SYSTEMTWIN)) !< bcc twin systems real(pReal), dimension(3+3,BCC_NCLEAVAGE), parameter :: & BCC_SYSTEMCLEAVAGE = reshape(real([& - ! Cleavage direction Plane normal + ! <001>{001} systems 0, 1, 0, 1, 0, 0, & 0, 0, 1, 0, 1, 0, & 1, 0, 0, 0, 0, 1 & - ],pReal),shape(BCC_SYSTEMCLEAVAGE)) + ],pReal),shape(BCC_SYSTEMCLEAVAGE)) !< bcc cleavage systems !-------------------------------------------------------------------------------------------------- ! hexagonal (hP) @@ -194,37 +206,31 @@ module lattice HEX_NTWINSYSTEM = [6, 6, 6, 6] !< # of slip systems per family for hex integer, parameter :: & -#ifndef __PGI HEX_NSLIP = sum(HEX_NSLIPSYSTEM), & !< total # of slip systems for hex HEX_NTWIN = sum(HEX_NTWINSYSTEM) !< total # of twin systems for hex -#else - HEX_NSLIP = 33, & - HEX_NTWIN = 24 -#endif real(pReal), dimension(4+4,HEX_NSLIP), parameter :: & HEX_SYSTEMSLIP = reshape(real([& - ! Slip direction Plane normal - ! Basal systems <-1-1.0>{00.1} (independent of c/a-ratio, Bravais notation (4 coordinate base)) + ! <-1-1.0>{00.1}/basal systems (independent of c/a-ratio) 2, -1, -1, 0, 0, 0, 0, 1, & -1, 2, -1, 0, 0, 0, 0, 1, & -1, -1, 2, 0, 0, 0, 0, 1, & - ! 1st type prismatic systems <-1-1.0>{1-1.0} (independent of c/a-ratio) + ! <-1-1.0>{1-1.0}/prismatic systems (independent of c/a-ratio) 2, -1, -1, 0, 0, 1, -1, 0, & -1, 2, -1, 0, -1, 0, 1, 0, & -1, -1, 2, 0, 1, -1, 0, 0, & - ! 2nd type prismatic systems <-11.0>{11.0} -- a slip; plane normals independent of c/a-ratio + ! <-11.0>{11.0}/2nd order prismatic compound systems (plane normal independent of c/a-ratio) -1, 1, 0, 0, 1, 1, -2, 0, & 0, -1, 1, 0, -2, 1, 1, 0, & 1, 0, -1, 0, 1, -2, 1, 0, & - ! 1st type 1st order pyramidal systems <-1-1.0>{-11.1} -- plane normals depend on the c/a-ratio + ! <-1-1.0>{-11.1}/1st order pyramidal systems (direction independent of c/a-ratio) -1, 2, -1, 0, 1, 0, -1, 1, & -2, 1, 1, 0, 0, 1, -1, 1, & -1, -1, 2, 0, -1, 1, 0, 1, & 1, -2, 1, 0, -1, 0, 1, 1, & 2, -1, -1, 0, 0, -1, 1, 1, & 1, 1, -2, 0, 1, -1, 0, 1, & - ! pyramidal system: c+a slip <11.3>{-10.1} -- plane normals depend on the c/a-ratio + ! <11.3>{-10.1}/1st order pyramidal systems (direction independent of c/a-ratio) -2, 1, 1, 3, 1, 0, -1, 1, & -1, -1, 2, 3, 1, 0, -1, 1, & -1, -1, 2, 3, 0, 1, -1, 1, & @@ -237,96 +243,96 @@ module lattice -1, 2, -1, 3, 0, -1, 1, 1, & -1, 2, -1, 3, 1, -1, 0, 1, & -2, 1, 1, 3, 1, -1, 0, 1, & - ! pyramidal system: c+a slip <11.3>{-1-1.2} -- as for hexagonal ice (Castelnau et al. 1996, similar to twin system found below) - -1, -1, 2, 3, 1, 1, -2, 2, & ! <11.3>{-1-1.2} shear = 2((c/a)^2-2)/(3 c/a) + ! <11.3>{-1-1.2}/2nd order pyramidal systems + -1, -1, 2, 3, 1, 1, -2, 2, & 1, -2, 1, 3, -1, 2, -1, 2, & 2, -1, -1, 3, -2, 1, 1, 2, & 1, 1, -2, 3, -1, -1, 2, 2, & -1, 2, -1, 3, 1, -2, 1, 2, & -2, 1, 1, 3, 2, -1, -1, 2 & - ],pReal),shape(HEX_SYSTEMSLIP)) !< slip systems for hex, sorted by P. Eisenlohr CCW around starting next to a_1 axis + ],pReal),shape(HEX_SYSTEMSLIP)) !< hex slip systems, sorted by P. Eisenlohr CCW around starting next to a_1 axis real(pReal), dimension(4+4,HEX_NTWIN), parameter :: & HEX_SYSTEMTWIN = reshape(real([& - ! Compression or Tension = f(twinning shear=f(c/a)) for each metal ! (according to Yoo 1981) - -1, 0, 1, 1, 1, 0, -1, 2, & ! <-10.1>{10.2} shear = (3-(c/a)^2)/(sqrt(3) c/a) + ! <-10.1>{10.2} systems, shear = (3-(c/a)^2)/(sqrt(3) c/a) + ! tension in Co, Mg, Zr, Ti, and Be; compression in Cd and Zn + -1, 0, 1, 1, 1, 0, -1, 2, & ! 0, -1, 1, 1, 0, 1, -1, 2, & 1, -1, 0, 1, -1, 1, 0, 2, & 1, 0, -1, 1, -1, 0, 1, 2, & 0, 1, -1, 1, 0, -1, 1, 2, & -1, 1, 0, 1, 1, -1, 0, 2, & -! - -1, -1, 2, 6, 1, 1, -2, 1, & ! <11.6>{-1-1.1} shear = 1/(c/a) + ! <11.6>{-1-1.1} systems, shear = 1/(c/a) + ! tension in Co, Re, and Zr + -1, -1, 2, 6, 1, 1, -2, 1, & 1, -2, 1, 6, -1, 2, -1, 1, & 2, -1, -1, 6, -2, 1, 1, 1, & 1, 1, -2, 6, -1, -1, 2, 1, & -1, 2, -1, 6, 1, -2, 1, 1, & -2, 1, 1, 6, 2, -1, -1, 1, & -! - 1, 0, -1, -2, 1, 0, -1, 1, & ! <10.-2>{10.1} shear = (4(c/a)^2-9)/(4 sqrt(3) c/a) + ! <10.-2>{10.1} systems, shear = (4(c/a)^2-9)/(4 sqrt(3) c/a) + ! compression in Mg + 1, 0, -1, -2, 1, 0, -1, 1, & 0, 1, -1, -2, 0, 1, -1, 1, & -1, 1, 0, -2, -1, 1, 0, 1, & -1, 0, 1, -2, -1, 0, 1, 1, & 0, -1, 1, -2, 0, -1, 1, 1, & 1, -1, 0, -2, 1, -1, 0, 1, & -! - 1, 1, -2, -3, 1, 1, -2, 2, & ! <11.-3>{11.2} shear = 2((c/a)^2-2)/(3 c/a) + ! <11.-3>{11.2} systems, shear = 2((c/a)^2-2)/(3 c/a) + ! compression in Ti and Zr + 1, 1, -2, -3, 1, 1, -2, 2, & -1, 2, -1, -3, -1, 2, -1, 2, & -2, 1, 1, -3, -2, 1, 1, 2, & -1, -1, 2, -3, -1, -1, 2, 2, & 1, -2, 1, -3, 1, -2, 1, 2, & 2, -1, -1, -3, 2, -1, -1, 2 & - ],pReal),shape(HEX_SYSTEMTWIN)) !< twin systems for hex, sorted by P. Eisenlohr CCW around starting next to a_1 axis + ],pReal),shape(HEX_SYSTEMTWIN)) !< hex twin systems, sorted by P. Eisenlohr CCW around starting next to a_1 axis !-------------------------------------------------------------------------------------------------- ! body centered tetragonal (tI) integer, dimension(*), parameter :: & - BCT_NSLIPSYSTEM = [2, 2, 2, 4, 2, 4, 2, 2, 4, 8, 4, 8, 8 ] !< # of slip systems per family for bct (Sn) Bieler J. Electr Mater 2009 + BCT_NSLIPSYSTEM = [2, 2, 2, 4, 2, 4, 2, 2, 4, 8, 4, 8, 8 ] !< # of slip systems per family for bct integer, parameter :: & -#ifndef __PGI BCT_NSLIP = sum(BCT_NSLIPSYSTEM) !< total # of slip systems for bct -#else - BCT_NSLIP = 52 -#endif + real(pReal), dimension(3+3,BCT_NSLIP), parameter :: & BCT_SYSTEMSLIP = reshape(real([& - ! Slip direction Plane normal - ! Slip family 1 {100)<001] (Bravais notation {hkl) acting - 2, 1, 2, 6, 4, 5, 5, 4, 6, 5, 3, 5, 9,10,11,12, 9,10, & ! | - 2, 2, 1, 5, 5, 3, 5, 6, 4, 6, 5, 4, 11,12, 9,10, 9,10, & ! | - 4, 6, 5, 1, 2, 2, 4, 5, 6, 3, 5, 5, 9,10,10, 9,12,11, & ! v - 6, 4, 5, 2, 1, 2, 5, 3, 5, 5, 4, 6, 9,10,12,11,10, 9, & ! reacting - 5, 5, 3, 2, 2, 1, 6, 5, 4, 5, 6, 4, 11,12,10, 9,10, 9, & - 3, 5, 5, 4, 5, 6, 1, 2, 2, 4, 6, 5, 10, 9,10, 9,11,12, & - 5, 4, 6, 5, 3, 5, 2, 1, 2, 6, 4, 5, 10, 9,12,11, 9,10, & - 5, 6, 4, 6, 5, 4, 2, 2, 1, 5, 5, 3, 12,11,10, 9, 9,10, & - 4, 5, 6, 3, 5, 5, 4, 6, 5, 1, 2, 2, 10, 9, 9,10,12,11, & - 5, 3, 5, 5, 4, 6, 6, 4, 5, 2, 1, 2, 10, 9,11,12,10, 9, & - 6, 5, 4, 5, 6, 4, 5, 5, 3, 2, 2, 1, 12,11, 9,10,10, 9, & + 1, 2, 2, 4, 7, 5, 3, 5, 5, 4, 6, 7, 10,11,10,11,12,13, & ! -----> acting (forest) + 2, 1, 2, 7, 4, 5, 6, 4, 7, 5, 3, 5, 10,11,12,13,10,11, & ! | + 2, 2, 1, 5, 5, 3, 6, 7, 4, 7, 6, 4, 12,13,10,11,10,11, & ! | + 4, 7, 6, 1, 2, 2, 4, 6, 7, 3, 5, 5, 10,11,11,10,13,12, & ! v + 7, 4, 6, 2, 1, 2, 5, 3, 5, 6, 4, 7, 10,11,13,12,11,10, & ! reacting (primary) + 5, 5, 3, 2, 2, 1, 7, 6, 4, 6, 7, 4, 12,13,11,10,11,10, & + 3, 5, 5, 4, 6, 7, 1, 2, 2, 4, 7, 6, 11,10,11,10,12,13, & + 6, 4, 7, 5, 3, 5, 2, 1, 2, 7, 4, 6, 11,10,13,12,10,11, & + 6, 7, 4, 7, 6, 4, 2, 2, 1, 5, 5, 3, 13,12,11,10,10,11, & + 4, 6, 7, 3, 5, 5, 4, 7, 6, 1, 2, 2, 11,10,10,11,13,12, & + 5, 3, 5, 6, 4, 7, 7, 4, 6, 2, 1, 2, 11,10,12,13,11,10, & + 7, 6, 4, 6, 7, 4, 5, 5, 3, 2, 2, 1, 13,12,10,11,11,10, & - 9, 9,11, 9, 9,11,10,10,12,10,10,12, 1, 7, 8, 8, 8, 8, & - 10,10,12,10,10,12, 9, 9,11, 9, 9,11, 7, 1, 8, 8, 8, 8, & - 9,11, 9,10,12,10,10,12,10, 9,11, 9, 8, 8, 1, 7, 8, 8, & - 10,12,10, 9,11, 9, 9,11, 9,10,12,10, 8, 8, 7, 1, 8, 8, & - 11, 9, 9,12,10,10,11, 9, 9,12,10,10, 8, 8, 8, 8, 1, 7, & - 12,10,10,11, 9, 9,12,10,10,11, 9, 9, 8, 8, 8, 8, 7, 1 & - ],shape(FCC_INTERACTIONSLIPSLIP)) !< Slip--slip interaction types for fcc - !< 1: self interaction - !< 2: coplanar interaction - !< 3: collinear interaction - !< 4: Hirth locks - !< 5: glissile junctions - !< 6: Lomer locks - !< 7: crossing (similar to Hirth locks in <110>{111} for two {110} planes) - !< 8: similar to Lomer locks in <110>{111} for two {110} planes - !< 9: similar to Lomer locks in <110>{111} btw one {110} and one {111} plane - !<10: similar to glissile junctions in <110>{111} btw one {110} and one {111} plane - !<11: crossing btw one {110} and one {111} plane - !<12: collinear btw one {110} and one {111} plane + 10,10,12,10,10,12,11,11,13,11,11,13, 1, 8, 9, 9, 9, 9, & + 11,11,13,11,11,13,10,10,12,10,10,12, 8, 1, 9, 9, 9, 9, & + 10,12,10,11,13,11,11,13,11,10,12,10, 9, 9, 1, 8, 9, 9, & + 11,13,11,10,12,10,10,12,10,11,13,11, 9, 9, 8, 1, 9, 9, & + 12,10,10,13,11,11,12,10,10,13,11,11, 9, 9, 9, 9, 1, 8, & + 13,11,11,12,10,10,13,11,11,12,10,10, 9, 9, 9, 9, 8, 1 & + ],shape(FCC_INTERACTIONSLIPSLIP)) !< Slip-slip interaction types for fcc / Madec 2017 (https://doi.org/10.1016/j.actamat.2016.12.040) + !< 1: self interaction --> alpha 0 + !< 2: coplanar interaction --> alpha copla + !< 3: collinear interaction --> alpha coli + !< 4: Hirth locks --> alpha 1 + !< 5: glissile junctions I --> alpha 2 + !< 6: glissile junctions II --> alpha 2* + !< 7: Lomer locks --> alpha 3 + !< 8: crossing (similar to Hirth locks in <110>{111} for two {110} planes) + !< 9: similar to Lomer locks in <110>{111} for two {110} planes + !<10: similar to Lomer locks in <110>{111} btw one {110} and one {111} plane + !<11: similar to glissile junctions in <110>{111} btw one {110} and one {111} plane + !<12: crossing btw one {110} and one {111} plane + !<13: collinear btw one {110} and one {111} plane integer, dimension(BCC_NSLIP,BCC_NSLIP), parameter :: & BCC_INTERACTIONSLIPSLIP = reshape( [& - 1,2,6,6,5,4,4,3,4,3,5,4, 6,6,4,3,3,4,6,6,4,3,6,6, & ! -----> acting - 2,1,6,6,4,3,5,4,5,4,4,3, 6,6,3,4,4,3,6,6,3,4,6,6, & ! | - 6,6,1,2,4,5,3,4,4,5,3,4, 4,3,6,6,6,6,3,4,6,6,4,3, & ! | - 6,6,2,1,3,4,4,5,3,4,4,5, 3,4,6,6,6,6,4,3,6,6,3,4, & ! v - 5,4,4,3,1,2,6,6,3,4,5,4, 3,6,4,6,6,4,6,3,4,6,3,6, & ! reacting - 4,3,5,4,2,1,6,6,4,5,4,3, 4,6,3,6,6,3,6,4,3,6,4,6, & - 4,5,3,4,6,6,1,2,5,4,3,4, 6,3,6,4,4,6,3,6,6,4,6,3, & - 3,4,4,5,6,6,2,1,4,3,4,5, 6,4,6,3,3,6,4,6,6,3,6,4, & - 4,5,4,3,3,4,5,4,1,2,6,6, 3,6,6,4,4,6,6,3,6,4,3,6, & - 3,4,5,4,4,5,4,3,2,1,6,6, 4,6,6,3,3,6,6,4,6,3,4,6, & - 5,4,3,4,5,4,3,4,6,6,1,2, 6,3,4,6,6,4,3,6,4,6,6,3, & - 4,3,4,5,4,3,4,5,6,6,2,1, 6,4,3,6,6,3,4,6,3,6,6,4, & - ! - 6,6,4,3,3,4,6,6,3,4,6,6, 1,5,6,6,5,6,6,3,5,6,3,6, & - 6,6,3,4,6,6,3,4,6,6,3,4, 5,1,6,6,6,5,3,6,6,5,6,3, & - 4,3,6,6,4,3,6,6,6,6,4,3, 6,6,1,5,6,3,5,6,3,6,5,6, & - 3,4,6,6,6,6,4,3,4,3,6,6, 6,6,5,1,3,6,6,5,6,3,6,5, & - 3,4,6,6,6,6,4,3,4,3,6,6, 5,6,6,3,1,6,5,6,5,3,6,6, & - 4,3,6,6,4,3,6,6,6,6,4,3, 6,5,3,6,6,1,6,5,3,5,6,6, & - 6,6,3,4,6,6,3,4,6,6,3,4, 6,3,5,6,5,6,1,6,6,6,5,3, & - 6,6,4,3,3,4,6,6,3,4,6,6, 3,6,6,5,6,5,6,1,6,6,3,5, & - 4,3,6,6,4,3,6,6,6,6,4,3, 5,6,3,6,5,3,6,6,1,6,6,5, & - 3,4,6,6,6,6,4,3,4,3,6,6, 6,5,6,3,3,5,6,6,6,1,5,6, & - 6,6,4,3,3,4,6,6,3,4,6,6, 3,6,5,6,6,6,5,3,6,5,1,6, & - 6,6,3,4,6,6,3,4,6,6,3,4, 6,3,6,5,6,6,3,5,5,6,6,1 & - ],shape(BCC_INTERACTIONSLIPSLIP)) !< Slip--slip interaction types for bcc from Queyreau et al. Int J Plast 25 (2009) 361–377 - !< 1: self interaction - !< 2: coplanar interaction - !< 3: collinear interaction - !< 4: mixed-asymmetrical junction - !< 5: mixed-symmetrical junction - !< 6: edge junction + 1, 3, 6, 6, 7, 5, 4, 2, 4, 2, 7, 5, 18, 18, 11, 8, 9, 13, 17, 14, 13, 9, 17, 14, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 28, 28, 25, 28, 28, 28, 25, &! -----> acting (forest) + 3, 1, 6, 6, 4, 2, 7, 5, 7, 5, 4, 2, 18, 18, 8, 11, 13, 9, 14, 17, 9, 13, 14, 17, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 28, 25, 28, &! | + 6, 6, 1, 3, 5, 7, 2, 4, 5, 7, 2, 4, 11, 8, 18, 18, 17, 14, 9, 13, 17, 14, 13, 9, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 25, 28, 28, 28, 25, 28, 28, &! | + 6, 6, 3, 1, 2, 4, 5, 7, 2, 4, 5, 7, 8, 11, 18, 18, 14, 17, 13, 9, 14, 17, 9, 13, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 25, 28, 28, 28, 25, 28, 28, 28, &! v + 7, 5, 4, 2, 1, 3, 6, 6, 2, 4, 7, 5, 9, 17, 13, 14, 18, 11, 18, 8, 13, 17, 9, 14, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 25, 28, 28, 28, 25, 28, 28, 28, &! reacting (primary) + 4, 2, 7, 5, 3, 1, 6, 6, 5, 7, 4, 2, 13, 14, 9, 17, 18, 8, 18, 11, 9, 14, 13, 17, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 28, 25, 28, & + 5, 7, 2, 4, 6, 6, 1, 3, 7, 5, 2, 4, 17, 9, 14, 13, 11, 18, 8, 18, 17, 13, 14, 9, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 25, 28, 28, 28, 25, 28, 28, & + 2, 4, 5, 7, 6, 6, 3, 1, 4, 2, 5, 7, 14, 13, 17, 9, 8, 18, 11, 18, 14, 9, 17, 13, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 28, 28, 25, 28, 28, 28, 25, & + 5, 7, 4, 2, 2, 4, 7, 5, 1, 3, 6, 6, 9, 17, 14, 13, 13, 17, 14, 9, 18, 11, 8, 18, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 25, 28, 28, 28, 25, 28, 28, 28, & + 2, 4, 7, 5, 5, 7, 4, 2, 3, 1, 6, 6, 13, 14, 17, 9, 9, 14, 17, 13, 18, 8, 11, 18, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 28, 28, 25, 28, 28, 28, 25, & + 7, 5, 2, 4, 7, 5, 2, 4, 6, 6, 1, 3, 17, 9, 13, 14, 17, 13, 9, 14, 11, 18, 18, 8, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 25, 28, 28, 28, 25, 28, 28, & + 4, 2, 5, 7, 4, 2, 5, 7, 6, 6, 3, 1, 14, 13, 9, 17, 14, 9, 13, 17, 8, 18, 18, 11, 25, 28, 28, 28, 28, 25, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 28, 28, 25, 28, 28, 28, 25, 28, & + ! + 19, 19, 10, 8, 9, 12, 16, 15, 9, 12, 16, 15, 1, 20, 24, 24, 23, 22, 21, 2, 23, 22, 2, 21, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 26, 28, 28, 28, 26, 28, 28, 28, & + 19, 19, 8, 10, 16, 15, 9, 12, 16, 15, 9, 12, 20, 1, 24, 24, 22, 23, 2, 21, 22, 23, 21, 2, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 26, 28, 28, 28, 26, 28, 28, & + 10, 8, 19, 19, 12, 9, 15, 16, 15, 16, 12, 9, 24, 24, 1, 20, 21, 2, 23, 22, 2, 21, 23, 22, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 26, 28, & + 8, 10, 19, 19, 15, 16, 12, 9, 12, 9, 15, 16, 24, 24, 20, 1, 2, 21, 22, 23, 21, 2, 22, 23, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 28, 28, 26, 28, 28, 28, 26, & + 9, 12, 16, 15, 19, 19, 10, 8, 12, 9, 16, 15, 23, 21, 22, 2, 1, 24, 20, 24, 23, 2, 22, 21, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 28, 28, 26, 28, 28, 28, 26, & + 12, 9, 15, 16, 10, 8, 19, 19, 16, 15, 12, 9, 21, 23, 2, 21, 24, 1, 24, 20, 2, 23, 21, 22, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 26, 28, & + 16, 15, 9, 12, 19, 19, 8, 10, 15, 16, 9, 12, 22, 2, 23, 22, 20, 24, 1, 24, 22, 21, 23, 2, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 26, 28, 28, 28, 26, 28, 28, & + 15, 16, 12, 9, 8, 10, 19, 19, 9, 12, 15, 16, 2, 22, 21, 23, 24, 20, 24, 1, 21, 22, 2, 23, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 26, 28, 28, 28, 26, 28, 28, 28, & + 12, 9, 16, 15, 12, 9, 16, 15, 19, 19, 10, 8, 23, 21, 2, 22, 23, 2, 21, 22, 1, 24, 24, 20, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 26, 28, & + 9, 12, 15, 16, 16, 15, 12, 9, 10, 8, 19, 19, 21, 23, 22, 2, 2, 23, 22, 21, 24, 1, 20, 24, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 28, 28, 26, 28, 28, 28, 26, & + 16, 15, 12, 9, 9, 12, 15, 16, 8, 10, 19, 19, 2, 22, 23, 21, 21, 22, 23, 2, 24, 20, 1, 24, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 26, 28, 28, 28, 26, 28, 28, 28, & + 15, 16, 9, 12, 15, 16, 9, 12, 19, 19, 8, 10, 22, 2, 21, 23, 22, 21, 2, 23, 20, 24, 24, 1, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 26, 28, 28, 28, 26, 28, 28, & + ! + 28, 25, 28, 28, 28, 25, 28, 28, 28, 28, 28, 25, 28, 28, 26, 28, 28, 26, 28, 28, 26, 28, 28, 28, 1, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 28, 27, 28, & + 25, 28, 28, 28, 28, 28, 28, 25, 28, 25, 28, 28, 28, 28, 28, 26, 26, 28, 28, 28, 28, 26, 28, 28, 28, 1, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 28, 28, 27, 28, 28, 28, 27, & + 28, 28, 28, 25, 25, 28, 28, 28, 25, 28, 28, 28, 26, 28, 28, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 1, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 27, 28, 28, 28, 27, 28, 28, 28, & + 28, 28, 25, 28, 28, 28, 25, 28, 28, 28, 25, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 1, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 27, 28, 28, 28, 27, 28, 28, & + 25, 28, 28, 28, 28, 28, 28, 25, 28, 25, 28, 28, 28, 28, 28, 26, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 1, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 28, 28, 27, 28, 28, 28, 27, & + 28, 25, 28, 28, 28, 25, 28, 28, 28, 28, 28, 25, 28, 28, 26, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 1, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 28, 27, 28, & + 28, 28, 25, 28, 28, 28, 25, 28, 28, 28, 25, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 1, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 27, 28, 28, 28, 27, 28, 28, & + 28, 28, 28, 25, 25, 28, 28, 28, 25, 28, 28, 28, 26, 28, 28, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 1, 28, 28, 27, 28, 28, 28, 28, 27, 27, 28, 28, 28, 27, 28, 28, 28, & + 28, 25, 28, 28, 28, 25, 28, 28, 28, 28, 28, 25, 28, 28, 26, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 1, 28, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 28, 27, 28, & + 25, 28, 28, 28, 28, 28, 28, 25, 28, 25, 28, 28, 28, 28, 28, 26, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 1, 28, 28, 27, 28, 28, 28, 28, 28, 28, 27, 28, 28, 28, 27, & + 28, 28, 28, 25, 25, 28, 28, 28, 25, 28, 28, 28, 26, 28, 28, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 1, 28, 28, 28, 28, 27, 27, 28, 28, 28, 27, 28, 28, 28, & + 28, 28, 25, 28, 28, 28, 25, 28, 28, 28, 25, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 1, 28, 28, 27, 28, 28, 27, 28, 28, 28, 27, 28, 28, & + 25, 28, 28, 28, 28, 28, 28, 25, 28, 25, 28, 28, 28, 28, 28, 26, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 1, 28, 28, 28, 28, 28, 28, 27, 28, 28, 28, 27, & + 28, 25, 28, 28, 28, 25, 28, 28, 28, 28, 28, 25, 28, 28, 26, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 1, 28, 28, 28, 28, 27, 28, 28, 28, 27, 28, & + 28, 28, 25, 28, 28, 28, 25, 28, 28, 28, 25, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 1, 28, 28, 27, 28, 28, 28, 27, 28, 28, & + 28, 28, 28, 25, 25, 28, 28, 28, 25, 28, 28, 28, 26, 28, 28, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 1, 27, 28, 28, 28, 27, 28, 28, 28, & + 28, 28, 28, 25, 25, 28, 28, 28, 25, 28, 28, 28, 26, 28, 28, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 1, 28, 28, 28, 27, 28, 28, 28, & + 28, 28, 25, 28, 28, 28, 25, 28, 28, 28, 25, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 1, 28, 28, 28, 27, 28, 28, & + 28, 25, 28, 28, 28, 25, 28, 28, 28, 28, 28, 25, 28, 28, 26, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 28, 28, 1, 28, 28, 28, 27, 28, & + 25, 28, 28, 28, 28, 28, 28, 25, 28, 25, 28, 28, 28, 28, 28, 26, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 28, 28, 1, 28, 28, 28, 27, & + 28, 28, 28, 25, 25, 28, 28, 28, 25, 28, 28, 28, 26, 28, 28, 28, 28, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 27, 28, 28, 28, 1, 28, 28, 28, & + 28, 28, 25, 28, 28, 28, 25, 28, 28, 28, 25, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 27, 28, 28, 28, 1, 28, 28, & + 28, 25, 28, 28, 28, 25, 28, 28, 28, 28, 28, 25, 28, 28, 26, 28, 28, 26, 28, 28, 26, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 28, 1, 28, & + 25, 28, 28, 28, 28, 28, 28, 25, 28, 25, 28, 28, 28, 28, 28, 26, 26, 28, 28, 28, 28, 26, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 27, 28, 28, 27, 28, 28, 28, 28, 28, 28, 27, 28, 28, 28, 1 & + ],shape(BCC_INTERACTIONSLIPSLIP)) !< Slip-slip interaction types for bcc / Madec 2017 (https://doi.org/10.1016/j.actamat.2016.12.040) + !< 1: self interaction --> alpha 0 + !< 2: collinear interaction --> alpha 1 + !< 3: coplanar interaction --> alpha 2 + !< 4-7: other coefficients + !< 8: {110}-{112}, collinear and perpendicular planes --> alpha 6 + !< 9: {110}-{112}, just collinear --> alpha 7 + !< 10-24: other coefficients + !< 25: {110}-{123}, just collinear + !< 26: {112}-{123}, just collinear + !< 27: {123}-{123}, just collinear + !< 28: other interaction integer, dimension(HEX_NSLIP,HEX_NSLIP), parameter :: & HEX_INTERACTIONSLIPSLIP = reshape( [& - 1, 2, 2, 3, 3, 3, 7, 7, 7, 13,13,13,13,13,13, 21,21,21,21,21,21,21,21,21,21,21,21, 31,31,31,31,31,31, & ! -----> acting + 1, 2, 2, 3, 3, 3, 7, 7, 7, 13,13,13,13,13,13, 21,21,21,21,21,21,21,21,21,21,21,21, 31,31,31,31,31,31, & ! -----> acting (forest) 2, 1, 2, 3, 3, 3, 7, 7, 7, 13,13,13,13,13,13, 21,21,21,21,21,21,21,21,21,21,21,21, 31,31,31,31,31,31, & ! | 2, 2, 1, 3, 3, 3, 7, 7, 7, 13,13,13,13,13,13, 21,21,21,21,21,21,21,21,21,21,21,21, 31,31,31,31,31,31, & ! | ! ! v - 6, 6, 6, 4, 5, 5, 8, 8, 8, 14,14,14,14,14,14, 22,22,22,22,22,22,22,22,22,22,22,22, 32,32,32,32,32,32, & ! reacting + 6, 6, 6, 4, 5, 5, 8, 8, 8, 14,14,14,14,14,14, 22,22,22,22,22,22,22,22,22,22,22,22, 32,32,32,32,32,32, & ! reacting (primary) 6, 6, 6, 5, 4, 5, 8, 8, 8, 14,14,14,14,14,14, 22,22,22,22,22,22,22,22,22,22,22,22, 32,32,32,32,32,32, & 6, 6, 6, 5, 5, 4, 8, 8, 8, 14,14,14,14,14,14, 22,22,22,22,22,22,22,22,22,22,22,22, 32,32,32,32,32,32, & ! @@ -867,7 +898,7 @@ function lattice_interaction_SlipBySlip(Nslip,interactionValues,structure) resul 42,42,42, 41,41,41, 40,40,40, 39,39,39,39,39,39, 38,38,38,38,38,38,38,38,38,38,38,38, 37,37,37,36,37,37, & 42,42,42, 41,41,41, 40,40,40, 39,39,39,39,39,39, 38,38,38,38,38,38,38,38,38,38,38,38, 37,37,37,37,36,37, & 42,42,42, 41,41,41, 40,40,40, 39,39,39,39,39,39, 38,38,38,38,38,38,38,38,38,38,38,38, 37,37,37,37,37,36 & - ],shape(HEX_INTERACTIONSLIPSLIP)) !< Slip--slip interaction types for hex (onion peel naming scheme) + ],shape(HEX_INTERACTIONSLIPSLIP)) !< Slip-slip interaction types for hex (onion peel naming scheme) integer, dimension(BCT_NSLIP,BCT_NSLIP), parameter :: & BCT_INTERACTIONSLIPSLIP = reshape( [& @@ -1167,11 +1198,38 @@ function lattice_interaction_SlipByTwin(Nslip,Ntwin,interactionValues,structure) 3,3,2,3,3,2,3,3,1,3,3,3, & 3,3,3,2,2,3,3,3,3,1,3,3, & 2,3,3,3,3,3,3,2,3,3,1,3, & - 3,2,3,3,3,3,2,3,3,3,3,1 & + 3,2,3,3,3,3,2,3,3,3,3,1, & + ! + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4, & + 4,4,4,4,4,4,4,4,4,4,4,4 & ],shape(BCC_INTERACTIONSLIPTWIN)) !< Slip-twin interaction types for bcc !< 1: coplanar interaction !< 2: screw trace between slip system and twin habit plane (easy cross slip) !< 3: other interaction + !< 4: other interaction with slip family {123} + integer, dimension(HEX_NTWIN,HEX_NSLIP), parameter :: & HEX_INTERACTIONSLIPTWIN = reshape( [& 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, & ! ----> twin (acting) @@ -1574,7 +1632,7 @@ end function lattice_slip_normal !-------------------------------------------------------------------------------------------------- -!> @brief Transverse direction of slip systems ( || t = b x n) +!> @brief Transverse direction of slip systems (|| t = b x n) !-------------------------------------------------------------------------------------------------- function lattice_slip_transverse(Nslip,structure,cOverA) result(t) diff --git a/src/phase.f90 b/src/phase.f90 index 93743d64a..4c9f170ee 100644 --- a/src/phase.f90 +++ b/src/phase.f90 @@ -58,10 +58,6 @@ module phase grain end type tDebugOptions - integer, dimension(:), allocatable, public :: & !< ToDo: should be protected (bug in Intel compiler) - phase_elasticityInstance, & - phase_NstiffnessDegradations - logical, dimension(:), allocatable, public :: & ! ToDo: should be protected (bug in Intel Compiler) phase_localPlasticity !< flags phases with local constitutive law @@ -71,10 +67,6 @@ module phase damageState - integer, public, protected :: & - phase_plasticity_maxSizeDotState, & - phase_source_maxSizeDotState - interface ! == cleaned:begin ================================================================================= @@ -239,11 +231,16 @@ module phase logical :: converged_ end function crystallite_stress + !ToDo: Try to merge the all stiffness functions module function phase_homogenizedC(ph,en) result(C) integer, intent(in) :: ph, en real(pReal), dimension(6,6) :: C end function phase_homogenizedC - + module function phase_damage_C(C_homogenized,ph,en) result(C) + real(pReal), dimension(3,3,3,3), intent(in) :: C_homogenized + integer, intent(in) :: ph,en + real(pReal), dimension(3,3,3,3) :: C + end function phase_damage_C module function phase_f_phi(phi,co,ce) result(f) integer, intent(in) :: ce,co @@ -298,7 +295,6 @@ module phase end interface - type(tDebugOptions) :: debugConstitutive #if __INTEL_COMPILER >= 1900 public :: & @@ -377,19 +373,6 @@ subroutine phase_init call damage_init call thermal_init(phases) - - phase_source_maxSizeDotState = 0 - PhaseLoop2:do ph = 1,phases%length -!-------------------------------------------------------------------------------------------------- -! partition and initialize state - plasticState(ph)%state = plasticState(ph)%state0 - if(damageState(ph)%sizeState > 0) & - damageState(ph)%state = damageState(ph)%state0 - enddo PhaseLoop2 - - phase_source_maxSizeDotState = maxval(damageState%sizeDotState) - phase_plasticity_maxSizeDotState = maxval(plasticState%sizeDotState) - end subroutine phase_init @@ -545,8 +528,7 @@ subroutine crystallite_init() phases => config_material%get('phase') do ph = 1, phases%length - if (damageState(ph)%sizeState > 0) & - allocate(damageState(ph)%subState0,source=damageState(ph)%state0) ! ToDo: hack + if (damageState(ph)%sizeState > 0) allocate(damageState(ph)%subState0,source=damageState(ph)%state0) ! ToDo: hack enddo print'(a42,1x,i10)', ' # of elements: ', eMax @@ -560,7 +542,7 @@ subroutine crystallite_init() do ip = 1, size(material_phaseMemberAt,2) do co = 1,homogenization_Nconstituents(material_homogenizationAt(el)) call crystallite_orientations(co,ip,el) - call plastic_dependentState(co,ip,el) ! update dependent state variables to be consistent with basic states + call plastic_dependentState(co,ip,el) ! update dependent state variables to be consistent with basic states enddo enddo enddo @@ -576,9 +558,9 @@ end subroutine crystallite_init subroutine crystallite_orientations(co,ip,el) integer, intent(in) :: & - co, & !< counter in integration point component loop - ip, & !< counter in integration point loop - el !< counter in element loop + co, & !< counter in integration point component loop + ip, & !< counter in integration point loop + el !< counter in element loop call crystallite_orientation(co,ip,el)%fromMatrix(transpose(math_rotationalPart(& diff --git a/src/phase_damage.f90 b/src/phase_damage.f90 index 5a1d5970b..45898f54d 100644 --- a/src/phase_damage.f90 +++ b/src/phase_damage.f90 @@ -15,13 +15,15 @@ submodule(phase) damage DAMAGE_ANISOBRITTLE_ID end enum + integer :: phase_damage_maxSizeDotState + type :: tDataContainer real(pReal), dimension(:), allocatable :: phi, d_phi_d_dot_phi end type tDataContainer integer(kind(DAMAGE_UNDEFINED_ID)), dimension(:), allocatable :: & - phase_source !< active sources mechanisms of each phase + phase_damage !< active sources mechanisms of each phase type(tDataContainer), dimension(:), allocatable :: current @@ -126,17 +128,38 @@ module subroutine damage_init enddo - allocate(phase_source(phases%length), source = DAMAGE_UNDEFINED_ID) + allocate(phase_damage(phases%length), source = DAMAGE_UNDEFINED_ID) if (damage_active) then - where(isobrittle_init() ) phase_source = DAMAGE_ISOBRITTLE_ID - where(isoductile_init() ) phase_source = DAMAGE_ISODUCTILE_ID - where(anisobrittle_init()) phase_source = DAMAGE_ANISOBRITTLE_ID + where(isobrittle_init() ) phase_damage = DAMAGE_ISOBRITTLE_ID + where(isoductile_init() ) phase_damage = DAMAGE_ISODUCTILE_ID + where(anisobrittle_init()) phase_damage = DAMAGE_ANISOBRITTLE_ID endif + phase_damage_maxSizeDotState = maxval(damageState%sizeDotState) + end subroutine damage_init +!-------------------------------------------------------------------------------------------------- +!> @brief returns the degraded/modified elasticity matrix +!-------------------------------------------------------------------------------------------------- +module function phase_damage_C(C_homogenized,ph,en) result(C) + + real(pReal), dimension(3,3,3,3), intent(in) :: C_homogenized + integer, intent(in) :: ph,en + real(pReal), dimension(3,3,3,3) :: C + + damageType: select case (phase_damage(ph)) + case (DAMAGE_ISOBRITTLE_ID) damageType + C = C_homogenized * damage_phi(ph,en)**2 + case default damageType + C = C_homogenized + end select damageType + +end function phase_damage_C + + !---------------------------------------------------------------------------------------------- !< @brief returns local part of nonlocal damage driving force !---------------------------------------------------------------------------------------------- @@ -155,7 +178,7 @@ module function phase_f_phi(phi,co,ce) result(f) ph = material_phaseID(co,ce) en = material_phaseEntry(co,ce) - select case(phase_source(ph)) + select case(phase_damage(ph)) case(DAMAGE_ISOBRITTLE_ID,DAMAGE_ISODUCTILE_ID,DAMAGE_ANISOBRITTLE_ID) f = 1.0_pReal & - phi*damageState(ph)%state(1,en) @@ -187,9 +210,9 @@ module function integrateDamageState(dt,co,ip,el) result(broken) size_so real(pReal) :: & zeta - real(pReal), dimension(phase_source_maxSizeDotState) :: & + real(pReal), dimension(phase_damage_maxSizeDotState) :: & r ! state residuum - real(pReal), dimension(phase_source_maxSizeDotState,2) :: source_dotState + real(pReal), dimension(phase_damage_maxSizeDotState,2) :: source_dotState logical :: & converged_ @@ -275,10 +298,10 @@ module subroutine damage_results(group,ph) integer, intent(in) :: ph - if (phase_source(ph) /= DAMAGE_UNDEFINED_ID) & + if (phase_damage(ph) /= DAMAGE_UNDEFINED_ID) & call results_closeGroup(results_addGroup(group//'damage')) - sourceType: select case (phase_source(ph)) + sourceType: select case (phase_damage(ph)) case (DAMAGE_ISOBRITTLE_ID) sourceType call isobrittle_results(ph,group//'damage/') @@ -309,7 +332,7 @@ function phase_damage_collectDotState(ph,me) result(broken) if (damageState(ph)%sizeState > 0) then - sourceType: select case (phase_source(ph)) + sourceType: select case (phase_damage(ph)) case (DAMAGE_ISODUCTILE_ID) sourceType call isoductile_dotState(ph,me) @@ -376,7 +399,7 @@ function phase_damage_deltaState(Fe, ph, me) result(broken) if (damageState(ph)%sizeState == 0) return - sourceType: select case (phase_source(ph)) + sourceType: select case (phase_damage(ph)) case (DAMAGE_ISOBRITTLE_ID) sourceType call isobrittle_deltaState(phase_homogenizedC(ph,me), Fe, ph,me) diff --git a/src/phase_damage_isobrittle.f90 b/src/phase_damage_isobrittle.f90 index e88a87352..64b49b1bc 100644 --- a/src/phase_damage_isobrittle.f90 +++ b/src/phase_damage_isobrittle.f90 @@ -66,14 +66,14 @@ module function isobrittle_init() result(mySources) Nmembers = count(material_phaseID==ph) call phase_allocateState(damageState(ph),Nmembers,1,1,1) - damageState(ph)%atol = src%get_asFloat('isoBrittle_atol',defaultVal=1.0e-3_pReal) + damageState(ph)%atol = src%get_asFloat('isobrittle_atol',defaultVal=1.0e-3_pReal) if(any(damageState(ph)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isobrittle_atol' end associate !-------------------------------------------------------------------------------------------------- ! exit if any parameter is out of range - if (extmsg /= '') call IO_error(211,ext_msg=trim(extmsg)//'(damage_isoBrittle)') + if (extmsg /= '') call IO_error(211,ext_msg=trim(extmsg)//'(damage_isobrittle)') endif enddo diff --git a/src/phase_damage_isoductile.f90 b/src/phase_damage_isoductile.f90 index 997b948fe..42b7db434 100644 --- a/src/phase_damage_isoductile.f90 +++ b/src/phase_damage_isoductile.f90 @@ -70,14 +70,14 @@ module function isoductile_init() result(mySources) Nmembers=count(material_phaseID==ph) call phase_allocateState(damageState(ph),Nmembers,1,1,0) - damageState(ph)%atol = src%get_asFloat('isoDuctile_atol',defaultVal=1.0e-3_pReal) + damageState(ph)%atol = src%get_asFloat('isoductile_atol',defaultVal=1.0e-3_pReal) if(any(damageState(ph)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isoductile_atol' end associate !-------------------------------------------------------------------------------------------------- ! exit if any parameter is out of range - if (extmsg /= '') call IO_error(211,ext_msg=trim(extmsg)//'(damage_isoDuctile)') + if (extmsg /= '') call IO_error(211,ext_msg=trim(extmsg)//'(damage_isoductile)') endif enddo diff --git a/src/phase_mechanical.f90 b/src/phase_mechanical.f90 index 298ff57f3..be0bca892 100644 --- a/src/phase_mechanical.f90 +++ b/src/phase_mechanical.f90 @@ -5,10 +5,6 @@ submodule(phase) mechanical enum, bind(c); enumerator :: & - ELASTICITY_UNDEFINED_ID, & - ELASTICITY_HOOKE_ID, & - STIFFNESS_DEGRADATION_UNDEFINED_ID, & - STIFFNESS_DEGRADATION_DAMAGE_ID, & PLASTICITY_UNDEFINED_ID, & PLASTICITY_NONE_ID, & PLASTICITY_ISOTROPIC_ID, & @@ -23,11 +19,6 @@ submodule(phase) mechanical KINEMATICS_THERMAL_EXPANSION_ID end enum - integer(kind(ELASTICITY_UNDEFINED_ID)), dimension(:), allocatable :: & - phase_elasticity !< elasticity of each phase - integer(kind(STIFFNESS_DEGRADATION_UNDEFINED_ID)), dimension(:,:), allocatable :: & - phase_stiffnessDegradation !< active stiffness degradation mechanisms of each phase - type(tTensorContainer), dimension(:), allocatable :: & ! current value phase_mechanical_Fe, & @@ -50,6 +41,7 @@ submodule(phase) mechanical integer(kind(PLASTICITY_undefined_ID)), dimension(:), allocatable :: & phase_plasticity !< plasticity of each phase + integer :: phase_plasticity_maxSizeDotState interface @@ -57,9 +49,27 @@ submodule(phase) mechanical class(tNode), pointer :: phases end subroutine eigendeformation_init + module subroutine elastic_init(phases) + class(tNode), pointer :: phases + end subroutine elastic_init + module subroutine plastic_init end subroutine plastic_init + module subroutine phase_hooke_SandItsTangents(S,dS_dFe,dS_dFi,Fe,Fi,ph,en) + integer, intent(in) :: & + ph, & + en + real(pReal), intent(in), dimension(3,3) :: & + Fe, & !< elastic deformation gradient + Fi !< intermediate deformation gradient + real(pReal), intent(out), dimension(3,3) :: & + S !< 2nd Piola-Kirchhoff stress tensor in lattice configuration + real(pReal), intent(out), dimension(3,3,3,3) :: & + dS_dFe, & !< derivative of 2nd P-K stress with respect to elastic deformation gradient + dS_dFi !< derivative of 2nd P-K stress with respect to intermediate deformation gradient + end subroutine phase_hooke_SandItsTangents + module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,ph,en) real(pReal), dimension(3,3), intent(out) :: & Li !< inleastic velocity gradient @@ -73,7 +83,6 @@ submodule(phase) mechanical end subroutine plastic_isotropic_LiAndItsTangent module function plastic_dotState(subdt,co,ip,el,ph,en) result(broken) - integer, intent(in) :: & co, & !< component-ID of integration point ip, & !< integration point @@ -198,17 +207,11 @@ module subroutine mechanical_init(materials,phases) constituents, & constituent, & phase, & - mech, & - elastic, & - stiffDegradation + mech print'(/,a)', ' <<<+- phase:mechanical init -+>>>' !------------------------------------------------------------------------------------------------- -! initialize elasticity (hooke) !ToDO: Maybe move to elastic submodule along with function homogenizedC? - allocate(phase_elasticity(phases%length), source = ELASTICITY_undefined_ID) - allocate(phase_elasticityInstance(phases%length), source = 0) - allocate(phase_NstiffnessDegradations(phases%length),source=0) allocate(output_constituent(phases%length)) allocate(phase_mechanical_Fe(phases%length)) @@ -253,32 +256,8 @@ module subroutine mechanical_init(materials,phases) #else output_constituent(ph)%label = mech%get_as1dString('output',defaultVal=emptyStringArray) #endif - elastic => mech%get('elastic') - if (IO_lc(elastic%get_asString('type')) == 'hooke') then ! accept small letter h for the moment - phase_elasticity(ph) = ELASTICITY_HOOKE_ID - else - call IO_error(200,ext_msg=elastic%get_asString('type')) - endif - stiffDegradation => mech%get('stiffness_degradation',defaultVal=emptyList) ! check for stiffness degradation mechanisms - phase_NstiffnessDegradations(ph) = stiffDegradation%length enddo - allocate(phase_stiffnessDegradation(maxval(phase_NstiffnessDegradations),phases%length), & - source=STIFFNESS_DEGRADATION_undefined_ID) - - if(maxVal(phase_NstiffnessDegradations)/=0) then - do ph = 1, phases%length - phase => phases%get(ph) - mech => phase%get('mechanical') - stiffDegradation => mech%get('stiffness_degradation',defaultVal=emptyList) - do stiffDegradationCtr = 1, stiffDegradation%length - if(stiffDegradation%get_asString(stiffDegradationCtr) == 'damage') & - phase_stiffnessDegradation(stiffDegradationCtr,ph) = STIFFNESS_DEGRADATION_damage_ID - enddo - enddo - endif - - do el = 1, size(material_phaseMemberAt,3); do ip = 1, size(material_phaseMemberAt,2) do co = 1, homogenization_Nconstituents(material_homogenizationAt(el)) material => materials%get(discretization_materialAt(el)) @@ -306,6 +285,9 @@ module subroutine mechanical_init(materials,phases) enddo; enddo +! initialize elasticity + call elastic_init(phases) + ! initialize plasticity allocate(plasticState(phases%length)) allocate(phase_plasticity(phases%length),source = PLASTICITY_undefined_ID) @@ -313,9 +295,11 @@ module subroutine mechanical_init(materials,phases) call plastic_init() - do ph = 1, phases%length - phase_elasticityInstance(ph) = count(phase_elasticity(1:ph) == phase_elasticity(ph)) + do ph = 1,phases%length + plasticState(ph)%state0 = plasticState(ph)%state enddo + phase_plasticity_maxSizeDotState = maxval(plasticState%sizeDotState) + num_crystallite => config_numerics%get('crystallite',defaultVal=emptyDict) @@ -348,51 +332,6 @@ module subroutine mechanical_init(materials,phases) end subroutine mechanical_init -!-------------------------------------------------------------------------------------------------- -!> @brief returns the 2nd Piola-Kirchhoff stress tensor and its tangent with respect to -!> the elastic and intermediate deformation gradients using Hooke's law -!-------------------------------------------------------------------------------------------------- -subroutine phase_hooke_SandItsTangents(S, dS_dFe, dS_dFi, & - Fe, Fi, ph, en) - - integer, intent(in) :: & - ph, & - en - real(pReal), intent(in), dimension(3,3) :: & - Fe, & !< elastic deformation gradient - Fi !< intermediate deformation gradient - real(pReal), intent(out), dimension(3,3) :: & - S !< 2nd Piola-Kirchhoff stress tensor in lattice configuration - real(pReal), intent(out), dimension(3,3,3,3) :: & - dS_dFe, & !< derivative of 2nd P-K stress with respect to elastic deformation gradient - dS_dFi !< derivative of 2nd P-K stress with respect to intermediate deformation gradient - - real(pReal), dimension(3,3) :: E - real(pReal), dimension(3,3,3,3) :: C - integer :: & - d, & !< counter in degradation loop - i, j - - C = math_66toSym3333(phase_homogenizedC(ph,en)) - - DegradationLoop: do d = 1, phase_NstiffnessDegradations(ph) - degradationType: select case(phase_stiffnessDegradation(d,ph)) - case (STIFFNESS_DEGRADATION_damage_ID) degradationType - C = C * damage_phi(ph,en)**2 - end select degradationType - enddo DegradationLoop - - E = 0.5_pReal*(matmul(transpose(Fe),Fe)-math_I3) !< Green-Lagrange strain in unloaded configuration - S = math_mul3333xx33(C,matmul(matmul(transpose(Fi),E),Fi)) !< 2PK stress in lattice configuration in work conjugate with GL strain pulled back to lattice configuration - - do i =1, 3;do j=1,3 - dS_dFe(i,j,1:3,1:3) = matmul(Fe,matmul(matmul(Fi,C(i,j,1:3,1:3)),transpose(Fi))) !< dS_ij/dFe_kl = C_ijmn * Fi_lm * Fi_on * Fe_ko - dS_dFi(i,j,1:3,1:3) = 2.0_pReal*matmul(matmul(E,Fi),C(i,j,1:3,1:3)) !< dS_ij/dFi_kl = C_ijln * E_km * Fe_mn - enddo; enddo - -end subroutine phase_hooke_SandItsTangents - - module subroutine mechanical_results(group,ph) character(len=*), intent(in) :: group @@ -1082,26 +1021,6 @@ module subroutine mechanical_forward() end subroutine mechanical_forward - -!-------------------------------------------------------------------------------------------------- -!> @brief returns the homogenize elasticity matrix -!> ToDo: homogenizedC66 would be more consistent -!-------------------------------------------------------------------------------------------------- -module function phase_homogenizedC(ph,en) result(C) - - real(pReal), dimension(6,6) :: C - integer, intent(in) :: ph, en - - plasticType: select case (phase_plasticity(ph)) - case (PLASTICITY_DISLOTWIN_ID) plasticType - C = plastic_dislotwin_homogenizedC(ph,en) - case default plasticType - C = lattice_C66(1:6,1:6,ph) - end select plasticType - -end function phase_homogenizedC - - !-------------------------------------------------------------------------------------------------- !> @brief calculate stress (P) !-------------------------------------------------------------------------------------------------- diff --git a/src/phase_mechanical_eigen_thermalexpansion.f90 b/src/phase_mechanical_eigen_thermalexpansion.f90 index cee60d675..94d15cfb1 100644 --- a/src/phase_mechanical_eigen_thermalexpansion.f90 +++ b/src/phase_mechanical_eigen_thermalexpansion.f90 @@ -40,7 +40,6 @@ module function thermalexpansion_init(kinematics_length) result(myKinematics) print'(/,a)', ' <<<+- phase:mechanical:eigen:thermalexpansion init -+>>>' myKinematics = kinematics_active('thermalexpansion',kinematics_length) - print*, myKinematics Ninstances = count(myKinematics) print'(a,i2)', ' # phases: ',Ninstances; flush(IO_STDOUT) if(Ninstances == 0) return @@ -101,16 +100,16 @@ module subroutine thermalexpansion_LiAndItsTangent(Li, dLi_dTstar, ph,me) associate(prm => param(kinematics_thermal_expansion_instance(ph))) Li = dot_T * ( & - prm%A(1:3,1:3,1)*(T - prm%T_ref)**0 & ! constant coefficient - + prm%A(1:3,1:3,2)*(T - prm%T_ref)**1 & ! linear coefficient - + prm%A(1:3,1:3,3)*(T - prm%T_ref)**2 & ! quadratic coefficient + prm%A(1:3,1:3,1)*(T - prm%T_ref)**0 & ! constant coefficient + + prm%A(1:3,1:3,2)*(T - prm%T_ref)**1 & ! linear coefficient + + prm%A(1:3,1:3,3)*(T - prm%T_ref)**2 & ! quadratic coefficient ) / & (1.0_pReal & + prm%A(1:3,1:3,1)*(T - prm%T_ref)**1 / 1. & + prm%A(1:3,1:3,2)*(T - prm%T_ref)**2 / 2. & + prm%A(1:3,1:3,3)*(T - prm%T_ref)**3 / 3. & ) - end associate + end associate dLi_dTstar = 0.0_pReal end subroutine thermalexpansion_LiAndItsTangent diff --git a/src/phase_mechanical_elastic.f90 b/src/phase_mechanical_elastic.f90 new file mode 100644 index 000000000..e25e5e8d2 --- /dev/null +++ b/src/phase_mechanical_elastic.f90 @@ -0,0 +1,102 @@ +submodule(phase:mechanical) elastic + + enum, bind(c); enumerator :: & + ELASTICITY_UNDEFINED_ID, & + ELASTICITY_HOOKE_ID + end enum + + integer(kind(ELASTICITY_UNDEFINED_ID)), dimension(:), allocatable :: & + phase_elasticity !< elasticity of each phase + +contains + + +module subroutine elastic_init(phases) + + class(tNode), pointer :: & + phases + + integer :: & + ph + class(tNode), pointer :: & + phase, & + mech, & + elastic + + print'(/,a)', ' <<<+- phase:mechanical:elastic init -+>>>' + + allocate(phase_elasticity(phases%length), source = ELASTICITY_undefined_ID) + + do ph = 1, phases%length + phase => phases%get(ph) + mech => phase%get('mechanical') + elastic => mech%get('elastic') + if(IO_lc(elastic%get_asString('type')) == 'hooke') then ! accept small letter h for the moment + phase_elasticity(ph) = ELASTICITY_HOOKE_ID + else + call IO_error(200,ext_msg=elastic%get_asString('type')) + endif + enddo + +end subroutine elastic_init + + +!-------------------------------------------------------------------------------------------------- +!> @brief returns the 2nd Piola-Kirchhoff stress tensor and its tangent with respect to +!> the elastic and intermediate deformation gradients using Hooke's law +!-------------------------------------------------------------------------------------------------- +module subroutine phase_hooke_SandItsTangents(S, dS_dFe, dS_dFi, & + Fe, Fi, ph, en) + + integer, intent(in) :: & + ph, & + en + real(pReal), intent(in), dimension(3,3) :: & + Fe, & !< elastic deformation gradient + Fi !< intermediate deformation gradient + real(pReal), intent(out), dimension(3,3) :: & + S !< 2nd Piola-Kirchhoff stress tensor in lattice configuration + real(pReal), intent(out), dimension(3,3,3,3) :: & + dS_dFe, & !< derivative of 2nd P-K stress with respect to elastic deformation gradient + dS_dFi !< derivative of 2nd P-K stress with respect to intermediate deformation gradient + + real(pReal), dimension(3,3) :: E + real(pReal), dimension(3,3,3,3) :: C + integer :: & + d, & !< counter in degradation loop + i, j + + C = math_66toSym3333(phase_homogenizedC(ph,en)) + C = phase_damage_C(C,ph,en) + + E = 0.5_pReal*(matmul(transpose(Fe),Fe)-math_I3) !< Green-Lagrange strain in unloaded configuration + S = math_mul3333xx33(C,matmul(matmul(transpose(Fi),E),Fi)) !< 2PK stress in lattice configuration in work conjugate with GL strain pulled back to lattice configuration + + do i =1, 3;do j=1,3 + dS_dFe(i,j,1:3,1:3) = matmul(Fe,matmul(matmul(Fi,C(i,j,1:3,1:3)),transpose(Fi))) !< dS_ij/dFe_kl = C_ijmn * Fi_lm * Fi_on * Fe_ko + dS_dFi(i,j,1:3,1:3) = 2.0_pReal*matmul(matmul(E,Fi),C(i,j,1:3,1:3)) !< dS_ij/dFi_kl = C_ijln * E_km * Fe_mn + enddo; enddo + +end subroutine phase_hooke_SandItsTangents + + +!-------------------------------------------------------------------------------------------------- +!> @brief returns the homogenized elasticity matrix +!> ToDo: homogenizedC66 would be more consistent +!-------------------------------------------------------------------------------------------------- +module function phase_homogenizedC(ph,en) result(C) + + real(pReal), dimension(6,6) :: C + integer, intent(in) :: ph, en + + plasticType: select case (phase_plasticity(ph)) + case (PLASTICITY_DISLOTWIN_ID) plasticType + C = plastic_dislotwin_homogenizedC(ph,en) + case default plasticType + C = lattice_C66(1:6,1:6,ph) + end select plasticType + +end function phase_homogenizedC + + +end submodule elastic diff --git a/src/phase_mechanical_plastic_dislotungsten.f90 b/src/phase_mechanical_plastic_dislotungsten.f90 index bb53684e1..0372d7f63 100644 --- a/src/phase_mechanical_plastic_dislotungsten.f90 +++ b/src/phase_mechanical_plastic_dislotungsten.f90 @@ -254,8 +254,6 @@ module function plastic_dislotungsten_init() result(myPlasticity) allocate(dst%Lambda_sl(prm%sum_N_sl,Nmembers), source=0.0_pReal) allocate(dst%threshold_stress(prm%sum_N_sl,Nmembers), source=0.0_pReal) - plasticState(ph)%state0 = plasticState(ph)%state ! ToDo: this could be done centrally - end associate !-------------------------------------------------------------------------------------------------- @@ -414,7 +412,7 @@ module subroutine plastic_dislotungsten_results(ph,group) 'mobile dislocation density','1/m²') case('rho_dip') if(prm%sum_N_sl>0) call results_writeDataset(group,stt%rho_dip,trim(prm%output(o)), & - 'dislocation dipole density''1/m²') + 'dislocation dipole density','1/m²') case('gamma_sl') if(prm%sum_N_sl>0) call results_writeDataset(group,stt%gamma_sl,trim(prm%output(o)), & 'plastic shear','1') diff --git a/src/phase_mechanical_plastic_dislotwin.f90 b/src/phase_mechanical_plastic_dislotwin.f90 index 151250f14..fdfd630c6 100644 --- a/src/phase_mechanical_plastic_dislotwin.f90 +++ b/src/phase_mechanical_plastic_dislotwin.f90 @@ -467,8 +467,6 @@ module function plastic_dislotwin_init() result(myPlasticity) allocate(dst%tau_r_tr (prm%sum_N_tr,Nmembers),source=0.0_pReal) allocate(dst%V_tr (prm%sum_N_tr,Nmembers),source=0.0_pReal) - plasticState(ph)%state0 = plasticState(ph)%state ! ToDo: this could be done centrally - end associate !-------------------------------------------------------------------------------------------------- diff --git a/src/phase_mechanical_plastic_isotropic.f90 b/src/phase_mechanical_plastic_isotropic.f90 index 5ab73895f..0e7920f1d 100644 --- a/src/phase_mechanical_plastic_isotropic.f90 +++ b/src/phase_mechanical_plastic_isotropic.f90 @@ -140,8 +140,6 @@ module function plastic_isotropic_init() result(myPlasticity) ! global alias plasticState(ph)%slipRate => plasticState(ph)%dotState(2:2,:) - plasticState(ph)%state0 = plasticState(ph)%state ! ToDo: this could be done centrally - end associate !-------------------------------------------------------------------------------------------------- diff --git a/src/phase_mechanical_plastic_kinehardening.f90 b/src/phase_mechanical_plastic_kinehardening.f90 index 18f467617..2e9d713f2 100644 --- a/src/phase_mechanical_plastic_kinehardening.f90 +++ b/src/phase_mechanical_plastic_kinehardening.f90 @@ -213,8 +213,6 @@ module function plastic_kinehardening_init() result(myPlasticity) stt%gamma0 => plasticState(ph)%state (startIndex :endIndex ,:) dlt%gamma0 => plasticState(ph)%deltaState(startIndex-o:endIndex-o,:) - plasticState(ph)%state0 = plasticState(ph)%state ! ToDo: this could be done centrally - end associate !-------------------------------------------------------------------------------------------------- diff --git a/src/phase_mechanical_plastic_nonlocal.f90 b/src/phase_mechanical_plastic_nonlocal.f90 index 05ba371ce..9250ccc83 100644 --- a/src/phase_mechanical_plastic_nonlocal.f90 +++ b/src/phase_mechanical_plastic_nonlocal.f90 @@ -1767,22 +1767,22 @@ end subroutine kinetics !> @brief returns copy of current dislocation densities from state !> @details raw values is rectified !-------------------------------------------------------------------------------------------------- -pure function getRho(ph,en) +pure function getRho(ph,en) result(rho) integer, intent(in) :: ph, en - real(pReal), dimension(param(ph)%sum_N_sl,10) :: getRho + real(pReal), dimension(param(ph)%sum_N_sl,10) :: rho associate(prm => param(ph)) - getRho = reshape(state(ph)%rho(:,en),[prm%sum_N_sl,10]) + rho = reshape(state(ph)%rho(:,en),[prm%sum_N_sl,10]) ! ensure positive densities (not for imm, they have a sign) - getRho(:,mob) = max(getRho(:,mob),0.0_pReal) - getRho(:,dip) = max(getRho(:,dip),0.0_pReal) + rho(:,mob) = max(rho(:,mob),0.0_pReal) + rho(:,dip) = max(rho(:,dip),0.0_pReal) - where(abs(getRho) < max(prm%rho_min/geom(ph)%V_0(en)**(2.0_pReal/3.0_pReal),prm%rho_significant)) & - getRho = 0.0_pReal + where(abs(rho) < max(prm%rho_min/geom(ph)%V_0(en)**(2.0_pReal/3.0_pReal),prm%rho_significant)) & + rho = 0.0_pReal end associate @@ -1793,22 +1793,22 @@ end function getRho !> @brief returns copy of current dislocation densities from state !> @details raw values is rectified !-------------------------------------------------------------------------------------------------- -pure function getRho0(ph,en) +pure function getRho0(ph,en) result(rho0) integer, intent(in) :: ph, en - real(pReal), dimension(param(ph)%sum_N_sl,10) :: getRho0 + real(pReal), dimension(param(ph)%sum_N_sl,10) :: rho0 associate(prm => param(ph)) - getRho0 = reshape(state0(ph)%rho(:,en),[prm%sum_N_sl,10]) + rho0 = reshape(state0(ph)%rho(:,en),[prm%sum_N_sl,10]) ! ensure positive densities (not for imm, they have a sign) - getRho0(:,mob) = max(getRho0(:,mob),0.0_pReal) - getRho0(:,dip) = max(getRho0(:,dip),0.0_pReal) + rho0(:,mob) = max(rho0(:,mob),0.0_pReal) + rho0(:,dip) = max(rho0(:,dip),0.0_pReal) - where (abs(getRho0) < max(prm%rho_min/geom(ph)%V_0(en)**(2.0_pReal/3.0_pReal),prm%rho_significant)) & - getRho0 = 0.0_pReal + where (abs(rho0) < max(prm%rho_min/geom(ph)%V_0(en)**(2.0_pReal/3.0_pReal),prm%rho_significant)) & + rho0 = 0.0_pReal end associate diff --git a/src/phase_mechanical_plastic_phenopowerlaw.f90 b/src/phase_mechanical_plastic_phenopowerlaw.f90 index b2db9e346..13335df4c 100644 --- a/src/phase_mechanical_plastic_phenopowerlaw.f90 +++ b/src/phase_mechanical_plastic_phenopowerlaw.f90 @@ -265,8 +265,6 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) plasticState(ph)%atol(startIndex:endIndex) = pl%get_asFloat('atol_gamma',defaultVal=1.0e-6_pReal) if(any(plasticState(ph)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_gamma' - plasticState(ph)%state0 = plasticState(ph)%state ! ToDo: this could be done centrally - end associate !-------------------------------------------------------------------------------------------------- diff --git a/src/rotations.f90 b/src/rotations.f90 index 171db2466..8d7b842fc 100644 --- a/src/rotations.f90 +++ b/src/rotations.f90 @@ -178,8 +178,7 @@ subroutine fromQuaternion(self,qu) class(rotation), intent(out) :: self real(pReal), dimension(4), intent(in) :: qu - if (dNeq(norm2(qu),1.0_pReal)) & - call IO_error(402,ext_msg='fromQuaternion') + if (dNeq(norm2(qu),1.0_pReal,1.0e-8_pReal)) call IO_error(402,ext_msg='fromQuaternion') self%q = qu @@ -639,7 +638,7 @@ function om2ax(om) result(ax) else call dgeev('N','V',3,om_,3,Wr,Wi,devNull,3,VR,3,work,size(work,1),ierr) if (ierr /= 0) error stop 'LAPACK error' -#if defined(__GFORTRAN__) && __GNUC__<9 || defined(__INTEL_COMPILER) && INTEL_COMPILER<1800 || defined(__PGI) +#if defined(__GFORTRAN__) && __GNUC__<9 i = maxloc(merge(1,0,cEq(cmplx(Wr,Wi,pReal),cmplx(1.0_pReal,0.0_pReal,pReal),tol=1.0e-14_pReal)),dim=1) #else i = findloc(cEq(cmplx(Wr,Wi,pReal),cmplx(1.0_pReal,0.0_pReal,pReal),tol=1.0e-14_pReal),.true.,dim=1) !find eigenvalue (1,0) @@ -1394,7 +1393,7 @@ end function conjugate_quaternion !-------------------------------------------------------------------------------------------------- !> @brief Check correctness of some rotations functions. !-------------------------------------------------------------------------------------------------- -subroutine selfTest +subroutine selfTest() type(rotation) :: R real(pReal), dimension(4) :: qu, ax, ro @@ -1405,7 +1404,7 @@ subroutine selfTest integer :: i - do i = 1, 10 + do i = 1, 20 #if defined(__GFORTRAN__) && __GNUC__<9 if(i<7) cycle @@ -1433,6 +1432,7 @@ subroutine selfTest sin(2.0_pReal*PI*x(1))*A] if(qu(1)<0.0_pReal) qu = qu * (-1.0_pReal) endif + if(.not. quaternion_equal(om2qu(qu2om(qu)),qu)) error stop 'om2qu/qu2om' if(.not. quaternion_equal(eu2qu(qu2eu(qu)),qu)) error stop 'eu2qu/qu2eu' if(.not. quaternion_equal(ax2qu(qu2ax(qu)),qu)) error stop 'ax2qu/qu2ax' @@ -1479,6 +1479,8 @@ subroutine selfTest if(all(dNeq(R%rotTensor4(R%rotTensor4(t3333),active=.true.),t3333,1.0e-12_pReal))) & error stop 'rotTensor4' + call R%fromQuaternion(qu * (1.0_pReal + merge(+5.e-9_pReal,-5.e-9_pReal, mod(i,2) == 0))) ! allow reasonable tolerance for ASCII/YAML + enddo contains