Format the workflow examples (model pages, first simulation, Gmsh import, web worker, llms.txt) in the same five numbered steps.

Author: nikoscham

docs/api/gmsh-import.md

@@ -56,18 +56,22 @@ Node.js has no `File` object for files on disk, so pass an object with a `text()
 import fs from "fs";
 import { FEAScriptModel, importGmshMesh } from "feascript";
 
-const meshContent = fs.readFileSync("rect_quad.msh", "utf8");
-const parsedMesh = await importGmshMesh({ text: async () => meshContent });
-
+// 1. Create the model and select the physics
 const model = new FEAScriptModel();
 model.setModelConfig("heatConductionScript");
+
+// 2. Load the Gmsh mesh
+const meshContent = fs.readFileSync("rect_quad.msh", "utf8");
+const parsedMesh = await importGmshMesh({ text: async () => meshContent });
 model.setMeshConfig({ parsedMesh });
 
+// 3. Apply the boundary conditions
 model.addBoundaryCondition("0", ["constantTemperature", 200]); // Physical curve "bottom", tag 1
 model.addBoundaryCondition("1", ["symmetry"]); // Physical curve "left", tag 2
 model.addBoundaryCondition("2", ["convection", 1, 20]); // Physical curve "top", tag 3
 model.addBoundaryCondition("3", ["constantTemperature", 200]); // Physical curve "right", tag 4
 
+// 4. Solve
 const { solutionVector, nodesCoordinates } = model.solve();
 ```
 

docs/api/web-worker.md

@@ -16,10 +16,12 @@ method is asynchronous and returns a `Promise`.
 ```javascript
 import { FEAScriptWorker, plotSolution } from "https://core.feascript.com/dist/feascript.esm.js";
 
+// 1. Create the model in a web worker and select the physics
 const model = new FEAScriptWorker();
 await model.ping(); // Wait until the worker is ready
-
 await model.setModelConfig("frontPropagationScript");
+
+// 2. Define the mesh
 await model.setMeshConfig({
   meshDimension: "2D",
   elementOrder: "quadratic",
@@ -28,18 +30,22 @@ await model.setMeshConfig({
   maxX: 4,
   maxY: 2,
 });
-await model.addBoundaryCondition("0", ["constantValue", 0]);
-await model.addBoundaryCondition("1", ["constantValue", 0]);
-await model.addBoundaryCondition("2", ["zeroGradient"]);
-await model.addBoundaryCondition("3", ["constantValue", 0]);
-await model.setSolverMethod("lusolve");
 
+// 3. Apply the boundary conditions
+await model.addBoundaryCondition("0", ["constantValue", 0]); // Bottom: initial front
+await model.addBoundaryCondition("1", ["constantValue", 0]); // Left: initial front
+await model.addBoundaryCondition("2", ["zeroGradient"]); // Top
+await model.addBoundaryCondition("3", ["constantValue", 0]); // Right: initial front
+
+// 4. Solve
+await model.setSolverMethod("lusolve");
 const result = await model.solve();
 
-// Plot in the main thread
+// 5. Plot the results in the main thread
 const modelInfo = await model.getModelInfo();
 plotSolution(modelInfo, result, "contour", "resultsCanvas");
 
+// 6. Terminate the worker
 model.terminate();
 ```
 

docs/getting-started/first-simulation.md

