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

merge evaluates both expression, can lead to division by zero

This commit is contained in:
Martin Diehl 2021-04-13 21:30:56 +02:00
parent d47c116ebd
commit 285889b48e
2 changed files with 14 additions and 12 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

14
src/phase_damage.f90
View File

@ -257,20 +257,20 @@ module function integrateDamageState(dt,co,ip,el) result(broken)
contains contains
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief calculate the damping for correction of state and dot state !> @brief calculate the damping for correction of state and dot state
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
real(pReal) pure function damper(current,previous,previous2) real(pReal) pure function damper(omega_0,omega_1,omega_2)
real(pReal), dimension(:), intent(in) ::& real(pReal), dimension(:), intent(in) :: &
current, previous, previous2 omega_0, omega_1, omega_2
real(pReal) :: dot_prod12, dot_prod22 real(pReal) :: dot_prod12, dot_prod22
dot_prod12 = dot_product(current - previous, previous - previous2) dot_prod12 = dot_product(omega_0-omega_1, omega_1-omega_2)
dot_prod22 = dot_product(previous - previous2, previous - previous2) dot_prod22 = dot_product(omega_1-omega_2, omega_1-omega_2)
if ((dot_product(current,previous) < 0.0_pReal .or. dot_prod12 < 0.0_pReal) .and. dot_prod22 > 0.0_pReal) then
if (min(dot_product(omega_0,omega_1),dot_prod12) < 0.0_pReal .and. dot_prod22 > 0.0_pReal) then
damper = 0.75_pReal + 0.25_pReal * tanh(2.0_pReal + 4.0_pReal * dot_prod12 / dot_prod22) damper = 0.75_pReal + 0.25_pReal * tanh(2.0_pReal + 4.0_pReal * dot_prod12 / dot_prod22)
else else
damper = 1.0_pReal damper = 1.0_pReal

12
src/phase_mechanical.f90
View File

@ -720,12 +720,14 @@ function integrateStateFPI(F_0,F,subFp0,subFi0,subState0,Delta_t,co,ip,el) resul
real(pReal) :: dot_prod12, dot_prod22 real(pReal) :: dot_prod12, dot_prod22
dot_prod12 = dot_product(omega_0 - omega_1, omega_1 - omega_2) dot_prod12 = dot_product(omega_0-omega_1, omega_1-omega_2)
dot_prod22 = dot_product(omega_1 - omega_2, omega_1 - omega_2) dot_prod22 = dot_product(omega_1-omega_2, omega_1-omega_2)
damper = merge(0.75_pReal + 0.25_pReal * tanh(2.0_pReal + 4.0_pReal * dot_prod12 / dot_prod22), & if (min(dot_product(omega_0,omega_1),dot_prod12) < 0.0_pReal .and. dot_prod22 > 0.0_pReal) then
1.0_pReal, & damper = 0.75_pReal + 0.25_pReal * tanh(2.0_pReal + 4.0_pReal * dot_prod12 / dot_prod22)
(min(dot_product(omega_0,omega_1), dot_prod12) < 0.0_pReal) .and. dot_prod22 > 0.0_pReal) else
damper = 1.0_pReal
endif
end function damper end function damper