From bf7efad250ed74e6c8301d9d514d2e55ecc68b2d Mon Sep 17 00:00:00 2001 From: Mingxuan Lin Date: Wed, 23 May 2018 19:51:11 +0200 Subject: [PATCH 01/27] Fix errors in CMake script --- CMakeLists.txt | 6 ++++-- src/CMakeLists.txt | 4 +++- src/homogenization_RGC.f90 | 0 3 files changed, 7 insertions(+), 3 deletions(-) mode change 100755 => 100644 src/homogenization_RGC.f90 diff --git a/CMakeLists.txt b/CMakeLists.txt index cd33a4b31..2507777df 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -113,6 +113,8 @@ elseif ("${DAMASK_SOLVER}" STREQUAL "FEM") project (DAMASK_FEM Fortran C) add_definitions (-DFEM) message ("Building FEM Solver\n") +else () + message (FATAL_ERROR "DAMASK_SOLVER is not defined!") endif () # set linker commands (needs to be done after defining the project) @@ -188,7 +190,7 @@ set (DAMASK_INCLUDE_FLAGS "${DAMASK_INCLUDE_FLAGS} ${PETSC_INCLUDES}") ################################################################################################### # Intel Compiler ################################################################################################### -if (${CMAKE_Fortran_COMPILER_ID} STREQUAL "Intel") +if ( CMAKE_Fortran_COMPILER_ID STREQUAL "Intel") if (OPENMP) set (OPENMP_FLAGS "-qopenmp -parallel") @@ -308,7 +310,7 @@ if (${CMAKE_Fortran_COMPILER_ID} STREQUAL "Intel") ################################################################################################### # GNU Compiler ################################################################################################### -elseif(${CMAKE_Fortran_COMPILER_ID} STREQUAL "GNU") +elseif( CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") if (OPENMP) set (OPENMP_FLAGS "-fopenmp") diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 9789ec67d..e5865b829 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -1,5 +1,5 @@ # special flags for some files -if (${CMAKE_Fortran_COMPILER_ID} STREQUAL "GNU") +if (CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") SET_SOURCE_FILES_PROPERTIES( "lattice.f90" PROPERTIES COMPILE_FLAGS "-ffree-line-length-240") # long lines for interaction matrix @@ -22,6 +22,8 @@ if ("${PROJECT_NAME}" STREQUAL "DAMASK_spectral") add_library(DAMASK_INTERFACE OBJECT "spectral_interface.f90") elseif ("${PROJECT_NAME}" STREQUAL "DAMASK_FEM") add_library(DAMASK_INTERFACE OBJECT "FEM_interface.f90") +else () + message ( FATAL_ERROR "PROJECT_NAME is not defined!") endif() add_dependencies(DAMASK_INTERFACE PREC) list(APPEND OBJECTFILES $) diff --git a/src/homogenization_RGC.f90 b/src/homogenization_RGC.f90 old mode 100755 new mode 100644 From 0ee34d608cde6f783441693d69425e4bb743f57a Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 24 May 2018 00:06:03 +0200 Subject: [PATCH 02/27] fixing all appaerances --- CMakeLists.txt | 4 ++-- src/CMakeLists.txt | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 2507777df..e3714b6e6 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -190,7 +190,7 @@ set (DAMASK_INCLUDE_FLAGS "${DAMASK_INCLUDE_FLAGS} ${PETSC_INCLUDES}") ################################################################################################### # Intel Compiler ################################################################################################### -if ( CMAKE_Fortran_COMPILER_ID STREQUAL "Intel") +if ("${CMAKE_Fortran_COMPILER_ID}" STREQUAL "Intel") if (OPENMP) set (OPENMP_FLAGS "-qopenmp -parallel") @@ -310,7 +310,7 @@ if ( CMAKE_Fortran_COMPILER_ID STREQUAL "Intel") ################################################################################################### # GNU Compiler ################################################################################################### -elseif( CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") +elseif("${CMAKE_Fortran_COMPILER_ID}" STREQUAL "GNU") if (OPENMP) set (OPENMP_FLAGS "-fopenmp") diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index e5865b829..656f3edd4 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -1,5 +1,5 @@ # special flags for some files -if (CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") +if ("${CMAKE_Fortran_COMPILER_ID}" STREQUAL "GNU") SET_SOURCE_FILES_PROPERTIES( "lattice.f90" PROPERTIES COMPILE_FLAGS "-ffree-line-length-240") # long lines for interaction matrix From 1c75a2e9cd770981eca48c24655437e52de3c363 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 24 May 2018 07:13:12 +0200 Subject: [PATCH 03/27] using syntax with better error handling --- CMakeLists.txt | 44 ++++++++++++++++++++++---------------------- src/CMakeLists.txt | 16 ++++++++-------- 2 files changed, 30 insertions(+), 30 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index e3714b6e6..5f0ebea41 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -5,7 +5,7 @@ cmake_minimum_required (VERSION 2.8.8 FATAL_ERROR) #--------------------------------------------------------------------------------------- # Find PETSc from system environment set(PETSC_DIR $ENV{PETSC_DIR}) -if ("${PETSC_DIR}" STREQUAL "") +if (PETSC_DIR STREQUAL "") message (FATAL_ERROR "PETSC_DIR is not defined") endif () @@ -105,11 +105,11 @@ set (CMAKE_C_COMPILER "${PETSC_MPICC}") # Now start to care about DAMASK # DAMASK solver defines project to build -if ("${DAMASK_SOLVER}" STREQUAL "SPECTRAL") +if (DAMASK_SOLVER STREQUAL "SPECTRAL") project (DAMASK_spectral Fortran C) add_definitions (-DSpectral) message ("Building Spectral Solver\n") -elseif ("${DAMASK_SOLVER}" STREQUAL "FEM") +elseif (DAMASK_SOLVER STREQUAL "FEM") project (DAMASK_FEM Fortran C) add_definitions (-DFEM) message ("Building FEM Solver\n") @@ -120,39 +120,39 @@ endif () # set linker commands (needs to be done after defining the project) set (CMAKE_LINKER "${PETSC_LINKER}") -if ("${CMAKE_BUILD_TYPE}" STREQUAL "") +if (CMAKE_BUILD_TYPE STREQUAL "") set (CMAKE_BUILD_TYPE "RELEASE") endif () # Predefined sets for OPTIMIZATION/OPENMP based on BUILD_TYPE -if ("${CMAKE_BUILD_TYPE}" STREQUAL "DEBUG" OR "${CMAKE_BUILD_TYPE}" STREQUAL "SYNTAXONLY" ) +if (CMAKE_BUILD_TYPE STREQUAL "DEBUG" OR CMAKE_BUILD_TYPE STREQUAL "SYNTAXONLY") set (DEBUG_FLAGS "${DEBUG_FLAGS} -DDEBUG") set (PARALLEL "OFF") set (OPTI "OFF") -elseif ("${CMAKE_BUILD_TYPE}" STREQUAL "RELEASE") +elseif (CMAKE_BUILD_TYPE STREQUAL "RELEASE") set (PARALLEL "ON") set (OPTI "DEFENSIVE") -elseif ("${CMAKE_BUILD_TYPE}" STREQUAL "PERFORMANCE") +elseif (CMAKE_BUILD_TYPE STREQUAL "PERFORMANCE") set (PARALLEL "ON") set (OPTI "AGGRESSIVE") endif () # $OPTIMIZATION takes precedence over $BUILD_TYPE defaults -if ("${OPTIMIZATION}" STREQUAL "") +if (OPTIMIZATION STREQUAL "") set (OPTIMIZATION "${OPTI}") else () set (OPTIMIZATION "${OPTIMIZATION}") endif () # $OPENMP takes precedence over $BUILD_TYPE defaults -if ("${OPENMP}" STREQUAL "") +if (OPENMP STREQUAL "") set (OPENMP "${PARALLEL}") else () set(OPENMP "${OPENMP}") endif () # syntax check only (mainly for pre-receive hook, works only with gfortran) -if ("${CMAKE_BUILD_TYPE}" STREQUAL "SYNTAXONLY" ) +if (CMAKE_BUILD_TYPE STREQUAL "SYNTAXONLY") set (BUILDCMD_POST "${BUILDCMD_POST} -fsyntax-only") endif () @@ -190,17 +190,17 @@ set (DAMASK_INCLUDE_FLAGS "${DAMASK_INCLUDE_FLAGS} ${PETSC_INCLUDES}") ################################################################################################### # Intel Compiler ################################################################################################### -if ("${CMAKE_Fortran_COMPILER_ID}" STREQUAL "Intel") +if (CMAKE_Fortran_COMPILER_ID STREQUAL "Intel") if (OPENMP) set (OPENMP_FLAGS "-qopenmp -parallel") endif () - if ("${OPTIMIZATION}" STREQUAL "OFF") + if (OPTIMIZATION STREQUAL "OFF") set (OPTIMIZATION_FLAGS "-O0 -no-ip") - elseif ("${OPTIMIZATION}" STREQUAL "DEFENSIVE") + elseif (OPTIMIZATION STREQUAL "DEFENSIVE") set (OPTIMIZATION_FLAGS "-O2") - elseif ("${OPTIMIZATION}" STREQUAL "AGGRESSIVE") + elseif (OPTIMIZATION STREQUAL "AGGRESSIVE") set (OPTIMIZATION_FLAGS "-ipo -O3 -no-prec-div -fp-model fast=2 -xHost") # -fast = -ipo, -O3, -no-prec-div, -static, -fp-model fast=2, and -xHost" endif () @@ -310,17 +310,17 @@ if ("${CMAKE_Fortran_COMPILER_ID}" STREQUAL "Intel") ################################################################################################### # GNU Compiler ################################################################################################### -elseif("${CMAKE_Fortran_COMPILER_ID}" STREQUAL "GNU") +elseif(CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") if (OPENMP) set (OPENMP_FLAGS "-fopenmp") endif () - if ("${OPTIMIZATION}" STREQUAL "OFF") + if (OPTIMIZATION STREQUAL "OFF") set (OPTIMIZATION_FLAGS "-O0" ) - elseif ("${OPTIMIZATION}" STREQUAL "DEFENSIVE") + elseif (OPTIMIZATION STREQUAL "DEFENSIVE") set (OPTIMIZATION_FLAGS "-O2") - elseif ("${OPTIMIZATION}" STREQUAL "AGGRESSIVE") + elseif (OPTIMIZATION STREQUAL "AGGRESSIVE") set (OPTIMIZATION_FLAGS "-O3 -ffast-math -funroll-loops -ftree-vectorize") endif () @@ -450,7 +450,7 @@ endif () 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_LINKER} ${OPENMP_FLAGS} ${OPTIMIZATION_FLAGS} ${LINKER_FLAGS}") -if ("${CMAKE_BUILD_TYPE}" STREQUAL "DEBUG") +if (CMAKE_BUILD_TYPE STREQUAL "DEBUG") set (CMAKE_Fortran_FLAGS_${CMAKE_BUILD_TYPE} "${CMAKE_Fortran_FLAGS_${CMAKE_BUILD_TYPE}} ${DEBUG_FLAGS}") set (CMAKE_Fortran_LINK_EXECUTABLE "${CMAKE_Fortran_LINK_EXECUTABLE} ${DEBUG_FLAGS}") endif () @@ -466,15 +466,15 @@ message ("Fortran Linker Command:\n${CMAKE_Fortran_LINK_EXECUTABLE}\n") add_subdirectory (src) # INSTALL BUILT BINARIES -if ("${CMAKE_BUILD_TYPE}" STREQUAL "SYNTAXONLY") +if (CMAKE_BUILD_TYPE STREQUAL "SYNTAXONLY") exec_program (mktemp ARGS -d OUTPUT_VARIABLE BLACK_HOLE) install (PROGRAMS ${PROJECT_BINARY_DIR}/src/prec.mod DESTINATION ${BLACK_HOLE}) else () - if ("${PROJECT_NAME}" STREQUAL "DAMASK_spectral") + if (PROJECT_NAME STREQUAL "DAMASK_spectral") install (PROGRAMS ${PROJECT_BINARY_DIR}/src/DAMASK_spectral DESTINATION ${CMAKE_INSTALL_PREFIX}) - elseif ("${PROJECT_NAME}" STREQUAL "DAMASK_FEM") + elseif (PROJECT_NAME STREQUAL "DAMASK_FEM") install (PROGRAMS ${PROJECT_BINARY_DIR}/src/DAMASK_FEM DESTINATION ${CMAKE_INSTALL_PREFIX}) endif () diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 656f3edd4..f0abd9707 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -1,5 +1,5 @@ # special flags for some files -if ("${CMAKE_Fortran_COMPILER_ID}" STREQUAL "GNU") +if (CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") SET_SOURCE_FILES_PROPERTIES( "lattice.f90" PROPERTIES COMPILE_FLAGS "-ffree-line-length-240") # long lines for interaction matrix @@ -18,9 +18,9 @@ add_library(PREC OBJECT "prec.f90") add_dependencies(PREC SYSTEM_ROUTINES) list(APPEND OBJECTFILES $) -if ("${PROJECT_NAME}" STREQUAL "DAMASK_spectral") +if (PROJECT_NAME STREQUAL "DAMASK_spectral") add_library(DAMASK_INTERFACE OBJECT "spectral_interface.f90") -elseif ("${PROJECT_NAME}" STREQUAL "DAMASK_FEM") +elseif (PROJECT_NAME STREQUAL "DAMASK_FEM") add_library(DAMASK_INTERFACE OBJECT "FEM_interface.f90") else () message ( FATAL_ERROR "PROJECT_NAME is not defined!") @@ -49,11 +49,11 @@ add_dependencies(DAMASK_MATH FEsolving) list(APPEND OBJECTFILES $) # SPECTRAL solver and FEM solver use different mesh files -if ("${PROJECT_NAME}" STREQUAL "DAMASK_spectral") +if (PROJECT_NAME STREQUAL "DAMASK_spectral") add_library(MESH OBJECT "mesh.f90") add_dependencies(MESH DAMASK_MATH) list(APPEND OBJECTFILES $) -elseif ("${PROJECT_NAME}" STREQUAL "DAMASK_FEM") +elseif (PROJECT_NAME STREQUAL "DAMASK_FEM") add_library(FEZoo OBJECT "FEZoo.f90") add_dependencies(FEZoo DAMASK_MATH) list(APPEND OBJECTFILES $) @@ -160,7 +160,7 @@ add_library(DAMASK_CPFE OBJECT "CPFEM2.f90") add_dependencies(DAMASK_CPFE DAMASK_ENGINE) list(APPEND OBJECTFILES $) -if ("${PROJECT_NAME}" STREQUAL "DAMASK_spectral") +if (PROJECT_NAME STREQUAL "DAMASK_spectral") add_library(SPECTRAL_UTILITIES OBJECT "spectral_utilities.f90") add_dependencies(SPECTRAL_UTILITIES DAMASK_CPFE) list(APPEND OBJECTFILES $) @@ -172,13 +172,13 @@ if ("${PROJECT_NAME}" STREQUAL "DAMASK_spectral") "spectral_mech_Basic.f90") add_dependencies(SPECTRAL_SOLVER SPECTRAL_UTILITIES) list(APPEND OBJECTFILES $) - if(NOT "${CMAKE_BUILD_TYPE}" STREQUAL "SYNTAXONLY") + if(NOT CMAKE_BUILD_TYPE STREQUAL "SYNTAXONLY") add_executable(DAMASK_spectral "DAMASK_spectral.f90" ${OBJECTFILES}) else() add_library(DAMASK_spectral OBJECT "DAMASK_spectral.f90") endif() add_dependencies(DAMASK_spectral SPECTRAL_SOLVER) -elseif ("${PROJECT_NAME}" STREQUAL "DAMASK_FEM") +elseif (PROJECT_NAME STREQUAL "DAMASK_FEM") add_library(FEM_UTILITIES OBJECT "FEM_utilities.f90") add_dependencies(FEM_UTILITIES DAMASK_CPFE) From f0309ac4ddd3332e61bac84445cf9e05adfd89a0 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 24 May 2018 07:51:16 +0200 Subject: [PATCH 04/27] Bessel functions not used --- src/math.f90 | 131 --------------------------------------------------- 1 file changed, 131 deletions(-) diff --git a/src/math.f90 b/src/math.f90 index 8f228582b..ea3754869 100644 --- a/src/math.f90 +++ b/src/math.f90 @@ -2632,135 +2632,4 @@ real(pReal) pure function math_limit(a, left, right) end function math_limit - -!-------------------------------------------------------------------------------------------------- -!> @brief Modified Bessel I function of order 0 -!> @author John Burkardt -!> @details original version available on https://people.sc.fsu.edu/~jburkardt/f_src/toms715/toms715.html -!-------------------------------------------------------------------------------------------------- -real(pReal) function bessel_i0 (x) - use, intrinsic :: IEEE_ARITHMETIC - - implicit none - real(pReal), intent(in) :: x - integer(pInt) :: i - real(pReal) :: sump_p, sump_q, xAbs, xx - real(pReal), parameter, dimension(15) :: p_small = real( & - [-5.2487866627945699800e-18, -1.5982226675653184646e-14, -2.6843448573468483278e-11, & - -3.0517226450451067446e-08, -2.5172644670688975051e-05, -1.5453977791786851041e-02, & - -7.0935347449210549190e+00, -2.4125195876041896775e+03, -5.9545626019847898221e+05, & - -1.0313066708737980747e+08, -1.1912746104985237192e+10, -8.4925101247114157499e+11, & - -3.2940087627407749166e+13, -5.5050369673018427753e+14, -2.2335582639474375249e+15], pReal) - real(pReal), parameter, dimension(5) :: q_small = real( & - [-3.7277560179962773046e+03, 6.5158506418655165707e+06, -6.5626560740833869295e+09, & - 3.7604188704092954661e+12, -9.7087946179594019126e+14], pReal) - real(pReal), parameter, dimension(8) :: p_large = real( & - [-3.9843750000000000000e-01, 2.9205384596336793945e+00, -2.4708469169133954315e+00, & - 4.7914889422856814203e-01, -3.7384991926068969150e-03, -2.6801520353328635310e-03, & - 9.9168777670983678974e-05, -2.1877128189032726730e-06], pReal) - real(pReal), parameter, dimension(7) :: q_large = real( & - [-3.1446690275135491500e+01, 8.5539563258012929600e+01, -6.0228002066743340583e+01, & - 1.3982595353892851542e+01, -1.1151759188741312645e+00, 3.2547697594819615062e-02, & - -5.5194330231005480228e-04], pReal) - - - xAbs = abs(x) - - argRange: if (xAbs < 5.55e-17_pReal) then - bessel_i0 = 1.0_pReal - else if (xAbs < 15.0_pReal) then argRange - xx = xAbs**2.0_pReal - sump_p = p_small(1) - do i = 2, 15 - sump_p = sump_p * xx + p_small(i) - end do - xx = xx - 225.0_pReal - sump_q = ((((xx+q_small(1))*xx+q_small(2))*xx+q_small(3))*xx+q_small(4))*xx+q_small(5) - bessel_i0 = sump_p / sump_q - else if (xAbs <= 713.986_pReal) then argRange - xx = 1.0_pReal / xAbs - 2.0_pReal/30.0_pReal - sump_p = ((((((p_large(1)*xx+p_large(2))*xx+p_large(3))*xx+p_large(4))*xx+ & - p_large(5))*xx+p_large(6))*xx+p_large(7))*xx+p_large(8) - sump_q = ((((((xx+q_large(1))*xx+q_large(2))*xx+q_large(3))*xx+ & - q_large(4))*xx+q_large(5))*xx+q_large(6))*xx+q_large(7) - bessel_i0 = sump_p / sump_q - - avoidOverflow: if (xAbs > 698.986_pReal) then - bessel_i0 = ((bessel_i0*exp(xAbs-40.0_pReal)-p_large(1)*exp(xAbs-40.0_pReal))/sqrt(xAbs))*exp(40.0) - else avoidOverflow - bessel_i0 = ((bessel_i0*exp(xAbs)-p_large(1)*exp(xAbs))/sqrt(xAbs)) - endif avoidOverflow - - else argRange - bessel_i0 = IEEE_value(bessel_i0,IEEE_positive_inf) - end if argRange - -end function bessel_i0 - - -!-------------------------------------------------------------------------------------------------- -!> @brief Modified Bessel I function of order 1 -!> @author John Burkardt -!> @details original version available on https://people.sc.fsu.edu/~jburkardt/f_src/toms715/toms715.html -!-------------------------------------------------------------------------------------------------- -real(pReal) function bessel_i1 (x) - use, intrinsic :: IEEE_ARITHMETIC - - implicit none - real(pReal), intent(in) :: x - integer(pInt) :: i - real(pReal) :: sump_p, sump_q, xAbs, xx - real(pReal), dimension(15), parameter :: p_small = real( & - [-1.9705291802535139930e-19, -6.5245515583151902910e-16, -1.1928788903603238754e-12, & - -1.4831904935994647675e-09, -1.3466829827635152875e-06, -9.1746443287817501309e-04, & - -4.7207090827310162436e-01, -1.8225946631657315931e+02, -5.1894091982308017540e+04, & - -1.0588550724769347106e+07, -1.4828267606612366099e+09, -1.3357437682275493024e+11, & - -6.9876779648010090070e+12, -1.7732037840791591320e+14, -1.4577180278143463643e+15], pReal) - real(pReal), dimension(5), parameter :: q_small = real( & - [-4.0076864679904189921e+03, 7.4810580356655069138e+06, -8.0059518998619764991e+09, & - 4.8544714258273622913e+12, -1.3218168307321442305e+15], pReal) - real(pReal), dimension(8), parameter :: p_large = real( & - [-6.0437159056137600000e-02, 4.5748122901933459000e-01, -4.2843766903304806403e-01, & - 9.7356000150886612134e-02, -3.2457723974465568321e-03, -3.6395264712121795296e-04, & - 1.6258661867440836395e-05, -3.6347578404608223492e-07], pReal) - real(pReal), dimension(6), parameter :: q_large = real( & - [-3.8806586721556593450e+00, 3.2593714889036996297e+00, -8.5017476463217924408e-01, & - 7.4212010813186530069e-02, -2.2835624489492512649e-03, 3.7510433111922824643e-05], pReal) - real(pReal), parameter :: pbar = 3.98437500e-01 - - - xAbs = abs(x) - - argRange: if (xAbs < 5.55e-17_pReal) then - bessel_i1 = 0.5_pReal * xAbs - else if (xAbs < 15.0_pReal) then argRange - xx = xAbs**2.0_pReal - sump_p = p_small(1) - do i = 2, 15 - sump_p = sump_p * xx + p_small(i) - end do - xx = xx - 225.0_pReal - sump_q = ((((xx+q_small(1))*xx+q_small(2))*xx+q_small(3))*xx+q_small(4)) * xx + q_small(5) - bessel_i1 = (sump_p / sump_q) * xAbs - else if (xAbs <= 713.986_pReal) then argRange - xx = 1.0_pReal / xAbs - 2.0_pReal/30.0_pReal - sump_p = ((((((p_large(1)*xx+p_large(2))*xx+p_large(3))*xx+p_large(4))*xx+& - p_large(5))*xx+p_large(6))*xx+p_large(7))*xx+p_large(8) - sump_q = (((((xx+q_large(1))*xx+q_large(2))*xx+q_large(3))*xx+ q_large(4))*xx+q_large(5))*xx+q_large(6) - bessel_i1 = sump_p / sump_q - - avoidOverflow: if (xAbs > 698.986_pReal) then - bessel_i1 = ((bessel_i1 * exp(xAbs-40.0_pReal) + pbar * exp(xAbs-40.0_pReal)) / sqrt(xAbs)) * exp(40.0_pReal) - else avoidOverflow - bessel_i1 = ((bessel_i1 * exp(xAbs) + pbar * exp(xAbs)) / sqrt(xAbs)) - endif avoidOverflow - - else argRange - bessel_i1 = IEEE_value(bessel_i1,IEEE_positive_inf) - end if argRange - - if (x < 0.0_pReal) bessel_i1 = -bessel_i1 - -end function bessel_i1 - end module math From 1bad719abec3e9ef053a9a80c608f58589e19f02 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 24 May 2018 07:51:43 +0200 Subject: [PATCH 05/27] undefined variables can be removed, style as in other plastic modules --- src/plastic_kinematichardening.f90 | 45 +++++++++++++----------------- 1 file changed, 20 insertions(+), 25 deletions(-) diff --git a/src/plastic_kinematichardening.f90 b/src/plastic_kinematichardening.f90 index c33a14db6..553d8b6cb 100644 --- a/src/plastic_kinematichardening.f90 +++ b/src/plastic_kinematichardening.f90 @@ -177,8 +177,6 @@ subroutine plastic_kinehardening_init(fileUnit) tag = '', & line = '', & extmsg = '' - character(len=64) :: & - outputtag = '' write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_KINEHARDENING_label//' init -+>>>' write(6,'(a15,a)') ' Current time: ',IO_timeStamp() @@ -241,34 +239,42 @@ subroutine plastic_kinehardening_init(fileUnit) tag = IO_lc(IO_stringValue(line,chunkPos,1_pInt)) ! extract key select case(tag) case ('(output)') - outputtag = IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - output_ID = undefined_ID - select case(outputtag) + plastic_kinehardening_Noutput(instance) = plastic_kinehardening_Noutput(instance) + 1_pInt + plastic_kinehardening_output(plastic_kinehardening_Noutput(instance),instance) = & + IO_lc(IO_stringValue(line,chunkPos,2_pInt)) + select case(IO_lc(IO_stringValue(line,chunkPos,2_pInt))) case ('resistance') output_ID = crss_ID + case ('backstress') output_ID = crss_back_ID + case ('sense') output_ID = sense_ID + case ('chi0') output_ID = chi0_ID + case ('gamma0') output_ID = gamma0_ID + case ('accumulatedshear') output_ID = accshear_ID + case ('totalshear') output_ID = sumGamma_ID + case ('shearrate') output_ID = shearrate_ID + case ('resolvedstress') output_ID = resolvedstress_ID + + case default + plastic_kinehardening_Noutput(instance) = plastic_kinehardening_Noutput(instance) - 1_pInt ! correct for invalid + end select - if (output_ID /= undefined_ID) then - plastic_kinehardening_Noutput(instance) = plastic_kinehardening_Noutput(instance) + 1_pInt - plastic_kinehardening_output(plastic_kinehardening_Noutput(instance),instance) = outputtag - param(instance)%outputID (plastic_kinehardening_Noutput(instance)) = output_ID - endif !-------------------------------------------------------------------------------------------------- ! parameters depending on number of slip families case ('nslip') @@ -619,7 +625,6 @@ subroutine plastic_kinehardening_LpAndItsTangent(Lp,dLp_dTstar99, & math_transpose33 use lattice, only: & lattice_Sslip, & !< schmid matrix - lattice_Sslip_v, & lattice_maxNslipFamily, & lattice_NslipSystem, & lattice_NnonSchmid @@ -739,8 +744,6 @@ subroutine plastic_kinehardening_deltaState(Tstar_v,ipc,ip,el) ipc, & !< component-ID of integration point ip, & !< integration point el !< element - real(pReal), dimension(6) :: & - Tstar_dev_v !< deviatoric 2nd Piola Kirchhoff stress tensor in Mandel notation real(pReal), dimension(plastic_kinehardening_totalNslip(phase_plasticityInstance(phaseAt(ipc,ip,el)))) :: & gdot_pos,gdot_neg, & tau_pos,tau_neg, & @@ -799,14 +802,10 @@ end subroutine plastic_kinehardening_deltaState !-------------------------------------------------------------------------------------------------- subroutine plastic_kinehardening_dotState(Tstar_v,ipc,ip,el) use lattice, only: & - lattice_Sslip_v, & - lattice_maxNslipFamily, & - lattice_NslipSystem, & - lattice_NnonSchmid + lattice_maxNslipFamily use material, only: & material_phase, & phaseAt, phasememberAt, & - plasticState, & phase_plasticityInstance implicit none @@ -819,10 +818,8 @@ subroutine plastic_kinehardening_dotState(Tstar_v,ipc,ip,el) integer(pInt) :: & instance,ph, & - f,i,j,k, & - index_myFamily,index_otherFamily, & + f,i,j, & nSlip, & - offset_accshear, & of real(pReal), dimension(plastic_kinehardening_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))) :: & @@ -873,14 +870,12 @@ end subroutine plastic_kinehardening_dotState function plastic_kinehardening_postResults(Tstar_v,ipc,ip,el) use material, only: & material_phase, & - plasticState, & phaseAt, phasememberAt, & phase_plasticityInstance use lattice, only: & lattice_Sslip_v, & lattice_maxNslipFamily, & - lattice_NslipSystem, & - lattice_NnonSchmid + lattice_NslipSystem implicit none real(pReal), dimension(6), intent(in) :: & @@ -896,7 +891,7 @@ function plastic_kinehardening_postResults(Tstar_v,ipc,ip,el) integer(pInt) :: & instance,ph, of, & nSlip,& - o,f,i,c,j,k, & + o,f,i,c,j,& index_myFamily real(pReal), dimension(plastic_kinehardening_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))) :: & From 5bf6ede661ab275d2f7887208a0eadbe645703c4 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 24 May 2018 07:57:10 +0200 Subject: [PATCH 06/27] polishing --- CMakeLists.txt | 21 ++++++++++++--------- src/CMakeLists.txt | 2 +- 2 files changed, 13 insertions(+), 10 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 5f0ebea41..ee376ef02 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -6,7 +6,7 @@ cmake_minimum_required (VERSION 2.8.8 FATAL_ERROR) # Find PETSc from system environment set(PETSC_DIR $ENV{PETSC_DIR}) if (PETSC_DIR STREQUAL "") - message (FATAL_ERROR "PETSC_DIR is not defined") + message (FATAL_ERROR "PETSc location (PETSC_DIR) is not defined") endif () set (petsc_conf_variables "${PETSC_DIR}/lib/petsc/conf/variables") @@ -114,7 +114,7 @@ elseif (DAMASK_SOLVER STREQUAL "FEM") add_definitions (-DFEM) message ("Building FEM Solver\n") else () - message (FATAL_ERROR "DAMASK_SOLVER is not defined!") + message (FATAL_ERROR "Build target (DAMASK_SOLVER) is not defined") endif () # set linker commands (needs to be done after defining the project) @@ -197,12 +197,12 @@ if (CMAKE_Fortran_COMPILER_ID STREQUAL "Intel") endif () if (OPTIMIZATION STREQUAL "OFF") - set (OPTIMIZATION_FLAGS "-O0 -no-ip") + set (OPTIMIZATION_FLAGS "-O0 -no-ip") elseif (OPTIMIZATION STREQUAL "DEFENSIVE") - set (OPTIMIZATION_FLAGS "-O2") + set (OPTIMIZATION_FLAGS "-O2") elseif (OPTIMIZATION STREQUAL "AGGRESSIVE") - set (OPTIMIZATION_FLAGS "-ipo -O3 -no-prec-div -fp-model fast=2 -xHost") - # -fast = -ipo, -O3, -no-prec-div, -static, -fp-model fast=2, and -xHost" + set (OPTIMIZATION_FLAGS "-ipo -O3 -no-prec-div -fp-model fast=2 -xHost") + # -fast = -ipo, -O3, -no-prec-div, -static, -fp-model fast=2, and -xHost" endif () # -assume std_mod_proc_name (included in -standard-semantics) causes problems if other modules @@ -317,11 +317,11 @@ elseif(CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") endif () if (OPTIMIZATION STREQUAL "OFF") - set (OPTIMIZATION_FLAGS "-O0" ) + set (OPTIMIZATION_FLAGS "-O0" ) elseif (OPTIMIZATION STREQUAL "DEFENSIVE") - set (OPTIMIZATION_FLAGS "-O2") + set (OPTIMIZATION_FLAGS "-O2") elseif (OPTIMIZATION STREQUAL "AGGRESSIVE") - set (OPTIMIZATION_FLAGS "-O3 -ffast-math -funroll-loops -ftree-vectorize") + set (OPTIMIZATION_FLAGS "-O3 -ffast-math -funroll-loops -ftree-vectorize") endif () set (STANDARD_CHECK "-std=f2008ts -pedantic-errors" ) @@ -445,8 +445,11 @@ elseif(CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") # Additional options # -fdefault-integer-8: Use it to set precision to 8 bytes for integer, don't use it for the standard case of pInt=4 (there is no -fdefault-integer-4) +else () + message (FATAL_ERROR "Compiler type (CMAKE_Fortran_COMPILER_ID) not recognized") endif () + 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_LINKER} ${OPENMP_FLAGS} ${OPTIMIZATION_FLAGS} ${LINKER_FLAGS}") diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index f0abd9707..9a828669b 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -23,7 +23,7 @@ if (PROJECT_NAME STREQUAL "DAMASK_spectral") elseif (PROJECT_NAME STREQUAL "DAMASK_FEM") add_library(DAMASK_INTERFACE OBJECT "FEM_interface.f90") else () - message ( FATAL_ERROR "PROJECT_NAME is not defined!") + message (FATAL_ERROR "Build target (PROJECT_NAME) is not defined") endif() add_dependencies(DAMASK_INTERFACE PREC) list(APPEND OBJECTFILES $) From a183fb5c471376956b3a6bcd580c36ae4226bdc8 Mon Sep 17 00:00:00 2001 From: Test User Date: Thu, 24 May 2018 21:56:03 +0200 Subject: [PATCH 07/27] [skip ci] updated version information after successful test of v2.0.2-4-g35ba260 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index e3e95eadf..affc2508b 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v2.0.1-1227-g7f0275e +v2.0.2-4-g35ba260 From f97800658f768e5a26b9235e31a0cc58ff434fd0 Mon Sep 17 00:00:00 2001 From: Eureka Pai Date: Thu, 24 May 2018 17:56:09 -0400 Subject: [PATCH 08/27] replaced param(instance) with p => pointer --- src/plastic_isotropic.f90 | 147 +++++++++++++++++++++----------------- 1 file changed, 81 insertions(+), 66 deletions(-) diff --git a/src/plastic_isotropic.f90 b/src/plastic_isotropic.f90 index a3a1d5caf..631e9661a 100644 --- a/src/plastic_isotropic.f90 +++ b/src/plastic_isotropic.f90 @@ -53,7 +53,7 @@ module plastic_isotropic dilatation = .false. end type - type(tParameters), dimension(:), allocatable, private :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable, target, private :: param !< containers of constitutive parameters (len Ninstance) type, private :: tIsotropicState !< internal state aliases real(pReal), pointer, dimension(:) :: & ! scalars along NipcMyInstance @@ -130,6 +130,7 @@ subroutine plastic_isotropic_init(fileUnit) implicit none integer(pInt), intent(in) :: fileUnit + type tParameters, pointer :: p integer(pInt), allocatable, dimension(:) :: chunkPos integer(pInt) :: & @@ -185,13 +186,14 @@ subroutine plastic_isotropic_init(fileUnit) if (IO_getTag(line,'[',']') /= '') then ! next section phase = phase + 1_pInt ! advance section counter if (phase_plasticity(phase) == PLASTICITY_ISOTROPIC_ID) then - instance = phase_plasticityInstance(phase) ! count instances of my constitutive law - allocate(param(instance)%outputID(phase_Noutput(phase))) ! allocate space for IDs of every requested output + p => param(phase_plasticityInstance(phase)) ! shorthand pointer to parameter object of my constitutive law + allocate(p%outputID(phase_Noutput(phase))) ! allocate space for IDs of every requested output endif cycle ! skip to next line endif if (phase > 0_pInt) then; if (phase_plasticity(phase) == PLASTICITY_ISOTROPIC_ID) then ! one of my phases. Do not short-circuit here (.and. between if-statements), it's not safe in Fortran instance = phase_plasticityInstance(phase) ! which instance of my plasticity is present phase + p => param(instance) chunkPos = IO_stringPos(line) tag = IO_lc(IO_stringValue(line,chunkPos,1_pInt)) ! extract key @@ -201,58 +203,58 @@ subroutine plastic_isotropic_init(fileUnit) select case(outputtag) case ('flowstress') plastic_isotropic_Noutput(instance) = plastic_isotropic_Noutput(instance) + 1_pInt - param(instance)%outputID (plastic_isotropic_Noutput(instance)) = flowstress_ID + p%outputID (plastic_isotropic_Noutput(instance)) = flowstress_ID plastic_isotropic_output(plastic_isotropic_Noutput(instance),instance) = outputtag case ('strainrate') plastic_isotropic_Noutput(instance) = plastic_isotropic_Noutput(instance) + 1_pInt - param(instance)%outputID (plastic_isotropic_Noutput(instance)) = strainrate_ID + p%outputID (plastic_isotropic_Noutput(instance)) = strainrate_ID plastic_isotropic_output(plastic_isotropic_Noutput(instance),instance) = outputtag end select case ('/dilatation/') - param(instance)%dilatation = .true. + p%dilatation = .true. case ('tau0') - param(instance)%tau0 = IO_floatValue(line,chunkPos,2_pInt) + p%tau0 = IO_floatValue(line,chunkPos,2_pInt) case ('gdot0') - param(instance)%gdot0 = IO_floatValue(line,chunkPos,2_pInt) + p%gdot0 = IO_floatValue(line,chunkPos,2_pInt) case ('n') - param(instance)%n = IO_floatValue(line,chunkPos,2_pInt) + p%n = IO_floatValue(line,chunkPos,2_pInt) case ('h0') - param(instance)%h0 = IO_floatValue(line,chunkPos,2_pInt) + p%h0 = IO_floatValue(line,chunkPos,2_pInt) case ('h0_slope','slopelnrate') - param(instance)%h0_slopeLnRate = IO_floatValue(line,chunkPos,2_pInt) + p%h0_slopeLnRate = IO_floatValue(line,chunkPos,2_pInt) case ('tausat') - param(instance)%tausat = IO_floatValue(line,chunkPos,2_pInt) + p%tausat = IO_floatValue(line,chunkPos,2_pInt) case ('tausat_sinhfita') - param(instance)%tausat_SinhFitA = IO_floatValue(line,chunkPos,2_pInt) + p%tausat_SinhFitA = IO_floatValue(line,chunkPos,2_pInt) case ('tausat_sinhfitb') - param(instance)%tausat_SinhFitB = IO_floatValue(line,chunkPos,2_pInt) + p%tausat_SinhFitB = IO_floatValue(line,chunkPos,2_pInt) case ('tausat_sinhfitc') - param(instance)%tausat_SinhFitC = IO_floatValue(line,chunkPos,2_pInt) + p%tausat_SinhFitC = IO_floatValue(line,chunkPos,2_pInt) case ('tausat_sinhfitd') - param(instance)%tausat_SinhFitD = IO_floatValue(line,chunkPos,2_pInt) + p%tausat_SinhFitD = IO_floatValue(line,chunkPos,2_pInt) case ('a', 'w0') - param(instance)%a = IO_floatValue(line,chunkPos,2_pInt) + p%a = IO_floatValue(line,chunkPos,2_pInt) case ('taylorfactor') - param(instance)%fTaylor = IO_floatValue(line,chunkPos,2_pInt) + p%fTaylor = IO_floatValue(line,chunkPos,2_pInt) case ('atol_flowstress') - param(instance)%aTolFlowstress = IO_floatValue(line,chunkPos,2_pInt) + p%aTolFlowstress = IO_floatValue(line,chunkPos,2_pInt) case ('atol_shear') - param(instance)%aTolShear = IO_floatValue(line,chunkPos,2_pInt) + p%aTolShear = IO_floatValue(line,chunkPos,2_pInt) case default @@ -269,17 +271,18 @@ subroutine plastic_isotropic_init(fileUnit) myPhase: if (phase_plasticity(phase) == PLASTICITY_isotropic_ID) then ! isolate instances of own constitutive description NipcMyPhase = count(material_phase == phase) ! number of own material points (including point components ipc) instance = phase_plasticityInstance(phase) + p => param(instance) extmsg = '' !-------------------------------------------------------------------------------------------------- ! sanity checks - if (param(instance)%aTolShear <= 0.0_pReal) param(instance)%aTolShear = 1.0e-6_pReal ! default absolute tolerance 1e-6 - if (param(instance)%tau0 < 0.0_pReal) extmsg = trim(extmsg)//' tau0' - if (param(instance)%gdot0 <= 0.0_pReal) extmsg = trim(extmsg)//' gdot0' - if (param(instance)%n <= 0.0_pReal) extmsg = trim(extmsg)//' n' - if (param(instance)%tausat <= 0.0_pReal) extmsg = trim(extmsg)//' tausat' - if (param(instance)%a <= 0.0_pReal) extmsg = trim(extmsg)//' a' - if (param(instance)%fTaylor <= 0.0_pReal) extmsg = trim(extmsg)//' taylorfactor' - if (param(instance)%aTolFlowstress <= 0.0_pReal) extmsg = trim(extmsg)//' atol_flowstress' + if (p%aTolShear <= 0.0_pReal) p%aTolShear = 1.0e-6_pReal ! default absolute tolerance 1e-6 + if (p%tau0 < 0.0_pReal) extmsg = trim(extmsg)//' tau0' + if (p%gdot0 <= 0.0_pReal) extmsg = trim(extmsg)//' gdot0' + if (p%n <= 0.0_pReal) extmsg = trim(extmsg)//' n' + if (p%tausat <= 0.0_pReal) extmsg = trim(extmsg)//' tausat' + if (p%a <= 0.0_pReal) extmsg = trim(extmsg)//' a' + if (p%fTaylor <= 0.0_pReal) extmsg = trim(extmsg)//' taylorfactor' + if (p%aTolFlowstress <= 0.0_pReal) extmsg = trim(extmsg)//' atol_flowstress' if (extmsg /= '') then extmsg = trim(extmsg)//' ('//PLASTICITY_ISOTROPIC_label//')' ! prepare error message identifier call IO_error(211_pInt,ip=instance,ext_msg=extmsg) @@ -287,7 +290,7 @@ subroutine plastic_isotropic_init(fileUnit) !-------------------------------------------------------------------------------------------------- ! Determine size of postResults array outputsLoop: do o = 1_pInt,plastic_isotropic_Noutput(instance) - select case(param(instance)%outputID(o)) + select case(p%outputID(o)) case(flowstress_ID,strainrate_ID) mySize = 1_pInt case default @@ -349,13 +352,13 @@ subroutine plastic_isotropic_init(fileUnit) !-------------------------------------------------------------------------------------------------- ! init state - state0(instance)%flowstress = param(instance)%tau0 + state0(instance)%flowstress = p%tau0 state0(instance)%accumulatedShear = 0.0_pReal !-------------------------------------------------------------------------------------------------- ! init absolute state tolerances - stateAbsTol(instance)%flowstress = param(instance)%aTolFlowstress - stateAbsTol(instance)%accumulatedShear = param(instance)%aTolShear + stateAbsTol(instance)%flowstress = p%aTolFlowstress + stateAbsTol(instance)%accumulatedShear = p%aTolShear endif myPhase enddo initializeInstances @@ -400,6 +403,8 @@ subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip,el) ip, & !< integration point el !< element + type tParameters, pointer :: p + real(pReal), dimension(3,3) :: & Tstar_dev_33 !< deviatoric part of the 2nd Piola Kirchhoff stress tensor as 2nd order tensor real(pReal), dimension(3,3,3,3) :: & @@ -414,7 +419,8 @@ subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip,el) of = phasememberAt(ipc,ip,el) ! phasememberAt should be tackled by material and be renamed to material_phasemember instance = phase_plasticityInstance(material_phase(ipc,ip,el)) - + p => param(instance) + Tstar_dev_33 = math_deviatoric33(math_Mandel6to33(Tstar_v)) ! deviatoric part of 2nd Piola-Kirchhoff stress squarenorm_Tstar_dev = math_mul33xx33(Tstar_dev_33,Tstar_dev_33) norm_Tstar_dev = sqrt(squarenorm_Tstar_dev) @@ -423,11 +429,11 @@ subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip,el) Lp = 0.0_pReal dLp_dTstar99 = 0.0_pReal else - gamma_dot = param(instance)%gdot0 & - * ( sqrt(1.5_pReal) * norm_Tstar_dev / param(instance)%fTaylor / state(instance)%flowstress(of) ) & - **param(instance)%n + gamma_dot = p%gdot0 & + * ( sqrt(1.5_pReal) * norm_Tstar_dev / p%fTaylor / state(instance)%flowstress(of) ) & + **p%n - Lp = Tstar_dev_33/norm_Tstar_dev * gamma_dot/param(instance)%fTaylor + Lp = Tstar_dev_33/norm_Tstar_dev * gamma_dot/p%fTaylor if (iand(debug_level(debug_constitutive), debug_levelExtensive) /= 0_pInt & .and. ((el == debug_e .and. ip == debug_i .and. ipc == debug_g) & @@ -441,13 +447,13 @@ subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip,el) !-------------------------------------------------------------------------------------------------- ! Calculation of the tangent of Lp forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) & - dLp_dTstar_3333(k,l,m,n) = (param(instance)%n-1.0_pReal) * & + dLp_dTstar_3333(k,l,m,n) = (p%n-1.0_pReal) * & Tstar_dev_33(k,l)*Tstar_dev_33(m,n) / squarenorm_Tstar_dev forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt) & dLp_dTstar_3333(k,l,k,l) = dLp_dTstar_3333(k,l,k,l) + 1.0_pReal forall (k=1_pInt:3_pInt,m=1_pInt:3_pInt) & dLp_dTstar_3333(k,k,m,m) = dLp_dTstar_3333(k,k,m,m) - 1.0_pReal/3.0_pReal - dLp_dTstar99 = math_Plain3333to99(gamma_dot / param(instance)%fTaylor * & + dLp_dTstar99 = math_Plain3333to99(gamma_dot / p%fTaylor * & dLp_dTstar_3333 / norm_Tstar_dev) end if end subroutine plastic_isotropic_LpAndItsTangent @@ -479,9 +485,11 @@ subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dTstar_3333,Tstar_v,ipc,ip,e ip, & !< integration point el !< element + type tParameters, pointer :: p + real(pReal), dimension(3,3) :: & Tstar_sph_33 !< sphiatoric part of the 2nd Piola Kirchhoff stress tensor as 2nd order tensor -real(pReal) :: & + real(pReal) :: & gamma_dot, & !< strainrate norm_Tstar_sph, & !< euclidean norm of Tstar_sph squarenorm_Tstar_sph !< square of the euclidean norm of Tstar_sph @@ -491,27 +499,28 @@ real(pReal) :: & of = phasememberAt(ipc,ip,el) ! phasememberAt should be tackled by material and be renamed to material_phasemember instance = phase_plasticityInstance(material_phase(ipc,ip,el)) - + p => param(instance) + Tstar_sph_33 = math_spherical33(math_Mandel6to33(Tstar_v)) ! spherical part of 2nd Piola-Kirchhoff stress squarenorm_Tstar_sph = math_mul33xx33(Tstar_sph_33,Tstar_sph_33) norm_Tstar_sph = sqrt(squarenorm_Tstar_sph) - if (param(instance)%dilatation .and. norm_Tstar_sph > 0.0_pReal) then ! Tstar == 0 or J2 plascitiy --> both Li and dLi_dTstar are zero - gamma_dot = param(instance)%gdot0 & - * (sqrt(1.5_pReal) * norm_Tstar_sph / param(instance)%fTaylor / state(instance)%flowstress(of) ) & - **param(instance)%n + if (p%dilatation .and. norm_Tstar_sph > 0.0_pReal) then ! Tstar == 0 or J2 plascitiy --> both Li and dLi_dTstar are zero + gamma_dot = p%gdot0 & + * (sqrt(1.5_pReal) * norm_Tstar_sph / p%fTaylor / state(instance)%flowstress(of) ) & + **p%n - Li = Tstar_sph_33/norm_Tstar_sph * gamma_dot/param(instance)%fTaylor + Li = Tstar_sph_33/norm_Tstar_sph * gamma_dot/p%fTaylor !-------------------------------------------------------------------------------------------------- ! Calculation of the tangent of Li forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) & - dLi_dTstar_3333(k,l,m,n) = (param(instance)%n-1.0_pReal) * & + dLi_dTstar_3333(k,l,m,n) = (p%n-1.0_pReal) * & Tstar_sph_33(k,l)*Tstar_sph_33(m,n) / squarenorm_Tstar_sph forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt) & dLi_dTstar_3333(k,l,k,l) = dLi_dTstar_3333(k,l,k,l) + 1.0_pReal - dLi_dTstar_3333 = gamma_dot / param(instance)%fTaylor * & + dLi_dTstar_3333 = gamma_dot / p%fTaylor * & dLi_dTstar_3333 / norm_Tstar_sph else Li = 0.0_pReal @@ -541,6 +550,7 @@ subroutine plastic_isotropic_dotState(Tstar_v,ipc,ip,el) ipc, & !< component-ID of integration point ip, & !< integration point el !< element + type tParameters, pointer :: p real(pReal), dimension(6) :: & Tstar_dev_v !< deviatoric 2nd Piola Kirchhoff stress tensor in Mandel notation real(pReal) :: & @@ -554,10 +564,11 @@ subroutine plastic_isotropic_dotState(Tstar_v,ipc,ip,el) of = phasememberAt(ipc,ip,el) ! phasememberAt should be tackled by material and be renamed to material_phasemember instance = phase_plasticityInstance(material_phase(ipc,ip,el)) - + p => param(instance) + !-------------------------------------------------------------------------------------------------- ! norm of (deviatoric) 2nd Piola-Kirchhoff stress - if (param(instance)%dilatation) then + if (p%dilatation) then norm_Tstar_v = sqrt(math_mul6x6(Tstar_v,Tstar_v)) else Tstar_dev_v(1:3) = Tstar_v(1:3) - sum(Tstar_v(1:3))/3.0_pReal @@ -566,26 +577,26 @@ subroutine plastic_isotropic_dotState(Tstar_v,ipc,ip,el) end if !-------------------------------------------------------------------------------------------------- ! strain rate - gamma_dot = param(instance)%gdot0 * ( sqrt(1.5_pReal) * norm_Tstar_v & + gamma_dot = p%gdot0 * ( sqrt(1.5_pReal) * norm_Tstar_v & / &!----------------------------------------------------------------------------------- - (param(instance)%fTaylor*state(instance)%flowstress(of) ))**param(instance)%n + (p%fTaylor*state(instance)%flowstress(of) ))**p%n !-------------------------------------------------------------------------------------------------- ! hardening coefficient if (abs(gamma_dot) > 1e-12_pReal) then - if (dEq0(param(instance)%tausat_SinhFitA)) then - saturation = param(instance)%tausat + if (dEq0(p%tausat_SinhFitA)) then + saturation = p%tausat else - saturation = param(instance)%tausat & - + asinh( (gamma_dot / param(instance)%tausat_SinhFitA& - )**(1.0_pReal / param(instance)%tausat_SinhFitD)& - )**(1.0_pReal / param(instance)%tausat_SinhFitC) & - / ( param(instance)%tausat_SinhFitB & - * (gamma_dot / param(instance)%gdot0)**(1.0_pReal / param(instance)%n) & + saturation = p%tausat & + + asinh( (gamma_dot / p%tausat_SinhFitA& + )**(1.0_pReal / p%tausat_SinhFitD)& + )**(1.0_pReal / p%tausat_SinhFitC) & + / ( p%tausat_SinhFitB & + * (gamma_dot / p%gdot0)**(1.0_pReal / p%n) & ) endif - hardening = ( param(instance)%h0 + param(instance)%h0_slopeLnRate * log(gamma_dot) ) & - * abs( 1.0_pReal - state(instance)%flowstress(of)/saturation )**param(instance)%a & + hardening = ( p%h0 + p%h0_slopeLnRate * log(gamma_dot) ) & + * abs( 1.0_pReal - state(instance)%flowstress(of)/saturation )**p%a & * sign(1.0_pReal, 1.0_pReal - state(instance)%flowstress(of)/saturation) else hardening = 0.0_pReal @@ -614,6 +625,9 @@ function plastic_isotropic_postResults(Tstar_v,ipc,ip,el) ipc, & !< component-ID of integration point ip, & !< integration point el !< element + + type tParameters, pointer :: p + real(pReal), dimension(plastic_isotropic_sizePostResults(phase_plasticityInstance(material_phase(ipc,ip,el)))) :: & plastic_isotropic_postResults @@ -629,10 +643,11 @@ function plastic_isotropic_postResults(Tstar_v,ipc,ip,el) of = phasememberAt(ipc,ip,el) ! phasememberAt should be tackled by material and be renamed to material_phasemember instance = phase_plasticityInstance(material_phase(ipc,ip,el)) + p => param(instance) !-------------------------------------------------------------------------------------------------- ! norm of (deviatoric) 2nd Piola-Kirchhoff stress - if (param(instance)%dilatation) then + if (p%dilatation) then norm_Tstar_v = sqrt(math_mul6x6(Tstar_v,Tstar_v)) else Tstar_dev_v(1:3) = Tstar_v(1:3) - sum(Tstar_v(1:3))/3.0_pReal @@ -644,15 +659,15 @@ function plastic_isotropic_postResults(Tstar_v,ipc,ip,el) plastic_isotropic_postResults = 0.0_pReal outputsLoop: do o = 1_pInt,plastic_isotropic_Noutput(instance) - select case(param(instance)%outputID(o)) + select case(p%outputID(o)) case (flowstress_ID) plastic_isotropic_postResults(c+1_pInt) = state(instance)%flowstress(of) c = c + 1_pInt case (strainrate_ID) plastic_isotropic_postResults(c+1_pInt) = & - param(instance)%gdot0 * ( sqrt(1.5_pReal) * norm_Tstar_v & + p%gdot0 * ( sqrt(1.5_pReal) * norm_Tstar_v & / &!---------------------------------------------------------------------------------- - (param(instance)%fTaylor * state(instance)%flowstress(of)) ) ** param(instance)%n + (p%fTaylor * state(instance)%flowstress(of)) ) ** p%n c = c + 1_pInt end select enddo outputsLoop From 982c0fb90ae27665f7d6c403c9b2751b894bf66a Mon Sep 17 00:00:00 2001 From: Eureka Pai Date: Thu, 24 May 2018 18:24:58 -0400 Subject: [PATCH 09/27] replaced param(instance) with p => pointer, corrected error --- src/plastic_isotropic.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/plastic_isotropic.f90 b/src/plastic_isotropic.f90 index 631e9661a..fe415f479 100644 --- a/src/plastic_isotropic.f90 +++ b/src/plastic_isotropic.f90 @@ -130,7 +130,7 @@ subroutine plastic_isotropic_init(fileUnit) implicit none integer(pInt), intent(in) :: fileUnit - type tParameters, pointer :: p + type(tParameters), pointer :: p integer(pInt), allocatable, dimension(:) :: chunkPos integer(pInt) :: & From c09a7fd157af44441a899b79401a687d8a36d298 Mon Sep 17 00:00:00 2001 From: Eureka Pai Date: Thu, 24 May 2018 18:31:32 -0400 Subject: [PATCH 10/27] replaced param(instance) with p => pointer, corrected errors --- src/plastic_isotropic.f90 | 11 +++++------ 1 file changed, 5 insertions(+), 6 deletions(-) diff --git a/src/plastic_isotropic.f90 b/src/plastic_isotropic.f90 index fe415f479..da6d80e60 100644 --- a/src/plastic_isotropic.f90 +++ b/src/plastic_isotropic.f90 @@ -403,7 +403,7 @@ subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip,el) ip, & !< integration point el !< element - type tParameters, pointer :: p + type(tParameters), pointer :: p real(pReal), dimension(3,3) :: & Tstar_dev_33 !< deviatoric part of the 2nd Piola Kirchhoff stress tensor as 2nd order tensor @@ -485,7 +485,7 @@ subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dTstar_3333,Tstar_v,ipc,ip,e ip, & !< integration point el !< element - type tParameters, pointer :: p + type(tParameters), pointer :: p real(pReal), dimension(3,3) :: & Tstar_sph_33 !< sphiatoric part of the 2nd Piola Kirchhoff stress tensor as 2nd order tensor @@ -526,8 +526,7 @@ subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dTstar_3333,Tstar_v,ipc,ip,e Li = 0.0_pReal dLi_dTstar_3333 = 0.0_pReal endif - -end subroutine plastic_isotropic_LiAndItsTangent + end subroutine plastic_isotropic_LiAndItsTangent !-------------------------------------------------------------------------------------------------- @@ -550,7 +549,7 @@ subroutine plastic_isotropic_dotState(Tstar_v,ipc,ip,el) ipc, & !< component-ID of integration point ip, & !< integration point el !< element - type tParameters, pointer :: p + type(tParameters), pointer :: p real(pReal), dimension(6) :: & Tstar_dev_v !< deviatoric 2nd Piola Kirchhoff stress tensor in Mandel notation real(pReal) :: & @@ -626,7 +625,7 @@ function plastic_isotropic_postResults(Tstar_v,ipc,ip,el) ip, & !< integration point el !< element - type tParameters, pointer :: p + type(tParameters), pointer :: p real(pReal), dimension(plastic_isotropic_sizePostResults(phase_plasticityInstance(material_phase(ipc,ip,el)))) :: & plastic_isotropic_postResults From bdfd3974d5ac06d21ab57db4d54d4bed79a589f8 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 25 May 2018 09:03:07 +0200 Subject: [PATCH 11/27] PRIVATE was cleaned --- PRIVATE | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PRIVATE b/PRIVATE index 74bae0779..8a2e3ccce 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 74bae0779137ad7ecbab44624b602e32f48a899d +Subproject commit 8a2e3ccce34be2a19cdcaecb8456a5dc6c2564f0 From 8184d51a99ad69ee595889fa50ca8f543a9708c0 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 26 May 2018 09:31:36 +0200 Subject: [PATCH 12/27] following style of more complex constitutive laws offset for different states needs to be computed, so it makes sense to define global and local aliases together. no need to introduce variables for state0 and aTolstate, they are only used once --- src/plastic_isotropic.f90 | 39 ++++++++------------------------------- 1 file changed, 8 insertions(+), 31 deletions(-) diff --git a/src/plastic_isotropic.f90 b/src/plastic_isotropic.f90 index da6d80e60..f26eb1ca4 100644 --- a/src/plastic_isotropic.f90 +++ b/src/plastic_isotropic.f90 @@ -61,20 +61,10 @@ module plastic_isotropic accumulatedShear end type - type, private :: tIsotropicAbsTol !< internal alias for abs tolerance in state - real(pReal), pointer :: & ! scalars - flowstress, & - accumulatedShear - end type - type(tIsotropicState), allocatable, dimension(:), private :: & !< state aliases per instance state, & - state0, & dotState - type(tIsotropicAbsTol), allocatable, dimension(:), private :: & !< state aliases per instance - stateAbsTol - public :: & plastic_isotropic_init, & plastic_isotropic_LpAndItsTangent, & @@ -263,9 +253,7 @@ subroutine plastic_isotropic_init(fileUnit) enddo parsingFile allocate(state(maxNinstance)) ! internal state aliases - allocate(state0(maxNinstance)) allocate(dotState(maxNinstance)) - allocate(stateAbsTol(maxNinstance)) initializeInstances: do phase = 1_pInt, size(phase_plasticity) ! loop over every plasticity myPhase: if (phase_plasticity(phase) == PLASTICITY_isotropic_ID) then ! isolate instances of own constitutive description @@ -334,31 +322,20 @@ subroutine plastic_isotropic_init(fileUnit) allocate(plasticState(phase)%RKCK45dotState (6,sizeDotState,NipcMyPhase),source=0.0_pReal) !-------------------------------------------------------------------------------------------------- -! globally required state aliases - plasticState(phase)%slipRate => plasticState(phase)%dotState(2:2,1:NipcMyPhase) - plasticState(phase)%accumulatedSlip => plasticState(phase)%state (2:2,1:NipcMyPhase) +! locally defined state aliases and initialization of state0 and aTolState -!-------------------------------------------------------------------------------------------------- -! locally defined state aliases state(instance)%flowstress => plasticState(phase)%state (1,1:NipcMyPhase) - state0(instance)%flowstress => plasticState(phase)%state0 (1,1:NipcMyPhase) dotState(instance)%flowstress => plasticState(phase)%dotState (1,1:NipcMyPhase) - stateAbsTol(instance)%flowstress => plasticState(phase)%aTolState(1) + plasticState(phase)%state0(1,1:NipcMyPhase) = p%tau0 + plasticState(phase)%aTolState(1) = p%aTolFlowstress state(instance)%accumulatedShear => plasticState(phase)%state (2,1:NipcMyPhase) - state0(instance)%accumulatedShear => plasticState(phase)%state0 (2,1:NipcMyPhase) dotState(instance)%accumulatedShear => plasticState(phase)%dotState (2,1:NipcMyPhase) - stateAbsTol(instance)%accumulatedShear => plasticState(phase)%aTolState(2) - -!-------------------------------------------------------------------------------------------------- -! init state - state0(instance)%flowstress = p%tau0 - state0(instance)%accumulatedShear = 0.0_pReal - -!-------------------------------------------------------------------------------------------------- -! init absolute state tolerances - stateAbsTol(instance)%flowstress = p%aTolFlowstress - stateAbsTol(instance)%accumulatedShear = p%aTolShear + plasticState(phase)%state0 (2,1:NipcMyPhase) = 0.0_pReal + plasticState(phase)%aTolState(2) = p%aTolShear + ! global alias + plasticState(phase)%slipRate => plasticState(phase)%dotState(2:2,1:NipcMyPhase) + plasticState(phase)%accumulatedSlip => plasticState(phase)%state (2:2,1:NipcMyPhase) endif myPhase enddo initializeInstances From c7c39922f071f4a31e495b1ba92921d4ad389816 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 26 May 2018 09:44:14 +0200 Subject: [PATCH 13/27] pointer assignment was done twice pointer is re-assigned automatically, but I found it confusing. Also using automatic left hand side reallocation to simplify handling of outputID --- src/plastic_isotropic.f90 | 10 +++------- 1 file changed, 3 insertions(+), 7 deletions(-) diff --git a/src/plastic_isotropic.f90 b/src/plastic_isotropic.f90 index f26eb1ca4..b30dc441e 100644 --- a/src/plastic_isotropic.f90 +++ b/src/plastic_isotropic.f90 @@ -175,15 +175,11 @@ subroutine plastic_isotropic_init(fileUnit) endif if (IO_getTag(line,'[',']') /= '') then ! next section phase = phase + 1_pInt ! advance section counter - if (phase_plasticity(phase) == PLASTICITY_ISOTROPIC_ID) then - p => param(phase_plasticityInstance(phase)) ! shorthand pointer to parameter object of my constitutive law - allocate(p%outputID(phase_Noutput(phase))) ! allocate space for IDs of every requested output - endif cycle ! skip to next line endif if (phase > 0_pInt) then; if (phase_plasticity(phase) == PLASTICITY_ISOTROPIC_ID) then ! one of my phases. Do not short-circuit here (.and. between if-statements), it's not safe in Fortran instance = phase_plasticityInstance(phase) ! which instance of my plasticity is present phase - p => param(instance) + p => param(instance) ! shorthand pointer to parameter object of my constitutive law chunkPos = IO_stringPos(line) tag = IO_lc(IO_stringValue(line,chunkPos,1_pInt)) ! extract key @@ -193,12 +189,12 @@ subroutine plastic_isotropic_init(fileUnit) select case(outputtag) case ('flowstress') plastic_isotropic_Noutput(instance) = plastic_isotropic_Noutput(instance) + 1_pInt - p%outputID (plastic_isotropic_Noutput(instance)) = flowstress_ID plastic_isotropic_output(plastic_isotropic_Noutput(instance),instance) = outputtag + p%outputID = [p%outputID,flowstress_ID] case ('strainrate') plastic_isotropic_Noutput(instance) = plastic_isotropic_Noutput(instance) + 1_pInt - p%outputID (plastic_isotropic_Noutput(instance)) = strainrate_ID plastic_isotropic_output(plastic_isotropic_Noutput(instance),instance) = outputtag + p%outputID = [p%outputID,strainrate_ID] end select case ('/dilatation/') From 2fbe60b949ff7b7fc902a398c9255e16cfee565a Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 26 May 2018 09:54:50 +0200 Subject: [PATCH 14/27] anticipate (proper) change in 23_BasticPETSc_2_PETSc as PRIVATE repository is ahead --- PRIVATE | 2 +- src/spectral_mech_Basic.f90 | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/PRIVATE b/PRIVATE index 8a2e3ccce..120a0251d 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 8a2e3ccce34be2a19cdcaecb8456a5dc6c2564f0 +Subproject commit 120a0251d368ee68b35f48e34a249b7eb8bf54a7 diff --git a/src/spectral_mech_Basic.f90 b/src/spectral_mech_Basic.f90 index 8fd2b7d96..7b6931848 100644 --- a/src/spectral_mech_Basic.f90 +++ b/src/spectral_mech_Basic.f90 @@ -22,7 +22,7 @@ module spectral_mech_basic private character (len=*), parameter, public :: & - DAMASK_spectral_SolverBasicPETSC_label = 'basicpetsc' + DAMASK_spectral_SolverBasicPETSC_label = 'petsc' !-------------------------------------------------------------------------------------------------- ! derived types From bcaab1c0688524082cea21662a5236ccbeb0cc08 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 26 May 2018 10:19:27 +0200 Subject: [PATCH 15/27] restored working behavior --- src/plastic_kinematichardening.f90 | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/src/plastic_kinematichardening.f90 b/src/plastic_kinematichardening.f90 index 553d8b6cb..c402ccd44 100644 --- a/src/plastic_kinematichardening.f90 +++ b/src/plastic_kinematichardening.f90 @@ -218,7 +218,6 @@ subroutine plastic_kinehardening_init(fileUnit) Nchunks_SlipFamilies = count(lattice_NslipSystem(:,phase) > 0_pInt) ! maximum number of slip families according to lattice type of current phase Nchunks_SlipSlip = maxval(lattice_interactionSlipSlip(:,:,phase)) Nchunks_nonSchmid = lattice_NnonSchmid(phase) - allocate(param(instance)%outputID(phase_Noutput(phase)), source=undefined_ID) ! allocate space for IDs of every requested output allocate(param(instance)%crss0 (Nchunks_SlipFamilies), source=0.0_pReal) allocate(param(instance)%tau1 (Nchunks_SlipFamilies), source=0.0_pReal) allocate(param(instance)%tau1_b (Nchunks_SlipFamilies), source=0.0_pReal) @@ -234,14 +233,11 @@ subroutine plastic_kinehardening_init(fileUnit) cycle ! skip to next line endif if (phase > 0_pInt ) then; if (phase_plasticity(phase) == PLASTICITY_KINEHARDENING_ID) then ! one of my phases. Do not short-circuit here (.and. between if-statements), it's not safe in Fortran - instance = phase_plasticityInstance(phase) ! which instance of my plasticity is present phase chunkPos = IO_stringPos(line) tag = IO_lc(IO_stringValue(line,chunkPos,1_pInt)) ! extract key select case(tag) case ('(output)') - plastic_kinehardening_Noutput(instance) = plastic_kinehardening_Noutput(instance) + 1_pInt - plastic_kinehardening_output(plastic_kinehardening_Noutput(instance),instance) = & - IO_lc(IO_stringValue(line,chunkPos,2_pInt)) + output_ID = undefined_ID select case(IO_lc(IO_stringValue(line,chunkPos,2_pInt))) case ('resistance') output_ID = crss_ID @@ -270,11 +266,15 @@ subroutine plastic_kinehardening_init(fileUnit) case ('resolvedstress') output_ID = resolvedstress_ID - case default - plastic_kinehardening_Noutput(instance) = plastic_kinehardening_Noutput(instance) - 1_pInt ! correct for invalid - end select + if (output_ID /= undefined_ID) then + plastic_kinehardening_Noutput(instance) = plastic_kinehardening_Noutput(instance) + 1_pInt + plastic_kinehardening_output(plastic_kinehardening_Noutput(instance),instance) = & + IO_lc(IO_stringValue(line,chunkPos,2_pInt)) + param(instance)%outputID = [param(instance)%outputID, output_ID] + endif + !-------------------------------------------------------------------------------------------------- ! parameters depending on number of slip families case ('nslip') From 0172987115ea4930777d94d92ee19c32c6be88e8 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 26 May 2018 12:08:28 +0200 Subject: [PATCH 16/27] forgot to update default --- src/numerics.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/numerics.f90 b/src/numerics.f90 index 3a91a36bd..27b04cd67 100644 --- a/src/numerics.f90 +++ b/src/numerics.f90 @@ -108,7 +108,7 @@ module numerics character(len=64), private :: & fftw_plan_mode = 'FFTW_PATIENT' !< reads the planing-rigor flag, see manual on www.fftw.org, Default FFTW_PATIENT: use patient planner flag character(len=64), protected, public :: & - spectral_solver = 'basicpetsc' , & !< spectral solution method + spectral_solver = 'basic', & !< spectral solution method spectral_derivative = 'continuous' !< spectral spatial derivative method character(len=1024), protected, public :: & petsc_defaultOptions = '-mech_snes_type ngmres & From 6891a48c823f43eefd866f4c56e1d6308d889178 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 26 May 2018 12:35:02 +0200 Subject: [PATCH 17/27] wrong label --- src/spectral_mech_Basic.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spectral_mech_Basic.f90 b/src/spectral_mech_Basic.f90 index 7b6931848..1c949bb7b 100644 --- a/src/spectral_mech_Basic.f90 +++ b/src/spectral_mech_Basic.f90 @@ -22,7 +22,7 @@ module spectral_mech_basic private character (len=*), parameter, public :: & - DAMASK_spectral_SolverBasicPETSC_label = 'petsc' + DAMASK_spectral_SolverBasicPETSC_label = 'basic' !-------------------------------------------------------------------------------------------------- ! derived types From 6f1c4ecdfd2d221ad4ff492cca330af89ed3cc18 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 26 May 2018 16:16:45 +0200 Subject: [PATCH 18/27] tests are still failing --- PRIVATE | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PRIVATE b/PRIVATE index 120a0251d..cd02f6c1a 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 120a0251d368ee68b35f48e34a249b7eb8bf54a7 +Subproject commit cd02f6c1a481491eb4517651516b8311348b4777 From 60cdffacd62957bdfab9b7aa3188002d828d0209 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 27 May 2018 08:09:32 +0200 Subject: [PATCH 19/27] integer kind need to match for array constructor hence, a=[a,b] requires a and b of the same kind --- src/plastic_kinematichardening.f90 | 65 +++++++++++++++--------------- 1 file changed, 33 insertions(+), 32 deletions(-) diff --git a/src/plastic_kinematichardening.f90 b/src/plastic_kinematichardening.f90 index c402ccd44..4c5bd2ea5 100644 --- a/src/plastic_kinematichardening.f90 +++ b/src/plastic_kinematichardening.f90 @@ -12,32 +12,32 @@ module plastic_kinehardening implicit none private integer(pInt), dimension(:), allocatable, public, protected :: & - plastic_kinehardening_sizePostResults !< cumulative size of post results + plastic_kinehardening_sizePostResults !< cumulative size of post results integer(pInt), dimension(:,:), allocatable, target, public :: & - plastic_kinehardening_sizePostResult !< size of each post result output + plastic_kinehardening_sizePostResult !< size of each post result output character(len=64), dimension(:,:), allocatable, target, public :: & - plastic_kinehardening_output !< name of each post result output + plastic_kinehardening_output !< name of each post result output integer(pInt), dimension(:), allocatable, target, public :: & - plastic_kinehardening_Noutput !< number of outputs per instance + plastic_kinehardening_Noutput !< number of outputs per instance integer(pInt), dimension(:), allocatable, public, protected :: & - plastic_kinehardening_totalNslip !< no. of slip system used in simulation + plastic_kinehardening_totalNslip !< no. of slip system used in simulation integer(pInt), dimension(:,:), allocatable, private :: & - plastic_kinehardening_Nslip !< active number of slip systems per family (input parameter, per family) + plastic_kinehardening_Nslip !< active number of slip systems per family (input parameter, per family) enum, bind(c) enumerator :: undefined_ID, & - crss_ID, & !< critical resolved stress - crss_back_ID, & !< critical resolved back stress - sense_ID, & !< sense of acting shear stress (-1 or +1) - chi0_ID, & !< backstress at last switch of stress sense (positive?) - gamma0_ID, & !< accumulated shear at last switch of stress sense (at current switch?) + crss_ID, & !< critical resolved stress + crss_back_ID, & !< critical resolved back stress + sense_ID, & !< sense of acting shear stress (-1 or +1) + chi0_ID, & !< backstress at last switch of stress sense (positive?) + gamma0_ID, & !< accumulated shear at last switch of stress sense (at current switch?) accshear_ID, & sumGamma_ID, & shearrate_ID, & @@ -46,26 +46,26 @@ module plastic_kinehardening end enum - type, private :: tParameters !< container type for internal constitutive parameters + type, private :: tParameters !< container type for internal constitutive parameters integer(kind(undefined_ID)), dimension(:), allocatable, private :: & - outputID !< ID of each post result output + outputID !< ID of each post result output real(pReal) :: & - gdot0, & !< reference shear strain rate for slip (input parameter) - n_slip, & !< stress exponent for slip (input parameter) + gdot0, & !< reference shear strain rate for slip (input parameter) + n_slip, & !< stress exponent for slip (input parameter) aTolResistance, & aTolShear real(pReal), dimension(:), allocatable, private :: & - crss0, & !< initial critical shear stress for slip (input parameter, per family) - theta0, & !< initial hardening rate of forward stress for each slip - theta1, & !< asymptotic hardening rate of forward stress for each slip > - theta0_b, & !< initial hardening rate of back stress for each slip > - theta1_b, & !< asymptotic hardening rate of back stress for each slip > + crss0, & !< initial critical shear stress for slip (input parameter, per family) + theta0, & !< initial hardening rate of forward stress for each slip + theta1, & !< asymptotic hardening rate of forward stress for each slip > + theta0_b, & !< initial hardening rate of back stress for each slip > + theta1_b, & !< asymptotic hardening rate of back stress for each slip > tau1, & tau1_b, & - interaction_slipslip, & !< latent hardening matrix + interaction_slipslip, & !< latent hardening matrix nonSchmidCoeff real(pReal), dimension(:,:), allocatable, private :: & @@ -73,20 +73,20 @@ module plastic_kinehardening end type type, private :: tKinehardeningState - real(pReal), pointer, dimension(:,:) :: & !< vectors along NipcMyInstance - crss, & !< critical resolved stress - crss_back, & !< critical resolved back stress - sense, & !< sense of acting shear stress (-1 or +1) - chi0, & !< backstress at last switch of stress sense - gamma0, & !< accumulated shear at last switch of stress sense - accshear !< accumulated (absolute) shear + real(pReal), pointer, dimension(:,:) :: & !< vectors along NipcMyInstance + crss, & !< critical resolved stress + crss_back, & !< critical resolved back stress + sense, & !< sense of acting shear stress (-1 or +1) + chi0, & !< backstress at last switch of stress sense + gamma0, & !< accumulated shear at last switch of stress sense + accshear !< accumulated (absolute) shear - real(pReal), pointer, dimension(:) :: & !< scalars along NipcMyInstance - sumGamma !< accumulated shear across all systems + real(pReal), pointer, dimension(:) :: & !< scalars along NipcMyInstance + sumGamma !< accumulated shear across all systems end type type(tParameters), dimension(:), allocatable, private :: & - param !< containers of constitutive parameters (len Ninstance) + param !< containers of constitutive parameters (len Ninstance) type(tKinehardeningState), allocatable, dimension(:), private :: & dotState, & @@ -155,9 +155,10 @@ subroutine plastic_kinehardening_init(fileUnit) integer(pInt), intent(in) :: fileUnit integer(pInt), allocatable, dimension(:) :: chunkPos + integer(kind(undefined_ID)) :: & + output_ID integer(pInt) :: & o, j, k, f, & - output_ID, & phase, & instance, & maxNinstance, & From e0e97abda463e1d52ebaad8e4c6a1d17dd6c7063 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 27 May 2018 08:16:05 +0200 Subject: [PATCH 20/27] citation does not match --- src/plastic_isotropic.f90 | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/src/plastic_isotropic.f90 b/src/plastic_isotropic.f90 index b30dc441e..162e0e51c 100644 --- a/src/plastic_isotropic.f90 +++ b/src/plastic_isotropic.f90 @@ -141,8 +141,6 @@ subroutine plastic_isotropic_init(fileUnit) integer(pInt) :: NipcMyPhase write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_ISOTROPIC_label//' init -+>>>' - write(6,'(/,a)') ' Ma et al., Computational Materials Science, 109:323–329, 2015' - write(6,'(/,a)') ' https://doi.org/10.1016/j.commatsci.2015.07.041' write(6,'(a15,a)') ' Current time: ',IO_timeStamp() #include "compilation_info.f90" @@ -289,7 +287,7 @@ subroutine plastic_isotropic_init(fileUnit) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - sizeDotState = 2_pInt ! flowstress, accumulated_shear + sizeDotState = len(["flowstress ","accumulated_shear"]) sizeDeltaState = 0_pInt ! no sudden jumps in state sizeState = sizeDotState + sizeDeltaState plasticState(phase)%sizeState = sizeState From 39e2e8a305e604a8496be1f70660266b4af0162a Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 27 May 2018 14:07:34 +0200 Subject: [PATCH 21/27] len is length of string, here size is required --- src/plastic_isotropic.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/plastic_isotropic.f90 b/src/plastic_isotropic.f90 index 162e0e51c..ad62ed398 100644 --- a/src/plastic_isotropic.f90 +++ b/src/plastic_isotropic.f90 @@ -287,7 +287,7 @@ subroutine plastic_isotropic_init(fileUnit) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - sizeDotState = len(["flowstress ","accumulated_shear"]) + sizeDotState = size(["flowstress ","accumulated_shear"]) sizeDeltaState = 0_pInt ! no sudden jumps in state sizeState = sizeDotState + sizeDeltaState plasticState(phase)%sizeState = sizeState From 87a16b775e01cad2d4c80e60b0795a7a01c1a5f4 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 25 May 2018 23:22:32 +0200 Subject: [PATCH 22/27] function to change working directory --- src/C_routines.c | 13 +++++++++---- src/system_routines.f90 | 35 +++++++++++++++++++++++++++++++++-- 2 files changed, 42 insertions(+), 6 deletions(-) diff --git a/src/C_routines.c b/src/C_routines.c index 5bc09745f..e3891765a 100644 --- a/src/C_routines.c +++ b/src/C_routines.c @@ -11,9 +11,9 @@ int isdirectory_c(const char *dir){ struct stat statbuf; - if(stat(dir, &statbuf) != 0) - return 0; - return S_ISDIR(statbuf.st_mode); + if(stat(dir, &statbuf) != 0) /* error */ + return 0; /* return "NO, this is not a directory" */ + return S_ISDIR(statbuf.st_mode); /* 1 => is directory, 0 => this is NOT a directory */ } @@ -29,7 +29,7 @@ void getcurrentworkdir_c(char cwd[], int *stat ){ } -void gethostname_c(char hostname[], int *stat ){ +void gethostname_c(char hostname[], int *stat){ char hostname_tmp[1024]; if(gethostname(hostname_tmp, sizeof(hostname_tmp)) == 0){ strcpy(hostname,hostname_tmp); @@ -39,3 +39,8 @@ void gethostname_c(char hostname[], int *stat ){ *stat = 1; } } + + +int chdir_c(const char *dir){ + return chdir(dir); +} diff --git a/src/system_routines.f90 b/src/system_routines.f90 index 07e12a20b..2740011b4 100644 --- a/src/system_routines.f90 +++ b/src/system_routines.f90 @@ -10,11 +10,12 @@ module system_routines public :: & isDirectory, & getCWD, & - getHostName + getHostName, & + setCWD interface - function isDirectory_C(path) BIND(C) + function isDirectory_C(path) bind(C) use, intrinsic :: ISO_C_Binding, only: & C_INT, & C_CHAR @@ -38,6 +39,14 @@ interface integer(C_INT),intent(out) :: stat end subroutine getHostName_C + function chdir_C(path) bind(C) + use, intrinsic :: ISO_C_Binding, only: & + C_INT, & + C_CHAR + integer(C_INT) :: chdir_C + character(kind=C_CHAR), dimension(1024), intent(in) :: path ! C string is an array + end function chdir_C + end interface @@ -123,5 +132,27 @@ logical function getHostName(str) end function getHostName +!-------------------------------------------------------------------------------------------------- +!> @brief changes the current working directory +!-------------------------------------------------------------------------------------------------- +logical function setCWD(path) + use, intrinsic :: ISO_C_Binding, only: & + C_INT, & + C_CHAR, & + C_NULL_CHAR + + implicit none + character(len=*), intent(in) :: path + character(kind=C_CHAR), dimension(1024) :: strFixedLength ! C string is an array + integer :: i + + strFixedLength = repeat(C_NULL_CHAR,len(strFixedLength)) + do i=1,len(path) ! copy array components + strFixedLength(i)=path(i:i) + enddo + setCWD=merge(.True.,.False.,chdir_C(strFixedLength) /= 0_C_INT) + +end function setCWD + end module system_routines From e661a42a4c2f4115a48f744fad6e09837203449c Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 28 May 2018 07:01:50 +0200 Subject: [PATCH 23/27] citation block should be separated by empty lines --- src/DAMASK_spectral.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/DAMASK_spectral.f90 b/src/DAMASK_spectral.f90 index ba7c93144..f340b683d 100644 --- a/src/DAMASK_spectral.f90 +++ b/src/DAMASK_spectral.f90 @@ -151,7 +151,7 @@ program DAMASK_spectral ! init DAMASK (all modules) call CPFEM_initAll(el = 1_pInt, ip = 1_pInt) write(6,'(/,a)') ' <<<+- DAMASK_spectral init -+>>>' - write(6,'(/,a)') ' Roters et al., Computational Materials Science, 2018' + write(6,'(/,a,/)') ' Roters et al., Computational Materials Science, 2018' write(6,'(a15,a)') ' Current time: ',IO_timeStamp() #include "compilation_info.f90" From 4bc89ad5aaddcc8945cc6f16cc259b234e5b1a3b Mon Sep 17 00:00:00 2001 From: Test User Date: Mon, 28 May 2018 19:17:52 +0200 Subject: [PATCH 24/27] [skip ci] updated version information after successful test of v2.0.2-22-g60e30e4 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index affc2508b..88f0b9112 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v2.0.2-4-g35ba260 +v2.0.2-22-g60e30e4 From 1ad700f7a073581e45245d54eb07001f99e6e1eb Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 29 May 2018 08:36:32 +0200 Subject: [PATCH 25/27] symlink code does not exist anymore --- Makefile | 1 - 1 file changed, 1 deletion(-) diff --git a/Makefile b/Makefile index 0ba67963e..4149078d1 100644 --- a/Makefile +++ b/Makefile @@ -25,7 +25,6 @@ build/FEM: .PHONY: marc marc: - @./installation/symLink_Code.sh @./installation/mods_MarcMentat/apply_DAMASK_modifications.sh ${MAKEFLAGS} .PHONY: clean From 80bfd39846743ec8dd9b2b7d88cb585398dc93c3 Mon Sep 17 00:00:00 2001 From: Test User Date: Sat, 2 Jun 2018 10:38:43 +0200 Subject: [PATCH 26/27] [skip ci] updated version information after successful test of v2.0.2-35-g4b5401e --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index 88f0b9112..388484f2d 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v2.0.2-22-g60e30e4 +v2.0.2-35-g4b5401e From c903880d196fb891f431e0b54ddd6663a7bec29f Mon Sep 17 00:00:00 2001 From: Vitesh Shah Date: Mon, 4 Jun 2018 12:36:52 +0200 Subject: [PATCH 27/27] using on unitialized error caused segfault in some cases --- src/math.f90 | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/math.f90 b/src/math.f90 index ea3754869..4b8a00646 100644 --- a/src/math.f90 +++ b/src/math.f90 @@ -1820,6 +1820,8 @@ function math_sampleFiberOri(alpha,beta,FWHM) integer(pInt):: j,& !< index of smallest component i + allocate(a(0)) + allocate(idx(0)) fInC = [sin(alpha(1))*cos(alpha(2)), sin(alpha(1))*sin(alpha(2)), cos(alpha(1))] fInS = [sin(beta(1))*cos(beta(2)), sin(beta(1))*sin(beta(2)), cos(beta(1))]