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

sitting together and finding open issues.. ;-)

This commit is contained in:
Philip Eisenlohr 2007-03-26 15:32:58 +00:00
parent 86d7f222ae
commit 5752fe2ca7
1 changed files with 19 additions and 30 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

41
trunk/CPFEM.f90
View File

@ -504,13 +504,12 @@
integer(pInt), parameter :: ninner = 2000_pInt integer(pInt), parameter :: ninner = 2000_pInt
real(pReal), parameter :: tol_inner = 1.0e-3_pReal real(pReal), parameter :: tol_inner = 1.0e-3_pReal
real(pReal), parameter :: eta = 13.7_pReal real(pReal), parameter :: eta = 13.7_pReal
integer(pInt), parameter :: numerical = 0_pInt
real(pReal), parameter :: pert_nr = 1.0e-8_pReal
crite=eta*constitutive_s0_slip/constitutive_n_slip !ÄÄÄ crite=eta*constitutive_s0_slip/constitutive_n_slip !ÄÄÄ
! !
! *** Tolerances *** ! *** Tolerances ***
tol_in = tol_inner*s0_slip tol_in = tol_inner*s0_slip !ÄÄÄ
tol_out = tol_outer*s0_slip tol_out = tol_outer*s0_slip !ÄÄÄ
! !
! *** Error treatment *** ! *** Error treatment ***
iconv = 0 iconv = 0
@ -546,30 +545,26 @@
! !
! *** Jacobi Calculation *** ! *** Jacobi Calculation ***
help=matmul(A, I3tLp) help=matmul(A, I3tLp)
! MISSING 1..6 outer loop, inner sum required
do i=1,3 do i=1,3
do j=1,3 do j=1,3
do k=1,3 do k=1,6
! dol=1,3
help1(k)=dLp(j,i,k)*help(i,j) help1(k)=dLp(j,i,k)*help(i,j)
! enddo
enddo enddo
enddo enddo
enddo enddo
! help=help1+transpose(help1) ! help=help1+transpose(help1)
Jacobi= matmul(C66, help1) + mat_identity(6) Jacobi= matmul(C66, help1) + math_identity(6)
call math_invert6x6(Jacobi, invJacobi, dummy, err) !ÄÄÄ call math_invert6x6(Jacobi, invJacobi, dummy, err) !ÄÄÄ
if (err==1_pInt) then if (err==1_pInt) then
do i=1,6 forall (i=1:6) Jacobi(i,i)=1.05d0*maxval(Jacobi(i,:)) ! regularization
Jacobi(i,i)=1.05d0*maxval(Jacobi(i,:))
enddo
invJacobi=Jacobi
call math_invert6x6(Jacobi, invJacobi, dummy, err) call math_invert6x6(Jacobi, invJacobi, dummy, err)
if (err==1_pInt) then if (err==1_pInt) then ! sorry, can't help here!!
ising=1 ising=1
return return
endif endif
endif endif
dTstar_v=matmul(invJacobi,R1) dTstar_v=matmul(invJacobi,R1) ! correction to Tstar
! *** Correction (see Kalidindi) *** ! *** Correction (see Kalidindi) ***
do i=1,6 do i=1,6
@ -584,26 +579,19 @@
return return
! *** End of the first level of iterative procedure *** ! *** End of the first level of iterative procedure ***
100 continue 100 dstate=dt*constitutive_dotState(Tstar_v, iori, CPFEM_in, cp_en)
! call hardening(tauc_slip_new,gdot_slip,dtauc_slip)
dstate=constitutive_dotState(Tstar_v, iori, CPFEM_in, cp_en)
! *** Arrays of residuals *** ! *** Arrays of residuals ***
R2=state_new-state_old-dt*dstate R2=state_new-state_old-dstate
norm2=maxval(abs(R2)) norm2=maxval(abs(R2))
if (norm2<tol_out) goto 200 if (norm2<tol_out) goto 200
state_new=state_old+dt*dstate state_new=state_old+dstate
enddo enddo
iconv=2 iconv=2
return return
! *** End of the second level of iterative procedure *** ! *** End of the second level of iterative procedure ***
200 continue
!
! call plastic_vel_grad(dt,tau_slip,tauc_slip_new,Lp)
!
! *** Calculation of Fp(t+dt) (see Kalidindi) *** ! *** Calculation of Fp(t+dt) (see Kalidindi) ***
invFp_new=matmul(Fp_old, I3tLp) 200 invFp_new=matmul(Fp_old, I3tLp)
call math_invert3x3(invFp_new, Fp_new, det, err) !ÄÄÄ call math_invert3x3(invFp_new, Fp_new, det, err) !ÄÄÄ
if (err==1_pInt) then if (err==1_pInt) then
ising=1 ising=1
@ -619,7 +607,8 @@
call math_conv33to6(Estar,Estar_v) call math_conv33to6(Estar,Estar_v)
! !
! *** Calculation of the Cauchy stress *** ! *** Calculation of the Cauchy stress ***
call CPFEM_cauchy_stress(Estar_v,Fe,cs) ! QUESTION seems to need Tstar, not Estar..??
cs = CPFEM_cauchy_stress(Estar_v,Fe)
! !
return return
end subroutine end subroutine