From 72b30729bbba8adb028bc2f3c529252d31a37163 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 5 Apr 2019 21:42:07 +0200 Subject: [PATCH] submodule for homogenization first draft, RGC not included because of name clash with isostrain --- CMakeLists.txt | 3 -- src/commercialFEM_fileList.f90 | 4 +-- src/homogenization.f90 | 48 ++++++++++++++++++--------- src/homogenization_mech_isostrain.f90 | 34 ++++++++----------- src/homogenization_mech_none.f90 | 12 +++---- 5 files changed, 52 insertions(+), 49 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 8be07198a..0cfe47248 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -337,9 +337,6 @@ elseif(CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") set (COMPILE_FLAGS "${COMPILE_FLAGS} -fimplicit-none") # assume "implicit none" even if not present in source - set (COMPILE_FLAGS "${COMPILE_FLAGS} -fmodule-private") - # assume "private" even if not present in source - set (COMPILE_FLAGS "${COMPILE_FLAGS} -Wall") # sets the following Fortran options: # -Waliasing: warn about possible aliasing of dummy arguments. Specifically, it warns if the same actual argument is associated with a dummy argument with "INTENT(IN)" and a dummy argument with "INTENT(OUT)" in a call with an explicit interface. diff --git a/src/commercialFEM_fileList.f90 b/src/commercialFEM_fileList.f90 index 5108fe853..5131eeaa9 100644 --- a/src/commercialFEM_fileList.f90 +++ b/src/commercialFEM_fileList.f90 @@ -46,8 +46,6 @@ #include "plastic_nonlocal.f90" #include "constitutive.f90" #include "crystallite.f90" -#include "homogenization_mech_none.f90" -#include "homogenization_mech_isostrain.f90" #include "homogenization_mech_RGC.f90" #include "thermal_isothermal.f90" #include "thermal_adiabatic.f90" @@ -56,4 +54,6 @@ #include "damage_local.f90" #include "damage_nonlocal.f90" #include "homogenization.f90" +#include "homogenization_mech_none.f90" +#include "homogenization_mech_isostrain.f90" #include "CPFEM.f90" diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 919680d6e..06da6ab2e 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -38,6 +38,30 @@ module homogenization materialpoint_converged logical, dimension(:,:,:), allocatable, private :: & materialpoint_doneAndHappy + + interface + + module subroutine mech_none_init + end subroutine mech_none_init + + module subroutine mech_isostrain_init + end subroutine mech_isostrain_init + + module subroutine mech_isostrain_partitionDeformation(F,avgF) + real(pReal), dimension (:,:,:), intent(out) :: F !< partitioned deformation gradient + real(pReal), dimension (3,3), intent(in) :: avgF !< average deformation gradient at material point + end subroutine mech_isostrain_partitionDeformation + + module subroutine mech_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) + real(pReal), dimension (3,3), intent(out) :: avgP !< average stress at material point + real(pReal), dimension (3,3,3,3), intent(out) :: dAvgPdAvgF !< average stiffness at material point + + real(pReal), dimension (:,:,:), intent(in) :: P !< partitioned stresses + real(pReal), dimension (:,:,:,:,:), intent(in) :: dPdF !< partitioned stiffnesses + integer, intent(in) :: instance + end subroutine mech_isostrain_averageStressAndItsTangent + + end interface public :: & homogenization_init, & @@ -77,9 +101,7 @@ subroutine homogenization_init config_homogenization, & homogenization_name use material - use homogenization_none - use homogenization_isostrain - use homogenization_RGC + use homogenization_mech_RGC use thermal_isothermal use thermal_adiabatic use thermal_conduction @@ -100,8 +122,8 @@ subroutine homogenization_init logical :: valid - if (any(homogenization_type == HOMOGENIZATION_NONE_ID)) call homogenization_none_init - if (any(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID)) call homogenization_isostrain_init + if (any(homogenization_type == HOMOGENIZATION_NONE_ID)) call mech_none_init + if (any(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID)) call mech_isostrain_init if (any(homogenization_type == HOMOGENIZATION_RGC_ID)) call homogenization_RGC_init if (any(thermal_type == THERMAL_isothermal_ID)) call thermal_isothermal_init @@ -709,9 +731,7 @@ subroutine partitionDeformation(ip,el) HOMOGENIZATION_RGC_ID use crystallite, only: & crystallite_partionedF - use homogenization_isostrain, only: & - homogenization_isostrain_partitionDeformation - use homogenization_RGC, only: & + use homogenization_mech_RGC, only: & homogenization_RGC_partitionDeformation implicit none @@ -725,7 +745,7 @@ subroutine partitionDeformation(ip,el) crystallite_partionedF(1:3,1:3,1,ip,el) = materialpoint_subF(1:3,1:3,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization - call homogenization_isostrain_partitionDeformation(& + call mech_isostrain_partitionDeformation(& crystallite_partionedF(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & materialpoint_subF(1:3,1:3,ip,el)) @@ -760,7 +780,7 @@ function updateState(ip,el) crystallite_dPdF, & crystallite_partionedF,& crystallite_partionedF0 - use homogenization_RGC, only: & + use homogenization_mech_RGC, only: & homogenization_RGC_updateState use thermal_adiabatic, only: & thermal_adiabatic_updateState @@ -824,9 +844,7 @@ subroutine averageStressAndItsTangent(ip,el) HOMOGENIZATION_RGC_ID use crystallite, only: & crystallite_P,crystallite_dPdF - use homogenization_isostrain, only: & - homogenization_isostrain_averageStressAndItsTangent - use homogenization_RGC, only: & + use homogenization_mech_RGC, only: & homogenization_RGC_averageStressAndItsTangent implicit none @@ -840,7 +858,7 @@ subroutine averageStressAndItsTangent(ip,el) materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el) = crystallite_dPdF(1:3,1:3,1:3,1:3,1,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization - call homogenization_isostrain_averageStressAndItsTangent(& + call mech_isostrain_averageStressAndItsTangent(& materialpoint_P(1:3,1:3,ip,el), & materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el),& crystallite_P(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), & @@ -887,7 +905,7 @@ function postResults(ip,el) DAMAGE_none_ID, & DAMAGE_local_ID, & DAMAGE_nonlocal_ID - use homogenization_RGC, only: & + use homogenization_mech_RGC, only: & homogenization_RGC_postResults use thermal_adiabatic, only: & thermal_adiabatic_postResults diff --git a/src/homogenization_mech_isostrain.f90 b/src/homogenization_mech_isostrain.f90 index 071420f6e..1fdf5435c 100644 --- a/src/homogenization_mech_isostrain.f90 +++ b/src/homogenization_mech_isostrain.f90 @@ -4,36 +4,32 @@ !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH !> @brief Isostrain (full constraint Taylor assuption) homogenization scheme !-------------------------------------------------------------------------------------------------- -module homogenization_mech_isostrain +submodule(homogenization) homogenization_mech_isostrain implicit none - private + enum, bind(c) enumerator :: & parallel_ID, & average_ID end enum - type, private :: tParameters !< container type for internal constitutive parameters + type :: tParameters !< container type for internal constitutive parameters integer :: & Nconstituents integer(kind(average_ID)) :: & mapping end type - type(tParameters), dimension(:), allocatable, private :: param !< containers of constitutive parameters (len Ninstance) - - public :: & - homogenization_isostrain_init, & - homogenization_isostrain_partitionDeformation, & - homogenization_isostrain_averageStressAndItsTangent + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all neccessary fields, reads information from material configuration file !-------------------------------------------------------------------------------------------------- -subroutine homogenization_isostrain_init +module subroutine mech_isostrain_init use debug, only: & debug_HOMOGENIZATION, & debug_level, & @@ -94,15 +90,13 @@ subroutine homogenization_isostrain_init enddo -end subroutine homogenization_isostrain_init +end subroutine mech_isostrain_init !-------------------------------------------------------------------------------------------------- !> @brief partitions the deformation gradient onto the constituents !-------------------------------------------------------------------------------------------------- -subroutine homogenization_isostrain_partitionDeformation(F,avgF) - use prec, only: & - pReal +module subroutine mech_isostrain_partitionDeformation(F,avgF) implicit none real(pReal), dimension (:,:,:), intent(out) :: F !< partitioned deformation gradient @@ -111,15 +105,13 @@ subroutine homogenization_isostrain_partitionDeformation(F,avgF) F = spread(avgF,3,size(F,3)) -end subroutine homogenization_isostrain_partitionDeformation +end subroutine mech_isostrain_partitionDeformation !-------------------------------------------------------------------------------------------------- !> @brief derive average stress and stiffness from constituent quantities !-------------------------------------------------------------------------------------------------- -subroutine homogenization_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) - use prec, only: & - pReal +module subroutine mech_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance) implicit none real(pReal), dimension (3,3), intent(out) :: avgP !< average stress at material point @@ -128,7 +120,7 @@ subroutine homogenization_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P real(pReal), dimension (:,:,:), intent(in) :: P !< partitioned stresses real(pReal), dimension (:,:,:,:,:), intent(in) :: dPdF !< partitioned stiffnesses integer, intent(in) :: instance - + associate(prm => param(instance)) select case (prm%mapping) @@ -142,6 +134,6 @@ subroutine homogenization_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P end associate -end subroutine homogenization_isostrain_averageStressAndItsTangent +end subroutine mech_isostrain_averageStressAndItsTangent -end module homogenization_mech_isostrain +end submodule homogenization_mech_isostrain diff --git a/src/homogenization_mech_none.f90 b/src/homogenization_mech_none.f90 index a528203b1..4ac509363 100644 --- a/src/homogenization_mech_none.f90 +++ b/src/homogenization_mech_none.f90 @@ -4,20 +4,16 @@ !> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH !> @brief dummy homogenization homogenization scheme for 1 constituent per material point !-------------------------------------------------------------------------------------------------- -module homogenization_mech_none +submodule(homogenization) homogenization_mech_none implicit none - private - - public :: & - homogenization_none_init contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all neccessary fields, reads information from material configuration file !-------------------------------------------------------------------------------------------------- -subroutine homogenization_none_init() +module subroutine mech_none_init use debug, only: & debug_HOMOGENIZATION, & debug_level, & @@ -55,6 +51,6 @@ subroutine homogenization_none_init() enddo -end subroutine homogenization_none_init +end subroutine mech_none_init -end module homogenization_mech_none +end submodule homogenization_mech_none