RGC as submodule

submodules inherit use-associated entities and implicit none/private
statements
This commit is contained in:
Martin Diehl 2019-05-18 07:23:46 +02:00
parent ed8af98d69
commit 2258bfb221
4 changed files with 83 additions and 69 deletions

View File

@ -60,11 +60,24 @@ module homogenization
module subroutine mech_isostrain_init
end subroutine mech_isostrain_init
module subroutine mech_RGC_init
end subroutine mech_RGC_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_RGC_partitionDeformation(F,avgF,instance,of)
real(pReal), dimension (:,:,:), intent(out) :: F !< partitioned deformation gradient
real(pReal), dimension (3,3), intent(in) :: avgF !< average deformation gradient at material point
integer, intent(in) :: &
instance, &
of
end subroutine mech_RGC_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
@ -73,7 +86,37 @@ module homogenization
real(pReal), dimension (:,:,:,:,:), intent(in) :: dPdF !< partitioned stiffnesses
integer, intent(in) :: instance
end subroutine mech_isostrain_averageStressAndItsTangent
module subroutine mech_RGC_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_RGC_averageStressAndItsTangent
module function mech_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el)
logical, dimension(2) :: mech_RGC_updateState
real(pReal), dimension(:,:,:), intent(in) :: &
P,& !< array of P
F,& !< array of F
F0 !< array of initial F
real(pReal), dimension(:,:,:,:,:), intent(in) :: dPdF !< array of current grain stiffness
real(pReal), dimension(3,3), intent(in) :: avgF !< average F
real(pReal), intent(in) :: dt !< time increment
integer, intent(in) :: &
ip, & !< integration point number
el !< element number
end function
module subroutine mech_RGC_results(instance,group)
integer, intent(in) :: instance
character(len=*), intent(in) :: group
end subroutine
end interface
public :: &
@ -112,7 +155,7 @@ subroutine homogenization_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(homogenization_type == HOMOGENIZATION_RGC_ID)) call mech_RGC_init
if (any(thermal_type == THERMAL_isothermal_ID)) call thermal_isothermal_init
if (any(thermal_type == THERMAL_adiabatic_ID)) call thermal_adiabatic_init
@ -610,8 +653,6 @@ end subroutine materialpoint_postResults
!> @brief partition material point def grad onto constituents
!--------------------------------------------------------------------------------------------------
subroutine partitionDeformation(ip,el)
use homogenization_mech_RGC, only: &
homogenization_RGC_partitionDeformation
integer, intent(in) :: &
ip, & !< integration point
@ -628,7 +669,7 @@ subroutine partitionDeformation(ip,el)
materialpoint_subF(1:3,1:3,ip,el))
case (HOMOGENIZATION_RGC_ID) chosenHomogenization
call homogenization_RGC_partitionDeformation(&
call mech_RGC_partitionDeformation(&
crystallite_partionedF(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), &
materialpoint_subF(1:3,1:3,ip,el),&
ip, &
@ -643,8 +684,6 @@ end subroutine partitionDeformation
!> "happy" with result
!--------------------------------------------------------------------------------------------------
function updateState(ip,el)
use homogenization_mech_RGC, only: &
homogenization_RGC_updateState
use thermal_adiabatic, only: &
thermal_adiabatic_updateState
use damage_local, only: &
@ -660,14 +699,14 @@ function updateState(ip,el)
case (HOMOGENIZATION_RGC_ID) chosenHomogenization
updateState = &
updateState .and. &
homogenization_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), &
crystallite_partionedF(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), &
crystallite_partionedF0(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el),&
materialpoint_subF(1:3,1:3,ip,el),&
materialpoint_subdt(ip,el), &
crystallite_dPdF(1:3,1:3,1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), &
ip, &
el)
mech_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), &
crystallite_partionedF(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), &
crystallite_partionedF0(1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el),&
materialpoint_subF(1:3,1:3,ip,el),&
materialpoint_subdt(ip,el), &
crystallite_dPdF(1:3,1:3,1:3,1:3,1:homogenization_Ngrains(mesh_element(3,el)),ip,el), &
ip, &
el)
end select chosenHomogenization
chosenThermal: select case (thermal_type(mesh_element(3,el)))
@ -695,8 +734,6 @@ end function updateState
!> @brief derive average stress and stiffness from constituent quantities
!--------------------------------------------------------------------------------------------------
subroutine averageStressAndItsTangent(ip,el)
use homogenization_mech_RGC, only: &
homogenization_RGC_averageStressAndItsTangent
integer, intent(in) :: &
ip, & !< integration point
@ -716,7 +753,7 @@ subroutine averageStressAndItsTangent(ip,el)
homogenization_typeInstance(mesh_element(3,el)))
case (HOMOGENIZATION_RGC_ID) chosenHomogenization
call homogenization_RGC_averageStressAndItsTangent(&
call mech_RGC_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), &

