compile c++ gcc online

Language: Layout:

// This code is in the public domain. Where that dedication is not
// recognized, you are granted a perpetual, irrevocable license to copy
// and modify it as you see fit.

// WARNING 1: the code relies on rounding mode being round-to-nearest!
// WARNING 2: the error bounds below are not rigorously tested.

#include <cstdio>
#include <cmath>
#include <cstdint>
#include <cstring>
#include <ctime>
#include <emmintrin.h>

//==============================================================================

// SSE2 assembly version (scalar).

// |x|<65536:       rel. error < 1.9*10^-7 [ 2 ulp ]
// |x|<262144:      abs. error < 4*10^-6
// Garbage for |x|>8388608.
static inline float sin_f32_sse2(float x)
{
    static const float inv_pi[1]={0.318309886f};                //~ 1/(pi)
    static const float bias[1]={12582912.0f};                   // 1.5*2^23
    static const float pi0[1]={3.140625f};                      // bits -2..+9
    static const float pi1[1]={9.670257568359375e-4f};          // bits +10..+19
    static const float pi2[1]={6.2771141529083251953125e-7f};   // bits +20..+28
    static const float pi3[1]={1.2154201013012385202950e-10f};  // bits +29..+51
    static const float p0[1]={2.6050554015241175e-6f};          //
    static const float p1[1]={-1.9808833894785493e-4f};         // -13612527*2^-36
    static const float p2[1]={+8.3330469205975532e-3f};         // +8947541*2^-30
    static const float p3[1]={-1.6666658222675323e-1f};         // -11184805*2^-26

asm(R"(
        movss %[x],%%xmm0

movss %[inv_pi],%%xmm1
            movss %[pi0],%%xmm3
        mulss %%xmm0,%%xmm1
        addss %[bias],%%xmm1
        movss %%xmm1,%%xmm2
        subss %[bias],%%xmm1
            movss %[pi1],%%xmm4

mulss %%xmm1,%%xmm3
            pslld $31,%%xmm2
        subss %%xmm3,%%xmm0

mulss %%xmm1,%%xmm4
            movss %[p0],%%xmm6
        subss %%xmm4,%%xmm0

movss %[pi2],%%xmm5
        mulss %%xmm1,%%xmm5
        subss %%xmm5,%%xmm0

mulss %[pi3],%%xmm1
        subss %%xmm1,%%xmm0
        movss %%xmm0,%%xmm7

mulss %%xmm7,%%xmm7
        mulss %%xmm7,%%xmm6
            pxor %%xmm2,%%xmm0
        addss %[p1],%%xmm6
        mulss %%xmm7,%%xmm6
        addss %[p2],%%xmm6
        mulss %%xmm7,%%xmm6
        addss %[p3],%%xmm6
        mulss %%xmm7,%%xmm6
        mulss %%xmm0,%%xmm6
        addss %%xmm6,%%xmm0

movss %%xmm0,%[x]
    )"
        :[x]"+m"(x)

:[inv_pi]"m"(inv_pi[0])
        ,[bias]"m"(bias[0])
        ,[pi0]"m"(pi0[0])
        ,[pi1]"m"(pi1[0])
        ,[pi2]"m"(pi2[0])
        ,[pi3]"m"(pi3[0])
        ,[p0]"m"(p0[0])
        ,[p1]"m"(p1[0])
        ,[p3]"m"(p3[0])
        ,[p2]"m"(p2[0])

:"xmm0","xmm1","xmm2","xmm3","xmm4","xmm5","xmm6","xmm7","cc");
    return x;
}

