I am having some doubts about publishing a code that uses a routine based on the globally convergent Newton method from Numerical Recipes. I do have the book and the routines. However, it wasn’t possible to incorporate their routine in the original form in my code. I had to heavily adapt it. Instead of a general routine, now it is a problem-specific one written in my coding style, using the globals of my main code, some features needed for the parallelization, a routine to compute Jacobian analytically, etc. Still there is some percentage of their original code spread over my version and it is certainly the same algorithm.
I know that the NR authors are particularly strict about their rights and that they do not allow sharing their source code. However, I have trouble to understand if my version of their code still falls under the license: http://numerical.recipes/licenses/redistribute.html
It’s a bit awkward situation. On one hand, it is, of course, perfectly fair to respect their rights. On the other hand, a large part of their own routine is very basic matrix algebra (LU decomposition) that is part of many numerical math textbooks and should not be monopolized. The part that I’ve learned from their book/routine is how to use the linear search in the quasi-Newton method. As they do not provide any references in that section (9.7), I am not sure if this is their original addition to the method or not.
It is a relatively simple algorithm and I could program it from the scratch (I guess that should be perfectly legal?), but even in that case my implementation would likely be similar or at least influenced by their code.
Does anyone have experience with similar situation?
I think you are safe if you learn the algorithm from Numerical Recipes and then program the algorithm from scratch. For solving sets of nonlinear equations, can you use open-source codes, for example the version of minpack that @certik et al. are maintaining, or nonlin? There are also open-source Fortran codes for LU decomposition.
Thanks for your reply. As I said, if I do it from the scratch, I’d likely be biased towards their solution. It feels somehow awkward to look for a different way of implementation only to avoid legal issues. In this part (LU decomposition), their implementation is not very original.
In the meantime I found in the library the book that they cite at the end of S.9, by Dennis and Schnabel. The linear search algorithm in NR comes mostly from there (I still need to compare the derivation, but it seems to be the same thing). So, the question is really how much one can claim rights on implementation of well known methods described in literature?
Thanks for mentioning minpack. I’ll check carefully what is in there.
NR has copyrighted its code but has not patented the underlying algorithms, which they do not claim to have discovered and which you are free to use implementations of.
tensolve.zip This contains both TOMS 739 (UNCMIN) for unconstrained minimization, and TOMS 768 (TENSOLVE) for the solution of sets of non-linear equations.
I just tried compiling the code with gfortran blas_prt.f90 uncmin.f90 tensolve.f90 problem.f90. First I got an error
tensolve.f90:7:27:
3 | USE unconstrained_min
| 2
......
7 | INTEGER, PARAMETER :: dp = SELECTED_REAL_KIND(14, 60)
| 1
Error: Symbol 'dp' at (1) conflicts with symbol from module 'unconstrained_min', use-associated at (2)
problem.f90:2:5:
which is resolved by commenting out line 7. Doing this, I still get an error
The algorithm for what I need is actually rather simple, so I guess I’ll go for the version written from the scratch (or, better, to rewrite from the scratch the remaining parts that are not rewritten yet).
As others said, best is not to use or look at any code in numerical recipes at all to avoid these copyright issues. I believe it is ok to learn the algorithm (from any source?), and then implement it yourself.
The following in particular shall not be considered inventions within the meaning of the first paragraph of this Article
(a) discoveries as well as scientific theories and mathematical methods ;
A numerical method can not be patented. The text of the source code implementing it can be copyrighted (licensed). But if you implement the same method from scratch, the source code will not be identical.
The nice thing with NR is that there is detailed mathematical explanation accompanying the code. Some alternative libraries have very good code, but there is no detailed description of what the code is doing. So, perhaps the effort of creating alternative numerical libraries should also include writing alternative “book” (or even better some wiki corresponding to the routines).
The late Dr. Alan Miller used a modified version of Mike Metcalf’s convert.f90 and edited the output to get it running with the free (at that time) ELF90 compiler from Lahey. He converted many of the ACM TOMS codes this way, and most of the converted codes work fine with modern Fortran compilers. However, there are many instances where a modern compiler will not accept Miller’s converted codes, and Tensolve is one of them. Do not use compiler options that will check whether subprogram arguments have correct INTENT specified; Miller specifies INTENT(OUT) in many cases where he should have not specified any INTENT or should have specified INTENT(IN OUT).
Two simple fixes for the converted TENSOLVE are as follows.
In UNCMIN.F90, add the attribute PRIVATE to the declaration of dp on line 17.
In TENSOLVE.F90, in subroutine TSNESI, comment out or delete the lines with the interface to JAC. In the same subroutine, in the call to TSNESV, replace the 22nd argument, JAC, by TSDUMJ.
An especially egregious example from numerical recipes is the MEDFIT code for L1 fitting. It is built on an inherently incorrect mathematical premise: that the minimum of a function with discontinuous derivatives occurs where the derivative is zero. In the case of L1 fitting, the minimum occurs precisely where the derivative does not exist (or it’s along a horizontal line). Starting from that flawed mathematical premise, they built a flawed algorithm. Then they coded the algorithm carelessly or incorrectly. Here’s what Charles Lawson wrote in 1991, in Why not use NR:
“One data set which causes looping is [x = 1, 2, 3; y = 1, 1, 1]. Another which causes looping in a different part of the code is [x = 2, 3, 4; y = 1, 3, 2]. A data set on which the code terminates, but with a significantly wrong result is [x = 3, 4, 5, 6, 7; y = 1, 3, 2, 4, 3]. Because of the faulty theoretical foundation, there is no reason to believe any particular result obtained by this code is correct, although by chance it will sometimes get a correct result…”
I had several cycles of correspondence with the authors. They never exhibited any evidence of understanding the problem.
Around 2016 I attended a university course on numerical methods by Simon Širca, a Slovenian physicist and author of the Springer textbook – Computational Methods in Physics (2018). Concerning Numerical Recipes his opinion was something like “good book, bad code”.
