fixed plastic work calculation

This commit is contained in:
Fengbo Han 2017-08-10 15:40:18 +02:00
parent 9bbc0d4803
commit 0750f7fd01
2 changed files with 13 additions and 4 deletions

View File

@ -317,6 +317,14 @@ program DAMASK_spectral
enddo; enddo enddo; enddo
close(FILEUNIT) close(FILEUNIT)
if(yieldStop) then ! initialize variables related to yield stop
yieldStressNew = 0.0_pReal
plasticStrainNew = 0.0_pReal
eqStressNew = 0.0_pReal
eqTotalStrainNew = 0.0_pReal
eqPlasticStrainNew = 0.0_pReal
plasticWorkNew = 0.0_pReal
endif
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! consistency checks and output of load case ! consistency checks and output of load case
loadCases(1)%followFormerTrajectory = .false. ! cannot guess along trajectory for first inc of first currentLoadCase loadCases(1)%followFormerTrajectory = .false. ! cannot guess along trajectory for first inc of first currentLoadCase

View File

@ -1080,9 +1080,12 @@ subroutine utilities_calcPlasticity(yieldStress, plasticStrain, eqStress, eqTota
real(pReal), dimension(15) :: WORK !< previous deformation gradient real(pReal), dimension(15) :: WORK !< previous deformation gradient
integer(pInt) :: INFO, i, j, k, l, ierr integer(pInt) :: INFO, i, j, k, l, ierr
real(pReal) :: stress, stress_av, strain_total, strain_total_av, strain_plastic, strain_plastic_av, wgtm real(pReal) :: stress, stress_av, strain_total, strain_total_av, strain_plastic, strain_plastic_av, wgtm
real(pReal) :: eqStressOld, eqPlasticStrainOld
external :: dgesvd external :: dgesvd
eqStressOld = eqStress
eqPlasticStrainOld = eqPlasticStrain
wgtm = 1.0/real(mesh_NcpElems*mesh_maxNips*homogenization_maxNgrains,pReal) wgtm = 1.0/real(mesh_NcpElems*mesh_maxNips*homogenization_maxNgrains,pReal)
Vp_av = 0.0_pReal Vp_av = 0.0_pReal
V_total_av = 0.0_pReal V_total_av = 0.0_pReal
@ -1142,20 +1145,18 @@ subroutine utilities_calcPlasticity(yieldStress, plasticStrain, eqStress, eqTota
enddo; enddo; enddo enddo; enddo; enddo
yieldStress = cauchy_av * wgtm yieldStress = cauchy_av * wgtm
call MPI_Allreduce(MPI_IN_PLACE,yieldStress,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,yieldStress,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
plasticStrain = Vp_av * wgtm plasticStrain = Vp_av * wgtm
call MPI_Allreduce(MPI_IN_PLACE,plasticStrain,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,plasticStrain,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
plasticWork = plasticWork + 0.5*(eqStress + stress_av * wgtm) * (strain_total_av * wgtm - eqTotalStrain)
call MPI_Allreduce(MPI_IN_PLACE,plasticWork,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
eqStress = stress_av * wgtm eqStress = stress_av * wgtm
call MPI_Allreduce(MPI_IN_PLACE,eqStress,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,eqStress,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
eqTotalStrain = strain_total_av * wgtm eqTotalStrain = strain_total_av * wgtm
call MPI_Allreduce(MPI_IN_PLACE,eqTotalStrain,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,eqTotalStrain,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
eqPlasticStrain = strain_plastic_av * wgtm eqPlasticStrain = strain_plastic_av * wgtm
call MPI_Allreduce(MPI_IN_PLACE,eqPlasticStrain,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,eqPlasticStrain,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
plasticWork = plasticWork + 0.5*(eqStressOld + eqStress) * (eqPlasticStrain - eqPlasticStrainOld)
end subroutine utilities_calcPlasticity end subroutine utilities_calcPlasticity