你好同学,如此即可:
关于x、tao的
clc;close;clear;
x=0:1:200;
D=0.25;
Q=0.50;
lambad=0.1;% 0<lambad<1;
alpha = 0.01;
beta = -0.1;
tao_arr = -1:0.01:-0.01;
Pst_Mat = zeros(length(tao_arr), numel(x));
for i = 1:1:numel(tao_arr)
tao = tao_arr(i);
f= @(x)((beta+tao*beta^2)*x.^2/4*D)-(beta+tao*beta^2)*(lambad*sqrt(D*Q)+2*D)*x/D^2;
C=(2*D*(alpha+tao*alpha*beta)-(beta+tao*beta^2)*(Q-4*Q*(lambad^2)-4*lambad*sqrt(D*Q)-3))/4*D^2;
E=-(2*D*lambad*(alpha+tao*alpha*beta)+(beta+tao*beta^2)*(lambad*Q*sqrt(D*Q)*(4*(lambad^2)-3)+D*Q*(8*(lambad^2)-4)+3*lambad*D*sqrt(D*Q)))/2*D^2;
a=@(x)D*x.^2+2*lambad*x*sqrt(D*Q)+Q;
b=E/(((1-lambad^2)*D*Q)^0.5);
m=@(x)(D*x+lambad*sqrt(D*Q))/(((1-lambad^2)*D*Q)^0.5);
P=@(x)exp((f(x)+b.*atan(m(x)))).*a(x).^C;
M = integral(P,0,inf);
N=1/M;
Pst=N*P(x);
Pst_Mat(i,:) = Pst;
end
mesh( x,tao_arr, Pst_Mat)
xlabel('x');ylabel('tao');zlabel('Pst')
关于x、alpha的
clc;close;clear;
x=0:1:200;
D=0.25;
Q=0.50;
lambad=0.1;% 0<lambad<1;
alpha_arr = 0:0.01:10;
tao = 0.1;
beta = -0.1;
Pst_Mat = zeros(length(alpha_arr), numel(x));
for i = 1:1:numel(alpha_arr)
alpha = alpha_arr(i);
f= @(x)((beta+tao*beta^2)*x.^2/4*D)-(beta+tao*beta^2)*(lambad*sqrt(D*Q)+2*D)*x/D^2;
C=(2*D*(alpha+tao*alpha*beta)-(beta+tao*beta^2)*(Q-4*Q*(lambad^2)-4*lambad*sqrt(D*Q)-3))/4*D^2;
E=-(2*D*lambad*(alpha+tao*alpha*beta)+(beta+tao*beta^2)*(lambad*Q*sqrt(D*Q)*(4*(lambad^2)-3)+D*Q*(8*(lambad^2)-4)+3*lambad*D*sqrt(D*Q)))/2*D^2;
a=@(x)D*x.^2+2*lambad*x*sqrt(D*Q)+Q;
b=E/(((1-lambad^2)*D*Q)^0.5);
m=@(x)(D*x+lambad*sqrt(D*Q))/(((1-lambad^2)*D*Q)^0.5);
P=@(x)exp((f(x)+b.*atan(m(x)))).*a(x).^C;
M = integral(P,0,inf);
N=1/M;
Pst=N*P(x);
Pst_Mat(i,:) = Pst;
end
mesh( x,alpha_arr, Pst_Mat)
xlabel('x');ylabel('alpha');zlabel('Pst')
关于x、beta的
clc;close;clear;
x=0:1:200;
D=0.25;
Q=0.50;
lambad=0.1;% 0<lambad<1;
alpha = 0.01;
tao = 0.1;
beta_arr = -1:0.01:-0.01;
Pst_Mat = zeros(length(beta_arr), numel(x));
for i = 1:1:numel(beta_arr)
beta = beta_arr(i);
f= @(x)((beta+tao*beta^2)*x.^2/4*D)-(beta+tao*beta^2)*(lambad*sqrt(D*Q)+2*D)*x/D^2;
C=(2*D*(alpha+tao*alpha*beta)-(beta+tao*beta^2)*(Q-4*Q*(lambad^2)-4*lambad*sqrt(D*Q)-3))/4*D^2;
E=-(2*D*lambad*(alpha+tao*alpha*beta)+(beta+tao*beta^2)*(lambad*Q*sqrt(D*Q)*(4*(lambad^2)-3)+D*Q*(8*(lambad^2)-4)+3*lambad*D*sqrt(D*Q)))/2*D^2;
a=@(x)D*x.^2+2*lambad*x*sqrt(D*Q)+Q;
b=E/(((1-lambad^2)*D*Q)^0.5);
m=@(x)(D*x+lambad*sqrt(D*Q))/(((1-lambad^2)*D*Q)^0.5);
P=@(x)exp((f(x)+b.*atan(m(x)))).*a(x).^C;
M = integral(P,0,inf);
N=1/M;
Pst=N*P(x);
Pst_Mat(i,:) = Pst;
end
mesh( x,beta_arr, Pst_Mat)
xlabel('x');ylabel('beta');zlabel('Pst')