diff --git a/CMakeLists.txt b/CMakeLists.txt
index d8c292ea..fe4f39c9 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -168,15 +168,12 @@ target_link_libraries(bm_cum_sum benchmark)
 add_executable(bm_write benchmarks/BM_writeInterleaved.cpp sources/SIMDSSE.cpp)
 target_link_libraries(bm_write benchmark absl::span)
 
+add_executable(bm_read benchmarks/BM_readInterleaved.cpp sources/SIMDSSE.cpp)
+target_link_libraries(bm_read benchmark absl::span)
+
 add_executable(bm_fill benchmarks/BM_fill.cpp sources/SIMDSSE.cpp)
 target_link_libraries(bm_fill benchmark absl::span)
 
-add_executable(bm_math_functions benchmarks/Math_functions.cpp)
-target_link_libraries(bm_math_functions benchmark)
-if (UNIX)
-    target_compile_options(bm_math_functions PRIVATE -funsafe-math-optimizations -msse4.2)
-endif()
-
 add_executable(bm_mathfuns benchmarks/BM_mathfuns.cpp sources/SIMDSSE.cpp)
 target_link_libraries(bm_mathfuns benchmark absl::span absl::span)
 if (UNIX)
diff --git a/benchmarks/BM_mathfuns.cpp b/benchmarks/BM_mathfuns.cpp
index 97bcd987..28e8fe39 100644
--- a/benchmarks/BM_mathfuns.cpp
+++ b/benchmarks/BM_mathfuns.cpp
@@ -6,13 +6,26 @@
 #include <iostream>
 #include "../sources/SIMDHelpers.h"
 
