Clear Filters
Clear Filters

Save Variable from ODE45

8 views (last 30 days)
Holden Tranquillo
Holden Tranquillo on 2 May 2022
Edited: Torsten on 3 May 2022
Hello All,
I am trying to use an outputFcn to save a variable that I calculate in ODE45. So in the end, I would like the value of u1 and u2 for all time steps taken.
I cannot figure out why this isn't saving it to workspace, can anyone see what I am missing?
Thanks

Answers (2)

Torsten
Torsten on 2 May 2022
can anyone see what I am missing?
You do not pass "options" to ode45.
  1 Comment
Holden Tranquillo
Holden Tranquillo on 2 May 2022
Thanks, not sure how I missed that. Anyway, with that I still get some errors. I have never used an output function so I'm pretty unfamiliar with debugging it. But in each ODE45 time step, an input u1 and u2 are calculated. I then want to store the timeseries of these

Sign in to comment.


Torsten
Torsten on 3 May 2022
Edited: Torsten on 3 May 2022
Try this:
figNum = 1;
%% System Parameters
% Aircraft
M = 5e4; % mass [kg]
J = 1.25e4; % moment of inertia [kg.m^2]
g = 9.8; % gravitational acceleration [m/s^2]
l = 5; % 1/2 the wingspan [m]
params = [M J g l]';
alpha = pi()/6; % tilt angle of applied couple [rad]
%% Set Point Controlled Response (New Outputs)
% Simulation Time
tmax = 30; % simulation time [s]
tsim = linspace(0,tmax,100)'; % time vector [s]
colors = ['r' 'g' 'b' 'm' 'c'];
% Loop through varying initial conditions
for i = 1:1
% Initial Conditons
x0 = [i 0 5*i 10 0.5 0]';
% Desired Set Point
setPoint = [0 0]';
% ODE45
odeFcn = @(t,x) aircraftSetPointB(t,x,alpha,params,setPoint);
[tSet,xSet] = ode45(odeFcn,tsim,x0);
U = postprocess(tSet,xSet,alpha,params,setPoint);
figure(1)
plot(tSet,U(:,1))
figure(2)
plot(tSet,U(:,2))
end
function U = postprocess(tSet,xSet,alpha,params,setPoint)
U = zeros(numel(tSet),2);
% Read parameters
M = params(1);
J = params(2);
g = params(3);
l = params(4);
x = zeros(6,1);
for i = 1:numel(tSet)
x = xSet(i,:);
aff = [-g 0]';
E = [1/M*cos(x(5)) 2/M*sin(x(5))*sin(alpha) ; 0 2*l/J*cos(alpha)];
V = [-(x(4)-0) - 0.9*(x(3)-setPoint(1)) ; -(x(6)-0) - (x(5)-setPoint(2))];
U(i,:) = -E^-1*aff + E^-1*V;
end
end
function xDot = aircraftSetPointB(~,x,alpha,params,setPoint)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FUNCTION NAME: aircraftSetPoint
%
% PURPOSE: ODE45 function for control of aircraft dynamics to fixed set points (outputs y and theta)
%
% INPUTS:
% t:
% x: current state vector (nx1)
% alpha: value of alpha
% params: System parameters as vector [M J g l]^T
% setPoint: vector of desired set points [x1_des ; x3_des]
%
% OUTPUTS:
% xDot: state vector derivatives
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% AUTHOR: Holden Tranquillo
% DATE: 5-2-22
%
% DESCRIPTION OF LOCAL VARIABLES
%
% FUNCTIONS CALLED
%
% START OF EXECUTABLE CODE
%
% Read parameters
M = params(1);
J = params(2);
g = params(3);
l = params(4);
% Calculate Control Input
aff = [-g 0]';
E = [1/M*cos(x(5)) 2/M*sin(x(5))*sin(alpha) ; 0 2*l/J*cos(alpha)];
V = [-(x(4)-0) - 0.9*(x(3)-setPoint(1)) ; -(x(6)-0) - (x(5)-setPoint(2))];
U = -E^-1*aff + E^-1*V;
u1 = U(1);
u2 = U(2);
% Calculate Derivatives
xDot = zeros(6,1); % initialize
xDot(1) = x(2);
xDot(2) = ( -1/M*sin(x(5)) )*u1 + ( 2/M*cos(x(5))*sin(alpha) )*u2;
xDot(3) = x(4);
xDot(4) = -g + ( 1/M*cos(x(5)) )*u1 + ( 2/M*sin(x(5))*sin(alpha) )*u2;
xDot(5) = x(6);
xDot(6) = ( 2*l/J*cos(alpha) )*u2;
end

Community Treasure Hunt

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

Start Hunting!