// |x|<32768:       rel. error < 1.9*10^-7 [ 2 ulp ]
// |x|<131072:      abs. error < 1*10^-6
// Garbage for |x|>8388608.
static inline float cos_f32_sse2(float x)
{
    static const float inv_pi[1]={0.318309886f};                //~ 1/(pi)
    static const float bias[1]={12582912.0f};                   // 1.5*2^23
    static const float half[1]={0.5f};                          // 1/2
    static const float pi0[1]={3.140625f};                      // bits -2..+9
    static const float pi1[1]={9.670257568359375e-4f};          // bits +10..+19
    static const float pi2[1]={6.2771141529083251953125e-7f};   // bits +20..+28
    static const float pi3[1]={1.2154201013012385202950e-10f};  // bits +29..+51
    static const float p0[1]={2.6050554015241175e-6f};          //
    static const float p1[1]={-1.9808833894785493e-4f};         // -13612527*2^-36
    static const float p2[1]={+8.3330469205975532e-3f};         // +8947541*2^-30
    static const float p3[1]={-1.6666658222675323e-1f};         // -11184805*2^-26

asm(R"(
        movss %[x],%%xmm0

movss %[inv_pi],%%xmm1
            movss %[pi0],%%xmm3
        mulss %%xmm0,%%xmm1
        addss %[half],%%xmm1
        addss %[bias],%%xmm1
        movss %%xmm1,%%xmm2
        subss %[bias],%%xmm1
            movss %[pi1],%%xmm4
        subss %[half],%%xmm1

mulss %%xmm1,%%xmm3
            pslld $31,%%xmm2
        subss %%xmm3,%%xmm0

mulss %%xmm1,%%xmm4
            movss %[p0],%%xmm6
        subss %%xmm4,%%xmm0

movss %[pi2],%%xmm5
        mulss %%xmm1,%%xmm5
        subss %%xmm5,%%xmm0

mulss %[pi3],%%xmm1
        subss %%xmm1,%%xmm0
        movss %%xmm0,%%xmm7

movss %%xmm0,%[x]
    )"
        :[x]"+m"(x)

:[inv_pi]"m"(inv_pi[0])
        ,[bias]"m"(bias[0])
        ,[half]"m"(half[0])
        ,[pi0]"m"(pi0[0])
        ,[pi1]"m"(pi1[0])
        ,[pi2]"m"(pi2[0])
        ,[pi3]"m"(pi3[0])
        ,[p0]"m"(p0[0])
        ,[p1]"m"(p1[0])
        ,[p3]"m"(p3[0])
        ,[p2]"m"(p2[0])

:"xmm0","xmm1","xmm2","xmm3","xmm4","xmm5","xmm6","xmm7","cc");
    return x;
}

//==============================================================================

// SSE2 assembly version (stream).

