使用“绘图”命令时矢量长度不一致

Question

我在下面创建了一个函数。该函数是使用 Euler-Maruyama 方法的随机微分方程的通用求解器。

function [Y,t]=eulermaruyama_solver(Y0,T,a,b,c,F,G,deltaW)
% deltaW is called on in the test script

dt=T/b; t=[0:dt:T]; sd=sqrt(dt);
dw = zeros(a,b+1);
dw(:,2:b+1) = sd * deltaW;  
Y=zeros(c,b+1);
Y(:,1) = Y0; Yn = Y0;
for n = 1:b
    Y(:,n+1) = Yn + F(t(n),Yn)*dt + G(t(n),Yn) * dw(:,n+1);
    Yn = Y(:,n+1);
end

为了在 a=c=1 的简单情况下测试我的函数，我定义了函数 F 和 G。下面是我的测试脚本。

a=1; c=1; Y0=1; T=1; b=5;
F=@(t,y) y^(1/2); G=@(t,y) 2*y;



% Euler-Maruyama solution
[deltaW,W]=randn_numbers(a,b);
[Y,t]=eulermaruyama_solver(Y0,T,a,b,c,F,G,deltaW);

% Exact solution
ExactY=zeros(1,b+1);
ExactY=(2^(1/5)+(1/2)*W).^3;

% Plot 
subplot(1,2,1)
plot(t, Y, 'r', t, ExactY,'b')

但显示以下错误消息：“使用绘图时出错，向量的长度必须相同。

我弄错了吗？

Answer 1

您的 ExactY 是一个 1x5 数组，而 t 是一个 1x6 数组（Y 也是）。对于 plot，每个 x/y 对的输入数组必须具有相同的长度，如错误消息所示。

这一行创建了一个 1x6 ExactY 但它是多余的，因为 ExactY 在下一行被覆盖

ExactY = zeros(1, b+1);

我不清楚正确的索引是什么，但像这样的事情可以解决错误，也许只有你才能决定正确的方法是什么，但同样的原则将适用：

ExactY=zeros(1,b+1);
ExactY(2:end)=(2^(1/3)+(1/3)*W).^3;

Answer 2

从数学上讲，问题在于您的确切解在两个方面是错误的。

最严重的是您没有使用 sqrt(dt) 应用缩放，给出了错误的量级。
不太严重，但仍然给出一个重大错误是您的初始值不匹配。您的精确解 Y(t)=(2^(1/3)+W(t)/3)^3 有 Y(0)=2，但是您使用 Y0=1。也许你应该设置ExactY(2:b+1,:)=(Y0^(1/3)+(1/3)*W*sqrt(dt)).^3;，别忘了设置ExactY(1,:)=Y0;。或者用一个元素构造 W 使得 W(1,:)==0.

顺便说一句，如果 a 和 c 的值不同，您认为 eulermaruyama_solver 中应该发生什么？

保持结构，您可以将方法更改为

function [dW,W]=randn_numbers(a,b,T)
  W=zeros(a,b+1);
  dW=randn(a,b)*sqrt(T/b);
  W(:,2:end)=cumsum(dW,2);
end%function

function [Y,t]=eulermaruyama_solver(Y0,T,a,b,c,F,G,dW)
  dt=T/b; t=[0:dt:T]; 
  Y=zeros(a,b+1);
  Y(:,1) = Y0; Yn = Y0;
  for n = 1:b
    Y(:,n+1) = Yn + F(t(n),Yn)*dt + G(t(n),Yn) .* dW(:,n);
    Yn = Y(:,n+1);
  end%for
end%function

a=4; c=a; Y0=1; T=1; b=50;
F=@(t,y) (1/3)*y.^(1/3); G=@(t,y) y.^(2/3);

[dW,W]=randn_numbers(a,b,T);
[Y,t]=eulermaruyama_solver(Y0,T,a,b,c,F,G,dW);

% Exact solution
ExactY=(Y0^(1/3)+(1/3)*W).^3;

% Plot 
plot(t, Y, 'r', t, ExactY,'b')

使用“绘图”命令时矢量长度不一致

Inconsistent vector length when using the `plot` command

function

matlab

plot

numerical-methods