// Copyright (C) 2026 Kiyotsugu Arai
// SPDX-License-Identifier: LGPL-3.0-or-later
//
// example_doc_vector_samples.cpp
// Verification harness for samples in api/Vector.html (ja+en).

#include <math/core/vector.hpp>
#include <cmath>
#include <iomanip>
#include <iostream>

using namespace sangi;

template <typename V>
void print_vec(const V& v) {
    std::cout << '{';
    for (std::size_t i = 0; i < v.size(); ++i) {
        if (i) std::cout << ", ";
        std::cout << v[i];
    }
    std::cout << '}';
}

int main() {
    std::cout << std::setprecision(15);

    // --- v1, v2, v3, v4 constructors ----------------------------------
    {
        std::cout << "[constructors]\n";
        Vector<double> v1;                  // empty
        Vector<double> v2(5);               // size 5, zero-initialized
        Vector<double> v3(5, 1.0);          // size 5, all 1.0
        Vector<double> v4{1.0, 2.0, 3.0};   // initializer list
        std::cout << "  v1.size=" << v1.size() << "  v2="; print_vec(v2);
        std::cout << "  v3="; print_vec(v3);
        std::cout << "  v4="; print_vec(v4); std::cout << "\n";
    }

    // --- element access ----------------------------------------------
    {
        std::cout << "[element access]\n";
        Vector<double> v{10, 20, 30};
        double a = v[0];
        double b = v.at(1);
        v[2] = 99;
        std::cout << "  v[0]=" << a << "  v.at(1)=" << b << "  v="; print_vec(v); std::cout << "\n";
    }

    // --- arithmetic + expression templates ---------------------------
    {
        std::cout << "[arithmetic]\n";
        Vector<double> a{1, 2, 3}, b{4, 5, 6};
        Vector<double> c = a + b;           // {5, 7, 9}
        Vector<double> d = a - b;           // {-3, -3, -3}
        Vector<double> e = 2.0 * a;         // {2, 4, 6}
        Vector<double> f = a / 2.0;         // {0.5, 1, 1.5}
        std::cout << "  a+b="; print_vec(c);
        std::cout << "  a-b="; print_vec(d);
        std::cout << "  2*a="; print_vec(e);
        std::cout << "  a/2="; print_vec(f);
        std::cout << "\n";

        Vector<double> result = 2.0 * a + b - 0.5 * Vector<double>{7, 8, 9};
        std::cout << "  2a+b-0.5*{7,8,9}="; print_vec(result); std::cout << "\n";
    }

    // --- dot / norms --------------------------------------------------
    {
        std::cout << "[dot / norms]\n";
        Vector<double> a{1, 2, 3}, b{4, 5, 6};
        std::cout << "  dot(a,b)=" << dot(a, b);
        std::cout << "  norm(a)=" << norm(a);
        std::cout << "  norm_l1(a)=" << norm_l1(a);
        std::cout << "  norm_linf(a)=" << norm_linf(a);
        std::cout << "  norm_lp(a,3)=" << norm_lp(a, 3.0);
        std::cout << "\n";
    }

    // --- normalize ---------------------------------------------------
    {
        std::cout << "[normalize]\n";
        Vector<double> v{3, 4};
        Vector<double> u = normalized(v);
        std::cout << "  normalized({3,4})="; print_vec(u);
        std::cout << "  norm(u)=" << norm(u) << "\n";
    }

    // --- block operations --------------------------------------------
    {
        std::cout << "[block ops]\n";
        Vector<double> v{10, 20, 30, 40, 50};
        auto h = v.head(2);
        auto t = v.tail(2);
        auto s = v.segment(1, 3);
        std::cout << "  head(2)="; print_vec(h);
        std::cout << "  tail(2)="; print_vec(t);
        std::cout << "  segment(1,3)="; print_vec(s);
        std::cout << "\n";
    }

    // --- cross (StaticVector<T,3>) ------------------------------------
    {
        std::cout << "[cross]\n";
        StaticVector<double, 3> i{1, 0, 0}, j{0, 1, 0};
        auto k = cross(i, j);
        std::cout << "  i x j="; print_vec(k); std::cout << "\n";
    }

    // --- axpy / axpby -------------------------------------------------
    {
        std::cout << "[axpy / axpby]\n";
        Vector<double> x{1, 2, 3}, y{10, 20, 30};
        axpy(2.0, x, y);  // y += 2*x  =>  {12, 24, 36}
        std::cout << "  y after axpy(2,x)="; print_vec(y); std::cout << "\n";

        Vector<double> a{1, 2, 3}, b{4, 5, 6}, r(3);
        axpby(2.0, a, 3.0, b, r);  // r = 2a + 3b
        std::cout << "  axpby(2,a,3,b)="; print_vec(r); std::cout << "\n";
    }

    // --- VectorMap (zero-copy view) -----------------------------------
    {
        std::cout << "[VectorMap]\n";
        double raw[] = {1.0, 2.0, 3.0};
        VectorMap<double> v(raw, 3);
        v[1] = 5.0;
        std::cout << "  raw[1] after v[1]=5: " << raw[1] << "\n";

        double data[] = {1, 2, 3, 4, 5};
        ConstVectorMap<double> view(data, 5);
        // dot has no overload that accepts ConstVectorMap directly -- convert
        // to Vector<T> first (the implicit conversion operator is provided).
        Vector<double> view_vec = view;
        std::cout << "  dot(view,view)=" << dot(view_vec, view_vec) << "\n";
    }

    return 0;
}