Feat/laasonen

.gitignore

@@ -15,7 +15,6 @@
 .cache/
 
 # Outputs
-/result/*
 /results/*
 !.gitkeep
 
@@ -26,5 +25,4 @@
 .temp/
 
 # Logs
-*.log
-
+*.log

include/method.hpp

@@ -1,68 +1,67 @@
-#pragma once
-/**
- * @file method.hpp
- * @brief Numerical scheme interface for the 1D heat equation.
- */
+#ifndef METHOD_HPP
+#define METHOD_HPP
 
 #include "grid.hpp"
 #include <memory>
-#include <string>
 #include <vector>
 
 /**
- * @brief Supported numerical schemes.
+ * @brief Enumeration of the numerical schemes available in the solver.
  */
 enum class SchemeKind
 {
-    DuFortFrankel,
-    Richardson
-    // Laasonen, CrankNicolson (to be added later)
+    Richardson,     ///< Central time, central space (explicit, unstable).
+    DuFortFrankel,  ///< Modified Richardson (explicit, stable).
+    Laasonen,       ///< Simple Implicit (forward time, central space).
+    CrankNicolson   ///< Trapezoidal Implicit (second order accuracy).
 };
 
 /**
- * @brief Abstract base class for time-stepping methods.
+ * @class Method
+ * @brief Abstract interface for time-integration schemes solving the 1D heat equation.
+ *
+ * This interface abstracts away the difference between explicit/implicit schemes
+ * and 2-level/3-level time stepping schemes.
  */
 class Method
 {
-  public:
+public:
     virtual ~Method() = default;
 
-    /// @return human-readable scheme name.
-    virtual std::string name() const = 0;
-
     /**
-     * @brief Advance one time step, computing T^{n+1} from T^n (and possibly T^{n-1}).
-     * @param g Grid metadata (spacing and node count)
-     * @param D Diffusivity [cm^2/h]
-     * @param dt Time step [h]
-     * @param Tprev Read-only temperature at the previous time layer T^{n-1}
-     *              (may be empty when the scheme is two-level)
-     * @param Tcurr Read-only temperature at the current time layer T^{n}
-     * @param Tnext Output container that must be filled with T^{n+1}
+     * @brief Computes the temperature field at the next time step (T^{n+1}).
+     *
+     * @param g     Grid metadata (spacing dx and node count Nx).
+     * @param D     Thermal diffusivity [cm^2/h].
+     * @param dt    Time step size [h].
+     * @param Tprev Temperature at time step n-1 (Read-only).
+     * Required only for 3-level schemes (Richardson, DuFort-Frankel).
+     * For 2-level schemes, this vector may be ignored.
+     * @param Tcurr Temperature at time step n (Read-only).
+     * @param Tnext Output vector to be filled with temperature at time step n+1.
      */
     virtual void step(const Grid& g,
                       double D,
                       double dt,
                       const std::vector<double>& Tprev,
                       const std::vector<double>& Tcurr,
-                      std::vector<double>& Tnext) = 0;
+                      std::vector<double>& Tnext) const = 0;
 
     /**
-     * @brief Indicate whether the scheme needs the previous time layer.
+     * @brief Indicates if the scheme requires data from time step n-1.
      *
-     * Three-level schemes (e.g. DuFort-Frankel) require both T^{n} and T^{n-1}
-     * to produce T^{n+1}. Two-level schemes can ignore the returned value.
+     * @return true  For 3-level schemes (e.g., Richardson, DuFort-Frankel).
+     * @return false For 2-level schemes (e.g., Laasonen, Crank-Nicolson).
      */
-    virtual bool uses_previous_step() const noexcept
-    {
-        return false;
-    }
+    virtual bool uses_previous_step() const = 0;
 };
 
 /**
- * @brief Factory helper that builds a Method instance for the requested scheme.
- * @param scheme Scheme identifier (e.g. SchemeKind::DuFortFrankel)
- * @return Unique pointer owning the concrete method implementation
- * @throws std::invalid_argument if the scheme is not supported
+ * @brief Factory function to create a concrete method instance.
+ *
+ * @param scheme The type of numerical scheme requested.
+ * @return std::unique_ptr<Method> Pointer to the created solver method.
  */
