Scipy.optimize.fsolve() 在我的复杂方程上

Question

我尝试用 scipy 求解复杂函数。我读到 fsolve 只适用于实数。所以我对我的公式做了一些修改，现在它看起来像这样：

import numpy as np
from scipy.optimize import fsolve

delta = 84.37228318652858 * np.pi / 180
psi = 55.2217535040673 * np.pi / 180
n_S = 2.6726 + 3.0375j
phi_i = 70 * np.pi / 180
d_L = 300  # thickness of layer in nm
n_air = 1  # refractive index of air


def snell(phi, n1, n2):
    phi_ref = np.arcsin((n1 / n2) * np.sin(phi))
    return phi_ref


def fresnel(n1, phi1, n2, phi2):
    rs = (n1 * np.cos(phi1) - n2 * np.cos(phi2)) / (
        n1 * np.cos(phi1) + n2 * np.cos(phi2))
    rp = (n2 * np.cos(phi1) - n1 * np.cos(phi2)) / (
        n2 * np.cos(phi1) + n1 * np.cos(phi2))
    return rs, rp


def calc_rho(n_k):
    n = n_k[0]
    k = n_k[1]
    n_L = n + 1j * k
    phi_L = snell(phi_i, n_air, n_L)
    phi_S = snell(phi_L, n_L, n_S)
    rs_al, rp_al = fresnel(n_air, phi_i, n_L, phi_L)
    rs_ls, rp_ls = fresnel(n_L, phi_L, n_S, phi_S)
    beta = 2 * np.pi * d_L * n_L * np.cos(phi_L) / lambda_vac
    rp_L = (rp_al + rp_ls * np.exp(-2 * 1j * beta)) / (
        1 + rp_al * rp_ls * np.exp(-2 * 1j * beta))
    rs_L = (rs_al + rs_ls * np.exp(-2 * 1j * beta)) / (
        1 + rs_al * rs_ls * np.exp(-2 * 1j * beta))
    rho_L = rp_L / rs_L
    return abs(rho_L - rho_given)


lambda_vac = 300
n_S = 2.6726 + 3.0375j
rho_given = np.tan(psi) * np.exp(
    1j * delta)  # should be 1.4399295435287844+0.011780279522394433j
initialGuess = [1.5, 0.1]
nsolve = fsolve(calc_rho, initialGuess)
print(nsolve)

这是我第一次使用求解器（scipy），但我认为这是正确的使用方法。

以上代码出现以下错误：

TypeError: fsolve: there is a mismatch between the input and output shape of the 'func' argument 'calc_rho'.Shape should be (2,) but it is (1,).

我不知道如何解决这个错误，也不知道如何得到一个复数 n 来求解我的方程。

Answer 1

正如评论中讨论的那样，minimize似乎更适合这种问题。以下工作并成功终止：

import numpy as np
from scipy.optimize import minimize

delta = 84.37228318652858 * np.pi / 180
psi = 55.2217535040673 * np.pi / 180
n_S = 2.6726 + 3.0375j
phi_i = 70 * np.pi / 180
d_L = 300  # thickness of layer in nm
n_air = 1  # refractive index of air


def snell(phi, n1, n2):
    phi_ref = np.arcsin((n1 / n2) * np.sin(phi))
    return phi_ref


def fresnel(n1, phi1, n2, phi2):
    rs = (n1 * np.cos(phi1) - n2 * np.cos(phi2)) / (
        n1 * np.cos(phi1) + n2 * np.cos(phi2))
    rp = (n2 * np.cos(phi1) - n1 * np.cos(phi2)) / (
        n2 * np.cos(phi1) + n1 * np.cos(phi2))
    return rs, rp


def calc_rho(n_k):
    n = n_k[0]
    k = n_k[1]
    n_L = n + 1j * k
    phi_L = snell(phi_i, n_air, n_L)
    phi_S = snell(phi_L, n_L, n_S)
    rs_al, rp_al = fresnel(n_air, phi_i, n_L, phi_L)
    rs_ls, rp_ls = fresnel(n_L, phi_L, n_S, phi_S)
    beta = 2 * np.pi * d_L * n_L * np.cos(phi_L) / lambda_vac
    rp_L = (rp_al + rp_ls * np.exp(-2 * 1j * beta)) / (
        1 + rp_al * rp_ls * np.exp(-2 * 1j * beta))
    rs_L = (rs_al + rs_ls * np.exp(-2 * 1j * beta)) / (
        1 + rs_al * rs_ls * np.exp(-2 * 1j * beta))
    rho_L = rp_L / rs_L
    return abs(rho_L - rho_given)


lambda_vac = 300
n_S = 2.6726 + 3.0375j
rho_given = np.tan(psi) * np.exp(
    1j * delta)  # should be 1.4399295435287844+0.011780279522394433j
initialGuess = [1.5, 0.1]
nsolve = minimize(calc_rho, initialGuess, method='Nelder-Mead')
print(nsolve)

这将打印：

 final_simplex: (array([[1.50419259, 0.10077416],
       [1.50419591, 0.10076434],
       [1.50421181, 0.1007712 ]]), array([0.00015298, 0.00015507, 0.00022935]))
           fun: 0.00015297827023618208
       message: 'Optimization terminated successfully.'
          nfev: 48
           nit: 25
        status: 0
       success: True
             x: array([1.50419259, 0.10077416])

我们使用method='Nelder-Mead'因为这是一种无梯度优化方法，似乎需要解决这类问题。

Scipy.optimize.fsolve() 在我的复杂方程上

Scipy.optimize.fsolve() on my complex equation

python

scientific-computing

solver

scipy