function [ n_data, k, fx ] = bessel_j1_zero_values ( n_data )

%*****************************************************************************80
%
%% bessel_j1_values() returns some zeros of the J1 Bessel function.
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    06 June 2025
%
%  Author:
%
%    John Burkardt
%
%  Input:
%
%    integer N_DATA.  The user sets N_DATA to 0 before the first call.  
%    Thereafter, it should simply be the value returned by the previous call.
%
%  Output:
%
%    integer N_DATA.  On each call, the routine increments N_DATA by 1, and
%    returns the corresponding data; when there is no more data, the
%    output value of N_DATA will be 0 again.
%
%    integer K, the index of the zero.
%
%    real FX, the value of the function.
%
  n_max = 30;

  fx_vec = [ ...
     3.831705970207512, ...
     7.015586669815619, ...
    10.17346813506272, ...
    13.32369193631422, ...
    16.47063005087763, ...
    19.61585851046824, ...
    22.76008438059277, ...
    25.90367208761838, ...
    29.04682853491686, ...
    32.18967991097440, ...
    35.33230755008387, ...
    38.47476623477162, ...
    41.61709421281445, ...
    44.75931899765282, ...
    47.90146088718545, ...
    51.04353518357151, ...
    54.18555364106132, ...
    57.32752543790101, ...
    60.46945784534749, ...
    63.61135669848123, ...
    66.75322673409849, ...
    69.89507183749577, ...
    73.03689522557383, ...
    76.17869958464146, ...
    79.32048717547630, ...
    82.46225991437356, ...
    85.60401943635023, ...
    88.74576714492631, ...
    91.88750425169499, ...
    95.02923180804470 ];

  k_vec = [ ...
     1, ...
     2, ...
     3, ...
     4, ...
     5, ...
     6, ...
     7, ...
     8, ...
     9, ...
    10, ...
    11, ...
    12, ...
    13, ...
    14, ...
    15, ...
    16, ...
    17, ...
    18, ...
    19, ...
    20, ...
    21, ...
    22, ...
    23, ...
    24, ...
    25, ...
    26, ...
    27, ...
    28, ...
    29, ...
    30  ];

  if ( n_data < 0 )
    n_data = 0;
  end

  n_data = n_data + 1;

  if ( n_max < n_data )
    n_data = 0;
    k = 0;
    fx = 0.0;
  else
    k = k_vec(n_data);
    fx = fx_vec(n_data);
  end

  return
end

