From 491e2ec0b2b963353e09347c38383c78710477fa Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 20 May 2022 06:30:07 +0200 Subject: [PATCH] avoid negative zero when not needed --- src/material.f90 | 2 +- src/rotations.f90 | 10 +++++----- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/src/material.f90 b/src/material.f90 index ce4cabea9..eafc411cd 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -164,7 +164,7 @@ subroutine parse() call material_O_0(ma)%data(co)%fromQuaternion(constituent%get_as1dFloat('O',requiredSize=4)) material_V_e_0(ma)%data(1:3,1:3,co) = constituent%get_as2dFloat('V_e',defaultVal=math_I3,requiredShape=[3,3]) if (any(dNeq(material_V_e_0(ma)%data(1:3,1:3,co),transpose(material_V_e_0(ma)%data(1:3,1:3,co))))) & - call IO_error(147) + call IO_error(147) end do if (dNeq(sum(v_of(ma,:)),1.0_pReal,1.e-9_pReal)) call IO_error(153,ext_msg='constituent') diff --git a/src/rotations.f90 b/src/rotations.f90 index 3340bf2a0..6257fe7e4 100644 --- a/src/rotations.f90 +++ b/src/rotations.f90 @@ -574,7 +574,7 @@ end function qu2cu !--------------------------------------------------------------------------------------------------- !> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH -!> @brief convert rotation matrix to cubochoric +!> @brief convert rotation matrix to unit quaternion !> @details the original formulation (direct conversion) had (numerical?) issues !--------------------------------------------------------------------------------------------------- pure function om2qu(om) result(qu) @@ -601,14 +601,14 @@ pure function om2qu(om) result(qu) endif endif if(sign(1.0_pReal,qu(1))<0.0_pReal) qu =-1.0_pReal * qu - qu = qu*[1.0_pReal,P,P,P] + qu(2:4) = merge(qu(2:4),qu(2:4)*P,dEq0(qu(2:4))) end function om2qu !--------------------------------------------------------------------------------------------------- !> @author Marc De Graef, Carnegie Mellon University -!> @brief orientation matrix to Euler angles +!> @brief convert orientation matrix to Euler angles !> @details Two step check for special cases to avoid invalid operations (not needed for python) !--------------------------------------------------------------------------------------------------- pure function om2eu(om) result(eu) @@ -1333,8 +1333,8 @@ pure function cu2ho(cu) result(ho) ! transform to sphere grid (inverse Lambert) ! [note that there is no need to worry about dividing by zero, since XYZ(3) can not become zero] c = sum(T**2) - s = PI * c/(24.0*XYZ(3)**2) - c = sqrt(PI) * c / sqrt(24.0_pReal) / XYZ(3) + s = c * PI/(24.0*XYZ(3)**2) + c = c * sqrt(PI/24.0_pReal) / XYZ(3) q = sqrt( 1.0 - s ) LamXYZ = [ T(order(2)) * q, T(order(1)) * q, PREF * XYZ(3) - c ] end if special