把求根公式换成chareqn2/chareqn3 跑的时候会显示计算溢出 搞不明白 有没有老哥救一下
底下放源码
Clear;
\[Lambda] = 0.633;
nSiO2 = Sqrt[
1 + x1*\[Lambda]^2/(\[Lambda]^2 - x2^2) +
x3*\[Lambda]^2/(\[Lambda]^2 - x4^2) +
x5*\[Lambda]^2/(\[Lambda]^2 - x6^2) /. {x1 -> 0.6961663,
x2 -> 0.0684043, x3 -> 0.4079426, x4 -> 0.1162414,
x5 -> 0.8974794, x6 -> 9.896161}];
nSi = Sqrt[
11.6858 + 0.939816/\[Lambda]^2 +
0.000993358/(\[Lambda]^2 - 1.22567)];
n1 = nSiO2;
n2 = 1;
k0 = 2*\[Pi]/\[Lambda]
c = 2.99979*10^14;
d0 = 0.12; (*0.12um=120nm*)
\[CapitalDelta]d = 0.001;
\[CapitalDelta] = (n1^2 - n2^2)/(2*n1^2); (*折射率差*)
Num = 6000;
b = Table[0, {n, 1, Num}];
diameter = Table[0, {n, 1, Num}];
\[Eta] = Table[0, {n, 1, Num}];
vg = Table[0, {n, 1, Num}];
i = 1;
beta0 = k0*(1.0 + 10^-14);
v = 1;
R = r/a;
d = d0;
Do[
a = d/2;
U = d/2*Sqrt[(k0^2*n1^2 - \[Beta]^2)];
W = d/2*Sqrt[(\[Beta]^2 - k0^2*n2^2)];
V = k0*a*Sqrt[(n1^2 - n2^2)];
chareqn1 = {(
D[BesselJ[v, u], u]/(u*BesselJ[v, u]) + D[BesselK[v, w], w]/(
w*BesselK[v, w]))*(D[BesselJ[v, u], u]/(
u*BesselJ[v,
u]) + (n2^2*D[BesselK[v, w], w])/(n1^2*w*
BesselK[v, w])) == ((v*\[Beta])/(k0*n1))^2 (V/(
u*w))^4 /. {u -> U, w -> W}};
chareqn2 = {BesselJ[1, u]/(u*BesselJ[0, u]) + BesselK[1, w]/(
w*BesselK[0, w]) == 0 /. {u -> U, w -> W}}; (*TEmodel*)
chareqn3 = {(n1^2*BesselJ[1, u])/(u*BesselJ[0, u]) + (
n2^2*BesselK[1, w])/(w*BesselK[0, w]) == 0 /. {u -> U,
w -> W}}; (*TMmodel*)
beta0 = \[Beta] /.
FindRoot[chareqn1, {\[Beta], beta0}, WorkingPrecision -> 44,
PrecisionGoal -> 43, AccuracyGoal -> 44];
b[[i]] = beta0;
diameter[[i]] = d;
d = d + \[CapitalDelta]d, {i, Num}]
dip\[Beta] = Table[{diameter[[i]], b[[i]]}, {i, 1, Num}];
ListPlot[dip\[Beta]]
