Pythonでの範囲指定がうまくいかない

Question

前提・実現したいこと

現在、Pythonでバーガース方程式についての数値解析を行っています。
下記のURLにあるバーガース方程式のソースコードを元に、発展させた内容を解こうとしているのですがどうにもうまくいきません。

https://qiita.com/moootoko/items/3f89efc0d6aff7c1c066

解こうとしている計算式がこちらになります。

この式を解こうとするとルート内が負になる場合があるので、0≦ｋ≦√(1/nu)　と　√(1/nu)＜k≦max_k　で式を分割し、負になる場合のみを実部と虚部に分けて計算する方法を教えていただけると幸いです。

発生している問題・エラーメッセージ

発生しているエラーは、次になります。

Traceback (most recent call last):
  File "　", line 111, in <module>
    main()
  File "　", line 100, in main
    diffusion(u_F)
  File "　", line 83, in diffusion
    for k in range(0,math.sqrt(1/nu)+1):
TypeError: 'float' object cannot be interpreted as an integer

該当のソースコード

変更したところは　def diffusion(data_F): のところになります。

import math
import cmath
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.animation as animation

# 円周率
pi = math.pi
# t = Tmax まで計算
Tmax = 3.0
# 時間ステップ
dt = 0.03
# ステップ数
step = int(Tmax/dt)
# x方向の分点数
nx = 2**6
# 領域幅
L = 2*pi
# x方向のメッシュ幅
dx = L / nx
# x座標
x = np.arange(0, L, dx)
# u = u(x)の初期条件
u = np.sin(x)
# 拡散係数
nu = 0.05
# グラフの格納場所
ims = []


# 波数を返す
def get_k(i):
    if i <= nx/2:
        return (2*pi/L)*i
    else:
        return (2*pi/L)*(i-nx)


# 1階導関数を返す
def differentiate(in_data_F):
    # 出力データ
    out_data_F = np.zeros(int(nx/2+1), dtype=np.complex)
    # 微分
    for i in range(int(nx/2+1)):
        k = get_k(i)
        Re = in_data_F[i].real 
        Im = in_data_F[i].imag
        out_data_F[i] = np.complex(-k*Im, k*Re)
    # return
    return out_data_F


# 2/3ルール
def two_thirds_rule(nl_F):
    max_k = int(nx/2)
    for i in range(max_k+1):
        if i >= max_k*2/3:
            nl_F[i] = 0


def advection(data_F):
    # 導関数(Fourier)
    drv_data_F = differentiate(data_F)
    # 導関数(実)
    drv_data = np.fft.irfft(drv_data_F)
    # もとの関数
    data = np.fft.irfft(data_F)
    # 非線形項
    nonlinear = data * drv_data
    # 非線形項をフーリエ変換
    nonlinear_F = np.fft.rfft(nonlinear)
    # 2/3ルール
    two_thirds_rule(nonlinear_F)
    # time marching
    data_F -= dt * nonlinear_F





＃問題の箇所
def diffusion(data_F):
    for i in range(int(nx/2+1)):
        k = get_k(i)
        Re = data_F[i].real
        Im = data_F[i].imag
    for k in range(0,math.sqrt(1/nu)+1):
        f1 = math.exp(math.sqrt(1-nu*k**2)*dt)
    for k in range(math.sqrt(1/nu)+1,max_k+1):
        f2 = cmath.exp(cmath.sqrt(1-nu*k**2)*dt).real
        f3 = cmath.exp(cmath.sqrt(1-nu*k**2)*dt).imag

        Re *= f1+f2
        Im *= f1+f3
        data_F[i] = np.complex(Re, Im)







# メイン
def main():
    fig = plt.figure()
    u_F = np.fft.rfft(u)

    for n in range(step):
        advection(u_F)
        diffusion(u_F)

        u_plotted = np.fft.irfft(u_F)
        im = plt.plot(x, u_plotted, 'r')
        ims.append(im)

    ani = animation.ArtistAnimation(fig, ims, interval=1)
    ani.save('i.gif', writer='pillow')


if __name__ == '__main__':
    main()

試したこと

if文による場合分けなど

Comments

[hayataka2049]

エラーをTracebackの先頭から示すよう編集してください。

[huruhi]

編集いたしました

Answer 1

rangeの引数は整数型しか取らないということで、まあこの場合は丸めてあげれば間違いはないでしょう。

range(0, round(math.sqrt(1/nu))+1)とでもしてみてください。

Comments

[huruhi]

迅速な指摘と回答ありがとうございます！
回答内容の通りにしたところ発生していたエラーは解決しましたが、Re *= f1+f2 の箇所でRuntimeWarning: overflow encountered in double_scalarsのエラーが発生しました。
このエラーはどういう意味を持つのでしょうか？

追記
kの刻み幅が0.01単位なのに対して、Roundでまとめると境界が1単位でまとまることになるので場合分けがうまくいかないようです。

[hayataka2049]

＞0.01単位なのに対して
であれば最初からnumpy.arageかlinspaceでデータを作るなどでしょうか。