Re: libcerf function voigt_hwhm not suitable for external use

[Ethan A M. <me...@uw...>]

Joachim:

It is good to hear that you are looking into this.
It was pretty easy to trigger the non-convergence error from gnuplot;
did I not send test cases?  
Testing these functions in gnuplot is really easy because you can
sample any desired domain in the complex plane (or in this case [sigma,gamma])
with a single command and plot the result.  Glitches or non-convergence or
outright errors and crashes show up immediately.
Sample plots using functions from libcerf are shown here:
    http://gnuplot.sourceforge.net/demo_5.5/cerf.html

Some time after I sent you that original note, I went ahead and re-implemented
the function directly in gnuplot.  I attach the code to this mail.
Note that it uses a different optimization algorithm (regula falsi)
that unlike the code in libcerf at that time is guaranteed to converge.

This code is from the file libcerf.c in the gnuplot project git repository.
You are welcome to adapt it for inclusion in libcerf, or I would be
willing to do so myself and contribute it.

	best regards,

		Ethan


On Thursday, 10 June 2021 07:18:03 PDT Joachim Wuttke wrote:
> Dear Ethan,
> 
> thank you very much for the kind mail you wrote me more than a year ago.
> 
> Today, finally, I looked into the issue.  I wrote a new test that computes
> voigt_hwhm for 100'000 different ratios gamma/sigma.  As expected, the Newton
> algorithm did not fail a single time.  The deviation of the final result from
> the initial estimate is always smaller then 1E-2.  The cases I had guarded
> with messages to stderr and an exit statement will never happen.  Therefore
> I converted them to assert statements.  These improvements are available
> in the freshly released v1.15.
> 
> Sorry for not answering earlier.
> 
> Kind greetings, Joachim
> 
> 
> 
> On 25/02/2020 20.03, Ethan A Merritt wrote:
> > In response to a gnuplot feature request, I added a wrapper in
> > gnuplot to find the width of the Voight profile by calling libcerf function
> > voigt_hwhm().
> > 
> > However on further evaluation this has turned out to be problematic.
> > The current libcerf source for this function deals badly with out-of-range
> > input.  It prints to stderr (problematic by itself) and then either
> > - returns an incorrect value
> > - calls exit(-1)
> > 
> > A call to exit() from a library routine is almost never a good thing to 
> do,
> > since it prevents error-handling and possible recovery by the calling
> > program.  Limiting error reporting to a message on stderr also prevents
> > detection and recovery by the calling program.
> > 
> > I suggest that it would be preferable to return NaN whenever the
> > input is out of range or the function does not converge.
> > 
> > Alternatively libcerf could provide an error-check routine or status variable,
> > but then you would run into issues of thread-safety and so on.
> > 
> > I am hoping that the libcerf function can be improved in a subsequent
> > version, but for now I will consider implementing equivalent code
> > directly in gnuplot.
> > 
> > 	Ethan
> > 
> > 
> 
> 


-- 
Ethan A Merritt
Biomolecular Structure Center,  K-428 Health Sciences Bldg
MS 357742,   University of Washington, Seattle 98195-7742