Bessel函數的零點計算 MATLAB

由於MATLAB自己沒有附帶貝塞爾函數零點，因此使用起來很不方便，特別是在繪制仿真場量時。

下面給出0～9階的貝塞爾函數零點的計算公式，其中理論上計算零點個數N在50以內時較為精確：

function rootBessel = besal_pi0(n,N)    

% 求解n階貝塞爾函數的零點(0-9)
% n為貝塞爾函數階數 
% N為要求的零點數量

j = zeros(n+1, N);    % 貝塞爾函數的根
incr = 4.0;
for v = 0 : n
   h = v + 1.9*v^(1/3)+1;
   if (v == 0)             % 0階貝塞爾函數的第一個零點
       j(v+1,1) = fzero(@(x)besselj(v,x),2);
   else                    % 1階及以上階貝塞爾函數的第一個零點
       j(v+1,1) = fzero(@(x)besselj(v,x),h);
   end
   for s = 2 : N           % 貝塞爾函數的第2個及后面的零點
       j(v+1,s) = fzero(@(x)besselj(v,x),j(v+1,s-1)+incr);
   end    
end

rootBessel = j;