首页 文章

求解python中的非线性方程

提问于
浏览
0

我试图找到介质波导的基本TE模式 . 我尝试解决它的方法是计算两个函数并尝试在图上找到它们的交集 . 但是,我在绘制情节时遇到了交叉点 . 我的代码:

def LHS(w):
    theta = 2*np.pi*1.455*10*10**(-6)*np.cos(np.deg2rad(w))/(900*10**(-9))
    if(theta>(np.pi/2) or theta < 0):
        y1 = 0
    else:
        y1 = np.tan(theta)
    return y1

def RHS(w):
    y = ((np.sin(np.deg2rad(w)))**2-(1.440/1.455)**2)**0.5/np.cos(np.deg2rad(w))
    return y

x = np.linspace(80, 90, 500)

LHS_vals = [LHS(v) for v in x]
RHS_vals = [RHS(v) for v in x]


# plot
fig, ax = plt.subplots(1, 1, figsize=(6,3))
ax.plot(x,LHS_vals)
ax.plot(x,RHS_vals)
ax.legend(['LHS','RHS'],loc='center left', bbox_to_anchor=(1, 0.5))

plt.ylim(0,20)
plt.xlabel('degree')
plt.ylabel('magnitude')
plt.show()

我有这样的情节:
graph

交叉点大约是89度,但是,我无法计算x的精确值 . 我试过fsolve,求解找到解决方案,但仍然徒劳无功 . 如果它不是唯一的解决方案,它似乎无法打印出解决方案 . 是否有可能只找到x在一定范围内的解?有人可以在这里给我任何建议吗?谢谢!

编辑:方程式是这样的(m = 0):
enter image description here

我试图通过找到交叉点来解决这个问题

编辑:我尝试的方式之一是这样的:

from scipy.optimize import fsolve
def f(wy):
   w, y = wy
   z = np.array([y - LHS(w), y - RHS(w)])
   return z

fsolve(f,[85, 90])

但它给了我错误的答案 . 我也尝试过这样的事情:

import matplotlib.pyplot as plt

x = np.arange(85, 90, 0.1)
f = LHS(x)
g = RHS(x)

plt.plot(x, f, '-')
plt.plot(x, g, '-')

idx = np.argwhere(np.diff(np.sign(f - g)) != 0).reshape(-1) + 0
plt.plot(x[idx], f[idx], 'ro')
print(x[idx])
plt.show()

但它显示:ValueError:具有多个元素的数组的真值是不明确的 . 使用a.any()或a.all()

2 回答

  • 1

    快速和(非常)脏的东西似乎有效(至少它为你的参数提供了〜89的theta值) - 在图之前的代码中添加以下代码,在 RHS_vals = [RHS(v) for v in x] 行之后:

    # build a list of differences between the values of RHS and LHS
    diff = [abs(r_val- l_val) for r_val, l_val in zip(RHS_vals, LHS_vals)]
    
    # find the minimum of absolute values of the differences
    # find the index of this minimum difference, then find at which angle it occured
    min_diff = min(diff)
    print "Minimum difference {}".format(min_diff)
    print "Theta = {}".format(x[diff.index(min_diff)])
    

    我看了85-90的范围:

    x = np.linspace(85, 90, 500)
    
  • 1

    首先,您需要确保该函数可以实际处理numpy数组 . 用于定义分段函数的几个选项显示在Plot Discrete Distribution using np.linspace()中 . 例如 .

    def LHS(w):
        theta = 2*np.pi*1.455*10e-6*np.cos(np.deg2rad(w))/(900e-9)
        y1 = np.select([theta < 0, theta <= np.pi/2, theta>np.pi/2], [np.nan, np.tan(theta), np.nan])
        return y1
    

    这已经允许使用第二种方法,在索引处绘制一个点,该点最接近两条曲线之间的差值的最小值 .

    import numpy as np
    import matplotlib.pyplot as plt
    
    def LHS(w):
        theta = 2*np.pi*1.455*10e-6*np.cos(np.deg2rad(w))/(900e-9)
        y1 = np.select([theta < 0, theta <= np.pi/2, theta>np.pi/2], [np.nan, np.tan(theta), np.nan])
        return y1
    
    def RHS(w):
        y = ((np.sin(np.deg2rad(w)))**2-(1.440/1.455)**2)**0.5/np.cos(np.deg2rad(w))
        return y
    
    x = np.linspace(82.1, 89.8, 5000)
    f = LHS(x)
    g = RHS(x)
    
    plt.plot(x, f, '-')
    plt.plot(x, g, '-')
    
    idx = np.argwhere(np.diff(np.sign(f - g)) != 0).reshape(-1) + 0
    plt.plot(x[idx], f[idx], 'ro')
    plt.ylim(0,40)
    plt.show()
    

    enter image description here

    然后,人们也可以使用 scipy.optimize.fsolve 来获得实际的解决方案 .

    idx = np.argwhere(np.diff(np.sign(f - g)) != 0)[-1]
    
    from scipy.optimize import fsolve
    
    h = lambda x: LHS(x)-RHS(x)
    sol = fsolve(h,x[idx])
    
    plt.plot(sol, LHS(sol), 'ro')
    plt.ylim(0,40)
    plt.show()
    

相关问题