Solving ODE45 with input vectors

1 visualizzazione (ultimi 30 giorni)
Jestoni Orcejola
Jestoni Orcejola il 31 Lug 2020
Commentato: Alan Stevens il 10 Set 2020
I would like to elvaluate the speed of a float as it travels through the water column with depth/time dependent seawater density. The vector rho is the density profile of the sea water from 1000m to 0m. I would like o pass it through the function to evaluate speed.
clear;
close ;
%% m files for obtaining sea water denstiy profile
Sea_Water_Density_s_t_p;
Float_Velocity_Calculation;
tspan=[0 1 1000];
y0 =[1000,0,1030];
rho=descending_rho_profile;
[t,y]=ode45(@(t,y) f(t,y,rho),tspan,y0);
hold off
plot(t,y(:,2))
grid on
Here is my function:
function speed=f(t,v,rho)
Cd=.16;
A=.0613;
m=82.25;
V=0.0809815;
%rho=1030;
g=9.81;
speed=[v(2); (rho*V*g/(m+V*rho))-(m*g/(m+V*rho))-((1/(2*m+V*rho)*Cd*A.*rho.*v(2).^2))];
end

Risposta accettata

Alan Stevens
Alan Stevens il 31 Lug 2020
I assume this represents some object sinking through changing desity sea water. The following shows how you might do it, though you will need to replace my simple linear density change with whatever you think is appropriate.
tspan=0:1000;
y0 =[0 0];
[t,y]=ode45(@f,tspan,y0);
subplot(2,1,1)
plot(t,y(:,1))
grid on
xlabel('time'),ylabel('depth')
subplot(2,1,2)
plot(t,y(:,2))
grid
xlabel('time'),ylabel('speed')
function speed=f(~,y)
Cd=.16;
A=.0613;
m=82.25;
V=0.0809815;
depth = y(1);
rho = interp1([0 10000],[1030 1070], depth); %This is just a linear interpolation between two values
% Needs to be replaced with
% your own relationship.
g=9.81;
speed=[y(2); (rho*V*g/(m+V*rho))-(m*g/(m+V*rho))-((1/(2*m+V*rho)*Cd*A.*rho.*y(2).^2))];
end
  6 Commenti
Jestoni Orcejola
Jestoni Orcejola il 9 Set 2020
What if I want to count down whith the interpolation. Is that possible?
Alan Stevens
Alan Stevens il 10 Set 2020
Interpolation is just that - it doesn't matter about direction! The interpolation routine will just take the depths you specify and spit out the corresponding densities.

Accedi per commentare.

Più risposte (0)

Categorie

Scopri di più su Splines in Help Center e File Exchange

Tag

Community Treasure Hunt

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

Start Hunting!

Translated by