From c74e071f5e67425edef5819a4be281940a884496 Mon Sep 17 00:00:00 2001 From: Philip Eisenlohr Date: Thu, 21 Feb 2008 12:50:37 +0000 Subject: [PATCH] updated constitutive_xxx to reflect new 9x9 dLp_dTstar - used to be 3x3x3x3 --- trunk/constitutive_dislo.f90 | 17 +++++++++-------- trunk/constitutive_pheno.f90 | 9 +++++---- trunk/prec.f90 | 7 ++++--- 3 files changed, 18 insertions(+), 15 deletions(-) diff --git a/trunk/constitutive_dislo.f90 b/trunk/constitutive_dislo.f90 index 1b595b280..6196a69b6 100644 --- a/trunk/constitutive_dislo.f90 +++ b/trunk/constitutive_dislo.f90 @@ -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 diff --git a/trunk/constitutive_pheno.f90 b/trunk/constitutive_pheno.f90 index b3f3224b6..cb8204df0 100644 --- a/trunk/constitutive_pheno.f90 +++ b/trunk/constitutive_pheno.f90 @@ -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 diff --git a/trunk/prec.f90 b/trunk/prec.f90 index 371ab13b3..b41ff7a80 100644 --- a/trunk/prec.f90 +++ b/trunk/prec.f90 @@ -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)