Use a 2×2 Gauss rule for linear quadrilaterals (#105)

examples/creepingFlowScript/README.md

@@ -19,6 +19,13 @@ This example solves a 2D lid-driven cavity flow using the creeping flow solver w
 [tutorial](https://feascript.com/tutorials/lid-driven-cavity-2d-creeping-flow.html) on the
 FEAScript website.
 
+### 2. Couette Flow in a 2D Channel (`couetteFlow2DCreepingFlow.js`)
+
+A channel whose top wall moves with velocity U, with a fixed bottom wall and stress-free (open)
+ends. The exact solution is the linear velocity profile u = U y / H with v = 0 and zero pressure,
+which the Taylor-Hood elements reproduce exactly. The example prints the velocity along the left
+end and the largest deviation from the exact solution.
+
 ## Running the Node.js examples
 
 ### 1. Create `package.json` with ES module support
@@ -39,4 +46,6 @@ npm install feascript
 node lidDrivenCavity2DCreepingFlow.js
 ```
 
+(or `node couetteFlow2DCreepingFlow/couetteFlow2DCreepingFlow.js`)
+
 (for the "Creeping Flow on a Lid-Driven 2D Cavity" example)

examples/creepingFlowScript/couetteFlow2DCreepingFlow/couetteFlow2DCreepingFlow.js

@@ -0,0 +1,65 @@
+/**
+ * ════════════════════════════════════════════════════════════════
+ *  FEAScript Core Library
+ *  Finite Element Simulation in JavaScript
+ *  Version: 0.3.0 | https://feascript.com
+ *  MIT License © 2023–2026 FEAScript
+ * ════════════════════════════════════════════════════════════════
+ */
+
+// Couette flow: a 4 × 2 channel whose top wall moves with velocity U, with stress-free (open) ends.
+// The exact solution is the linear velocity profile u = U y / H, v = 0, with zero pressure, which
+// Taylor-Hood elements reproduce exactly
+
+// Import FEAScript library
+import { FEAScriptModel, printVersion } from "feascript";
+
+console.log("FEAScript Version:", printVersion);
+
+const wallVelocity = 1; // U
+const channelHeight = 2; // H
+
+// 1. Create the model and select the physics
+const model = new FEAScriptModel();
+model.setModelConfig("creepingFlowScript", {
+  coefficientFunctions: {
+    viscosity: 1, // μ
+  },
+});
+
+// 2. Define the mesh (creeping flow needs quadratic elements)
+model.setMeshConfig({
+  meshDimension: "2D",
+  elementOrder: "quadratic",
+  numElementsX: 4,
+  numElementsY: 2,
+  maxX: 4,
+  maxY: channelHeight,
+});
+
+// 3. Apply the boundary conditions
+model.addBoundaryCondition("0", ["constantVelocity", 0, 0]); // Bottom: fixed wall
+model.addBoundaryCondition("1", ["stressFree"]); // Left: open end
+model.addBoundaryCondition("2", ["constantVelocity", wallVelocity, 0]); // Top: moving wall
+model.addBoundaryCondition("3", ["stressFree"]); // Right: open end
+
+// 4. Solve
+model.setSolverMethod("lusolve");
+const { nodesCoordinates, fields } = model.solve();
+
+// 5. Print the results
+const { nodesYCoordinates } = nodesCoordinates;
+const largestDeviation = Math.max(
+  ...nodesYCoordinates.map((y, i) => Math.abs(fields.u[i] - (wallVelocity * y) / channelHeight)),
+);
+console.log("Velocity u along the left end (x = 0):");
+nodesYCoordinates.forEach((y, i) => {
+  if (nodesCoordinates.nodesXCoordinates[i] === 0) {
+    console.log(
+      `  y = ${y.toFixed(2)}: u = ${fields.u[i].toFixed(4)} (exact ${(wallVelocity * y) / channelHeight})`,
+    );
+  }
+});
+console.log(`Largest deviation from u = U y / H: ${largestDeviation.toExponential(2)}`);
+console.log(`Largest |v|: ${Math.max(...fields.v.map(Math.abs)).toExponential(2)} (exact 0)`);
+console.log(`Largest |p|: ${Math.max(...fields.p.map(Math.abs)).toExponential(2)} (exact 0)`);

examples/creepingFlowScript/lidDrivenCavity2DCreepingFlow/lidDrivenCavity2DCreepingFlow.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/eulerBernoulliBeamScript/README.md

@@ -20,6 +20,14 @@ one end, supported by a roller at midspan, and connected to a linear elastic spr
 end. It carries a uniformly distributed load over the clamped span, a concentrated moment at the
 roller, and a point load at the free end.
 
+### 2. Cantilever Beam with a Tip Load (`cantileverBeam1D.js`)
+
+A 4 m cantilever, clamped at x = 0, with a point load P at the free end. The example prints the tip
+deflection and rotation next to the closed-form values w(L) = P L³ / (3 EI) and
+θ(L) = P L² / (2 EI), which the cubic Hermite elements reproduce exactly. For a plotted version, see
+the corresponding [tutorial](https://feascript.com/tutorials/cantilever-beam-1d.html) on the
+FEAScript website.
+
 ## Model overview
 
 The 1D Euler-Bernoulli beam solver (`eulerBernoulliBeamScript`) assembles the fourth-order beam
@@ -108,3 +116,5 @@ npm install feascript
 ```bash
 node clampedSpringSupportedBeam1D/clampedSpringSupportedBeam1D.js
 ```
+
+(or `node cantileverBeam1D/cantileverBeam1D.js`)

examples/eulerBernoulliBeamScript/cantileverBeam1D/cantileverBeam1D.js

@@ -0,0 +1,56 @@
+/**
+ * ════════════════════════════════════════════════════════════════
+ *  FEAScript Core Library
+ *  Finite Element Simulation in JavaScript
+ *  Version: 0.3.0 | https://feascript.com
+ *  MIT License © 2023–2026 FEAScript
+ * ════════════════════════════════════════════════════════════════
+ */
+
+// Cantilever beam of length L, clamped at x = 0, with a point load P at the free end. The exact
+// tip deflection and rotation are w(L) = P L³ / (3 EI) and θ(L) = P L² / (2 EI)
+
+// Import FEAScript library
+import { FEAScriptModel, printVersion } from "feascript";
+
+console.log("FEAScript Version:", printVersion);
+
+const beamLength = 4; // L (m)
+const bendingStiffness = 1.0e5; // EI (N·m²)
+const tipLoad = -1000; // P (N), downward
+
+// 1. Create the model and select the physics
+const model = new FEAScriptModel();
+model.setModelConfig("eulerBernoulliBeamScript", {
+  coefficientFunctions: {
+    EI: (x) => bendingStiffness, // Bending stiffness
+  },
+});
+
+// 2. Define the mesh
+model.setMeshConfig({
+  meshDimension: "1D",
+  elementOrder: "linear",
+  numElementsX: 4,
+  maxX: beamLength,
+});
+
+// 3. Apply the boundary conditions (keyed by node number)
+model.addBoundaryCondition("1", [["fixed"]]); // Clamped support at x = 0
+model.addBoundaryCondition("5", [["force", tipLoad]]); // Tip load at x = L (last node)
+
+// 4. Solve
+model.setSolverMethod("lusolve");
+const { fields } = model.solve();
+
+// 5. Print the results
+const tipDeflection = fields.w.at(-1);
+const tipRotation = fields.theta.at(-1);
+const exactTipDeflection = (tipLoad * beamLength ** 3) / (3 * bendingStiffness);
+const exactTipRotation = (tipLoad * beamLength ** 2) / (2 * bendingStiffness);
+console.log(
+  `Tip deflection w(L): ${tipDeflection.toExponential(6)} m (exact ${exactTipDeflection.toExponential(6)} m)`,
+);
+console.log(
+  `Tip rotation θ(L): ${tipRotation.toExponential(6)} rad (exact ${exactTipRotation.toExponential(6)} rad)`,
+);

examples/eulerBernoulliBeamScript/clampedSpringSupportedBeam1D/clampedSpringSupportedBeam1D.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/frontPropagationScript/solidificationFront2D/solidificationFront2D.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/generalFormPDEScript/advectionDiffusion1D/advectionDiffusion1D.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/generalFormPDEScript/nonlinearReactionDiffusion1D/nonlinearReactionDiffusion1D.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/heatConductionScript/README.md

@@ -52,6 +52,15 @@ Demonstrates both `thermalConductivity(x, y)` and `heatSource(x, y)` as function
 2D structured mesh. The domain is split into a high-conductivity metal half and a low-conductivity
 ceramic half, with a localised volumetric heat source in the upper strip.
 
+### 7. NAFEMS T4 Benchmark (`heatConduction2DNafemsT4.js`)
+
+Reproduces the NAFEMS T4 benchmark for 2D heat conduction with convection: a 0.6 m × 1.0 m plate with
+a fixed temperature of 100 °C on the bottom edge, an insulated left edge, and convection
+(h = 750 W/(m²·K), ambient 0 °C) on the top and right edges. It prints the temperature at point E
+(x = 0.6 m, y = 0.2 m) next to the published target of 18.3 °C. The 12 × 20 quadratic mesh gives
+18.26 °C, within 0.01 °C of the converged value of 18.25 °C. This is an independent reproduction
+of the benchmark, not sponsored or endorsed by NAFEMS.
+
 ## Spatially varying coefficients
 
 Both `thermalConductivity` and `heatSource` can be provided either as constants (scalars) or as
@@ -99,4 +108,5 @@ npm install feascript
 node heatConduction2DFin.js
 ```
 
-(or `node heatConduction2DFinGmsh.js` / `node heatConduction2DRhomFinGmsh.js`)
+(or `node heatConduction2DFinGmsh.js` / `node heatConduction2DRhomFinGmsh.js` /
+`node heatConduction2DNafemsT4/heatConduction2DNafemsT4.js`)

examples/heatConductionScript/heatConduction1DWall/heatConduction1DVaryingK.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/heatConductionScript/heatConduction1DWall/heatConduction1DWall.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/heatConductionScript/heatConduction2DFin/heatConduction2DFin.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/heatConductionScript/heatConduction2DFin/heatConduction2DFinGmsh.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/heatConductionScript/heatConduction2DFin/heatConduction2DRhomFinGmsh.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/heatConductionScript/heatConduction2DFin/heatConduction2DVaryingKQ.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

examples/heatConductionScript/heatConduction2DNafemsT4/heatConduction2DNafemsT4.js

@@ -0,0 +1,53 @@
+/**
+ * ════════════════════════════════════════════════════════════════
+ *  FEAScript Core Library
+ *  Finite Element Simulation in JavaScript
+ *  Version: 0.3.0 | https://feascript.com
+ *  MIT License © 2023–2026 FEAScript
+ * ════════════════════════════════════════════════════════════════
+ */
+
+// NAFEMS T4 benchmark: a 0.6 m × 1.0 m plate with a fixed temperature on the bottom edge, an
+// insulated left edge, and convection to 0 °C on the top and right edges. The published target
+// temperature at point E (x = 0.6 m, y = 0.2 m) is 18.3 °C (18.25 °C before rounding)
+
+// Import FEAScript library
+import { FEAScriptModel, printVersion } from "feascript";
+
+console.log("FEAScript Version:", printVersion);
+
+// 1. Create the model and select the physics
+const model = new FEAScriptModel();
+model.setModelConfig("heatConductionScript", {
+  coefficientFunctions: {
+    thermalConductivity: 52, // k = 52 W/(m·K)
+  },
+});
+
+// 2. Define the mesh
+model.setMeshConfig({
+  meshDimension: "2D",
+  elementOrder: "quadratic",
+  numElementsX: 12,
+  numElementsY: 20,
+  maxX: 0.6,
+  maxY: 1.0,
+});
+
+// 3. Apply the boundary conditions
+model.addBoundaryCondition("0", ["constantTemperature", 100]); // Bottom: T = 100 °C
+model.addBoundaryCondition("1", ["symmetry"]); // Left: insulated
+model.addBoundaryCondition("2", ["convection", 750, 0]); // Top: h = 750 W/(m²·K), T∞ = 0 °C
+model.addBoundaryCondition("3", ["convection", 750, 0]); // Right: h = 750 W/(m²·K), T∞ = 0 °C
+
+// 4. Solve
+model.setSolverMethod("lusolve");
+const { solutionVector, nodesCoordinates } = model.solve();
+
+// 5. Print the results
+const { nodesXCoordinates, nodesYCoordinates } = nodesCoordinates;
+const pointE = nodesXCoordinates.findIndex(
+  (x, i) => Math.abs(x - 0.6) < 1e-10 && Math.abs(nodesYCoordinates[i] - 0.2) < 1e-10,
+);
+console.log(`Temperature at point E (x = 0.6 m, y = 0.2 m): ${solutionVector[pointE].toFixed(2)} °C`);
+console.log("NAFEMS target temperature at point E: 18.3 °C (converged value 18.25 °C)");

package.json

@@ -1,6 +1,6 @@
 {
   "name": "feascript",
-  "version": "0.3.0-RC",
+  "version": "0.3.0",
   "description": "Finite element simulation library built in JavaScript",
   "main": "dist/feascript.cjs.js",
   "module": "dist/feascript.esm.js",

src/FEAScript.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

src/elements/basisFunctions.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */

src/elements/elementGeometry.js

@@ -2,7 +2,7 @@
  * ════════════════════════════════════════════════════════════════
  *  FEAScript Core Library
  *  Finite Element Simulation in JavaScript
- *  Version: 0.3.0 (RC) | https://feascript.com
+ *  Version: 0.3.0 | https://feascript.com
  *  MIT License © 2023–2026 FEAScript
  * ════════════════════════════════════════════════════════════════
  */