-
Notifications
You must be signed in to change notification settings - Fork 4
/
Bifurcation_Turing2.m
33 lines (27 loc) · 1.13 KB
/
Bifurcation_Turing2.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
function out = Bifurcation_Turing2(k1,k2,k3,k4,k5,k6,k7,nc,Dx,Dy,qspace)
% Bifurcation_Turing2
% Copyright (c) Microsoft Corporation. All rights reserved.
% Licensed under the MIT License.
Xss=max([(-k1*k4^2*k6-k1*k4*k5*k7+sqrt(k4)*sqrt(k4*k6+k5*k7)*sqrt(k1^2*k4^2*k6+...
8*k3*k4*k5^2*k6+k1^2*k4*k5*k7+8*k3*k5^3*k7+8*k3*k5^2*k6*k7*nc))/(4*(k3*...
k4*k5*k6+k3*k5^2*k7)),(-k1*k4^2*k6-k1*k4*k5*k7-sqrt(k4)*sqrt(k4*k6+k5*k7)*sqrt(k1^2*k4^2*k6+...
8*k3*k4*k5^2*k6+k1^2*k4*k5*k7+8*k3*k5^3*k7+8*k3*k5^2*k6*k7*nc))/(4*(k3*...
k4*k5*k6+k3*k5^2*k7))]);
Yss=k4/k5;Wss=k3*Xss^2/k2+k1*Xss*Yss/k2;Css=k7*nc/(k7+k6*Yss);
J=[-k1*Yss-4*k3*Xss,-k1*Xss-k6^2*k7*nc*Yss/(k7+k6*Yss)^2+k6*k7*nc/(k7+k6*Yss);
k1*Yss+4*k3*Xss,k1*Xss-k5+k6^2*k7*nc*Yss/(k7+k6*Yss)^2-k6*k7*nc/(k7+k6*Yss)];
D=[Dx,0;0,Dy];
DispersionRelation=arrayfun(@(q)max(real(eig(J-q^2*D))),qspace);
swd=(max(real(eig(J)))<0); %StableWithoutDiffusion
na=false; %NoiseAmplifying
pfsk=(max(DispersionRelation)>0);%PositiveForSomeK
if swd&&(~na)&&pfsk
out=0;%patterns
elseif swd&&(~na)&&(~pfsk)
out=1;%always stable
elseif swd&&na
out=2; %noise amplifying
elseif ~swd
out=3; %unstable
end
end