// src and dst may not overlap.
// |x|<65536:       rel. error < 1.9*10^-7 [ 2 ulp ]
// |x|<262144:      abs. error < 4*10^-6
// Garbage for |x|>8388608.
static inline void sin_f32_sse2_stream(const float *src,float *dst,uint32_t n)
{
    static const float inv_pi[4] __attribute__ ((aligned (16))) = // ~ 1/(pi)
        {0.318309886f,0.318309886f,0.318309886f,0.318309886f};
    static const float bias[4] __attribute__ ((aligned (16))) = // 1.5*2^23
        {12582912.0f,12582912.0f,12582912.0f,12582912.0f};
    static const float pi0[4] __attribute__ ((aligned (16))) = // bits -2..+9
        {3.140625f,3.140625f,3.140625f,3.140625f};
    static const float pi1[4] __attribute__ ((aligned (16))) = // bits +10..+19
        {9.670257568359375e-4f,9.670257568359375e-4f,9.670257568359375e-4f,9.670257568359375e-4f};
    static const float pi2[4] __attribute__ ((aligned (16))) = // bits +20..+28
        {6.2771141529083251953125e-7f,6.2771141529083251953125e-7f,6.2771141529083251953125e-7f,6.2771141529083251953125e-7f};
    static const float pi3[4] __attribute__ ((aligned (16))) = // bits +29..+51
        {1.2154201013012385202950e-10f,1.2154201013012385202950e-10f,1.2154201013012385202950e-10f,1.2154201013012385202950e-10f};
    static const float p0[4] __attribute__ ((aligned (16))) = //
        {2.6050554015241175e-6f,2.6050554015241175e-6f,2.6050554015241175e-6f,2.6050554015241175e-6f};
    static const float p1[4] __attribute__ ((aligned (16))) = // -13612527*2^-36
        {-1.9808833894785493e-4f,-1.9808833894785493e-4f,-1.9808833894785493e-4f,-1.9808833894785493e-4f};
    static const float p2[4] __attribute__ ((aligned (16))) = // +8947541*2^-30
        {+8.3330469205975532e-3f,+8.3330469205975532e-3f,+8.3330469205975532e-3f,+8.3330469205975532e-3f};
    static const float p3[4] __attribute__ ((aligned (16))) = // -11184805*2^-26
        {-1.6666658222675323e-1f,-1.6666658222675323e-1f,-1.6666658222675323e-1f,-1.6666658222675323e-1f};

if(n==0) return;
    unsigned int n4=n/4;
    unsigned int nd=n-4*n4;
	if(n4>0)
	{
		asm(R"(
			movl %[n],%%eax
			movl $0,%%ecx
			movl %[x],%%ebx
			movl %[y],%%edx
			shll $4,%%eax
			movaps %[inv_pi],%%xmm4
			movaps %[bias],%%xmm5
		0:
			movups (%%ebx,%%ecx),%%xmm0
			movaps %%xmm0,%%xmm1
			mulps %%xmm4,%%xmm1
			addps %%xmm5,%%xmm1
			movaps %%xmm1,%%xmm2
			subps %%xmm5,%%xmm1
			pslld $31,%%xmm2
			movaps %[pi0],%%xmm3
			mulps %%xmm1,%%xmm3
			subps %%xmm3,%%xmm0
			movaps %[pi1],%%xmm3
			mulps %%xmm1,%%xmm3
			subps %%xmm3,%%xmm0
			movaps %[pi2],%%xmm3
			mulps %%xmm1,%%xmm3
			subps %%xmm3,%%xmm0
			mulps %[pi3],%%xmm1
			subps %%xmm1,%%xmm0
			xorps %%xmm2,%%xmm0
			movups %%xmm0,(%%edx,%%ecx)
			addl $16,%%ecx
			cmpl %%eax,%%ecx
			jnz  0b

movl $0,%%ecx
			movaps %[p0],%%xmm4
			movaps %[p1],%%xmm5
			movaps %[p2],%%xmm6
			movaps %[p3],%%xmm7
		1:
			movups (%%edx,%%ecx),%%xmm0
			movaps %%xmm0,%%xmm1
			mulps %%xmm1,%%xmm1
			movaps %%xmm1,%%xmm2
			mulps %%xmm4,%%xmm2
			addps %%xmm5,%%xmm2
			mulps %%xmm1,%%xmm2
			addps %%xmm6,%%xmm2
			mulps %%xmm1,%%xmm2
			addps %%xmm7,%%xmm2
			mulps %%xmm1,%%xmm2
			mulps %%xmm0,%%xmm2
			addps %%xmm2,%%xmm0
			movups %%xmm0,(%%edx,%%ecx)
			addl $16,%%ecx
			cmpl %%eax,%%ecx
			jnz  1b
		)"  :
			:[n]"m"(n4)
			,[x]"m"(src)
			,[y]"m"(dst)
			,[inv_pi]"m"(inv_pi[0])
			,[bias]"m"(bias[0])
			,[pi0]"m"(pi0[0])
			,[pi1]"m"(pi1[0])
			,[pi2]"m"(pi2[0])
			,[pi3]"m"(pi3[0])
			,[p0]"m"(p0[0])
			,[p1]"m"(p1[0])
			,[p3]"m"(p3[0])
			,[p2]"m"(p2[0])
		   :"eax","ebx","ecx","edx","xmm0","xmm1","xmm2","xmm3","xmm4","xmm5","xmm6","xmm7","memory","cc");
	}
	for(uint32_t i=0;i<nd;++i)
        dst[4*n4+i]=sin_f32_sse2(src[4*n4+i]);
}

