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

1 回表示 (過去 30 日間)
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

採用された回答

Walter Roberson
Walter Roberson 2013 年 3 月 17 日
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
  1 件のコメント
Andrew
Andrew 2013 年 3 月 18 日
編集済み: Andrew 2013 年 3 月 18 日
Thanks for your help! I can't believe I missed that.

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeProgramming についてさらに検索

製品

Community Treasure Hunt

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

Start Hunting!

Translated by