Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Matlab Passing variables out of ODE45

  1. Jan 18, 2017 #1
    Hi all,

    I have the following ode code, and I want to know how I could output the variable "k_mb" within my main script so that I can plot it to see if it's actually active when the solver is integrating.
    odescript
    Code (Matlab M):

    function zdot = odenesttest(t,z)
    global r_p r_g I_p I_g m_p m_g V k_mb
    A_or_B = [5,10];
    V = [0,5,15,20,30,35,45,50,60];
    %% state space equations are defined below:
    zdot(2) = ((r_p*k_mb)/I_p)*(-r_p*z(1)+r_g*z(3)+z(5)-z(7));  % Eq 1ss
    zdot(6) = (k_mb/m_p)*(r_p*z(1)-r_g*z(3)-z(5)+z(7));         % Eq 2ss
    zdot(4) = ((r_g*k_mb)/I_g)*(-r_g*z(3)+r_p*z(1)+z(7)-z(5));  % Eq 3ss
    zdot(8) = (k_mb/m_g)*(-r_p*z(1)+r_g*z(3)-z(7)+z(5));        % Eq 4ss

    zdot(1) = z(2); % Eq 5ss
    zdot(3) = z(4); % Eq 6ss
    zdot(5) = z(6); % Eq 7ss
    zdot(7) = z(8); % Eq 8ss

    zdot(9) = z(1); % Pinion angular displacement (control angle).
    torsional_stiffness_changer
       function torsional_stiffness_changer
           if 0 <= zdot(9)
               lower_index = find(V < zdot(9), 1, 'last');
               even_or_odd = mod(lower_index,2);
               k_mb = A_or_B(2 - even_or_odd);
           end
       end
    % Outputs
    zdot = zdot';
    end
     
    Mainscript
    Code (Matlab M):

    % Mainscript for vibration analysis (19/01/2017).
    clear all;
    clc;
    %% Defining global variables
    global r_p r_g I_p I_g m_p m_g k_mb
    k_mb = 50;  % torsional stiffness.
    r_p = 10;   % pinion pitch(?) radius
    r_g = 15;   % gear pitch(?) radius
    I_p = 100;  % pinion inertia.
    I_g = 200;  % gear inertia.
    m_p = 50;   % pinion mass
    m_g = 100;  % gear mass
    %==========================================================================
    %% initial coditions for ode solver
    %==========================================================================
    z0 = [1,0,0,0,0,0,0,0,0]; % zero initial conditions.
    t_step = [0:1:10];        % timestep 0 to 5 seconds.
    opts = odeset('Reltol', 1e-3, 'AbsTol', 1e-3);  % conservative tolerances.
    [t,z] = ode45(@odenesttest, t_step, z0, opts);
    %==========================================================================
    %% Screen Output and Plotting
    %==========================================================================
    [t,z(:,:)]
     
    Thanks a lot!
     
  2. jcsd
  3. Jan 19, 2017 #2

    Dr Transport

    User Avatar
    Science Advisor
    Gold Member

    if it is a global, you should be able to see it all the time, write it to a file or vector and update the plot....
     
  4. Jan 19, 2017 #3
    Sorry, in the mainscript code that I listed I meant to take out the "k_mb = 50;". I have edited it and put it below.

    Updated mainscript
    Code (Text):


    % Mainscript for vibration analysis (19/01/2017).
    clear all;
    clc;
    %% Defining global variables
    global r_p r_g I_p I_g m_p m_g k_mb
    r_p = 10;   % pinion pitch(?) radius
    r_g = 15;   % gear pitch(?) radius
    I_p = 100;  % pinion inertia.
    I_g = 200;  % gear inertia.
    m_p = 50;   % pinion mass
    m_g = 100;  % gear mass
    %==========================================================================
    %% initial coditions for ode solver
    %==========================================================================
    z0 = [1,0,0,0,0,0,0,0,0]; % zero initial conditions.
    t_step = [0:1:10];        % timestep 0 to 5 seconds.
    opts = odeset('Reltol', 1e-3, 'AbsTol', 1e-3);  % conservative tolerances.
    [t,z] = ode45(@odenesttest, t_step, z0, opts);
    %==========================================================================
    %% Screen Output and Plotting
    %==========================================================================
    [t,z(:,:)]
     
     
    Last edited: Jan 19, 2017
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Have something to add?
Draft saved Draft deleted