Вызов Fortran DLL из python с использованием CFFI с многомерными массивами - PullRequest
Новая
       68

Вызов Fortran DLL из python с использованием CFFI с многомерными массивами

0 голосов
/ 17 января 2020

Я использую DLL, которая содержит решатели дифференциальных уравнений среди других полезных математических инструментов. К сожалению, эта DLL написана на Фортране. Моя программа написана на python 3.7, и я использую spyder в качестве IDE.

Я успешно вызвал легкие функции из dll. Однако я не могу заставить работать функции, требующие многомерных массивов.

Это онлайн-документация к функции, которую я пытаюсь вызвать: https://www.nag.co.uk/numeric/fl/nagdoc_fl26/html/f01/f01adf.html

Ядро умирает без сообщения об ошибке, если я выполню следующий код: 

import numpy as np
import cffi as cf

ffi=cf.FFI()
lib=ffi.dlopen("C:\Windows\SysWOW64\DLL20DDS")
ffi.cdef("""void F01ADF (const int *n, double** a, const int *lda, int *ifail);""")

#Integer
nx = 4
n = ffi.new('const int*', nx)
lda = nx + 1
lda = ffi.new('const int*', lda)
ifail = 0
ifail = ffi.new('int*', ifail)

#matrix to be inversed

ax1 = np.array([5,7,6,5],dtype = float, order = 'F')
ax2 = np.array([7,10,8,7],dtype = float, order = 'F')
ax3 = np.array([6,8,10,9],dtype = float, order = 'F')
ax4 = np.array([5,7,9,10], dtype = float, order = 'F')
ax5 = np.array([0,0,0,0], dtype = float, order = 'F')

ax = (ax1,ax2,ax3,ax4,ax5)

#Array
zx = np.zeros(nx, dtype = float, order = 'F')

a = ffi.cast("double** ", zx.__array_interface__['data'][0])
for i in range(lda[0]):
    a[i] = ffi.cast("double* ", ax[i].__array_interface__['data'][0])

lib.F01ADF(n, a, lda, ifail)

Поскольку функция с одномерными массивами работает, я предполагаю, что проблема заключается в многомерном массиве.

Любая помощь очень ценится, Тило

1 Ответ

1 голос
/ 24 января 2020

Отсутствие доступа к dll, на который вы ссылаетесь, затрудняет точный ответ, однако документации по dll и предоставленному сценарию Python может быть достаточно для диагностики проблемы. В вашем примере есть как минимум две проблемы:

  • Интерфейс заголовка C:

    В вашей ссылке на документацию четко указано, что функция * Интерфейс заголовка 1051 * должен выглядеть так. Я не очень хорошо разбираюсь в C, Python в cffi или cdef, но объявление параметра для a в интерфейсе вашей функции кажется неправильным. double** a (указатель на двойной указатель) в вашем интерфейсе функции, скорее всего, должен быть double a[] или double* a (указатель на двойной), как указано в документации.

  • Определение массива 2d Numpy с порядком Фортрана:

    1. Обратите внимание, что Numpy массивы ax1..5 являются одномерными массивами, поскольку массивы имеют только одно измерение order='F' и order='C' эквивалентны с точки зрения расположения памяти и доступа. Таким образом, указание order='F' здесь, вероятно, не даст ожидаемого эффекта (Fortran использует упорядочение по главному столбцу для многомерных массивов).
    2. Переменная ax является кортежем Numpy массивов, а не массив 2d Numpy и, следовательно, будет иметь совсем другое представление в памяти (что крайне важно при передаче данных в dll Фортран), чем массив 2d.

На пути к решению

Моим первым шагом будет исправление интерфейса заголовка C. Затем я объявил бы ax как правильный массив Numpy с двумя измерениями, используя порядок Фортрана, а затем привел бы его к соответствующему типу данных, как в этом примере: 

#file: test.py
import numpy as np
import cffi as cf

ffi=cf.FFI()
lib=ffi.dlopen("./f01adf.dll")
ffi.cdef("""void f01adf_ (const int *n, double a[], const int *lda, int *ifail);""")

# integers
nx = 4
n = ffi.new('const int*', nx)
lda = nx + 1
lda = ffi.new('const int*', lda)
ifail = 0
ifail = ffi.new('int*', ifail)

# matrix to be inversed
ax = np.array([[5,  7,  6,  5],
               [7, 10,  8,  7],
               [6,  8, 10,  9],
               [5,  7,  9, 10],
               [0,  0,  0,  0]], dtype=float, order='F')

#  operation on matrix using dll
print("BEFORE:")
print(ax.astype(int))

a = ffi.cast("double* ", ax.__array_interface__['data'][0])
lib.f01adf_(n, a, lda, ifail)

print("\nAFTER:")
print(ax.astype(int))

Для целей тестирования рассмотрите следующую подпрограмму Fortran, которая имеет тот же интерфейс, что и ваша фактическая dll , в качестве замены вашей dll. Он просто добавит 10 ** (i-1) в i-й столбец входного массива a. Это позволит проверить, что интерфейс между Python и Fortran работает так, как предполагалось, и что предполагаемые элементы массива a работают: 

!file: f01adf.f90
Subroutine f01adf(n, a, lda, ifail)
  Integer, Intent (In) :: n, lda
  Integer, Intent (Inout) :: ifail
  Real(Kind(1.d0)), Intent (Inout) :: a(lda,*)

  Integer :: i

  print *, "Fortran DLL says: Hello world!"

  If ((n < 1) .or. (lda < n+1)) Then
    ! Input variables not conforming to requirements
    ifail = 2

  Else
    ! Input variables acceptable
    ifail = 0

    ! add 10**(i-1) to the i'th column of 2d array 'a'
    Do i = 1, n
      a(:, i) = a(:, i) + 10**(i-1)
    End Do

  End If

End Subroutine

Компиляция кода Fortran, а затем запуск предлагаемого Python скрипт, дает мне следующий вывод: 

> gfortran -O3 -shared -fPIC -fcheck=all -Wall -Wextra -std=f2008 -o f01adf.dll f01adf.f90
> python test.py
  BEFORE:
  [[ 5  7  6  5]
  [ 7 10  8  7]
  [ 6  8 10  9]
  [ 5  7  9 10]
  [ 0  0  0  0]]

  Fortran DLL says: Hello world!

  AFTER:
  [[   6   17  106 1005]
  [   8   20  108 1007]
  [   7   18  110 1009]
  [   6   17  109 1010]
  [   1   10  100 1000]]
...