223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity

nonSchmid behavior now ready to use for bcc with phenopowerlaw and nonlocal constitutive model. 

nonSchmid tensors according to Koester,Ma,Hartmaier,2012.
extended lattice_Sslip with non-Schmid tensors (the full non-symmetric tensors are required for the tangent and cannot be retrieved from the symmetrized Mandel notation)
This commit is contained in:
Christoph Kords 2013-08-05 09:23:21 +00:00
parent 45c26ee471
commit 813ed130b0
5 changed files with 144 additions and 99 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
code/constitutive_dislotwin.f90
View File

@ -1139,14 +1139,14 @@ do f = 1_pInt,lattice_maxNslipFamily ! loop over
StressRatio_pminus1*(1-StressRatio_p)**(constitutive_dislotwin_q(myInstance)-1.0_pReal) StressRatio_pminus1*(1-StressRatio_p)**(constitutive_dislotwin_q(myInstance)-1.0_pReal)
!* Plastic velocity gradient for dislocation glide !* Plastic velocity gradient for dislocation glide
Lp = Lp + (1.0_pReal - sumf)*gdot_slip(j)*lattice_Sslip(:,:,index_myFamily+i,myStructure) Lp = Lp + (1.0_pReal - sumf)*gdot_slip(j)*lattice_Sslip(:,:,1,index_myFamily+i,myStructure)
!* Calculation of the tangent of Lp !* Calculation of the tangent of Lp
forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) & forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) &
dLp_dTstar3333(k,l,m,n) = & dLp_dTstar3333(k,l,m,n) = &
dLp_dTstar3333(k,l,m,n) + dgdot_dtauslip(j)*& dLp_dTstar3333(k,l,m,n) + dgdot_dtauslip(j)*&
lattice_Sslip(k,l,index_myFamily+i,myStructure)*& lattice_Sslip(k,l,1,index_myFamily+i,myStructure)*&
lattice_Sslip(m,n,index_myFamily+i,myStructure) lattice_Sslip(m,n,1,index_myFamily+i,myStructure)
enddo enddo
enddo enddo

92
code/constitutive_nonlocal.f90
View File