-[[nodiscard]] std::unique_ptr<Method> make_method(SchemeKind scheme);
+std::unique_ptr<Method> make_method(SchemeKind scheme);
+
+#endif // METHOD_HPP

include/methods/laasonen.hpp

@@ -0,0 +1,30 @@
+#pragma once
+
+#include "method.hpp"
+#include <string>
+#include <vector>
+
+/**
+ * @brief Laasonen implicit method (Backward Euler in time + Central in space).
+ *        Unconditionally stable for the 1D heat equation.
+ */
+class Laasonen : public Method
+{
+  public:
+    std::string name() const override
+    {
+        return "Laasonen";
+    }
+
+    void step(const Grid& g,
+              double D,
+              double dt,
+              const std::vector<double>& Tprev,
+              const std::vector<double>& Tcurr,
+              std::vector<double>& Tnext) const override;
+
+    bool uses_previous_step() const noexcept override
+    {
+        return false;
+    }
+};

src/io.cpp

@@ -10,21 +10,23 @@ namespace fs = std::filesystem;
 /**
  * @brief Save a temperature profile (x, T) to a CSV file.
  *
- * This function automatically creates the "result/" directory if missing.
+ * This function automatically creates the directory structure if missing.
  * Example output path: "result/dufort_frankel/t_0.10.csv"
  *
- * @param path  Relative path to the file (inside result/)
- * @param x     Spatial coordinates [cm]
- * @param T     Temperature values [°C]
+ * @param path  Relative or absolute path to the file.
+ * @param x     Spatial coordinates [cm].
+ * @param T     Temperature values [°C].
  */
 void save_profile_csv(const std::string& path,
                       const std::vector<double>& x,
                       const std::vector<double>& T)
 {
+    // Validate input sizes
     if (x.size() != T.size())
         throw std::runtime_error("save_profile_csv: vector size mismatch");
 
-    // --- Ensure the directory exists ---
+    // --- Ensure the directory exists (Resolved Conflict) ---
+    // We use the safer check from the 'dev' branch.
     fs::path fullPath(path);
     if (fullPath.has_parent_path())
     {
@@ -41,4 +43,4 @@ void save_profile_csv(const std::string& path,
     file << std::fixed << std::setprecision(6);
     for (std::size_t i = 0; i < x.size(); ++i)
         file << x[i] << "," << T[i] << "\n";
-}
+}

src/main.cpp

@@ -1,7 +1,7 @@
 /**
  * @file main.cpp
  * @brief Test driver for the 1D heat conduction solver.
- *        Supports DuFort-Frankel and Richardson schemes.
+ * Supports DuFort-Frankel, Richardson, and Laasonen schemes.
  */
 
 #include "grid.hpp"
@@ -15,16 +15,21 @@
 #include <iostream>
 #include <sstream>
 #include <string>
+#include <stdexcept>
+#include <vector>
 
 namespace fs = std::filesystem;
 
 /**
  * @brief Run a single simulation with the specified scheme.
  *
- * @param phys Physical parameters
- * @param num Numerical parameters
- * @param scheme The numerical scheme to use
- * @param schemeName String name for output directory
+ * This helper function avoids code duplication by handling the solver setup,
+ * directory creation, and output callback for any given scheme.
+ *
+ * @param phys       Physical parameters.
+ * @param num        Numerical parameters.
+ * @param scheme     The numerical scheme to use.
+ * @param schemeName String name for output directory (e.g. "Laasonen").
  */
 void run_simulation(const PhysParams& phys,
                     const NumParams& num,
@@ -37,21 +42,23 @@ void run_simulation(const PhysParams& phys,
     HeatSolver solver(phys, num, scheme);
 
     // --- Prepare output directory ---
+    // Results will be saved in: results/<schemeName>/
     fs::path outDir = fs::path("results") / schemeName;
     if (!fs::exists(outDir))
     {
         fs::create_directories(outDir);
     }
 
     // --- Output callback ---
+    // Called by the solver at t=0 and every 'outEvery' hours.
     auto onDump = [&](double t, const Grid& g, const std::vector<double>& T)
     {
         std::ostringstream name;
+        // Filename format: profile_t_0.10h.csv
         name << "profile_t_" << std::fixed << std::setprecision(2) << t << "h.csv";
+
         fs::path fullPath = outDir / name.str();
         save_profile_csv(fullPath.string(), g.x, T);
-        // Optional: reduce verbosity for batch runs
-        // std::cout << "Saved " << fullPath.string() << "\n";
     };
 
     // --- Run ---
@@ -67,23 +74,25 @@ int main()
 
         // --- Physical parameters ---
         PhysParams phys;
-        phys.L_cm = 31.0;
-        phys.Tin = 38.0;
-        phys.Tsur = 149.0;
+        phys.L_cm   = 31.0;
+        phys.Tin    = 38.0;
+        phys.Tsur   = 149.0;
         phys.D_cm2h = 93.0;
 
         // --- Numerical parameters ---
         NumParams num;
+
         // Ask for dx once
         num.dx = askForDouble(
                 "Enter Δx [cm] (suggested 0.05): ",
                 [](double v) { return v > 0.0 && v <= 1.0; },
                 "Δx must be positive and <= 1.0 cm.");
-        num.dt = 0.01;
-        num.tEnd = 0.5;
-        num.outEvery = 0.1;
+
+        num.dt       = 0.01;  // Time step [h]
+        num.tEnd     = 0.5;   // Total duration [h]
+        num.outEvery = 0.1;   // Output interval [h]
 
-        // --- Grid validation ---
+        // --- Grid validation (Early check) ---
         Grid g(phys.L_cm, num.dx);
         validate_grid(g);
         std::cout << "Grid created ✅  Nx = " << g.Nx << ", dx = " << g.dx << " cm\n";
@@ -92,23 +101,35 @@ int main()
         std::cout << "\nSelect numerical method:\n";
         std::cout << "1. DuFort-Frankel\n";
         std::cout << "2. Richardson\n";
-        std::cout << "3. Run All (Batch)\n";
+        std::cout << "3. Laasonen\n";
+        std::cout << "4. Run All (Batch)\n";
 
+        // We use askForDouble (or askForInteger) to get a valid choice between 1 and 4
         int choice = static_cast<int>(askForDouble(
-                "Enter choice (1-3): ",
-                [](double v) { return v >= 1.0 && v <= 3.0; },
-                "Please enter 1, 2, or 3."));
+                "Enter choice (1-4): ",
+                [](double v) { return v >= 1.0 && v <= 4.0; },
+                "Please enter a number between 1 and 4."));
+
+        // --- Execution Logic ---
 
-        if (choice == 1 || choice == 3)
+        // 1. DuFort-Frankel
+        if (choice == 1 || choice == 4)
         {
             run_simulation(phys, num, SchemeKind::DuFortFrankel, "DuFort-Frankel");
         }
 
-        if (choice == 2 || choice == 3)
+        // 2. Richardson
+        if (choice == 2 || choice == 4)
         {
             run_simulation(phys, num, SchemeKind::Richardson, "Richardson");
         }
 
+        // 3. Laasonen
+        if (choice == 3 || choice == 4)
+        {
+            run_simulation(phys, num, SchemeKind::Laasonen, "Laasonen");
+        }
+
         std::cout << "\nAll requested simulations completed successfully ✅\n";
         return 0;
     }
@@ -117,4 +138,4 @@ int main()
         std::cerr << "Fatal error: " << e.what() << "\n";
         return 1;
     }
-}
+}