DAMASK_EICMD/src/homogenization_isostrain.f90

150 lines
5.5 KiB
Fortran
Raw Normal View History

!--------------------------------------------------------------------------------------------------
!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH
!> @author Franz Roters, Max-Planck-Institut für Eisenforschung GmbH
!> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH
!> @brief Isostrain (full constraint Taylor assuption) homogenization scheme
!--------------------------------------------------------------------------------------------------
module homogenization_isostrain
use prec, only: &
pInt
2019-01-13 02:34:03 +05:30
implicit none
private
enum, bind(c)
2019-01-13 02:34:03 +05:30
enumerator :: &
parallel_ID, &
average_ID
end enum
type, private :: tParameters !< container type for internal constitutive parameters
integer(pInt) :: &
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
contains
!--------------------------------------------------------------------------------------------------
!> @brief allocates all neccessary fields, reads information from material configuration file
!--------------------------------------------------------------------------------------------------
subroutine homogenization_isostrain_init()
use debug, only: &
debug_HOMOGENIZATION, &
debug_level, &
debug_levelBasic
use IO, only: &
2018-11-04 13:19:40 +05:30
IO_error
use material, only: &
homogenization_type, &
2019-03-10 15:32:32 +05:30
material_homogenizationAt, &
homogState, &
HOMOGENIZATION_ISOSTRAIN_ID, &
HOMOGENIZATION_ISOSTRAIN_LABEL, &
homogenization_typeInstance
use config, only: &
config_homogenization
implicit none
integer(pInt) :: &
2019-01-13 02:34:03 +05:30
Ninstance, &
h, &
NofMyHomog
character(len=65536) :: &
tag = ''
2019-01-13 02:34:03 +05:30
write(6,'(/,a)') ' <<<+- homogenization_'//HOMOGENIZATION_ISOSTRAIN_label//' init -+>>>'
2019-01-13 02:34:03 +05:30
Ninstance = int(count(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID),pInt)
if (iand(debug_level(debug_HOMOGENIZATION),debug_levelBasic) /= 0_pInt) &
2018-11-03 21:20:43 +05:30
write(6,'(a16,1x,i5,/)') '# instances:',Ninstance
2014-09-26 16:04:36 +05:30
2018-11-03 21:20:43 +05:30
allocate(param(Ninstance)) ! one container of parameters per instance
do h = 1_pInt, size(homogenization_type)
if (homogenization_type(h) /= HOMOGENIZATION_ISOSTRAIN_ID) cycle
2019-01-13 02:34:03 +05:30
associate(prm => param(homogenization_typeInstance(h)),&
config => config_homogenization(h))
prm%Nconstituents = config_homogenization(h)%getInt('nconstituents')
tag = 'sum'
2019-01-13 02:34:03 +05:30
select case(trim(config%getString('mapping',defaultVal = tag)))
case ('sum')
prm%mapping = parallel_ID
case ('avg')
prm%mapping = average_ID
case default
call IO_error(211_pInt,ext_msg=trim(tag)//' ('//HOMOGENIZATION_isostrain_label//')')
end select
2019-03-10 15:32:32 +05:30
NofMyHomog = count(material_homogenizationAt == h)
homogState(h)%sizeState = 0_pInt
homogState(h)%sizePostResults = 0_pInt
2018-11-04 13:19:40 +05:30
allocate(homogState(h)%state0 (0_pInt,NofMyHomog))
allocate(homogState(h)%subState0(0_pInt,NofMyHomog))
allocate(homogState(h)%state (0_pInt,NofMyHomog))
2019-01-13 02:34:03 +05:30
end associate
2014-09-26 16:04:36 +05:30
enddo
end subroutine homogenization_isostrain_init
!--------------------------------------------------------------------------------------------------
!> @brief partitions the deformation gradient onto the constituents
!--------------------------------------------------------------------------------------------------
subroutine homogenization_isostrain_partitionDeformation(F,avgF)
use prec, only: &
pReal
implicit none
real(pReal), dimension (:,:,:), intent(out) :: F !< partitioned deformation gradient
2019-01-13 02:34:03 +05:30
real(pReal), dimension (3,3), intent(in) :: avgF !< average deformation gradient at material point
F = spread(avgF,3,size(F,3))
end subroutine homogenization_isostrain_partitionDeformation
!--------------------------------------------------------------------------------------------------
!> @brief derive average stress and stiffness from constituent quantities
!--------------------------------------------------------------------------------------------------
subroutine homogenization_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF,instance)
use prec, only: &
pReal
implicit none
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
2019-01-13 02:34:03 +05:30
real(pReal), dimension (:,:,:), intent(in) :: P !< partitioned stresses
real(pReal), dimension (:,:,:,:,:), intent(in) :: dPdF !< partitioned stiffnesses
integer(pInt), intent(in) :: instance
associate(prm => param(instance))
2019-01-13 02:34:03 +05:30
select case (prm%mapping)
case (parallel_ID)
avgP = sum(P,3)
dAvgPdAvgF = sum(dPdF,5)
case (average_ID)
avgP = sum(P,3) /real(prm%Nconstituents,pReal)
dAvgPdAvgF = sum(dPdF,5)/real(prm%Nconstituents,pReal)
end select
2019-01-13 02:34:03 +05:30
end associate
end subroutine homogenization_isostrain_averageStressAndItsTangent
end module homogenization_isostrain