// E06: fixed 2x2 toy, THREE steps; not the official gamma experiment.
#include <array>
#include <cmath>
#include <filesystem>
#include <fstream>
#include <iomanip>
#include <iostream>
#include <limits>
#include <stdexcept>
using Vec = std::array<double, 2>;

Vec matvec(const std::array<double, 4>& T, const Vec& x) {
    Vec y{};
    for (int i = 0; i < 2; ++i) {
        for (int j = 0; j < 2; ++j) y[i] += T[i*2+j] * x[j];
    }
    return y;
}
int main(int argc, char* argv[]) {
    try {
        if (argc > 2) throw std::runtime_error("Use: program [output.csv]");
        const std::filesystem::path path =
            (argc == 2) ? argv[1] : "data/toy_run01.csv";
        if (std::filesystem::exists(path))
            throw std::runtime_error("Output exists; choose a new filename.");
        std::ofstream csv(path);
        if (!csv) throw std::runtime_error("Cannot open output; check parent.");
        const std::array<double, 4> T = {2, 1, 1, 2};
        Vec x = {1, 1};
        csv << "iteration,norm2\n";
        csv << std::setprecision(std::numeric_limits<double>::max_digits10);
        std::cout << std::setprecision(12);
        for (int k = 1; k <= 3; ++k) {
            const Vec next = matvec(T, x);  // The whole old x is read first.
            x = next;
            const double norm = std::sqrt(x[0]*x[0] + x[1]*x[1]);
            const double expected = std::pow(3.0, k) * std::sqrt(2.0);
            if (!std::isfinite(norm) || std::abs(norm-expected) > 1e-12*expected)
                throw std::runtime_error("Toy calculation check failed.");
            csv << k << ',' << norm << '\n';
            std::cout << k << ": x=" << x[0] << ',' << x[1]
                      << " norm=" << norm << '\n';
        }
        csv.close();
        if (!csv) throw std::runtime_error("Output write/close failed.");
        std::cout << "Saved: " << path.string() << '\n';
        return 0;
    } catch (const std::exception& e) {
        std::cerr << "Error: " << e.what() << '\n';
        return 1;
    }
}