sorting according to physics

This commit is contained in:
Martin Diehl 2021-01-27 07:10:53 +01:00 committed by Martin Diehl
parent 6c99f1a234
commit 4f4adf7d68
7 changed files with 209 additions and 198 deletions

View File

@ -145,8 +145,6 @@ module homogenization
integer, intent(in) :: &
ip, & !< integration point number
el !< element number
integer :: &
co
real(pReal) :: M
end function damage_nonlocal_getMobility

View File

@ -15,6 +15,7 @@ module phase
use discretization
use parallelization
use HDF5_utilities
implicit none
private
@ -255,8 +256,6 @@ module phase
TDot
end subroutine constitutive_thermal_getRate
module subroutine plastic_nonlocal_updateCompatibility(orientation,instance,i,e)
integer, intent(in) :: &
instance, &
@ -266,54 +265,6 @@ module phase
orientation !< crystal orientation
end subroutine plastic_nonlocal_updateCompatibility
module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,me)
real(pReal), dimension(3,3), intent(out) :: &
Li !< inleastic velocity gradient
real(pReal), dimension(3,3,3,3), intent(out) :: &
dLi_dMi !< derivative of Li with respect to Mandel stress
real(pReal), dimension(3,3), intent(in) :: &
Mi !< Mandel stress
integer, intent(in) :: &
instance, &
me
end subroutine plastic_isotropic_LiAndItsTangent
module subroutine kinematics_cleavage_opening_LiAndItsTangent(Ld, dLd_dTstar, S, co, ip, el)
integer, intent(in) :: &
co, & !< grain number
ip, & !< integration point number
el !< element number
real(pReal), intent(in), dimension(3,3) :: &
S
real(pReal), intent(out), dimension(3,3) :: &
Ld !< damage velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLd_dTstar !< derivative of Ld with respect to Tstar (4th-order tensor)
end subroutine kinematics_cleavage_opening_LiAndItsTangent
module subroutine kinematics_slipplane_opening_LiAndItsTangent(Ld, dLd_dTstar, S, co, ip, el)
integer, intent(in) :: &
co, & !< grain number
ip, & !< integration point number
el !< element number
real(pReal), intent(in), dimension(3,3) :: &
S
real(pReal), intent(out), dimension(3,3) :: &
Ld !< damage velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLd_dTstar !< derivative of Ld with respect to Tstar (4th-order tensor)
end subroutine kinematics_slipplane_opening_LiAndItsTangent
module subroutine thermalexpansion_LiAndItsTangent(Li, dLi_dTstar, ph,me)
integer, intent(in) :: ph, me
!< element number
real(pReal), intent(out), dimension(3,3) :: &
Li !< thermal velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLi_dTstar !< derivative of Li with respect to Tstar (4th-order tensor defined to be zero)
end subroutine thermalexpansion_LiAndItsTangent
module subroutine plastic_dependentState(co,ip,el)
integer, intent(in) :: &
co, & !< component-ID of integration point
@ -350,7 +301,6 @@ module phase
constitutive_forward, &
constitutive_restore, &
plastic_nonlocal_updateCompatibility, &
kinematics_active, &
converged, &
crystallite_init, &
crystallite_stress, &
@ -424,38 +374,6 @@ subroutine constitutive_init
end subroutine constitutive_init
!--------------------------------------------------------------------------------------------------
!> @brief checks if a kinematic mechanism is active or not
!--------------------------------------------------------------------------------------------------
function kinematics_active(kinematics_label,kinematics_length) result(active_kinematics)
character(len=*), intent(in) :: kinematics_label !< name of kinematic mechanism
integer, intent(in) :: kinematics_length !< max. number of kinematics in system
logical, dimension(:,:), allocatable :: active_kinematics
class(tNode), pointer :: &
phases, &
phase, &
kinematics, &
kinematics_type
integer :: p,k
phases => config_material%get('phase')
allocate(active_kinematics(kinematics_length,phases%length), source = .false. )
do p = 1, phases%length
phase => phases%get(p)
kinematics => phase%get('kinematics',defaultVal=emptyList)
do k = 1, kinematics%length
kinematics_type => kinematics%get(k)
if(kinematics_type%get_asString('type') == kinematics_label) active_kinematics(k,p) = .true.
enddo
enddo
end function kinematics_active
!--------------------------------------------------------------------------------------------------
!> @brief Allocate the components of the state structure for a given phase
!--------------------------------------------------------------------------------------------------

View File