@ -183,7 +183,8 @@ rhoDotFluxOutput, &
rhoDotMultiplicationOutput, & rhoDotMultiplicationOutput, &
rhoDotSingle2DipoleGlideOutput, & rhoDotSingle2DipoleGlideOutput, &
rhoDotAthermalAnnihilationOutput, & rhoDotAthermalAnnihilationOutput, &
rhoDotThermalAnnihilationOutput rhoDotThermalAnnihilationOutput, &
screwStressProjection !< combined projection of Schmid and non-Schmid stresses for resolved shear stress acting on screws
real(pReal), dimension(:,:,:,:,:,:), allocatable, private :: & real(pReal), dimension(:,:,:,:,:,:), allocatable, private :: &
compatibility !< slip system compatibility between me and my neighbors compatibility !< slip system compatibility between me and my neighbors
@ -414,7 +415,7 @@ allocate(peierlsStressPerSlipFamily(lattice_maxNslipFamily,2,maxNinstance))
minDipoleHeightPerSlipFamily = -1.0_pReal minDipoleHeightPerSlipFamily = -1.0_pReal
peierlsStressPerSlipFamily = 0.0_pReal peierlsStressPerSlipFamily = 0.0_pReal
allocate(nonSchmidCoeff(lattice_maxNonSchmid,maxNinstance)) allocate(nonSchmidCoeff(lattice_maxNnonSchmid,maxNinstance))
nonSchmidCoeff = 0.0_pReal nonSchmidCoeff = 0.0_pReal
!*** readout data from material.config file !*** readout data from material.config file
@ -577,7 +578,7 @@ do while (trim(line) /= '#EOF#')
case('shortrangestresscorrection') case('shortrangestresscorrection')
shortRangeStressCorrection(i) = IO_floatValue(line,positions,2_pInt) > 0.0_pReal shortRangeStressCorrection(i) = IO_floatValue(line,positions,2_pInt) > 0.0_pReal
case ('nonschmid_coefficients') case ('nonschmid_coefficients')
do f = 1_pInt, lattice_maxNonSchmid do f = 1_pInt, lattice_maxNnonSchmid
nonSchmidCoeff(f,i) = IO_floatValue(line,positions,1_pInt+f) nonSchmidCoeff(f,i) = IO_floatValue(line,positions,1_pInt+f)
enddo enddo
case('deadzonescaling','deadzone','deadscaling') case('deadzonescaling','deadzone','deadscaling')
@ -769,6 +770,10 @@ peierlsStress = 0.0_pReal
allocate(colinearSystem(maxTotalNslip,maxNinstance)) allocate(colinearSystem(maxTotalNslip,maxNinstance))
colinearSystem = 0_pInt colinearSystem = 0_pInt
allocate(screwStressProjection(3,3,4,maxTotalNslip,maxNinstance))
screwStressProjection = 0.0_pReal
do i = 1,maxNinstance do i = 1,maxNinstance
myStructure = constitutive_nonlocal_structure(i) ! lattice structure of this instance myStructure = constitutive_nonlocal_structure(i) ! lattice structure of this instance
@ -1004,6 +1009,29 @@ do i = 1,maxNinstance
lattice_sn(1:3, slipSystemLattice(s1,i), myStructure)], [3,3])) lattice_sn(1:3, slipSystemLattice(s1,i), myStructure)], [3,3]))
enddo enddo
!*** combined projection of Schmid and non-Schmid stress contributions to resolved shear stress for screws
!* four types t:
!* 1) positive screw at positive resolved stress
!* 2) positive screw at negative resolved stress
!* 3) negative screw at positive resolved stress
!* 4) negative screw at negative resolved stress
screwStressProjection = 0.0_pReal
do s = 1_pInt,ns
do l = 1_pInt,lattice_NnonSchmid(myStructure)
screwStressProjection(1:3,1:3,1,s,i) = screwStressProjection(1:3,1:3,1,s,i) &
+ nonSchmidCoeff(l,i) * lattice_Sslip(1:3,1:3,2*l,slipSystemLattice(s,i),myStructure)
screwStressProjection(1:3,1:3,2,s,i) = screwStressProjection(1:3,1:3,2,s,i) &
+ nonSchmidCoeff(l,i) * lattice_Sslip(1:3,1:3,2*l+1,slipSystemLattice(s,i),myStructure)
enddo
screwStressProjection(1:3,1:3,3,s,i) = -screwStressProjection(1:3,1:3,2,s,i)
screwStressProjection(1:3,1:3,4,s,i) = -screwStressProjection(1:3,1:3,1,s,i)
forall (t = 1:4) &
screwStressProjection(1:3,1:3,t,s,i) = screwStressProjection(1:3,1:3,t,s,i) &
+ lattice_Sslip(1:3,1:3,1,slipSystemLattice(s,i),myStructure)
enddo
enddo enddo
endsubroutine endsubroutine
@ -1708,6 +1736,7 @@ subroutine constitutive_nonlocal_LpAndItsTangent(Lp, dLp_dTstar99, Tstar_v, Temp
use math, only: math_Plain3333to99, & use math, only: math_Plain3333to99, &
math_mul6x6, & math_mul6x6, &
math_mul33xx33, &
math_Mandel6to33 math_Mandel6to33
use debug, only: debug_level, & use debug, only: debug_level, &
debug_constitutive, & debug_constitutive, &
@ -1722,7 +1751,7 @@ use material, only: homogenization_maxNgrains, &
phase_plasticityInstance phase_plasticityInstance
use lattice, only: lattice_Sslip, & use lattice, only: lattice_Sslip, &
lattice_Sslip_v, & lattice_Sslip_v, &
NnonSchmid lattice_NnonSchmid
use mesh, only: mesh_ipVolume use mesh, only: mesh_ipVolume
implicit none implicit none
@ -1754,8 +1783,6 @@ integer(pInt) myInstance, & ! curren
s, & ! index of my current slip system s, & ! index of my current slip system
sLattice ! index of my current slip system according to lattice order 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,3,3) :: dLp_dTstar3333 ! derivative of Lp with respect to Tstar (3x3x3x3 matrix)
real(pReal), dimension(3,3,2,totalNslip(phase_plasticityInstance(material_phase(g,ip,el)))) :: &
nonSchmidTensor
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(g,ip,el))),8) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(g,ip,el))),8) :: &
rhoSgl ! single dislocation densities (including blocked) rhoSgl ! single dislocation densities (including blocked)
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(g,ip,el))),4) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(g,ip,el))),4) :: &
@ -1773,7 +1800,6 @@ real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(g,ip,e
Lp = 0.0_pReal Lp = 0.0_pReal
dLp_dTstar3333 = 0.0_pReal dLp_dTstar3333 = 0.0_pReal
nonSchmidTensor = 0.0_pReal
myInstance = phase_plasticityInstance(material_phase(g,ip,el)) myInstance = phase_plasticityInstance(material_phase(g,ip,el))
myStructure = constitutive_nonlocal_structure(myInstance) myStructure = constitutive_nonlocal_structure(myInstance)
@ -1795,31 +1821,26 @@ tauBack = state%p(iTauB(1:ns,myInstance))
!*** get effective resolved shear stress !*** get effective resolved shear stress
!*** add non schmid contributions to ONLY screw components if present (i.e. if NnonSchmid(myStructure) > 0) !*** for screws possible non-schmid contributions are also taken into account
do s = 1_pInt,ns do s = 1_pInt,ns
sLattice = slipSystemLattice(s,myInstance) sLattice = slipSystemLattice(s,myInstance)
tau(s,1:4) = math_mul6x6(Tstar_v, lattice_Sslip_v(1:6,1,sLattice,myStructure)) + tauBack(s) tau(s,1:2) = math_mul6x6(Tstar_v, lattice_Sslip_v(1:6,1,sLattice,myStructure))
nonSchmidTensor(1:3,1:3,1,s) = lattice_Sslip(1:3,1:3,sLattice,myStructure) if (tau(s,1) > 0.0_pReal) then
nonSchmidTensor(1:3,1:3,2,s) = nonSchmidTensor(1:3,1:3,1,s) tau(s,3) = math_mul33xx33(math_Mandel6to33(Tstar_v), screwStressProjection(1:3,1:3,1,s,myInstance))
do k = 1_pInt, NnonSchmid(myStructure) tau(s,4) = math_mul33xx33(math_Mandel6to33(Tstar_v), screwStressProjection(1:3,1:3,3,s,myInstance))
tau(s,3) = tau(s,3) + nonSchmidCoeff(k,myInstance) & else
* math_mul6x6(Tstar_v, lattice_Sslip_v(1:6,2*k,sLattice,myStructure)) tau(s,3) = math_mul33xx33(math_Mandel6to33(Tstar_v), screwStressProjection(1:3,1:3,2,s,myInstance))
tau(s,4) = tau(s,4) + nonSchmidCoeff(k,myInstance) & tau(s,4) = math_mul33xx33(math_Mandel6to33(Tstar_v), screwStressProjection(1:3,1:3,4,s,myInstance))
* math_mul6x6(Tstar_v, lattice_Sslip_v(1:6,2*k+1,sLattice,myStructure)) endif
nonSchmidTensor(1:3,1:3,1,s) = nonSchmidTensor(1:3,1:3,1,s) & forall (t = 1_pInt:4_pInt) &
+ nonSchmidCoeff(k,myInstance) & tau(s,t) = tau(s,t) + tauBack(s) ! add backstress
* math_Mandel6to33(lattice_Sslip_v(1:6,2*k,sLattice,myStructure))
nonSchmidTensor(1:3,1:3,2,s) = nonSchmidTensor(1:3,1:3,2,s) &
+ nonSchmidCoeff(k,myInstance) &
* math_Mandel6to33(lattice_Sslip_v(1:6,2*k+1,sLattice,myStructure))
enddo
enddo enddo
!*** get dislocation velocity and its tangent and store the velocity in the state array !*** get dislocation velocity and its tangent and store the velocity in the state array
if (myStructure == 1_pInt .and. NnonSchmid(myStructure) == 0_pInt) then ! for fcc all velcities are equal if (myStructure == 1_pInt .and. lattice_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, & call constitutive_nonlocal_kinetics(v(1:ns,1), tau(1:ns,1), 1_pInt, Temperature, state, &
g, ip, el, dv_dtau(1:ns,1)) g, ip, el, dv_dtau(1:ns,1))
do t = 1_pInt,4_pInt do t = 1_pInt,4_pInt
@ -1860,13 +1881,22 @@ enddo
do s = 1_pInt,ns do s = 1_pInt,ns
sLattice = slipSystemLattice(s,myInstance) sLattice = slipSystemLattice(s,myInstance)
Lp = Lp + gdotTotal(s) * lattice_Sslip(1:3,1:3,sLattice,myStructure) Lp = Lp + gdotTotal(s) * lattice_Sslip(1:3,1:3,1,sLattice,myStructure)
forall (i=1_pInt:3_pInt,j=1_pInt:3_pInt,k=1_pInt:3_pInt,l=1_pInt:3_pInt) & if (tau(s,1) > 0.0_pReal) then
dLp_dTstar3333(i,j,k,l) = dLp_dTstar3333(i,j,k,l) & forall (i=1_pInt:3_pInt,j=1_pInt:3_pInt,k=1_pInt:3_pInt,l=1_pInt:3_pInt) &
+ dgdot_dtau(s,1) * 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,2) * lattice_Sslip(i,j,sLattice,myStructure) * lattice_Sslip(k,l,sLattice,myStructure) & + dgdot_dtau(s,1) * lattice_Sslip(i,j,1,sLattice,myStructure) * lattice_Sslip(k,l,1,sLattice,myStructure) &
+ dgdot_dtau(s,3) * lattice_Sslip(i,j,sLattice,myStructure) * nonSchmidTensor(k,l,1,s) & + dgdot_dtau(s,2) * lattice_Sslip(i,j,1,sLattice,myStructure) * lattice_Sslip(k,l,1,sLattice,myStructure) &
+ dgdot_dtau(s,4) * lattice_Sslip(i,j,sLattice,myStructure) * nonSchmidTensor(k,l,2,s) + dgdot_dtau(s,3) * lattice_Sslip(i,j,1,sLattice,myStructure) * screwStressProjection(k,l,1,s,myInstance) &
+ dgdot_dtau(s,4) * lattice_Sslip(i,j,1,sLattice,myStructure) * screwStressProjection(k,l,3,s,myInstance)
else
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) &
+ dgdot_dtau(s,1) * lattice_Sslip(i,j,1,sLattice,myStructure) * lattice_Sslip(k,l,1,sLattice,myStructure) &
+ dgdot_dtau(s,2) * lattice_Sslip(i,j,1,sLattice,myStructure) * lattice_Sslip(k,l,1,sLattice,myStructure) &
+ dgdot_dtau(s,3) * lattice_Sslip(i,j,1,sLattice,myStructure) * screwStressProjection(k,l,2,s,myInstance) &
+ dgdot_dtau(s,4) * lattice_Sslip(i,j,1,sLattice,myStructure) * screwStressProjection(k,l,4,s,myInstance)
endif
enddo enddo
dLp_dTstar99 = math_Plain3333to99(dLp_dTstar3333) dLp_dTstar99 = math_Plain3333to99(dLp_dTstar3333)

