基于mex文件的Matlab和Fortran混合编程问题

AsIlll · 发表于 2023-11-3 12:14:04

最近正在学习基于mex文件的Matlab和Fortran的混合编程，编写了下面这样一个f90文件，在matlab中编译成mex_function也没有报错，但是一调用生成的function，Matlab就直接崩溃，我是初学者不知道是什么问题，求助大家看看代码有没有什么问题，感激不尽。

[Fortran] 纯文本查看 复制代码

001#include "fintrf.h"
002 
003 
004subroutine mexFunction(OutSum,OutVar,InSum,InVar)         !函数接口名称必须为mexFunction,
005 
006    !OutSum:输出参数个数
007    !OutVar:输出参数数组指针
008    !InSum:输入参数个数
009    !InVar:输入参数数组指针
010 
011    !参数顺序不能随意更改
012 
013    Integer InSum,OutSum
014     
015    mwPointer InVar(*),OutVar(*)                           !mwPointer专门用于表示指针变量,这个不能随意用Integer代替
016     
017    mwPointer mxGetPr, mxGetDimensions, mxCreateNumericArray  !这个对返回指针函数的再次申明,
018    integer, parameter :: dp = 8
019    Integer , parameter :: myINT  = SELECTED_INT_KIND(8)
020 
021    Real(dp),Allocatable :: P(:,:,:,:), F1(:,:,:,:), F2(:,:,:,:), Gamma(:,:,:,:), rc(:,:,:,:), uind(:,:,:,:) ! input and output uind
022    Integer,Allocatable :: Dim_P(:), Dim_F(:), Dim_G(:), Dim_rc(:) ! input array dimension
023    Real(dp) :: n,co   ! input n and cut_off distance
024    Real(dp) :: uindx,uindy,uindz,temuindx,temuindy,temuindz
025 
026    Integer I1,J1,K1,I2,J2,K2,NBP,NTP,NSP,NBF,NTF,NSF,PD,FD,GD,RD
027 
028    If(InSum/=7)Then
029     
030        call mexErrMsgIdAndTxt('MATLAB:InputTooBig','输入参数个数必须为7个')
031     
032        Return
033     
034    EndIf
035 
036    PD = mxGetNumberOfDimensions(P)
037    FD = mxGetNumberOfDimensions(F1)
038    GD = mxGetNumberOfDimensions(Gamma)
039    RD = mxGetNumberOfDimensions(rc)
040 
041 
042    Allocate(Dim_P(PD), Dim_F(FD), Dim_G(GD), Dim_rc(RD))
043 
044    Dim_P = mxGetDimensions(InVar(1)) !获取第1个输入参数P的维度
045    Dim_F = mxGetDimensions(InVar(2)) !获取第2 (3)个输入参数F1 (F2)的维度
046    Dim_G = mxGetDimensions(InVar(4)) !获取第4个输入参数Gamma的维度
047    Dim_rc = mxGetDimensions(InVar(5)) !获取第5个输入参数rc的维度  
048     
049     
050    Allocate(P(Dim_P(1),Dim_P(2),Dim_P(3),Dim_P(4)))
051    Allocate(F1(Dim_F(1),Dim_F(2),Dim_F(3),Dim_F(4)))
052    Allocate(F2(Dim_F(1),Dim_F(2),Dim_F(3),Dim_F(4)))
053    Allocate(Gamma(Dim_G(1),Dim_G(2),Dim_G(3),Dim_G(4)))
054    Allocate(rc(Dim_rc(1),Dim_rc(2),Dim_rc(3),Dim_rc(4)))
055 
056    Allocate(Uind(Dim_P(1),Dim_P(2),Dim_P(3),Dim_P(4)))
057     
058 
059    Call mxCopyPtrToReal8(mxGetPr(InVar(1)),P,Dim_P(1)*Dim_P(2)*Dim_P(3)*Dim_P(4)) !将第1个参数数组赋值给P变量
060    Call mxCopyPtrToReal8(mxGetPr(InVar(2)),F1,Dim_F(1)*Dim_F(2)*Dim_F(3)*Dim_F(4)) !将第2个参数数组赋值给F1变量
061    Call mxCopyPtrToReal8(mxGetPr(InVar(3)),F2,Dim_F(1)*Dim_F(2)*Dim_F(3)*Dim_F(4)) !将第3个参数数组赋值给F2变量
062    Call mxCopyPtrToReal8(mxGetPr(InVar(4)),Gamma,Dim_G(1)*Dim_G(2)*Dim_G(3)*Dim_G(4)) !将第4个参数数组赋值给Gamma变量
063    Call mxCopyPtrToReal8(mxGetPr(InVar(5)),rc,Dim_rc(1)*Dim_rc(2)*Dim_rc(3)*Dim_rc(4)) !将第5个参数数组赋值给rc变量
064    Call mxCopyPtrToReal8(mxGetPr(InVar(6)),n,1)!将第6个整数变量赋值给n
065    Call mxCopyPtrToReal8(mxGetPr(InVar(7)),co,1)!将第7个整数变量赋值给co cut_off_distance
066     
067    NBP = Dim_P(1)
068    NTP = Dim_P(3)
069    NSP = Dim_P(4)
070    NBF = Dim_F(1)
071    NTF = Dim_F(3)
072    NSF = Dim_F(4)
073 
074    DO K1 = 1,NBP
075        DO J1 = 1,NTP
076            DO I1 = 1,NSP
077 
078                temuindx = 0.0
079                temuindy = 0.0
080                temuindz = 0.0
081 
082                DO K2 = 1,NBF
083                    DO J2 = 1,NTF
084                        DO I2 = 1,NSF
085 
086                            call Biot_p2p(P(I1,1,J1,K1),P(I1,2,J1,K1),P(I1,3,J1,K1),F1(I2,1,J2,K2),F1(I2,2,J2,K2),F1(I2,3,J2,K2),F2(I2,1,J2,K2),F2(I2,2,J2,K2),F2(I2,3,J2,K2),&
087                                    Gamma(I1,1,J1,K1),rc(I1,1,J1,K1),co,n,temuindx,temuindy,temuindz)
088 
089                             
090                            temuindx = temuindx + temuindx
091                            temuindy = temuindy + temuindy
092                            temuindy = temuindy + temuindy
093 
094 
095 
096                        ENDDO
097                    ENDDO
098                ENDDO
099 
100 
101                Uind(I1,1,J1,K1) = temuindx
102                Uind(I1,2,J1,K1) = temuindy
103                Uind(I1,3,J1,K1) = temuindz
104 
105 
106 
107 
108            ENDDO
109        ENDDO
110    ENDDO
111 
112    OutVar(1)=mxCreateNumericArray(PD,Dim_P,REAL*8,0)!给返回参数分配内存
113 
114    Call mxCopyReal8ToPtr(Uind,mxGetPr(OutVar(1)),Dim_P(1)*Dim_P(2)*Dim_P(3)*Dim_P(4))!将返回参数赋值给分配的内存
115     
116    DeAllocate(P,F1,F2,Gamma,rc)!释放临时分配的内存
117 
118    Return
119 
120End SubRoutine
121 
122Subroutine Biot_p2p(px,py,pz,f1x,f1y,f1z,f2x,f2y,f2z,gamma,rc,co,n,temuindx,temuindy,temuindz)
123 
124    implicit none
125 
126    integer, parameter :: dp = 8
127    Real(dp),Intent(In):: px,py,pz,f1x,f1y,f1z,f2x,f2y,f2z,gamma
128    Real(dp),Intent(Out):: temuindx,temuindy,temuindz
129    Real(dp) :: ldx, ldy, ldz, pxx1, pxx2, pyy1, pyy2, pzz1, pzz2, r1, r2, ldr, len, r1dr2, r1tr2, cnu, den, ubar,cmn,rc
130    Real ::  co,cutoff,cored,n
131    Real, parameter :: pi = acos(-1.0)
132    Real, parameter :: infinity = 999999.0
133 
134    cutoff = co
135    cored = n
136 
137     
138    temuindx = 0.0
139    temuindy = 0.0
140    temuindz = 0.0
141 
142    ldx = f2x-f1x
143    ldy = f2y-f1y
144    ldz = f2z-f1z
145 
146    pxx1 = px-f1x
147    pxx2 = px-f2x
148    pyy1 = py-f1y
149    pyy2 = py-f2y
150    pzz1 = pz-f1z
151    pzz2 = pz-f2z
152 
153    cmn = cored
154 
155    r1 = sqrt(pxx1*pxx1+pyy1*pyy1+pzz1*pzz1)
156    r2 = sqrt(pxx2*pxx2+pyy2*pyy2+pzz2*pzz2)
157 
158 
159    if((r1<cutoff).AND.(r2<cutoff)) then
160     
161        ldr = (ldx*pxx1 + ldy*pxx2 + ldz*pzz1)
162        ldr = ldr*ldr
163        len = ldx*ldx + ldy*ldy + ldz*ldz !!L^2
164        r1dr2 = pxx1*pxx2+pyy1*pyy2+pzz1*pzz2
165        r1tr2 = r1*r2
166 
167        cnu = (r1*r1)-(ldr/len)
168        cnu = cnu*((rc**cmn+cnu**cmn)**(-1.0/cmn))
169 
170        den = 4.0*pi*(r1tr2*(r1tr2 + r1dr2))
171 
172        ubar = cnu*gamma*(r1+r2)
173        ubar = ubar/den
174 
175        if (isnan(ubar)) then
176           temuindx = 0.0
177           temuindy = 0.0
178           temuindz = 0.0
179        else if (ubar.GE.infinity) then
180           temuindx = 0.0
181           temuindy = 0.0
182           temuindz = 0.0        
183        else
184            temuindx = ubar*(pyy1*pzz2-pzz1*pyy2)
185            temuindy = ubar*(pzz1*pxx2-pxx1*pzz2)
186            temuindz = ubar*(pxx1*pyy2-pyy1*pxx2)
187        endif
188 
189    endif
190 
191 
192End SubRoutine Biot_p2p

