1 #include "cs_mex.h"
2 /* x(p) = b */
mexFunction(int nargout,mxArray * pargout[],int nargin,const mxArray * pargin[])3 void mexFunction
4 (
5 int nargout,
6 mxArray *pargout [ ],
7 int nargin,
8 const mxArray *pargin [ ]
9 )
10 {
11 CS_INT n, k, *p ;
12 double *xx ;
13 if (nargout > 1 || nargin != 2)
14 {
15 mexErrMsgTxt ("Usage: x = cs_ipvec(b,p)") ;
16 }
17 n = mxGetNumberOfElements (pargin [0]) ;
18 if (n != mxGetNumberOfElements (pargin [1]))
19 {
20 mexErrMsgTxt ("b or p wrong size") ;
21 }
22
23 xx = mxGetPr (pargin [1]) ;
24 p = cs_dl_malloc (n, sizeof (CS_INT)) ;
25 for (k = 0 ; k < n ; k++) p [k] = xx [k] - 1 ;
26
27 if (mxIsComplex (pargin [0]))
28 {
29 #ifdef NCOMPLEX
30 mexErrMsgTxt ("complex case not supported") ;
31 #else
32 cs_complex_t *x, *b ;
33 b = cs_cl_mex_get_double (n, pargin [0]) ;
34 x = cs_dl_malloc (n, sizeof (cs_complex_t)) ;
35 cs_cl_ipvec (p, b, x, n) ;
36 pargout [0] = cs_cl_mex_put_double (n, x) ;
37 cs_free (b) ; /* free copy of complex values */
38 #endif
39 }
40 else
41 {
42 double *x, *b ;
43 b = cs_dl_mex_get_double (n, pargin [0]) ;
44 pargout [0] = mxCreateDoubleMatrix (n, 1, mxREAL) ;
45 x = mxGetPr (pargout [0]) ;
46 cs_dl_ipvec (p, b, x, n) ;
47 }
48 cs_free (p) ;
49 }
50