@ -57,7 +57,17 @@ submodule(phase) mechanics
end subroutine plastic_init
module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,me)
real(pReal), dimension(3,3), intent(out) :: &
Li !< inleastic velocity gradient
real(pReal), dimension(3,3,3,3), intent(out) :: &
dLi_dMi !< derivative of Li with respect to Mandel stress
real(pReal), dimension(3,3), intent(in) :: &
Mi !< Mandel stress
integer, intent(in) :: &
instance, &
me
end subroutine plastic_isotropic_LiAndItsTangent
module function plastic_dotState(subdt,co,ip,el,ph,me) result(broken)
@ -85,6 +95,25 @@ submodule(phase) mechanics
end function plastic_deltaState
module subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, &
S, Fi, co, ip, el)
integer, intent(in) :: &
co, & !< component-ID of integration point
ip, & !< integration point
el !< element
real(pReal), intent(in), dimension(3,3) :: &
S !< 2nd Piola-Kirchhoff stress
real(pReal), intent(in), dimension(3,3) :: &
Fi !< intermediate deformation gradient
real(pReal), intent(out), dimension(3,3) :: &
Li !< intermediate velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLi_dS, & !< derivative of Li with respect to S
dLi_dFi
end subroutine constitutive_LiAndItsTangents
module subroutine plastic_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, &
S, Fi, co, ip, el)
@ -173,8 +202,7 @@ module subroutine mech_init(phases)
phase, &
mech, &
elastic, &
stiffDegradation, &
kinematics
stiffDegradation
print'(/,a)', ' <<<+- constitutive_mech init -+>>>'
@ -318,31 +346,7 @@ module subroutine mech_init(phases)
end subroutine mech_init
!--------------------------------------------------------------------------------------------------
!> @brief checks if a plastic module is active or not
!--------------------------------------------------------------------------------------------------
function plastic_active(plastic_label) result(active_plastic)
character(len=*), intent(in) :: plastic_label !< type of plasticity model
logical, dimension(:), allocatable :: active_plastic
class(tNode), pointer :: &
phases, &
phase, &
mech, &
pl
integer :: ph
phases => config_material%get('phase')
allocate(active_plastic(phases%length), source = .false. )
do ph = 1, phases%length
phase => phases%get(ph)
mech => phase%get('mechanics')
pl => mech%get('plasticity')
if(pl%get_asString('type') == plastic_label) active_plastic(ph) = .true.
enddo
end function plastic_active
!--------------------------------------------------------------------------------------------------
@ -486,7 +490,6 @@ function integrateStress(F,subFp0,subFi0,Delta_t,co,ip,el) result(broken)
ierr, & ! error indicator for LAPACK
o, &
p, &
m, &
ph, &
me, &
jacoCounterLp, &
@ -1503,86 +1506,4 @@ module subroutine constitutive_mech_setF(F,co,ip,el)
end subroutine constitutive_mech_setF
!--------------------------------------------------------------------------------------------------
!> @brief contains the constitutive equation for calculating the velocity gradient
! ToDo: MD: S is Mi?
!--------------------------------------------------------------------------------------------------
subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, &
S, Fi, co, ip, el)
integer, intent(in) :: &
co, & !< component-ID of integration point
ip, & !< integration point
el !< element
real(pReal), intent(in), dimension(3,3) :: &
S !< 2nd Piola-Kirchhoff stress
real(pReal), intent(in), dimension(3,3) :: &
Fi !< intermediate deformation gradient
real(pReal), intent(out), dimension(3,3) :: &
Li !< intermediate velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLi_dS, & !< derivative of Li with respect to S
dLi_dFi
real(pReal), dimension(3,3) :: &
my_Li, & !< intermediate velocity gradient
FiInv, &
temp_33
real(pReal), dimension(3,3,3,3) :: &
my_dLi_dS
real(pReal) :: &
detFi
integer :: &
k, i, j, &
instance, of, me, ph
Li = 0.0_pReal
dLi_dS = 0.0_pReal
dLi_dFi = 0.0_pReal
plasticityType: select case (phase_plasticity(material_phaseAt(co,el)))
case (PLASTICITY_isotropic_ID) plasticityType
of = material_phasememberAt(co,ip,el)
instance = phase_plasticityInstance(material_phaseAt(co,el))
call plastic_isotropic_LiAndItsTangent(my_Li, my_dLi_dS, S ,instance,of)
case default plasticityType
my_Li = 0.0_pReal
my_dLi_dS = 0.0_pReal
end select plasticityType
Li = Li + my_Li
dLi_dS = dLi_dS + my_dLi_dS
KinematicsLoop: do k = 1, phase_Nkinematics(material_phaseAt(co,el))
kinematicsType: select case (phase_kinematics(k,material_phaseAt(co,el)))
case (KINEMATICS_cleavage_opening_ID) kinematicsType
call kinematics_cleavage_opening_LiAndItsTangent(my_Li, my_dLi_dS, S, co, ip, el)
case (KINEMATICS_slipplane_opening_ID) kinematicsType
call kinematics_slipplane_opening_LiAndItsTangent(my_Li, my_dLi_dS, S, co, ip, el)
case (KINEMATICS_thermal_expansion_ID) kinematicsType
me = material_phaseMemberAt(co,ip,el)
ph = material_phaseAt(co,el)
call thermalexpansion_LiAndItsTangent(my_Li, my_dLi_dS, ph,me)
case default kinematicsType
my_Li = 0.0_pReal
my_dLi_dS = 0.0_pReal
end select kinematicsType
Li = Li + my_Li
dLi_dS = dLi_dS + my_dLi_dS
enddo KinematicsLoop
FiInv = math_inv33(Fi)
detFi = math_det33(Fi)
Li = matmul(matmul(Fi,Li),FiInv)*detFi !< push forward to intermediate configuration
temp_33 = matmul(FiInv,Li)
do i = 1,3; do j = 1,3
dLi_dS(1:3,1:3,i,j) = matmul(matmul(Fi,dLi_dS(1:3,1:3,i,j)),FiInv)*detFi
dLi_dFi(1:3,1:3,i,j) = dLi_dFi(1:3,1:3,i,j) + Li*FiInv(j,i)
dLi_dFi(1:3,i,1:3,j) = dLi_dFi(1:3,i,1:3,j) + math_I3*temp_33(j,i) + Li*FiInv(j,i)
enddo; enddo
end subroutine constitutive_LiAndItsTangents
end submodule mechanics

