質問をすることでしか得られない、回答やアドバイスがある。

15分調べてもわからないことは、質問しよう!

新規登録して質問してみよう
ただいま回答率
85.50%
Python 3.x

Python 3はPythonプログラミング言語の最新バージョンであり、2008年12月3日にリリースされました。

Python

Pythonは、コードの読みやすさが特徴的なプログラミング言語の1つです。 強い型付け、動的型付けに対応しており、後方互換性がないバージョン2系とバージョン3系が使用されています。 商用製品の開発にも無料で使用でき、OSだけでなく仮想環境にも対応。Unicodeによる文字列操作をサポートしているため、日本語処理も標準で可能です。

Q&A

解決済

1回答

544閲覧

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

huruhi

総合スコア5

Python 3.x

Python 3はPythonプログラミング言語の最新バージョンであり、2008年12月3日にリリースされました。

Python

Pythonは、コードの読みやすさが特徴的なプログラミング言語の1つです。 強い型付け、動的型付けに対応しており、後方互換性がないバージョン2系とバージョン3系が使用されています。 商用製品の開発にも無料で使用でき、OSだけでなく仮想環境にも対応。Unicodeによる文字列操作をサポートしているため、日本語処理も標準で可能です。

0グッド

0クリップ

投稿2019/11/28 09:28

編集2019/11/28 09:37

前提・実現したいこと

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

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

解こうとしている計算式がこちらになります。
イメージ説明

この式を解こうとするとルート内が負になる場合があるので、0≦k≦√(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): のところになります。

Python

1import math 2import cmath 3import numpy as np 4import matplotlib.pyplot as plt 5import matplotlib.animation as animation 6 7# 円周率 8pi = math.pi 9# t = Tmax まで計算 10Tmax = 3.0 11# 時間ステップ 12dt = 0.03 13# ステップ数 14step = int(Tmax/dt) 15# x方向の分点数 16nx = 2**6 17# 領域幅 18L = 2*pi 19# x方向のメッシュ幅 20dx = L / nx 21# x座標 22x = np.arange(0, L, dx) 23# u = u(x)の初期条件 24u = np.sin(x) 25# 拡散係数 26nu = 0.05 27# グラフの格納場所 28ims = [] 29 30 31# 波数を返す 32def get_k(i): 33 if i <= nx/2: 34 return (2*pi/L)*i 35 else: 36 return (2*pi/L)*(i-nx) 37 38 39# 1階導関数を返す 40def differentiate(in_data_F): 41 # 出力データ 42 out_data_F = np.zeros(int(nx/2+1), dtype=np.complex) 43 # 微分 44 for i in range(int(nx/2+1)): 45 k = get_k(i) 46 Re = in_data_F[i].real 47 Im = in_data_F[i].imag 48 out_data_F[i] = np.complex(-k*Im, k*Re) 49 # return 50 return out_data_F 51 52 53# 2/3ルール 54def two_thirds_rule(nl_F): 55 max_k = int(nx/2) 56 for i in range(max_k+1): 57 if i >= max_k*2/3: 58 nl_F[i] = 0 59 60 61def advection(data_F): 62 # 導関数(Fourier) 63 drv_data_F = differentiate(data_F) 64 # 導関数(実) 65 drv_data = np.fft.irfft(drv_data_F) 66 # もとの関数 67 data = np.fft.irfft(data_F) 68 # 非線形項 69 nonlinear = data * drv_data 70 # 非線形項をフーリエ変換 71 nonlinear_F = np.fft.rfft(nonlinear) 72 # 2/3ルール 73 two_thirds_rule(nonlinear_F) 74 # time marching 75 data_F -= dt * nonlinear_F 76 77 78 79 80 81#問題の箇所 82def diffusion(data_F): 83 for i in range(int(nx/2+1)): 84 k = get_k(i) 85 Re = data_F[i].real 86 Im = data_F[i].imag 87 for k in range(0,math.sqrt(1/nu)+1): 88 f1 = math.exp(math.sqrt(1-nu*k**2)*dt) 89 for k in range(math.sqrt(1/nu)+1,max_k+1): 90 f2 = cmath.exp(cmath.sqrt(1-nu*k**2)*dt).real 91 f3 = cmath.exp(cmath.sqrt(1-nu*k**2)*dt).imag 92 93 Re *= f1+f2 94 Im *= f1+f3 95 data_F[i] = np.complex(Re, Im) 96 97 98 99 100 101 102 103# メイン 104def main(): 105 fig = plt.figure() 106 u_F = np.fft.rfft(u) 107 108 for n in range(step): 109 advection(u_F) 110 diffusion(u_F) 111 112 u_plotted = np.fft.irfft(u_F) 113 im = plt.plot(x, u_plotted, 'r') 114 ims.append(im) 115 116 ani = animation.ArtistAnimation(fig, ims, interval=1) 117 ani.save('i.gif', writer='pillow') 118 119 120if __name__ == '__main__': 121 main() 122

試したこと

if文による場合分けなど

気になる質問をクリップする

クリップした質問は、後からいつでもMYページで確認できます。

またクリップした質問に回答があった際、通知やメールを受け取ることができます。

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

hayataka2049

2019/11/28 09:32

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

2019/11/28 09:37

編集いたしました
guest

回答1

0

ベストアンサー

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

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

投稿2019/11/28 09:40

hayataka2049

総合スコア30933

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

huruhi

2019/11/29 14:25 編集

迅速な指摘と回答ありがとうございます! 回答内容の通りにしたところ発生していたエラーは解決しましたが、Re *= f1+f2 の箇所でRuntimeWarning: overflow encountered in double_scalarsのエラーが発生しました。 このエラーはどういう意味を持つのでしょうか? 追記 kの刻み幅が0.01単位なのに対して、Roundでまとめると境界が1単位でまとまることになるので場合分けがうまくいかないようです。
hayataka2049

2019/11/30 07:54

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

あなたの回答

tips

太字

斜体

打ち消し線

見出し

引用テキストの挿入

コードの挿入

リンクの挿入

リストの挿入

番号リストの挿入

表の挿入

水平線の挿入

プレビュー

15分調べてもわからないことは
teratailで質問しよう!

ただいまの回答率
85.50%

質問をまとめることで
思考を整理して素早く解決

テンプレート機能で
簡単に質問をまとめる

質問する

関連した質問