How do I output the second derivative from an ODE solver for further use?
4 views (last 30 days)
Show older comments
I am currently working on a coupled ode problem, and using ode45 solver to solve this.
My code is something like this:
bubble.m
function rdot = f(t, r)
%Some mathematical expressions
....
...
...
rdot(1) = r(2); % first derivative of r1
rdot(2) = (X11 - X21 - (X31*(X41 - X51)))/X81; % second derivative of r1
rdot(3) = r(4); % first derivative of r2
rdot(4) = (X12 - X22 - (X32*(X42 - X52)))/X82; % second derivative of r2
rdot = rdot';
And,
bubble_plotter.m
clc;
clear all;
%close all;
time_range = [0 1d-4];
r1_eq = 10d-6;
r2_eq = 1d-6;
f_s = 20000;
T_s = 1/f_s;
initial_conditions = [r1_eq 0.d0 r2_eq 0.d0];
[t, r] = ode45('bubble', time_range, initial_conditions);
r1_us=1000000*r(:,1);
r1dot=r(:,2);
r2_us=1000000*r(:,3);
r2dot=r(:,4);
normalized_time = t./T_s;
figure(101);
h1 = plot(normalized_time, r1_us, 'b-', normalized_time, r2_us, 'k-.');
set(h1,'linewidth',3);
txt = text(0.03, 21.5, (['d = ' num2str(d_close), '\mum']));
txt.FontSize = 25;
legend('Bubble1', 'Bubble2')
legend('Location','Northwest')
set(gca,'xscale','linear','FontSize',26)
set(gca,'yscale','linear','FontSize',26)
set(gca,'XMinorTick','on')
set(gca,'YMinorTick','on')
You can see, I am pulling out r(:,1), r(:,2), r(:,3) and r(:,4) as r1_us, r1dot and r2_us, r2dot for plotting and further uses. But I also need the values of rdot(2) and rdot(4) from the main function file. How can I pull those double derivative values of r1 and r2 into my plotting file for further use?
0 Comments
Accepted Answer
Torsten
on 14 Nov 2018
After the line
[t, r] = ode45('bubble', time_range, initial_conditions);
insert
for i=1:numel(t)
t_actual = t(i);
r_actual = r(i,:);
rdot(i,:) = bubble(t_actual,r_actual);
end
Best wishes
Torsten.
7 Comments
More Answers (0)
See Also
Categories
Find more on Interactive Control and Callbacks in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!