diff --git a/src/layer/x86/sse_mathfun.h b/src/layer/x86/sse_mathfun.h index 4df439323..c6c5e19f6 100644 --- a/src/layer/x86/sse_mathfun.h +++ b/src/layer/x86/sse_mathfun.h @@ -836,4 +836,89 @@ static NCNN_FORCEINLINE __m128 floor_ps(__m128 x) _mm_andnot_ps(no_fraction, fixed_result)); } +static NCNN_FORCEINLINE __m128 asin_ps(__m128 x) +{ + const __m128 magic_negative_zero = _mm_set_ps1(-0.0f); + const __m128 magic_half_one = _mm_set_ps1(0.5f); + const __m128 magic_one = _mm_set_ps1(1.0f); + const __m128 magic_a4 = _mm_set_ps1(0.023994016f); + const __m128 magic_a5 = _mm_set_ps1(0.042417344f); + const __m128 magic_a2 = _mm_set_ps1(0.07494697f); + const __m128 magic_a3 = _mm_set_ps1(0.045520633f); + const __m128 magic_a0 = _mm_set_ps1(1.0f); + const __m128 magic_a1 = _mm_set_ps1(0.166667819f); + const __m128 magic_half_pi = _mm_set_ps1(1.5707964f); + const __m128 magic_three = _mm_set_ps1(3.0f); + + // negative_mask = magic_negative_zero && x; + __m128 negative_mask = _mm_and_ps(magic_negative_zero, x); + + // absolute = abs(x); + __m128 absolute = _mm_andnot_ps(magic_negative_zero, x); + + // Reference: https://en.wikipedia.org/wiki/Small-angle_approximation + + // is_small_input = (absolute <= 0.5f); + __m128 is_small_input = _mm_cmple_ps(absolute, magic_half_one); + + // is_big_input = (!is_small_input); + __m128 is_big_input = _mm_andnot_ps(is_small_input, magic_one); + + // big_input_approx = sqrt(0.5f * (1 - absolute)); + __m128 big_input_approx = _mm_sqrt_ps(_mm_mul_ps( + magic_half_one, + _mm_sub_ps(magic_one, absolute))); + + // input_approx = (is_small_input ? absolute : big_input_approx); + __m128 input_approx = _mm_or_ps( + _mm_and_ps(is_small_input, absolute), + _mm_andnot_ps(is_small_input, big_input_approx)); + + // square_of_input_approx = input_approx * input_approx; + __m128 square_of_input_approx = _mm_mul_ps(input_approx, input_approx); + + // fourth_power_of_input_approx = + // square_of_input_approx * square_of_input_approx; + __m128 fourth_power_of_input_approx = _mm_mul_ps( + square_of_input_approx, square_of_input_approx); + + // TODO: Need more explanations. + // x1 = ((magic_a4 * fourth_power_of_input_approx) + magic_a2); + // x2 = ((magic_a5 * fourth_power_of_input_approx) + magic_a3); + // x3 = ((x1 * fourth_power_of_input_approx) + magic_a0); + // x4 = ((fourth_power_of_input_approx * x2) + magic_a1); + // output_approx = (x3 + (square_of_input_approx * x4)); + __m128 output_approx = _mm_add_ps( + _mm_add_ps( + _mm_mul_ps( + _mm_add_ps( + _mm_mul_ps(magic_a4, fourth_power_of_input_approx), + magic_a2), + fourth_power_of_input_approx), + magic_a0), + _mm_mul_ps( + square_of_input_approx, + _mm_add_ps( + _mm_mul_ps( + fourth_power_of_input_approx, + _mm_add_ps( + _mm_mul_ps(magic_a5, fourth_power_of_input_approx), + magic_a3)), + magic_a1))); + + // TODO: Need more explanations. + // x1 = ((0.5 * PI) * is_big_input); + // x2 = (output_approx * input_approx); + // x3 = (1.0f - (3.0f * is_big_input)); + // final_approx = (x1 + (x2 * x3)); + __m128 final_approx = _mm_add_ps( + _mm_mul_ps(magic_half_pi, is_big_input), + _mm_mul_ps( + _mm_mul_ps(output_approx, input_approx), + _mm_sub_ps(magic_one, _mm_mul_ps(magic_three, is_big_input)))); + + // return (final_approx || negative_mask); + return _mm_or_ps(final_approx, negative_mask); +} + #endif // SSE_MATHFUN_H diff --git a/src/layer/x86/unaryop_x86.cpp b/src/layer/x86/unaryop_x86.cpp index 987eaffcf..ecb361f40 100644 --- a/src/layer/x86/unaryop_x86.cpp +++ b/src/layer/x86/unaryop_x86.cpp @@ -422,14 +422,7 @@ struct unary_op_asin #if __SSE2__ __m128 func_pack4(const __m128& x) const { - //TODO sse optimize - float tmp[4]; - _mm_storeu_ps(tmp, x); - tmp[0] = asin(tmp[0]); - tmp[1] = asin(tmp[1]); - tmp[2] = asin(tmp[2]); - tmp[3] = asin(tmp[3]); - return _mm_loadu_ps(tmp); + return asin_ps(x); } #if __AVX__ __m256 func_pack8(const __m256& x) const