updated constitutive_xxx to reflect new 9x9 dLp_dTstar - used to be 3x3x3x3

This commit is contained in:
Philip Eisenlohr 2008-02-21 12:50:37 +00:00
parent 9c1e8a7944
commit c74e071f5e
3 changed files with 18 additions and 15 deletions

View File

@ -982,7 +982,8 @@ integer(pInt) matID,startIdxTwin,i,j,k,l,m,n
real(pReal) Tp,Ftwin
real(pReal), dimension(6) :: Tstar_v
real(pReal), dimension(3,3) :: Lp,Sslip,Stwin
real(pReal), dimension(3,3,3,3) :: dLp_dTstar
real(pReal), dimension(3,3,3,3) :: dLp_dTstar3333
real(pReal), dimension(9,9) :: dLp_dTstar
real(pReal), dimension(constitutive_Nstatevars(ipc,ip,el)) :: state
!* Get the material-ID from the triplet(ipc,ip,el)
@ -1040,27 +1041,27 @@ dLp_dTstar=0.0_pReal
do i=1,material_Nslip(matID)
Sslip = crystal_Sslip(:,:,i,material_CrystalStructure(matID))
forall (k=1:3,l=1:3,m=1:3,n=1:3)
dLp_dTstar(k,l,m,n) = dLp_dTstar(k,l,m,n)+ &
(1-Ftwin)*constitutive_dgdot_dtauslip(i)*Sslip(k,l)*(Sslip(m,n)+Sslip(n,m))/2.0_pReal !force m,n symmetry
dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n)+ &
(1-Ftwin)*constitutive_dgdot_dtauslip(i)*Sslip(k,l)*(Sslip(m,n)+Sslip(n,m))/2.0_pReal !force m,n symmetry
endforall
enddo
do i=1,material_Ntwin(matID)
Stwin = crystal_Stwin(:,:,i,material_CrystalStructure(matID))
forall (k=1:3,l=1:3,m=1:3,n=1:3)
dLp_dTstar(k,l,m,n) = dLp_dTstar(k,l,m,n)+ &
state(material_Nslip(matID)+i)*crystal_TwinShear(material_CrystalStructure(matID))*&
constitutive_dfdot_dtautwin(i)*Stwin(k,l)*(Stwin(m,n)+Stwin(n,m))/2.0_pReal !force m,n symmetry
dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n)+ &
state(material_Nslip(matID)+i)*crystal_TwinShear(material_CrystalStructure(matID))*&
constitutive_dfdot_dtautwin(i)*Stwin(k,l)*(Stwin(m,n)+Stwin(n,m))/2.0_pReal !force m,n symmetry
endforall
do j=1,material_Nslip(matID)
Sslip = crystal_Sslip(:,:,j,material_CrystalStructure(matID))
forall (k=1:3,l=1:3,m=1:3,n=1:3)
dLp_dTstar(k,l,m,n) = dLp_dTstar(k,l,m,n)+ &
dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n)+ &
state(material_Nslip(matID)+i)*crystal_TwinShear(material_CrystalStructure(matID))*&
constitutive_dfdot_dtauslip(i,j)*Stwin(k,l)*(Sslip(m,n)+Sslip(n,m))/2.0_pReal !force m,n symmetry
endforall
enddo
enddo
dLp_dTstar = math_Plain3333to99(dLp_dTstar3333)
return
end subroutine

View File

@ -752,6 +752,7 @@ subroutine constitutive_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,state,Temperature,
use prec, only: pReal,pInt
use crystal, only: crystal_Sslip,crystal_Sslip_v
use math, only: math_Plain3333to99
use debug
implicit none
@ -761,7 +762,7 @@ integer(pInt) matID,i,k,l,m,n
real(pReal) Temperature
real(pReal), dimension(6) :: Tstar_v
real(pReal), dimension(3,3) :: Lp
real(pReal), dimension(3,3,3,3) :: dLp
real(pReal), dimension(3,3,3,3) :: dLp_dTstar3333
real(pReal), dimension(9,9) :: dLp_dTstar
real(pReal), dimension(constitutive_Nstatevars(ipc,ip,el)) :: state
real(pReal), dimension(material_Nslip(constitutive_matID(ipc,ip,el))) :: gdot_slip,dgdot_dtauslip,tau_slip
@ -779,17 +780,17 @@ do i=1,material_Nslip(matID)
enddo
!* Calculation of the tangent of Lp
dLp = 0.0_pReal
dLp_dTstar3333 = 0.0_pReal
dLp_dTstar = 0.0_pReal
do i=1,material_Nslip(matID)
dgdot_dtauslip(i) = material_gdot0_slip(matID)*(abs(tau_slip(i))/state(i))**&
(material_n_slip(matID)-1.0_pReal)*material_n_slip(matID)/state(i)
forall (k=1:3,l=1:3,m=1:3,n=1:3) &
dLp(k,l,m,n) = dLp(k,l,m,n) + &
dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n) + &
dgdot_dtauslip(i)*crystal_Sslip(k,l,i,material_CrystalStructure(matID))* &
crystal_Sslip(m,n,i,material_CrystalStructure(matID))
enddo
dLp_dTstar = math_Plain3333to99(dLp)
dLp_dTstar = math_Plain3333to99(dLp_dTstar3333)
return
end subroutine

View File

@ -4,20 +4,21 @@
!##############################################################
implicit none
! *** Precision of real and integer variables ***
integer, parameter :: pReal = 8
integer, parameter :: pInt = 4
! *** Strain increment considered significant ***
real(pReal), parameter :: relevantStrain = 1.0e-7_pReal
! *** Numerical parameters ***
! *** How frequently the jacobian is recalculated ***
integer(pInt), parameter :: ijaco = 1_pInt ! frequency of FEM Jacobi update
integer(pInt), parameter :: nCutback = 10_pInt ! cutbacks in time-step integration
integer(pInt), parameter :: nReg = 1_pInt ! regularization attempts for Jacobi inversion
real(pReal), parameter :: pert_Fg = 1.0e-5_pReal ! strain perturbation for FEM Jacobi
integer(pInt), parameter :: nOuter = 10_pInt ! outer loop limit
integer(pInt), parameter :: nInner = 200_pInt ! inner loop limit
integer(pInt), parameter :: nInner = 1000_pInt ! inner loop limit
real(pReal), parameter :: reltol_Outer = 1.0e-4_pReal ! relative tolerance in outer loop (state)
real(pReal), parameter :: reltol_Inner = 1.0e-6_pReal ! relative tolerance in inner loop (Lp)
real(pReal), parameter :: abstol_Inner = 1.0e-8_pReal ! absolute tolerance in inner loop (Lp)