Filter löschen
Filter löschen

Polyfit error warning for SDE

6 Ansichten (letzte 30 Tage)
Mel Hernandez
Mel Hernandez am 5 Apr. 2019
Kommentiert: Mel Hernandez am 5 Apr. 2019
I'm getting this Warning everytime I use polyfit, and I can't seem to find the error behind my code. Here's the warning:
Warning: Polynomial is not unique; degree >= number of data points.
I'm trying to find the slope of 4 points from using the Euler-Maruyama method applied to my SDE. Here's the code:
% Euler-Maruyama method on linear SDE
% SDE is dX = -cos^3(X)*sin(X)*dt + cos^2(X)*dW
% X(0) = Xzero, Xzero = 0
% For 0 <= t <= 1
% Time steps: 1/50, 1/100, 1/150, 1/200
% Exact Solution: X(t) = tan^-1(W(t))
% K = 50, 100, 150, 200; % number of steps
% T = 1; % end of interval
% tau = T/K; % time step
% A for loop for the 4 timesteps (1/50, 1/100, 1/150, 1/200)
for f = 1:4
K = 50*f; % So we can get 50, 100, 150, 200
T = 1;
tau = T/K;
N_s = 2000; % number of trajectories
W = zeros(N_s,K+1); % setting up discrete Brownian path
dW = zeros(N_s,K); % setting up Brownian increments
X_app = zeros(N_s,K+1); % setting up Approx solution
X_exact = zeros(N_s,K+1); % setting up Exact solution
Xerr = zeros(1,K+1); % setting up Global error
% Compute 2000 independent Brownian Motion trajectories
for q = 1:N_s
for p = 1:K
W(q,p+1) = W(q,p) + sqrt(tau)*normrnd(0,1);
end
end
% Use these Brownian Motion trajectories to define the increments
for c = 1:N_s
for n = 1:K
dW(c,n) = W(c,n+1) - W(c,n);
end
end
% Generate 2000 solution trajectories using the EM update & these increments
for j = 1:N_s
for k = 1:K
Xzero = 0;
X_app(j,1) = Xzero;
X_app(j,k+1) = X_app(j,k) + (-cos(X_app(j,k))^3)*(sin(X_app(j,k)))*tau + (cos(X_app(j,k))^2)*dW(j,k);
end
end
% Estimate the global error of tau for the timestep using the exact solution
for tk = 1:K+1
for g = 1:N_s
X_exact(g,tk) = atan(W(g,tk));
Xerr(g,tk) = abs(X_exact(g,tk) - X_app(g,tk));
end
Err(tk) = (1/N_s)*sum(Xerr(:,tk));
end
GErr = max(Err)
% Make a log-log plot of the data through a scatter plot
scatter(log(tau),log(GErr),'filled'), hold on
end
hold off
% Finding the slop of the line is the strong order of convergence
di = polyfit(log(tau),log(GErr),1);
slope = di(1)
I don't understand why I'm getting the error. It should be a 1 term polynomial (y = mx+b), and I only need the slope for my strong order of convergence. I get an answer with the warning, but I don't want to feel contempt with my answer when the warning is happening.

Akzeptierte Antwort

Jan
Jan am 5 Apr. 2019
Bearbeitet: Jan am 5 Apr. 2019
tau and GErr are scalars, so you have one single point for the polynomial only. Then you cannot determine a slope.
I guess you want to create vectors instead. Then replace all tau and GErr to tau(f) and GErr(f):
...
tau(f) = T/K;
...
W(q,p+1) = W(q,p) + sqrt(tau(f))*normrnd(0,1);
...
GErr(f) = max(Err);
...
Then this will be working:
di = polyfit(log(tau),log(GErr),1);
By the way, this is not useful:
for k = 1:K
Xzero = 0; % ??? Set this once outside the loops
X_app(j,1) = Xzero; % ??? Why setting this in each iteration again
X_app(j,k+1) = X_app(j,k) + (-cos(X_app(j,k))^3)*(sin(X_app(j,k)))*tau + (cos(X_app(j,k))^2)*dW(j,k);
end
This:
for c = 1:N_s
for n = 1:K
dW(c,n) = W(c,n+1) - W(c,n);
end
end
can be simpified to:
dW = diff(W, 1, 2)
and
for tk = 1:K+1
for g = 1:N_s
X_exact(g,tk) = atan(W(g,tk));
Xerr(g,tk) = abs(X_exact(g,tk) - X_app(g,tk));
end
Err(tk) = (1/N_s)*sum(Xerr(:,tk));
end
to
X_exact = atan(W);
Xerr = abs(X_exact - X_app);
Err = sum(Xerr, 1) / N_s;
  1 Kommentar
Mel Hernandez
Mel Hernandez am 5 Apr. 2019
Thank you so much! I completely forgot about tau(f) and GErr(f). It fixed it! And you didn't have to clean up my code. I do appreciate it though. Realized I did have a few too many lines going on. Thank you again!

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (0)

Kategorien

Mehr zu Time Series Objects finden Sie in Help Center und File Exchange

Produkte

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by