/* NOTE simon (28/04/25 17:37:11): cl simd_test.c -Fesimd_test.exe -nologo -GR- -EHa- -Oi -W4 -wd4100 -wd4189 -wd4201 -FC -Zi -Zo -diagnostics:caret -diagnostics:column -O2 -MT -GS- -WX -std:c11 -link -INCREMENTAL:NO -opt:ref,icf */ #include #include #include #include #define _USE_MATH_DEFINES #include typedef float f32; typedef double f64; typedef uintptr_t umm; typedef uint8_t u8; typedef uint32_t u32; typedef uint64_t u64; typedef int32_t s32; typedef uint32_t b32; #define cast( type, expr ) ( type ) ( expr ) #define array_count( array ) ( ( umm ) ( sizeof( array ) / sizeof( ( array )[ 0 ] ) ) ) #define _assert( expr ) if ( !( expr ) ) { __debugbreak( ); } static f32 interpolation_linear( f32 start, f32 end, f32 t ) { f32 result = start * ( 1 - t ) + end * t; return result; } static f32 interpolation_cubic( f32 start, f32 cp_1, f32 cp_2, f32 end, f32 t ) { f32 u = 1 - t; f32 result = ( u * u * u * start ) + ( 3 * u * u * t * cp_1 ) + ( 3 * u * t * t * cp_2 ) + ( t * t * t * end ); return result; } #include #include #include _Alignas( 64 ) u8 shuffle_lut[ 16 * 16 ] = { 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b0000 0, 1, 2, 3, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b0001 4, 5, 6, 7, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b0010 0, 1, 2, 3, 4, 5, 6, 7, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b0011 8, 9, 10, 11, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b0100 0, 1, 2, 3, 8, 9, 10, 11, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b0101 4, 5, 6, 7, 8, 9, 10, 11, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b0110 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 0x80, 0x80, 0x80, 0x80, // 0b0111 12, 13, 14, 15, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b1000 0, 1, 2, 3, 12, 13, 14, 15, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b1001 4, 5, 6, 7, 12, 13, 14, 15, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b1010 0, 1, 2, 3, 4, 5, 6, 7, 12, 13, 14, 15, 0x80, 0x80, 0x80, 0x80, // 0b1011 8, 9, 10, 11, 12, 13, 14, 15, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, 0x80, // 0b1100 0, 1, 2, 3, 8, 9, 10, 11, 12, 13, 14, 15, 0x80, 0x80, 0x80, 0x80, // 0b1101 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 0x80, 0x80, 0x80, 0x80, // 0b1110 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15 // 0b1111 }; _Alignas( 64 ) u8 shuffle_lut_2[ 16 * 16 ] = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, // 0b0000 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, // 0b0001 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, // 0b0010 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, // 0b0011 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, // 0b0100 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, // 0b0101 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, // 0b0110 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, // 0b0111 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, // 0b1000 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, // 0b1001 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, // 0b1010 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, // 0b1011 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, // 0b1100 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, // 0b1101 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, // 0b1110 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, // 0b1111 }; _Alignas( 64 ) u8 shuffle_lut_3[ 16 * 16 ] = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, // 0b0000 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, // 0b0001 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, // 0b0010 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, // 0b0011 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, // 0b0100 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, // 0b0101 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, // 0b0110 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, // 0b0111 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, 12, 13, 14, 15, // 0b1000 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, // 0b1001 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, // 0b1010 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, // 0b1011 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, 8, 9, 10, 11, // 0b1100 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, // 0b1101 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, 4, 5, 6, 7, // 0b1110 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, // 0b1111 }; static f32 interpolation_cubic_t_from_x_bisection_simd_1( f32 start, f32 cp_1, f32 cp_2, f32 end, f32 x, f32 guess, f32 epsilon ) { f32 a = ( -start + 3 * cp_1 - 3 * cp_2 + end ); f32 b = ( 3 * start - 6 * cp_1 + 3 * cp_2 ); f32 c = ( -3 * start + 3 * cp_1 ); f32 d = ( start - x ); f32 result = guess; /* NOTE simon (24/04/25 18:41:29): The order is inverted to my mental model. */ __m128 factors_1 = _mm_set_ps( 0.2f, 0.4f, 0.6f, 0.8f ); __m128 factors_2 = _mm_set_ps( 0.8f, 0.6f, 0.4f, 0.2f ); __m128 ts = factors_2; __m128 as = _mm_set1_ps( a ); __m128 bs = _mm_set1_ps( b ); __m128 cs = _mm_set1_ps( c ); __m128 ds = _mm_set1_ps( d ); __m128 epsilons = _mm_set_ps1( epsilon ); __m128 xs = _mm_set1_ps( x ); __m128 ls = _mm_set1_ps( 0 ); __m128 rs = _mm_set1_ps( 1 ); while ( 1 ) { /* NOTE simon (26/04/25 13:57:04): Compute cubic. */ __m128 tt = _mm_mul_ps( ts, ts ); __m128 ttt = _mm_mul_ps( tt, ts ); __m128 attt = _mm_mul_ps( as, ttt ); __m128 btt = _mm_mul_ps( bs, tt ); __m128 ct = _mm_mul_ps( cs, ts ); __m128 f1 = _mm_add_ps( attt, btt ); __m128 f2 = _mm_add_ps( ct, ds ); __m128 fr = _mm_add_ps( f1, f2 ); /* NOTE simon (26/04/25 13:57:12): Absolute value */ __m128 n_zero = _mm_set1_ps( -0.0f ); __m128 abs = _mm_andnot_ps( n_zero, fr ); __m128 cmp = _mm_cmple_ps( abs, epsilons ); /* NOTE simon (26/04/25 13:57:37): Have we found the result ? */ u32 mask = _mm_movemask_ps( cmp ); if ( mask ) { /* NOTE simon (26/04/25 13:58:08): Extract the result. */ /* NOTE simon (28/04/25 15:37:10): Those two line produce the same code, except one tends to add security cookie. => /GS-*/ __m128i shuffle_mask = *( __m128i* ) ( shuffle_lut + ( mask * 16 ) ); // __m128i shuffle_mask = _mm_load_si128( ( __m128i* ) ( shuffle_lut + ( mask * 16 ) ) ); __m128i shuffled = _mm_shuffle_epi8( _mm_castps_si128( ts ), shuffle_mask ); result = _mm_cvtss_f32( _mm_castsi128_ps( shuffled ) ); break; } /* NOTE simon (26/04/25 13:58:20): Update left and right limits. */ __m128 fxs = _mm_sub_ps( fr, ds ); __m128 setl = _mm_cmplt_ps( fxs, xs ); // __m128 setr = _mm_cmpgt_ps( fxs, xs ); u32 lmask = _mm_movemask_ps( setl ); // u32 rmask = _mm_movemask_ps( setr ); u32 rmask = 15 ^ lmask; #if 1 if ( lmask ) #endif { __m128i shuffle_mask = *( __m128i* ) ( shuffle_lut_2 + ( lmask * 16 ) ); /* NOTE simon (29/04/25 16:20:07): Storing in a variable isn't great. The compiler can still do the right thing and use the suffle with the memory operand, but it's not guaranteed. */ // __m128i shuffle_mask = _mm_load_si128( ( __m128i* ) ( shuffle_lut_2 + ( lmask * 16 ) ) ); __m128i shuffled = _mm_shuffle_epi8( _mm_castps_si128( ts ), shuffle_mask ); // __m128i l_blend_mask = _mm_set1_epi8( lmask ? 0xff : 0 ); // ls = _mm_blendv_ps( ls, _mm_castsi128_ps( shuffled ), _mm_castsi128_ps( l_blend_mask ) ); ls = _mm_castsi128_ps( shuffled ); } #if 1 if ( rmask ) #endif { __m128i shuffle_mask = *( __m128i* ) ( shuffle_lut_3 + ( rmask * 16 ) ); // __m128i shuffle_mask = _mm_load_si128( ( __m128i* ) ( shuffle_lut_3 + ( rmask * 16 ) ) ); __m128i shuffled = _mm_shuffle_epi8( _mm_castps_si128( ts ), shuffle_mask ); // __m128i r_blend_mask = _mm_set1_epi8( rmask ? 0xff : 0 ); // rs = _mm_blendv_ps( rs, _mm_castsi128_ps( shuffled ), _mm_castsi128_ps( r_blend_mask ) ); rs = _mm_castsi128_ps( shuffled ); } /* NOTE simon (24/04/25 18:12:28): Lerp */ __m128 left = _mm_mul_ps( ls, factors_1 ); __m128 right = _mm_mul_ps( rs, factors_2 ); ts = _mm_add_ps( left, right ); } return result; } static f32 interpolation_cubic_t_from_x_bisection_simd_2( f32 start, f32 cp_1, f32 cp_2, f32 end, f32 x, f32 guess, f32 epsilon ) { f32 a = ( -start + 3 * cp_1 - 3 * cp_2 + end ); f32 b = ( 3 * start - 6 * cp_1 + 3 * cp_2 ); f32 c = ( -3 * start + 3 * cp_1 ); f32 d = ( start - x ); f32 result = guess; /* NOTE simon (24/04/25 18:41:29): The order is inverted to my mental model. */ __m128 factors_1 = _mm_set_ps( 0.2f, 0.4f, 0.6f, 0.8f ); __m128 factors_2 = _mm_set_ps( 0.8f, 0.6f, 0.4f, 0.2f ); __m128 ts = factors_2; __m128 as = _mm_set1_ps( a ); __m128 bs = _mm_set1_ps( b ); __m128 cs = _mm_set1_ps( c ); __m128 ds = _mm_set1_ps( d ); __m128 epsilons = _mm_set_ps1( epsilon ); __m128 xs = _mm_set1_ps( x ); __m128 ls = _mm_set1_ps( 0 ); __m128 rs = _mm_set1_ps( 1 ); while ( 1 ) { /* NOTE simon (26/04/25 13:57:04): Compute cubic. */ __m128 tt = _mm_mul_ps( ts, ts ); __m128 ttt = _mm_mul_ps( tt, ts ); __m128 attt = _mm_mul_ps( as, ttt ); __m128 btt = _mm_mul_ps( bs, tt ); __m128 ct = _mm_mul_ps( cs, ts ); __m128 f1 = _mm_add_ps( attt, btt ); __m128 f2 = _mm_add_ps( ct, ds ); __m128 fr = _mm_add_ps( f1, f2 ); /* NOTE simon (26/04/25 13:57:12): Absolute value */ __m128 n_zero = _mm_set1_ps( -0.0f ); __m128 abs = _mm_andnot_ps( n_zero, fr ); /* NOTE simon (26/04/25 13:57:37): Have we found the result ? */ __m128 cmp = _mm_cmple_ps( abs, epsilons ); u32 mask = _mm_movemask_ps( cmp ); if ( mask ) { /* NOTE simon (26/04/25 13:58:08): Extract the result. */ __m128 shuffled; if ( mask & 0x1 ) { shuffled = _mm_shuffle_ps( ts, ts, 0 ); } else if ( mask & 0x2 ) { shuffled = _mm_shuffle_ps( ts, ts, 1 ); } else if ( mask & 0x4 ) { shuffled = _mm_shuffle_ps( ts, ts, 2 ); } else { // _assert( mask & 0x8 ); shuffled = _mm_shuffle_ps( ts, ts, 3 ); } result = _mm_cvtss_f32( shuffled ); break; } /* NOTE simon (26/04/25 13:58:20): Update left and right limits. */ __m128 fxs = _mm_sub_ps( fr, ds ); __m128 setl = _mm_cmplt_ps( fxs, xs ); u32 lmask = _mm_movemask_ps( setl ); u32 rmask = ~lmask; if ( lmask & 0x8 ) { // 0b11111111 ls = _mm_shuffle_ps( ts, ts, 0xff ); } else if ( lmask & 0x4 ) { // 0b10101010 ls = _mm_shuffle_ps( ts, ts, 0xaa ); } else if ( lmask & 0x2 ) { // 0b01010101 ls = _mm_shuffle_ps( ts, ts, 0x55 ); } else if ( lmask & 0x1 ) { // 0b00000000 ls = _mm_shuffle_ps( ts, ts, 0 ); } if ( rmask & 0x1 ) { rs = _mm_shuffle_ps( ts, ts, 0 ); } else if ( rmask & 0x2 ) { rs = _mm_shuffle_ps( ts, ts, 0x55 ); } else if ( rmask & 0x4 ) { rs = _mm_shuffle_ps( ts, ts, 0xaa ); } else if ( rmask & 0x8 ) { rs = _mm_shuffle_ps( ts, ts, 0xff ); } /* NOTE simon (24/04/25 18:12:28): Lerp */ __m128 left = _mm_mul_ps( ls, factors_1 ); __m128 right = _mm_mul_ps( rs, factors_2 ); ts = _mm_add_ps( left, right ); } return result; } static f32 interpolation_cubic_t_from_x_bisection_simd_3( f32 start, f32 cp_1, f32 cp_2, f32 end, f32 x, f32 guess, f32 epsilon ) { f32 a = ( -start + 3 * cp_1 - 3 * cp_2 + end ); f32 b = ( 3 * start - 6 * cp_1 + 3 * cp_2 ); f32 c = ( -3 * start + 3 * cp_1 ); f32 d = ( start - x ); f32 result = guess; /* NOTE simon (24/04/25 18:41:29): The order is inverted to my mental model. */ __m128 factors_1 = _mm_set_ps( 0.2f, 0.4f, 0.6f, 0.8f ); __m128 factors_2 = _mm_set_ps( 0.8f, 0.6f, 0.4f, 0.2f ); __m128 ts = factors_2; __m128 as = _mm_set1_ps( a ); __m128 bs = _mm_set1_ps( b ); __m128 cs = _mm_set1_ps( c ); __m128 ds = _mm_set1_ps( d ); __m128 epsilons = _mm_set_ps1( epsilon ); __m128 xs = _mm_set1_ps( x ); __m128 ls = _mm_set1_ps( 0 ); __m128 rs = _mm_set1_ps( 1 ); while ( 1 ) { /* NOTE simon (26/04/25 13:57:04): Compute cubic. */ __m128 tt = _mm_mul_ps( ts, ts ); __m128 ttt = _mm_mul_ps( tt, ts ); __m128 attt = _mm_mul_ps( as, ttt ); __m128 btt = _mm_mul_ps( bs, tt ); __m128 ct = _mm_mul_ps( cs, ts ); __m128 f1 = _mm_add_ps( attt, btt ); __m128 f2 = _mm_add_ps( f1, ct ); __m128 fr = _mm_add_ps( f2, ds ); /* NOTE simon (26/04/25 13:57:12): Absolute value */ __m128 n_zero = _mm_set1_ps( -0.0f ); __m128 abs = _mm_andnot_ps( n_zero, fr ); __m128 cmp = _mm_cmple_ps( abs, epsilons ); /* NOTE simon (26/04/25 13:57:37): Have we found the result ? */ u32 mask = _mm_movemask_ps( cmp ); if ( mask ) { /* NOTE simon (26/04/25 13:58:08): Extract the result. */ __m128 shuffled; if ( mask & 0x1 ) { shuffled = _mm_shuffle_ps( ts, ts, 0 ); } else if ( mask & 0x2 ) { shuffled = _mm_shuffle_ps( ts, ts, 1 ); } else if ( mask & 0x4 ) { shuffled = _mm_shuffle_ps( ts, ts, 2 ); } else { // _assert( mask & 0x8 ); shuffled = _mm_shuffle_ps( ts, ts, 3 ); } result = _mm_cvtss_f32( shuffled ); break; } /* NOTE simon (26/04/25 13:58:20): Update left and right limits. */ __m128 setl = _mm_cmplt_ps( f2, xs ); u32 lmask = _mm_movemask_ps( setl ); u32 rmask = ~lmask; if ( lmask & 0x8 ) { // 0b11111111 ls = _mm_shuffle_ps( ts, ts, 0xff ); } else if ( lmask & 0x4 ) { // 0b10101010 ls = _mm_shuffle_ps( ts, ts, 0xaa ); } else if ( lmask & 0x2 ) { // 0b01010101 ls = _mm_shuffle_ps( ts, ts, 0x55 ); } else if ( lmask & 0x1 ) { // 0b00000000 ls = _mm_shuffle_ps( ts, ts, 0 ); } if ( rmask & 0x1 ) { rs = _mm_shuffle_ps( ts, ts, 0 ); } else if ( rmask & 0x2 ) { rs = _mm_shuffle_ps( ts, ts, 0x55 ); } else if ( rmask & 0x4 ) { rs = _mm_shuffle_ps( ts, ts, 0xaa ); } else if ( rmask & 0x8 ) { rs = _mm_shuffle_ps( ts, ts, 0xff ); } /* NOTE simon (24/04/25 18:12:28): Lerp */ __m128 left = _mm_mul_ps( ls, factors_1 ); __m128 right = _mm_mul_ps( rs, factors_2 ); ts = _mm_add_ps( left, right ); } return result; } static f32 interpolation_cubic_t_from_x_bisection( f32 start, f32 cp_1, f32 cp_2, f32 end, f32 x, f32 guess, f32 epsilon ) { f32 a = ( -start + 3 * cp_1 - 3 * cp_2 + end ); f32 b = ( 3 * start - 6 * cp_1 + 3 * cp_2 ); f32 c = ( -3 * start + 3 * cp_1 ); f32 d = ( start - x ); f32 l = 0; f32 w = 1; f32 t = 0; while ( 1 ) { w *= 0.5f; t = l + w; f32 f = a * t * t * t + b * t * t + c * t; f32 f_d = f + d; if ( fabsf( f_d ) < epsilon ) { break; } l = ( f > x ) ? l : t; } return t; } static f32 interpolation_cubic_t_from_x_bisection_simd_4( f32 start, f32 cp_1, f32 cp_2, f32 end, f32 x, f32 guess, f32 epsilon ) { f32 t = 0; f32 a = ( -start + 3 * cp_1 - 3 * cp_2 + end ); f32 b = ( 3 * start - 6 * cp_1 + 3 * cp_2 ); f32 c = ( -3 * start + 3 * cp_1 ); f32 d = ( start - x ); __m128 as = _mm_set1_ps( a ); __m128 bs = _mm_set1_ps( b ); __m128 cs = _mm_set1_ps( c ); __m128 ds = _mm_set1_ps( d ); __m128 epsilons = _mm_set_ps1( epsilon ); __m128 xs = _mm_set_ps1( x ); __m128 ls = _mm_set_ps1( 0 ); __m128 factors = _mm_set_ps( 0.8f, 0.6f, 0.4f, 0.2f ); __m128 distance = _mm_set_ps1( 1.0f ); __m128 scale = _mm_set_ps1( 0.2f ); while ( 1 ) { __m128 ws = _mm_mul_ps( distance, factors ); __m128 ts = _mm_add_ps( ls, ws ); /* NOTE simon (26/04/25 13:57:04): Compute cubic. */ __m128 tt = _mm_mul_ps( ts, ts ); __m128 ttt = _mm_mul_ps( tt, ts ); __m128 attt = _mm_mul_ps( as, ttt ); __m128 btt = _mm_mul_ps( bs, tt ); __m128 ct = _mm_mul_ps( cs, ts ); __m128 f1 = _mm_add_ps( attt, btt ); __m128 f2 = _mm_add_ps( ct, ds ); __m128 fr = _mm_add_ps( f1, f2 ); /* NOTE simon (26/04/25 13:57:12): Absolute value */ __m128 n_zero = _mm_set1_ps( -0.0f ); __m128 abs = _mm_andnot_ps( n_zero, fr ); /* NOTE simon (26/04/25 13:57:37): Have we found the result ? */ __m128 cmp = _mm_cmple_ps( abs, epsilons ); u32 mask = _mm_movemask_ps( cmp ); if ( mask ) { /* NOTE simon (26/04/25 13:58:08): Extract the result. */ __m128 shuffled; if ( mask & 0x1 ) { shuffled = _mm_shuffle_ps( ts, ts, 0 ); } else if ( mask & 0x2 ) { shuffled = _mm_shuffle_ps( ts, ts, 1 ); } else if ( mask & 0x4 ) { shuffled = _mm_shuffle_ps( ts, ts, 2 ); } else { // _assert( mask & 0x8 ); shuffled = _mm_shuffle_ps( ts, ts, 3 ); } t = _mm_cvtss_f32( shuffled ); break; } /* NOTE simon (26/04/25 13:58:20): Update left and distance. */ __m128 fxs = _mm_sub_ps( fr, ds ); __m128 setl = _mm_cmplt_ps( fxs, xs ); u32 lmask = _mm_movemask_ps( setl ); if ( lmask & 0x8 ) { // 0b11111111 ls = _mm_shuffle_ps( ts, ts, 0xff ); } else if ( lmask & 0x4 ) { // 0b10101010 ls = _mm_shuffle_ps( ts, ts, 0xaa ); } else if ( lmask & 0x2 ) { // 0b01010101 ls = _mm_shuffle_ps( ts, ts, 0x55 ); } else if ( lmask & 0x1 ) { // 0b00000000 ls = _mm_shuffle_ps( ts, ts, 0 ); } distance = _mm_mul_ps( distance, scale ); } return t; } int main( int argc, char** argv ) { #if 1 { u64 total = 0; u64 hit = 0; for ( u32 i = 0; i < 101; i++ ) { for ( u32 j = 0; j < 101; j++ ) { for ( u32 k = 0; k < 101; k++ ) { for ( u32 l = 0; l < 101; l++ ) { f32 cp_1 = i * 0.01f; f32 cp_2 = j * 0.01f; f32 x = k * 0.01f; f32 guess = l * 0.01f; u64 start = __rdtsc( ); f32 t = interpolation_cubic_t_from_x_bisection_simd_1( 0, cp_1, cp_2, 1, x, guess, 0.00001f * 0.75f ); u64 end = __rdtsc( ); f32 check = interpolation_cubic( 0, cp_1, cp_2, 1, t ); _assert( fabsf( check - x ) < 0.00001f ); total += ( end - start ); hit++; } } } } printf( "left pack total: %lld | hit: %lld | avg: %lld\n", total, hit, total / hit ); } #endif #if 1 { u64 total = 0; u64 hit = 0; for ( u32 i = 0; i < 101; i++ ) { for ( u32 j = 0; j < 101; j++ ) { for ( u32 k = 0; k < 101; k++ ) { for ( u32 l = 0; l < 101; l++ ) { f32 cp_1 = i * 0.01f; f32 cp_2 = j * 0.01f; f32 x = k * 0.01f; f32 guess = l * 0.01f; u64 start = __rdtsc( ); f32 t = interpolation_cubic_t_from_x_bisection_simd_2( 0, cp_1, cp_2, 1, x, guess, 0.00001f * 0.75f ); u64 end = __rdtsc( ); f32 check = interpolation_cubic( 0, cp_1, cp_2, 1, t ); _assert( fabsf( check - x ) < 0.00001f ); total += ( end - start ); hit++; } } } } printf( "ifs total: %lld | hit: %lld | avg: %lld\n", total, hit, total / hit ); } #endif #if 1 { u64 total = 0; u64 hit = 0; for ( u32 i = 0; i < 101; i++ ) { for ( u32 j = 0; j < 101; j++ ) { for ( u32 k = 0; k < 101; k++ ) { for ( u32 l = 0; l < 101; l++ ) { f32 cp_1 = i * 0.01f; f32 cp_2 = j * 0.01f; f32 x = k * 0.01f; f32 guess = l * 0.01f; u64 start = __rdtsc( ); f32 t = interpolation_cubic_t_from_x_bisection_simd_3( 0, cp_1, cp_2, 1, x, guess, 0.00001f * 0.75f ); u64 end = __rdtsc( ); f32 check = interpolation_cubic( 0, cp_1, cp_2, 1, t ); _assert( fabsf( check - x ) < 0.00001f ); total += ( end - start ); hit++; } } } } printf( "ifs 2 total: %lld | hit: %lld | avg: %lld\n", total, hit, total / hit ); } #endif #if 1 { u64 total = 0; u64 hit = 0; for ( u32 i = 0; i < 101; i++ ) { for ( u32 j = 0; j < 101; j++ ) { for ( u32 k = 0; k < 101; k++ ) { for ( u32 l = 0; l < 101; l++ ) { f32 cp_1 = i * 0.01f; f32 cp_2 = j * 0.01f; f32 x = k * 0.01f; f32 guess = l * 0.01f; u64 start = __rdtsc( ); f32 t = interpolation_cubic_t_from_x_bisection( 0, cp_1, cp_2, 1, x, guess, 0.00001f * 0.75f ); u64 end = __rdtsc( ); f32 check = interpolation_cubic( 0, cp_1, cp_2, 1, t ); _assert( fabsf( check - x ) < 0.00001f ); total += ( end - start ); hit++; } } } } printf( "left dist total: %lld | hit: %lld | avg: %lld\n", total, hit, total / hit ); } #endif { u64 total = 0; u64 hit = 0; for ( u32 i = 0; i < 101; i++ ) { for ( u32 j = 0; j < 101; j++ ) { for ( u32 k = 0; k < 101; k++ ) { for ( u32 l = 0; l < 101; l++ ) { f32 cp_1 = i * 0.01f; f32 cp_2 = j * 0.01f; f32 x = k * 0.01f; f32 guess = l * 0.01f; u64 start = __rdtsc( ); f32 t = interpolation_cubic_t_from_x_bisection_simd_4( 0, cp_1, cp_2, 1, x, guess, 0.00001f * 0.75f ); u64 end = __rdtsc( ); f32 check = interpolation_cubic( 0, cp_1, cp_2, 1, t ); _assert( fabsf( check - x ) < 0.00001f ); total += ( end - start ); hit++; } } } } printf( "simd4 total: %lld | hit: %lld | avg: %lld\n", total, hit, total / hit ); } return 0; }