31
code/constitutive_phenopowerlaw.f90
View File

@ -240,7 +240,7 @@ subroutine constitutive_phenopowerlaw_init(myFile)
constitutive_phenopowerlaw_aTolShear = 0.0_pReal constitutive_phenopowerlaw_aTolShear = 0.0_pReal
allocate(constitutive_phenopowerlaw_aTolTwinfrac(maxNinstance)) allocate(constitutive_phenopowerlaw_aTolTwinfrac(maxNinstance))
constitutive_phenopowerlaw_aTolTwinfrac = 0.0_pReal constitutive_phenopowerlaw_aTolTwinfrac = 0.0_pReal
allocate(constitutive_phenopowerlaw_nonSchmidCoeff(lattice_maxNonSchmid,maxNinstance)) allocate(constitutive_phenopowerlaw_nonSchmidCoeff(lattice_maxNnonSchmid,maxNinstance))
constitutive_phenopowerlaw_nonSchmidCoeff = 0.0_pReal constitutive_phenopowerlaw_nonSchmidCoeff = 0.0_pReal
rewind(myFile) rewind(myFile)
@ -368,7 +368,7 @@ subroutine constitutive_phenopowerlaw_init(myFile)
constitutive_phenopowerlaw_interaction_TwinTwin(j,i) = IO_floatValue(line,positions,1_pInt+j) constitutive_phenopowerlaw_interaction_TwinTwin(j,i) = IO_floatValue(line,positions,1_pInt+j)
enddo enddo
case ('nonschmid_coefficients') case ('nonschmid_coefficients')
do j = 1_pInt, lattice_maxNonSchmid do j = 1_pInt, lattice_maxNnonSchmid
constitutive_phenopowerlaw_nonSchmidCoeff(j,i) = IO_floatValue(line,positions,1_pInt+j) constitutive_phenopowerlaw_nonSchmidCoeff(j,i) = IO_floatValue(line,positions,1_pInt+j)
enddo enddo
case default case default
@ -682,7 +682,7 @@ pure subroutine constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar
lattice_maxNtwinFamily, & lattice_maxNtwinFamily, &
lattice_NslipSystem, & lattice_NslipSystem, &
lattice_NtwinSystem, & lattice_NtwinSystem, &
NnonSchmid lattice_NnonSchmid
use mesh, only: & use mesh, only: &
mesh_NcpElems, & mesh_NcpElems, &
mesh_maxNips mesh_maxNips
@ -745,17 +745,17 @@ pure subroutine constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar
! Calculation of Lp ! Calculation of Lp
tau_slip_pos(j) = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,structID)) tau_slip_pos(j) = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,structID))
tau_slip_neg(j) = tau_slip_pos(j) tau_slip_neg(j) = tau_slip_pos(j)
nonSchmid_tensor(1:3,1:3,1) = math_Mandel6to33(lattice_Sslip_v(1:6,1,index_myFamily+i,structID)) nonSchmid_tensor(1:3,1:3,1) = lattice_Sslip(1:3,1:3,1,index_myFamily+i,structID)
nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,1) nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,1)
do k = 1, NnonSchmid(structID) do k = 1,lattice_NnonSchmid(structID)
tau_slip_pos(j) = tau_slip_pos(j) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* & tau_slip_pos(j) = tau_slip_pos(j) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* &
dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k,index_myFamily+i,structID)) dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k,index_myFamily+i,structID))
tau_slip_neg(j) = tau_slip_neg(j) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* & tau_slip_neg(j) = tau_slip_neg(j) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* &
dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k+1,index_myFamily+i,structID)) dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k+1,index_myFamily+i,structID))
nonSchmid_tensor(1:3,1:3,1) = nonSchmid_tensor(1:3,1:3,1) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)*& nonSchmid_tensor(1:3,1:3,1) = nonSchmid_tensor(1:3,1:3,1) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)*&
math_Mandel6to33(lattice_Sslip_v(1:6,2*k,index_myFamily+i,structID)) lattice_Sslip(1:3,1:3,2*k,index_myFamily+i,structID)
nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,2) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)*& nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,2) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)*&
math_Mandel6to33(lattice_Sslip_v(1:6,2*k+1,index_myFamily+i,structID)) lattice_Sslip(1:3,1:3,2*k+1,index_myFamily+i,structID)
enddo enddo
gdot_slip_pos(j) = 0.5_pReal*constitutive_phenopowerlaw_gdot0_slip(matID)* & gdot_slip_pos(j) = 0.5_pReal*constitutive_phenopowerlaw_gdot0_slip(matID)* &
((abs(tau_slip_pos(j))/state(ipc,ip,el)%p(j))**constitutive_phenopowerlaw_n_slip(matID))*& ((abs(tau_slip_pos(j))/state(ipc,ip,el)%p(j))**constitutive_phenopowerlaw_n_slip(matID))*&
@ -764,7 +764,7 @@ pure subroutine constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar
((abs(tau_slip_neg(j))/state(ipc,ip,el)%p(j))**constitutive_phenopowerlaw_n_slip(matID))*& ((abs(tau_slip_neg(j))/state(ipc,ip,el)%p(j))**constitutive_phenopowerlaw_n_slip(matID))*&
sign(1.0_pReal,tau_slip_neg(j)) sign(1.0_pReal,tau_slip_neg(j))
Lp = Lp + (1.0_pReal-state(ipc,ip,el)%p(index_F))*& ! 1-F Lp = Lp + (1.0_pReal-state(ipc,ip,el)%p(index_F))*& ! 1-F
(gdot_slip_pos(j)+gdot_slip_neg(j))*lattice_Sslip(1:3,1:3,index_myFamily+i,structID) (gdot_slip_pos(j)+gdot_slip_neg(j))*lattice_Sslip(1:3,1:3,1,index_myFamily+i,structID)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! Calculation of the tangent of Lp ! Calculation of the tangent of Lp
@ -772,7 +772,7 @@ pure subroutine constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar
dgdot_dtauslip_pos(j) = gdot_slip_pos(j)*constitutive_phenopowerlaw_n_slip(matID)/tau_slip_pos(j) dgdot_dtauslip_pos(j) = gdot_slip_pos(j)*constitutive_phenopowerlaw_n_slip(matID)/tau_slip_pos(j)
forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) & forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) &
dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n) + & dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n) + &
dgdot_dtauslip_pos(j)*lattice_Sslip(k,l,index_myFamily+i,structID)* & dgdot_dtauslip_pos(j)*lattice_Sslip(k,l,1,index_myFamily+i,structID)* &
nonSchmid_tensor(m,n,1) nonSchmid_tensor(m,n,1)
endif endif
@ -780,7 +780,7 @@ pure subroutine constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar
dgdot_dtauslip_neg(j) = gdot_slip_neg(j)*constitutive_phenopowerlaw_n_slip(matID)/tau_slip_neg(j) dgdot_dtauslip_neg(j) = gdot_slip_neg(j)*constitutive_phenopowerlaw_n_slip(matID)/tau_slip_neg(j)
forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) & forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) &
dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n) + & dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n) + &
dgdot_dtauslip_neg(j)*lattice_Sslip(k,l,index_myFamily+i,structID)* & dgdot_dtauslip_neg(j)*lattice_Sslip(k,l,1,index_myFamily+i,structID)* &
nonSchmid_tensor(m,n,2) nonSchmid_tensor(m,n,2)
endif endif
enddo enddo
@ -832,7 +832,7 @@ function constitutive_phenopowerlaw_dotState(Tstar_v,temperature,state,ipc,ip,el
lattice_NslipSystem, & lattice_NslipSystem, &
lattice_NtwinSystem, & lattice_NtwinSystem, &
lattice_shearTwin, & lattice_shearTwin, &
NnonSchmid lattice_NnonSchmid
use mesh, only: & use mesh, only: &
mesh_NcpElems,& mesh_NcpElems,&
mesh_maxNips mesh_maxNips
@ -911,7 +911,7 @@ function constitutive_phenopowerlaw_dotState(Tstar_v,temperature,state,ipc,ip,el
! Calculation of dot gamma ! Calculation of dot gamma
tau_slip_pos(j) = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,structID)) tau_slip_pos(j) = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,structID))
tau_slip_neg(j) = tau_slip_pos(j) tau_slip_neg(j) = tau_slip_pos(j)
do k = 1, NnonSchmid(structID) do k = 1,lattice_NnonSchmid(structID)
tau_slip_pos(j) = tau_slip_pos(j) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* & tau_slip_pos(j) = tau_slip_pos(j) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* &
dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k,index_myFamily+i,structID)) dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k,index_myFamily+i,structID))
tau_slip_neg(j) = tau_slip_neg(j) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* & tau_slip_neg(j) = tau_slip_neg(j) + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* &
@ -1069,7 +1069,8 @@ pure function constitutive_phenopowerlaw_postResults(Tstar_v,temperature,dt,stat
lattice_maxNslipFamily, & lattice_maxNslipFamily, &
lattice_maxNtwinFamily, & lattice_maxNtwinFamily, &
lattice_NslipSystem, & lattice_NslipSystem, &
lattice_NtwinSystem,NnonSchmid lattice_NtwinSystem, &
lattice_NnonSchmid
use mesh, only: & use mesh, only: &
mesh_NcpElems, & mesh_NcpElems, &
mesh_maxNips mesh_maxNips
@ -1132,7 +1133,7 @@ pure function constitutive_phenopowerlaw_postResults(Tstar_v,temperature,dt,stat
j = j + 1_pInt j = j + 1_pInt
tau_slip_pos = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,structID)) tau_slip_pos = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,structID))
tau_slip_neg = tau_slip_pos tau_slip_neg = tau_slip_pos
do k = 1, NnonSchmid(structID) do k = 1,lattice_NnonSchmid(structID)
tau_slip_pos = tau_slip_pos + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* & tau_slip_pos = tau_slip_pos + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* &
dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k,index_myFamily+i,structID)) dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k,index_myFamily+i,structID))
tau_slip_neg = tau_slip_neg + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* & tau_slip_neg = tau_slip_neg + constitutive_phenopowerlaw_nonSchmidCoeff(k,matID)* &