崩溃信息：This error was detected while a MEX-file was running. If the MEX-file

is not an official MathWorks function, please examine its source code

for errors. Please consult the External Interfaces Guide for information

on debugging MEX-files.

Alohafancy · 发表于 2024-3-8 16:16:41

AsIlll 发表于 2023-11-29 22:02
已经解决了。谢谢大家

楼主怎么解决的？就是改了变量声明的精度吗？

AsIlll · 发表于 2023-11-29 22:02:04

已经解决了。谢谢大家

AsIlll · 发表于 2023-11-12 11:22:04

li913 发表于 2023-11-5 19:25
matlab 默认浮点数是双精度，整数是8字节。

抱歉，我没太明白您的意思，是声明的问题吗，我又声明了一下变量的精度，浮点数为双精度，整数为8字节，但是仍然出现了同样的问题

li913 · 发表于 2023-11-5 19:25:19

matlab 默认浮点数是双精度，整数是8字节。

AsIlll · 发表于 2023-11-4 18:53:34

本帖最后由 AsIlll 于 2023-11-4 18:54 编辑

我将计算语句全部删除之后，仍然出现了程序崩溃的问题，应该是接口的问题，但是我没看出来是什么问题，请问有大佬帮我看看吗

[Fortran] 纯文本查看 复制代码

