# Feedback wanted on library design

**URL:** <https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085>\
**Category:** Uncategorized\
**Created:** [July 29, 2022, 1:00am UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085 "2022-07-29T01:00:37Z")\
**Posts on this page:** 9\
**Page:** 1

<div class="post-metadata">

**Author:** ![KjellJorner](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/kjelljorner/32/590_2.png) [@KjellJorner](https://fortran-lang.discourse.group/u/KjellJorner)\
**Post date:** [July 29, 2022, 1:00am UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085/1 "2022-07-29T01:00:38Z")

</div>

I’ve written a [Python C API](https://github.com/kjelljorner/libdlfind) to the molecular geometry optimization package DL-FIND. The design of DL-FIND requires the user of the package to supply a number of functions at compile time that are later called by the library during runs: `dlf_error`, `dlf_get_gradient`, `dlf_get_hessian` and so on. This design is obviously not compatible with a Python library where the user could dynamically provide different versions of these functions at runtime. My solution to this was to construct dummy functions on the Fortran side that just pass on their arguments to a similarly named procedure imported from the module `mod_globals`.

```fortran
subroutine dlf_get_gradient(nvar, coords, energy, gradient, iimage, kiter, status)
  use mod_globals, only: dlf_get_gradient_
  use dlf_parameter_module, only: rk

  implicit none
  integer, intent(in) :: nvar ! number of xyz variables (3*nat)
  real(rk), intent(in) :: coords(nvar) ! coordinates
  real(rk), intent(out) :: energy ! energy
  real(rk), intent(out) :: gradient(nvar) ! gradient
  integer, intent(in) :: iimage ! current image (for NEB)
  integer, intent(in) :: kiter ! flag related to microiterations
  integer, intent(out) :: status ! return code

  call dlf_get_gradient_(nvar, coords, energy, gradient, iimage, kiter, status)
end subroutine

```

`mod_globals` stores all the procedures as procedure pointers to the functions obtained from the C (Python) side.

```fortran
...
procedure(c_dlf_get_gradient), pointer :: dlf_get_gradient_ => null()
...

```

The abstract interfaces for these functions are defined in the `mod_api` module

```fortran
  abstract interface
    subroutine c_dlf_get_gradient(nvar, coords, energy, gradient, iimage, kiter, status)
      import c_double, c_int

      implicit none
      integer(c_int), intent(in), value :: nvar
      real(c_double), intent(in) :: coords(nvar)
      real(c_double), intent(out) :: energy
      real(c_double), intent(out) :: gradient(nvar)
      integer(c_int), intent(in), value :: iimage
      integer(c_int), intent(in), value :: kiter
      integer(c_int), intent(out) :: status
    end subroutine
  end interface

```

Finally, the way the procedure pointers in `mod_globals` get associated to the functions from the C side is by arguments as `c_funptr` to the function `api_dl_find` that is called from the C/Python side. `c_f_procpointer` is used to convert from C pointers to Fortran pointers.

```fortran
subroutine api_dl_find(... , c_dlf_get_gradient_, ...) bind(c)
...
  type(c_funptr), intent(in), value :: ..., c_dlf_get_gradient_, .... ! Functions received from C side
...
 call c_f_procpointer(c_dlf_get_gradient_, dlf_get_gradient_)
...

```

This design for the API seems to work well in my tests with different energy and gradient calculators, but also feels rather hacky. If anyone has any ideas for how to improve or simplify it would be greatly appreciated.

I’m [building](https://github.com/kjelljorner/libdlfind/blob/main/.github/workflows/build_wheels.yml) Python wheels on GitHub Actions, linking statically on Windows and “repairing” the wheels on Linux and MacOS to include the required libraries. Once I make an initial release I will continue with conda-forge packages as well.

Some references:

- [DL-FIND paper](https://10.1021/jp9028968)
- [ChemShell information pages](https://www.chemshell.org/dl-find)
- [Johannes Kästner’s information page](https://www.itheoc.uni-stuttgart.de/research/kaestner/research/dlfind/)

DL-FIND is also available via Py-ChemShell, but I wanted something that is:

- Leaner and more modular
- Compatible with arbitrary energy and gradient codes
- Available to install as a dependency via package managers

---

<div class="post-metadata">

**Author:** ![awvwgk](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/awvwgk/32/154_2.png) [@awvwgk](https://fortran-lang.discourse.group/u/awvwgk)\
**Post date:** [July 29, 2022, 1:04pm UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085/2 "2022-07-29T13:04:36Z")

</div>

There are certainly some discussions around about designing interfaces to optimization libraries, best checkout the thread at

> [@Modern Fortran interface for NLopt](https://fortran-lang.discourse.group/t/2685/):
>
> @jacobwilliams and others, nlopt looks like a great optimization tool, but its fortran interface currently is Fortran 77. Could be updated to use iso\_c\_binding.

A central topic is the design and separation of the objective function callback and the solver object. I collected a couple of links at [Callback function design · Issue #16 · grimme-lab/nlopt-f · GitHub](https://github.com/grimme-lab/nlopt-f/issues/16).

---

<div class="post-metadata">

**Author:** ![FortranFan](https://avatars.discourse-cdn.com/v4/letter/f/96bed5/32.png) [@FortranFan](https://fortran-lang.discourse.group/u/FortranFan)\
**Post date:** [July 29, 2022, 4:50pm UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085/3 "2022-07-29T16:50:25Z")

</div>

> [@KjellJorner](#):
>
> If anyone has any ideas for how to improve or simplify it would be greatly appreciated.

A quick question: does your actual code include `bind(C)` clause for the interfaces? e.g., for `c_dlf_get_gradient`?

---

<div class="post-metadata">

**Author:** ![ivanpribec](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/ivanpribec/32/3290_2.png) [@ivanpribec](https://fortran-lang.discourse.group/u/ivanpribec)\
**Post date:** [July 29, 2022, 6:23pm UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085/4 "2022-07-29T18:23:30Z")

</div>

> [@KjellJorner](#):
>
> `mod_globals` stores all the procedures as procedure pointers to the functions obtained from the C (Python) side.

I’d be very cautious placing the pointers in a global scope. This means your Python library can’t be used in a parallel context (assuming that DL-FIND is thread-safe to begin with…). You should document the “parallel status” explicitly in your README.

---

<div class="post-metadata">

**Author:** ![FortranFan](https://avatars.discourse-cdn.com/v4/letter/f/96bed5/32.png) [@FortranFan](https://fortran-lang.discourse.group/u/FortranFan)\
**Post date:** [July 30, 2022, 1:03am UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085/5 "2022-07-30T01:03:32Z")

</div>

> [@ivanpribec](#):
>
> … very cautious placing the pointers in a global scope. This means your Python library can’t be used in a parallel context (assuming that DL-FIND is thread-safe to begin with …

As opposed to data pointers, how do procedure pointers affect thread safety and/or using the library in a parallel context?

---

<div class="post-metadata">

**Author:** ![ivanpribec](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/ivanpribec/32/3290_2.png) [@ivanpribec](https://fortran-lang.discourse.group/u/ivanpribec)\
**Post date:** [July 30, 2022, 10:20am UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085/6 "2022-07-30T10:20:39Z")

</div>

It appears my understanding of how procedure pointers work is incomplete. Wouldn’t two or more threads trying to change the same procedure pointer potentially lead to a race?

The example I built to test this, shows no wrong behavior using a global pointer. (I only inspected small `n` sizes, and didn’t attempt to verify that multiple threads are actually being used.)

> **(No) Race**
>
> Compile with: `$ gfortran -Wall -fopenmp -O2 polyrace.f90`
> 
> ```fortran
> module poly
> 
> implicit none
> private
> 
> public :: the_poly, poly1, poly2, poly3, poly4, wp
> 
> integer, parameter :: wp = kind(1.0d0)
> 
> abstract interface
> function poly_interface(x) result(y)
> import wp
> real(wp), intent(in) :: x
> real(wp) :: y
> end function
> end interface
> 
> procedure(poly_interface), pointer :: the_poly => null()
> 
> contains
> 
> function poly1(x) result(y)
> real(wp), intent(in) :: x
> real(wp) :: y
> print *, "calling poly1 with x =", x
> y = 1.0_wp + x
> end function
> 
> function poly2(x) result(y)
> real(wp), intent(in) :: x
> real(wp) :: y
> print *, "calling poly2 with x =", x
> y = 1.0_wp + x*(1.0_wp + x)
> end function
> 
> function poly3(x) result(y)
> real(wp), intent(in) :: x
> real(wp) :: y
> print *, "calling poly3 with x =", x
> y = 1.0_wp + x*(1.0_wp + x*(1.0_wp + x))
> end function
> 
> function poly4(x) result(y)
> real(wp), intent(in) :: x
> real(wp) :: y
> print *, "calling poly4 with x =", x
> y = 1.0_wp + x*(1.0_wp + x*(1.0_wp + x*(1.0_wp + x)))
> end function
> 
> end module
> 
> program polyrace
> 
> use poly 
> implicit none
> 
> integer, parameter :: n = 1000
> 
> real(wp) :: y(n), r
> integer :: i 
>   
> !$omp parallel do schedule(static) shared(the_poly,y)
> do i = 1, size(y)
> 
> call random_number(r)
> if (r < 0.25_wp) then
> the_poly => poly1
> else if (r < 0.5_wp) then
> the_poly => poly2
> else if (r < 0.75_wp) then
> the_poly => poly3
> else
> the_poly => poly4
> end if 
> 
> y(i) = the_poly(r)
> 
> end do
> !$omp end parallel do
> 
> end program
> 
> ```

---

<div class="post-metadata">

**Author:** ![KjellJorner](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/kjelljorner/32/590_2.png) [@KjellJorner](https://fortran-lang.discourse.group/u/KjellJorner)\
**Post date:** [July 30, 2022, 12:15pm UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085/7 "2022-07-30T12:15:54Z")

</div>

> [@FortranFan](#):
>
> A quick question: does your actual code include `bind(C)` clause for the interfaces? e.g., for `c_dlf_get_gradient`?

[It](https://github.com/kjelljorner/libdlfind/blob/main/src/mod_api.f90) does not. After some experimentation, compilation seemed to crash when I included it, if I recall correctly. I’m not sure if it’s because I am using the `c_f_procpointer` method. The GCC [documentation](https://gcc.gnu.org/onlinedocs/gcc-9.3.0/gfortran/Working-with-Pointers.html) that I followed omits `bind(c)`.

---

<div class="post-metadata">

**Author:** ![KjellJorner](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/kjelljorner/32/590_2.png) [@KjellJorner](https://fortran-lang.discourse.group/u/KjellJorner)\
**Post date:** [July 30, 2022, 12:37pm UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085/8 "2022-07-30T12:37:50Z")

</div>

> [@ivanpribec](#):
>
> I’d be very cautious placing the pointers in a global scope. This means your Python library can’t be used in a parallel context (assuming that DL-FIND is thread-safe to begin with…). You should document the “parallel status” explicitly in your README.

The problem is that DL-FIND expects to find the callback functions in the global scope. I don’t really want to change the design of DL-FIND itself (to be able to upgrade to the latest version with minimal changes) and storing the pointers in a module scope was the only solution that I could come up with where they would be accessible from both (a) the function called from the C side and (b) the dummy compiled callback. I will investigate the thread safety, but it could potentially work well with multiprocessing where each Python process would get its own copy of the library. DL-FIND has internal parallelization via MPI which I have not investigated.

---

<div class="post-metadata">

**Author:** ![KjellJorner](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/kjelljorner/32/590_2.png) [@KjellJorner](https://fortran-lang.discourse.group/u/KjellJorner)\
**Post date:** [August 3, 2022, 3:32pm UTC](https://fortran-lang.discourse.group/t/feedback-wanted-on-library-design/4085/9 "2022-08-03T15:32:12Z")

</div>

Thanks everyone who gave me feedback on this project. As a result I made a number of changes that I think improved the library a lot:

- Added `bind(c)` to the abstract interfaces on the suggestion of @FortranFan. This seems to be an omission in the example in the GCC [documentation](https://gcc.gnu.org/onlinedocs/gcc-12.1.0/gfortran/Working-with-C-Pointers.html#Working-with-C-Pointers).
- Reworked most of the factory functions into decorators on the suggestion of @ivanpribec
- Clarified the thread safety in the README.md on the suggestion of @ivanpribec. Seems to crash with multiple threads but runs fine with multiple processes.

I opted not to make any greater changes to the design of the API to be as close as possible to DL-FIND itself. This saves me from having to write a detailed documentation as I can just refer to the documentation of DL-FIND.
