我在 Matlab 中有一段代碼(見下文)運行良好:
clc; clear all;
f = 8500;
c0 = 343;
rho = 1.225;
omega = 2*pi*f;
k = omega/c0;
Z = -426;
lx = 0.1;
ly = 0.1;
nx = 50;
ny = nx/2;
integrand1 = @(x,y,kx) real(((exp(1i*(kx*x sqrt(k.^2 - kx.^2).*abs(y))))./sqrt(k.^2 - kx.^2)));
integrand2 = @(x,y,kx) imag(((exp(1i*(kx*x sqrt(k.^2 - kx.^2).*abs(y))))./sqrt(k.^2 - kx.^2)));
Gz = @(x,y) integral(@(kx)integrand1(x,y,kx), -100*k, 100*k, "AbsTol", 1e-4) 1i*...
integral(@(kx)integrand2(x,y,kx), -100*k, 100*k, "AbsTol", 1e-4);
Test = Gz(lx,ly)
但我需要在 Python 中實作它。我嘗試使用與 lambda 運算式相同的方法將我的函式定義為函式句柄,但不幸的是,它不起作用并且我遇到了錯誤。這是我的 Python 代碼:
import numpy as np
from scipy import integrate
import matplotlib.pyplot as plt
f = 8500
rho = 1.225
c0 = 343
omega = 2*np.pi*f
k = omega/c0
Z = -426
Lx = 0.1
Ly = 0.1
nx = 50
ny = int(nx/2)
integrandReal = lambda kx, x, y: np.real(((2*np.sqrt(k**2 - kx**2)*Z)/(np.sqrt(k**2 - kx**2)*Z omega*rho))*((np.exp(1j*(kx*x np.sqrt(k**2 - kx**2)*y)))/(np.sqrt(k**2 - kx**2))))
integrandImag = lambda kx, x, y: np.imag(((2*np.sqrt(k**2 - kx**2)*Z)/(np.sqrt(k**2 - kx**2)*Z omega*rho))*((np.exp(1j*(kx*x np.sqrt(k**2 - kx**2)*y)))/(np.sqrt(k**2 - kx**2))))
integral = lambda x, y: integrate.quad(integrandReal, -100*k, 100*k, args=(x,y)) 1j*integrate.quad(integrandImag, -100*k, 100*k, args=(x,y))
G = integral(Lx,Ly)
這些是錯誤:
runfile(XYZ/'untitled0.py', wdir='/Users/XYZ ')
/Users/XYZ/untitled0.py:29: RuntimeWarning: invalid value encountered in sqrt
integrandReal = lambda kx, x, y: np.real(((2*np.sqrt(k**2 - kx**2)*Z)/(np.sqrt(k**2 - kx**2)*Z omega*rho))*((np.exp(1j*(kx*x np.sqrt(k**2 - kx**2)*y)))/(np.sqrt(k**2 - kx**2))))
/Users/ /untitled0.py:32: IntegrationWarning: The occurrence of roundoff error is detected, which prevents
the requested tolerance from being achieved. The error may be
underestimated.
integral = lambda x, y: integrate.quad(integrandReal, -100*k, 100*k, args=(x,y)) 1j*integrate.quad(integrandImag, -100*k, 100*k, args=(x,y))
/Users/XYZ/untitled0.py:30: RuntimeWarning: invalid value encountered in sqrt
integrandImag = lambda kx, x, y: np.imag(((2*np.sqrt(k**2 - kx**2)*Z)/(np.sqrt(k**2 - kx**2)*Z omega*rho))*((np.exp(1j*(kx*x np.sqrt(k**2 - kx**2)*y)))/(np.sqrt(k**2 - kx**2))))
Traceback (most recent call last):
File "/Users/XYZ/untitled0.py", line 34, in <module>
G = integral(0.1,0.1)
File "/Users/XYZ/untitled0.py", line 32, in <lambda>
integral = lambda x, y: integrate.quad(integrandReal, -100*k, 100*k, args=(x,y)) 1j*integrate.quad(integrandImag, -100*k, 100*k, args=(x,y))
TypeError: can't multiply sequence by non-int of type 'complex'
任何幫助都非常受歡迎。先感謝您!
uj5u.com熱心網友回復:
scipy.integrate.quad回傳一個元組,而不是一個數字,因此您應該從中選擇第一個回傳值,即數字結果。
并且使用 np.sqrt 將在負輸入時回傳 nan,您應該改用numpy.lib.scimath.sqrt因為它可以處理負輸入并回傳復數,如下所示。
from numpy.lib.scimath import sqrt
integrandReal = lambda kx, x, y: np.real(((2*sqrt(k**2 - kx**2)*Z)/(sqrt(k**2 - kx**2)*Z omega*rho))*((np.exp(1j*(kx*x sqrt(k**2 - kx**2)*y)))/(sqrt(k**2 - kx**2))))
integrandImag = lambda kx, x, y: np.imag(((2*sqrt(k**2 - kx**2)*Z)/(sqrt(k**2 - kx**2)*Z omega*rho))*((np.exp(1j*(kx*x sqrt(k**2 - kx**2)*y)))/(sqrt(k**2 - kx**2))))
integral = lambda x, y: integrate.quad(integrandReal, -100*k, 100*k, args=(x,y))[0] 1j*integrate.quad(integrandImag, -100*k, 100*k, args=(x,y))[0]
轉載請註明出處,本文鏈接:https://www.uj5u.com/yidong/436237.html
標籤:Python 麻木的 matlab scipy 四边形
上一篇:通過std::pair<std::array<std::array<u_int16_t,2>,1>,std::string>>直接讀取參考值
