필터 지우기
필터 지우기

Interp1 returning NaN for ode45

조회 수: 4 (최근 30일)
Rafael Félix Soriano
Rafael Félix Soriano 2018년 11월 10일
답변: Bruno Luong 2018년 11월 10일
The interp1 function returns a NaN value for the time-dependent functions that are later introduced in the ode45 function. Using pchip or spline instead of interp1 returns a value other than NaN, but the solutions are not valid.
My code is the one that follows, thank you in advance.
function [T_V , W_V , x_V , y_V , h_V , Vfin] = viraje(x0 , y0 , h0 , W0 ,...
V0 , Vf , gamma0 , gammaf , chi0 , chif , chip , CD0 , k , S , cE)
clc
g = 9.80665; %m/s2
tfin = (chif - chi0)/chip;
gammap = (gammaf-gamma0)/tfin;
Vp = (Vf-V0)/tfin;
t_interp = linspace(0 , tfin , 10);
V = V0 + Vp*t_interp;
gamma = gamma0 + gammap*t_interp;
mu = atan(chip/g*V.*cos(gamma)./(V/g*gammap + cos(gamma)));
h = Hviraje(t_interp , h0 , V0 , Vp , gamma0 , gammap);
rho = densidad(h);
function dWdt = dif(t,W,t_interp,V,gamma,mu,rho)
if Vp == 0
V = V0;
else
V = interp1(V,t_interp,t,'pchip');
end
if gammap == 0
gamma = gamma0;
else
gamma = interp1(gamma,t_interp,t,'pchip');
end
if Vp == 0 && gammap == 0
mu = mu(1);
else
mu = interp1(mu,t_interp,t,'pchip');
end
if gamma0 == 0 && gammap == 0
rho = densidad(h0);
else
rho = interp1(rho,t_interp,t,'pchip');
end
dWdt = -cE*(0.5*rho.*(V.^2)*S*CD0 + 2*k./(rho*S)*(chip*W.*cos(gamma)./...
(g*sin(mu))).^2 + W.*sin(gamma) + Vp/g*W)
end
tspan = [0 tfin];
ic = W0;
[t,W] = ode45(@(t,W) dif(t,W,t_interp,V,gamma,mu,rho) , tspan , ic);
W_V = W;
h_V = Hviraje(t , h0 , V0 , Vp , gamma0 , gammap);
x_V = Xviraje(t , x0 , V0 , Vp , gamma0 , gammap , chi0 , chip);
y_V = Yviraje(t , y0 , V0 , Vp , gamma0 , gammap , chi0 , chip);
rho = densidad(h_V)';
gamma = gamma0 + gammap*t;
V = V0 + Vp*t;
mu = atan(chip/g*V.*cos(gamma)./(V/g*gammap + cos(gamma)));
T_V = 0.5*rho.*V.^2*S*CD0 + 2*k./(rho*S).*(chip*W.*cos(gamma)./...
(g*sin(mu))).^2 + W.*sin(gamma) + Vp/g*W;
Vfin = Vf;
end
function H = Hviraje(T , h0 , V0 , Vp , gamma0 , gammap)
syms h(t)
eqn = diff(h,t) == (V0 + Vp*t)*sin(gamma0 + gammap*t);
cond = h(0) == h0;
ySol(t) = dsolve(eqn,cond);
H = ySol(T);
end
function X = Xviraje(T , x0 , V0 , Vp , gamma0 , gammap , chi0 , chip)
syms x(t)
eqn = diff(x,t) == (V0 + Vp*t)*cos(gamma0 + gammap*t)*cos(chi0 + chip*t);
cond = x(0) == x0;
ySol(t) = dsolve(eqn,cond);
X = ySol(T);
end
function Y = Yviraje(T , y0 , V0 , Vp , gamma0 , gammap , chi0 , chip)
syms y(t)
eqn = diff(y,t) == (V0 + Vp*t)*cos(gamma0 + gammap*t)*sin(chi0 + chip*t);
cond = y(0) == y0;
ySol(t) = dsolve(eqn,cond);
Y = ySol(T);
end
  댓글 수: 1
Walter Roberson
Walter Roberson 2018년 11월 10일
interp1() returns nan for most interpolation methods when you ask it to project at locations that are outside the range of the first parameter, unless you pass it the 'extrap' parameter.

댓글을 달려면 로그인하십시오.

채택된 답변

Bruno Luong
Bruno Luong 2018년 11월 10일
V = interp1(V,t_interp,t,'pchip');
It looks arguments are swapped; I think it should be
V = interp1(t_interp,V,t,'pchip');
and similar for other INTERP1 calls.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Dimensionality Reduction and Feature Extraction에 대해 자세히 알아보기

제품


릴리스

R2018a

Community Treasure Hunt

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

Start Hunting!

Translated by