View File

@ -15,6 +15,42 @@ submodule(phase:mechanics) eigendeformation
integer, intent(in) :: kinematics_length
logical, dimension(:,:), allocatable :: myKinematics
end function kinematics_thermal_expansion_init
module subroutine kinematics_cleavage_opening_LiAndItsTangent(Ld, dLd_dTstar, S, co, ip, el)
integer, intent(in) :: &
co, & !< grain number
ip, & !< integration point number
el !< element number
real(pReal), intent(in), dimension(3,3) :: &
S
real(pReal), intent(out), dimension(3,3) :: &
Ld !< damage velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLd_dTstar !< derivative of Ld with respect to Tstar (4th-order tensor)
end subroutine kinematics_cleavage_opening_LiAndItsTangent
module subroutine kinematics_slipplane_opening_LiAndItsTangent(Ld, dLd_dTstar, S, co, ip, el)
integer, intent(in) :: &
co, & !< grain number
ip, & !< integration point number
el !< element number
real(pReal), intent(in), dimension(3,3) :: &
S
real(pReal), intent(out), dimension(3,3) :: &
Ld !< damage velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLd_dTstar !< derivative of Ld with respect to Tstar (4th-order tensor)
end subroutine kinematics_slipplane_opening_LiAndItsTangent
module subroutine thermalexpansion_LiAndItsTangent(Li, dLi_dTstar, ph,me)
integer, intent(in) :: ph, me
!< element number
real(pReal), intent(out), dimension(3,3) :: &
Li !< thermal velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLi_dTstar !< derivative of Li with respect to Tstar (4th-order tensor defined to be zero)
end subroutine thermalexpansion_LiAndItsTangent
end interface
@ -54,4 +90,118 @@ module subroutine eigendeformation_init(phases)
end subroutine eigendeformation_init
!--------------------------------------------------------------------------------------------------
!> @brief checks if a kinematic mechanism is active or not
!--------------------------------------------------------------------------------------------------
function kinematics_active(kinematics_label,kinematics_length) result(active_kinematics)
character(len=*), intent(in) :: kinematics_label !< name of kinematic mechanism
integer, intent(in) :: kinematics_length !< max. number of kinematics in system
logical, dimension(:,:), allocatable :: active_kinematics
class(tNode), pointer :: &
phases, &
phase, &
kinematics, &
kinematics_type
integer :: p,k
phases => config_material%get('phase')
allocate(active_kinematics(kinematics_length,phases%length), source = .false. )
do p = 1, phases%length
phase => phases%get(p)
kinematics => phase%get('kinematics',defaultVal=emptyList)
do k = 1, kinematics%length
kinematics_type => kinematics%get(k)
if(kinematics_type%get_asString('type') == kinematics_label) active_kinematics(k,p) = .true.
enddo
enddo
end function kinematics_active
!--------------------------------------------------------------------------------------------------
!> @brief contains the constitutive equation for calculating the velocity gradient
! ToDo: MD: S is Mi?
!--------------------------------------------------------------------------------------------------
module subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, &
S, Fi, co, ip, el)
integer, intent(in) :: &
co, & !< component-ID of integration point
ip, & !< integration point
el !< element
real(pReal), intent(in), dimension(3,3) :: &
S !< 2nd Piola-Kirchhoff stress
real(pReal), intent(in), dimension(3,3) :: &
Fi !< intermediate deformation gradient
real(pReal), intent(out), dimension(3,3) :: &
Li !< intermediate velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLi_dS, & !< derivative of Li with respect to S
dLi_dFi
real(pReal), dimension(3,3) :: &
my_Li, & !< intermediate velocity gradient
FiInv, &
temp_33
real(pReal), dimension(3,3,3,3) :: &
my_dLi_dS
real(pReal) :: &
detFi
integer :: &
k, i, j, &
instance, of, me, ph
Li = 0.0_pReal
dLi_dS = 0.0_pReal
dLi_dFi = 0.0_pReal
plasticityType: select case (phase_plasticity(material_phaseAt(co,el)))
case (PLASTICITY_isotropic_ID) plasticityType
of = material_phasememberAt(co,ip,el)
instance = phase_plasticityInstance(material_phaseAt(co,el))
call plastic_isotropic_LiAndItsTangent(my_Li, my_dLi_dS, S ,instance,of)
case default plasticityType
my_Li = 0.0_pReal
my_dLi_dS = 0.0_pReal
end select plasticityType
Li = Li + my_Li
dLi_dS = dLi_dS + my_dLi_dS
KinematicsLoop: do k = 1, phase_Nkinematics(material_phaseAt(co,el))
kinematicsType: select case (phase_kinematics(k,material_phaseAt(co,el)))
case (KINEMATICS_cleavage_opening_ID) kinematicsType
call kinematics_cleavage_opening_LiAndItsTangent(my_Li, my_dLi_dS, S, co, ip, el)
case (KINEMATICS_slipplane_opening_ID) kinematicsType
call kinematics_slipplane_opening_LiAndItsTangent(my_Li, my_dLi_dS, S, co, ip, el)
case (KINEMATICS_thermal_expansion_ID) kinematicsType
me = material_phaseMemberAt(co,ip,el)
ph = material_phaseAt(co,el)
call thermalexpansion_LiAndItsTangent(my_Li, my_dLi_dS, ph,me)
case default kinematicsType
my_Li = 0.0_pReal
my_dLi_dS = 0.0_pReal
end select kinematicsType
Li = Li + my_Li
dLi_dS = dLi_dS + my_dLi_dS
enddo KinematicsLoop
FiInv = math_inv33(Fi)
detFi = math_det33(Fi)
Li = matmul(matmul(Fi,Li),FiInv)*detFi !< push forward to intermediate configuration
temp_33 = matmul(FiInv,Li)
do i = 1,3; do j = 1,3
dLi_dS(1:3,1:3,i,j) = matmul(matmul(Fi,dLi_dS(1:3,1:3,i,j)),FiInv)*detFi
dLi_dFi(1:3,1:3,i,j) = dLi_dFi(1:3,1:3,i,j) + Li*FiInv(j,i)
dLi_dFi(1:3,i,1:3,j) = dLi_dFi(1:3,i,1:3,j) + math_I3*temp_33(j,i) + Li*FiInv(j,i)
enddo; enddo
end subroutine constitutive_LiAndItsTangents
end submodule eigendeformation

