Variable is in a function file, how do I use a loop to change the value of this variable?

4 visualizzazioni (ultimi 30 giorni)
If the variable is in the main program, I can compile a loop to change the value of a variable and get the result with the variable, but what if the variable is in the function file?
it means i want to use loop to get the plot that zmax changes as the cd or cm changes.
main program are as this
clc;
clear all;
close all;
step=0.1;
t0=50;
length1=t0/step;
opts = odeset('RelTol',1e-4,'AbsTol',1e-4);
Z=zeros(length1,2);
[T,Z]=ode45('vdp', [0:step:t0], [0,1], opts);
figure(1)
plot(T,Z(:,1),'r.');
zend=length(Z(:,1));
zmax=max(Z(20:zend,1))%get the plot that zmax changes as the cd or cm changes
function file is as this
function dz = vdp(t,z)
g=9.8;
gamma=10094;
d=20;
J0=15;
J1=5;
c1=7680;
r=0.21;
m=13500;
G=m*g;
c2=32000;
m1=13500;
m2=0;
%************************************************************
%************************************************************
cd=0.8;%which I want to use the loop to change
%************************************************************
%************************************************************
k0=2000;
p0=1000;
area=15.9;
%************************************************************
%************************************************************
ca=1;%which I want to use the loop to change
%************************************************************
%************************************************************
h=0.5;
t0=3.5;
d1=0.8283479061205706;
height=1;
rho=1025;
dz=zeros(2,1);
dz(1)=z(2);
if z(2)>=0
a1x=1;
else a1x=0;
end
if z(2)>=0
Jx=J0;
else Jx=J1;
end
if z(2)>=0
b1x=0;
else b1x=0;
end
dz(2)=...
...
-(p0/3 + g*m + (k0*z(1))/9 - area*g*rho*(d1 - z(1) + (h*cos((2*pi*t)/t0))/2) + (a1x*c2*z(2))/r + (area*cd*rho*abs(z(2) + (pi*h*sin((2*pi*t)/t0))/t0)*(z(2) + (pi*h*sin((2*pi*t)/t0))/t0))/2 + (2*pi^2*area*h*rho*cos((2*pi*t)/t0)*(d1 - z(1) + (h*cos((2*pi*t)/t0))/2))/t0^2 + (2*pi^2*area*ca*h*rho*cos((2*pi*t)/t0)*(d1 - z(1) + (h*cos((2*pi*t)/t0))/2))/t0^2)/(m + Jx/r^2 + area*ca*rho*(d1 - z(1) + (h*cos((2*pi*t)/t0))/2));
end

Risposte (2)

Jan
Jan il 11 Giu 2019
Do not redefine the important Matlab function cd as a variable, because this can cause serious troubles during debugging.
[T,Z]=ode45('vdp', [0:step:t0], [0,1], opts);
Providing a function to be integrated as a char vector is outdated for over 15 years now. Although this is still supported for backwar compatibility, it is less powerful and secure than a function handle. In addition a function handle allows to define a parameter:
p = 3.1415; % Your parameter
[T,Z] = ode45(@(t,y) vdp(t, y, p), 0:step:t0, [0,1], opts);
function dy = vdp(t, y, p)
...
end
By the way, o:step:to is a vector already. Using the Matlab operator for array concatenation [] is a waste of time here only.
Another problem can cause invalid outputs: Your function to be integrated is not smooth:
if z(2)>=0
Jx=J0;
else Jx=J1;
end
Matlab's ODE integrators are designed to process smooth functions only. Otherwise the step size controller can drive completely mad, which can cause different effects: A stop with an error message (minimal step size exceeded), or even worse a result, which is dominated by rounding errors due to a massive increase of the evaluated steps. Don't do this. Do not use a numerical tool for a computation, which does not match the specified limitations.
Use event functions instead, which stop the integration when the state of the function changes.
  1 Commento
dcydhb dcydhb
dcydhb dcydhb il 12 Giu 2019
why we can't use the
if z(2)>=0
Jx=J0;
else Jx=J1;
end
and can you clarify the event functions usage by an example?
thanks a lot!!!

Accedi per commentare.


Rik
Rik il 11 Giu 2019
Modificato: Rik il 11 Giu 2019
According to the doc you can pass extra parameters by setting an anonymous function. I see no indication in the R2014a doc that this wouldn't work in that release.
If you make the changes indicated below, you should be able to have cd and ca as inputs, so you can change them in a loop.
cd=0.7;
ca=1.1;
[T,Z]=ode45(@(t,z) vdp(t,z,cd,ca), [0:step:t0], [0,1], opts);
function dz = vdp(t,z,cd,ca)
if nargin<4
%set defaults
ca=1;
if nargin<3
cd=0.8;
end
end
%%
%rest of your function
%%
end

Categorie

Scopri di più su Loops and Conditional Statements in Help Center e File Exchange

Prodotti


Release

R2014a

Community Treasure Hunt

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

Start Hunting!

Translated by