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

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

新規登録して質問してみよう
ただいま回答率
85.45%
Matplotlib

MatplotlibはPythonのおよび、NumPy用のグラフ描画ライブラリです。多くの場合、IPythonと連携して使われます。

Python

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

Q&A

解決済

1回答

1735閲覧

【Python】極座標上へのベクトル図作図

Masa06

総合スコア2

Matplotlib

MatplotlibはPythonのおよび、NumPy用のグラフ描画ライブラリです。多くの場合、IPythonと連携して使われます。

Python

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

0グッド

0クリップ

投稿2020/08/28 11:51

編集2020/08/28 12:27

質問内容

下記作成中のコードにて,極座標(r-θ平面)上にベクトル図を作図しようと考えています.
r,θを変数とし,各座標の値を代入して,その場の速度を算出してベクトル図を描ければと考えましたが上手くいきませんでした.

下記エラーコードが出たため,その対策も調べ確認してみましたが,上手くいかずご質問させていただきました.
アドバイスなど頂けるとありがたいです.

作成コード

import cmath import numpy as np import scipy from scipy.integrate import quad import matplotlib.pyplot as plt #積分1 def M(func, a, b, **kwargs): def real_func(r): return scipy.real(func(r)) def imag_func(r): return scipy.imag(func(r)) real_integral = quad(real_func, a, b, **kwargs) imag_integral = quad(imag_func, a, b, **kwargs) return (real_integral[0] + 1j*imag_integral[0], real_integral[1:], imag_integral[1:]) #複素積分計算1 def F_rc(func, a, b, **kwargs): def real_func(r): return scipy.real(func(r)) def imag_func(r): return scipy.imag(func(r)) real_integral = quad(real_func, a, b, **kwargs) imag_integral = quad(imag_func, a, b, **kwargs) return (real_integral[0] + 1j*imag_integral[0], real_integral[1:], imag_integral[1:]) #複素積分計算2 def G_rc(func, a, b, **kwargs): def real_func(r): return scipy.real(func(r)) def imag_func(r): return scipy.imag(func(r)) real_integral = quad(real_func, a, b, **kwargs) imag_integral = quad(imag_func, a, b, **kwargs) return (real_integral[0] + 1j*imag_integral[0], real_integral[1:], imag_integral[1:]) #複素積分計算3 def F_r2(func, a, b, **kwargs): def real_func(r): return scipy.real(func(r)) def imag_func(r): return scipy.imag(func(r)) real_integral = quad(real_func, a, b, **kwargs) imag_integral = quad(imag_func, a, b, **kwargs) return (real_integral[0] + 1j*imag_integral[0], real_integral[1:], imag_integral[1:]) #複素積分計算4 def G_r2(func, a, b, **kwargs): def real_func(r): return scipy.real(func(r)) def imag_func(r): return scipy.imag(func(r)) real_integral = quad(real_func, a, b, **kwargs) imag_integral = quad(imag_func, a, b, **kwargs) return (real_integral[0] + 1j*imag_integral[0], real_integral[1:], imag_integral[1:]) #Vr,Vtheta変換 def pol2cart(th,v_th,v_r): return v_r*np.cos(th) - v_th*np.sin(th), v_r*np.sin(th) + v_th*np.cos(th) RadiusColumn = 0.12 NumberRadii = 10 NumberThetas = 20 theta = np.linspace(0,2*np.pi,NumberThetas,endpoint=False) radius = np.linspace(RadiusColumn, 10 * RadiusColumn, NumberRadii)[:,None] f = plt.figure() ax = f.add_subplot(111, polar=True) r0 = 0.14 r1 = 0.11 r2 = 0.096 r3 = 0.0384 vr2 = 1.65 a_a = 0.14 +0.06*1j beta2 = cmath.pi / 9 time = 360 n = 1 Q = 1 ganma = - Q / cmath.tan(beta2) omega = 1.1e-06 +1.6e-05 * 1j A_M = M(lambda r: r2 / (r * cmath.sin(beta2)), r2, r1) A_F_rc = F_rc(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (r / radius)**(n + 1), radius, r3) A_G_rc = G_rc(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (radius / r)**(n - 1), r3, radius) A_F_r2 = F_r2(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (r / r2)**(n + 1), r3, r2) A_G_r2 = G_r2(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (r2 / r)**(n - 1), r2, r2) vr_tilde = (-1j *((0.5 * A_F_rc[0]) + (0.5 * A_G_rc[0]))) - (r2 / radius)**(n + 1) * (-1j *((0.5 * A_F_r2[0]) + (0.5 * A_G_r2[0]))) - a_a *((radius / r2)**(n-1) - (r2 / radius)**(n+1) * cmath.exp(-2 * 1j * beta2)) vtheta_tilde = ((0.5 * A_F_rc[0]) - (0.5 * A_G_rc[0])) - 1j * (r2/radius)**(n + 1) * (-1j *((0.5 * A_F_r2[0]) + (0.5 * A_G_r2[0]))) - 1j * a_a * ((radius / r2)**(n - 1) - (r2 / radius)**(n + 1) * cmath.exp(-2 * 1j * beta2)) VelocityRadius = Q/(2* cmath.pi * radius) + vr_tilde * cmath.exp(1j * (omega * time - n * theta)) VelocityTheta = ganma/(2* cmath.pi * radius) + vtheta_tilde * cmath.exp(1j * (omega * time - n * theta)) TotalVelocity = np.linalg.norm([VelocityRadius.real, VelocityTheta.real],axis=0) VelocityX,VelocityY = pol2cart(theta,VelocityTheta.real ,VelocityRadius.real) ax.quiver(theta,radius, VelocityX/TotalVelocity, VelocityY/TotalVelocity) plt.show()

