7019
|
1 /* |
|
2 |
|
3 Copyright (C) 2006, 2007 David Bateman |
|
4 |
|
5 This file is part of Octave. |
|
6 |
|
7 Octave is free software; you can redistribute it and/or modify it |
|
8 under the terms of the GNU General Public License as published by the |
|
9 Free Software Foundation; either version 3 of the License, or (at your |
|
10 option) any later version. |
|
11 |
|
12 Octave is distributed in the hope that it will be useful, but WITHOUT |
|
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
|
14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
|
15 for more details. |
|
16 |
|
17 You should have received a copy of the GNU General Public License |
|
18 along with Octave; see the file COPYING. If not, see |
|
19 <http://www.gnu.org/licenses/>. |
|
20 |
|
21 */ |
|
22 |
5903
|
23 #include "mex.h" |
|
24 |
|
25 void |
|
26 mexFunction (int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) |
|
27 { |
6686
|
28 mwSize n, m, nz; |
5903
|
29 mxArray *v; |
6686
|
30 mwIndex i; |
5903
|
31 double *pr, *pi; |
|
32 double *pr2, *pi2; |
6686
|
33 mwIndex *ir, *jc; |
|
34 mwIndex *ir2, *jc2; |
5903
|
35 |
|
36 if (nrhs != 1 || ! mxIsSparse (prhs[0])) |
|
37 mexErrMsgTxt ("expects sparse matrix"); |
|
38 |
|
39 m = mxGetM (prhs [0]); |
|
40 n = mxGetN (prhs [0]); |
|
41 nz = mxGetNzmax (prhs [0]); |
|
42 |
|
43 if (mxIsComplex (prhs[0])) |
|
44 { |
6580
|
45 mexPrintf ("Matrix is %d-by-%d complex sparse matrix", m, n); |
|
46 mexPrintf (" with %d elements\n", nz); |
5903
|
47 |
|
48 pr = mxGetPr (prhs[0]); |
|
49 pi = mxGetPi (prhs[0]); |
|
50 ir = mxGetIr (prhs[0]); |
|
51 jc = mxGetJc (prhs[0]); |
|
52 |
|
53 i = n; |
|
54 while (jc[i] == jc[i-1] && i != 0) i--; |
6580
|
55 mexPrintf ("last non-zero element (%d, %d) = (%g, %g)\n", |
|
56 ir[nz-1]+ 1, i, pr[nz-1], pi[nz-1]); |
5903
|
57 |
|
58 v = mxCreateSparse (m, n, nz, mxCOMPLEX); |
|
59 pr2 = mxGetPr (v); |
|
60 pi2 = mxGetPi (v); |
|
61 ir2 = mxGetIr (v); |
|
62 jc2 = mxGetJc (v); |
|
63 |
|
64 for (i = 0; i < nz; i++) |
6580
|
65 { |
|
66 pr2[i] = 2 * pr[i]; |
|
67 pi2[i] = 2 * pi[i]; |
|
68 ir2[i] = ir[i]; |
|
69 } |
5903
|
70 for (i = 0; i < n + 1; i++) |
6580
|
71 jc2[i] = jc[i]; |
5903
|
72 |
|
73 if (nlhs > 0) |
6580
|
74 plhs[0] = v; |
5903
|
75 } |
|
76 else if (mxIsLogical (prhs[0])) |
|
77 { |
|
78 bool *pbr, *pbr2; |
6580
|
79 mexPrintf ("Matrix is %d-by-%d logical sparse matrix", m, n); |
|
80 mexPrintf (" with %d elements\n", nz); |
5903
|
81 |
|
82 pbr = mxGetLogicals (prhs[0]); |
|
83 ir = mxGetIr (prhs[0]); |
|
84 jc = mxGetJc (prhs[0]); |
|
85 |
|
86 i = n; |
|
87 while (jc[i] == jc[i-1] && i != 0) i--; |
|
88 mexPrintf ("last non-zero element (%d, %d) = %d\n", ir[nz-1]+ 1, |
6580
|
89 i, pbr[nz-1]); |
5903
|
90 |
|
91 v = mxCreateSparseLogicalMatrix (m, n, nz); |
|
92 pbr2 = mxGetLogicals (v); |
|
93 ir2 = mxGetIr (v); |
|
94 jc2 = mxGetJc (v); |
|
95 |
|
96 for (i = 0; i < nz; i++) |
6580
|
97 { |
|
98 pbr2[i] = pbr[i]; |
|
99 ir2[i] = ir[i]; |
|
100 } |
5903
|
101 for (i = 0; i < n + 1; i++) |
6580
|
102 jc2[i] = jc[i]; |
5903
|
103 |
|
104 if (nlhs > 0) |
6580
|
105 plhs[0] = v; |
5903
|
106 } |
|
107 else |
|
108 { |
6580
|
109 mexPrintf ("Matrix is %d-by-%d real sparse matrix", m, n); |
|
110 mexPrintf (" with %d elements\n", nz); |
5903
|
111 |
|
112 pr = mxGetPr (prhs[0]); |
|
113 ir = mxGetIr (prhs[0]); |
|
114 jc = mxGetJc (prhs[0]); |
|
115 |
|
116 i = n; |
|
117 while (jc[i] == jc[i-1] && i != 0) i--; |
|
118 mexPrintf ("last non-zero element (%d, %d) = %g\n", ir[nz-1]+ 1, |
6580
|
119 i, pr[nz-1]); |
5903
|
120 |
|
121 v = mxCreateSparse (m, n, nz, mxREAL); |
|
122 pr2 = mxGetPr (v); |
|
123 ir2 = mxGetIr (v); |
|
124 jc2 = mxGetJc (v); |
|
125 |
|
126 for (i = 0; i < nz; i++) |
6580
|
127 { |
|
128 pr2[i] = 2 * pr[i]; |
|
129 ir2[i] = ir[i]; |
|
130 } |
5903
|
131 for (i = 0; i < n + 1; i++) |
6580
|
132 jc2[i] = jc[i]; |
5903
|
133 |
|
134 if (nlhs > 0) |
6580
|
135 plhs[0] = v; |
5903
|
136 } |
|
137 } |