Goal: Get this minpack used by SciPy

[certik]

SciPy contains minpack: https://github.com/scipy/scipy/tree/main/scipy/optimize/minpack; we should develop this fortran-lang/minpack in a way so that SciPy could use it. We should add some improvements (perhaps some new minimization algorithms), and then include them in SciPy. If we can pull this off, that would be a huge win, for both Fortran and SciPy. SciPy gets a maintainable version, we can all contribute to it in Fortran, and modern Fortran gets a "user" in a high profile Python library.

Fortran Forum discussion thread: https://fortran-lang.discourse.group/t/lets-get-fortran-lang-minpack-into-scipy/2722

[jacobwilliams]

One thing TODO: check to see what changes they have made to the original code and port those to this version. See #23.

[milancurcic]

I agree, this would be a huge boon for Fortran and SciPy. Great idea!

[awvwgk]

Here is a quick outline of the things I found in SciPy on this:

• minpack itself is compiled as static library (this is the part we actually care about here):
  • scipy/optimize/meson.build#L9-L35
• a Python extension module (scipy.optimize._minpack) is defined directly via the Python and numpy headers
  • scipy/optimize/_minpackmodule.c
• two headers are defined for this extension module:
  • scipy/optimize/minpack.h
  • scipy/optimize/__minpack.h
• the extension module for minpack is compiled in meson
  • scipy/optimize/meson.build#L37-L46
• a Python wrapper module (scipy.optimize._minpack_py) is build around this setup
  • scipy/optimize/_minpack_py.py
• this module is reexported as scipy.optimize.minpack, however direct access seems deprecated
  • scipy/optimize/minpack.py
• instead of actually using scipy.optimize.minpack for anything, scipy.optimize._minpack_py is reexported
  • scipy/optimize/__init__.py#L403

---

Note that there is a scipy.optimize.minpack2 module, this is interfaced via f2py

• interface definition in scipy/optimize/minpack2/minpack2.pyf
• generation of extension module at scipy/optimize/meson.build#L136-L151
• reexport in as Python module in scipy/optimize/minpack2.py
• not sure whether this module is later reexported somewhere else

[sakamoti]

How about considering adding optional argument class(*) for fortran user?

In the f77 version, the user defined subroutine, which is the argument like hybrd1, is defined as func(n, x, fvec, iflag), etc. The user defined derived type cannot be passed.
So, I use module variables or the internal subroutine to refer to the type.

By adding class(*) as an optional argument to func, for example, I think it is straightforward way to pass user-defined types. Personally, Such an additional argument class(*),optional will provide well userbility for fortran user.

Is it practically possible to add an optional argument like class(*),optional when considering the integration with Scipy?

[jacobwilliams]

See my comments in #4. I think that's the way to go, rather than having an optional argument.

[awvwgk]

There is also a relevant discussion on the callback design in https://fortran-lang.discourse.group/t/modern-fortran-interface-for-nlopt/2685/6 and grimme-lab/nlopt-f#16. I think this is a quite fundamental issue, and I'm hitting against the limitation of the current callback when building the C API already.

[sakamoti]

@jacobwilliams
Thank you for the link. I had missed that discussion. The code in the link is enough to serve the purpose.

@awvwgk
I really appreciate it, because discussions like the one you linked to are hard to find outside of fortran-lang.

[awvwgk]

Is it practically possible to add an optional argument like class(*),optional when considering the integration with Scipy?

I would go so far to say that is necessary to pass additional data with the procedure pointer to allow integration with SciPy, Right now we have to rely on C and Fortran using similar enough calling conventions to use the same callback function interface. Once we switch to a more Fortran-friendly interface (assumed shape arrays instead of explicit shape arrays) this is no longer true.

[Nicholaswogan]

One other thing to add that I think is relevant. Currently released GNU compilers (e.g. GNU Fortran (Homebrew GCC 11.2.0_3) 11.2.0) do not support nested functions on Apple Silicon. This will change someday but it isn't clear when (see here). If Scipy relies on GNU compilers to build binaries, then passing data via nested functions will prevent building Apple Silicon binaries.

[ilayn]

Hi folks big fan of this initiative!

Let me know if you need any help from SciPy side please, @mdhaber would also be interested I presume, as the scipy.optimize resident and may involve others.

There is quite some legacy luggage with respect to xPACK routines in scipy. Hence if you think there is a more modern and "sane" way to wrap it, we are all ears. That is to say, please don't feel obliged to fit a circle to our strange square build mechanisms.

[mdhaber]

Yes, this sounds like it could be good. Have you found any issues that this would fix?

[certik]

@ilayn, @mdhaber thanks! Here is list of ideas that @ivan-pi collected that we could fix in minpack that would fix several SciPy issues:

