Processingでシミュレーション～熱伝導

はじめに

Processingを使って熱伝導をシミュレーションします。

熱伝導方程式

時間とともに変化する温度の分布を教えてくれる方程式です。

温度$T(x,t)$は時間と位置の関数になります。

\frac{\partial T}{\partial t} = k \frac{\partial ^2 T}{\partial x^2}

差分化・離散化

差分に変換します。

\frac{1}{Δt}(T(x,t+Δt) - T(x,t)) = k\frac{1}{Δx^2}(T(x+Δx,t)-2T(x,t)+T(x-Δx,t))

離散化
Δt=1, t=0,1,2,...
Δx=1, x=0,1,2,...,len
とします。

T(x,t+1) - T(x,t) = k(T(x+1,t)-2T(x,t)+T(x-1,t))

t=0のとき

T(x,1) = T(x,0) + k(T(x+1,0)-2T(x,0)+T(x-1,0))

$T(x,0)$は初期条件です。時刻0での温度分布です。
$T(x,1)$は1時刻進んだ温度分布です。
xが0~lenの範囲をとると、$x(-1,t),x(len+1,t)$のデータが必要になります。これは境界条件になります。

Processingで可視化

サーモグラフィーのような色付け

HSBカラーモードを使ってHue値を0%~68%まで使うと0%が赤、68%が青色になります。
高温を赤、低温を青で色づけします。

void setup() 
{
  size(680, 120);
  int step = 100;
  PImage img = createImage(step, 1, RGB);
  colorMode(HSB, step*100/68, 100, 100);
  noSmooth();
  for (int i=0; i<step; i++) {
    img.pixels[i] = color(i+1, 100, 100);
  }
  img.updatePixels();
  image(img, 0, 0, width, height);
}

シミュレーション

t1は現在の温度分布、t0は1時刻前の温度分布です。
xminは左端の温度、xmaxは右端の温度です。初期条件はすべて温度0です。
右が高温、左が低温になっており、温度が右から左に移動します。

float k = 0.5;
float[] t0;
float[] t1;
float xmin;
float xmax;
PImage img;
int step = 100;

void setup() {
  size(640, 100);
  init();
  noSmooth();
  frameRate(60);
}

void draw() {
  for (int i=0; i<step; i++) {
    img.pixels[i] = color(255-t1[i], 100, 100);
  }
  img.updatePixels();
  image(img, 0, 0, width, height);
  for (int i=0; i<10; i++) // speed up
    updateTime();
}

void init()
{ 
  t0 = new float[step];
  t1 = new float[step];
  xmin = 0;
  xmax = 255;
  img = createImage(step, 1, RGB);
  colorMode(HSB, 255*100/68, 100, 100);
}

void updateTime()
{
  for (int i=0; i<step; i++) {
    t0[i] = t1[i];
  }

  t1[0] = t0[0] + k * (xmin - 2 * t0[0] + t0[1]);
  for (int i=1; i<step-1; i++) {
    t1[i] = t0[i] + k * (t0[i-1] - 2 * t0[i] + t0[i+1]);
  }
  t1[step-1] = t0[step-1] + k * (t0[step-2] - 2 * t0[step-1] + xmax);
}