偏微分方程式、差分法
Fitzhugh-Nagumo方程式
dU/dt = d/dx(dU/dx) + (a-U)(U-1)U-V (d/dx(dU/dx)はUのxに対する二回微分)
dV/dt =e(bU-gV)
a=0.1, b=0.5, g=1.0, e=0.01
初期条件
(U,V)=(1,0) if x=0
(U,V)=(0,0) if x>0
境界条件
dU/dx=0 at x=0
dV/dx=0 at x=0
を差分して陽解法で解くと添付の答えになるらしいのですが、自分で解いたところ添付の結果のようになりました。答えと一致しないため、プログラム上で何を間違えているのかご指摘頂けると助かります。
#include <stdio.h>
#include <stdlib.h> //exsit()で必要
#include <math.h>
int main(){
double a,b,g,e;
double dt,dx;
int x,t;
double **u,**v;
int i,j,k;
a=0.1;b=0.5;g=1.0;e=0.01;
x=100;t=5000.;
dt=0.001;dx=0.1;
FILE *fp;
if((fp = fopen("ResultNagumo.txt","w"))==NULL){
printf("Can't open file\n");
exit(2);
}
u = (double**)malloc(sizeof(double*)*t);
v = (double**)malloc(sizeof(double*)*t);
for(i=0;i<t;i++){
u[i]=(double*)malloc(sizeof(double)*x);
v[i]=(double*)malloc(sizeof(double)*x);
}
//初期値
u[0][0]=1.0;
v[0][0]=0.0;
for(i=1;i<x;i++){
u[0][i]=0.0;
v[0][i]=0.0;
}
//差分計算
for(i=0;i<t-1;i++){
for(j=1;j<x-1;j++){
u[i+1][j] = u[i][j] + dt*((u[i][j+1]-2*u[i][j]+u[i][j-1])/pow(dx,2)+(a-u[i][j])*(u[i][j]-1)*u[i][j]-v[i][j]);
v[i+1][j] = v[i][j] + dt*(e*(b*u[i][j]-g*v[i][j]));
}
//境界条件
u[i+1][0]=u[i+1][1];
v[i+1][0]=v[i+1][1];
u[i+1][x-1]=u[i+1][x-2];
v[i+1][x-1]=v[i+1][x-2];
}
//結果の出力
for(i=0;i<t;i++){
if((i%100)==0){
fprintf(fp,"%d\n",i);
for(j=0;j<x;j++){
fprintf(fp,"%2.4e,",u[i][j]);
}
}
}
for(i=0;i<t;i++){
free(u[i]);
free(v[i]);
}
free(u);free(v);
fclose(fp);
getchar();
return 0;
}
お礼
返事が遅れて申し訳ありません。 最初からWikipediaで検索するべきでしたね 大変ありがとうございました。