Schnakenberg kinetics + diffusion
1次元空間
'13 4/22

反応項のみ （空間一様解）

微分方程式の定義

deq ={u'[t]==k1 - k2 u[t] + k3 u[t]^2 v[t],
	v'[t]==k4 - k3 u[t]^2 v[t],u[0]==u0,v[0]==v0}

パラメータと初期条件の設定

para = {k1->0.2, k2->1, k3->1, k4->0.5}

deqN = deq/.{u0->0.1,v0->0.2}/.para

微分方程式の数値解

tEnd = 200;
sol = NDSolve[deqN,{u[t],v[t]},{t,0,tEnd}]

gDyna=ParametricPlot[ Evaluate[{u[t],v[t]}/.sol],{t,0,tEnd},PlotRange->{{0,3},{0,3}},
AxesLabel->{"u","v"},AspectRatio->1]

Plot[ Evaluate[{u[t],v[t]}/.sol],{t,0,tEnd},PlotRange->All,AxesLabel->{"t","u, v"}]

相平面解析

fn1 = k1 - k2 u + k3 u^2 v;
fn2 = k4 - k3 u^2 v;

para

{fn1N, fn2N} = {fn1, fn2}/.para

g1 = ContourPlot[fn1N==0,{u,0,5},{v,0,5},ContourStyle->RGBColor[1,0,0]]

g2 =ContourPlot[fn2N==0,{u,0,5},{v,0,5},ContourStyle->RGBColor[0,0,1]]

Show[g1,g2]

fn1N

fn2N

FindRoot[{fn1N==0,fn2N==0},{u,0.5},{v,0.9}]

Show[gDyna,g1,g2]

反応項 + 拡散（Cプログラムによる数値計算の視覚化）

SetDirectory["/Users/takasu/Library/Developer/Xcode/DerivedData/reaction_diffusion-drcxqrxbavdoubdxqelreislhgkm/Build/Products/Debug/"]

drawFigure[x1_List,x2_List,opt___]:=Block[{g1,g2},
g1=ListPlot[x1,Joined->True,PlotStyle->RGBColor[1,0,0],opt,DisplayFunction->Identity];g2=ListPlot[x2,Joined->True,PlotStyle->RGBColor[0,0,1],opt,DisplayFunction->Identity];Return[Show[g1,g2,DisplayFunction->$DisplayFunction]]
]

データの読み込みと視覚化

size = 101; (* size の値に注意 *)
datau = ReadList["RD-1dim-u",{Real,Real}];
datau = Partition[datau,size];
Length[datau]

(* size の値に注意 *)
datav = ReadList["RD-1dim-v",{Real,Real}];
datav = Partition[datav,size];
Length[datav]

drawFigure[datau[[1]], datav[[1]],PlotRange->All,AxesLabel->{"x","u, v"}]

drawFigure[Last[datau], Last[datav],PlotRange->All,AxesLabel->{"x","u, v"}]

glist = {};
Do[
g = drawFigure[ datau[[i]],datav[[i]] ,PlotRange->{0,1.2},Joined->True,AxesLabel->{"x","u, v"}];
AppendTo[glist,g],{i,1,Length[datau],10}
]

Show[glist,PlotRange->All]

ListAnimate[glist]

<<Graphics`Graphics3D`

Show[StackGraphics[glist],BoxRatios->{2,2,1},PlotRange->{0,1.5}]

1