• Port features from SciPy's minpack version #23 (comment)

If you have any other ideas, please let us know. We'll try to fix some of these, to "motivate" the transition to the new version, so that SciPy gains something concrete from it. Of those issues, which do you think would be the highest priority for SciPy?

[mdhaber]

Of those issues, which do you think would be the highest priority for SciPy?

Personally I'd suggest bug fixes first, then new features. Otherwise, I don't have a sense of priority, so I'd just go from easiest to hardest.

[awvwgk]

Is anyone opposed to adding meson build files to minpack? This way SciPy could easily give our minpack a try via a meson subproject now that the C API is merged.

[certik]

Yes, go ahead and add Meson!

[ivan-pi]

Given that SciPy moved to Meson it makes a lot of sense to support a build-system which makes it easier for them. 👍🏻

[awvwgk]

Sounds good, we now have meson support.

[ilayn]

@rgommers Do we need anything Meson specific with the modules? I think this is going to fit like a glove otherwise.

[rgommers]

No, things should work out of the box easily build-system wise. The work will be to do the integration step (differences in Fortran sources, f2py wrapper changes if needed, etc.).

Here is a quick outline of the things I found in SciPy on this:
....

Note that the header files, minpack.py, minpack2.py etc. are all private. The only thing we care about is the public functionality in the scipy.optimize namespace and that all the tests pass. Anything else is an implementation detail that can be changed as needed.

[certik]

Awesome, thanks @rgommers for the feedback. I think the move to Meson was a good move for SciPy, it should make this a lot easier.

[awvwgk]

I started working on Python bindings for Minpack in #49 using CFFI. To keep this thread on topic I opened a separate issue to discuss the general design of the Python API in #43, feel free to add links and resources to the description.

So far I only provide wrappers for the high-level functions (hybrd1, hybrj1, ...), with #49 the following snippet should work:

import numpy as np
from math import sqrt
from minpack import hybrd1

def fcn(x, fvec) -> None:
    """
    The problem is to determine the values of x(1), x(2), ..., x(9),
    which solve the system of tridiagonal equations.::

        (3-2*x(1))*x(1)           -2*x(2)                   = -1
                -x(i-1) + (3-2*x(i))*x(i)         -2*x(i+1) = -1, i=2-8
                                    -x(8) + (3-2*x(9))*x(9) = -1
    """

    for k in range(x.size):
        tmp = (3.0 - 2.0 * x[k]) * x[k]
        tmp1 = x[k - 1] if k > 0 else 0.0
        tmp2 = x[k + 1] if k < len(x) - 1 else 0.0
        fvec[k] = tmp - tmp1 - 2.0 * tmp2 + 1.0

x = np.array(9 * [-1.0])
fvec = np.zeros(9, dtype=np.float64)
tol = sqrt(np.finfo(np.float64).eps)

hybrd1(fcn, x, fvec, tol)

print(x)
# [-0.57065451 -0.68162834 -0.70173245
#  -0.70421294 -0.70136905 -0.69186564
#  -0.66579201 -0.5960342  -0.41641206]

Feedback is welcome. I have some experience in interop, but since I'm mostly working on the compiled language side, my Python skill is somewhat lacking behind ;).

[awvwgk]

Any thoughts on releasing the current version of minpack and start packaging it for conda-forge and other?

• see Release version 2.0.0 #44

We could create a release candidate, but I doubt that there will be significant changes in the API in the next time.

[certik]

Sounds good to me.

[xecej4]