@@ -38,13 +38,11 @@ Save the following as an `.html` file and open it through a local web server (fo
     <script type="module">
       import { FEAScriptModel, plotSolution } from "https://core.feascript.com/dist/feascript.esm.js";
 
-      // 1. Create a model
+      // 1. Create the model and select the physics
       const model = new FEAScriptModel();
-
-      // 2. Select the physical model
       model.setModelConfig("heatConductionScript");
 
-      // 3. Define the mesh: 8 × 4 quadratic quadrilaterals on 0 ≤ x ≤ 4, 0 ≤ y ≤ 2
+      // 2. Define the mesh: 8 × 4 quadratic quadrilaterals on 0 ≤ x ≤ 4, 0 ≤ y ≤ 2
       model.setMeshConfig({
         meshDimension: "2D",
         elementOrder: "quadratic",
@@ -54,17 +52,17 @@ Save the following as an `.html` file and open it through a local web server (fo
         maxY: 2,
       });
 
-      // 4. Apply the boundary conditions
+      // 3. Apply the boundary conditions
       model.addBoundaryCondition("0", ["constantTemperature", 200]); // Bottom
       model.addBoundaryCondition("1", ["symmetry"]); // Left
       model.addBoundaryCondition("2", ["convection", 1, 20]); // Top
       model.addBoundaryCondition("3", ["constantTemperature", 200]); // Right
 
-      // 5. Select the solver and solve
+      // 4. Solve
       model.setSolverMethod("lusolve");
       const result = model.solve();
 
-      // 6. Plot the temperature
+      // 5. Plot the results
       plotSolution(model, result, "contour", "resultsCanvas");
     </script>
   </body>
@@ -73,18 +71,17 @@ Save the following as an `.html` file and open it through a local web server (fo
 
 ## Step by step
 
-1. **Create a model.** [`FEAScriptModel`](../api/feascript-model.md) holds the whole configuration of the
-   simulation.
-2. **Select the physical model.** `setModelConfig` takes the model name and, optionally, its coefficients. Here
-   the defaults of the [heat conduction](../models/heat-conduction.md) model are used ($k = 1$, $Q = 0$).
-3. **Define the mesh.** `setMeshConfig` generates a structured mesh from the size of the domain and the number of
+1. **Create the model and select the physics.** [`FEAScriptModel`](../api/feascript-model.md) holds the whole
+   configuration of the simulation. `setModelConfig` takes the model name and, optionally, its coefficients. Here the
+   defaults of the [heat conduction](../models/heat-conduction.md) model are used ($k = 1$, $Q = 0$).
+2. **Define the mesh.** `setMeshConfig` generates a structured mesh from the size of the domain and the number of
    elements (see [Mesh](../fundamentals/mesh.md)), or uses a mesh imported from Gmsh.
-4. **Apply the boundary conditions.** The sides of a structured 2D mesh have the keys `"0"` (bottom), `"1"`
+3. **Apply the boundary conditions.** The sides of a structured 2D mesh have the keys `"0"` (bottom), `"1"`
    (left), `"2"` (top) and `"3"` (right). See [Boundary conditions](../fundamentals/boundary-conditions.md) for the
    condition types of every model.
-5. **Solve.** `setSolverMethod` is optional (the default is `"lusolve"`). `solve()` returns the
+4. **Solve.** `setSolverMethod` is optional (the default is `"lusolve"`). `solve()` returns the
    `solutionVector` and the `nodesCoordinates` of the mesh (see [Solvers](../fundamentals/solvers.md)).
-6. **Post-process.** [`plotSolution`](../api/visualization.md) draws the result in the `resultsCanvas` element.
+5. **Plot the results.** [`plotSolution`](../api/visualization.md) draws the result in the `resultsCanvas` element.
 
 ## The result
 
@@ -105,8 +102,11 @@ for (let i = 0; i < solutionVector.length; i++) {
 ```javascript
 import { FEAScriptModel } from "feascript";
 
+// 1. Create the model and select the physics
 const model = new FEAScriptModel();
 model.setModelConfig("heatConductionScript");
+
+// 2. Define the mesh
 model.setMeshConfig({
   meshDimension: "2D",
   elementOrder: "quadratic",
@@ -115,12 +115,17 @@ model.setMeshConfig({
   maxX: 4,
   maxY: 2,
 });
-model.addBoundaryCondition("0", ["constantTemperature", 200]);
-model.addBoundaryCondition("1", ["symmetry"]);
-model.addBoundaryCondition("2", ["convection", 1, 20]);
-model.addBoundaryCondition("3", ["constantTemperature", 200]);
 
+// 3. Apply the boundary conditions
+model.addBoundaryCondition("0", ["constantTemperature", 200]); // Bottom
+model.addBoundaryCondition("1", ["symmetry"]); // Left
+model.addBoundaryCondition("2", ["convection", 1, 20]); // Top
+model.addBoundaryCondition("3", ["constantTemperature", 200]); // Right
+
+// 4. Solve
 const { solutionVector, nodesCoordinates } = model.solve();
+
+// 5. Print the results
 console.log(`Number of nodes: ${nodesCoordinates.nodesXCoordinates.length}`);
 console.log("Temperature:", solutionVector);
 ```

docs/models/creeping-flow.md

@@ -72,28 +72,28 @@ The lid-driven cavity: a $4 \times 2$ cavity whose top wall moves with unit velo
 ```javascript
 import { FEAScriptModel } from "feascript";
 
+// 1. Create the model and select the physics
 const model = new FEAScriptModel();
-
 model.setModelConfig("creepingFlowScript");
 
-const numElementsX = 12;
-const numElementsY = 6;
+// 2. Define the mesh
 model.setMeshConfig({
   meshDimension: "2D",
   elementOrder: "quadratic", // Required (Taylor–Hood elements)
-  numElementsX,
-  numElementsY,
+  numElementsX: 12,
+  numElementsY: 6,
   maxX: 4,
   maxY: 2,
 });
 
+// 3. Apply the boundary conditions
 model.addBoundaryCondition("0", ["constantVelocity", 0, 0]); // Bottom wall
 model.addBoundaryCondition("1", ["constantVelocity", 0, 0]); // Left wall
 model.addBoundaryCondition("2", ["constantVelocity", 1, 0]); // Moving lid
 model.addBoundaryCondition("3", ["constantVelocity", 0, 0]); // Right wall
 
+// 4. Solve
 model.setSolverMethod("lusolve");
-
 const { solutionVector, nodesCoordinates } = model.solve();
 ```
 

docs/models/euler-bernoulli-beam.md

@@ -61,31 +61,33 @@ from J. N. Reddy, _An Introduction to the Finite Element Method_, 3rd ed., McGra
 ```javascript
 import { FEAScriptModel } from "feascript";
 
+// 1. Create the model and select the physics
 const model = new FEAScriptModel();
-
 model.setModelConfig("eulerBernoulliBeamScript", {
   coefficientFunctions: {
     EI: (x) => 2.0e6, // Bending stiffness
     q: (x) => (x <= 5 ? -1000 : 0), // Distributed load over the first span
   },
 });
 
+// 2. Define the mesh
 model.setMeshConfig({
   meshDimension: "1D",
   elementOrder: "linear", // Required
   numElementsX: 2, // Nodes "1", "2", "3" at x = 0, 5, 10
   maxX: 10,
 });
 
+// 3. Apply the boundary conditions
 model.addBoundaryCondition("1", [["fixed"]]); // Clamped end
 model.addBoundaryCondition("2", [["pinned"], ["moment", 1250]]); // Roller with an applied moment
 model.addBoundaryCondition("3", [
   ["spring", 200],
   ["force", -2500],
 ]); // Spring support with an applied force
 
+// 4. Solve
 model.setSolverMethod("lusolve");
-
 const { solutionVector, nodesCoordinates } = model.solve();
 ```
 

docs/models/front-propagation.md

@@ -78,10 +78,11 @@ The solidification front in a $4 \times 2$ domain, starting from the bottom, lef
 ```javascript
 import { FEAScriptModel } from "feascript";
 
+// 1. Create the model and select the physics
 const model = new FEAScriptModel();
-
 model.setModelConfig("frontPropagationScript");
 
+// 2. Define the mesh
 model.setMeshConfig({
   meshDimension: "2D",
   elementOrder: "quadratic",
@@ -91,13 +92,14 @@ model.setMeshConfig({
   maxY: 2,
 });
 
+// 3. Apply the boundary conditions
 model.addBoundaryCondition("0", ["constantValue", 0]); // Bottom: initial front
 model.addBoundaryCondition("1", ["constantValue", 0]); // Left: initial front
 model.addBoundaryCondition("2", ["zeroGradient"]); // Top
 model.addBoundaryCondition("3", ["constantValue", 0]); // Right: initial front
 
+// 4. Solve
 model.setSolverMethod("lusolve");
-
 const { solutionVector, nodesCoordinates } = model.solve();
 ```
 

docs/models/general-form-pde.md

@@ -76,8 +76,8 @@ The equation $\frac{d^2u}{dx^2} - 10 \frac{du}{dx} = 10 \, e^{-200 (x - 0.5)^2}$
 ```javascript
 import { FEAScriptModel } from "feascript";
 
+// 1. Create the model and select the physics
 const model = new FEAScriptModel();
-
 model.setModelConfig("generalFormPDEScript", {
   coefficientFunctions: {
     A: (x) => 1, // Diffusion coefficient
@@ -87,18 +87,20 @@ model.setModelConfig("generalFormPDEScript", {
   },
 });
 
+// 2. Define the mesh
 model.setMeshConfig({
   meshDimension: "1D",
   elementOrder: "quadratic",
   numElementsX: 20,
   maxX: 1.0,
 });
 
+// 3. Apply the boundary conditions
 model.addBoundaryCondition("0", ["constantValue", 1]); // Left, x = 0
 model.addBoundaryCondition("1", ["zeroGradient"]); // Right, x = 1
 
+// 4. Solve
 model.setSolverMethod("lusolve");
-
 const { solutionVector, nodesCoordinates } = model.solve();
 ```
 

docs/models/heat-conduction.md

@@ -52,15 +52,16 @@ The model is linear, so the system is solved once.
 ```javascript
 import { FEAScriptModel } from "feascript";
 
+// 1. Create the model and select the physics
 const model = new FEAScriptModel();
-
 model.setModelConfig("heatConductionScript", {
   coefficientFunctions: {
     thermalConductivity: 1, // k
     heatSource: 0, // Q
   },
 });
 
+// 2. Define the mesh
 model.setMeshConfig({
   meshDimension: "2D",
   elementOrder: "quadratic",
@@ -70,13 +71,14 @@ model.setMeshConfig({
   maxY: 2,
 });
 
+// 3. Apply the boundary conditions
 model.addBoundaryCondition("0", ["constantTemperature", 200]); // Bottom
 model.addBoundaryCondition("1", ["symmetry"]); // Left
 model.addBoundaryCondition("2", ["convection", 1, 20]); // Top
 model.addBoundaryCondition("3", ["constantTemperature", 200]); // Right
 
+// 4. Solve
 model.setSolverMethod("lusolve");
-
 const { solutionVector, nodesCoordinates } = model.solve();
 ```
 

static/llms.txt

@@ -22,12 +22,13 @@ Exports: `FEAScriptModel`, `importGmshMesh`, `logSystem`, `printVersion` (a stri
 Every simulation follows the same steps:
 
 ```javascript
+// 1. Create the model and select the physics
 const model = new FEAScriptModel();
-
 model.setModelConfig("heatConductionScript", {
   coefficientFunctions: { thermalConductivity: 1, heatSource: 0 },
 });
 
+// 2. Define the mesh
 model.setMeshConfig({
   meshDimension: "2D",
   elementOrder: "quadratic",
@@ -37,13 +38,14 @@ model.setMeshConfig({
   maxY: 2,
 });
 
+// 3. Apply the boundary conditions
 model.addBoundaryCondition("0", ["constantTemperature", 200]); // Bottom
 model.addBoundaryCondition("1", ["symmetry"]); // Left
 model.addBoundaryCondition("2", ["convection", 1, 20]); // Top: h = 1, T∞ = 20
 model.addBoundaryCondition("3", ["constantTemperature", 200]); // Right
 
+// 4. Solve
 model.setSolverMethod("lusolve"); // Optional, "lusolve" is the default
-
 const { solutionVector, nodesCoordinates } = model.solve();
 // solutionVector: Float64Array
 // nodesCoordinates: { nodesXCoordinates: number[], nodesYCoordinates: number[] }