From f3bd920c23055f5fcbdd80d4075e547e21879236 Mon Sep 17 00:00:00 2001
From: Pratheek Shanthraj
Date: Mon, 21 Jan 2013 23:50:28 +0000
Subject: [PATCH] added non-schmid structure
---
code/constitutive_nonlocal.f90 | 58 ++++++++++++++++++++++++++--------
1 file changed, 44 insertions(+), 14 deletions(-)
diff --git a/code/constitutive_nonlocal.f90 b/code/constitutive_nonlocal.f90
index d185c3bdd..31c07292c 100644
--- a/code/constitutive_nonlocal.f90
+++ b/code/constitutive_nonlocal.f90
@@ -172,6 +172,9 @@ constitutive_nonlocal_rhoDotThermalAnnihilation
real(pReal), dimension(:,:,:,:,:,:), allocatable, private :: &
constitutive_nonlocal_compatibility ! slip system compatibility between me and my neighbors
+real(pReal), dimension(:,:), allocatable, private :: &
+constitutive_nonlocal_nonSchmidCoeff
+
logical, dimension(:), allocatable, private :: &
constitutive_nonlocal_shortRangeStressCorrection, & ! flag indicating the use of the short range stress correction by a excess density gradient term
constitutive_nonlocal_deadZoneScaling, &
@@ -234,7 +237,8 @@ use lattice, only: lattice_maxNslipFamily, &
lattice_sd, &
lattice_sn, &
lattice_st, &
- lattice_interactionSlipSlip
+ lattice_interactionSlipSlip, &
+ lattice_maxNonSchmid
!*** output variables
@@ -398,6 +402,9 @@ allocate(constitutive_nonlocal_peierlsStressPerSlipFamily(lattice_maxNslipFamily
constitutive_nonlocal_minimumDipoleHeightPerSlipFamily = -1.0_pReal
constitutive_nonlocal_peierlsStressPerSlipFamily = 0.0_pReal
+allocate(constitutive_nonlocal_nonSchmidCoeff(lattice_maxNonSchmid,maxNinstance))
+constitutive_nonlocal_nonSchmidCoeff = 0.0_pReal
+
!*** readout data from material.config file
rewind(myFile)
@@ -542,6 +549,9 @@ do
constitutive_nonlocal_fEdgeMultiplication(i) = IO_floatValue(line,positions,2_pInt)
case('shortrangestresscorrection')
constitutive_nonlocal_shortRangeStressCorrection(i) = IO_floatValue(line,positions,2_pInt) > 0.0_pReal
+ case ('nonschmid_coefficients')
+ forall (f = 1_pInt:lattice_maxNonSchmid) &
+ constitutive_nonlocal_nonSchmidCoeff(f,i) = IO_floatValue(line,positions,1_pInt+f)
case('deadzonescaling','deadzone','deadscaling')
constitutive_nonlocal_deadZoneScaling(i) = IO_floatValue(line,positions,2_pInt) > 0.0_pReal
case('probabilisticmultiplication','randomsources','randommultiplication','discretesources')
@@ -1475,7 +1485,7 @@ real(pReal), dimension(constitutive_nonlocal_totalNslip(phase_plasticityInstance
!*** local variables
integer(pInt) instance, & ! current instance of this plasticity
ns, & ! short notation for the total number of active slip systems
- s ! index of my current slip system
+ s, t ! index of my current slip system
real(pReal), dimension(constitutive_nonlocal_totalNslip(phase_plasticityInstance(material_phase(g,ip,el)))) :: &
tauThreshold, & ! threshold shear stress
tauEff ! effective shear stress
@@ -1631,7 +1641,8 @@ use prec, only: pReal, &
pInt, &
p_vec
use math, only: math_Plain3333to99, &
- math_mul6x6
+ math_mul6x6, &
+ math_Mandel6to33
use debug, only: debug_level, &
debug_constitutive, &
debug_levelBasic, &
@@ -1644,7 +1655,8 @@ use material, only: homogenization_maxNgrains, &
material_phase, &
phase_plasticityInstance
use lattice, only: lattice_Sslip, &
- lattice_Sslip_v
+ lattice_Sslip_v, &
+ NnonSchmid
use mesh, only: mesh_ipVolume
implicit none
@@ -1676,15 +1688,16 @@ integer(pInt) myInstance, & ! curren
s, & ! index of my current slip system
sLattice ! index of my current slip system according to lattice order
real(pReal), dimension(3,3,3,3) :: dLp_dTstar3333 ! derivative of Lp with respect to Tstar (3x3x3x3 matrix)
+real(pReal), dimension(3,3,2) :: nonSchmid_tensor
real(pReal), dimension(constitutive_nonlocal_totalNslip(phase_plasticityInstance(material_phase(g,ip,el))),8) :: &
rhoSgl ! single dislocation densities (including used)
real(pReal), dimension(constitutive_nonlocal_totalNslip(phase_plasticityInstance(material_phase(g,ip,el))),4) :: &
v, & ! velocity
+ tau, & ! resolved shear stress including non Schmid and backstress terms
+ dgdot_dtau, & ! derivative of the shear rate with respect to the shear stress
dv_dtau ! velocity derivative with respect to the shear stress
real(pReal), dimension(constitutive_nonlocal_totalNslip(phase_plasticityInstance(material_phase(g,ip,el)))) :: &
- tau, & ! resolved shear stress including non Schmid and backstress terms
gdotTotal, & ! shear rate
- dgdotTotal_dtau, & ! derivative of the shear rate with respect to the shear stress
tauBack, & ! back stress from dislocation gradients on same slip system
deadZoneSize
@@ -1714,15 +1727,28 @@ where (abs(rhoSgl) * mesh_ipVolume(ip,el) ** 0.667_pReal < constitutive_nonlocal
!*** get effective resolved shear stress
do s = 1_pInt,ns
- tau(s) = math_mul6x6(Tstar_v, lattice_Sslip_v(:,1,constitutive_nonlocal_slipSystemLattice(s,myInstance),myStructure)) &
+ tau(s,1:4) = math_mul6x6(Tstar_v, lattice_Sslip_v(:,1,constitutive_nonlocal_slipSystemLattice(s,myInstance),myStructure)) &
+ tauBack(s)
+!*** adding non schmid contributions to ONLY screw components if present (i.e. if NnonSchmid(myStructure) > 0)
+ nonSchmid_tensor(1:3,1:3,1) = math_Mandel6to33(lattice_Sslip_v(:,1,constitutive_nonlocal_slipSystemLattice(s,myInstance),myStructure))
+ nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,1)
+ do k = 1_pInt, NnonSchmid(myStructure)
+ tau(s,3) = tau(s,3) + constitutive_nonlocal_nonSchmidCoeff(k,myInstance)* &
+ math_mul6x6(Tstar_v, lattice_Sslip_v(:,2*k,constitutive_nonlocal_slipSystemLattice(s,myInstance),myStructure))
+ tau(s,4) = tau(s,4) + constitutive_nonlocal_nonSchmidCoeff(k,myInstance)* &
+ math_mul6x6(Tstar_v, lattice_Sslip_v(:,2*k+1,constitutive_nonlocal_slipSystemLattice(s,myInstance),myStructure))
+ nonSchmid_tensor(1:3,1:3,1) = nonSchmid_tensor(1:3,1:3,1) + constitutive_nonlocal_nonSchmidCoeff(k,myInstance)*&
+ math_Mandel6to33(lattice_Sslip_v(:,2*k,constitutive_nonlocal_slipSystemLattice(s,myInstance),myStructure))
+ nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,2) + constitutive_nonlocal_nonSchmidCoeff(k,myInstance)*&
+ math_Mandel6to33(lattice_Sslip_v(:,2*k+1,constitutive_nonlocal_slipSystemLattice(s,myInstance),myStructure))
+ enddo
enddo
!*** get dislocation velocity and its tangent and store the velocity in the state array
-if (myStructure == 1_pInt) then ! for fcc all velcities are equal
- call constitutive_nonlocal_kinetics(v(1:ns,1), tau, 1_pInt, Temperature, state, g, ip, el, dv_dtau(1:ns,1))
+if (myStructure == 1_pInt .and. NnonSchmid(myStructure) == 0_pInt) then ! for fcc all velcities are equal
+ call constitutive_nonlocal_kinetics(v(1:ns,1), tau(1:ns,1), 1_pInt, Temperature, state, g, ip, el, dv_dtau(1:ns,1))
do t = 1_pInt,4_pInt
v(1:ns,t) = v(1:ns,1)
dv_dtau(1:ns,t) = dv_dtau(1:ns,1)
@@ -1731,7 +1757,7 @@ if (myStructure == 1_pInt) then ! for fcc all velcities are equal
else ! for all other lattice structures the velcities may vary with character and sign
do t = 1_pInt,4_pInt
c = (t-1_pInt)/2_pInt+1_pInt
- call constitutive_nonlocal_kinetics(v(1:ns,t), tau, c, Temperature, state, g, ip, el, dv_dtau(1:ns,t))
+ call constitutive_nonlocal_kinetics(v(1:ns,t), tau(1:ns,t), c, Temperature, state, g, ip, el, dv_dtau(1:ns,t))
state%p((12+t)*ns+1:(13+t)*ns) = v(1:ns,t)
enddo
endif
@@ -1751,8 +1777,9 @@ if (constitutive_nonlocal_deadZoneScaling(myInstance)) then
deadZoneSize(s) = maxval(abs(rhoSgl(s,5:8)) / (rhoSgl(s,1:4) + abs(rhoSgl(s,5:8))))
endif
gdotTotal = sum(rhoSgl(1:ns,1:4) * v, 2) * constitutive_nonlocal_burgers(1:ns,myInstance) * (1.0_pReal - deadZoneSize)
-dgdotTotal_dtau = sum(rhoSgl(1:ns,1:4) * dv_dtau, 2) * constitutive_nonlocal_burgers(1:ns,myInstance) * (1.0_pReal - deadZoneSize)
-
+do t = 1_pInt,4_pInt
+ dgdot_dtau(:,t) = rhoSgl(1:ns,t) * dv_dtau(1:ns,t) * constitutive_nonlocal_burgers(1:ns,myInstance) * (1.0_pReal - deadZoneSize)
+enddo
!*** Calculation of Lp and its tangent
@@ -1760,8 +1787,11 @@ do s = 1_pInt,ns
sLattice = constitutive_nonlocal_slipSystemLattice(s,myInstance)
Lp = Lp + gdotTotal(s) * lattice_Sslip(1:3,1:3,sLattice,myStructure)
forall (i=1_pInt:3_pInt,j=1_pInt:3_pInt,k=1_pInt:3_pInt,l=1_pInt:3_pInt) &
- dLp_dTstar3333(i,j,k,l) = dLp_dTstar3333(i,j,k,l) + dgdotTotal_dtau(s) * lattice_Sslip(i,j, sLattice,myStructure) &
- * lattice_Sslip(k,l, sLattice,myStructure)
+ dLp_dTstar3333(i,j,k,l) = dLp_dTstar3333(i,j,k,l) + &
+ dgdot_dtau(s,1)*lattice_Sslip(i,j,sLattice,myStructure)*lattice_Sslip(k,l,sLattice,myStructure) +&
+ dgdot_dtau(s,2)*lattice_Sslip(i,j,sLattice,myStructure)*lattice_Sslip(k,l,sLattice,myStructure) +&
+ dgdot_dtau(s,3)*lattice_Sslip(i,j,sLattice,myStructure)*nonSchmid_tensor(k,l,1) +&
+ dgdot_dtau(s,4)*lattice_Sslip(i,j,sLattice,myStructure)*nonSchmid_tensor(k,l,2)
enddo
dLp_dTstar99 = math_Plain3333to99(dLp_dTstar3333)