エラーコード

ValueError Traceback (most recent call last) <ipython-input-1-ee1d8f65e98c> in <module> 86 87 A_M = M(lambda r: r2 / (r * cmath.sin(beta2)), r2, r1) ---> 88 A_F_rc = F_rc(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (r / radius)**(n + 1), radius, r3) 89 A_G_rc = G_rc(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (radius / r)**(n - 1), r3, radius) 90 A_F_r2 = F_r2(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (r / r2)**(n + 1), r3, r2) <ipython-input-1-ee1d8f65e98c> in F_rc(func, a, b, **kwargs) 22 def imag_func(r): 23 return scipy.imag(func(r)) ---> 24 real_integral = quad(real_func, a, b, **kwargs) 25 imag_integral = quad(imag_func, a, b, **kwargs) 26 return (real_integral[0] + 1j*imag_integral[0], real_integral[1:], imag_integral[1:]) ~/opt/anaconda3/lib/python3.7/site-packages/scipy/integrate/quadpack.py in quad(func, a, b, args, full_output, epsabs, epsrel, limit, points, weight, wvar, wopts, maxp1, limlst) 336 337 # check the limits of integration: \int_a^b, expect a < b --> 338 flip, a, b = b < a, min(a, b), max(a, b) 339 340 if weight is None: ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()

出力結果

出力結果

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

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

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

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

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

can110

2020/08/28 11:53

具体的になにどううまくいかないのかを記載ください。
Masa06

2020/08/28 12:29

質問文に追記したためご確認ください.
guest

回答1

0

自己解決

for文に書き換えたら上手く動きました.

for ra in range(len(radius)): for th in range(len(theta)): A_M = M(lambda r: r2 / (r * (cmath.sin(beta2))**2), r2, r1) A_F_rc = F_rc(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (r / radius[ra])**(n + 1), r3, radius[ra]) A_G_rc = G_rc(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (radius[ra] / r)**(n - 1), radius[ra], r2) A_F_r2 = F_r2(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (r / r2)**(n + 1), r3, r2) A_G_r2 = G_r2(lambda r: cmath.exp((-1j* (omega * cmath.pi * (r**2 - r2**2) / Q - n * ganma * cmath.log(r / r2) / Q ))) * (r2 / r)**(n - 1), r2, r2) vr_tilde = (-1j *((0.5 * A_F_rc[0]) + (0.5 * A_G_rc[0]))) - (r2 / radius[ra])**(n + 1) * (-1j *((0.5 * A_F_r2[0]) + (0.5 * A_G_r2[0]))) - a31_a32 *((radius[ra] / r2)**(n-1) - (r2 / radius[ra])**(n+1) * cmath.exp(-2 * 1j * beta2)) vtheta_tilde = ((0.5 * A_F_rc[0]) - (0.5 * A_G_rc[0])) - 1j * (r2/radius[ra])**(n + 1) * (-1j *((0.5 * A_F_r2[0]) + (0.5 * A_G_r2[0]))) - 1j * a31_a32 * ((radius[ra] / r2)**(n - 1) - (r2 / radius[ra])**(n + 1) * cmath.exp(-2 * 1j * beta2)) VelocityRadius =Q/(2* cmath.pi * radius[ra]) + vr_tilde * cmath.exp(1j * (omega * time - n * theta[th])) VelocityTheta = ganma/(2* cmath.pi * radius[ra]) + vtheta_tilde * cmath.exp(1j * (omega * time - n * theta[th])) TotalVelocity = np.linalg.norm([VelocityRadius.real, VelocityTheta.real],axis=0) ax.quiver(theta[th], radius[ra], VelocityRadius.real * cmath.cos(theta[th]) - VelocityTheta.real * cmath.sin(theta[th]), VelocityRadius.real * cmath.sin(theta[th]) + VelocityTheta.real * cmath.cos(theta[th])) print(math.degrees(math.atan2(VelocityRadius.real, VelocityTheta.real))) plt.show()

投稿2020/08/29 14:54

Masa06

総合スコア2

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

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

あなたの回答

tips

太字

斜体

打ち消し線

見出し

引用テキストの挿入

コードの挿入

リンクの挿入

リストの挿入

番号リストの挿入

表の挿入

水平線の挿入

プレビュー

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

ただいまの回答率
85.45%

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

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

質問する

関連した質問