The formation of a rational S matrix using phase shifts from elastic scattering (Q1110988)

In scattering theory as currently applied to nuclear physics, bound states, antibound states and resonances for given angular momentum appear in the S matrix as simple poles. It is therefore useful to approximate \(D(k)=\exp [2ikR]S(k)\), where R is the potential range (possibly 0), as a rational function of the momentum k. For the s-wave, \[ D(k)=\frac{G(k^ 2)+ik}{G(k^ 2)-ik}\quad where\quad G(k^ 2)=k\cot (\delta_ 0(k)+kR). \] The meromorphic continuation of G(x) can be characterized statistically by minimizing \[ \chi^ 2=\sum^{N}_{j=1}[G_ e(x_ j)-G_ t(x_ j,{\underline \alpha})]^ 2/\sigma^ 2(x_ j). \] Here \(\{G_ e(x_ j)\}\) are the data, \(\{\sigma (x_ j)\}\) their associated standard errors and \(\{G_ t(x,\alpha)\}\) an approximation depending on a parameter vector \({\underline \alpha}\). Assume \[ G_ t(x,{\underline \alpha})=P_ L(x)/Q_ M(x)=[L/M] \] where \(P_ L(x)=a_ 0+a_ 1x+...+a_ Lx^ L\), \(Q_ M(x)=1+b_ 1x+...+b_ Mx^ M\) and \({\underline \alpha}=(a_ 0,...,a_ L,b_ 1,...,b_ M).\) In class I problems the data are smooth and the errors uniformly small - often \(G_ 2\) is an analytical expression. In class II problems the data contain significant random variability. This method of approximation is a statistical generalization of Padé approximants (PAs), and defective PAs, in which P and Q have common or nearly common zeros, must be avoided (with care!). No method examined is superior for both classes of problem. Five methods are compared: two old and three new. The Newton-Raphson method (NR) requires the time-consuming evaluation of the Hessian matrix of \(\chi^ 2\). The Levenberg-Marquardt method (LM) is tailored to minimizing a sum of squares and only requires the first derivative of \(\chi^ 2\) making each iteration faster, and using finite differences makes it even faster. However, the convergence of LM is slow, and both it and NR may fail to converge. A method due to \textit{K. Miller} [SIAM J. Appl. Math. 18, 346-363 (1970; Zbl 0211.193)] minimizes \[ \chi^ 2_ 0=\sum^{N}_{j=1}[G_ e(x_ j)Q_ M(x_ j)-P_ L(x_ j)]^ 2/\sigma^ 2(x_ j) \] by solving the linear equations \[ (*)\quad \partial \chi^ 2_ 0/\partial a_ j=0,\quad 0\leq j\leq L;\quad \partial \chi^ 2_ 0/\partial b_ k=0,\quad 0\leq k\leq M. \] The KKH algorithm [\textit{V. M. Krasnopolsky}, \textit{V. I. Kukulin} and \textit{J. Horáček}, Czech. J. Phys. B35, 805 (1985)] starts from the solution of (*) defining \(P_ L^{(0)}\) and \(Q_ M^{(0)}\), defines \(\sigma^{(1)}(x)= \sigma(x)Q_ M^{(0)}(x)\) and iteratively re-solves (*). The HY method of [\textit{K. Hartt} and \textit{P. V. A. Yidana}, Phys. Rev. C31, 1105 (1985)] iteratively minimizes \(\chi^ 2\) with respect to \(\{a_ j\}\), and then minimizes \[ \chi_{H^ 2}=\sum^{N}_{j=1}[1/G_ e(x_ j)-Q_ M(x_ j)/P_ L(x_ j)]^ 2/\sigma_{H^ 2}(x_ j),\quad \sigma_ H(x)=\sigma (x)/| G_ e^ 2(x)-\sigma^ 2(x)|, \] with respect to \(\{b_ j\}\). Sometimes the Miller method gives good starting values. Both KKH and HY tend to converge quickly and then diverge. The computations reported used quadruple precision (32 digit) arithmetic because of the ill-conditioning of the root-finding problem. To reconstruct a linear rational [1/1] function with no noise requires 3 data points, and Miller's method is best. For a more realistic [3/2] fit of phase-shift data on 40 points KKH is best. With experimental error included HY is best, and HY can be easily adapted to include linear constraints. NR and LM are generally inferior, but ideally all 5 methods should be available.