我不熟悉使用现代 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 最佳实践 所推荐的内容.
期待学习如何正确执行此操作。
这纯粹是Fortran不正确的问题,与f2py或Python无关。
你的种类常量被命名为
REAL_KIND
,因此文字必须看起来像
0.0_REAL_KIND
没有理由重复 REAL_ 这个词。该常数称为
REAL_KIND
,而不是 REAL_REAL_KIND
。