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

mitigate rounding errors

This commit is contained in:
Martin Diehl 2022-05-20 06:43:32 +02:00
parent 491e2ec0b2
commit f42bb5d175
2 changed files with 2 additions and 2 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
src/phase_mechanical.f90
View File

@ -264,8 +264,6 @@ module subroutine mechanical_init(phases)
en = material_phaseEntry(co,ce)
phase_mechanical_F(ph)%data(1:3,1:3,en) = math_I3
phase_mechanical_Fp(ph)%data(1:3,1:3,en) = material_O_0(ma)%data(co)%asMatrix() ! Fp reflects initial orientation (see 10.1016/j.actamat.2006.01.005)
phase_mechanical_Fp(ph)%data(1:3,1:3,en) = phase_mechanical_Fp(ph)%data(1:3,1:3,en) &
/ math_det33(phase_mechanical_Fp(ph)%data(1:3,1:3,en))**(1.0_pReal/3.0_pReal)
phase_mechanical_Fe(ph)%data(1:3,1:3,en) = matmul(material_V_e_0(ma)%data(1:3,1:3,co), &
transpose(phase_mechanical_Fp(ph)%data(1:3,1:3,en)))
phase_mechanical_Fi(ph)%data(1:3,1:3,en) = material_O_0(ma)%data(co)%rotate(math_inv33(material_V_e_0(ma)%data(1:3,1:3,co)))

2
src/rotations.f90
View File

@ -449,6 +449,7 @@ pure function qu2om(qu) result(om)
om(1,3) = 2.0_pReal*(qu(2)*qu(4)+qu(1)*qu(3))
if (sign(1.0_pReal,P) < 0.0_pReal) om = transpose(om)
om = om/math_det33(om)**(1.0_pReal/3.0_pReal)
end function qu2om
@ -602,6 +603,7 @@ pure function om2qu(om) result(qu)
endif
if(sign(1.0_pReal,qu(1))<0.0_pReal) qu =-1.0_pReal * qu
qu(2:4) = merge(qu(2:4),qu(2:4)*P,dEq0(qu(2:4)))
qu = qu/norm2(qu)
end function om2qu