In dislo:

- model for stacking fault energy computation w.r.t temperature and orientation
- classical thermodynamical approach to twin nucleation (delta_E)
- rename some parameters

In crystallite:
- add msg='ok', even if components in consistent tangent are not re-computed
- allow for pseudo-convergence in the outer-loop (case of slow convergence of dislocation densities)

in mattex:
- add new parameters
This commit is contained in:
Luc Hantcherli 2008-06-17 14:40:46 +00:00
parent 2472de77c2
commit 434ac7f06e
3 changed files with 304 additions and 178 deletions

View File

@ -24,6 +24,10 @@ character(len=300), parameter :: mattexFile = 'mattex.mpie'
real(pReal), parameter :: attack_frequency = 1.0e10_pReal real(pReal), parameter :: attack_frequency = 1.0e10_pReal
!* Physical parameter, Boltzmann constant in J/Kelvin !* Physical parameter, Boltzmann constant in J/Kelvin
real(pReal), parameter :: kB = 1.38e-23_pReal real(pReal), parameter :: kB = 1.38e-23_pReal
!* Physical parameter, Avogadro number in 1/mol
real(pReal), parameter :: avogadro = 6.022e23_pReal
!* Physical parameter, Gaz constant in J.mol/Kelvin
real(pReal), parameter :: Rgaz = 8.314_pReal
!************************************* !*************************************
!* Definition of material properties * !* Definition of material properties *
@ -31,7 +35,7 @@ real(pReal), parameter :: kB = 1.38e-23_pReal
!* Number of materials !* Number of materials
integer(pInt) material_maxN integer(pInt) material_maxN
!* Lattice structure and number of selected slip or twin systems per material !* Lattice structure and number of selected slip or twin systems per material
integer(pInt), dimension(:) , allocatable :: material_LatticeStructure integer(pInt), dimension(:) , allocatable :: material_CrystalStructure
integer(pInt), dimension(:) , allocatable :: material_Nslip integer(pInt), dimension(:) , allocatable :: material_Nslip
integer(pInt), dimension(:) , allocatable :: material_Ntwin integer(pInt), dimension(:) , allocatable :: material_Ntwin
!* Maximum number of selected slip or twin systems over materials !* Maximum number of selected slip or twin systems over materials
@ -47,15 +51,15 @@ real(pReal), dimension(:) , allocatable :: material_Gmod
real(pReal), dimension(:,:,:) , allocatable :: material_Cslip_66 real(pReal), dimension(:,:,:) , allocatable :: material_Cslip_66
real(preal), dimension(:,:,:,:) , allocatable :: material_Ctwin_66 real(preal), dimension(:,:,:,:) , allocatable :: material_Ctwin_66
!* Visco-plastic material parameters !* Visco-plastic material parameters
real(pReal) material_sfe
real(pReal), dimension(:) , allocatable :: material_rho0 real(pReal), dimension(:) , allocatable :: material_rho0
real(pReal), dimension(:) , allocatable :: material_bg real(pReal), dimension(:) , allocatable :: material_bg
real(pReal), dimension(:) , allocatable :: material_Qedge real(pReal), dimension(:) , allocatable :: material_Qedge
real(pReal), dimension(:) , allocatable :: material_tau0
real(pReal), dimension(:) , allocatable :: material_GrainSize real(pReal), dimension(:) , allocatable :: material_GrainSize
real(pReal), dimension(:) , allocatable :: material_StackSize real(pReal), dimension(:) , allocatable :: material_StackSize
real(pReal), dimension(:) , allocatable :: material_ActivationLength real(pReal), dimension(:) , allocatable :: material_ActivationLength
real(pReal), dimension(:) , allocatable :: material_TwinSaturation real(pReal), dimension(:) , allocatable :: material_TwinSaturation
real(pReal), dimension(:) , allocatable :: material_twin_res real(pReal), dimension(:) , allocatable :: material_SiteScaling
real(pReal), dimension(:) , allocatable :: material_c1 real(pReal), dimension(:) , allocatable :: material_c1
real(pReal), dimension(:) , allocatable :: material_c2 real(pReal), dimension(:) , allocatable :: material_c2
real(pReal), dimension(:) , allocatable :: material_c3 real(pReal), dimension(:) , allocatable :: material_c3
@ -63,8 +67,8 @@ real(pReal), dimension(:) , allocatable :: material_c4
real(pReal), dimension(:) , allocatable :: material_c5 real(pReal), dimension(:) , allocatable :: material_c5
real(pReal), dimension(:) , allocatable :: material_c6 real(pReal), dimension(:) , allocatable :: material_c6
real(pReal), dimension(:) , allocatable :: material_c7 real(pReal), dimension(:) , allocatable :: material_c7
real(pReal), dimension(:) , allocatable :: material_c8 real(pReal), dimension(:) , allocatable :: material_q1
real(pReal), dimension(:) , allocatable :: material_c9 real(pReal), dimension(:) , allocatable :: material_q2
real(pReal), dimension(:,:) , allocatable :: material_SlipIntCoeff real(pReal), dimension(:,:) , allocatable :: material_SlipIntCoeff
!************************************ !************************************
@ -126,7 +130,8 @@ real(pReal), dimension(:) , allocatable :: constitutive_rho_f
real(pReal), dimension(:) , allocatable :: constitutive_rho_p real(pReal), dimension(:) , allocatable :: constitutive_rho_p
real(pReal), dimension(:) , allocatable :: constitutive_g0_slip real(pReal), dimension(:) , allocatable :: constitutive_g0_slip
real(pReal), dimension(:) , allocatable :: constitutive_twin_volume real(pReal), dimension(:) , allocatable :: constitutive_twin_volume
real(pReal), dimension(:) , allocatable :: constitutive_inv_intertwin_len real(pReal), dimension(:) , allocatable :: constitutive_inv_intertwin_len_s
real(pReal), dimension(:) , allocatable :: constitutive_inv_intertwin_len_t
real(pReal), dimension(:) , allocatable :: constitutive_twin_mfp real(pReal), dimension(:) , allocatable :: constitutive_twin_mfp
!************************************ !************************************
@ -318,8 +323,8 @@ do while(.true.)
!$OMP CRITICAL (write2out) !$OMP CRITICAL (write2out)
select case(tag) select case(tag)
case ('lattice_structure') case ('lattice_structure')
material_LatticeStructure(section)=IO_intValue(line,positions,2) material_CrystalStructure(section)=IO_intValue(line,positions,2)
write(6,*) 'lattice_structure', material_LatticeStructure(section) write(6,*) 'lattice_structure', material_CrystalStructure(section)
case ('nslip') case ('nslip')
material_Nslip(section)=IO_intValue(line,positions,2) material_Nslip(section)=IO_intValue(line,positions,2)
write(6,*) 'nslip', material_Nslip(section) write(6,*) 'nslip', material_Nslip(section)
@ -355,24 +360,18 @@ do while(.true.)
case ('qedge') case ('qedge')
material_Qedge(section)=IO_floatValue(line,positions,2) material_Qedge(section)=IO_floatValue(line,positions,2)
write(6,*) 'Qedge', material_Qedge(section) write(6,*) 'Qedge', material_Qedge(section)
case ('tau0')
material_tau0(section)=IO_floatValue(line,positions,2)
write(6,*) 'tau0', material_tau0(section)
case ('grain_size') case ('grain_size')
material_GrainSize(section)=IO_floatValue(line,positions,2) material_GrainSize(section)=IO_floatValue(line,positions,2)
write(6,*) 'grain_size', material_GrainSize(section) write(6,*) 'grain_size', material_GrainSize(section)
case ('stack_size') case ('stack_size')
material_StackSize(section)=IO_floatValue(line,positions,2) material_StackSize(section)=IO_floatValue(line,positions,2)
write(6,*) 'stack_size', material_StackSize(section) write(6,*) 'stack_size', material_StackSize(section)
case ('d_star')
material_ActivationLength(section)=IO_floatValue(line,positions,2)
write(6,*) 'activation length', material_ActivationLength(section)
case ('f_sat') case ('f_sat')
material_TwinSaturation(section)=IO_floatValue(line,positions,2) material_TwinSaturation(section)=IO_floatValue(line,positions,2)
write(6,*) 'twin saturation', material_TwinSaturation(section) write(6,*) 'twin saturation', material_TwinSaturation(section)
case ('twin_resistance') case ('site_scaling')
material_twin_res(section)=IO_floatValue(line,positions,2) material_SiteScaling(section)=IO_floatValue(line,positions,2)
write(6,*) 'twin_resistance', material_twin_res(section) write(6,*) 'site_scaling', material_SiteScaling(section)
case ('c1') case ('c1')
material_c1(section)=IO_floatValue(line,positions,2) material_c1(section)=IO_floatValue(line,positions,2)
write(6,*) 'c1', material_c1(section) write(6,*) 'c1', material_c1(section)
@ -394,12 +393,12 @@ do while(.true.)
case ('c7') case ('c7')
material_c7(section)=IO_floatValue(line,positions,2) material_c7(section)=IO_floatValue(line,positions,2)
write(6,*) 'c7', material_c7(section) write(6,*) 'c7', material_c7(section)
case ('c8') case ('q1')
material_c8(section)=IO_floatValue(line,positions,2) material_q1(section)=IO_floatValue(line,positions,2)
write(6,*) 'c8', material_c8(section) write(6,*) 'q1', material_q1(section)
case ('c9') case ('q2')
material_c9(section)=IO_floatValue(line,positions,2) material_q2(section)=IO_floatValue(line,positions,2)
write(6,*) 'c9', material_c9(section) write(6,*) 'q2', material_q2(section)
end select end select
!$OMP END CRITICAL (write2out) !$OMP END CRITICAL (write2out)
endif endif
@ -536,7 +535,7 @@ do while (part/='')
end select end select
enddo enddo
!* Array allocation !* Array allocation
allocate(material_LatticeStructure(material_maxN)) ; material_LatticeStructure=0_pInt allocate(material_CrystalStructure(material_maxN)) ; material_CrystalStructure=0_pInt
allocate(material_Nslip(material_maxN)) ; material_Nslip=0_pInt allocate(material_Nslip(material_maxN)) ; material_Nslip=0_pInt
allocate(material_Ntwin(material_maxN)) ; material_Ntwin=0_pInt allocate(material_Ntwin(material_maxN)) ; material_Ntwin=0_pInt
allocate(material_C11(material_maxN)) ; material_C11=0.0_pReal allocate(material_C11(material_maxN)) ; material_C11=0.0_pReal
@ -551,12 +550,11 @@ allocate(material_rho0(material_maxN)) ; mate
allocate(material_SlipIntCoeff(lattice_MaxMaxNslipOfStructure,material_maxN)) ; material_SlipIntCoeff=0.0_pReal allocate(material_SlipIntCoeff(lattice_MaxMaxNslipOfStructure,material_maxN)) ; material_SlipIntCoeff=0.0_pReal
allocate(material_bg(material_maxN)) ; material_bg=0.0_pReal allocate(material_bg(material_maxN)) ; material_bg=0.0_pReal
allocate(material_Qedge(material_maxN)) ; material_Qedge=0.0_pReal allocate(material_Qedge(material_maxN)) ; material_Qedge=0.0_pReal
allocate(material_tau0(material_maxN)) ; material_tau0=0.0_pReal
allocate(material_GrainSize(material_maxN)) ; material_GrainSize=0.0_pReal allocate(material_GrainSize(material_maxN)) ; material_GrainSize=0.0_pReal
allocate(material_StackSize(material_maxN)) ; material_StackSize=0.0_pReal allocate(material_StackSize(material_maxN)) ; material_StackSize=0.0_pReal
allocate(material_ActivationLength(material_maxN)) ; material_ActivationLength=0.0_pReal allocate(material_ActivationLength(material_maxN)) ; material_ActivationLength=0.0_pReal
allocate(material_TwinSaturation(material_maxN)) ; material_TwinSaturation=0.0_pReal allocate(material_TwinSaturation(material_maxN)) ; material_TwinSaturation=0.0_pReal
allocate(material_twin_res(material_maxN)) ; material_twin_res=0.0_pReal allocate(material_SiteScaling(material_maxN)) ; material_SiteScaling=0.0_pReal
allocate(material_c1(material_maxN)) ; material_c1=0.0_pReal allocate(material_c1(material_maxN)) ; material_c1=0.0_pReal
allocate(material_c2(material_maxN)) ; material_c2=0.0_pReal allocate(material_c2(material_maxN)) ; material_c2=0.0_pReal
allocate(material_c3(material_maxN)) ; material_c3=0.0_pReal allocate(material_c3(material_maxN)) ; material_c3=0.0_pReal
@ -564,8 +562,8 @@ allocate(material_c4(material_maxN))
allocate(material_c5(material_maxN)) ; material_c5=0.0_pReal allocate(material_c5(material_maxN)) ; material_c5=0.0_pReal
allocate(material_c6(material_maxN)) ; material_c6=0.0_pReal allocate(material_c6(material_maxN)) ; material_c6=0.0_pReal
allocate(material_c7(material_maxN)) ; material_c7=0.0_pReal allocate(material_c7(material_maxN)) ; material_c7=0.0_pReal
allocate(material_c8(material_maxN)) ; material_c8=0.0_pReal allocate(material_q1(material_maxN)) ; material_q1=0.0_pReal
allocate(material_c9(material_maxN)) ; material_c9=0.0_pReal allocate(material_q2(material_maxN)) ; material_q2=0.0_pReal
allocate(texture_ODFfile(texture_maxN)) ; texture_ODFfile='' allocate(texture_ODFfile(texture_maxN)) ; texture_ODFfile=''
allocate(texture_Ngrains(texture_maxN)) ; texture_Ngrains=0_pInt allocate(texture_Ngrains(texture_maxN)) ; texture_Ngrains=0_pInt
allocate(texture_symmetry(texture_maxN)) ; texture_symmetry='' allocate(texture_symmetry(texture_maxN)) ; texture_symmetry=''
@ -612,7 +610,7 @@ close(fileunit)
!* Construction of the elasticity matrices !* Construction of the elasticity matrices
do i=1,material_maxN do i=1,material_maxN
select case (material_LatticeStructure(i)) select case (material_CrystalStructure(i))
case(1:2) ! cubic(s) case(1:2) ! cubic(s)
material_Gmod(i)=material_C44(i) material_Gmod(i)=material_C44(i)
forall(k=1:3) forall(k=1:3)
@ -711,7 +709,7 @@ material_maxNslip = maxval(material_Nslip) ! max of slip systems among mat
material_maxNtwin = maxval(material_Ntwin) ! max of twin systems among materials present material_maxNtwin = maxval(material_Ntwin) ! max of twin systems among materials present
constitutive_maxNstatevars = maxval(material_Nslip) + maxval(material_Ntwin) constitutive_maxNstatevars = maxval(material_Nslip) + maxval(material_Ntwin)
! ----------------------------------------------------------------------------------------------------------------------- ! -----------------------------------------------------------------------------------------------------------------------
constitutive_maxNresults = 24_pInt constitutive_maxNresults = 2_pInt
! ----------------------------------------------------------------------------------------------------------------------- ! -----------------------------------------------------------------------------------------------------------------------
@ -754,7 +752,8 @@ allocate(constitutive_jump_width(material_maxNslip)) ; constitutive_jump_
allocate(constitutive_activation_volume(material_maxNslip)) ; constitutive_activation_volume=0.0_pReal allocate(constitutive_activation_volume(material_maxNslip)) ; constitutive_activation_volume=0.0_pReal
allocate(constitutive_g0_slip(material_maxNslip)) ; constitutive_g0_slip=0.0_pReal allocate(constitutive_g0_slip(material_maxNslip)) ; constitutive_g0_slip=0.0_pReal
allocate(constitutive_twin_volume(material_maxNtwin)) ; constitutive_twin_volume=0.0_pReal allocate(constitutive_twin_volume(material_maxNtwin)) ; constitutive_twin_volume=0.0_pReal
allocate(constitutive_inv_intertwin_len(material_maxNtwin)) ; constitutive_inv_intertwin_len=0.0_pReal allocate(constitutive_inv_intertwin_len_s(material_maxNslip)) ; constitutive_inv_intertwin_len_s=0.0_pReal
allocate(constitutive_inv_intertwin_len_t(material_maxNtwin)) ; constitutive_inv_intertwin_len_t=0.0_pReal
allocate(constitutive_twin_mfp(material_maxNtwin)) ; constitutive_twin_mfp=0.0_pReal allocate(constitutive_twin_mfp(material_maxNtwin)) ; constitutive_twin_mfp=0.0_pReal
allocate(constitutive_tau_slip(material_maxNslip)) ; constitutive_tau_slip=0.0_pReal allocate(constitutive_tau_slip(material_maxNslip)) ; constitutive_tau_slip=0.0_pReal
allocate(constitutive_tau_twin(material_maxNtwin)) ; constitutive_tau_twin=0.0_pReal allocate(constitutive_tau_twin(material_maxNtwin)) ; constitutive_tau_twin=0.0_pReal
@ -812,7 +811,7 @@ do e=1,mesh_NcpElems
constitutive_TexVolFrac(g,i,e) = texVolfrac(s)/multiplicity(texID)/Nsym(texID) constitutive_TexVolFrac(g,i,e) = texVolfrac(s)/multiplicity(texID)/Nsym(texID)
constitutive_Nstatevars(g,i,e) = material_Nslip(matID) + material_Ntwin(matID)! number of state variables (i.e. tau_c of each slip system) constitutive_Nstatevars(g,i,e) = material_Nslip(matID) + material_Ntwin(matID)! number of state variables (i.e. tau_c of each slip system)
! ----------------------------------------------------------------------------------------------------------------------- ! -----------------------------------------------------------------------------------------------------------------------
constitutive_Nresults(g,i,e) = 24 ! number of constitutive results output by constitutive_post_results constitutive_Nresults(g,i,e) = 2 ! number of constitutive results output by constitutive_post_results
! ----------------------------------------------------------------------------------------------------------------------- ! -----------------------------------------------------------------------------------------------------------------------
constitutive_EulerAngles(:,g,i,e) = Euler(:,s) ! store initial orientation constitutive_EulerAngles(:,g,i,e) = Euler(:,s) ! store initial orientation
forall (l=1:material_Nslip(matID)) ! initialize state variables forall (l=1:material_Nslip(matID)) ! initialize state variables
@ -829,7 +828,7 @@ enddo ! cp_element
do i=1,material_maxN do i=1,material_maxN
C_3333=math_Mandel66to3333(material_Cslip_66(:,:,i)) C_3333=math_Mandel66to3333(material_Cslip_66(:,:,i))
do j=1,material_Ntwin(i) do j=1,material_Ntwin(i)
Qtwin=lattice_Qtwin(:,:,j,material_LatticeStructure(i)) Qtwin=lattice_Qtwin(:,:,j,material_CrystalStructure(i))
do k=1,3 do k=1,3
do l=1,3 do l=1,3
do m=1,3 do m=1,3
@ -872,6 +871,7 @@ do i=1,material_maxN
enddo enddo
enddo enddo
end subroutine end subroutine
@ -920,7 +920,7 @@ subroutine constitutive_Microstructure(state,Tp,ipc,ip,el)
!********************************************************************* !*********************************************************************
use prec, only: pReal,pInt use prec, only: pReal,pInt
use math, only: pi use math, only: pi
use lattice, only: lattice_TwinIntType use lattice, only: lattice_TwinIntType,lattice_SlipTwinIntType
implicit none implicit none
!* Definition of variables !* Definition of variables
@ -928,6 +928,8 @@ integer(pInt) ipc,ip,el
integer(pInt) matID,i,j,startIdxTwin integer(pInt) matID,i,j,startIdxTwin
real(pReal) Tp,Ftwin real(pReal) Tp,Ftwin
real(pReal), dimension(constitutive_Nstatevars(ipc,ip,el)) :: state real(pReal), dimension(constitutive_Nstatevars(ipc,ip,el)) :: state
real(pReal) x_fe,x_Mn,x_C,beta_mart,Tp_mart,f_mart,beta_aust,Tp_aust,f_aust
real(pReal) deltaG1,deltaG2,deltaG3,deltaG4,deltaG5
!* Get the material-ID from the triplet(ipc,ip,el) !* Get the material-ID from the triplet(ipc,ip,el)
matID = constitutive_matID(ipc,ip,el) matID = constitutive_matID(ipc,ip,el)
@ -936,10 +938,10 @@ startIdxTwin = material_Nslip(matID)
!* Quantities derived from state - slip !* Quantities derived from state - slip
!$OMP CRITICAL (evilmatmul) !$OMP CRITICAL (evilmatmul)
constitutive_rho_f=matmul(constitutive_Pforest (1:material_Nslip(matID),1:material_Nslip(matID),matID),state) constitutive_rho_f=matmul(constitutive_Pforest (1:material_Nslip(matID),1:material_Nslip(matID),matID),state)
constitutive_rho_p=matmul(constitutive_Pparallel(1:material_Nslip(matID),1:material_Nslip(matID),matID),state) constitutive_rho_p=matmul(constitutive_Pparallel(1:material_Nslip(matID),1:material_Nslip(matID),matID),state)
!$OMP END CRITICAL (evilmatmul) !$OMP END CRITICAL (evilmatmul)
do i=1,material_Nslip(matID) do i=1,material_Nslip(matID)
constitutive_passing_stress(i) = material_tau0(matID)+material_c1(matID)*material_Gmod(matID)*material_bg(matID)*& constitutive_passing_stress(i) = material_c1(matID)*material_Gmod(matID)*material_bg(matID)*&
sqrt(constitutive_rho_p(i)) sqrt(constitutive_rho_p(i))
constitutive_jump_width(i) = material_c2(matID)/sqrt(constitutive_rho_f(i)) constitutive_jump_width(i) = material_c2(matID)/sqrt(constitutive_rho_f(i))
@ -957,22 +959,70 @@ enddo
!* Quantities derived from state - twin !* Quantities derived from state - twin
Ftwin = sum(state((startIdxTwin+1):(startIdxTwin+material_Ntwin(matID)))) Ftwin = sum(state((startIdxTwin+1):(startIdxTwin+material_Ntwin(matID))))
do i=1,material_Nslip(matID)
!* Inverse of the average distance between 2 twins of the same familly
constitutive_inv_intertwin_len_s(i)=0.0_pReal
do j=1,material_Ntwin(matID)
constitutive_inv_intertwin_len_s(i)=constitutive_inv_intertwin_len_s(i)+&
(lattice_SlipTwinIntType(i,j,material_CrystalStructure(matID))*state(startIdxTwin+j))/&
(2.0_pReal*material_StackSize(matID)*(1.0_pReal-Ftwin))
enddo
enddo
do i=1,material_Ntwin(matID) do i=1,material_Ntwin(matID)
!* Inverse of the average distance between 2 twins of the same familly !* Inverse of the average distance between 2 twins of the same familly
constitutive_inv_intertwin_len(i)=0.0_pReal constitutive_inv_intertwin_len_t(i)=0.0_pReal
do j=1,material_Ntwin(matID) do j=1,material_Ntwin(matID)
constitutive_inv_intertwin_len(i)=constitutive_inv_intertwin_len(i)+& constitutive_inv_intertwin_len_t(i)=constitutive_inv_intertwin_len_t(i)+&
(lattice_TwinIntType(i,j,material_LatticeStructure(matID))*state(startIdxTwin+j))/& (lattice_TwinIntType(i,j,material_CrystalStructure(matID))*state(startIdxTwin+j))/&
(2.0_pReal*material_StackSize(matID)*(1.0_pReal-Ftwin)) (2.0_pReal*material_StackSize(matID)*(1.0_pReal-Ftwin))
enddo enddo
constitutive_twin_mfp(i)=(1.0_pReal)/((1.0_pReal/material_GrainSize(matID))+constitutive_inv_intertwin_len(i)) constitutive_twin_mfp(i)=1.0_pReal/((1.0_pReal/material_GrainSize(matID))+constitutive_inv_intertwin_len_t(i))
constitutive_twin_volume(i)=(pi/6.0_pReal)*material_StackSize(matID)*constitutive_twin_mfp(i)**2.0_pReal constitutive_twin_volume(i)=((4.0_pReal*pi)/3.0_pReal)*material_StackSize(matID)*constitutive_twin_mfp(i)**2.0_pReal
enddo enddo
!* Stacking fault energy as function of temperature (see Allain PhD Thesis p40-42) *
x_Fe=0.774_pReal ! atomic %
x_Mn=0.218_pReal ! atomic %
x_C =0.027_pReal ! atomic %
!* Chemical contribution *
deltaG1=x_Fe*(4.309_pReal*Tp-2243.38_pReal)
deltaG2=x_Mn*(1.123_preal*Tp-1000.0_pReal)
deltaG3=x_Fe*x_Mn*(2873.0_pReal+717.0_pReal*(x_Fe-x_Mn))
deltaG4=1246.0_pReal*(1.0_pReal-exp(-24.29_pReal*x_C))-17175.0_pReal*x_Mn*x_C
!* Magnetical contribution *
beta_mart=0.62_pReal*x_Mn-4.0_pReal*x_C ! Magnetic spin in µB
Tp_mart=580.0_pReal*x_Mn ! Néel Temperature
beta_aust=0.7_pReal*x_Fe+0.62_pReal*x_Mn-0.64_pReal*x_Fe*x_Mn-4.0_pReal*x_C ! Magnetic spin in µB
Tp_aust=669.27_pReal*(1.0_pReal-exp(-5.46_pReal*x_Mn))-2408.0_pReal*x_C*x_Fe-109.0_pReal ! Néel Temperature
if (Tp<=Tp_mart) then
f_mart=1.0_pReal-(1.0_pReal/2.34_pReal)*((79.0_pReal*Tp_mart)/(140.0_pReal*0.28_pReal*Tp)+&
(474.0_pReal/497.0_pReal)*((1.0_pReal/0.28_pReal)-1.0_pReal)*(((Tp/Tp_mart)**3.0_pReal)/6.0_pReal+&
((Tp/Tp_mart)**9.0_pReal)/135.0_pReal+((Tp/Tp_mart)**15.0_pReal)/600.0_pReal))
else
f_mart=-(1.0_pReal/2.34_pReal)*(((Tp/Tp_mart)**-5.0_pReal)/10.0_pReal+((Tp/Tp_mart)**-15.0_pReal)/315.0_pReal+&
((Tp/Tp_mart)**-25.0_pReal)/1500.0_pReal)
endif
if (Tp<=Tp_aust) then
f_aust=1.0_pReal-(1.0_pReal/2.34_pReal)*((79.0_pReal*Tp_aust)/(140.0_pReal*0.28_pReal*Tp)+&
(474.0_pReal/497.0_pReal)*((1.0_pReal/0.28_pReal)-1.0_pReal)*(((Tp/Tp_aust)**3.0_pReal)/6.0_pReal+&
((Tp/Tp_aust)**9.0_pReal)/135.0_pReal+((Tp/Tp_aust)**15.0_pReal)/600.0_pReal))
else
f_aust=-(1.0_pReal/2.34_pReal)*(((Tp/Tp_aust)**-5.0_pReal)/10.0_pReal+((Tp/Tp_aust)**-15.0_pReal)/315.0_pReal+&
((Tp/Tp_aust)**-25.0_pReal)/1500.0_pReal)
endif
deltaG5=Rgaz*Tp*(log(beta_mart+1.0_pReal)*f_mart-log(beta_aust+1.0_pReal)*f_aust)
!* Final expression *
material_sfe=(4.0_pReal/(sqrt(3.0_pReal)*avogadro*material_bg(matID)**2.0_pReal))*&
(deltaG1+deltaG2+deltaG3+deltaG4+deltaG5)+2.0_pReal*0.0035_pReal
return return
end subroutine end subroutine
subroutine constitutive_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,state,Tp,ipc,ip,el) subroutine constitutive_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,state,Tp,ipc,ip,el)
!********************************************************************* !*********************************************************************
!* This subroutine contains the constitutive equation for * !* This subroutine contains the constitutive equation for *
@ -990,7 +1040,7 @@ subroutine constitutive_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,state,Tp,ipc,ip,el
!********************************************************************* !*********************************************************************
use prec, only: pReal,pInt use prec, only: pReal,pInt
use lattice, only: lattice_Sslip,lattice_Sslip_v,lattice_Stwin,lattice_Stwin_v,lattice_TwinShear use lattice, only: lattice_Sslip,lattice_Sslip_v,lattice_Stwin,lattice_Stwin_v,lattice_TwinShear
use math, only: math_Plain3333to99 use math, only: pi,math_Plain3333to99
implicit none implicit none
!* Definition of variables !* Definition of variables
@ -1002,92 +1052,155 @@ real(pReal), dimension(3,3) :: Lp,Sslip,Stwin
real(pReal), dimension(3,3,3,3) :: dLp_dTstar3333 real(pReal), dimension(3,3,3,3) :: dLp_dTstar3333
real(pReal), dimension(9,9) :: dLp_dTstar real(pReal), dimension(9,9) :: dLp_dTstar
real(pReal), dimension(constitutive_Nstatevars(ipc,ip,el)) :: state 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
real(pReal), dimension(material_Ntwin(constitutive_matID(ipc,ip,el))) :: sfe_eff,fdot_twin,dfdot_dtautwin,tau_twin,tauc_twin
real(pReal), dimension(material_Ntwin(constitutive_matID(ipc,ip,el)),material_Nslip(constitutive_matID(ipc,ip,el))) :: dfdot_dtauslip
!* Get the material-ID from the triplet(ipc,ip,el) !* Get the material-ID from the triplet(ipc,ip,el)
matID = constitutive_matID(ipc,ip,el) matID = constitutive_matID(ipc,ip,el)
startIdxTwin = material_Nslip(matID) startIdxTwin = material_Nslip(matID)
Ftwin = sum(state((startIdxTwin+1):(startIdxTwin+material_Ntwin(matID))))
!* Calculation of Lp - slip !* Calculation of Lp - slip
Ftwin = sum(state((startIdxTwin+1):(startIdxTwin+material_Ntwin(matID)))) gdot_slip = 0.0_pReal
Lp = 0.0_pReal dgdot_dtauslip = 0.0_pReal
Lp = 0.0_pReal
do i=1,material_Nslip(matID) do i=1,material_Nslip(matID)
constitutive_tau_slip(i)=dot_product(Tstar_v,lattice_Sslip_v(:,i,material_LatticeStructure(matID))) tau_slip(i)=dot_product(Tstar_v,lattice_Sslip_v(:,i,material_CrystalStructure(matID)))
if (abs(constitutive_tau_slip(i))>constitutive_passing_stress(i)) then if (abs(tau_slip(i))>constitutive_passing_stress(i)) then
constitutive_gdot_slip(i) = constitutive_g0_slip(i)*& gdot_slip(i) = constitutive_g0_slip(i)*sinh((tau_slip(i)*constitutive_activation_volume(i))/(kB*Tp))
sinh((constitutive_tau_slip(i)*constitutive_activation_volume(i))/(kB*Tp)) dgdot_dtauslip(i) = (constitutive_g0_slip(i)*constitutive_activation_volume(i))/(kB*Tp)*&
constitutive_dgdot_dtauslip(i) = (constitutive_g0_slip(i)*constitutive_activation_volume(i))/(kB*Tp)*& cosh((tau_slip(i)*constitutive_activation_volume(i))/(kB*Tp))
cosh((constitutive_tau_slip(i)*constitutive_activation_volume(i))/(kB*Tp)) endif
else Lp=Lp+(1.0_pReal-Ftwin)*gdot_slip(i)*lattice_Sslip(:,:,i,material_CrystalStructure(matID))
constitutive_gdot_slip(i) = 0.0_pReal
constitutive_dgdot_dtauslip(i) = 0.0_pReal
endif
Lp=Lp+(1.0_pReal-Ftwin)*constitutive_gdot_slip(i)*lattice_Sslip(:,:,i,material_LatticeStructure(matID))
enddo enddo
!write(6,*) '##############'
!write(6,*) '##############'
!write(6,*) 'Schmid_1', lattice_Sslip_v(:,1,material_CrystalStructure(matID))
!write(6,*) 'Schmid_2', lattice_Sslip_v(:,2,material_CrystalStructure(matID))
!write(6,*) 'Schmid_3', lattice_Sslip_v(:,3,material_CrystalStructure(matID))
!write(6,*) 'Schmid_4', lattice_Sslip_v(:,4,material_CrystalStructure(matID))
!write(6,*) 'Schmid_5', lattice_Sslip_v(:,5,material_CrystalStructure(matID))
!write(6,*) 'Schmid_6', lattice_Sslip_v(:,6,material_CrystalStructure(matID))
!write(6,*) 'Schmid_7', lattice_Sslip_v(:,7,material_CrystalStructure(matID))
!write(6,*) 'Schmid_8', lattice_Sslip_v(:,8,material_CrystalStructure(matID))
!write(6,*) 'Schmid_9', lattice_Sslip_v(:,9,material_CrystalStructure(matID))
!write(6,*) 'Schmid_10', lattice_Sslip_v(:,10,material_CrystalStructure(matID))
!write(6,*) 'Schmid_11', lattice_Sslip_v(:,11,material_CrystalStructure(matID))
!write(6,*) 'Schmid_12', lattice_Sslip_v(:,12,material_CrystalStructure(matID))
!write(6,*) 'Tstar_v',Tstar_v
!write(6,*) 'state',state
!write(6,*) 'Tp', Tp
!write(6,*) 'ssd_f', constitutive_rho_f
!write(6,*) 'ssd_p', constitutive_rho_p
!write(6,*) 'jump_width', constitutive_jump_width
!write(6,*) 'activation_volume', constitutive_activation_volume
!write(6,*) 'passing_stress', constitutive_passing_stress
!write(6,*) 'ssd_m', constitutive_rho_m
!write(6,*) 'g0_slip', constitutive_g0_slip
!write(6,*) 'tau_slip',tau_slip
!write(6,*) 'gdot_slip', gdot_slip
!write(6,*) 'dgdot_dtauslip', dgdot_dtauslip
!* Calculation of Lp - twin !* Calculation of Lp - twin
!do i=1,material_Ntwin(matID) sfe_eff = 0.0_pReal
! constitutive_tau_twin(i)=dot_product(Tstar_v,lattice_Stwin_v(:,i,material_LatticeStructure(matID))) fdot_twin = 0.0_pReal
! if (constitutive_tau_twin(i)>0.0_pReal) then dfdot_dtautwin = 0.0_pReal
! constitutive_fdot_twin(i) = (material_TwinSaturation(matID)-Ftwin)*constitutive_twin_volume(i)*& dfdot_dtauslip = 0.0_pReal
! material_c8(matID)*sum(state(1:material_Nslip(matID)))**(1.5_pReal)*& do i=1,material_Ntwin(matID)
! (material_ActivationLength(matID)/material_bg(matID))*sum(abs(constitutive_gdot_slip))*& tau_twin(i)=dot_product(Tstar_v,lattice_Stwin_v(:,i,material_CrystalStructure(matID)))
! exp(-((material_twin_res(matID)/constitutive_tau_twin(i))**material_c9(matID))) if ((tau_twin(i) > 0.0_pReal).AND.(material_TwinSaturation(matID)-Ftwin>=0)) then
! constitutive_dfdot_dtautwin(i) = (material_TwinSaturation(matID)-Ftwin)*constitutive_twin_volume(i)*& sfe_eff(i)=material_sfe-(sqrt(3.0_pReal)/3.0_pReal)*material_q1(matID)*material_q2(matID)*material_bg(matID)*tau_twin(i)
! material_c8(matID)*sum(state(1:material_Nslip(matID)))**(1.5_pReal)*& if (sfe_eff(i)<0.0_pReal) sfe_eff(i) = 0.0_pReal
! (material_ActivationLength(matID)/material_bg(matID))*sum(abs(constitutive_gdot_slip))*& fdot_twin(i) = (material_TwinSaturation(matID)-Ftwin)*&
! (material_c9(matID)/constitutive_tau_twin(i))*& constitutive_twin_volume(i)*&
! (material_twin_res(matID)/constitutive_tau_twin(i))**material_c9(matID)*& ((2.0_pReal*sqrt(6.0_pReal)*material_SiteScaling(matID)*sum(abs(gdot_slip))*&
! exp(-((material_twin_res(matID)/constitutive_tau_twin(i))**material_c9(matID))) sum(state(1:material_Nslip(matID)))**1.5_pReal)/3.0_pReal)*&
! do j=1,material_Nslip(matID) exp((-25.0_pReal*pi**3.0_pReal*material_Gmod(matID)**2.0_pReal*sfe_eff(i)**3.0_pReal)/&
! if (constitutive_gdot_slip(i)>0.0_pReal) then (3.0_pReal*Kb*Tp*(material_q2(matID)*tau_twin(i))**4.0_pReal))
! constitutive_dfdot_dtauslip(i,j) = (material_TwinSaturation(matID)-Ftwin)*constitutive_twin_volume(i)*&
! material_c8(matID)*sum(state(1:material_Nslip(matID)))**(1.5_pReal)*& dfdot_dtautwin(i) = (material_TwinSaturation(matID)-Ftwin)*&
! (material_ActivationLength(matID)/material_bg(matID))*constitutive_dgdot_dtauslip(j)*& constitutive_twin_volume(i)*&
! exp(-((material_twin_res(matID)/constitutive_tau_twin(i))**material_c9(matID))) ((2.0_pReal*sqrt(6.0_pReal)*material_SiteScaling(matID)*sum(abs(gdot_slip))*&
! else sum(state(1:material_Nslip(matID)))**1.5_pReal)/3.0_pReal)*&
! constitutive_dfdot_dtauslip(i,j) = (material_TwinSaturation(matID)-Ftwin)*constitutive_twin_volume(i)*& ((-25.0_pReal*pi**3.0_pReal*material_Gmod(matID)**2.0_pReal*sfe_eff(i)**2.0_pReal)/&
! material_c8(matID)*sum(state(1:material_Nslip(matID)))**(1.5_pReal)*& (3.0_pReal*Kb*Tp*material_q2(matID)**4.0_pReal*tau_twin(i)**5.0_pReal))*&
! (material_ActivationLength(matID)/material_bg(matID))*(-constitutive_dgdot_dtauslip(j))*& (-sqrt(3.0_pReal)*material_q1(matID)*material_q2(matID)*material_bg(matID)*tau_twin(i)-&
! exp(-((material_twin_res(matID)/constitutive_tau_twin(i))**material_c9(matID))) 4.0_pReal*sfe_eff(i))*&
! endif exp((-25.0_pReal*pi**3.0_pReal*material_Gmod(matID)**2.0_pReal*sfe_eff(i)**3.0_pReal)/&
! enddo (3.0_pReal*Kb*Tp*(material_q2(matID)*tau_twin(i))**4.0_pReal))
! else
! constitutive_fdot_twin(i)=0.0_pReal do j=1,material_Nslip(matID)
! constitutive_dfdot_dtautwin(i)=0.0_pReal if (gdot_slip(j)>0.0_pReal) then
! do j=1,material_Nslip(matID) dfdot_dtauslip(i,j) = (material_TwinSaturation(matID)-Ftwin)*&
! constitutive_dfdot_dtauslip(i,j)=0.0_pReal constitutive_twin_volume(i)*&
! enddo ((2.0_pReal*sqrt(6.0_pReal)*material_SiteScaling(matID)*dgdot_dtauslip(j)*&
! endif sum(state(1:material_Nslip(matID)))**1.5_pReal)/3.0_pReal)*&
! Lp=Lp+state(material_Nslip(matID)+i)*lattice_TwinShear(material_LatticeStructure(matID))*constitutive_fdot_twin(i)*& exp((-25.0_pReal*pi**3.0_pReal*material_Gmod(matID)**2.0_pReal*sfe_eff(i)**3.0_pReal)/&
! lattice_Stwin(:,:,i,material_LatticeStructure(matID)) (3.0_pReal*Kb*Tp*(material_q2(matID)*tau_twin(i))**4.0_pReal))
else
dfdot_dtauslip(i,j) = (material_TwinSaturation(matID)-Ftwin)*&
constitutive_twin_volume(i)*&
((2.0_pReal*sqrt(6.0_pReal)*material_SiteScaling(matID)*(-dgdot_dtauslip(j))*&
sum(state(1:material_Nslip(matID)))**1.5_pReal)/3.0_pReal)*&
exp((-25.0_pReal*pi**3.0_pReal*material_Gmod(matID)**2.0_pReal*sfe_eff(i)**3.0_pReal)/&
(3.0_pReal*Kb*Tp*(material_q2(matID)*tau_twin(i))**4.0_pReal))
endif
enddo
endif
Lp=Lp+state(material_Nslip(matID)+i)*lattice_TwinShear(material_CrystalStructure(matID))*constitutive_fdot_twin(i)*&
lattice_Stwin(:,:,i,material_CrystalStructure(matID))
enddo
!write(6,*) 'twin_mfp', constitutive_twin_mfp
!write(6,*) 'twin_volume', constitutive_twin_volume
!write(6,*) 'tau_twin',tau_twin
!write(6,*) 'part1:',material_TwinSaturation(matID)-Ftwin
!write(6,*) 'part2:',constitutive_twin_volume
!write(6,*) 'part3:',((2.0_pReal*sqrt(6.0_pReal)*sum(abs(gdot_slip))*&
! sum(state(1:material_Nslip(matID)))**1.5_pReal)/3.0_pReal)
!do i=1,12
!write(6,*) 'part4:',exp((-25.0_pReal*pi**3.0_pReal*material_Gmod(matID)**2.0_pReal*sfe_eff(i)**3.0_pReal)/&
! (3.0_pReal*Kb*Tp*(material_q2(matID)*tau_twin(i))**4.0_pReal))
!write(6,*) 'part5:',(-25.0_pReal*pi**3.0_pReal*material_Gmod(matID)**2.0_pReal*sfe_eff(i)**3.0_pReal)/&
! (3.0_pReal*Kb*Tp*(material_q2(matID)*tau_twin(i))**4.0_pReal)
!enddo !enddo
!write(6,*) 'sfe', material_sfe
!write(6,*) 'sfe_eff',sfe_eff
!write(6,*) 'fdot_twin', fdot_twin
!write(6,*) 'dfdot_dtautwin', dfdot_dtautwin
!write(6,*) 'dfdot_dtauslip', dfdot_dtauslip
!* Calculation of the tangent of Lp !* Calculation of the tangent of Lp
dLp_dTstar3333=0.0_pReal dLp_dTstar3333=0.0_pReal
do i=1,material_Nslip(matID) do i=1,material_Nslip(matID)
Sslip = lattice_Sslip(:,:,i,material_LatticeStructure(matID)) Sslip = lattice_Sslip(:,:,i,material_CrystalStructure(matID))
forall (k=1:3,l=1:3,m=1:3,n=1:3) forall (k=1:3,l=1:3,m=1:3,n=1:3)
dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n)+ & dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n)+ &
(1.0_pReal-Ftwin)*constitutive_dgdot_dtauslip(i)*Sslip(k,l)*Sslip(m,n) !force m,n symmetry (1.0_pReal-Ftwin)*dgdot_dtauslip(i)*Sslip(k,l)*Sslip(m,n) !force m,n symmetry
endforall endforall
enddo enddo
!do i=1,material_Ntwin(matID) do i=1,material_Ntwin(matID)
! Stwin = lattice_Stwin(:,:,i,material_LatticeStructure(matID)) Stwin = lattice_Stwin(:,:,i,material_CrystalStructure(matID))
! forall (k=1:3,l=1:3,m=1:3,n=1:3) forall (k=1:3,l=1:3,m=1:3,n=1:3)
! dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n)+ & dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n)+ &
! state(material_Nslip(matID)+i)*lattice_TwinShear(material_LatticeStructure(matID))*& state(material_Nslip(matID)+i)*lattice_TwinShear(material_CrystalStructure(matID))*&
! constitutive_dfdot_dtautwin(i)*Stwin(k,l)*Stwin(m,n) !force m,n symmetry dfdot_dtautwin(i)*Stwin(k,l)*Stwin(m,n) !force m,n symmetry
! endforall endforall
! do j=1,material_Nslip(matID) do j=1,material_Nslip(matID)
! Sslip = lattice_Sslip(:,:,j,material_LatticeStructure(matID)) Sslip = lattice_Sslip(:,:,j,material_CrystalStructure(matID))
! forall (k=1:3,l=1:3,m=1:3,n=1:3) forall (k=1:3,l=1:3,m=1:3,n=1:3)
! dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n)+ & dLp_dTstar3333(k,l,m,n) = dLp_dTstar3333(k,l,m,n)+ &
! state(material_Nslip(matID)+i)*lattice_TwinShear(material_LatticeStructure(matID))*& state(material_Nslip(matID)+i)*lattice_TwinShear(material_CrystalStructure(matID))*&
! constitutive_dfdot_dtauslip(i,j)*Stwin(k,l)*Sslip(m,n) !force m,n symmetry dfdot_dtauslip(i,j)*Stwin(k,l)*Sslip(m,n) !force m,n symmetry
! endforall endforall
! enddo enddo
!enddo enddo
dLp_dTstar = math_Plain3333to99(dLp_dTstar3333) dLp_dTstar = math_Plain3333to99(dLp_dTstar3333)
return return
@ -1109,6 +1222,7 @@ function constitutive_dotState(Tstar_v,state,Tp,ipc,ip,el)
!* - constitutive_DotState : evolution of state variable * !* - constitutive_DotState : evolution of state variable *
!********************************************************************* !*********************************************************************
use prec, only: pReal,pInt use prec, only: pReal,pInt
use math, only: pi
use lattice, only: lattice_Sslip_v,lattice_Stwin_v use lattice, only: lattice_Sslip_v,lattice_Stwin_v
implicit none implicit none
@ -1118,51 +1232,48 @@ integer(pInt) matID,i,j,startIdxTwin
real(pReal) Tp,Ftwin real(pReal) Tp,Ftwin
real(pReal), dimension(6) :: Tstar_v real(pReal), dimension(6) :: Tstar_v
real(pReal), dimension(constitutive_Nstatevars(ipc,ip,el)) :: constitutive_dotState,state real(pReal), dimension(constitutive_Nstatevars(ipc,ip,el)) :: constitutive_dotState,state
real(pReal), dimension(material_Nslip(constitutive_matID(ipc,ip,el))) :: gdot_slip,tau_slip
real(pReal), dimension(material_Ntwin(constitutive_matID(ipc,ip,el))) :: sfe_eff,fdot_twin,tau_twin,tauc_twin
real(pReal), dimension(material_Nslip(constitutive_matID(ipc,ip,el))) :: locks,grainboundaries,recovery
real(pReal), dimension(material_Nslip(constitutive_matID(ipc,ip,el))) :: twinboundaries
!* Get the material-ID from the triplet(ipc,ip,el) !* Get the material-ID from the triplet(ipc,ip,el)
matID = constitutive_matID(ipc,ip,el) matID = constitutive_matID(ipc,ip,el)
startIdxTwin = material_Nslip(matID) startIdxTwin = material_Nslip(matID)
Ftwin = sum(state((startIdxTwin+1):(startIdxTwin+material_Ntwin(matID))))
constitutive_dotState = 0.0_pReal constitutive_dotState = 0.0_pReal
!* Dislocation density evolution !* Dislocation density evolution
gdot_slip = 0.0_pReal
do i=1,material_Nslip(matID) do i=1,material_Nslip(matID)
constitutive_tau_slip(i)=dot_product(Tstar_v,lattice_Sslip_v(:,i,material_LatticeStructure(matID))) tau_slip(i)=dot_product(Tstar_v,lattice_Sslip_v(:,i,material_CrystalStructure(matID)))
if (abs(constitutive_tau_slip(i))>constitutive_passing_stress(i)) then if (abs(tau_slip(i))>constitutive_passing_stress(i)) &
constitutive_gdot_slip(i) = constitutive_g0_slip(i)*& gdot_slip(i) = constitutive_g0_slip(i)*sinh((tau_slip(i)*constitutive_activation_volume(i))/(kB*Tp))
sinh((constitutive_tau_slip(i)*constitutive_activation_volume(i))/(kB*Tp))
else locks(i) = (sqrt(constitutive_rho_f(i))*abs(gdot_slip(i)))/(material_c4(matID)*material_bg(matID))
constitutive_gdot_slip(i) = 0.0_pReal grainboundaries(i) = abs(gdot_slip(i))/(material_c5(matID)*material_bg(matID)*material_GrainSize(matID))
endif twinboundaries(i) = (abs(gdot_slip(i))*constitutive_inv_intertwin_len_s(i))/(material_c6(matID)*material_bg(matID))
constitutive_locks(i) = (sqrt(constitutive_rho_f(i))*abs(constitutive_gdot_slip(i)))/& recovery(i) = material_c7(matID)*state(i)*abs(gdot_slip(i))
(material_c4(matID)*material_bg(matID))
constitutive_dotState(i) = locks(i)+grainboundaries(i)+twinboundaries(i)-recovery(i)
constitutive_grainboundaries(i) = abs(constitutive_gdot_slip(i))/(material_c5(matID)*material_bg(matID)*material_GrainSize(matID))
! if (material_Ntwin(matID)>0) then
! constitutive_twinboundaries(i) = (abs(constitutive_gdot_slip(i))*constitutive_inv_intertwin_len(i))/&
! (material_c6(matID)*material_bg(matID))
! endif
constitutive_recovery(i) = material_c7(matID)*state(i)*abs(constitutive_gdot_slip(i))
constitutive_dotState(i) = constitutive_locks(i)+constitutive_grainboundaries(i)-constitutive_recovery(i)
! constitutive_dotState(i) = constitutive_locks(i)+constitutive_grainboundaries(i)+constitutive_twinboundaries(i)&
! -constitutive_recovery(i)
enddo enddo
!* Twin volume fraction evolution !* Twin volume fraction evolution
!Ftwin = sum(state((startIdxTwin+1):(startIdxTwin+material_Ntwin(matID)))) fdot_twin = 0.0_pReal
!do i=1,material_Ntwin(matID) do i=1,material_Ntwin(matID)
! constitutive_tau_twin(i)=dot_product(Tstar_v,lattice_Stwin_v(:,i,material_LatticeStructure(matID))) tau_twin(i)=dot_product(Tstar_v,lattice_Stwin_v(:,i,material_CrystalStructure(matID)))
! if (constitutive_tau_twin(i)>0.0_pReal) then if (tau_twin(i)>0.0_pReal) then
! constitutive_fdot_twin(i) = (material_TwinSaturation(matID)-Ftwin)*constitutive_twin_volume(i)*& sfe_eff(i)=material_sfe-(sqrt(3.0_pReal)/3.0_pReal)*material_q1(matID)*material_q2(matID)*material_bg(matID)*tau_twin(i)
! material_c8(matID)*sum(state(1:material_Nslip(matID)))**(1.5_pReal)*& if (sfe_eff(i)<0.0_pReal) sfe_eff(i) = 0.0_pReal
! (material_ActivationLength(matID)/material_bg(matID))*sum(abs(constitutive_gdot_slip))*& fdot_twin(i) = (material_TwinSaturation(matID)-Ftwin)*&
! exp(-((material_twin_res(matID)/constitutive_tau_twin(i))**material_c9(matID))) constitutive_twin_volume(i)*&
! else ((2.0_pReal*sqrt(6.0_pReal)*material_SiteScaling(matID)*sum(abs(gdot_slip))*&
! constitutive_fdot_twin(i) = 0.0_pReal sum(state(1:material_Nslip(matID)))**1.5_pReal)/3.0_pReal)*&
! endif exp((-25.0_pReal*pi**3.0_pReal*material_Gmod(matID)**2.0_pReal*sfe_eff(i)**3.0_pReal)/&
! constitutive_dotState(startIdxTwin+i)=constitutive_fdot_twin(i) (3.0_pReal*Kb*Tp*(material_q2(matID)*tau_twin(i))**4.0_pReal))
!enddo endif
constitutive_dotState(startIdxTwin+i) = fdot_twin(i)
!constitutive_dotState=0.0_pReal enddo
return return
end function end function
@ -1200,12 +1311,14 @@ startIdxTwin = material_Nslip(matID)
if(constitutive_Nresults(ipc,ip,el)==0) return if(constitutive_Nresults(ipc,ip,el)==0) return
constitutive_post_results=0.0_pReal constitutive_post_results=0.0_pReal
do i=1,material_Nslip(matID) !do i=1,material_Nslip(matID)
constitutive_post_results(i) = state(i) ! constitutive_post_results(i) = state(i)
enddo !enddo
do i=1,material_Ntwin(matID) !do i=1,material_Ntwin(matID)
constitutive_post_results(startIdxTwin+i) = state(startIdxTwin+i) ! constitutive_post_results(startIdxTwin+i) = state(startIdxTwin+i)
enddo !enddo
constitutive_post_results(1) = sum(state(1:material_Nslip(matID)))
constitutive_post_results(2) = sum(state((startIdxTwin+1):(startIdxTwin+material_Ntwin(matID))))
return return
end function end function