View File

@ -439,4 +439,31 @@ module function plastic_deltaState(co, ip, el, ph, of) result(broken)
end function plastic_deltaState
!--------------------------------------------------------------------------------------------------
!> @brief checks if a plastic module is active or not
!--------------------------------------------------------------------------------------------------
function plastic_active(plastic_label) result(active_plastic)
character(len=*), intent(in) :: plastic_label !< type of plasticity model
logical, dimension(:), allocatable :: active_plastic
class(tNode), pointer :: &
phases, &
phase, &
mech, &
pl
integer :: ph
phases => config_material%get('phase')
allocate(active_plastic(phases%length), source = .false. )
do ph = 1, phases%length
phase => phases%get(ph)
mech => phase%get('mechanics')
pl => mech%get('plasticity')
if(pl%get_asString('type') == plastic_label) active_plastic(ph) = .true.
enddo
end function plastic_active
end submodule plastic

View File

@ -988,7 +988,6 @@ module subroutine nonlocal_dotState(Mp, Temperature,timestep, &
el !< current element number
integer :: &
ph, &
ns, & !< short notation for the total number of active slip systems
c, & !< character of dislocation
t, & !< type of dislocation

View File

@ -187,8 +187,6 @@ module function thermal_stress(Delta_t,ph,me) result(converged_)
integer, intent(in) :: ph, me
logical :: converged_
integer :: so
converged_ = .not. integrateThermalState(Delta_t,ph,me)