0

Hi I am working on a script that will solve and plot an ODE using the Runge Kutta method. I want to have the script use different functions so I can expand upon it later. If I write it with out the function definition it works fine, but with the def() it will not open a plot window or print the results. Than you!

from numpy import *
import matplotlib.pyplot as plt

#H=p^2/2-cosq
#p=dp=-dH/dq
#q=dq=dH/dp

t = 0
h = 0.5
pfa = []                        #Create arrays that will hold pf,qf values
qfa = []

while t < 10:
    q = 1*t
    p = -sin(q*t)

    p1 = p
    q1 = q
    p2 = p + h/2*q1
    q2 = q + h/2*p1
    p3 = p+ h/2*q2
    q3 = q+ h/2*p2
    p4 = p+ h/2*q3
    q4 = q+ h/2*p4
    pf = (p +(h/6.0)*(p1+2*p2+3*p3+p4))
    qf = (q +(h/6.0)*(q1+2*q2+3*q3+q4))

    pfa.append(pf)                   #append arrays
    qfa.append(qf)
    t += h                           #increase time step                        

print("test")
plt.plot(pfa,qfa)
print("test1")
plt.show()
print("tes2t")
talonmies
  • 70,661
  • 34
  • 192
  • 269
Stripers247
  • 2,265
  • 11
  • 38
  • 40
  • You aren't calling `rk`, but even if you do, your example will only compute a single point, and the graph will appear, but be empty. – talonmies Feb 18 '12 at 23:48
  • @talonmies Can you explain how I can get it to plot? – Stripers247 Feb 19 '12 at 23:27
  • You are asking `plot()` to draw a line, but only supplying only one point. If you did something like `plot([p,pf],[q,qf])` it would try and plot two points. Even then, it will still be not plotting, because `p` and `pf` and `q` and `qf` will be equal if `h<6`, because you are doing integer calculations and `h/6 = 0` for all positive h less than 6. Also, if this is truly supposed to be a 4 stage explicit Runge Kutta method, it would appear you have the formula incorrect. – talonmies Feb 20 '12 at 06:23
  • @talonmies I edited the code but it is still not plotting – Stripers247 Feb 20 '12 at 06:29
  • Python is a white space sensitive language. You must use correct and consistent indentation in order for the code to work correctly. If the code you are running is as you posted, it contained an empty loop, and still only calculated 1 point. I have corrected the whitespace in an edit. Now it will do something, although the numerical method you have written is, I believe, still incorrect. – talonmies Feb 20 '12 at 06:57
  • @Talonmies Thank you so much!! Being kinda new to this I keep forgetting about the white space thing. I will look into the Runga-Kutta as you suggested. Thanks again! – Stripers247 Feb 20 '12 at 20:31

1 Answers1

1

If you have a function declared, you'll need to call it at some point, e.g.:

def rk(p,q,h):
    pass # your code here

if __name__ == '__main__':
    rk(1,2,1)

Putting the function call within the if __name__ == '__main__' block will ensure that the function is called only when you run the script directly, not when you import it from another script. ( More on that here in case you're interested: What does if __name__ == "__main__": do? )

And here's an even better option; to avoid hard-coding fn args (your real code should have some error-handling for unexpected command-line input):

def rk(p,q,h):
    pass # your code here

if __name__ == '__main__':
    import argparse
    the_parser = argparse.ArgumentParser()
    the_parser.add_argument('integers', type=int, nargs=3)
    args = the_parser.parse_args()
    p,q,h = args.integers
    rk(p,q,h)
Community
  • 1
  • 1
mechanical_meat
  • 163,903
  • 24
  • 228
  • 223