// src and dst may not overlap.
// |x|<32768:       rel. error < 1.9*10^-7 [ 2 ulp ]
// |x|<131072:      abs. error < 1*10^-6
// Garbage for |x|>8388608.
static inline void cos_f32_sse2_stream(const float *src,float *dst,uint32_t n)
{
    static const float inv_pi[4] __attribute__ ((aligned (16))) = // ~ 1/(pi)
        {0.318309886f,0.318309886f,0.318309886f,0.318309886f};
    static const float bias[4] __attribute__ ((aligned (16))) = // 1.5*2^23
        {12582912.0f,12582912.0f,12582912.0f,12582912.0f};
    static const float half[4] __attribute__ ((aligned (16))) = // 1/2
        {0.5f,0.5f,0.5f,0.5f};
    static const float pi0[4] __attribute__ ((aligned (16))) = // bits -2..+9
        {3.140625f,3.140625f,3.140625f,3.140625f};
    static const float pi1[4] __attribute__ ((aligned (16))) = // bits +10..+19
        {9.670257568359375e-4f,9.670257568359375e-4f,9.670257568359375e-4f,9.670257568359375e-4f};
    static const float pi2[4] __attribute__ ((aligned (16))) = // bits +20..+28
        {6.2771141529083251953125e-7f,6.2771141529083251953125e-7f,6.2771141529083251953125e-7f,6.2771141529083251953125e-7f};
    static const float pi3[4] __attribute__ ((aligned (16))) = // bits +29..+51
        {1.2154201013012385202950e-10f,1.2154201013012385202950e-10f,1.2154201013012385202950e-10f,1.2154201013012385202950e-10f};
    static const float p0[4] __attribute__ ((aligned (16))) = //
        {2.6050554015241175e-6f,2.6050554015241175e-6f,2.6050554015241175e-6f,2.6050554015241175e-6f};
    static const float p1[4] __attribute__ ((aligned (16))) = // -13612527*2^-36
        {-1.9808833894785493e-4f,-1.9808833894785493e-4f,-1.9808833894785493e-4f,-1.9808833894785493e-4f};
    static const float p2[4] __attribute__ ((aligned (16))) = // +8947541*2^-30
        {+8.3330469205975532e-3f,+8.3330469205975532e-3f,+8.3330469205975532e-3f,+8.3330469205975532e-3f};
    static const float p3[4] __attribute__ ((aligned (16))) = // -11184805*2^-26
        {-1.6666658222675323e-1f,-1.6666658222675323e-1f,-1.6666658222675323e-1f,-1.6666658222675323e-1f};

if(n==0) return;
    unsigned int n4=n/4;
    unsigned int nd=n-4*n4;
	if(n4>0)
	{
		asm(R"(
			movl %[n],%%eax
			movl $0,%%ecx
			movl %[x],%%ebx
			movl %[y],%%edx
			shll $4,%%eax
			movaps %[inv_pi],%%xmm4
			movaps %[bias],%%xmm5
            movaps %[half],%%xmm6
		0:
			movups (%%ebx,%%ecx),%%xmm0
			movaps %%xmm0,%%xmm1
			mulps %%xmm4,%%xmm1
            addps %%xmm6,%%xmm1
			addps %%xmm5,%%xmm1
			movaps %%xmm1,%%xmm2
			subps %%xmm5,%%xmm1
            subps %%xmm6,%%xmm1
			pslld $31,%%xmm2
			movaps %[pi0],%%xmm3
			mulps %%xmm1,%%xmm3
			subps %%xmm3,%%xmm0
			movaps %[pi1],%%xmm3
			mulps %%xmm1,%%xmm3
			subps %%xmm3,%%xmm0
			movaps %[pi2],%%xmm3
			mulps %%xmm1,%%xmm3
			subps %%xmm3,%%xmm0
			mulps %[pi3],%%xmm1
			subps %%xmm1,%%xmm0
			xorps %%xmm2,%%xmm0
			movups %%xmm0,(%%edx,%%ecx)
			addl $16,%%ecx
			cmpl %%eax,%%ecx
			jnz  0b

movl $0,%%ecx
			movaps %[p0],%%xmm4
			movaps %[p1],%%xmm5
			movaps %[p2],%%xmm6
			movaps %[p3],%%xmm7
		1:
			movups (%%edx,%%ecx),%%xmm0
			movaps %%xmm0,%%xmm1
			mulps %%xmm1,%%xmm1
			movaps %%xmm1,%%xmm2
			mulps %%xmm4,%%xmm2
			addps %%xmm5,%%xmm2
			mulps %%xmm1,%%xmm2
			addps %%xmm6,%%xmm2
			mulps %%xmm1,%%xmm2
			addps %%xmm7,%%xmm2
			mulps %%xmm1,%%xmm2
			mulps %%xmm0,%%xmm2
			addps %%xmm2,%%xmm0
			movups %%xmm0,(%%edx,%%ecx)
			addl $16,%%ecx
			cmpl %%eax,%%ecx
			jnz  1b
		)"  :
			:[n]"m"(n4)
			,[x]"m"(src)
			,[y]"m"(dst)
			,[inv_pi]"m"(inv_pi[0])
			,[bias]"m"(bias[0])
            ,[half]"m"(half[0])
			,[pi0]"m"(pi0[0])
			,[pi1]"m"(pi1[0])
			,[pi2]"m"(pi2[0])
			,[pi3]"m"(pi3[0])
			,[p0]"m"(p0[0])
			,[p1]"m"(p1[0])
			,[p3]"m"(p3[0])
			,[p2]"m"(p2[0])
		   :"eax","ebx","ecx","edx","xmm0","xmm1","xmm2","xmm3","xmm4","xmm5","xmm6","xmm7","memory","cc");
	}
	for(uint32_t i=0;i<nd;++i)
        dst[4*n4+i]=cos_f32_sse2(src[4*n4+i]);
}

