Q&A
plt.scatter (i, v, marker='.')
のすぐ上に
print(i)
print(v)
を(インデントを「plt.scatter」に合わせて)追加して実行してみて、値が正しく計算されてるのかを確認してみたらいかがでしょうか
python
1import numpy as np 2import matplotlib.pyplot as plt 3from numpy.lib.type_check import isrealobj 4 5#関数の定義 6def f1(fai2): 7 return fai2 8 9def f2(fai1,fai2): 10 return (i-np.sin(fai1)-fai2)/bc 11 12#初期値 13a = 0.0 14b = 2.0 15N = 100 16h = (b-a)/N 17j = 0 18 19bc = 2.0 20 21t = np.arange(a,b,h) 22fai1_0 = 0.0 23fai2_0 = 0.0 24 25#配列 26fai1 = np.empty(N) 27fai2 = np.empty(N) 28fai1[0] = fai1_0 29fai2[0] = fai2_0 30 31#iを増加させたとき 32for i in np.arange(0,2.0,0.01): 33 34 #ルンゲクッタ法 35 for j in range (N-2): 36 37 k1 = h * f1(fai2[j]) 38 d1 = h * f2(fai1[j] , fai2[j]) 39 40 k2 = h * f1(fai2[j] + d1/2) 41 d2 = h * f2(fai1[j] + k1/2, fai2[j] + d1/2) 42 43 k3 = h * f1(fai2[j] + d2/2) 44 d3 = h * f2(fai1[j] + k2/2 , fai2[j] +d2/2) 45 46 k4 = h * f1(fai2[j] + d3) 47 d4 = h * f2(fai1[j] + k3 , fai2[j] + d3) 48 49 fai1[j+1] = fai1[j] + 1/6 * (k1 + 2 * k2 + 2 * k3 + k4) 50 fai2[j+1] = fai2[j] + 1/6 * (d1 + 2 * d2 + 2 * d3 + d4) 51 52 if j == N-1: 53 break 54 55 j += 1 56 57 #足した値を平均化する 58 v = np.average(fai2[6000:]) 59 60 plt.scatter (i, v, marker='.') 61 62#iを減少させたとき 63for i in np.arange(2.0,0,-0.01): 64 65 #ルンゲクッタ法 66 for j in range (N-2): 67 68 k1 = h * f1(fai2[j]) 69 d1 = h * f2(fai1[j] , fai2[j]) 70 71 k2 = h * f1(fai2[j] + d1/2) 72 d2 = h * f2(fai1[j] + k1/2, fai2[j] + d1/2) 73 74 k3 = h * f1(fai2[j] + d2/2) 75 d3 = h * f2(fai1[j] + k2/2 , fai2[j] +d2/2) 76 77 k4 = h * f1(fai2[j] + d3) 78 d4 = h * f2(fai1[j] + k3 , fai2[j] + d3) 79 80 fai1[j+1] = fai1[j] + 1/6 * (k1 + 2 * k2 + 2 * k3 + k4) 81 fai2[j+1] = fai2[j] + 1/6 * (d1 + 2 * d2 + 2 * d3 + d4) 82 83 if j == N-1: 84 break 85 86 j += 1 87 88 #足した値を平均化する 89 v = np.average(fai2[6000:]) 90 91 plt.scatter (i, v, marker='.') 92 93plt.xlabel("v") 94plt.ylabel("i") 95plt.show() 96```ルンゲクッタ法を用いたコードなのですが、横軸をv、縦軸をiとして,ルンゲクッタ法を用いて計算した6000以降の値を平均した値を、iを０から2.0までの間で繰り返し、それぞれ得た値をプロットしていくプログラムなのですが、実行してもがvが発散してnanとなり、グラフにプロットされません。 97iをループせずに固定して、その値ひとつをプロットしたときはできたのですが、i をループさせた途端にできなくなってしまいました。 98また、iの範囲は０から2.0までと指定しているのに、結果として出てくるグラフは―0.04から0.04になってしまいます。どうすればよいでしょうか。
確認してみたところ、iは正しくループできていましたが、vが発散してnanとなってしまっていました。
iをループせずに固定して、vの点をひとつプロットしたときはできたのですが、なぜiをループさせたらできないのでしょうか、、
> vが発散してnanとなってしまっていました。
それでしたら、
> なぜか点がグラフにプロットされません。
という質問ではなく、「v」がうまく計算できません、という質問ですね
上記分かったことを、質問を編集して追記することをお勧めします
わかりました。
ありがとうございます！
N = 100 としていますが、これですと np.average(fai2[6000:]) で範囲外アクセスになってしまいます。
> np.average(fai2[6000:]) で範囲外アクセス
は、
> iをループせずに固定して、その値ひとつをプロットしたときはできた
の場合もダメなはずですよね
N ＝10000にしたらできました！
ありがとうございます！
回答1件
0
自己解決
python
1import numpy as np 2import matplotlib.pyplot as plt 3from numpy.lib.type_check import isrealobj 4 5#関数の定義 6def f1(fai2): 7 return fai2 8 9def f2(fai1,fai2): 10 return (i-np.sin(fai1)-fai2)/bc 11 12#初期値 13a = 0.0 14b = 2.0 15N = 10000 16h = (b-a)/N 17j = 0 18 19bc = 2.0 20 21t = np.arange(a,b,h) 22fai1_0 = 0.0 23fai2_0 = 0.0 24 25#配列 26fai1 = np.empty(N) 27fai2 = np.empty(N) 28fai1[0] = fai1_0 29fai2[0] = fai2_0 30 31#iを増加させたとき 32for i in np.arange(0,2.0,0.01): 33 34 #ルンゲクッタ法 35 for j in range (N-2): 36 k1 = h * f1(fai2[j]) 37 d1 = h * f2(fai1[j] , fai2[j]) 38 k2 = h * f1(fai2[j] + d1/2) 39 d2 = h * f2(fai1[j] + k1/2, fai2[j] + d1/2) 40 k3 = h * f1(fai2[j] + d2/2) 41 d3 = h * f2(fai1[j] + k2/2, fai2[j] +d2/2) 42 k4 = h * f1(fai2[j] + d3) 43 d4 = h * f2(fai1[j] + k3 , fai2[j] + d3) 44 fai1[j+1] = fai1[j] + 1/6 * (k1 + 2 * k2 + 2 * k3 + k4) 45 fai2[j+1] = fai2[j] + 1/6 * (d1 + 2 * d2 + 2 * d3 + d4) 46 47 if j == N-1: 48 break 49 50 j += 1 51 52 #足した値を平均化する 53 v = np.average(fai2[6000:]) 54 55 print(v) 56 print(i) 57 plt.scatter (v, i, marker='.') 58 59#プロット 60plt.xlim(0, 1.8) 61plt.ylim(0, 2.0) 62plt.xlabel("v") 63plt.ylabel("i") 64plt.show()
投稿2021/12/22 04:58
総合スコア27
あなたの回答
tips
プレビュー
