ode45 Error using odearguments (line 95)

14 visualizaciones (últimos 30 días)
Rick Sellers
Rick Sellers el 17 de Jul. de 2020
Comentada: Star Strider el 18 de Jul. de 2020
I know this one has been asked before, and I looked through several of those answers. I've been over it and over it and can't see where the solution vector is of lentgh 1. Thanks for the help. The first block of code is my main script. THe second is my ode function called by ode45. Here are the errors I get:
@(T,PRY)PRY_ODE(T,PRY,OMEGA,R,ME,G,IXX,IYY,IZZ) returns a vector of length 1, but the length of initial conditions vector is 6. The vector returned by
@(T,PRY)PRY_ODE(T,PRY,OMEGA,R,ME,G,IXX,IYY,IZZ) and the initial conditions vector must have the same number of elements.
Error in ode45 (line 115)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Error in PRY_main (line 35)
[t,prydot] = ode45(@(t,pry) pry_ODE(t,PRY,Omega,R,mE,G,Ixx,Iyy,Izz),tspan,PRY);
%% Main script
clc; clear;
% Set initial values
PRY = zeros(1,6);
phi = 0; tht = 0; psi = 0;
wx = 0; wy = 0; wz = 0;
RE = 6357; % km - radius of earth
RS = 500; % km - orbital altitude
R = RE+RS; % km
mE = 5.9722e24; % kg
G = 6.6743e-11/1000^3; % m^3/kg s^2
Omega = sqrt(G*mE/R^3); % rad/s
Ixx = 6; Iyy = 8; Izz = 4; % kg m^2
tend = (2*pi/Omega)/60;
tspan = [0 tend];
init = [1;5;10;15;20;25];
% for k= 1:6
% phi =ic(k);
% tht=ic(k);
% psi=ic(k);
phi = init(1);
tht = init(1);
psi = init(1);
ic = [phi;tht;psi;wx;wy;wz];
PRY = [phi,tht,psi,wx,wy,wz];
[t,prydot] = ode45(@(t,pry) pry_ODE(t,PRY,Omega,R,mE,G,Ixx,Iyy,Izz),tspan,PRY);
% end
%% pry_ODE function file
% pry_ODE.m
% myODE function for solver
% PRY = [phi,tht,psi,wx,wy,wz,Omega,R,mE,G,Ixx,Iyy,Izz];
function [t,prydot] = pry_ODE(t,pry,Omega,R,mE,G,Ixx,Iyy,Izz)
prydot = zeros(size(pry));
phi = pry(1);
tht = pry(2);
psi = pry(3);
wx = pry(4);
wy = pry(5);
wz = pry(6);
K = 3*G*mE/R^5;
X = -R*cos(tht)*cos(psi);
Y = -R*(cos(psi)*sin(phi)*sin(tht)-cos(phi)*sin(psi));
Z = -R*(cos(phi)*cos(psi)*sin(tht)+sin(phi)*sin(psi));
prydot(1) = (Omega+prydot(3))*sin(tht)+wx;
prydot(2) = (1/cos(phi))*(-(Omega+prydot(3))*cos(tht)*sin(phi)+wy);
prydot(3) = (1/(cos(phi)*cos(tht)))*(prydot(2)*sin(phi)-Omega*cos(phi)*cos(tht)+wz);
prydot(4) = (1/Ixx)*((K*Y*Z-wy*wz)*(Izz-Iyy));
prydot(5) = (1/Iyy)*((K*X*Z-wx*wz)*(Ixx-Izz));
prydot(6) = (1/Izz)*((K*X*Y-wx*wy)*(Iyy-Ixx));
end
Thank you,
Rick

Respuesta aceptada

Star Strider
Star Strider el 18 de Jul. de 2020
First, specify:
prydot = zeros(6,1);
to create the necessary column vector.
Then, since it is only neccessary for ‘pry_ODE’ to return the vector of derivatives, not the time as well, change the function declaration to:
prydot = pry_ODE(t,pry,Omega,R,mE,G,Ixx,Iyy,Izz)
and your system runs without error.
  2 comentarios
Rick Sellers
Rick Sellers el 18 de Jul. de 2020
Thank you, Star Strider.
Star Strider
Star Strider el 18 de Jul. de 2020
As always, my pleasure!

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Programming en Help Center y File Exchange.

Productos


Versión

R2019b

Community Treasure Hunt

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

Start Hunting!

Translated by