//==============================================================================

// Intrinsic versions.

static const __m128 c0=_mm_set_ps1(0.318309886f);
static const __m128 c1=_mm_set_ps1(12582912.0f);
static const __m128 c2=_mm_set_ps1(3.140625f);
static const __m128 c3=_mm_set_ps1(9.670257568359375e-4f);
static const __m128 c4=_mm_set_ps1(6.2771141529083251953125e-7f);
static const __m128 c5=_mm_set_ps1(1.2154201013012385202950e-10f);
static const __m128 c6=_mm_set_ps1(+2.6050554015241175e-6f);
static const __m128 c7=_mm_set_ps1(-1.9808833894785493e-4f);
static const __m128 c8=_mm_set_ps1(+8.3330469205975532e-3f);
static const __m128 c9=_mm_set_ps1(-1.6666658222675323e-1f);
static const __m128 c10=_mm_set_ps1(0.5f);

// |x|<65536:       rel. error < 1.9*10^-7 [ 2 ulp ]
// |x|<262144:      abs. error < 4*10^-6
// Garbage for |x|>8388608.
inline float sin_f32_sse2_intrin(float x)
{
	__m128 x_sse=_mm_load_ss(&x);
	__m128 y=_mm_mul_ss(c0,x_sse);

// Round y to closest integer.
	y=_mm_add_ss(y,c1);
	__m128 sign=y;
	y=_mm_sub_ss(y,c1);

sign=_mm_castsi128_ps(_mm_slli_epi32(_mm_castps_si128(sign),31));

x_sse=_mm_sub_ss(x_sse,_mm_mul_ss(c2,y));
	x_sse=_mm_sub_ss(x_sse,_mm_mul_ss(c3,y));
	x_sse=_mm_sub_ss(x_sse,_mm_mul_ss(c4,y));
	x_sse=_mm_sub_ss(x_sse,_mm_mul_ss(c5,y));
	
    __m128 x2=_mm_mul_ss(x_sse,x_sse);

__m128 ret_sse=_mm_mul_ss(c6,x2);
	ret_sse=_mm_mul_ss(_mm_add_ss(c7,ret_sse),x2);
	ret_sse=_mm_mul_ss(_mm_add_ss(c8,ret_sse),x2);
    x_sse=_mm_xor_ps(x_sse,sign);
	ret_sse=_mm_mul_ss(_mm_add_ss(c9,ret_sse),x2);
	ret_sse=_mm_add_ss(x_sse,_mm_mul_ss(ret_sse,x_sse));
	float ret;
	
	_mm_store_ss(&ret,ret_sse);
    return ret;
}

