2013-03-21 72 views
2

我是Python中的初級/中級。我已經編寫了一個4th-order Runge-Kutta method (RK4)到Python。它基本上解決了擺錘,但這不是重點。我希望能夠將函數f直接傳遞給RK4函數,即RK4(y_0,n,h)應該成爲RK4(f,y_0,n) ,H)。這可以帶來很大的好處,我可以將RK4用於描述其他系統的其他f函數,而不僅僅是這一個鐘擺。如何將函數傳遞給Python中的函數?

我曾經玩過簡單的功能到RK4,但我做錯了什麼。我如何在Python中做到這一點?

import numpy as np 

def RK4(y_0, n, h): 
    #4th order Runge-Kutta solver, takes as input 
    #initial value y_0, the number of steps n and stepsize h 
    #returns solution vector y and time vector t 
    #right now function f is defined below 

    t = np.linspace(0,n*h,n,endpoint = False) #create time vector t 
    y = np.zeros((n,len(y_0))) #create solution vector y 
    y[0] = y_0 #assign initial value to first position in y 
    for i in range(0,n-1): 
     #compute Runge-Kutta weights k_1 till k_4 
     k_1 = f(t[i],y[i]) 
     k_2 = f(t[i] + 0.5*h, y[i] + 0.5*h*k_1) 
     k_3 = f(t[i] + 0.5*h, y[i] + 0.5*h*k_2) 
     k_4 = f(t[i] + 0.5*h, y[i] + h*k_3) 
     #compute next y   
     y[i+1] = y[i] + h/6. * (k_1 + 2.*k_2 + 2.*k_3 + k_4) 
    return t,y 

def f(t,vec): 
    theta=vec[0] 
    omega = vec[1] 
    omegaDot = -np.sin(theta) - omega + np.cos(t) 
    result = np.array([omega,omegaDot])  
    return result 

test = np.array([0,0.5]) 
t,y = RK4(test,10,0.1) 

回答

4

Python函數也是對象。您可以通過他們周圍像任何其他對象:

>>> def foo(): print 'Hello world!' 
... 
>>> foo 
<function foo at 0x10c4685f0> 
>>> foo() 
Hello world! 
>>> bar = foo 
>>> bar() 
Hello world! 

簡單地傳遞一個函數作爲一個額外的參數你RK4功能和使用,作爲一個局部變量。

+0

感謝您的澄清。我來自MatLab(與一些C)。所以Python不停地驚訝於:-)。 – seb 2013-03-21 10:57:29

2

這很簡單。更改RK4函數的定義,像這樣:

def RK4(f, y_0, n, h): 

在這裏,我已經添加了一個額外的參數,函數。

然後,當你調用RK4,傳遞函數:

t, y = RK4(f, test, 10, 0.1) 

而現在,當然,您也可以替換不同的功能,而無需重新編寫集成代碼。

Python中的函數只是另一種對象。你可以將它們傳遞給你,就像你做更多的平行對象一樣。

+0

好吧,我把它放在裏面,現在它可以工作。我是一個血腥的新人!無論如何,將它們全部輸入並張貼在這裏是一個很好的練習。 – seb 2013-03-21 10:56:01

+0

是的,很多時候寫一個簡潔的例子有助於你的理解。順便說一下,這是一個寫得很好的問題,做得很好! – 2013-03-21 10:57:33

3

您可以傳遞一個函數的功能在Python,正如你所期望的:

def call_function(f): 
    f() 

def my_function(): 
    print "OK" 

call_function(my_function) # Prints OK 

也許你應該張貼發生故障的代碼?