01#include "fintrf.h"
02 
03 
04subroutine mexFunction(OutSum,OutVar,InSum,InVar)         !函数接口名称必须为mexFunction,
05 
06    !OutSum:输出参数个数
07    !OutVar:输出参数数组指针
08    !InSum:输入参数个数
09    !InVar:输入参数数组指针
10 
11    !参数顺序不能随意更改
12 
13    Integer InSum,OutSum
14     
15    mwPointer InVar(*),OutVar(*)                           !mwPointer专门用于表示指针变量,这个不能随意用Integer代替
16     
17    mwPointer mxGetPr, mxGetDimensions, mxCreateNumericArray  !这个对返回指针函数的再次申明,
18    integer, parameter :: dp = 8
19    Integer , parameter :: myINT  = SELECTED_INT_KIND(8)
20 
21    Real(dp),Allocatable :: P(:,:,:,:), F1(:,:,:,:), F2(:,:,:,:), Gamma(:,:,:,:), rc(:,:,:,:), uind(:,:,:,:) ! input and output uind
22    Integer,Allocatable :: Dim_P(:), Dim_F(:), Dim_G(:), Dim_rc(:) ! input array dimension
23    Real(dp) :: n,co   ! input n and cut_off distance
24    Real(dp) :: uindx,uindy,uindz,temuindx,temuindy,temuindz
25 
26    Integer I1,J1,K1,I2,J2,K2,NBP,NTP,NSP,NBF,NTF,NSF,PD,FD,GD,RD
27 
28    If(InSum/=7)Then
29     
30        call mexErrMsgIdAndTxt('MATLAB:InputTooBig','输入参数个数必须为7个')
31     
32        Return
33     
34    EndIf
35 
36    PD = mxGetNumberOfDimensions(P)
37    FD = mxGetNumberOfDimensions(F1)
38    GD = mxGetNumberOfDimensions(Gamma)
39    RD = mxGetNumberOfDimensions(rc)
40 
41 
42    Allocate(Dim_P(PD), Dim_F(FD), Dim_G(GD), Dim_rc(RD))
43 
44    Dim_P = mxGetDimensions(InVar(1)) !获取第1个输入参数P的维度
45    Dim_F = mxGetDimensions(InVar(2)) !获取第2 (3)个输入参数F1 (F2)的维度
46    Dim_G = mxGetDimensions(InVar(4)) !获取第4个输入参数Gamma的维度
47    Dim_rc = mxGetDimensions(InVar(5)) !获取第5个输入参数rc的维度  
48     
49     
50    Allocate(P(Dim_P(1),Dim_P(2),Dim_P(3),Dim_P(4)))
51    Allocate(F1(Dim_F(1),Dim_F(2),Dim_F(3),Dim_F(4)))
52    Allocate(F2(Dim_F(1),Dim_F(2),Dim_F(3),Dim_F(4)))
53    Allocate(Gamma(Dim_G(1),Dim_G(2),Dim_G(3),Dim_G(4)))
54    Allocate(rc(Dim_rc(1),Dim_rc(2),Dim_rc(3),Dim_rc(4)))
55 
56    Allocate(Uind(Dim_P(1),Dim_P(2),Dim_P(3),Dim_P(4)))
57     
58 
59    Call mxCopyPtrToReal8(mxGetPr(InVar(1)),P,Dim_P(1)*Dim_P(2)*Dim_P(3)*Dim_P(4)) !将第1个参数数组赋值给P变量
60    Call mxCopyPtrToReal8(mxGetPr(InVar(2)),F1,Dim_F(1)*Dim_F(2)*Dim_F(3)*Dim_F(4)) !将第2个参数数组赋值给F1变量
61    Call mxCopyPtrToReal8(mxGetPr(InVar(3)),F2,Dim_F(1)*Dim_F(2)*Dim_F(3)*Dim_F(4)) !将第3个参数数组赋值给F2变量
62    Call mxCopyPtrToReal8(mxGetPr(InVar(4)),Gamma,Dim_G(1)*Dim_G(2)*Dim_G(3)*Dim_G(4)) !将第4个参数数组赋值给Gamma变量
63    Call mxCopyPtrToReal8(mxGetPr(InVar(5)),rc,Dim_rc(1)*Dim_rc(2)*Dim_rc(3)*Dim_rc(4)) !将第5个参数数组赋值给rc变量
64    Call mxCopyPtrToReal8(mxGetPr(InVar(6)),n,1)!将第6个整数变量赋值给n
65    Call mxCopyPtrToReal8(mxGetPr(InVar(7)),co,1)!将第7个整数变量赋值给co cut_off_distance
66     
67    NBP = Dim_P(1)
68    NTP = Dim_P(3)
69    NSP = Dim_P(4)
70    NBF = Dim_F(1)
71    NTF = Dim_F(3)
72    NSF = Dim_F(4)
73 
74 
75    OutVar(1)=mxCreateNumericArray(PD,Dim_P,REAL*8,0)!给返回参数分配内存
76 
77    Call mxCopyReal8ToPtr(Uind,mxGetPr(OutVar(1)),Dim_P(1)*Dim_P(2)*Dim_P(3)*Dim_P(4))!将返回参数赋值给分配的内存
78     
79    DeAllocate(P,F1,F2,Gamma,rc)!释放临时分配的内存
80 
81    Return
82 
83End SubRoutine

AsIlll · 发表于 2023-11-4 12:51:51

li913 发表于 2023-11-3 20:13
你把fortran里面所有的计算语句都删除，看看有没有问题。这样可以确定到底是计算问题还是接口问题。 ...

好的，谢谢您回复，我试试

li913 · 发表于 2023-11-3 20:13:47

你把fortran里面所有的计算语句都删除，看看有没有问题。这样可以确定到底是计算问题还是接口问题。

[混编] 基于mex文件的Matlab和Fortran混合编程问题