View File

@ -176,6 +176,8 @@ CONTAINS
enddo enddo
enddo enddo
endif endif
!
msg = 'ok' ! a new consistent tangent was computed even if msg was not ok for all components
! !
return return
! !
@ -418,7 +420,13 @@ Inner: do ! inner iteration: Lp
dt*constitutive_dotState(Tstar_v,state,Temperature,& dt*constitutive_dotState(Tstar_v,state,Temperature,&
grain,ip,cp_en) ! residuum from evolution of microstructure grain,ip,cp_en) ! residuum from evolution of microstructure
!!$OMP END CRITICAL (stateupdate) !!$OMP END CRITICAL (stateupdate)
state = state - ROuter ! update of microstructure state = state - ROuter ! update of microstructure
if (iOuter==nOuter) then
!$OMP CRITICAL (write2out)
write (6,*) 'WARNING: Outer loop has not really converged'
!$OMP END CRITICAL (write2out)
exit Outer
endif
if (maxval(abs(Router/state),state /= 0.0_pReal) < reltol_Outer) exit Outer if (maxval(abs(Router/state),state /= 0.0_pReal) < reltol_Outer) exit Outer
enddo Outer enddo Outer
! !

View File

@ -21,52 +21,57 @@ w0 1.0
hardening_coefficients 1.0 1.4 hardening_coefficients 1.0 1.4
## Parameters for dislocation-based modeling ## Parameters for dislocation-based modeling
# Initial dislocation density [m]²
rho0 2.8e13
# Burgers vector [m] # Burgers vector [m]
burgers 2.56e-10 burgers 2.56e-10
# Activation energy for dislocation glide [J/K] # Activation energy for dislocation glide [J/K]
Qedge 3.0e-19 Qedge 3.0e-19
# Reference for passing stress [Pa] # Initial dislocation density [m]²
tau0 0.0 rho0 2.8e13
# Average grain size [m]
grain_size 2.0e-5
# Passing stress adjustment # Passing stress adjustment
c1 0.1 c1 0.1
# Jump width adjustment # Jump width adjustment
c2 2.0 c2 2.0
# Activation volume adjustment # Activation volume adjustment
c3 1.2 c3 1.2
# Dislocation storage adjustment # Average slip distance adjustment for lock formation
# = c4(Anxin)*c2(Anxin) !!!!!! # = c4(Anxin)*c2(Anxin) !!!!!!
c4 14.25 c4 14.25
# Grain boundaries storage adjustment # Average slip distance adjustment when grain boundaries
c5 1.0 c5 1.0
# Athermal annihilation adjustment # Average slip distance adjustment when twin boundaries
c6 0.1
# Athermal recovery adjustment
c7 23.5 c7 23.5
# Dislocation interaction coefficients # Dislocation interaction coefficients
interaction_coefficients 1.0 2.2 3.0 1.6 3.8 4.5 interaction_coefficients 1.0 2.2 3.0 1.6 3.8 4.5
# Twin parameters ## Parameters for mechanical twinning
# Grain size [m] # Average twin thickness (stacks) [m]
grain_size 2.0e-5
# Twin thickness (stacks) [m]
stack_size 5.0e-8 stack_size 5.0e-8
# Activation length for twin nucleation [m] # Total twin volume fraction saturation
d_star 5.0e-10 f_sat 0.2
# Twin saturation value # Scaling potential nucleation sites
f_sat 0.3 site_scaling 1.0e-7
# Twin boundaries storage adjustment # Scaling the P-K force on the twinning dislocation
c6 0.425 q1 0.75
# Scaling of really activated nucleation sites # Scaling the resolved shear stress
c8 2.0e-3 q2 1.0
# Selection of active twin systems
c9 10.0
# Twin resistance [Pa]
twin_resistance 1000.0e6
stacking_fault_energy 2.0e-2
<textures> <textures>
[cube SX] [cube SX]
symmetry no /monoclinic /orthorhombic symmetry no /monoclinic /orthorhombic
Ngrains 1 /2 /4 Ngrains 10 /2 /4
(gauss) phi1 0.0 phi 0.0 phi2 0.0 scatter 0.0 fraction 1.0 #(gauss) phi1 0.0 phi 29.21 phi2 -26.57 scatter 0.0 fraction 1.0
#(gauss) phi1 0.0 phi 54.74 phi2 -45.0 scatter 0.0 fraction 0.1
#(gauss) phi1 0.0 phi 45.0 phi2 0.0 scatter 0.0 fraction 0.1
#(gauss) phi1 0.0 phi 0.0 phi2 0.0 scatter 0.0 fraction 0.1
#(gauss) phi1 0.0 phi 35.26 phi2 -45.0 scatter 0.0 fraction 0.1
#(gauss) phi1 0.0 phi 48.19 phi2 -26.57 scatter 0.0 fraction 0.1
#(gauss) phi1 0.0 phi 26.57 phi2 0.0 scatter 0.0 fraction 0.1
#(gauss) phi1 0.0 phi 42.03 phi2 -33.69 scatter 0.0 fraction 0.1
#(gauss) phi1 0.0 phi 40.36 phi2 -11.31 scatter 0.0 fraction 0.1
#(gauss) phi1 0.0 phi 15.62 phi2 -26.57 scatter 0.0 fraction 0.1