Error using ode45 : Error using ==> horzcat CAT arguments dimensions are not consistent.

1 visualización (últimos 30 días)
Hello,
I am trying to use ode45 to integrate the nondimensional circular restricted three body problem equations of motion. In order to do this I have created the following function
function [ Yi,Y,mu,t_event,y_event ] = CR3BP( cb,pb,Y0,TOI )
%[Yi,Y,mu,t_event,y_event]=CR3BP(cb,pb,Y0,TOI) This function integrates the CR3BP EOMs and the state transition matrix
INPUTS: cb - a string containing the central body
pb - a string containing the perturbing body
Y0 - a (42,1) column vector contain the initial conditions
TOI- the end time of integration (note, the initial time is set
to 0)
OUTPUT: Yi - The nondimensionalized initial conditions
Y - The solution to the CR3BP and the State Transition Matrix
mu - The nondimensional gravitational parameter
t_event - The time that the integration is stopped (this occurs
when x=0)
y_event - The solution when the integration is stopped (this occurs
when x_sc=0)
Created by Andrew Liounis 3/17/13
%%mu
[~,~,mucb]=cbinfo(cb); %calling value for mu for central body
[~,rorbit,mupb]=cbinfo(pb); %calling value for mu for perturbing body
mustar=mucb+mupb; %calculating characteristic gravitaional parameter
mu=mupb/mustar; %caluculating non-dimensional gravitational parameter
%%lengths
Lstar=rorbit; %calling characteristic length
%%time
tstar=sqrt(Lstar.^3/mustar); %calculating characteristic time
%%Nondimensionalizing Y0
Yi(1:3)=Y0(1:3)./Lstar-mu.*[1,0,0];
Yi(4:6)=Y0(4:6)*tstar/Lstar-cross([0,0,1],Yi(1:3));
Yi(7:42)=Y0(7:42);
%%Nondimensionalizing Integration time
TOInd=TOI/tstar;
%%ODE defenition
d=@(t,IC) sqrt((IC(1)+mu)^2+IC(2)^2+IC(3)^2);
r=@(t,IC) sqrt((IC(1)-1+mu)^2+IC(2)^2+IC(3)^2);
uxx=@(t,IC) 1-(1-mu)/d(t,IC)^3-mu/r(t,IC)^3+3*(1-mu)*(IC(1)+mu)^2/d(t,IC)^5+3*mu*(IC(1)-1+mu)^2/r(t,IC)^5;
uyy=@(t,IC) 1-(1-mu)/d(t,IC)^3-mu/r(t,IC)^3+3*(1-mu)*IC(2)^2/d(t,IC)^5+3*mu*IC(2)^2/r(t,IC)^5;
uzz=@(t,IC) -(1-mu)/d(t,IC)^3-mu/r(t,IC)^3+3*(1-mu)*IC(3)^2/d(t,IC)^5+3*mu*IC(3^2)/r(t,IC)^5;
uxy=@(t,IC) 3*(1-mu)*(IC(1)+mu)*IC(2)/d(t,IC)^5+3*mu*(IC(1)-1+mu)*IC(2)/r(t,IC)^5;
uxz=@(t,IC) 3*(1-mu)*(IC(1)+mu)*IC(3)/d(t,IC)^5+3*mu*(IC(1)-1+mu)*IC(3)/r(t,IC)^5;
uyz=@(t,IC) 3*(1-mu)*IC(2)*IC(3)/d(t,IC)^5+3*mu*IC(2)*IC(3)/r(t,IC)^5;
A=@(t,IC) [zeros(3,3),eye(3);uxx(t,IC),uxy(t,IC),uxz(t,IC),0,2,0;uxy(t,IC),uyy(t,IC),uyz(t,IC),-2,0,0;uxz(t,IC),uyz(t,IC),uzz(t,IC),0,0,0];
ode=@(t,IC) [IC(4);IC(5);IC(6);2*IC(5)+IC(1)-(1-mu)*(IC(1)+mu)/sqrt((IC(1)+mu)^2+IC(2)^2+IC(3)^2)^3-mu*(IC(1)-1+mu)/sqrt((IC(1)-1+mu)^2+IC(2)^2+IC(3)^2)^3;-2*IC(4)+IC(2)-(1-mu)*IC(2)/sqrt((IC(1)+mu)^2+IC(2)^2+IC(3)^2)^3-mu*IC(2)/sqrt((IC(1)-1+mu)^2+IC(2)^2+IC(3)^2)^3;-(1-mu)*IC(3)/sqrt((IC(1)+mu)^2+IC(2)^2+IC(3)^2)^3-mu*IC(3)/sqrt((IC(1)-1+mu)^2+IC(2)^2+IC(3)^2)^3,...
A(t,IC)*(IC(7:12));A(t,IC)*(IC(13:18));A(t,IC)*(IC(19:24));A(t,IC)*(IC(25:30));A(t,IC)*(IC(31:36));A(t,IC)*(IC(37:42))];
%%ODE Calculation
options = odeset('RelTol',1e-12,'AbsTol',1e-15,'Events',@Event_function);
[~,Y,t_event,y_event]=ode45(ode,[0,TOInd],Yi,options);
end
The events function I am using is
function [Trigger,Terminate,Direction]=Event_function(t,y)
mu=.0121506;
Trigger=min(y(1)+mu,(1-mu)-y(1));
Terminate=1;
Direction=0;
end
I am using anonymous functions for the ode definition to avoid using global variables.
My problem is, when I make the following call:
clc, close all
Y0=[.83692-.029324,-0.029324,0,.085982,.039561,0,...
1,0,0,0,0,0,...
0,1,0,0,0,0,...
0,0,1,0,0,0,...
0,0,0,1,0,0,...
0,0,0,0,1,0,...
0,0,0,0,0,1];
global Mu
Mu=mu;
options = odeset('RelTol',1e-12,'AbsTol',1e-15,'Events',@Event_function);
[Y,Y1,mu1,t_event,y_event]=CR3BP('Earth','Moon',Y0,-10);
I get an error of
??? Error using ==> horzcat
CAT arguments dimensions are not consistent.
Error in ==>
CR3BP>@(t,IC)[IC(4);IC(5);IC(6);2*IC(5)+IC(1)-(1-mu)*(IC(1)+mu)/sqrt((IC(1)+mu)^2+IC(2)^2+IC(3)^2)^3-mu*(IC(1)-1+mu)/sqrt((IC(1)-1+mu)^2+IC(2)^2+IC(3)^2)^3;-2*IC(4)+IC(2)-(1-mu)*IC(2)/sqrt((IC(1)+mu)^2+IC(2)^2+IC(3)^2)^3-mu*IC(2)/sqrt((IC(1)-1+mu)^2+IC(2)^2+IC(3)^2)^3;-(1-mu)*IC(3)/sqrt((IC(1)+mu)^2+IC(2)^2+IC(3)^2)^3-mu*IC(3)/sqrt((IC(1)-1+mu)^2+IC(2)^2+IC(3)^2)^3,A(t,IC)*(IC(7:12));A(t,IC)*(IC(13:18));A(t,IC)*(IC(19:24));A(t,IC)*(IC(25:30));A(t,IC)*(IC(31:36));A(t,IC)*(IC(37:42))]
at 52
ode=@(t,IC)
[IC(4);IC(5);IC(6);2*IC(5)+IC(1)-(1-mu)*(IC(1)+mu)/sqrt((IC(1)+mu)^2+IC(2)^2+IC(3)^2)^3-mu*(IC(1)-1+mu)/sqrt((IC(1)-1+mu)^2+IC(2)^2+IC(3)^2)^3;-2*IC(4)+IC(2)-(1-mu)*IC(2)/sqrt((IC(1)+mu)^
Error in ==> odearguments at 98
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ==> ode45 at 172
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in ==> CR3BP at 59
[~,Y,t_event,y_event]=ode45(ode,[0,TOInd],Yi,options);
I have tried to figure out what the problem is but am at a loss. Can someone help me?
Thanks!
Andrew Liounis

Respuesta aceptada

Walter Roberson
Walter Roberson el 17 de Mzo. de 2013
This line
/sqrt((IC(1)-1+mu)^2+IC(2)^2+IC(3)^2)^3,A(t,IC)*(IC(7:12));A(t,IC)*(IC(13:18));
has a comma instead of a semi-colon

Más respuestas (0)

Categorías

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

Productos

Community Treasure Hunt

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

Start Hunting!

Translated by