Calling from Matlab a recusive function written in C

Question

Ameer Ahmed 2017-6-14

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/344686-calling-from-matlab-a-recusive-function-written-in-c

编辑： James Tursa 2017-6-15

Hi,

I have a written a recursive function in C to replace nchoosek() in Matlab. When i call this function in Matlab, it gives error if i use output parameter. If i don't use output parameter then it does not return the required answer. I have tried to debug it but of no use.

Code without output parameter is:

#include "math.h"
#include "mex.h"
double test(double n, double k)
{
      if (k ==0){
            return 1;
      }
      else
      {
       return(n*test(n-1,k-1))/k;
      }
  }
void mexFunction(int nlhs, mxArray *plhs[], int nrhs,
        const mxArray *prhs[])
{
    double y, n, k;
      /* Create matrix for the return argument. */
      plhs[0] = mxCreateDoubleMatrix(1,1, mxREAL);
      /* Assign pointers to each input and output. */
      n = mxGetScalar(prhs[0]);
      k = mxGetScalar(prhs[1]);
       /* Call the test subroutine. */
      test(n,k);
  }

Code with output parameter is:

#include "math.h"
#include "mex.h"
double test(double *y, double n, double k)
{
      if (k ==0){
            return 1;
      }
      else
      {
       *y=(n*test(n-1,k-1))/k;
      }
  }
void mexFunction(int nlhs, mxArray *plhs[], int nrhs,const mxArray *prhs[])
{
      double *y, n, k;
      /* Create matrix for the return argument. */
      plhs[0] = mxCreateDoubleMatrix(1,1, mxREAL);
      /* Assign pointers to each input and output. */
      n = mxGetScalar(prhs[0]);
      k = mxGetScalar(prhs[1]);
      y = mxGetPr(plhs[0]);
      /* Call the test subroutine. */
      test(y,n,k);
  }

I have to implement this function since nchoosek(n,k) gives warning if i compute it for higher values of n & k.

I will be highly obliged if any one is kind enough to help me out.

Ameer

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

James Tursa 2017-6-15

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/344686-calling-from-matlab-a-recusive-function-written-in-c#answer_270749

编辑：James Tursa 2017-6-15

在 MATLAB Online 中打开

In addition to what Geoff has already written, in your "code without output parameter" you are not setting the value of the return variable. E.g., these lines:

      plhs[0] = mxCreateDoubleMatrix(1,1, mxREAL);
      /* Assign pointers to each input and output. */
      n = mxGetScalar(prhs[0]);
      k = mxGetScalar(prhs[1]);
       /* Call the test subroutine. */
      test(n,k); // <-- Calculates a value then throws it away, doesn't get into plhs[0]

should look something like this instead:

      /* Assign pointers to each input and output. */
      n = mxGetScalar(prhs[0]);
      k = mxGetScalar(prhs[1]);
       /* Call the test subroutine. */
      plhs[0] = mxCreateDoubleScalar(test(n,k)); // <-- Calculates value and gets it into plhs[0]

As a side note, it would be much better to code this up in a loop inside mexFunction rather than chewing up a lot of stack space and overhead on all of those recursive function calls. The loop will be cleaner, use less resources, and certainly run faster.

Also, if you are getting warnings or errors from nchoosek because the numbers you are using are too large, then a mex routine isn't likely to fix that. You may need to rethink how you are solving your problem.

3 个评论
显示 1更早的评论隐藏 1更早的评论

Ameer Ahmed 2017-6-15

Thanks James

I just added

plhs[0] = mxCreateDoubleScalar(test(n,k)); into my code and i got the result. Thanks for your advice.

James Tursa 2017-6-15

编辑：James Tursa 2017-6-15

在 MATLAB Online 中打开

E.g., here is a bare bones loop implementation. This is for example purposes only. At the very least, much more argument checking would have to be inserted in order for it to be "production quality" code.

/* File:  nchoosekx.c                                                    */
/* Programmer: James Tursa                                               */
/* Bare bones C-mex replacement of nchoosek (for example purposes)       */
/* Includes ------------------------------------------------------------ */
#include "mex.h"
/* Gateway ------------------------------------------------------------- */
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
    double i, n, k, m, result;
    // Need to add more input checking here (complex, size, etc)
    if( nrhs == 2 && mxIsNumeric(prhs[0]) && mxIsNumeric(prhs[1]) ) {
        n = mxGetScalar(prhs[0]);
        k = mxGetScalar(prhs[1]);
        // Add checks here for valid n and k values (integer, k<=n, etc)
        result = 1;
        m = k < n-k ? k : n-k;
        for( i=1; i<=m; i++ ) {
            result = result * (n--) / i;
        }
        // Might need to add rounding here
        plhs[0] = mxCreateDoubleScalar(result);
    } else {
        mexErrMsgTxt("nchoosekx: Invalid inputs");
    }
}

Just spot checking results against a vpa implementation, this loop implementation does appear to give slightly more accurate results than the MATLAB nchoosek function when the inputs are too large to generate an exact double result. The mex loop implementation also runs faster. Maybe you can morph this code into what you are really wanting to do.

请先登录，再进行评论。

Answer 2

Geoff Hayes 2017-6-14

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/344686-calling-from-matlab-a-recusive-function-written-in-c#answer_270731

在 MATLAB Online 中打开

Ameer - you are getting an error with the output parameter code because of your recursive call to test

*y=(n*test(n-1,k-1))/k;

So two parameters are being passed in to test but three are required

double test(double *y, double n, double k)

so you will need to either pass in the third parameter to test or re-work this method in order to return a value instead of updating the input y

double test(double n, double k)
{
    if (k == 0)
    {
        return 1;
    }
    else
    {
        return (n*test(n-1,k-1))/k;
    }
}

Note also that you will want to re-think how you are calculating n choose k. Remember, that this is really

n!/((n-k)!*(k!))

rather than the

n!/k!

that you seem to have.

2 个评论
显示无隐藏无

James Tursa 2017-6-15

test(n,k) =

n! / ( (n-k)! * (k!) ) =

n * (n-1)! / ( (n-k)! * k * (k-1)! ) =

n * (n-1)! / ( ( (n-1) - (k-1) )! * (k-1)! ) / k =

n * test(n-1,k-1) / k

So the recursive call in the above code appears to be correct.

Geoff Hayes 2017-6-15

hmmm...I was getting a different answer when comparing it with the MATLAB nchoosek. But now I'm not - weird!

请先登录，再进行评论。

Answer 3

Ameer Ahmed 2017-6-15

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/344686-calling-from-matlab-a-recusive-function-written-in-c#answer_270760

在 MATLAB Online 中打开

Thanks James & Geoff.

I have tried doing it with two parameters as

double test(double n, double k) {

      if (k ==0){
            return 1;
      }
      else
      {
       return(n*test(n-1,k-1))/k;
      }
  }

But it does not show the result. Just returns 0. I have implemented it in C without Mex Function. This works fine and returns the required answer. I don't know where is the mistake.

Actually, i want to implement a complex expression in C using Mex function which could replace three nchoosek(n,k) functions

nchoosek(n,k)*nchoosek(n,k)/nchoosek(n,k)

where each nchoosek() has different n and k.

I will be really helping if some one help me in implementing it.