9
code/constitutive_titanmod.f90
View File

@ -1246,8 +1246,7 @@ subroutine constitutive_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,&
lattice_maxNslipFamily, & lattice_maxNslipFamily, &
lattice_maxNtwinFamily, & lattice_maxNtwinFamily, &
lattice_NslipSystem, & lattice_NslipSystem, &
lattice_NtwinSystem, & lattice_NtwinSystem
NnonSchmid
use mesh, only: & use mesh, only: &
mesh_NcpElems, & mesh_NcpElems, &
mesh_maxNips mesh_maxNips
@ -1441,14 +1440,14 @@ subroutine constitutive_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,&
!************************************************* !*************************************************
!sumf=0.0_pReal !sumf=0.0_pReal
!* Plastic velocity gradient for dislocation glide !* Plastic velocity gradient for dislocation glide
Lp = Lp + (1.0_pReal - sumf)*gdot_slip(j)*lattice_Sslip(:,:,index_myFamily+i,myStructure) Lp = Lp + (1.0_pReal - sumf)*gdot_slip(j)*lattice_Sslip(1:3,1:3,1,index_myFamily+i,myStructure)
!* Calculation of the tangent of Lp !* Calculation of the tangent of Lp
forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) & forall (k=1_pInt:3_pInt,l=1_pInt:3_pInt,m=1_pInt:3_pInt,n=1_pInt:3_pInt) &
dLp_dTstar3333(k,l,m,n) = & dLp_dTstar3333(k,l,m,n) = &
dLp_dTstar3333(k,l,m,n) + dgdot_dtauslip(j)*& dLp_dTstar3333(k,l,m,n) + dgdot_dtauslip(j)*&
lattice_Sslip(k,l,index_myFamily+i,myStructure)*& lattice_Sslip(k,l,1,index_myFamily+i,myStructure)*&
lattice_Sslip(m,n,index_myFamily+i,myStructure) lattice_Sslip(m,n,1,index_myFamily+i,myStructure)
enddo enddo
enddo enddo

