如何在 Fortran 子例程中定义意图输出的动态数组并使用 f2py 在 Python 中调用它们?

问题描述 投票:0回答:1

我不熟悉使用现代 Fortran 结构。我已成功应用动态分配数组作为输入,并在该 jupyter 笔记本单元中验证和确认了单个输出结果:

import numpy as np
import locateindex

tstarr=np.array([1,3,2,4,10])
chk=3
indxval=0

indxval=locateindex.spallocindx(tstarr,chk)
print('indxval: ',indxval)

单元格返回输入数组的正确 Fortran 索引位置(从 1 开始)。下面提供了 Fortran 子例程 spallocindx.f90:

SUBROUTINE spallocindx(indx, xarray, chkval) !RESULT(indx)

  IMPLICIT NONE
  INTEGER, PARAMETER :: INT64=8
  INTEGER(KIND=INT64), INTENT(IN) :: xarray(:)
  INTEGER(KIND=INT64), INTENT(IN) :: chkval
  INTEGER(KIND=INT64), INTENT(OUT) :: indx
  
  indx = findloc(xarray,chkval,dim=1)
 
  END SUBROUTINE spallocindx

该子例程已成功编译用于 python 中的 numpy.f2py 例程,如下所示:

(venv) MyPythonDirUbuntu$ python -m numpy.f2py -c spallocindx.f90 -m locateindex

所以,现在我尝试将这个概念扩展到求解下三角矩阵(本质上是 python scipy.sparse.linalg.spsolve_triangle 的功能)。

我的关键 Fortran 子例程是

SUBROUTINE fwdsubv2(xsol,ElIdArr,Lval,rhsvec,nsol) !RESULT(xsol)

  IMPLICIT NONE
  INTEGER, PARAMETER :: REAL_KIND=8
  INTEGER, PARAMETER :: INT64=8

  INTEGER(KIND=INT64), INTENT(IN) :: nsol
  INTEGER(KIND=INT64), INTENT(IN) :: ElIdArr(:)

  REAL(KIND=REAL_KIND), INTENT(IN) :: Lval(:), rhsvec(:)
  
  REAL(KIND=REAL_KIND), ALLOCATABLE, INTENT(OUT) :: xsol(:)

  REAL(KIND=REAL_KIND) :: rowsum, dgval 
  INTEGER(KIND=INT64) :: elidstart, elidend, dgelid, iz, offdgindx, dgindx
  INTEGER(KIND=INT64) :: ir, remraw, colid

  ALLOCATE(xsol(nsol), source=0.0_REAL_REAL_KIND)

  xsol(1)=rhsvec(1)/Lval(1)

  DO ir=2,nsol
    rowsum = 0.0
    dgelid=(ir-1)*nsol+ir
    CALL spallocindx(dgindx, ElIdArr, dgelid)

    IF (dgindx == 0) THEN
        print '(2g0)', 'WARNING: Diagonal value missing for Lower Triangle row ', ir
        dgval=0.000000001
     ELSE
        dgval=Lval(dgindx)
     END IF

    elidstart=(ir-1)*nsol+1
    elidend=(ir-1)*(nsol+1)

    DO iz = elidstart,elidend
        CALL spallocindx(offdgindx, ElIdArr, iz)
        IF (offdgindx > 0) THEN
           remraw=MOD(iz,nsol)
           IF (remraw > 0) THEN
              colid=remraw
           ELSE
              colid=nsol
           END IF
           rowsum = rowsum - Lval(iz) * xsol(colid)/dgval
        END IF
     END DO

     xsol(ir) = rhsvec(ir)/dgval + rowsum

  END DO
  
END SUBROUTINE fwdsubv2

内部子程序spalloc.f90可以工作,并已通过上面的jupyter cell验证。 当尝试使用 f2py 编译时,像以前一样使用:

(venv) MyPythonDirUbuntu$ python -m numpy.f2py -c fwdsubv2.f90 spallocindx.f90 -m tstfwdsub

我收到以下错误:

fwdsubv2.f90:19:53:

   19 |       ALLOCATE(xsol(nsol), source=0.0_REAL_REAL_KIND)
      |                                                     1
Error: Missing kind-parameter at (1)
error: Command "/usr/bin/gfortran -Wall -g -fno-second-underscore -fPIC -O3 -funroll-loops -I/tmp/tmpp07jlue4/src.linux-x86_64-3.10 -I/home/exouser/venv/lib/python3.10/site-packages/numpy/core/include -I/home/exouser/venv/include -I/usr/include/python3.10 -c -c fwdsubv2.f90 -o /tmp/tmpp07jlue4/fwdsubv2.o" failed with exit status 1

我尝试了源初始化的不同变体,例如 source=0.0 或 source=REAL_KIND(0.0) 或 source=REAL_REAL_KIND(0.0) 但显示了上面的示例,因为它似乎最接近此处链接的 Fortran 最佳实践 所推荐的内容.

期待学习如何正确执行此操作。

python arrays fortran gfortran f2py
1个回答
0
投票

这纯粹是Fortran不正确的问题,与f2py或Python无关。

你的种类常量被命名为

REAL_KIND
,因此文字必须看起来像

0.0_REAL_KIND

没有理由重复 REAL_ 这个词。该常数称为

REAL_KIND
,而不是
REAL_REAL_KIND

© www.soinside.com 2019 - 2024. All rights reserved.