Fzero function with Bessel functions

5 次查看（过去 30 天）
Diana Dawoud
Diana Dawoud 2021-2-5
回答： Walter Roberson 2021-2-5
I'm trying to excute the following function, bu the problem is that I can keep getting the error "functions values t interval endpoints mut be finite and real", how can I resolve this problem
% File name: HE11.m
% Function determines and plots universal relation for mode HE_11
clear all
format long
l=1;
N_max = 190;
for n = 1:N_max
%VV(n) = 0.5 + n*0.01;
VV(n) = 3 + n*0.01;
V = VV(n);
b_c(n) = fzero(@(b) func_HE11(b,V,l),[0.0000001 30]);
end
%
hold on
h1 = plot(VV,b_c); % plots b =b(V)
%
axis([0 2.5 0 1.5]);
% Redefine figure properties
xlabel('V','FontSize',22);
grid on
box on
set(h1,'LineWidth',1.5); % new thickness of plotting lines
set(gca,'FontSize',22); % new size of tick marks on both axes
function result = func_HE11(b,V,l)
% Function name: func_HE11.m
% Function which defines basic relation
if l==0
u= V*sqrt(1-b);
J1 = besselj(1,u);
J0 = besselj(0,u);
%
w = V*sqrt(b);
K1 = besselk(1,w);
K0 = besselk(0,w);
lhs = J1./J0;
rhs = (w./u).*(K1./K0);
result = lhs - rhs;
else
u= V*sqrt(1-b);
J1 = besselj(l-1,u);
J0 = besselj(l,u);
%
w = V*sqrt(b);
K1 = besselk(l-1,w);
K0 = besselk(l,w);
lhs = J1./J0;
rhs = (w./u).*(K1./K0);
result = lhs - rhs;
end
end

请先登录，再回答此问题。

回答（1 个）

Walter Roberson
Walter Roberson 2021-2-5
b_c(n) = fzero(@(b) func_HE11(b,V,l),[0.0000001 30]);
fzero() needs to evaluate at both endpoints, and the two endpoints have to give different sign of value, so that fzero can do whatever (bisection, etc..) So the function will be evaluated with b = 30 and l = 1
u= V*sqrt(1-b);
b > 1, so 1-b is negative, so sqrt(1-b) is complex, so u is complex valued.
J1 = besselj(l-1,u);
l is 1, so l-1 is 0, so that is a bessel0 computation of a complex number -- confusingly then stored in a variable named J1 which would tend to suggest bessel1 . Anyhow, the computation comes up real-valued for b = 30.
J0 = besselj(l,u);
l is 1, so this is a bessel1 computation of a complex number -- confusingly then stored in a variable named J0 whikch would tend to suggest bessel0. Anyhow, the computation comex up complex-valued for b = 30.
w = V*sqrt(b);
K1 = besselk(l-1,w);
K0 = besselk(l,w);
real, real, and real.
lhs = J1./J0;
Ratio of a real and a complex is complex, so lhs is complex.
rhs = (w./u).*(K1./K0);
w is real but u is complex, so their ratio is complex. K1 and K0 are real so their ratio is real. So the overall expression is complex.
result = lhs - rhs;
Complex minus complex. Do the imaginary components just happen to cancel out, giving something that is real-valued as the result? NO. And the imbalance is too much compared to the real component to believe that it would have balanced out if not for round-off error.
The easiest fix would appear to be to plot to no more than about 0.9999

类别

在 Help Center 和 File Exchange 中查找有关 Bessel functions 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by