103
code/lattice.f90
View File

@ -38,7 +38,7 @@ module lattice
lattice_maxNslip = 33_pInt, & !< max # of slip systems over lattice structures lattice_maxNslip = 33_pInt, & !< max # of slip systems over lattice structures
lattice_maxNtwin = 24_pInt, & !< max # of twin systems over lattice structures lattice_maxNtwin = 24_pInt, & !< max # of twin systems over lattice structures
lattice_maxNinteraction = 42_pInt, & !< max # of interaction types (in hardening matrix part) lattice_maxNinteraction = 42_pInt, & !< max # of interaction types (in hardening matrix part)
lattice_maxNonSchmid = 6_pInt !< max # of non schmid contributions over lattice structures lattice_maxNnonSchmid = 6_pInt !< max # of non schmid contributions over lattice structures
integer(pInt), allocatable, dimension(:,:), protected, public :: & integer(pInt), allocatable, dimension(:,:), protected, public :: &
lattice_NslipSystem, & !< total # of slip systems in each family lattice_NslipSystem, & !< total # of slip systems in each family
@ -51,14 +51,16 @@ module lattice
lattice_interactionTwinTwin !< Twin--twin interaction type lattice_interactionTwinTwin !< Twin--twin interaction type
real(pReal), allocatable, dimension(:,:,:,:,:), protected, public :: &
lattice_Sslip !< Schmid and non-Schmid matrices
real(pReal), allocatable, dimension(:,:,:,:), protected, public :: & real(pReal), allocatable, dimension(:,:,:,:), protected, public :: &
lattice_Sslip_v, & lattice_Sslip_v !< Mandel notation of lattice_Sslip
lattice_Sslip !< Schmid matrices, normal, shear direction and d x n of slip systems
real(pReal), allocatable, dimension(:,:,:), protected, public :: & real(pReal), allocatable, dimension(:,:,:), protected, public :: &
lattice_sn, & lattice_sn, & !< normal direction of slip system
lattice_sd, & lattice_sd, & !< slip direction of slip system
lattice_st lattice_st !< sd x sn
! rotation and Schmid matrices, normal, shear direction and d x n of twin systems ! rotation and Schmid matrices, normal, shear direction and d x n of twin systems
real(pReal), allocatable, dimension(:,:,:,:), protected, public :: & real(pReal), allocatable, dimension(:,:,:,:), protected, public :: &
@ -85,7 +87,7 @@ module lattice
interactionTwinTwin interactionTwinTwin
integer(pInt), allocatable, dimension(:), protected, public :: & integer(pInt), allocatable, dimension(:), protected, public :: &
NnonSchmid !< total # of non-Schmid contributions for each structure lattice_NnonSchmid !< total # of non-Schmid contributions for each structure
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! fcc (1) ! fcc (1)
@ -97,7 +99,8 @@ module lattice
integer(pInt), parameter, private :: & integer(pInt), parameter, private :: &
lattice_fcc_Nslip = 12_pInt, & ! sum(lattice_fcc_NslipSystem), & !< total # of slip systems for fcc lattice_fcc_Nslip = 12_pInt, & ! sum(lattice_fcc_NslipSystem), & !< total # of slip systems for fcc
lattice_fcc_Ntwin = 12_pInt ! sum(lattice_fcc_NtwinSystem) !< total # of twin systems for fcc lattice_fcc_Ntwin = 12_pInt, & ! sum(lattice_fcc_NtwinSystem) !< total # of twin systems for fcc
lattice_fcc_NnonSchmid = 0_pInt !< total # of non-Schmid contributions for fcc
integer(pInt), private :: & integer(pInt), private :: &
lattice_fcc_Nstructure = 0_pInt lattice_fcc_Nstructure = 0_pInt
@ -224,10 +227,6 @@ module lattice
2,2,2,2,2,2,2,2,2,1,1,1 & 2,2,2,2,2,2,2,2,2,1,1,1 &
],pInt),[lattice_fcc_Ntwin,lattice_fcc_Ntwin],order=[2,1]) !< Twin--twin interaction types for fcc ],pInt),[lattice_fcc_Ntwin,lattice_fcc_Ntwin],order=[2,1]) !< Twin--twin interaction types for fcc
integer(pInt), parameter, private :: NnonSchmid_fcc = 0_pInt !< total # of non-Schmid contributions for fcc
real(pReal), dimension(3,3,2,NnonSchmid_fcc,lattice_fcc_Nslip), parameter, private :: &
lattice_nonSchmid_fcc = 0.0_pReal ! reshape([],[3,3,2,NnonSchmid_fcc,lattice_fcc_Nslip]) !< Tensor for each non-Schmid contribution for fcc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -240,7 +239,8 @@ module lattice
integer(pInt), parameter, private :: & integer(pInt), parameter, private :: &
lattice_bcc_Nslip = 24_pInt, & ! sum(lattice_bcc_NslipSystem), & !< total # of slip systems for bcc lattice_bcc_Nslip = 24_pInt, & ! sum(lattice_bcc_NslipSystem), & !< total # of slip systems for bcc
lattice_bcc_Ntwin = 12_pInt ! sum(lattice_bcc_NtwinSystem) !< total # of twin systems for bcc lattice_bcc_Ntwin = 12_pInt, & ! sum(lattice_bcc_NtwinSystem) !< total # of twin systems for bcc
lattice_bcc_NnonSchmid = 6_pInt !< # of non-Schmid contributions for bcc. 6 known non schmid contributions for BCC (A. Koester, A. Ma, A. Hartmaier 2012)
integer(pInt), private :: & integer(pInt), private :: &
lattice_bcc_Nstructure = 0_pInt lattice_bcc_Nstructure = 0_pInt
@ -419,10 +419,7 @@ module lattice
!< 1: self interaction !< 1: self interaction
!< 2: collinear interaction !< 2: collinear interaction
!< 3: other interaction !< 3: other interaction
integer(pInt), parameter, private :: NnonSchmid_bcc = 0_pInt !< # of non-Schmid contributions for bcc
real(pReal), dimension(3,3,2,NnonSchmid_bcc,lattice_bcc_Nslip), parameter, private :: &
lattice_nonSchmid_bcc = 0.0_pReal ! reshape([],[3,3,2,NnonSchmid_bcc,lattice_bcc_Nslip]) !< Tensor for each non-Schmid contribution for bcc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! hex (3+) ! hex (3+)
@ -434,7 +431,8 @@ module lattice
integer(pInt), parameter , private :: & integer(pInt), parameter , private :: &
lattice_hex_Nslip = 33_pInt, & ! sum(lattice_hex_NslipSystem), !< total # of slip systems for hex lattice_hex_Nslip = 33_pInt, & ! sum(lattice_hex_NslipSystem), !< total # of slip systems for hex
lattice_hex_Ntwin = 24_pInt ! sum(lattice_hex_NtwinSystem) !< total # of twin systems for hex lattice_hex_Ntwin = 24_pInt, & ! sum(lattice_hex_NtwinSystem) !< total # of twin systems for hex
lattice_hex_NnonSchmid = 0_pInt !< # of non-Schmid contributions for hex
integer(pInt), private :: & integer(pInt), private :: &
lattice_hex_Nstructure = 0_pInt lattice_hex_Nstructure = 0_pInt
@ -689,10 +687,7 @@ module lattice
20,20,20,20,20,20, 19,19,19,19,19,19, 18,18,18,18,18,18, 17,17,17,17,17, 4 & 20,20,20,20,20,20, 19,19,19,19,19,19, 18,18,18,18,18,18, 17,17,17,17,17, 4 &
],pInt),[lattice_hex_Ntwin,lattice_hex_Ntwin],order=[2,1]) !< Twin--slip interaction types for hex (isotropic, 16 in total) ],pInt),[lattice_hex_Ntwin,lattice_hex_Ntwin],order=[2,1]) !< Twin--slip interaction types for hex (isotropic, 16 in total)
integer(pInt), parameter, private :: NnonSchmid_hex = 0_pInt !< # of non-Schmid contributions for hex
real(pReal), dimension(3,3,2,NnonSchmid_hex,lattice_hex_Nslip), parameter, private :: &
lattice_nonSchmid_hex = 0.0_pReal ! reshape([],[3,3,2,NnonSchmid_hex,lattice_hex_Nslip]) !< Tensor for each non-Schmid contribution for hex
public :: & public :: &
lattice_init, & lattice_init, &
@ -745,9 +740,9 @@ subroutine lattice_init
write(6,'(a16,1x,i5,/)') ' # structures:',lattice_Nstructure write(6,'(a16,1x,i5,/)') ' # structures:',lattice_Nstructure
endif endif
allocate(NnonSchmid(lattice_Nstructure)); NnonSchmid = 0_pInt allocate(lattice_NnonSchmid(lattice_Nstructure)); lattice_NnonSchmid = 0_pInt
allocate(lattice_Sslip(3,3,lattice_maxNslip,lattice_Nstructure)); lattice_Sslip = 0.0_pReal allocate(lattice_Sslip(3,3,1+2*lattice_maxNnonSchmid,lattice_maxNslip,lattice_Nstructure)); lattice_Sslip = 0.0_pReal
allocate(lattice_Sslip_v(6,1+2*lattice_maxNonSchmid,lattice_maxNslip,lattice_Nstructure)); lattice_Sslip_v = 0.0_pReal allocate(lattice_Sslip_v(6,1+2*lattice_maxNnonSchmid,lattice_maxNslip,lattice_Nstructure)); lattice_Sslip_v = 0.0_pReal
allocate(lattice_sd(3,lattice_maxNslip,lattice_Nstructure)); lattice_sd = 0.0_pReal allocate(lattice_sd(3,lattice_maxNslip,lattice_Nstructure)); lattice_sd = 0.0_pReal
allocate(lattice_st(3,lattice_maxNslip,lattice_Nstructure)); lattice_st = 0.0_pReal allocate(lattice_st(3,lattice_maxNslip,lattice_Nstructure)); lattice_st = 0.0_pReal
allocate(lattice_sn(3,lattice_maxNslip,lattice_Nstructure)); lattice_sn = 0.0_pReal allocate(lattice_sn(3,lattice_maxNslip,lattice_Nstructure)); lattice_sn = 0.0_pReal
@ -784,6 +779,7 @@ integer(pInt) function lattice_initializeStructure(struct,CoverA)
math_vectorproduct, & math_vectorproduct, &
math_tensorproduct, & math_tensorproduct, &
math_norm3, & math_norm3, &
math_mul33x3, &
math_trace33, & math_trace33, &
math_symmetric33, & math_symmetric33, &
math_Mandel33to6, & math_Mandel33to6, &
@ -795,9 +791,13 @@ integer(pInt) function lattice_initializeStructure(struct,CoverA)
implicit none implicit none
character(len=*) struct character(len=*) struct
real(pReal) CoverA real(pReal) CoverA
real(pReal), dimension(3) :: sdU = 0.0_pReal, &
snU = 0.0_pReal, &
np = 0.0_pReal, &
nn = 0.0_pReal
real(pReal), dimension(3,lattice_maxNslip) :: sd = 0.0_pReal, & real(pReal), dimension(3,lattice_maxNslip) :: sd = 0.0_pReal, &
sn = 0.0_pReal sn = 0.0_pReal
real(pReal), dimension(12,lattice_maxNonSchmid,lattice_maxNslip) :: sns = 0.0_pReal real(pReal), dimension(3,3,2,lattice_maxNnonSchmid,lattice_maxNslip) :: sns = 0.0_pReal
real(pReal), dimension(3,lattice_maxNtwin) :: td = 0.0_pReal, & real(pReal), dimension(3,lattice_maxNtwin) :: td = 0.0_pReal, &
tn = 0.0_pReal tn = 0.0_pReal
real(pReal), dimension(lattice_maxNtwin) :: ts = 0.0_pReal real(pReal), dimension(lattice_maxNtwin) :: ts = 0.0_pReal
@ -818,13 +818,13 @@ integer(pInt) function lattice_initializeStructure(struct,CoverA)
lattice_fcc_Nstructure = lattice_fcc_Nstructure + 1_pInt ! count fcc instances lattice_fcc_Nstructure = lattice_fcc_Nstructure + 1_pInt ! count fcc instances
if (lattice_fcc_Nstructure == 1_pInt) then ! me is first fcc structure if (lattice_fcc_Nstructure == 1_pInt) then ! me is first fcc structure
processMe = .true. processMe = .true.
NnonSchmid(myStructure) = NnonSchmid_fcc ! Currently no known non schmid contributions for FCC (to be changed later) lattice_NnonSchmid(myStructure) = lattice_fcc_NnonSchmid ! Currently no known non schmid contributions for FCC (to be changed later)
do i = 1_pInt,myNslip ! assign slip system vectors do i = 1_pInt,myNslip ! assign slip system vectors
sd(1:3,i) = lattice_fcc_systemSlip(1:3,i) sd(1:3,i) = lattice_fcc_systemSlip(1:3,i)
sn(1:3,i) = lattice_fcc_systemSlip(4:6,i) sn(1:3,i) = lattice_fcc_systemSlip(4:6,i)
do j = 1_pInt, NnonSchmid_fcc do j = 1_pInt,lattice_fcc_NnonSchmid
sns(1:6,j,i) = math_Mandel33to6(lattice_nonSchmid_fcc(1:3,1:3,1,j,i)) sns(1:3,1:3,1,j,i) = 0.0_pReal
sns(7:12,j,i) = math_Mandel33to6(lattice_nonSchmid_fcc(1:3,1:3,2,j,i)) sns(1:3,1:3,2,j,i) = 0.0_pReal
enddo enddo
enddo enddo
do i = 1_pInt,myNtwin ! assign twin system vectors and shears do i = 1_pInt,myNtwin ! assign twin system vectors and shears
@ -847,14 +847,26 @@ integer(pInt) function lattice_initializeStructure(struct,CoverA)
lattice_bcc_Nstructure = lattice_bcc_Nstructure + 1_pInt ! count bcc instances lattice_bcc_Nstructure = lattice_bcc_Nstructure + 1_pInt ! count bcc instances
if (lattice_bcc_Nstructure == 1_pInt) then ! me is first bcc structure if (lattice_bcc_Nstructure == 1_pInt) then ! me is first bcc structure
processMe = .true. processMe = .true.
NnonSchmid(myStructure) = NnonSchmid_BCC ! 5 known non schmid contributions for BCC (A. Koester, A. Ma, A. Hartmaier 2012) lattice_NnonSchmid(myStructure) = lattice_bcc_NnonSchmid
do i = 1_pInt,myNslip ! assign slip system vectors do i = 1_pInt,myNslip ! assign slip system vectors
sd(1:3,i) = lattice_bcc_systemSlip(1:3,i) sd(1:3,i) = lattice_bcc_systemSlip(1:3,i)
sn(1:3,i) = lattice_bcc_systemSlip(4:6,i) sn(1:3,i) = lattice_bcc_systemSlip(4:6,i)
do j = 1_pInt, NnonSchmid_bcc sdU = sd(1:3,i) / math_norm3(sd(1:3,i))
sns(1:6,j,i) = math_Mandel33to6(lattice_nonSchmid_bcc(1:3,1:3,1,j,i)) snU = sn(1:3,i) / math_norm3(sn(1:3,i))
sns(7:12,j,i) = math_Mandel33to6(lattice_nonSchmid_bcc(1:3,1:3,2,j,i)) np = math_mul33x3(math_axisAngleToR(sdU,60.0_pReal*INRAD), snU)
enddo nn = math_mul33x3(math_axisAngleToR(-sdU,60.0_pReal*INRAD), snU)
sns(1:3,1:3,1,1,i) = math_tensorproduct(sdU, np)
sns(1:3,1:3,2,1,i) = math_tensorproduct(-sdU, nn)
sns(1:3,1:3,1,2,i) = math_tensorproduct(math_vectorproduct(snU, sdU), snU)
sns(1:3,1:3,2,2,i) = math_tensorproduct(math_vectorproduct(snU, -sdU), snU)
sns(1:3,1:3,1,3,i) = math_tensorproduct(math_vectorproduct(np, sdU), np)
sns(1:3,1:3,2,3,i) = math_tensorproduct(math_vectorproduct(nn, -sdU), nn)
sns(1:3,1:3,1,4,i) = math_tensorproduct(snU, snU)
sns(1:3,1:3,2,4,i) = math_tensorproduct(snU, snU)
sns(1:3,1:3,1,5,i) = math_tensorproduct(math_vectorproduct(snU, sdU), math_vectorproduct(snU, sdU))
sns(1:3,1:3,2,5,i) = math_tensorproduct(math_vectorproduct(snU, -sdU), math_vectorproduct(snU, -sdU))
sns(1:3,1:3,1,6,i) = math_tensorproduct(sdU, sdU)
sns(1:3,1:3,2,6,i) = math_tensorproduct(-sdU, -sdU)
enddo enddo
do i = 1_pInt,myNtwin ! assign twin system vectors and shears do i = 1_pInt,myNtwin ! assign twin system vectors and shears
td(1:3,i) = lattice_bcc_systemTwin(1:3,i) td(1:3,i) = lattice_bcc_systemTwin(1:3,i)
@ -876,7 +888,7 @@ integer(pInt) function lattice_initializeStructure(struct,CoverA)
myNslip = lattice_hex_Nslip ! overall number of slip systems myNslip = lattice_hex_Nslip ! overall number of slip systems
myNtwin = lattice_hex_Ntwin ! overall number of twin systems myNtwin = lattice_hex_Ntwin ! overall number of twin systems
processMe = .true. processMe = .true.
NnonSchmid(myStructure) = NnonSchmid_hex ! Currently no known non schmid contributions for hex (to be changed later) lattice_NnonSchmid(myStructure) = lattice_hex_NnonSchmid ! Currently no known non schmid contributions for hex (to be changed later)
! converting from 4 axes coordinate system (a1=a2=a3=c) to ortho-hexgonal system (a, b, c) ! converting from 4 axes coordinate system (a1=a2=a3=c) to ortho-hexgonal system (a, b, c)
do i = 1_pInt,myNslip do i = 1_pInt,myNslip
sd(1,i) = lattice_hex_systemSlip(1,i)*1.5_pReal ! direction [uvtw]->[3u/2 (u+2v)*sqrt(3)/2 w*(c/a)] sd(1,i) = lattice_hex_systemSlip(1,i)*1.5_pReal ! direction [uvtw]->[3u/2 (u+2v)*sqrt(3)/2 w*(c/a)]
@ -885,9 +897,9 @@ integer(pInt) function lattice_initializeStructure(struct,CoverA)
sn(1,i) = lattice_hex_systemSlip(5,i) ! plane (hkil)->(h (h+2k)/sqrt(3) l/(c/a)) sn(1,i) = lattice_hex_systemSlip(5,i) ! plane (hkil)->(h (h+2k)/sqrt(3) l/(c/a))
sn(2,i) = (lattice_hex_systemSlip(5,i)+2.0_pReal*lattice_hex_systemSlip(6,i))/sqrt(3.0_pReal) sn(2,i) = (lattice_hex_systemSlip(5,i)+2.0_pReal*lattice_hex_systemSlip(6,i))/sqrt(3.0_pReal)
sn(3,i) = lattice_hex_systemSlip(8,i)/CoverA sn(3,i) = lattice_hex_systemSlip(8,i)/CoverA
do j = 1_pInt, NnonSchmid_hex do j = 1_pInt,lattice_hex_NnonSchmid
sns(1:6,j,i) = math_Mandel33to6(lattice_nonSchmid_hex(1:3,1:3,1,j,i)) sns(1:3,1:3,1,j,i) = 0.0_pReal
sns(7:12,j,i) = math_Mandel33to6(lattice_nonSchmid_hex(1:3,1:3,2,j,i)) sns(1:3,1:3,2,j,i) = 0.0_pReal
enddo enddo
enddo enddo
do i = 1_pInt,myNtwin do i = 1_pInt,myNtwin
@ -925,14 +937,17 @@ integer(pInt) function lattice_initializeStructure(struct,CoverA)
lattice_sn(1:3,i,myStructure) = sn(1:3,i)/math_norm3(sn(1:3,i)) ! make unit vector lattice_sn(1:3,i,myStructure) = sn(1:3,i)/math_norm3(sn(1:3,i)) ! make unit vector
lattice_st(1:3,i,myStructure) = math_vectorproduct(lattice_sd(1:3,i,myStructure), & lattice_st(1:3,i,myStructure) = math_vectorproduct(lattice_sd(1:3,i,myStructure), &
lattice_sn(1:3,i,myStructure)) lattice_sn(1:3,i,myStructure))
lattice_Sslip(1:3,1:3,i,myStructure) = math_tensorproduct(lattice_sd(1:3,i,myStructure), & lattice_Sslip(1:3,1:3,1,i,myStructure) = math_tensorproduct(lattice_sd(1:3,i,myStructure), &
lattice_sn(1:3,i,myStructure)) lattice_sn(1:3,i,myStructure))
lattice_Sslip_v(1:6,1,i,myStructure) = math_Mandel33to6(math_symmetric33(lattice_Sslip(1:3,1:3,i,myStructure))) do j = 1_pInt,lattice_NnonSchmid(myStructure)
do j = 1_pInt, NnonSchmid(myStructure) lattice_Sslip(1:3,1:3,2*j ,i,myStructure) = sns(1:3,1:3,1,j,i)
lattice_Sslip_v(1:6,2*j,i,myStructure) = sns(1:6,j,i) lattice_Sslip(1:3,1:3,2*j+1,i,myStructure) = sns(1:3,1:3,2,j,i)
lattice_Sslip_v(1:6,2*j+1,i,myStructure) = sns(7:12,j,i)
enddo enddo
if (abs(math_trace33(lattice_Sslip(1:3,1:3,i,myStructure))) > 1.0e-8_pReal) & do j = 1_pInt,1_pInt+2_pInt*lattice_NnonSchmid(myStructure)
lattice_Sslip_v(1:6,j,i,myStructure) = &
math_Mandel33to6(math_symmetric33(lattice_Sslip(1:3,1:3,j,i,myStructure)))
enddo
if (abs(math_trace33(lattice_Sslip(1:3,1:3,1,i,myStructure))) > 1.0e-8_pReal) &
call IO_error(0_pInt,myStructure,i,0_pInt,ext_msg = 'dilatational slip Schmid matrix') call IO_error(0_pInt,myStructure,i,0_pInt,ext_msg = 'dilatational slip Schmid matrix')
enddo enddo
do i = 1_pInt,myNtwin ! store twin system vectors and Schmid plus rotation matrix for my structure do i = 1_pInt,myNtwin ! store twin system vectors and Schmid plus rotation matrix for my structure