// |x|<32768:       rel. error < 1.9*10^-7 [ 2 ulp ]
// |x|<131072:      abs. error < 1*10^-6
// Garbage for |x|>8388608.
inline float cos_f32_sse2_intrin(float x)
{
	__m128 x_sse=_mm_load_ss(&x);
	__m128 y=_mm_mul_ss(c0,x_sse);

// Round y to closest integer.
    y=_mm_add_ss(y,c10);
	y=_mm_add_ss(y,c1);
	__m128 sign=y;
	y=_mm_sub_ss(y,c1);
    y=_mm_sub_ss(y,c10);

sign=_mm_castsi128_ps(_mm_slli_epi32(_mm_castps_si128(sign),31));

// src and dst may not overlap.
// |x|<65536:       rel. error < 1.9*10^-7 [ 2 ulp ]
// |x|<262144:      abs. error < 4*10^-6
// Garbage for |x|>8388608.
inline void sin_f32_sse2_intrin_stream(const float *src,float *dst,uint32_t n)
{
    if(n==0) return;
    uint32_t n4=n/4;
    uint32_t nd=n-4*n4;
    for(uint32_t i=0;i<n4;++i)
    {
        __m128 x=_mm_loadu_ps(src+4*i);
        __m128 y=_mm_mul_ps(c0,x);
        y=_mm_add_ps(y,c1);
        __m128 sign=y;
        y=_mm_sub_ps(y,c1);
        x=_mm_sub_ps(x,_mm_mul_ps(c2,y));
        x=_mm_sub_ps(x,_mm_mul_ps(c3,y));
        sign=_mm_castsi128_ps(_mm_slli_epi32(_mm_castps_si128(sign),31));
        x=_mm_sub_ps(x,_mm_mul_ps(c4,y));
        x=_mm_sub_ps(x,_mm_mul_ps(c5,y));
        x=_mm_xor_ps(x,sign);
        _mm_storeu_ps(dst+4*i,x);
    }
    for(uint32_t i=0;i<n4;++i)
    {
        __m128 x=_mm_loadu_ps(dst+4*i);
        __m128 x2=_mm_mul_ps(x,x);
        __m128 y=_mm_mul_ps(c6,x2);
        y=_mm_mul_ps(_mm_add_ps(c7,y),x2);
        y=_mm_mul_ps(_mm_add_ps(c8,y),x2);
    	y=_mm_mul_ps(_mm_add_ps(c9,y),x2);
        y=_mm_add_ps(x,_mm_mul_ps(y,x));
        _mm_storeu_ps(dst+4*i,y);
    }
    for(uint32_t i=0;i<nd;++i)
        dst[4*n4+i]=sin_f32_sse2_intrin(src[4*n4+i]);
}

// src and dst may not overlap.
// |x|<32768:       rel. error < 1.9*10^-7 [ 2 ulp ]
// |x|<131072:      abs. error < 1*10^-6
// Garbage for |x|>8388608.
inline void cos_f32_sse2_intrin_stream(const float *src,float *dst,uint32_t n)
{
    if(n==0) return;
    uint32_t n4=n/4;
    uint32_t nd=n-4*n4;
    for(uint32_t i=0;i<n4;++i)
    {
        __m128 x=_mm_loadu_ps(src+4*i);
        __m128 y=_mm_mul_ps(c0,x);
        y=_mm_add_ps(y,c10);
        y=_mm_add_ps(y,c1);
        __m128 sign=y;
        y=_mm_sub_ps(y,c1);
        y=_mm_sub_ps(y,c10);
        x=_mm_sub_ps(x,_mm_mul_ps(c2,y));
        x=_mm_sub_ps(x,_mm_mul_ps(c3,y));
        sign=_mm_castsi128_ps(_mm_slli_epi32(_mm_castps_si128(sign),31));
        x=_mm_sub_ps(x,_mm_mul_ps(c4,y));
        x=_mm_sub_ps(x,_mm_mul_ps(c5,y));
        x=_mm_xor_ps(x,sign);
        _mm_storeu_ps(dst+4*i,x);
    }
    for(uint32_t i=0;i<n4;++i)
    {
        __m128 x=_mm_loadu_ps(dst+4*i);
        __m128 x2=_mm_mul_ps(x,x);
        __m128 y=_mm_mul_ps(c6,x2);
        y=_mm_mul_ps(_mm_add_ps(c7,y),x2);
        y=_mm_mul_ps(_mm_add_ps(c8,y),x2);
    	y=_mm_mul_ps(_mm_add_ps(c9,y),x2);
        y=_mm_add_ps(x,_mm_mul_ps(y,x));
        _mm_storeu_ps(dst+4*i,y);
    }
    for(uint32_t i=0;i<nd;++i)
        dst[4*n4+i]=cos_f32_sse2_intrin(src[4*n4+i]);
}

//==============================================================================

static float test(void (*f)(const float*,float*,uint32_t),double (*ref)(double))
{
    static const int M=1024+3;
    static const int N=100*1024;
    static float x[M];
    static float y[M];
    for(int i=0;i<M;++i)
        x[i]=float(i);
    double t=double(clock());
    for(int i=0;i<N;++i)
    {
        f(x,y,M);
    }
    t=double(clock())-t;
    float b=0.0f;
    for(int i=0;i<M;++i)
    {
        if(fabs(ref(x[i])-y[i])>b) b=fabs(ref(x[i])-y[i]);
    }
    printf("Max. abs. error: %.7g\n",b);
    printf("%.1f ns/value\n",1.0e+9*t/(double(N)*double(M)*double(CLOCKS_PER_SEC)));

}

int main()
{
    printf("sin_f32_sse2_stream\n");
    test(sin_f32_sse2_stream,sin);
    printf("cos_f32_sse2_stream\n");
    test(cos_f32_sse2_stream,cos);
    printf("sin_f32_sse2_intrin_stream\n");
    test(sin_f32_sse2_intrin_stream,sin);
    printf("cos_f32_sse2_intrin_stream\n");
    test(cos_f32_sse2_intrin_stream,cos);
    return 0;
}

Show compiler warnings [ + ] Compiler args [ + ] Show input
Compilation time: 0.33 sec, absolute running time: 0.85 sec, cpu time: 0.7 sec, memory peak: 4 Mb, absolute service time: 1.19 sec