{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2efd7be2-5ea7-42cb-a2e2-9ddb409bc893",
   "metadata": {},
   "outputs": [],
   "source": [
    "#include <cstddef>\n",
    "#include <vector>\n",
    "#include <random>\n",
    "#include <iostream>\n",
    "\n",
    "#include \"xsimd/xsimd.hpp\""
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6f498bff-4366-4e1c-ab25-3d3589740e78",
   "metadata": {},
   "source": [
    "Comparing two implementation of element wise mean."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cf6b2c85-1427-49a3-8bef-0340601319e8",
   "metadata": {},
   "outputs": [],
   "source": [
    "std::vector<double> random_vector(std::size_t n)\n",
    "{\n",
    "    static std::mt19937 gen(42);\n",
    "    std::uniform_real_distribution<double> dist(0.0, 1.0);\n",
    "    std::vector<double> v(n);\n",
    "    for (auto& x : v)\n",
    "    {\n",
    "        x = dist(gen);\n",
    "    }\n",
    "    return v;\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a2ab6177-60e2-4d65-b323-51a6a21ad530",
   "metadata": {},
   "outputs": [],
   "source": [
    "std::vector<double> mean_scalar(const std::vector<double>& a, const std::vector<double>& b)\n",
    "{\n",
    "    std::size_t size = a.size();\n",
    "    std::vector<double> res(size);\n",
    "    for (std::size_t i = 0; i < size; ++i)\n",
    "    {\n",
    "        res[i] = (a[i] + b[i]) / 2;\n",
    "    }\n",
    "    return res;\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "082ed2cd-273a-4c3c-a206-03b1efe49482",
   "metadata": {},
   "outputs": [],
   "source": [
    "std::vector<double> mean_simd(const std::vector<double>& a, const std::vector<double>& b)\n",
    "{\n",
    "    using b_type = xsimd::batch<double>;\n",
    "    std::size_t inc = b_type::size;\n",
    "    std::size_t size = a.size();\n",
    "    std::vector<double> res(size);\n",
    "\n",
    "    // size for which the vectorization is possible\n",
    "    std::size_t vec_size = size - size % inc;\n",
    "    for (std::size_t i = 0; i < vec_size; i += inc)\n",
    "    {\n",
    "        b_type avec = b_type::load_unaligned(&a[i]);\n",
    "        b_type bvec = b_type::load_unaligned(&b[i]);\n",
    "        b_type rvec = (avec + bvec) / 2;\n",
    "        rvec.store_unaligned(&res[i]);\n",
    "    }\n",
    "    // Remaining part that cannot be vectorize\n",
    "    for (std::size_t i = vec_size; i < size; ++i)\n",
    "    {\n",
    "        res[i] = (a[i] + b[i]) / 2;\n",
    "    }\n",
    "    return res;\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5959864d",
   "metadata": {},
   "outputs": [],
   "source": [
    "// Print v[begin:end] as grayscale ASCII blocks (dark = 0, light = 1).\n",
    "void print_grayscale(const std::vector<double>& v, std::size_t begin, std::size_t end)\n",
    "{\n",
    "    for (std::size_t i = begin; i < end; ++i)\n",
    "    {\n",
    "        // map value in [0, 1] to a 24-bit truecolor gray (256 levels)\n",
    "        int gray = static_cast<int>(v[i] * 255);\n",
    "        std::cout << \"\\033[38;2;\" << gray << \";\" << gray << \";\" << gray << \"m█\";\n",
    "    }\n",
    "    std::cout << \"\\033[0m\\n\";\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2068e201",
   "metadata": {},
   "outputs": [],
   "source": [
    "std::size_t n = 1000;\n",
    "std::vector<double> const a = random_vector(n);\n",
    "std::vector<double> const b = random_vector(n);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9d5a7adc-5b55-4224-a4ca-5e438d1fc4d4",
   "metadata": {},
   "outputs": [],
   "source": [
    "auto const res_simd = mean_simd(a, b);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d279cd74-c397-4ea8-af15-3745f043f656",
   "metadata": {},
   "outputs": [],
   "source": [
    "auto const res_scalar = mean_scalar(a, b);"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a04ae145",
   "metadata": {},
   "source": [
    "Each block is a value in `[0, 1]` shown as a grayscale shade (dark = low, light = high). The vector is drawn in chunks of 100, each chunk stacking input `a`, the SIMD mean, the scalar mean, and input `b`. The two mean rows should look identical, each sitting between its inputs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "59ca84ad-8b97-4523-af0d-d3eca01b41d0",
   "metadata": {},
   "outputs": [],
   "source": [
    "std::size_t const width = 100;\n",
    "for (std::size_t begin = 0; begin < n; begin += width)\n",
    "{\n",
    "    std::size_t end = std::min(begin + width, n);\n",
    "    std::cout << \"a:      \";\n",
    "    print_grayscale(a, begin, end);\n",
    "    std::cout << \"simd:   \";\n",
    "    print_grayscale(res_simd, begin, end);\n",
    "    std::cout << \"scalar: \";\n",
    "    print_grayscale(res_scalar, begin, end);\n",
    "    std::cout << \"b:      \";\n",
    "    print_grayscale(b, begin, end);\n",
    "    std::cout << \"\\n\";\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "90b86fb4-0426-43fd-b5fc-8626ef8732a5",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "C++17",
   "language": "cpp",
   "name": "xcpp17"
  },
  "language_info": {
   "codemirror_mode": "text/x-c++src",
   "file_extension": ".cpp",
   "mimetype": "text/x-c++src",
   "name": "C++",
   "nbconvert_exporter": "",
   "pygments_lexer": "",
   "version": "cxx17"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}