After reading some of the bug reports on the use of Minpack in SciPy ( https://github.com/scipy/scipy/issues?page=2&q=minpack+bug ) , I came away with the feeling that bugs may be located in many places in the calling chain and the wrappers, or caused by compiling the Fortran part with incorrect compiler options. I suspect that most of the errors are not related to the stable Minpack Fortran code from Netlib, etc. However, if we replace the current Minpack sources with the new Fortran-Lang sources, and problems occur, the new Minpack code is likely to be suspected, even if the true culprit is the wrapper code or something else. Are there tests to verify the Fortran part independently, without the involvement of the wrapper code, Python, etc. ?

[certik]

@xecej4, good question. Yes, there are Fortran tests here: https://github.com/fortran-lang/minpack/tree/main/test and we can/will add more as needed. We also maintain Python wrappers in this repository here: https://github.com/fortran-lang/minpack/tree/main/python, we should add test for them as well. I suggest that we are responsible for the Fortran code and the wrappers and ensure everything works.

Then when a bug in SciPy happens, we simply reproduce it, test it, fix it. If things work in this repository, then we know the bug is caused by the way SciPy includes or compiles Minpack, either Fortran or the Python wrappers, so we can focus on just that and help them fix it.

That way there is a clean separation and it will workout great.

[ilayn]

Just to add some SciPy take, the first place we look at for faults are indeed our wrappers and other plumbing code that bridges over to compiled bits. Fortunately, there is not much in between and most bridges are quite battle tested. Thus if there is some trivial wrapping mistake we happen to find them quickly.

I, and maybe others would too, would like to help you folks with the tests and maybe carry them over to here from SciPy side to kickstart the testing suite.

Our biggest issue is that we keep track of patches of Fortran bits within SciPy with limited maintainer knowledge and not being able to judge fairly whether certain things are suitable for fixing-- say things like typed things with REAL*8 like stuff, threadsafe directives or other Fortran specific behaviors that we might not understand and/or evaluate accurately.

In that sense, while I don't foresee a big maintenance burden since it is stabilized by now for a long time already, there might be some bits that maybe modernized or replaced. It would be nice if we can follow your lead on that.

[awvwgk]

I just created a first release candidate at https://github.com/fortran-lang/minpack/releases/tag/v2.0.0-rc.1, including some Python wheels and sdists. Please let me know whether those work as expected.

[arjenmarkus]

Nice work, just a few remarks here: when I built the library and the test programs I got a number of warnings about array accesses out of bound from gfortran, While the message for "if ( k /= 1) temp1 = x(k - 1)" is clearly a false positive, the code for the test programs also contains constructions like: do icase = 1, ncases+1 if (icase == ncases+1) then write (nwrite, '(a,i3,a/)') '1SUMMARY OF ', lnp, ' TESTS OF CHKDER' write (nwrite, '(a/)') ' NPROB N STATUS ERRMIN ERRMAX' do i = 1, lnp write (nwrite, '(i4, i6, 6x, l1, 3x, 2d15.7)') np(i), na(i), a(i), errmin(i), errmax(i) end do stop else nprob = nprobs(icase) n = ns(icase) ldfjac = n Why do it that way? Why not simply put writing the summary AFTER the do-loop over the cases? There is also a warning about the function dfloat(i): Warning: 'dfloat' declared at (1) may shadow the intrinsic of the same name. In order to call the intrinsic, explicit INTRINSIC declarations may be required. [-Wintrinsic-shadow] I am not sure if any of these would be of influence of influence on the acceptance of the code, but I thought I'd bring it up. Op vr 20 mei 2022 om 21:53 schreef Sebastian Ehlert < ***@***.***>:

…

I just created a first release candidate at https://github.com/fortran-lang/minpack/releases/tag/v2.0.0-rc.1, including some Python wheels and sdists. Please let me know whether those work as expected. — Reply to this email directly, view it on GitHub <#14 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AAN6YR27O27IMO7FDIN4BYDVK7USXANCNFSM5NMJVQKQ> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

[awvwgk]

One thing where we could use some helping hands is the wheel builder, see .github/workflows/wheel.yml, which currently only creates Manylinux wheels (Linux/x86_64, CPython 3.7–3.10). Having MacOS and Windows wheels would be important for a PyPI upload once we release 2.0.0 (see #44). Building PyPy wheels on those platforms might be feasible as well.

[rgommers]

Maybe I'm missing some context @awvwgk, but I hope you are not working on wheels for SciPy's sake? We'd want to depend on this minpack as a git submodule, not as a separate dependency.

[awvwgk]

We'd want to depend on this minpack as a git submodule, not as a separate dependency.

A meson subproject seems preferable than. However, for channels like conda-forge I would prefer to not have SciPy vendor dependencies which can be packaged separately. Question is whether you want to use the provided Python bindings of this project or write your own bindings to the Minpack C-API.

[rgommers]

Sure, we can do a Meson subproject I think. That has been discussed for a few other git submodules as well - Boost in particular. Linux distros with strict rules may be happy with that. And I think we'd be happy to adopt the Python bindings from this project.

That said, I don't think we'd want a new runtime dependency for the scipy conda-forge package. But that's a thing that can be revisited later/separately.

[h-vetinari]

Hi, one of the conda-forge feedstock maintainers for scipy (not claiming to speak for the others).

A meson subproject seems preferable than. However, for channels like conda-forge I would prefer to not have SciPy vendor dependencies which can be packaged separately.

That's indeed the best practice for conda-forge.

Question is whether you want to use the provided Python bindings of this project or write your own bindings to the Minpack C-API.

I'm not 100% wure I understand the situation correctly, but from conda-forge side, it's not a problem to change build invocations and (less ideal but still possible) even patch things in or out, but I'd really like to keep the code-as-packaged as closely possible the same as what scipy tests and releases. IOW, switching the bindings should only happen if we can reasonably guarantee that they are using the exact same definition/release/commit