mirror/eigen
2
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2024-11-27 06:30:28 +08:00
Code Issues Packages Projects Releases Wiki Activity

Use vrsqrts for rsqrt Newton iterations.

Browse Source 
It's slightly faster and slightly more accurate, allowing our current
packetmath tests to pass for sqrt with a single iteration.
This commit is contained in:
Antonio Sanchez 2021-02-11 11:33:51 -08:00
parent 9fde9cce5d
commit 90ee821c56
1 changed files with 4 additions and 7 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

11
Eigen/src/Core/arch/NEON/PacketMath.h
View File

@ -3294,26 +3294,23 @@ template<> EIGEN_STRONG_INLINE Packet4ui psqrt(const Packet4ui& a) {
// effective latency. This is similar to Quake3's fast inverse square root. // effective latency. This is similar to Quake3's fast inverse square root.
// For more details see: http://www.beyond3d.com/content/articles/8/ // For more details see: http://www.beyond3d.com/content/articles/8/
template<> EIGEN_STRONG_INLINE Packet4f psqrt(const Packet4f& _x){ template<> EIGEN_STRONG_INLINE Packet4f psqrt(const Packet4f& _x){
Packet4f minus_half_x = vmulq_n_f32(_x, -0.5f);
Packet4ui denormal_mask = vandq_u32(vcgeq_f32(_x, vdupq_n_f32(0.0f)), Packet4ui denormal_mask = vandq_u32(vcgeq_f32(_x, vdupq_n_f32(0.0f)),
vcltq_f32(_x, pset1<Packet4f>((std::numeric_limits<float>::min)()))); vcltq_f32(_x, pset1<Packet4f>((std::numeric_limits<float>::min)())));
// Compute approximate reciprocal sqrt. // Compute approximate reciprocal sqrt.
Packet4f x = vrsqrteq_f32(_x); Packet4f x = vrsqrteq_f32(_x);
// Do a single step of Newton's iteration. // Do one Newton's iteration for 1/sqrt(x).
//the number 1.5f was set reference to Quake3's fast inverse square root x = vmulq_f32(vrsqrtsq_f32(vmulq_f32(_x, x), x), x);
x = pmul(x, pmadd(minus_half_x, pmul(x, x), pset1<Packet4f>(1.5f)));
// Flush results for denormals to zero. // Flush results for denormals to zero.
return vreinterpretq_f32_u32(vbicq_u32(vreinterpretq_u32_f32(pmul(_x, x)), denormal_mask)); return vreinterpretq_f32_u32(vbicq_u32(vreinterpretq_u32_f32(pmul(_x, x)), denormal_mask));
} }
template<> EIGEN_STRONG_INLINE Packet2f psqrt(const Packet2f& _x) { template<> EIGEN_STRONG_INLINE Packet2f psqrt(const Packet2f& _x) {
Packet2f minus_half_x = vmul_n_f32(_x, -0.5f);
Packet2ui denormal_mask = vand_u32(vcge_f32(_x, vdup_n_f32(0.0f)), Packet2ui denormal_mask = vand_u32(vcge_f32(_x, vdup_n_f32(0.0f)),
vclt_f32(_x, pset1<Packet2f>((std::numeric_limits<float>::min)()))); vclt_f32(_x, pset1<Packet2f>((std::numeric_limits<float>::min)())));
// Compute approximate reciprocal sqrt. // Compute approximate reciprocal sqrt.
Packet2f x = vrsqrte_f32(_x); Packet2f x = vrsqrte_f32(_x);
// Do a single step of Newton's iteration. // Do one Newton's iteration for 1/sqrt(x).
x = pmul(x, pmadd(minus_half_x, pmul(x, x), pset1<Packet2f>(1.5f))); x = vmul_f32(vrsqrts_f32(vmul_f32(_x, x), x), x);
// Flush results for denormals to zero. // Flush results for denormals to zero.
return vreinterpret_f32_u32(vbic_u32(vreinterpret_u32_f32(pmul(_x, x)), denormal_mask)); return vreinterpret_f32_u32(vbic_u32(vreinterpret_u32_f32(pmul(_x, x)), denormal_mask));
} }