3

オイラー法とレンジ クッタ法の両方を使用してばね質量の問題を解決し、プロットを比較しようとしています。Euler と Runge-Kutta の両方の関数を作成しましたが、問題に対して関数を呼び出した後、プロットにデータが表示されないようです。プロットを修正して、コードにエラーがないかどうかを確認してください。ありがとう

#function Euler
def euler ( y, t, dt, derivative):
    y_next = y + derivative(y, t) * dt
    return y_next

# function Runge-Kutta
# 2nd order Runge-Kutta method routine

def Runge_Kutta (y, time, dt, derivative):
    k0 = dt * derivative (y, time)
    k1 = dt * derivative (y + k0, time + dt)
    y_next = y + 0.5 * (k0 + k1)
    return y_next

ここに私が解決しようとしている問題があります

[![""" A spring and mass  system. the coefficient of friction \mu is not negligible.generate a position vs. time plot for the motion of the mass, given an initial displacement x = 0.2m , spring constant k = 42 N/m , mass m =0.25 Kg, coefficient of friction \mu = 0.15 and initial velocity v = 0

F = -kx +/-mu mg """


from pylab import *
from Runge_Kutta_routine import Runge_Kutta
from eulerODE import euler

N = 500      #input ("How many number of steps to take?")
x0 = 0.2
v0 = 0.0
tau = 3.0     #input ("What is the total time of the simulation in seconds?")
dt = tau /float ( N-1)

k = 41.0      #input (" what is the spring constant?")
m = 0.25      #input ("what is the mass of the bob?")
gravity = 9.8
mu = 0.15     #input ("what is the coefficient of friction?")

""" we create a Nx2 array for storing the results of our calculations. Each 2- element row will be used for the state of the system at one instant, and each instant is separated by time dt. the first element in each row will denote position, the second would be velocity"""

y = zeros (\[N,2\])

y \[0,0\] = x0
y \[0,1\] = v0

def SpringMass (state, time):
    """ we break this second order DE into two first order DE introducing dx/ dt = v & dv/dt = kx/ m +/- mu g....

Note that the direction of the frictional force changes depending on the sign of the velocity, we handle this with an if statement."""


    g0 = state\[1\]
    if g0 > 0:
        g1 = -k/m * state \[0\] - gravity * mu
    else:
        g1 = -k/m * state \[0\] + gravity * mu

    return array (\[g0, g1\])

# Now we do the calculations
# loop only N-1 so that we don;t run into a problem addresssing y\[N+1\] on the last point

    for j in range (N-1):
        #y \[j+1\] = euler ( y\[j\] , 0, dt, SpringMass)
        y \[j+1\] = Runge_Kutta ( y\[j\], 0 , dt, SpringMass)

# Now we plot the result

time = linspace ( 0 , tau, N)
plot ( time, y\[:,0\], 'b-', label ='position')

xlabel('time')
ylabel('position')


show()][1]][1] 
4

1 に答える 1

1

for j in range (N-1):配列を計算するで始まるループyがインデントされているように見えるため、Python はこれらの行を function の一部と見なしますSpringMass。これらの行はreturnステートメントの後に来るため、実行されることはありません。

これを修正するには、これらの行を移動して、そのfor行にはインデントがなく、他の行には 4 スペースのインデントしかないようにします。問題が解決するかどうかを確認してください。

ここに書かれたコードはまだ機能しないことに注意してください。角かっこの前に不要なバックスラッシュがあり、名前のない単一のモジュールに関数EulerRunge_Kutta関数を記述しますが、メインコードはそれらが2つの異なるモジュールにあることを期待しています。スタイルの悪い例もたくさんあります。これらの問題は、おそらく、(その他の) 回答がまだ得られていない理由です。ここに投稿する前にコードをクリーンアップして、スタイルをクリーンアップしてください

于 2016-05-24T16:48:36.597 に答える