View File

@ -6,20 +6,7 @@
!> @brief Relaxed grain cluster (RGC) homogenization scheme
!> Nconstituents is defined as p x q x r (cluster)
!--------------------------------------------------------------------------------------------------
module homogenization_mech_RGC
use prec
use IO
use config
use debug
use math
use material
use numerics
use constitutive
#if defined(PETSc) || defined(DAMASK_HDF5)
use results
#endif
implicit none
submodule(homogenization) homogenization_mech_RGC
enum, bind(c)
enumerator :: &
@ -79,7 +66,7 @@ contains
!--------------------------------------------------------------------------------------------------
!> @brief allocates all necessary fields, reads information from material configuration file
!--------------------------------------------------------------------------------------------------
subroutine homogenization_RGC_init
module subroutine mech_RGC_init
integer :: &
Ninstance, &
@ -197,17 +184,17 @@ subroutine homogenization_RGC_init
enddo
end subroutine homogenization_RGC_init
end subroutine mech_RGC_init
!--------------------------------------------------------------------------------------------------
!> @brief partitions the deformation gradient onto the constituents
!--------------------------------------------------------------------------------------------------
subroutine homogenization_RGC_partitionDeformation(F,avgF,instance,of)
module subroutine mech_RGC_partitionDeformation(F,avgF,instance,of)
real(pReal), dimension (:,:,:), intent(out) :: F !< partioned F per grain
real(pReal), dimension (:,:), intent(in) :: avgF !< averaged F
real(pReal), dimension (3,3), intent(in) :: avgF !< averaged F
integer, intent(in) :: &
instance, &
of
@ -247,14 +234,14 @@ subroutine homogenization_RGC_partitionDeformation(F,avgF,instance,of)
end associate
end subroutine homogenization_RGC_partitionDeformation
end subroutine mech_RGC_partitionDeformation
!--------------------------------------------------------------------------------------------------
!> @brief update the internal state of the homogenization scheme and tell whether "done" and
! "happy" with result
!--------------------------------------------------------------------------------------------------
function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el)
module function mech_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el)
real(pReal), dimension(:,:,:), intent(in) :: &
P,& !< array of P
@ -267,8 +254,6 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el)
ip, & !< integration point number
el !< element number
logical, dimension(2) :: homogenization_RGC_updateState
integer, dimension(4) :: intFaceN,intFaceP,faceID
integer, dimension(3) :: nGDim,iGr3N,iGr3P
integer :: instance,iNum,i,j,nIntFaceTot,iGrN,iGrP,iMun,iFace,k,l,ipert,iGrain,nGrain, of
@ -285,7 +270,7 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el)
#endif
zeroTimeStep: if(dEq0(dt)) then
homogenization_RGC_updateState = .true. ! pretend everything is fine and return
mech_RGC_updateState = .true. ! pretend everything is fine and return
return
endif zeroTimeStep
@ -404,12 +389,12 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el)
endif
#endif
homogenization_RGC_updateState = .false.
mech_RGC_updateState = .false.
!--------------------------------------------------------------------------------------------------
! If convergence reached => done and happy
if (residMax < relTol_RGC*stresMax .or. residMax < absTol_RGC) then
homogenization_RGC_updateState = .true.
mech_RGC_updateState = .true.
#ifdef DEBUG
if (iand(debug_level(debug_homogenization),debug_levelExtensive) /= 0 &
.and. prm%of_debug == of) write(6,'(1x,a55,/)')'... done and happy'
@ -451,7 +436,7 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el)
!--------------------------------------------------------------------------------------------------
! if residual blows-up => done but unhappy
elseif (residMax > relMax_RGC*stresMax .or. residMax > absMax_RGC) then ! try to restart when residual blows up exceeding maximum bound
homogenization_RGC_updateState = [.true.,.false.] ! with direct cut-back
mech_RGC_updateState = [.true.,.false.] ! with direct cut-back
#ifdef DEBUG
if (iand(debug_level(debug_homogenization),debug_levelExtensive) /= 0 &
@ -648,7 +633,7 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el)
enddo; enddo
stt%relaxationVector(:,of) = relax + drelax ! Updateing the state variable for the next iteration
if (any(abs(drelax) > maxdRelax_RGC)) then ! Forcing cutback when the incremental change of relaxation vector becomes too large
homogenization_RGC_updateState = [.true.,.false.]
mech_RGC_updateState = [.true.,.false.]
!$OMP CRITICAL (write2out)
write(6,'(1x,a,1x,i3,1x,a,1x,i3,1x,a)')'RGC_updateState: ip',ip,'| el',el,'enforces cutback'
write(6,'(1x,a,1x,e15.8)')'due to large relaxation change =',maxval(abs(drelax))
@ -935,13 +920,13 @@ function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el)
end subroutine grainDeformation
end function homogenization_RGC_updateState
end function mech_RGC_updateState
!--------------------------------------------------------------------------------------------------
!> @brief derive average stress and stiffness from constituent quantities
!--------------------------------------------------------------------------------------------------
subroutine homogenization_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance)
module subroutine mech_RGC_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
@ -953,7 +938,7 @@ subroutine homogenization_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,
avgP = sum(P,3) /real(product(param(instance)%Nconstituents),pReal)
dAvgPdAvgF = sum(dPdF,5)/real(product(param(instance)%Nconstituents),pReal)
end subroutine homogenization_RGC_averageStressAndItsTangent
end subroutine mech_RGC_averageStressAndItsTangent
!--------------------------------------------------------------------------------------------------
@ -963,8 +948,9 @@ end subroutine homogenization_RGC_averageStressAndItsTangent
subroutine mech_RGC_results(instance,group)
#if defined(PETSc) || defined(DAMASK_HDF5)
integer, intent(in) :: instance
character(len=*) :: group
integer, intent(in) :: instance
character(len=*), intent(in) :: group
integer :: o
associate(stt => state(instance), dst => dependentState(instance), prm => param(instance))
@ -1135,7 +1121,7 @@ integer pure function interface4to1(iFace4D, nGDim)
else
interface4to1 = iFace4D(4) + nGDim(3)*(iFace4D(2)-1) &
+ nGDim(3)*nGDim(1)*(iFace4D(3)-1) &
+ (nGDim(1)-1)*nGDim(2)*nGDim(3) ! total number of interfaces normal //e1
+ (nGDim(1)-1)*nGDim(2)*nGDim(3) ! total # of interfaces normal || e1
endif
case(3)
@ -1144,8 +1130,8 @@ integer pure function interface4to1(iFace4D, nGDim)
else
interface4to1 = iFace4D(2) + nGDim(1)*(iFace4D(3)-1) &
+ nGDim(1)*nGDim(2)*(iFace4D(4)-1) &
+ (nGDim(1)-1)*nGDim(2)*nGDim(3) & ! total number of interfaces normal //e1
+ nGDim(1)*(nGDim(2)-1)*nGDim(3) ! total number of interfaces normal //e2
+ (nGDim(1)-1)*nGDim(2)*nGDim(3) & ! total # of interfaces normal || e1
+ nGDim(1)*(nGDim(2)-1)*nGDim(3) ! total # of interfaces normal || e2
endif
case default
@ -1169,23 +1155,23 @@ pure function interface1to4(iFace1D, nGDim)
!--------------------------------------------------------------------------------------------------
! compute the total number of interfaces, which ...
nIntFace = [(nGDim(1)-1)*nGDim(2)*nGDim(3), & ! ... normal //e1
nGDim(1)*(nGDim(2)-1)*nGDim(3), & ! ... normal //e2
nGDim(1)*nGDim(2)*(nGDim(3)-1)] ! ... normal //e3
nIntFace = [(nGDim(1)-1)*nGDim(2)*nGDim(3), & ! ... normal || e1
nGDim(1)*(nGDim(2)-1)*nGDim(3), & ! ... normal || e2
nGDim(1)*nGDim(2)*(nGDim(3)-1)] ! ... normal || e3
!--------------------------------------------------------------------------------------------------
! get the corresponding interface ID in 4D (normal and local position)
if (iFace1D > 0 .and. iFace1D <= nIntFace(1)) then ! interface with normal //e1
if (iFace1D > 0 .and. iFace1D <= nIntFace(1)) then ! interface with normal || e1
interface1to4(1) = 1
interface1to4(3) = mod((iFace1D-1),nGDim(2))+1
interface1to4(4) = mod(int(real(iFace1D-1,pReal)/real(nGDim(2),pReal)),nGDim(3))+1
interface1to4(2) = int(real(iFace1D-1,pReal)/real(nGDim(2),pReal)/real(nGDim(3),pReal))+1
elseif (iFace1D > nIntFace(1) .and. iFace1D <= (nIntFace(2) + nIntFace(1))) then ! interface with normal //e2
elseif (iFace1D > nIntFace(1) .and. iFace1D <= (nIntFace(2) + nIntFace(1))) then ! interface with normal || e2
interface1to4(1) = 2
interface1to4(4) = mod((iFace1D-nIntFace(1)-1),nGDim(3))+1
interface1to4(2) = mod(int(real(iFace1D-nIntFace(1)-1,pReal)/real(nGDim(3),pReal)),nGDim(1))+1
interface1to4(3) = int(real(iFace1D-nIntFace(1)-1,pReal)/real(nGDim(3),pReal)/real(nGDim(1),pReal))+1
elseif (iFace1D > nIntFace(2) + nIntFace(1) .and. iFace1D <= (nIntFace(3) + nIntFace(2) + nIntFace(1))) then ! interface with normal //e3
elseif (iFace1D > nIntFace(2) + nIntFace(1) .and. iFace1D <= (nIntFace(3) + nIntFace(2) + nIntFace(1))) then ! interface with normal || e3
interface1to4(1) = 3
interface1to4(2) = mod((iFace1D-nIntFace(2)-nIntFace(1)-1),nGDim(1))+1
interface1to4(3) = mod(int(real(iFace1D-nIntFace(2)-nIntFace(1)-1,pReal)/real(nGDim(1),pReal)),nGDim(2))+1
@ -1195,4 +1181,4 @@ pure function interface1to4(iFace1D, nGDim)
end function interface1to4
end module homogenization_mech_RGC
end submodule homogenization_mech_RGC

View File

@ -5,11 +5,6 @@
!> @brief Isostrain (full constraint Taylor assuption) homogenization scheme
!--------------------------------------------------------------------------------------------------
submodule(homogenization) homogenization_mech_isostrain
use config
use debug
use IO
implicit none
enum, bind(c)
enumerator :: &

View File

@ -5,10 +5,6 @@
!> @brief dummy homogenization homogenization scheme for 1 constituent per material point
!--------------------------------------------------------------------------------------------------
submodule(homogenization) homogenization_mech_none
use config
use debug
implicit none
contains
@ -23,7 +19,7 @@ module subroutine mech_none_init
NofMyHomog
write(6,'(/,a)') ' <<<+- homogenization_'//HOMOGENIZATION_NONE_label//' init -+>>>'
Ninstance = count(homogenization_type == HOMOGENIZATION_NONE_ID)
if (iand(debug_level(debug_HOMOGENIZATION),debug_levelBasic) /= 0) &
write(6,'(a16,1x,i5,/)') '# instances:',Ninstance