bug#58929: 29.0.50; Calc: finding roots utpn doesn't work as advertised

bug#58929: 29.0.50; Calc: finding roots utpn doesn't work as advertised

[Matt Armstrong]

Quick facts:
- all of these examples work with emacs -Q
- reproduces in Emacs 29, 27, 25

I think I have found that the calc package's "root" finding command does
not work as documented with calc's probability distribution functions.
Or, at least, not with "utpn", the probability distribution function for
the normal (gaussian) distribution.  I have found similar issues with
"utpt" (the Student's t-distribution), but will focus on "utpn" here.

The Emacs Calc package info page for Probability Distribution Functions
(https://www.gnu.org/software/emacs/manual/html_node/calc/Probability-Distribution-Functions.html)
has this to say about the normal distribution function:

> The ‘utpn(x,m,s)’ function uses a normal (Gaussian) distribution with
> mean ‘m’ and standard deviation ‘s’. It is the probability that such a
> normal-distributed random variable would exceed ‘x’.

I am interested in solving `utpn(x,m,s) = p` for `p` given a fixed `m`
and `s`.  The same info page has this to say about that:

> While Calc does not provide inverses of the probability distribution
> functions, the a R command can be used to solve for the inverse. Since
> the distribution functions are monotonic, a R is guaranteed to be able
> to find a solution given any initial guess.

SUPPOSITION: the above is true, especially this: "a R is guaranteed to be
able to find a solution given any initial guess."

The Root Finding info page
(https://www.gnu.org/software/emacs/manual/html_node/calc/Root-Finding.html)
has this to say about the a R command:

> The a R (calc-find-root) [root] command finds a numerical solution (or
> root) of an equation. (This command treats inequalities the same as
> equations. If the input is any other kind of formula, it is interpreted
> as an equation of the form ‘X = 0’.)
>
> The a R command requires an initial guess on the top of the stack, and a
> formula in the second-to-top position. It prompts for a solution
> variable, which must appear in the formula. All other variables that
> appear in the formula must have assigned values, i.e., when a value is
> assigned to the solution variable and the formula is evaluated with =, it
> should evaluate to a number. Any assigned value for the solution variable
> itself is ignored and unaffected by this command.
>
> When the command completes, the initial guess is replaced on the stack by
> a vector of two numbers: The value of the solution variable that solves
> the equation, and the difference between the lefthand and righthand sides
> of the equation at that value. Ordinarily, the second number will be zero
> or very nearly zero. (Note that Calc uses a slightly higher precision
> while finding the root, and thus the second number may be slightly
> different from the value you would compute from the equation yourself.)

WHAT I HAVE TRIED

I have tried entering the formula `utpn(x, 10, 1) = 0.5` and invoking `a
R` with an initial guess of `10`.  I then solve for `x`, expecting
`x = 10`.

Note that the correct answer of `x = 10` is "obviously" correct.  I am
passing `m = 10` for the "mean" arg to `utpn`, and by definition a
random value from the normal distribution has an `0.5` chance of being
larger than the mean.

To do this I type "C-x * *" to enter the *Calculator* buffer and then
type:

    ` utpn(x, 10, 1) = 0.5 RET
    10 RET                      <-- initial guess
    a R RET                     <-- invoke "root" function
    x RET                       <-- answer to prompt; solve for "x"

and I end up with this in the `*Calculator*` buffer:

    2:  utpn(x, 10, 1) = 0.5
    1:  [10., 0.]

This is correct.  The 'a R' command has returned a vector whose first
element is `x` and the second is an error delta (not important here).

Alternatively, evaluate this lisp to do the same thing:

    ELISP> (calc-eval '("root(utpn(x, 10, 1) = 0.5, x, 10)"))
    "[10., 0.]"

While correct, here I have *given* the answer to calc in the initial
guess (the root function's third arg).  The problem is that the function
doesn't do well when given other guesses.

```
` utpn(x,10,1)=0.5 RET
1.0 RET                  <---- a different guess
a R RET
x RET
```

results in this:

```
2:  utpn(x, 10, 1) = 0.5
1:  root(utpn(x, 10, 1) = 0.5, x, 1)
```

With this in the minibuffer:

> *Newton's method failed to converge: 4.8639203627839e17*

I'm finding that if the initial guess is close enough to the right answer,
in this case `10`, then `a R` works, otherwise not.  E.g. guessing `8.6`
fails but `8.7` works.

This behavior makes `a R` much less useful.  Also, documentation
suggests that any guess is expected to work (re-quoting from above with
my emphasis):

> Since the distribution functions are monotonic, 'a R' is guaranteed to
> be able to find a solution *given any initial guess*.

Have I found a bug here, or am I doing something wrong?

ADDENDUM

This shows the function failing to find a root for guesses not all that
far away from the actual answer:

ELISP> (cl-loop for x from 8.0 to 12.0 by 0.1 collect (list (format "%.1f" x) (calc-eval '("root(utpn(x,10,1)=0.5,x,$)") nil (prin1-to-string x))))
(("8.0" "root(utpn(x, 10, 1) = 0.5, x, 8.)")
 ("8.1" "root(utpn(x, 10, 1) = 0.5, x, 8.1)")
 ("8.2" "root(utpn(x, 10, 1) = 0.5, x, 8.2)")
 ("8.3" "root(utpn(x, 10, 1) = 0.5, x, 8.3)")
 ("8.4" "root(utpn(x, 10, 1) = 0.5, x, 8.4)")
 ("8.5" "root(utpn(x, 10, 1) = 0.5, x, 8.5)")
 ("8.6" "root(utpn(x, 10, 1) = 0.5, x, 8.6)")
 ("8.7" "[10., 0.]")
 ("8.8" "[10., 0.]")
 ("8.9" "[10., 0.]")
 ("9.0" "[10., 0.]")
 ("9.1" "[10., 0.]")
 ("9.2" "[10., 0.]")
 ("9.3" "[10., 0.]")
 ("9.4" "[10., 0.]")
 ("9.5" "[10., 0.]")
 ("9.6" "[10., 0.]")
 ("9.7" "[10., 0.]")
 ("9.8" "[10., 0.]")
 ("9.9" "[10., 0.]")
 ("10.0" "[10., 0.]")
 ("10.1" "[10., 0.]")
 ("10.2" "[10., 0.]")
 ("10.3" "[10., 0.]")
 ("10.4" "[10., 0.]")
 ("10.5" "[10., 0.]")
 ("10.6" "[10., 0.]")
 ("10.7" "[10., 0.]")
 ("10.8" "[10., 0.]")
 ("10.9" "[10., 0.]")
 ("11.0" "[10., 0.]")
 ("11.1" "[10., 0.]")
 ("11.2" "[10., 0.]")
 ("11.3" "[10., 0.]")
 ("11.4" "root(utpn(x, 10, 1) = 0.5, x, 11.4)")
 ("11.5" "root(utpn(x, 10, 1) = 0.5, x, 11.5)")
 ("11.6" "root(utpn(x, 10, 1) = 0.5, x, 11.6)")
 ("11.7" "root(utpn(x, 10, 1) = 0.5, x, 11.7)")
 ("11.8" "root(utpn(x, 10, 1) = 0.5, x, 11.8)")
 ("11.9" "root(utpn(x, 10, 1) = 0.5, x, 11.9)")
 ("12.0" "root(utpn(x, 10, 1) = 0.5, x, 12.)"))



In GNU Emacs 29.0.50 (build 3, x86_64-pc-linux-gnu, GTK+ Version
 3.24.34, cairo version 1.16.0) of 2022-10-30 built on naz
Repository revision: 73953b23aad5af80275838cd3b1c88e045e1856e
Repository branch: master
Windowing system distributor 'The X.Org Foundation', version 11.0.12201003
System Description: Debian GNU/Linux bookworm/sid

Configured using:
 'configure 'CFLAGS=-O2 -g3''

Configured features:
ACL CAIRO DBUS FREETYPE GIF GLIB GMP GNUTLS GPM GSETTINGS HARFBUZZ JPEG
JSON LCMS2 LIBOTF LIBSELINUX LIBSYSTEMD LIBXML2 M17N_FLT MODULES NOTIFY
INOTIFY PDUMPER PNG RSVG SECCOMP SOUND SQLITE3 THREADS TIFF
TOOLKIT_SCROLL_BARS X11 XDBE XIM XINPUT2 XPM GTK3 ZLIB

Important settings:
  value of $LANG: en_US.UTF-8
  value of $XMODIFIERS: @im=ibus
  locale-coding-system: utf-8-unix

Major mode: Text

Minor modes in effect:
  msb-mode: t
  display-time-mode: t
  flyspell-mode: t
  global-tab-line-mode: t
  tab-line-mode: t
  envrc-global-mode: t
  envrc-mode: t
  shell-dirtrack-mode: t
  auto-insert-mode: t
  keyfreq-autosave-mode: t
  keyfreq-mode: t
  savehist-mode: t
  icomplete-vertical-mode: t
  icomplete-mode: t
  editorconfig-mode: t
  which-key-mode: t
  electric-pair-mode: t
  override-global-mode: t
  tooltip-mode: t
  global-eldoc-mode: t
  show-paren-mode: t
  electric-indent-mode: t
  mouse-wheel-mode: t
  tab-bar-mode: t
  menu-bar-mode: t
  file-name-shadow-mode: t
  global-font-lock-mode: t
  font-lock-mode: t
  blink-cursor-mode: t
  column-number-mode: t
  line-number-mode: t
  auto-fill-function: do-auto-fill
  transient-mark-mode: t
  auto-composition-mode: t
  auto-encryption-mode: t
  auto-compression-mode: t
  temp-buffer-resize-mode: t
  abbrev-mode: t

Load-path shadows:
~/env/elisp/ol-notmuch hides /home/matt/.config/emacs/elpa/ol-notmuch-20220428.1337/ol-notmuch
/home/matt/.config/emacs/elpa/transient-20220918.2101/transient hides /home/matt/git/e/daily-driver/lisp/transient
/home/matt/.config/emacs/elpa/eglot-20221005.1955/eglot hides /home/matt/git/e/daily-driver/lisp/progmodes/eglot

Features:
(shadow emacsbug ielm calc-yank shortdoc calc-undo edebug dired-aux
calcalg3 imenu cl-print calc-stuff calc-funcs dabbrev calcalg2 calc-poly
calc-vec calc-mtx calc-lang calc-cplx calccomp calc-menu calc-frac
calc-comb calc-aent calc-math calc-arith calc-alg calc-bin calc-tests
calc-misc calc-prog calc-forms calc-units calc-ext calc calc-loaddefs
rect calc-macs ert pp ewoc debug backtrace mule-util magit-base crm
misearch multi-isearch copyright help-fns radix-tree pcase smerge-mode
diff qp sort smiley gnus-cite mail-extr textsec uni-scripts idna-mapping
ucs-normalize uni-confusable textsec-check gnus-async gnus-bcklg
gnus-agent gnus-srvr gnus-score score-mode nnvirtual nntp gnus-ml
gnus-msg nndoc gnus-cache gnus-dup mm-archive network-stream url-cache
display-line-numbers debbugs-gnu add-log debbugs-compat debbugs
soap-client url-http url-auth url-gw nsm rng-xsd rng-dt rng-util
xsd-regexp debbugs-browse vc-git diff-mode vc-dispatcher bug-reference
editorconfig-core editorconfig-core-handle editorconfig-fnmatch protbuf
msb time flyspell ispell org-element avl-tree generator ol-w3m ol-rmail
ol-mhe ol-irc ol-info org-habit org-agenda org-refile ol-gnus nnselect
gnus-art mm-uu mml2015 mm-view mml-smime smime gnutls dig gnus-sum
gnus-group gnus-undo gnus-start gnus-dbus dbus gnus-cloud nnimap nnmail
mail-source utf7 nnoo parse-time gnus-spec gnus-int gnus-range message
sendmail yank-media rfc822 mml mml-sec epa derived epg rfc6068
epg-config mm-decode mm-bodies mm-encode mail-parse rfc2231 rfc2047
rfc2045 ietf-drums mailabbrev gmm-utils mailheader gnus-win ol-eww eww
xdg url-queue shr pixel-fill kinsoku url-file svg xml dom puny mm-url
gnus nnheader gnus-util text-property-search mail-utils range wid-edit
mm-util mail-prsvr ol-doi org-link-doi ol-docview doc-view filenotify
jka-compr image-mode exif dired dired-loaddefs ol-bibtex ol-bbdb
tab-line server envrc inheritenv web-mode disp-table nix-mode ffap
thingatpt smie nix-repl nix-shell nix-store magit-section dash compat-27
compat-26 nix-instantiate nix-shebang nix-format nix dirtrack ob-shell
shell ob-ruby ob-python python compat compat-macs ob-dot org-protocol
org ob ob-tangle ob-ref ob-lob ob-table ob-exp org-macro org-footnote
org-src ob-comint org-pcomplete pcomplete comint ansi-osc ansi-color
ring org-list org-faces org-entities noutline outline org-version
ob-emacs-lisp ob-core ob-eval org-table oc-basic bibtex iso8601
time-date org-keys oc org-loaddefs find-func cal-menu calendar
cal-loaddefs finder-inf ol-notmuch ol rx org-compat org-macs format-spec
skeleton autoinsert advice keyfreq project edmacro kmacro warnings icons
savehist icomplete editorconfig which-key package browse-url url
url-proxy url-privacy url-expand url-methods url-history url-cookie
generate-lisp-file url-domsuf url-util mailcap url-handlers url-parse
auth-source eieio eieio-core password-cache json subr-x map byte-opt
url-vars cl-extra help-mode cl-macs gv cl-seq elec-pair use-package
use-package-ensure use-package-delight use-package-diminish
use-package-bind-key bind-key easy-mmode use-package-core cl-loaddefs
cl-lib bytecomp byte-compile info bazel-autoloads
clang-format+-autoloads clang-format-autoloads cmake-mode-autoloads
d-mode-autoloads debbugs-autoloads editorconfig-autoloads
eglot-autoloads elpy-autoloads company-autoloads envrc-autoloads
exec-path-from-shell-autoloads flymake-ruby-autoloads
flymake-easy-autoloads flymake-yamllint-autoloads go-mode-autoloads
google-c-style-autoloads graphviz-dot-mode-autoloads
highlight-indentation-autoloads inheritenv-autoloads magit-autoloads
git-commit-autoloads markdown-mode-autoloads meson-mode-autoloads
nix-mode-autoloads magit-section-autoloads dash-autoloads
nixpkgs-fmt-autoloads ol-notmuch-autoloads notmuch-autoloads
orderless-autoloads org-drill-autoloads ox-hugo-autoloads
persist-autoloads pylint-autoloads pyvenv-autoloads s-autoloads
shfmt-autoloads reformatter-autoloads tomelr-autoloads
transient-autoloads use-package-autoloads bind-key-autoloads
vertico-autoloads web-mode-autoloads which-key-autoloads
with-editor-autoloads compat-autoloads yaml-mode-autoloads
yasnippet-autoloads rmc iso-transl tooltip cconv eldoc paren electric
uniquify ediff-hook vc-hooks lisp-float-type elisp-mode mwheel
term/x-win x-win term/common-win x-dnd tool-bar dnd fontset image
regexp-opt fringe tabulated-list replace newcomment text-mode lisp-mode
prog-mode register page tab-bar menu-bar rfn-eshadow isearch easymenu
timer select scroll-bar mouse jit-lock font-lock syntax font-core
term/tty-colors frame minibuffer nadvice seq simple cl-generic
indonesian philippine cham georgian utf-8-lang misc-lang vietnamese
tibetan thai tai-viet lao korean japanese eucjp-ms cp51932 hebrew greek
romanian slovak czech european ethiopic indian cyrillic chinese
composite emoji-zwj charscript charprop case-table epa-hook
jka-cmpr-hook help abbrev obarray oclosure cl-preloaded button loaddefs
theme-loaddefs faces cus-face macroexp files window text-properties
overlay sha1 md5 base64 format env code-pages mule custom widget keymap
hashtable-print-readable backquote threads dbusbind inotify lcms2
dynamic-setting system-font-setting font-render-setting cairo
move-toolbar gtk x-toolkit xinput2 x multi-tty make-network-process
emacs)

Memory information:
((conses 16 977488 266858)
 (symbols 48 43332 287)
 (strings 32 225511 12569)
 (string-bytes 1 5996788)
 (vectors 16 103426)
 (vector-slots 8 1994918 188498)
 (floats 8 517 940)
 (intervals 56 9190 882)
 (buffers 984 38))

bug#58929: 29.0.50; Calc: finding roots utpn doesn't work as advertised

[Mattias Engdegård]

Unfortunately straight Newton won't converge for this function (erf, essentially) if you stray too far from the origin (ie, the mean) so I suggest you always use that as a starting guess no matter what point you are trying to find the inverse for.

Another way is to use an interval as starting guess instead of a single point; that seems to make Calc switch tactics and persevere.

One thing that can happen is running out of stack, because many of the implementations like math-newton-root are recursive. This should be easy to remedy as they all seem to be tail calls.

Obviously the claim to guarantee a solution for "any initial guess" is fanciful.

bug#58929: 29.0.50; Calc: finding roots utpn doesn't work as advertised

[Matt Armstrong]

[-- Attachment #1: Type: text/plain, Size: 1130 bytes --]

Mattias Engdegård <mattiase@acm.org> writes:

> Obviously the claim to guarantee a solution for "any initial guess" is
> fanciful.

Yes.  See attached patch for what I wish the docs had said all along
(would have saved me a lot of time).


> Unfortunately straight Newton won't converge for this function (erf,
> essentially) if you stray too far from the origin (ie, the mean) so I
> suggest you always use that as a starting guess no matter what point
> you are trying to find the inverse for.
>
> Another way is to use an interval as starting guess instead of a
> single point; that seems to make Calc switch tactics and persevere.

I wonder if the implementation can be improved.
https://www.gnu.org/software/emacs/manual/html_node/calc/Root-Finding.html
suggests that there are some cases where calc will switch over to the
bisection method.  Perhaps the heuristics can be tweaked?


> One thing that can happen is running out of stack, because many of the
> implementations like math-newton-root are recursive. This should be
> easy to remedy as they all seem to be tail calls.

Good to know.


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #2: 0002-In-calc.texi-provide-guidance-for-taking-the-inverse.patch --]
[-- Type: text/x-diff, Size: 2134 bytes --]

From ba68ca5671f9772994bc5c0b2b0b9216a750b692 Mon Sep 17 00:00:00 2001
From: Matt Armstrong <matt@rfc20.org>
Date: Tue, 1 Nov 2022 14:51:36 -0700
Subject: [PATCH 2/2] In calc.texi provide guidance for taking the inverse of
 PDFs.

* doc/misc/calc.texi (Probability Distribution Functions): Stop
claiming that the 'a R' function will "always work" for any guess.
Explain the problem with the PDF functions being very flat.  Provide
recommended guess values for each function.  bug#58929
---
 doc/misc/calc.texi | 36 +++++++++++++++++++++++++++++++++---
 1 file changed, 33 insertions(+), 3 deletions(-)

diff --git a/doc/misc/calc.texi b/doc/misc/calc.texi
index 89a340e7343..e6b416a447a 100644
--- a/doc/misc/calc.texi
+++ b/doc/misc/calc.texi
@@ -19456,9 +19456,39 @@ Probability Distribution Functions
 
 While Calc does not provide inverses of the probability distribution
 functions, the @kbd{a R} command can be used to solve for the inverse.
-Since the distribution functions are monotonic, @kbd{a R} is guaranteed
-to be able to find a solution given any initial guess.
-@xref{Numerical Solutions}.
+@xref{Root Finding}.
+
+The distribution functions, while monotonic, are also be very flat,
+which can cause problems for the Newton's root-finding algorithm used
+by @kbd{a R}.  The problem can be avoided by providing a good starting
+guess.  Prefer guesses that occur where the slope of the plotted
+function is not flat.  Alternatively, use an interval range for the
+guess.  The table below provides a recommendation for each of the
+probability distribution functions.
+
+@multitable @columnfractions .25 .25 .25
+@headitem Function
+@tab When solving for
+@tab Recommended initial guess
+@item @samp{utpb(x,n,p)}
+@tab @expr{x}
+@tab @code{0.5}
+@item @samp{utpc(x,c)}
+@tab @expr{x}
+@tab @code{1.0}
+@item @samp{utpf(F,v1,v2)}
+@tab @expr{F}
+@tab @code{1.0}
+@item @samp{utpn(x,m,s)}
+@tab @expr{x}
+@tab @expr{m}
+@item @samp{utpp(n,x)}
+@tab @expr{n}
+@tab @expr{x}
+@item @samp{utpt(t,v)}
+@tab @expr{t}
+@tab @code{[0..1e4]}
+@end multitable
 
 @node Matrix Functions
 @chapter Vector/Matrix Functions
-- 
2.35.1

bug#58929: 29.0.50; Calc: finding roots utpn doesn't work as advertised

[Mattias Engdegård]

1 nov. 2022 kl. 23.00 skrev Matt Armstrong <matt@rfc20.org>:

> See attached patch for what I wish the docs had said all along
> (would have saved me a lot of time).

Thank you, I have not worked through your table to confirm all the details but it looks fine on the whole.

One detail: it's not just that the functions are "flat" (ie, all derivatives rapidly vanishing as we stray from the mean) but at least for the normal distribution, Newton actively diverges with a bad starting guess even if computed with infinite precision (draw the graph).

> I wonder if the implementation can be improved.
> https://www.gnu.org/software/emacs/manual/html_node/calc/Root-Finding.html
> suggests that there are some cases where calc will switch over to the
> bisection method.  Perhaps the heuristics can be tweaked?

From what I can tell, bisection is used if there is an interval bracketing the root and either no derivative could be computed or Newton fails to converge. Clearly Newton is preferable when it converges since it's much faster than bisection.

We could employ special tactics by pattern-matching the function expression (perhaps according to your table); not sure if it would be worth the trouble.

Of course an actual numerical analyst would know immediately what to do. Is there one in the audience?

bug#58929: 29.0.50; Calc: finding roots utpn doesn't work as advertised

[Eli Zaretskii]

> Cc: 58929@debbugs.gnu.org
> From: Mattias Engdegård <mattiase@acm.org>
> Date: Thu, 3 Nov 2022 09:53:02 +0100
> 
> Of course an actual numerical analyst would know immediately what to do. Is there one in the audience?

I have good experience from using this technique:

  https://www.researchgate.net/publication/226830173_On_the_Structure_of_Zero_Finders

It employs a hybrid method that doesn't need derivatives (but is
almost as fast as Newton).  However, its main advantage (which is a
huge one in some real-life situations) is that the structure of the
algorithm allows the caller to make arbitrary modifications to the
approximations that the algorithm produces, and in general have
complete control on the iterations' process.  In particular, one can
easily deal with situations where the root approximation goes out of
the expected range, or lands in the area where the function doesn't
behave well enough.

bug#58929: 29.0.50; Calc: finding roots utpn doesn't work as advertised

[Matt Armstrong]

Eli Zaretskii <eliz@gnu.org> writes:

>> Cc: 58929@debbugs.gnu.org
>> From: Mattias Engdegård <mattiase@acm.org>
>> Date: Thu, 3 Nov 2022 09:53:02 +0100
>> 
>> Of course an actual numerical analyst would know immediately what to do. Is there one in the audience?
>
> I have good experience from using this technique:
>
>   https://www.researchgate.net/publication/226830173_On_the_Structure_of_Zero_Finders
>
> It employs a hybrid method that doesn't need derivatives (but is
> almost as fast as Newton).  However, its main advantage (which is a
> huge one in some real-life situations) is that the structure of the
> algorithm allows the caller to make arbitrary modifications to the
> approximations that the algorithm produces, and in general have
> complete control on the iterations' process.  In particular, one can
> easily deal with situations where the root approximation goes out of
> the expected range, or lands in the area where the function doesn't
> behave well enough.

Interesting.  I like how that algorithm attempts to "do its best" when
given a single guess point, then adapts to bisection if it happens to
find an opportunity and need.  Calc's manual suggests that its "root"
function attempts to do the same, opening a window to use this other
approach.

I'm not sure how to expose the root finding "loop" in calc's UI, but
exposing it calc's lisp level interface would be a natural thing to do.

https://en.wikipedia.org/wiki/ITP_method also looks like an interesting
twist on the bisection/secant method.  I have no personal experience
with it, and on Wikipedia I always wonder if pages like this are
reflective of a true utility or merely the ambitions of the paper's
authors!  This function needs function values that are opposite in sign
before it is useful.  It might be possible to incorporate it into the
approach Eli linked to, but that method claims to be almost as fast as
Newton's method in most cases so it may not need the improvement.

I arrived here because I wanted to use calc to compute 95th and 99th
percentile confidence intervals of some benchmarking results using the
Student's t-distribution.  I have figured this out now, but it has been
more difficult than I anticipated, and the manual's over-promise with
respect to "find the root, any guess works" sent me on a long tangent
trying to figure out what I had done wrong.  That said, my background in
this area is weak, so I don't think I am the right person to improve
Emacs' root finding implementation.

The patch I sent previously, with my (hopefully) useful suggestions in
the user manual, is probably the end of my useful contribution here.

bug#58929: 29.0.50; Calc: finding roots utpn doesn't work as advertised

[Mattias Engdegård]

5 nov. 2022 kl. 20.29 skrev Matt Armstrong <matt@rfc20.org>:

> I arrived here because I wanted to use calc to compute 95th and 99th
> percentile confidence intervals of some benchmarking results using the
> Student's t-distribution.

Quantile computations are common and useful enough that they may merit convenience functions in Calc in the future.

> The patch I sent previously, with my (hopefully) useful suggestions in
> the user manual, is probably the end of my useful contribution here.

Your manual changes are fine as far as I'm concerned, but again I haven't worked through them in detail -- maybe you could add tests to confirm the convergence for your recommended starting guesses? That would both help validating your patch and protect against future regressions.