diff --git a/python/damask/_rotation.py b/python/damask/_rotation.py index 723664e18..42b32552b 100644 --- a/python/damask/_rotation.py +++ b/python/damask/_rotation.py @@ -1057,15 +1057,15 @@ class Rotation: """ if len(ho.shape) == 1: - ball_ = ho/np.linalg.norm(ho)*_R1 if np.isclose(np.linalg.norm(ho),_R1,atol=1e-6) \ + ho_ = ho/np.linalg.norm(ho)*_R1 if np.isclose(np.linalg.norm(ho),_R1,atol=1e-6) \ else ho - rs = np.linalg.norm(ball_) + rs = np.linalg.norm(ho_) - if np.allclose(ball_,0.0,rtol=0.0,atol=1.0e-16): - cube = np.zeros(3) + if np.allclose(ho_,0.0,rtol=0.0,atol=1.0e-16): + cu = np.zeros(3) else: - p = _get_order(ball_) - xyz3 = ball_[p[0]] + p = _get_order(ho_) + xyz3 = ho_[p[0]] # inverse M_3 xyz2 = xyz3[0:2] * np.sqrt( 2.0*rs/(rs+np.abs(xyz3[2])) ) @@ -1085,11 +1085,11 @@ class Rotation: Tinv = q * np.where(xyz2<0.0,-Tinv,Tinv) # inverse M_1 - cube = np.array([ Tinv[0], Tinv[1], (-1.0 if xyz3[2] < 0.0 else 1.0) * rs / np.sqrt(6.0/np.pi) ]) /_sc + cu = np.array([ Tinv[0], Tinv[1], (-1.0 if xyz3[2] < 0.0 else 1.0) * rs / np.sqrt(6.0/np.pi) ]) /_sc # reverse the coordinates back to the regular order according to the original pyramid number - cube = cube[p[1]] + cu = cu[p[1]] - return cube + return cu else: raise NotImplementedError @@ -1133,20 +1133,20 @@ class Rotation: """ if len(cu.shape) == 1: - cube_ = np.clip(cu,None,np.pi**(2./3.) * 0.5) if np.isclose(np.abs(np.max(cu)),np.pi**(2./3.) * 0.5,atol=1e-6) \ + cu_ = np.clip(cu,None,np.pi**(2./3.) * 0.5) if np.isclose(np.abs(np.max(cu)),np.pi**(2./3.) * 0.5,atol=1e-6) \ else cu # transform to the sphere grid via the curved square, and intercept the zero point - if np.allclose(cube_,0.0,rtol=0.0,atol=1.0e-16): - ball = np.zeros(3) + if np.allclose(cu_,0.0,rtol=0.0,atol=1.0e-16): + ho = np.zeros(3) else: # get pyramide and scale by grid parameter ratio - p = _get_order(cube_) - XYZ = cube_[p[0]] * _sc + p = _get_order(cu_) + XYZ = cu_[p[0]] * _sc # intercept all the points along the z-axis if np.allclose(XYZ[0:2],0.0,rtol=0.0,atol=1.0e-16): - ball = np.array([0.0, 0.0, np.sqrt(6.0/np.pi) * XYZ[2]]) + ho = np.array([0.0, 0.0, np.sqrt(6.0/np.pi) * XYZ[2]]) else: order = [1,0] if np.abs(XYZ[1]) <= np.abs(XYZ[0]) else [0,1] q = np.pi/12.0 * XYZ[order[0]]/XYZ[order[1]] @@ -1161,12 +1161,12 @@ class Rotation: s = c * np.pi/24.0 /XYZ[2]**2 c = c * np.sqrt(np.pi/24.0)/XYZ[2] q = np.sqrt( 1.0 - s ) - ball = np.array([ T[order[1]] * q, T[order[0]] * q, np.sqrt(6.0/np.pi) * XYZ[2] - c ]) + ho = np.array([ T[order[1]] * q, T[order[0]] * q, np.sqrt(6.0/np.pi) * XYZ[2] - c ]) # reverse the coordinates back to the regular order according to the original pyramid number - ball = ball[p[1]] + ho = ho[p[1]] - return ball + return ho else: raise NotImplementedError diff --git a/src/rotations.f90 b/src/rotations.f90 index 8ff5aca4e..f3110b7d4 100644 --- a/src/rotations.f90 +++ b/src/rotations.f90 @@ -1233,28 +1233,28 @@ end function cu2ho !> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH !> @brief map from 3D cubic grid to 3D ball !-------------------------------------------------------------------------- -pure function Lambert_CubeToBall(cube) result(ball) +pure function Lambert_CubeToBall(cu) result(ho) - real(pReal), intent(in), dimension(3) :: cube - real(pReal), dimension(3) :: ball, LamXYZ, XYZ + real(pReal), intent(in), dimension(3) :: cu + real(pReal), dimension(3) :: ho, LamXYZ, XYZ real(pReal), dimension(2) :: T real(pReal) :: c, s, q real(pReal), parameter :: eps = 1.0e-8_pReal integer, dimension(3,2) :: p integer, dimension(2) :: order - if (maxval(abs(cube)) > AP/2.0+eps) then - ball = IEEE_value(cube,IEEE_positive_inf) + if (maxval(abs(cu)) > AP/2.0+eps) then + ho = IEEE_value(cu,IEEE_positive_inf) return end if ! transform to the sphere grid via the curved square, and intercept the zero point - center: if (all(dEq0(cube))) then - ball = 0.0_pReal + center: if (all(dEq0(cu))) then + ho = 0.0_pReal else center ! get pyramide and scale by grid parameter ratio - p = GetPyramidOrder(cube) - XYZ = cube(p(:,1)) * sc + p = GetPyramidOrder(cu) + XYZ = cu(p(:,1)) * sc ! intercept all the points along the z-axis special: if (all(dEq0(XYZ(1:2)))) then @@ -1277,7 +1277,7 @@ pure function Lambert_CubeToBall(cube) result(ball) endif special ! reverse the coordinates back to order according to the original pyramid number - ball = LamXYZ(p(:,2)) + ho = LamXYZ(p(:,2)) endif center @@ -1289,25 +1289,25 @@ end function Lambert_CubeToBall !> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH !> @brief map from 3D ball to 3D cubic grid !-------------------------------------------------------------------------- -pure function Lambert_BallToCube(xyz) result(cube) +pure function Lambert_BallToCube(ho) result(cu) - real(pReal), intent(in), dimension(3) :: xyz - real(pReal), dimension(3) :: cube, xyz1, xyz3 + real(pReal), intent(in), dimension(3) :: ho + real(pReal), dimension(3) :: cu, xyz1, xyz3 real(pReal), dimension(2) :: Tinv, xyz2 real(pReal) :: rs, qxy, q2, sq2, q, tt integer, dimension(3,2) :: p - rs = norm2(xyz) + rs = norm2(ho) if (rs > R1+1.e-6_pReal) then - cube = IEEE_value(cube,IEEE_positive_inf) + cu = IEEE_value(cu,IEEE_positive_inf) return endif - center: if (all(dEq0(xyz))) then - cube = 0.0_pReal + center: if (all(dEq0(ho))) then + cu = 0.0_pReal else center - p = GetPyramidOrder(xyz) - xyz3 = xyz(p(:,1)) + p = GetPyramidOrder(ho) + xyz3 = ho(p(:,1)) ! inverse M_3 xyz2 = xyz3(1:2) * sqrt( 2.0*rs/(rs+abs(xyz3(3))) ) @@ -1331,7 +1331,7 @@ pure function Lambert_BallToCube(xyz) result(cube) xyz1 = [ Tinv(1), Tinv(2), sign(1.0_pReal,xyz3(3)) * rs / pref ] /sc ! reverse the coordinates back to order according to the original pyramid number - cube = xyz1(p(:,2)) + cu = xyz1(p(:,2)) endif center