Skip to content

Commit 2096259

Browse files
janbenderjanbender
committed
- added heat equation example
1 parent 96b956e commit 2096259

File tree

2 files changed

+353
-0
lines changed

2 files changed

+353
-0
lines changed

‎examples/heat_equation.html‎

Lines changed: 352 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,352 @@
1+
<!doctype html>
2+
<html class="no-js" lang="en">
3+
<head>
4+
<meta charset="utf-8">
5+
<style>
6+
body {font-family: Helvetica, sans-serif;}
7+
table {background-color:#CCDDEE;text-align:left}
8+
</style>
9+
<script type="text/x-mathjax-config">
10+
MathJax.Hub.Config({
11+
extensions: ["tex2jax.js"],
12+
jax: ["input/TeX", "output/HTML-CSS"],
13+
tex2jax: {
14+
inlineMath: [ ['$','$'], ["\\(","\\)"] ],
15+
displayMath: [ ['$$','$$'], ["\\[","\\]"] ],
16+
processEscapes: true
17+
},
18+
"HTML-CSS": { fonts: ["TeX"] }
19+
});
20+
</script>
21+
<script type="text/javascript" aync src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.4/MathJax.js"></script>
22+
<title>Heat Equation — 2D Diffusion Example</title>
23+
</head>
24+
<body>
25+
<main>
26+
<h1 style="text-align:center">Heat Equation — 2D Diffusion Example</h1>
27+
<table style="align_center;border-radius: 20px;padding: 20px;margin:auto">
28+
<col width="1100">
29+
<col width="400">
30+
<tr>
31+
<td>
32+
<canvas id="simCanvas" width="1024" height="768" style="border:2px solid #000000;border-radius: 20px;background-color:#EEEEEE">Your browser does not support the HTML5 canvas tag.</canvas>
33+
</td>
34+
<td>
35+
<table>
36+
<col width="180" style="padding-right:10px">
37+
<col width="100">
38+
<tr>
39+
<td><label>Current time</label></td>
40+
<td><span id="time">0.00</span> s</td>
41+
</tr>
42+
<tr>
43+
<td><label>Time per sim. step</label></td>
44+
<td><span id="timePerStep">0.00</span> ms</td>
45+
</tr>
46+
<tr>
47+
<td><label for="timeStepSizeInput">Time step size</label></td>
48+
<td><input onchange="gui.sim.timeStepSize=parseFloat(value)" id="timeStepSizeInput" type="number" value="0.005" step="0.001"></td>
49+
</tr>
50+
<tr>
51+
<td><label for="diffusivityInput">Diffusivity constant</label></td>
52+
<td><input onchange="gui.sim.diffusivity=parseFloat(value)" id="diffusivityInput" type="number" value="20" step="1"></td>
53+
</tr>
54+
<tr>
55+
<td></td>
56+
<td><button onclick="gui.restart()" type="button" id="restart">Restart</button></td>
57+
</tr>
58+
<tr>
59+
<td></td>
60+
<td><button onclick="gui.doPause()" type="button" id="Pause">Pause</button></td>
61+
</tr>
62+
</table>
63+
</td>
64+
</tr>
65+
<tr><td>
66+
<h2>Heat diffusion example</h2>
67+
<p>
68+
This example demonstrates a simple explicit finite-difference simulation of the 2D heat equation (diffusion).
69+
The simulation stores a temperature field on a lower-resolution grid for performance; the canvas displays the field upscaled.
70+
</p>
71+
72+
<h3>Initial condition</h3>
73+
<p>The temperature is initially zero everywhere except a centered hot box (set to 100). You can restart the simulation with the <em>Restart</em> button.</p>
74+
75+
<h3>Model (continuous)</h3>
76+
<p>The heat equation solved is
77+
\[ \frac{\partial T}{\partial t} = D \nabla^2 T, \]
78+
where \(D\) is the diffusivity.</p>
79+
80+
<h3>Numerical method</h3>
81+
<p>Space is discretized on a uniform Cartesian grid with spacing \(\Delta x\) (we assume \(\Delta x = \Delta y\) here). Let \(T_{i,j}^n\) denote the temperature at grid cell \((i,j)\) at time step \(n\). The continuous PDE is
82+
\[ \frac{\partial T}{\partial t} = D \nabla^2 T, \]
83+
which we discretize using a 5-point central difference for the Laplacian and forward-Euler in time (explicit):
84+
85+
\[ \nabla^2 T\big|_{i,j} \approx \frac{T_{i+1,j} + T_{i-1,j} + T_{i,j+1} + T_{i,j-1} - 4T_{i,j}}{\Delta x^2}. \]
86+
The time-stepping update is then
87+
\[ T_{i,j}^{n+1} = T_{i,j}^n + \Delta t \; D \; \frac{T_{i+1,j}^n + T_{i-1,j}^n + T_{i,j+1}^n + T_{i,j-1}^n - 4T_{i,j}^n}{\Delta x^2}. \]
88+
89+
<h3>Note on stability</h3>
90+
<p>Because the integrator is explicit, the time step must be small enough to satisfy the CFL-like condition for stability roughly: \(\Delta t \le \frac{\Delta x^2}{4D}\). If you increase \(D\) or \(\Delta t\) you may see the simulation blow up or become noisy.</p>
91+
92+
</td></tr>
93+
</table>
94+
95+
</main>
96+
97+
<script id="simulation_code" type="text/javascript">
98+
class Simulation
99+
{
100+
constructor(resX, resY, scaleFactor)
101+
{
102+
this.timeStepSize = 0.005;
103+
this.time = 0;
104+
// scaling factor between simulation resolution and canvas resolution
105+
// (simulation runs at lower resolution for performance reasons)
106+
this.scale = 1.0/scaleFactor
107+
108+
this.resX = Math.floor(resX * this.scale);
109+
this.resY = Math.floor(resY * this.scale);
110+
this.temperature = new Array(this.resX*this.resY).fill(0.0)
111+
// diffusivity constant
112+
// (higher values lead to faster heat diffusion)
113+
this.diffusivity = 20.0;
114+
115+
this.init()
116+
}
117+
118+
init()
119+
{
120+
// init temperature field: zero everywhere, non-zero inside a centered box
121+
const boxW = 50; // box width in simulation cells
122+
const boxH = 50; // box height in simulation cells
123+
const cx = Math.floor(this.resX / 2);
124+
const cy = Math.floor(this.resY / 2);
125+
const x0 = Math.max(0, cx - Math.floor(boxW / 2));
126+
const x1 = Math.min(this.resX - 1, cx + Math.floor((boxW - 1) / 2));
127+
const y0 = Math.max(0, cy - Math.floor(boxH / 2));
128+
const y1 = Math.min(this.resY - 1, cy + Math.floor((boxH - 1) / 2));
129+
for (let i = y0; i <= y1; i++)
130+
for (let j = x0; j <= x1; j++)
131+
this.temperature[i * this.resX + j] = 100.0;
132+
}
133+
134+
// simulation step
135+
simulationStep()
136+
{
137+
let temp = this.temperature.slice();
138+
for (let i = 1; i < this.resY-1; i++)
139+
{
140+
for (let j = 1; j < this.resX-1; j++)
141+
{
142+
// compute Laplace operator using 5-point stencil
143+
let right = temp[i*this.resX+(j+1)];
144+
let center = temp[i*this.resX+j];
145+
let left = temp[i*this.resX+(j-1)];
146+
let up = temp[(i-1)*this.resX+j];
147+
let down = temp[(i+1)*this.resX+j];
148+
149+
const dx = 1.0; // set physical spacing between grid points
150+
const invDx2 = 1.0 / (dx * dx);
151+
let laplace = (up + down + left + right - 4.0 * center) * invDx2;;
152+
this.temperature[i*this.resX+j] += this.diffusivity * laplace * this.timeStepSize;
153+
// enforce temperature limits
154+
this.temperature[i*this.resX+j] = Math.min(100.0, this.temperature[i*this.resX+j]);
155+
this.temperature[i*this.resX+j] = Math.max(0.0, this.temperature[i*this.resX+j]);
156+
}
157+
}
158+
// update simulation time
159+
this.time = this.time + this.timeStepSize;
160+
}
161+
}
162+
163+
class GUI
164+
{
165+
constructor()
166+
{
167+
this.canvas = document.getElementById("simCanvas");
168+
this.c = this.canvas.getContext("2d");
169+
this.requestID = -1;
170+
171+
this.origin = { x : this.canvas.width / 2, y : this.canvas.height/2};
172+
this.scaleFactor = 2;
173+
this.mouseButtonDown = false
174+
175+
// register mouse event listeners (zoom/selection)
176+
this.canvas.addEventListener("mousedown", this.mouseDown.bind(this), false);
177+
this.canvas.addEventListener("mousemove", this.mouseMove.bind(this), false);
178+
this.canvas.addEventListener("mouseup", this.mouseUp.bind(this), false);
179+
}
180+
181+
// set simulation parameters from GUI and start mainLoop
182+
restart()
183+
{
184+
window.cancelAnimationFrame(this.requestID);
185+
delete this.sim;
186+
this.sim = new Simulation(this.canvas.width, this.canvas.height, this.scaleFactor);
187+
188+
this.timeSum = 0.0;
189+
this.counter = 0;
190+
191+
this.sim.timeStepSize = parseFloat(document.getElementById('timeStepSizeInput').value);
192+
this.sim.diffusivity = parseFloat(document.getElementById('diffusivityInput').value);
193+
194+
this.mainLoop();
195+
}
196+
197+
hsvToRgb(h, s, v)
198+
{
199+
let i = Math.floor(h * 6);
200+
let f = h * 6 - i;
201+
let p = v * (1 - s);
202+
let q = v * (1 - f * s);
203+
let t = v * (1 - (1 - f) * s);
204+
205+
let r = 0;
206+
let g = 0;
207+
let b = 0;
208+
switch (i % 6)
209+
{
210+
case 0: r = v, g = t, b = p; break;
211+
case 1: r = q, g = v, b = p; break;
212+
case 2: r = p, g = v, b = t; break;
213+
case 3: r = p, g = q, b = v; break;
214+
case 4: r = t, g = p, b = v; break;
215+
case 5: r = v, g = p, b = q; break;
216+
}
217+
return [Math.round(r * 255), Math.round(g * 255), Math.round(b * 255)];
218+
}
219+
220+
draw()
221+
{
222+
this.c.clearRect(0, 0, this.canvas.width, this.canvas.height);
223+
224+
let imgData = this.c.getImageData(0, 0, this.canvas.width, this.canvas.height);
225+
for(let py = 0; py < this.sim.resY; py++)
226+
{
227+
for(let px = 0; px < this.sim.resX; px++)
228+
{
229+
let value = this.sim.temperature[py * this.sim.resX + px]/100.0
230+
value = Math.min(0.666, value);
231+
value = Math.max(0.0, value);
232+
// blue (cold) to red (hot)
233+
let color = this.hsvToRgb(0.666 - value, 1.0, 1.0);
234+
235+
for (let i=0; i < this.scaleFactor; i++)
236+
for (let j=0; j < this.scaleFactor; j++)
237+
{
238+
let ppx = this.scaleFactor*px + i;
239+
let ppy = this.scaleFactor*py + j;
240+
if (ppx < imgData.width && ppy < imgData.height)
241+
{
242+
imgData.data[4 * (ppy * imgData.width + ppx) ] = color[0];
243+
imgData.data[4 * (ppy * imgData.width + ppx) + 1] = color[1];
244+
imgData.data[4 * (ppy * imgData.width + ppx) + 2] = color[2];
245+
imgData.data[4 * (ppy * imgData.width + ppx) + 3] = 255;
246+
}
247+
}
248+
}
249+
}
250+
this.c.putImageData(imgData, 0, 0);
251+
}
252+
253+
mainLoop()
254+
{
255+
// perform multiple sim steps per render step
256+
for (let i=0; i < 8; i++)
257+
{
258+
let t0 = performance.now();
259+
260+
this.sim.simulationStep();
261+
262+
let t1 = performance.now();
263+
this.timeSum += t1 - t0;
264+
this.counter += 1;
265+
266+
if (this.counter % 50 == 0)
267+
{
268+
this.timeSum /= this.counter;
269+
document.getElementById("timePerStep").innerHTML = this.timeSum.toFixed(2);
270+
this.timeSum = 0.0;
271+
this.counter = 0;
272+
}
273+
document.getElementById("time").innerHTML = this.sim.time.toFixed(2);
274+
}
275+
276+
this.draw();
277+
if (!this.pause)
278+
this.requestID = window.requestAnimationFrame(this.mainLoop.bind(this));
279+
}
280+
281+
doPause()
282+
{
283+
this.pause = !this.pause;
284+
if (!this.pause)
285+
this.mainLoop();
286+
}
287+
288+
mouseDown(event)
289+
{
290+
// left mouse button down
291+
if (event.which == 1) {
292+
this.mouseButtonDown = true;
293+
let mousePos = this.getMousePos(this.canvas, event);
294+
let x = Math.floor(mousePos.x / this.scaleFactor);
295+
let y = Math.floor(mousePos.y / this.scaleFactor);
296+
this.lastMouse = { x: x, y: y };
297+
}
298+
}
299+
300+
getMousePos(canvas, event)
301+
{
302+
let rect = canvas.getBoundingClientRect();
303+
return {
304+
x: event.clientX - rect.left,
305+
y: event.clientY - rect.top
306+
};
307+
}
308+
309+
mouseMove(event)
310+
{
311+
let mousePos = this.getMousePos(this.canvas, event);
312+
if (this.mouseButtonDown)
313+
{
314+
let x = Math.floor(mousePos.x / this.scaleFactor);
315+
let y = Math.floor(mousePos.y / this.scaleFactor);
316+
if (!this.lastMouse) this.lastMouse = { x: x, y: y };
317+
let lx = this.lastMouse.x;
318+
let ly = this.lastMouse.y;
319+
let dx = x - lx;
320+
let dy = y - ly;
321+
let steps = Math.max(Math.abs(dx), Math.abs(dy), 1);
322+
for (let s = 0; s <= steps; s++) {
323+
let xi = Math.round(lx + dx * s / steps);
324+
let yi = Math.round(ly + dy * s / steps);
325+
for (let ddy = -1; ddy <= 1; ddy++) {
326+
for (let ddx = -1; ddx <= 1; ddx++) {
327+
let xj = xi + ddx;
328+
let yj = yi + ddy;
329+
if (xj > 0 && xj < this.sim.resX-1 && yj > 0 && yj < this.sim.resY-1)
330+
this.sim.temperature[yj * this.sim.resX + xj] = 100.0;
331+
}
332+
}
333+
}
334+
this.lastMouse.x = x;
335+
this.lastMouse.y = y;
336+
//this.requestID = window.requestAnimationFrame(this.mainLoop.bind(this));
337+
}
338+
}
339+
340+
mouseUp(event)
341+
{
342+
this.mouseButtonDown = false;
343+
this.lastMouse = null;
344+
}
345+
}
346+
347+
gui = new GUI();
348+
gui.restart();
349+
</script>
350+
351+
</body>
352+
</html>

‎index.md‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -11,6 +11,7 @@ I hope the examples and their documentation help you to learn something about th
1111
# General
1212

1313
* [Newton's method](examples/Newton_solver.html)
14+
* [Heat equation](examples/heat_equation.html)
1415

1516
# Particles
1617

0 commit comments

Comments
 (0)