-static void Dummy(benchmark::State& state) {
+class MyFixture : public benchmark::Fixture {
+public:
+  void SetUp(const ::benchmark::State& state) {
     std::random_device rd { };
     std::mt19937 gen { rd() };
-    std::vector<float> source(state.range(0));
-    std::vector<float> result(state.range(0));
-    std::normal_distribution<float> dist { };
+    std::uniform_real_distribution<float> dist { 0.1, 1 };
+    source = std::vector<float>(state.range(0));
+    result = std::vector<float>(state.range(0));
     std::generate(source.begin(), source.end(), [&]() { return dist(gen); });
+  }
+
+  void TearDown(const ::benchmark::State& state [[maybe_unused]]) {
+  }
+
+  std::vector<float> source;
+  std::vector<float> result;
+};
+
+
+BENCHMARK_DEFINE_F(MyFixture, Dummy)(benchmark::State& state) {
     for (auto _ : state)
     {
          for (int i = 0; i < state.range(0); ++i)
@@ -21,13 +34,7 @@ static void Dummy(benchmark::State& state) {
     }
 }
 
-static void StdExp(benchmark::State& state) {
-    std::random_device rd { };
-    std::mt19937 gen { rd() };
-    std::vector<float> source(state.range(0));
-    std::vector<float> result(state.range(0));
-    std::normal_distribution<float> dist { };
-    std::generate(source.begin(), source.end(), [&]() { return dist(gen); });
+BENCHMARK_DEFINE_F(MyFixture, StdExp)(benchmark::State& state) {
     for (auto _ : state)
     {
         for (int i = 0; i < state.range(0); ++i)
@@ -36,29 +43,7 @@ static void StdExp(benchmark::State& state) {
     }
 }
 
-// static void StdExpOMP(benchmark::State& state) {
-//     std::random_device rd { };
-//     std::mt19937 gen { rd() };
-//     std::vector<float> source(state.range(0));
-//     std::vector<float> result(state.range(0));
-//     std::normal_distribution<float> dist { };
-//     std::generate(source.begin(), source.end(), [&]() { return dist(gen); });
-//     for (auto _ : state)
-//     {
-//         #pragma omp simd
-//         for (int i = 0; i < state.range(0); ++i)
-//             result[i] = std::exp(source[i]);
-//         benchmark::DoNotOptimize(result);
-//     }
-// }
-
-static void Scalar(benchmark::State& state) {
-    std::random_device rd { };
-    std::mt19937 gen { rd() };
-    std::vector<float> source(state.range(0));
-    std::vector<float> result(state.range(0));
-    std::normal_distribution<float> dist { };
-    std::generate(source.begin(), source.end(), [&]() { return dist(gen); });
+BENCHMARK_DEFINE_F(MyFixture, ScalarExp)(benchmark::State& state) {
     for (auto _ : state)
     {
         exp<float, false>(source, absl::MakeSpan(result));
@@ -66,13 +51,7 @@ static void Scalar(benchmark::State& state) {
     }
 }
 
-static void SIMD(benchmark::State& state) {
-    std::random_device rd { };
-    std::mt19937 gen { rd() };
-    std::vector<float> source(state.range(0));
-    std::vector<float> result(state.range(0));
-    std::normal_distribution<float> dist { };
-    std::generate(source.begin(), source.end(), [&]() { return dist(gen); });
+BENCHMARK_DEFINE_F(MyFixture, SIMDExp)(benchmark::State& state) {
     for (auto _ : state)
     {
         exp<float, true>(source, absl::MakeSpan(result));
@@ -80,8 +59,93 @@ static void SIMD(benchmark::State& state) {
     }
 }
 
-BENCHMARK(Dummy)->RangeMultiplier(2)->Range(1 << 6, 1 << 10);
-BENCHMARK(StdExp)->RangeMultiplier(2)->Range(1 << 6, 1 << 10);
-BENCHMARK(Scalar)->RangeMultiplier(2)->Range(1 << 6, 1 << 10);
-BENCHMARK(SIMD)->RangeMultiplier(2)->Range(1 << 6, 1 << 10);
+BENCHMARK_DEFINE_F(MyFixture, StdLog)(benchmark::State& state) {
+    for (auto _ : state)
+    {
+        for (int i = 0; i < state.range(0); ++i)
+            result[i] = std::log(source[i]);
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK_DEFINE_F(MyFixture, ScalarLog)(benchmark::State& state) {
+    for (auto _ : state)
+    {
+        log<float, false>(source, absl::MakeSpan(result));
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK_DEFINE_F(MyFixture, SIMDLog)(benchmark::State& state) {
+    for (auto _ : state)
+    {
+        log<float, true>(source, absl::MakeSpan(result));
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK_DEFINE_F(MyFixture, StdSin)(benchmark::State& state) {
+    for (auto _ : state)
+    {
+        for (int i = 0; i < state.range(0); ++i)
+            result[i] = std::sin(source[i]);
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK_DEFINE_F(MyFixture, ScalarSin)(benchmark::State& state) {
+    for (auto _ : state)
+    {
+        sin<float, false>(source, absl::MakeSpan(result));
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK_DEFINE_F(MyFixture, SIMDSin)(benchmark::State& state) {
+    for (auto _ : state)
+    {
+        sin<float, true>(source, absl::MakeSpan(result));
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK_DEFINE_F(MyFixture, StdCos)(benchmark::State& state) {
+    for (auto _ : state)
+    {
+        for (int i = 0; i < state.range(0); ++i)
+            result[i] = std::cos(source[i]);
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK_DEFINE_F(MyFixture, ScalarCos)(benchmark::State& state) {
+    for (auto _ : state)
+    {
+        cos<float, false>(source, absl::MakeSpan(result));
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK_DEFINE_F(MyFixture, SIMDCos)(benchmark::State& state) {
+    for (auto _ : state)
+    {
+        cos<float, true>(source, absl::MakeSpan(result));
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK_REGISTER_F(MyFixture, Dummy)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, StdExp)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, ScalarExp)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, SIMDExp)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, StdLog)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, ScalarLog)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, SIMDLog)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, StdSin)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, ScalarSin)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, SIMDSin)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, StdCos)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, ScalarCos)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+BENCHMARK_REGISTER_F(MyFixture, SIMDCos)->RangeMultiplier(4)->Range(1 << 6, 1 << 10);
+
 BENCHMARK_MAIN();
\ No newline at end of file
diff --git a/benchmarks/BM_readInterleaved.cpp b/benchmarks/BM_readInterleaved.cpp
new file mode 100644
index 00000000..6d39d7ae
--- /dev/null
+++ b/benchmarks/BM_readInterleaved.cpp
@@ -0,0 +1,54 @@
+#include <benchmark/benchmark.h>
+#include "../sources/SIMDHelpers.h"
+#include "../sources/Buffer.h"
+#include <algorithm>
+#include <numeric>
+#include <absl/types/span.h>
+
+static void Scalar(benchmark::State& state) {
+  Buffer<float> input (state.range(0) * 2);
+  Buffer<float> outputLeft (state.range(0));
+  Buffer<float> outputRight (state.range(0));
+  std::iota(input.begin(), input.end(), 1.0f);
+
+  for (auto _ : state) {
+    readInterleaved<float, false>(input, absl::MakeSpan(outputLeft), absl::MakeSpan(outputRight));
+  }
+}
+
+static void SSE(benchmark::State& state) {
+  Buffer<float> input (state.range(0) * 2);
+  Buffer<float> outputLeft (state.range(0));
+  Buffer<float> outputRight (state.range(0));
+  std::iota(input.begin(), input.end(), 1.0f);
+
+  for (auto _ : state) {
+    readInterleaved<float, true>(input, absl::MakeSpan(outputLeft), absl::MakeSpan(outputRight));
+  }
+}
+
+static void Scalar_Unaligned(benchmark::State& state) {
+  Buffer<float> input (state.range(0) * 2);
+  Buffer<float> outputLeft (state.range(0));
+  Buffer<float> outputRight (state.range(0));
+  std::iota(input.begin(), input.end(), 1.0f);
+  for (auto _ : state) {
+    readInterleaved<float, false>(absl::MakeSpan(input).subspan(1), absl::MakeSpan(outputLeft).subspan(1), absl::MakeSpan(outputRight).subspan(1));
+  }
+}
+
+static void SSE_Unaligned(benchmark::State& state) {
+  Buffer<float> input (state.range(0) * 2);
+  Buffer<float> outputLeft (state.range(0));
+  Buffer<float> outputRight (state.range(0));
+  std::iota(input.begin(), input.end(), 1.0f);
+  for (auto _ : state) {
+    readInterleaved<float, true>(absl::MakeSpan(input).subspan(1), absl::MakeSpan(outputLeft).subspan(1), absl::MakeSpan(outputRight).subspan(1));
+  }
+}
+
+BENCHMARK(Scalar)->Range((8<<10), (8<<20));
+BENCHMARK(SSE)->Range((8<<10), (8<<20));
+BENCHMARK(Scalar_Unaligned)->Range((8<<10), (8<<20));
+BENCHMARK(SSE_Unaligned)->Range((8<<10), (8<<20));
+BENCHMARK_MAIN();
\ No newline at end of file
diff --git a/benchmarks/BM_writeInterleaved.cpp b/benchmarks/BM_writeInterleaved.cpp
index ee32a78f..d80cca04 100644
--- a/benchmarks/BM_writeInterleaved.cpp
+++ b/benchmarks/BM_writeInterleaved.cpp
@@ -3,6 +3,7 @@
 #include "../sources/Buffer.h"
 #include <algorithm>
 #include <numeric>
+#include <absl/types/span.h>
 
 static void Interleaved_Write(benchmark::State& state) {
   Buffer<float> inputLeft (state.range(0));
@@ -12,7 +13,7 @@ static void Interleaved_Write(benchmark::State& state) {
   std::iota(inputRight.begin(), inputRight.end(), 1.0f);
 
   for (auto _ : state) {
-    writeInterleaved<float, false>(inputLeft, inputRight, output);
+    writeInterleaved<float, false>(inputLeft, inputRight, absl::MakeSpan(output));
   }
 }
 
@@ -23,7 +24,7 @@ static void Interleaved_Write_SSE(benchmark::State& state) {
   std::iota(inputLeft.begin(), inputLeft.end(), 1.0f);
   std::iota(inputRight.begin(), inputRight.end(), 1.0f);
   for (auto _ : state) {
-    writeInterleaved<float, true>(inputLeft, inputRight, output);
+    writeInterleaved<float, true>(inputLeft, inputRight, absl::MakeSpan(output));
     benchmark::DoNotOptimize(output);
   }
 }
@@ -35,7 +36,7 @@ static void Unaligned_Interleaved_Write(benchmark::State& state) {
   std::iota(inputLeft.begin(), inputLeft.end(), 1.0f);
   std::iota(inputRight.begin(), inputRight.end(), 1.0f);
   for (auto _ : state) {
-    writeInterleaved<float, false>(gsl::span(inputLeft).subspan(1) , gsl::span(inputRight).subspan(1), gsl::span(output).subspan(1));
+    writeInterleaved<float, false>(absl::MakeSpan(inputLeft).subspan(1) , absl::MakeSpan(inputRight).subspan(1), absl::MakeSpan(output).subspan(1));
     benchmark::DoNotOptimize(output);
   }
 }
@@ -47,7 +48,7 @@ static void Unaligned_Interleaved_Write_SSE(benchmark::State& state) {
   std::iota(inputLeft.begin(), inputLeft.end(), 1.0f);
   std::iota(inputRight.begin(), inputRight.end(), 1.0f);
   for (auto _ : state) {
-    writeInterleaved<float, true>(gsl::span(inputLeft).subspan(1) , gsl::span(inputRight).subspan(1), gsl::span(output).subspan(1));
+    writeInterleaved<float, true>(absl::MakeSpan(inputLeft).subspan(1) , absl::MakeSpan(inputRight).subspan(1), absl::MakeSpan(output).subspan(1));
     benchmark::DoNotOptimize(output);
   }
 }
diff --git a/benchmarks/Math_functions.cpp b/benchmarks/Math_functions.cpp
deleted file mode 100644
index eabffa80..00000000
--- a/benchmarks/Math_functions.cpp
+++ /dev/null
@@ -1,111 +0,0 @@
-#include <benchmark/benchmark.h>
-#include <random>
-#include <cmath>
-#include <iostream>
-
-/*
-Cephes Math Library Release 2.2:  June, 1992
-Copyright 1984, 1987, 1989 by Stephen L. Moshier
-Direct inquiries to 30 Frost Street, Cambridge, MA 02140
-*/
-
-/* Single precision exponential function.
- * test interval: [-0.5, +0.5]
- * trials: 80000
- * peak relative error: 7.6e-8
- * rms relative error: 2.8e-8
- */
-
-static float MAXNUMF = 3.4028234663852885981170418348451692544e38;
-static float MAXLOGF = 88.72283905206835;
-static float MINLOGF = -103.278929903431851103; /* log(2^-149) */
-static float LOG2EF = 1.44269504088896341;
-static float C1 =   0.693359375;
-static float C2 =  -2.12194440e-4;
-
-float cephes_expf(float xx) {
-    float x, z;
-    int n;
-
-    x = xx;
-
-
-    if( x > MAXLOGF)
-        {
-        //mtherr( "expf", OVERFLOW );
-        return( MAXNUMF );
-        }
-
-    if( x < MINLOGF )
-        {
-        //mtherr( "expf", UNDERFLOW );
-        return(0.0);
-        }
-
-    /* Express e**x = e**g 2**n
-    *   = e**g e**( n loge(2) )
-    *   = e**( g + n loge(2) )
-    */
-    z = floorf( LOG2EF * x + 0.5 ); /* floor() truncates toward -infinity. */
-
-    x -= z * C1;
-    x -= z * C2;
-    n = z;
-
-    z = x * x;
-    /* Theoretical peak relative error in [-0.5, +0.5] is 4.2e-9. */
-    z =
-    ((((( 1.9875691500E-4f  * x
-    + 1.3981999507E-3f) * x
-    + 8.3334519073E-3f) * x
-    + 4.1665795894E-2f) * x
-    + 1.6666665459E-1f) * x
-    + 5.0000001201E-1f) * z
-    + x
-    + 1.0;
-
-    /* multiply by power of 2 */
-    x = ldexpf( z, n );
-
-    return( x );
-}
-
-static void Dummy(benchmark::State& state) {
-    std::random_device rd { };
-    std::mt19937 gen { rd() };
-    std::normal_distribution<float> dist { };
-
-    for (auto _ : state)
-    {
-        auto value = dist(gen);
-        benchmark::DoNotOptimize(value);
-    }
-}
-
-static void StdExp(benchmark::State& state) {
-    std::random_device rd { };
-    std::mt19937 gen { rd() };
-    std::normal_distribution<float> dist { };
-    for (auto _ : state)
-    {
-        auto value = std::exp(dist(gen));
-        benchmark::DoNotOptimize(value);
-    }
-}
-
-static void CephesExp(benchmark::State& state) {
-    std::random_device rd { };
-    std::mt19937 gen { rd() };
-    std::normal_distribution<float> dist { };
-    for (auto _ : state)
-    {
-        auto value = cephes_expf(dist(gen));
-        benchmark::DoNotOptimize(value);
-    }
-}
-
-
-BENCHMARK(Dummy);
-BENCHMARK(StdExp);
-BENCHMARK(CephesExp);
-BENCHMARK_MAIN();
\ No newline at end of file
diff --git a/sources/SIMDDummy.cpp b/sources/SIMDDummy.cpp
index 768add64..b63aef45 100644
--- a/sources/SIMDDummy.cpp
+++ b/sources/SIMDDummy.cpp
@@ -38,4 +38,22 @@ template<>
 void exp<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept
 {
     exp<float, false>(input, output);
+}
+
+template<>
+void log<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept
+{
+    log<float, false>(input, output);
+}
+
+template<>
+void sin<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept
+{
+    sin<float, false>(input, output);
+}
+
+template<>
+void cos<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept
+{
+    cos<float, false>(input, output);
 }
\ No newline at end of file
diff --git a/sources/SIMDHelpers.h b/sources/SIMDHelpers.h
index a823c74d..16c6505f 100644
--- a/sources/SIMDHelpers.h
+++ b/sources/SIMDHelpers.h
@@ -64,6 +64,42 @@ void exp(absl::Span<const Type> input, absl::Span<Type> output) noexcept
 template<>
 void exp<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept;
 
+template<class Type, bool SIMD=SIMDConfig::useSIMD>
+void log(absl::Span<const Type> input, absl::Span<Type> output) noexcept
+{
+    ASSERT(output.size() >= input.size());
+    auto sentinel = std::min(input.size(), output.size());
+    for (decltype(sentinel) i = 0; i < sentinel; ++i)
+        output[i] = std::log(input[i]);
+}
+
+template<>
+void log<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept;
+
+template<class Type, bool SIMD=SIMDConfig::useSIMD>
+void sin(absl::Span<const Type> input, absl::Span<Type> output) noexcept
+{
+    ASSERT(output.size() >= input.size());
+    auto sentinel = std::min(input.size(), output.size());
+    for (decltype(sentinel) i = 0; i < sentinel; ++i)
+        output[i] = std::sin(input[i]);
+}
+
+template<>
+void sin<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept;
+
+template<class Type, bool SIMD=SIMDConfig::useSIMD>
+void cos(absl::Span<const Type> input, absl::Span<Type> output) noexcept
+{
+    ASSERT(output.size() >= input.size());
+    auto sentinel = std::min(input.size(), output.size());
+    for (decltype(sentinel) i = 0; i < sentinel; ++i)
+        output[i] = std::cos(input[i]);
+}
+
+template<>
+void cos<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept;
+
 template<class T, bool SIMD=SIMDConfig::useSIMD>
 void loopingSFZIndex(absl::Span<const T> inputLeft, absl::Span<const T> inputRight, absl::Span<T> output);
 
diff --git a/sources/SIMDSSE.cpp b/sources/SIMDSSE.cpp
index 4dcc4581..850b2c57 100644
--- a/sources/SIMDSSE.cpp
+++ b/sources/SIMDSSE.cpp
@@ -120,4 +120,49 @@ void exp<float, true>(absl::Span<const float> input, absl::Span<float> output) n
         out += 4;
         in += 4;
     }
+}
+
+template<>
+void cos<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept
+{
+    ASSERT(output.size() >= input.size());
+    auto* in = input.begin();
+    auto* out = output.begin();
+    auto* sentinel = in + std::min(input.size(), output.size());
+    while (in < sentinel)
+    {
+        _mm_storeu_ps(out, cos_ps(_mm_loadu_ps(in)));
+        out += 4;
+        in += 4;
+    }
+}
+
+template<>
+void log<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept
+{
+    ASSERT(output.size() >= input.size());
+    auto* in = input.begin();
+    auto* out = output.begin();
+    auto* sentinel = in + std::min(input.size(), output.size());
+    while (in < sentinel)
+    {
+        _mm_storeu_ps(out, log_ps(_mm_loadu_ps(in)));
+        out += 4;
+        in += 4;
+    }
+}
+
+template<>
+void sin<float, true>(absl::Span<const float> input, absl::Span<float> output) noexcept
+{
+    ASSERT(output.size() >= input.size());
+    auto* in = input.begin();
+    auto* out = output.begin();
+    auto* sentinel = in + std::min(input.size(), output.size());
+    while (in < sentinel)
+    {
+        _mm_storeu_ps(out, sin_ps(_mm_loadu_ps(in)));
+        out += 4;
+        in += 4;
+    }
 }
\ No newline at end of file