1

I with difficult to understand how the function rhs(u, m, r) from the code below is receiving the parameters m and r. As it is possible to see from the code, the function rhs is called inside the function euler_step(u, rhs, dt), however neither the parameters m and r is passed as arguments to the function euler_step nor they are global variable. So, someone can explain to me how the parameters m and u are arriving go the function rhs.

# model parameters:
mpo = 100.   # initial mass of the rocket propellant in kg
ms = 50.     # mass of the rocket shell in kg
g = 9.81     # gravity in m s^{-2}
rho = 1.091  # average air density in kg/m^{3}
rad = 0.5    # radius of the maximum cross sectional area of the rocket in m
A = numpy.pi*(rad**2)# maximum cross sectional area of the rocket in m^{2}
v_e = 325.   # the exhaust speed in m/s
C_D = 0.15   # drag coefficient
rt = 20.0    # propellant burn rate in kg/s
dtp = 5.0    # time interval to empty the propellant in s

### set initial conditions ###
h0 = 0.0 # start at the zero height [m]
v0 = 0.0 # initial speed [m/s]

def rhs(u, m, r):
    """Returns the right-hand side of the phugoid system of equations.

    Parameters
    ----------
    u : array of float
        array containing the solution at time n.
    mp: float
        mass of the propellant at time t
    mp_rate: float
        propellant burn rate

    Returns
    -------
    dudt : array of float
        array containing the RHS given u.
    """
    print("[m,r]",[m,r])
    [h,v] = u.copy()

    return numpy.array( [ v, -g + pow((ms+m),-1)*(r*v_e - 0.5*rho*v*abs(v)*A*C_D) ] )

def euler_step(u, rhs, dt):
    """Returns the solution at the next time-step using Euler's method.

    Parameters
    ----------
    u : array of float
        solution at the previous time-step.
    rhs : function
        function to compute the right hand-side of the system of equation.
    dt : float
        time-increment.

    Returns
    -------
    u_n_plus_1 : array of float
        approximate solution at the next time step.
    """

    return u + dt * rhs(u, m, r)

if __name__ == "__main__":
    T = 17.0                     # final time
    dt = 0.1                     # time increment
    t = numpy.arange(0.0, T, dt) # time discretization
    N = len(t)                   # number of time-steps

    # initialize the array containing the solution for each time-step
    u = numpy.zeros((N, 2))
    u[0] = numpy.array([h0, v0]) # fill 1st element with initial values
    rate = numpy.zeros(N)
    mp = numpy.zeros(N)
    Np = int(((N)/(T))*dtp) # number of time-steps with propellant burn

    rate[0:Np] = rt # propellant burn rate in kg/s
    mp[0:Np] = mpo - rt*t[0:Np]

    # time loop - Euler method
    for n in range(1,N-1):

        r = rate[n]
        m = mp[n]
        print("[R,M]",[r,m])
        u[n+1] = euler_step(u[n], rhs, dt)

Thanks in advance.

1
  • euler_step is passing the global variable m and r to rhs function. Commented Mar 9, 2016 at 16:32

3 Answers 3

4

m and n are globals.

It can be confusing because it may seem that __main__ was a function but it is not. The if __name__ == "__main__" .... is running on a global scope.

Sign up to request clarification or add additional context in comments.

1 Comment

Other question: Why it is used [h,v] = u.copy(), to the get value of u and not simply a normal assignment [h,v] = u???
3

They are global variables. In Python, if, while and for do not create a separate variable scope, so they're still assigned values in global/module scope prior to first invocation of euler_step:

if __name__ == "__main__":  # does not start a new variable scope
    ...

    for n in range(1,N-1):  # does not start one either
        # thus these variables are set in global scope.
        r = rate[n]
        m = mp[n]

        # and euler_step is invoked only here, thus it will see
        # r and m being set.
        u[n+1] = euler_step(u[n], rhs, dt)

See also Short Description of the Scoping Rules?.

Comments

1

m and r are defined near the bottom of your script at the module level:

r = rate[n]
m = mp[n]

Therefore they are available to all functions within the module:

The following are blocks: a module, a function body, and a class definition. A scope defines the visibility of a name within a block. If a local variable is defined in a block, its scope includes that block. If the definition occurs in a function block, the scope extends to any blocks contained within the defining one, unless a contained block introduces a different binding for the name... When a name is used in a code block, it is resolved using the nearest enclosing scope.

Comments

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.