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