【python】scipyの特殊関数を用いた方程式を解きたい

Question

実現したいこと

下式に適当なβを代入してφを求めたいのですがうまくいきません。

jn(x)：球ベッセル関数
β：0~1の任意の数値

参考文献によるとβに対するφの値は下図のようになるそうです。
ご教授いただけますと助かります。

試したこと

以下コードの通り、optimize.fsolveとoptimize.newtonより計算しようとしたのですが解が出ず止まってしまいます。

#optimize.fsolveの場合
import numpy as np
import scipy as sc
from scipy import optimize
from scipy.special import spherical_jn
import cmath as cm

def j(n, z):
    return spherical_jn(n, z)

def f(x):
    return j(0, x[0])*j(2, cm.sqrt(b)*x[0])+j(2, x[0])*j(0, cm.sqrt(b)*x[0])

b = 1
x= optimize.fsolve(f, 0)

#optimize.newtonの場合
import numpy as np
import scipy as sc
from scipy import optimize
from scipy.special import spherical_jn, spherical_yn
import cmath as cm

def j(n, z):
    return spherical_jn(n, z)

def f(x):
    return j(0, x)*j(2, cm.sqrt(b)*x)+j(2, x)*j(0, cm.sqrt(b)*x)

b = 1
x = optimize.newton(f, 0)

補足資料

参考文献のURLはこちらです。
http://xrm.phys.northwestern.edu/research/pdf_papers/1996/kuno_physrevb_1996.pdf

Answer 1

※optimize.fsolveを使うのは初めてなので、手探りな回答です。
optimize.fsolve(f, 0)をoptimize.fsolve(f, 3.9)に変更：fsolveは開始時の推定値によって結果が大きく異なってしまうようなので、試しに測定値を振ってみて3.9が最も希望の形に違い答えを返しました。
理論の裏付けがないので、いいやり方ではないですね。
推定値を理論的に出せないなら、違うアルゴリズムを使ったほうが良さそうな。

#optimize.fsolveの場合
import numpy as np
import scipy as sc
from scipy import optimize
from scipy.special import spherical_jn
import cmath as cm
import matplotlib.pyplot as plt

def j(n, z):
    return spherical_jn(n, z)

def f(x):
    return j(0, x[0])*j(2, cm.sqrt(b)*x[0])+j(2, x[0])*j(0, cm.sqrt(b)*x[0])

B = np.linspace(0,1,11)
X = []
for b in B:
    X.append(optimize.fsolve(f, 3.9))

fig = plt.figure()
ax = fig.add_subplot(111, ylim=(2,7))
ax.plot(B,X)
fig.savefig('test.png')

Comments

[tallman0202]

ご回答いただき有難うございます。
計算がうまくいかない場合はfsolveの開始時の推定値を変更してみるという手があるのですね。
実際にグラフ化して説明いただきすごく助かりました。
理論の裏付けはない点に関してはもう少し勉強してみます。