这个c文件的mex函数怎么写

2024-04-14

函数是cyclic.c.

  void cyclic(float a[], float b[], float c[], float alpha, float beta,
    float r[], float x[], unsigned long n)

   // Solves for a vector x[1..n] the “cyclic” set of linear equations. a,
    //b, c, and r are input vectors, all dimensioned as [1..n], while alpha and beta are        //the corner
  //  entries in the matrix.

我对 Matlab 和 C 之间的接口很陌生。而且我已经好几年没有使用 C 了。 昨晚,我完成并编译。最后一件事是调用它。

#include "mex.h" 


#include "nrutil.h"

#define FREE_ARG char*
#define NR_END 1


#include <stdio.h>
#include <stddef.h>
#include <stdlib.h>
#define NR_END 1
#define FREE_ARG char*
void nrerror(char error_text[])
/* Numerical Recipes standard error handler */
{fprintf(stderr,"Numerical Recipes run-time error...\n");
fprintf(stderr,"%s\n",error_text);
fprintf(stderr,"...now exiting to system...\n");
exit(1);
}
float *vector(long nl, long nh)
/* allocate a float vector with subscript range v[nl..nh] */
{
float *v;
v=(float *)malloc((size_t) ((nh-nl+1+NR_END)*sizeof(float)));
if (!v) nrerror("allocation failure in vector()");
return v-nl+NR_END;
}

void free_vector(float *v, long nl, long nh)
/* free a float vector allocated with vector() */
{
free((FREE_ARG) (v+nl-NR_END));
}

void tridag(float a[], float b[], float c[], float r[], float u[],
unsigned long n)
{
unsigned long j;
float bet,*gam;
gam=vector(1,n); 
if (b[1] == 0.0) nrerror("Error 1 in tridag");

u[1]=r[1]/(bet=b[1]);
for (j=2;j<=n;j++) {
gam[j]=c[j-1]/bet;
bet=b[j]-a[j]*gam[j];
if (bet == 0.0) nrerror("Error 2 in tridag"); 
        u[j]=(r[j]-a[j]*u[j-1])/bet; 
}
for (j=(n-1);j>=1;j--)
u[j] -= gam[j+1]*u[j+1]; 
free_vector(gam,1,n);
}


void cyclic(float a[], float b[], float c[], float alpha, float beta,
float r[], float x[], unsigned long n)

{
void tridag(float a[], float b[], float c[], float r[], float u[],
unsigned long n);
unsigned long i;
float fact,gamma,*bb,*u,*z;
if (n <= 2) nrerror("n too small in cyclic");
bb=vector(1,n);
u=vector(1,n);
z=vector(1,n);
gamma = -b[1]; //Avoid subtraction error in forming bb[1].
bb[1]=b[1]-gamma; //Set up the diagonal of the modified tridiagonal
bb[n]=b[n]-alpha*beta/gamma; //system.
for (i=2;i<n;i++) bb[i]=b[i];
tridag(a,bb,c,r,x,n);// Solve A · x = r.
u[1]=gamma;// Set up the vector u.
u[n]=alpha;
for (i=2;i<n;i++) u[i]=0.0;
tridag(a,bb,c,u,z,n);// Solve A · z = u.
fact=(x[1]+beta*x[n]/gamma)/ //Form v · x/(1 + v · z).
(1.0+z[1]+beta*z[n]/gamma);
for (i=1;i<=n;i++) x[i] -= fact*z[i]; //Nowget the solution vector x.
free_vector(z,1,n);
free_vector(u,1,n);
free_vector(bb,1,n);
}


void mexFunction(int nlhs, mxArray *plhs[], 
        int nrhs, const mxArray *prhs[]) 
    { 

   float *a,*b,*c,*x,*r;
   float alpha,beta;
unsigned long n = (unsigned long) mxGetScalar(prhs[6]);
 //   a=mxGetPr(prhs[0]); 
  //  b=mxGetPr(prhs[1]); 
  //  c=mxGetPr(prhs[2]);
  //  r=mxGetPr(prhs[5]);
    a = (float*) mxGetData(prhs[0]);
    b = (float*) mxGetData(prhs[1]);
    c = (float*) mxGetData(prhs[2]);
    r = (float*) mxGetData(prhs[5]);
   // alpha=*(mxGetPr(prhs[3]));
   // beta=*(mxGetPr(prhs[4]));
    alpha = (float) mxGetScalar(prhs[3]);
 beta = (float) mxGetScalar(prhs[4]);

    plhs[0]= mxCreateDoubleMatrix(n, 1, mxREAL);
    x = mxGetPr(plhs[0]);
    mexPrintf("%f  ",alpha); 
    mexPrintf("\n"); 
    mexPrintf("%f  ",beta); 
    mexPrintf("\n"); 
    mexPrintf("%d  ",n);
    mexPrintf("\n"); 

    cyclic(a,b,c, alpha, beta,r,x,n) ;
    mexPrintf("%d  ",n);
    mexPrintf("\n"); 
    } 

最后我成功编译了cyclic(a,b,c, alpha, beta,r,x,n) ;。但答案是不对的。我认为这是因为r是一个虚向量。所以我的问题是我应该如何在之间转换 rC and Matlab?


C 函数cyclic期望数组floats, but mexFunction正在通过一个double*。不改变cyclic.c,您有两个选择:

  • 将数据转换为single在 MATLAB 中并得到float* with mxGetData.

    In mexFunction:

    float *a = (float*) mxGetData(prhs[0]);
    

    在 MATLAB 中:

    mexFunction(single(a),...)
    
  • 转换(复制,而不是强制转换!)中的数据mexFunction.

    In mexFunction,分配新的float数组,并复制双输入数组中的每个元素(mxGetPr(prhs[0]))进入临时float array.

    Call mexFunction与正常的doubleMATLAB 中的数组。

前者可能更容易做到。

在任何情况下,您都不应该简单地投射指针,而不是您打算这样做。


另外,标量alpha, beta and n需要从中读取prhs作为标量并传递给cyclic作为标量。在mexFunction, use:

float alpha = (float) mxGetScalar(prhs[...]);
float beta = (float) mxGetScalar(prhs[...]);
unsigned long n = (unsigned long) mxGetScalar(prhs[...]);

你已经完全忘记了c and r in mexFunction.

本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

这个c文件的mex函数怎么写 的相关文章

随机推荐