Add benchmark test for NAFEMS T4 (2D heat conduction with convection)

Author: nikoscham

tests/verification/benchmark/heatConduction2DNafemsT4/README.md

@@ -0,0 +1,56 @@
+# Benchmark Test — NAFEMS T4 (2D Heat Conduction with Convection)
+
+## Purpose
+
+This test replicates the NAFEMS T4 benchmark, also documented as a DIANA FEA tutorial
+([2DHeatTransferConvection.pdf](https://tutorials.dianafea.com/2DHeatTransferConvection.pdf))
+and by [Altair](https://help.altair.com/hwsolvers/os/topics/solvers/os/nafems_test_problem_t4_r.htm),
+and checks that `heatConductionScript` reproduces the published target temperature for a 2D
+steady-state conduction problem with a convective (Robin) boundary condition.
+
+This is an independent reproduction of the NAFEMS T4 problem and is not sponsored, endorsed, or
+affiliated with NAFEMS.
+
+## Problem setup
+
+A 0.6 m by 1.0 m rectangular plate `ABCD` (`A` at the origin, `B` at `(0.6, 0)`, `C` at
+`(0.6, 1.0)`, `D` at `(0, 1.0)`), thermal conductivity `k = 52 W/(m·K)`, no internal heat
+generation:
+
+| Edge        | Boundary condition                                |
+| ----------- | ------------------------------------------------- |
+| AB (bottom) | Constant temperature, T = 100 °C                  |
+| AD (left)   | Insulated (natural, zero-flux — left unspecified) |
+| CD (top)    | Convection, h = 750 W/(m²·K), ambient T = 0 °C    |
+| BC (right)  | Convection, h = 750 W/(m²·K), ambient T = 0 °C    |
+
+The published target is the temperature at point E, located on edge BC at `(x = 0.6, y = 0.2)`:
+**18.3 °C**.
+
+## Expected values
+
+| Case           | Temperature at E | Compared against        | Tolerance |
+| -------------- | ---------------- | ----------------------- | --------- |
+| Mesh (12 × 20) | ≈ 18.26 °C       | NAFEMS target (18.3 °C) | `0.1` °C  |
+
+## How to run
+
+From the repository root:
+
+```bash
+node tests/verification/benchmark/heatConduction2DNafemsT4/benchmark.test.js
+```
+
+The `test` script in `package.json` also runs this file, so `npm test` works too.
+
+A passing run prints two `PASS:` lines followed by `2 passed, 0 failed.`. A failing run prints
+one or more `FAIL:` lines, ends with the same summary line format, and exits with code 1.
+
+## After modifying the code
+
+| Situation                                                          | Action                                                                                                          |
+| ------------------------------------------------------------------ | --------------------------------------------------------------------------------------------------------------- |
+| Bug fix that should not change results                             | Run the test — it must still pass.                                                                              |
+| Intentional algorithm change (new integration rule, mesh ordering) | Re-check the result stays close to the NAFEMS target and adjust the tolerance if needed.                        |
+| Adding 8-node serendipity elements to FEAScript                    | Add a case using 8-node elements on the original 3 × 5 mesh for a tighter comparison against the NAFEMS target. |
+| New boundary condition API (e.g. convection parameter order)       | Update both this test and the reference JavaScript examples together.                                           |

tests/verification/benchmark/heatConduction2DNafemsT4/benchmark.test.js

@@ -0,0 +1,105 @@
+/**
+ * ════════════════════════════════════════════════════════════════
+ *  FEAScript Core Library
+ *  Lightweight Finite Element Simulation in JavaScript
+ *  Version: 0.3.0 (RC) | https://feascript.com
+ *  MIT License © 2023–2026 FEAScript
+ * ════════════════════════════════════════════════════════════════
+ */
+
+/**
+ * Benchmark test for 2D steady-state heat conduction with convection (NAFEMS T4)
+ *
+ * Replicates the NAFEMS T4 benchmark, also documented as a DIANA FEA tutorial
+ * (https://tutorials.dianafea.com/2DHeatTransferConvection.pdf) and by Altair
+ * (https://help.altair.com/hwsolvers/os/topics/solvers/os/nafems_test_problem_t4_r.htm):
+ * a 0.6 m x 1.0 m rectangular plate ABCD with a fixed 100 degC temperature on the bottom edge
+ * AB, an insulated left edge AD, and convection (h = 750 W/(m2.K), ambient = 0 degC) on the top
+ * edge CD and right edge BC. The published target temperature at point E (x = 0.6 m, y = 0.2 m,
+ * on the convective right edge) is 18.3 degC.
+ *
+ * This test is an independent reproduction of the NAFEMS T4 problem and is not sponsored,
+ * endorsed, or affiliated with NAFEMS.
+ *
+ * Run: node tests/verification/benchmark/heatConduction2DNafemsT4/benchmark.test.js (or npm test)
+ */
+
+import * as mathjs from "mathjs";
+import { FEAScriptModel } from "../../../../src/FEAScript.js";
+import { basicLog, errorLog } from "../../../../src/utilities/logging.js";
+
+// FEAScript.js references `math` as a global (loaded via CDN in browser).
+// Set it here before any solve() call.
+globalThis.math = mathjs;
+
+const POINT_E = { x: 0.6, y: 0.2 };
+const NAFEMS_TARGET_T = 18.3;
+
+function runSimulation(numElementsX, numElementsY) {
+  const model = new FEAScriptModel();
+
+  model.setModelConfig("heatConductionScript", { coefficientFunctions: { thermalConductivity: 52 } });
+  model.setMeshConfig({
+    meshDimension: "2D",
+    elementOrder: "quadratic",
+    numElementsX,
+    numElementsY,
+    maxX: 0.6,
+    maxY: 1.0,
+  });
+
+  model.addBoundaryCondition("0", ["constantTemperature", 100]); // Bottom edge (AB), fixed 100 degC
+  // Left edge (AD) is insulated and left unspecified (natural, zero-flux boundary)
+  model.addBoundaryCondition("2", ["convection", 750, 0]); // Top edge (CD), convection to 0 degC
+  model.addBoundaryCondition("3", ["convection", 750, 0]); // Right edge (BC), convection to 0 degC
+  model.setSolverMethod("lusolve");
+
+  const { solutionVector, nodesCoordinates } = model.solve();
+  const { nodesXCoordinates, nodesYCoordinates } = nodesCoordinates;
+
+  const nodeIndex = nodesXCoordinates.findIndex(
+    (x, i) => Math.abs(x - POINT_E.x) < 1e-10 && Math.abs(nodesYCoordinates[i] - POINT_E.y) < 1e-10,
+  );
+
+  // solutionVector from math.lusolve is a nested array: [[T0], [T1], ...]
+  const temperatureAtE = Array.isArray(solutionVector[nodeIndex])
+    ? solutionVector[nodeIndex][0]
+    : solutionVector[nodeIndex];
+
+  return { nodeIndex, temperatureAtE };
+}
+
+let passed = 0;
+let failed = 0;
+
+function assert(condition, message) {
+  if (!condition) {
+    errorLog(`FAIL: ${message}`);
+    failed++;
+  } else {
+    basicLog(`PASS: ${message}`);
+    passed++;
+  }
+}
+
+basicLog("");
+basicLog("================================");
+basicLog("Starting benchmark test for NAFEMS T4 (2D heat conduction with convection)...");
+
+const CONVERGENCE_TOLERANCE = 0.1; // degC
+const { nodeIndex, temperatureAtE } = runSimulation(12, 20);
+
+assert(nodeIndex !== -1, `Found node E at (x=${POINT_E.x}, y=${POINT_E.y})`);
+assert(
+  Math.abs(temperatureAtE - NAFEMS_TARGET_T) < CONVERGENCE_TOLERANCE,
+  `Mesh (12x20): temperature at E expected close to NAFEMS target ${NAFEMS_TARGET_T}, got ${temperatureAtE} (tolerance ${CONVERGENCE_TOLERANCE})`,
+);
+
+basicLog("");
+if (failed > 0) {
+  errorLog(`${passed} passed, ${failed} failed.`);
+} else {
+  basicLog(`${passed} passed, ${failed} failed.`);
+}
+basicLog("================================");
+if (failed > 0) process.exit(1);