Mercurial > hg > octave-nkf
annotate examples/code/mysparse.c @ 20811:a22d8a2eb0e5
fix adaptive strategy in ode solvers.
* script/ode/ode45.m: remove unused option OutputSave
* script/ode/private/integrate_adaptive.m: rewrite algorithm
to be more compatible.
* script/ode/private/runge_kutta_45_dorpri.m: use kahan summation
for time increment.
| author | Carlo de Falco <carlo.defalco@polimi.it> |
|---|---|
| date | Sun, 11 Oct 2015 18:44:58 +0200 |
| parents | c8240a60dd01 |
| children |
| rev | line source |
|---|---|
| 5903 | 1 #include "mex.h" |
| 2 | |
| 3 void | |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
4 mexFunction (int nlhs, mxArray *plhs[], |
|
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
5 int nrhs, const mxArray *prhs[]) |
| 5903 | 6 { |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
7 mwSize m, n, nz; |
| 5903 | 8 mxArray *v; |
| 6686 | 9 mwIndex i; |
| 5903 | 10 double *pr, *pi; |
| 11 double *pr2, *pi2; | |
| 6686 | 12 mwIndex *ir, *jc; |
| 13 mwIndex *ir2, *jc2; | |
|
17289
4d7f95eb8bfe
doc: Miscellaneous small tweaks to documentation for consistency.
Rik <rik@octave.org>
parents:
16867
diff
changeset
|
14 |
| 5903 | 15 if (nrhs != 1 || ! mxIsSparse (prhs[0])) |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
16 mexErrMsgTxt ("ARG1 must be a sparse matrix"); |
| 5903 | 17 |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
18 m = mxGetM (prhs[0]); |
|
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
19 n = mxGetN (prhs[0]); |
|
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
20 nz = mxGetNzmax (prhs[0]); |
|
17289
4d7f95eb8bfe
doc: Miscellaneous small tweaks to documentation for consistency.
Rik <rik@octave.org>
parents:
16867
diff
changeset
|
21 |
| 5903 | 22 if (mxIsComplex (prhs[0])) |
| 23 { | |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
24 mexPrintf ("Matrix is %d-by-%d complex sparse matrix", m, n); |
| 6580 | 25 mexPrintf (" with %d elements\n", nz); |
| 5903 | 26 |
| 27 pr = mxGetPr (prhs[0]); | |
| 28 pi = mxGetPi (prhs[0]); | |
| 29 ir = mxGetIr (prhs[0]); | |
| 30 jc = mxGetJc (prhs[0]); | |
| 31 | |
| 32 i = n; | |
| 33 while (jc[i] == jc[i-1] && i != 0) i--; | |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
34 |
|
19007
9ac2357f19bc
doc: Replace "non-zero" with "nonzero" to match existing usage.
Rik <rik@octave.org>
parents:
17289
diff
changeset
|
35 mexPrintf ("last nonzero element (%d, %d) = (%g, %g)\n", |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
36 ir[nz-1]+ 1, i, pr[nz-1], pi[nz-1]); |
| 5903 | 37 |
| 38 v = mxCreateSparse (m, n, nz, mxCOMPLEX); | |
| 39 pr2 = mxGetPr (v); | |
| 40 pi2 = mxGetPi (v); | |
| 41 ir2 = mxGetIr (v); | |
| 42 jc2 = mxGetJc (v); | |
|
17289
4d7f95eb8bfe
doc: Miscellaneous small tweaks to documentation for consistency.
Rik <rik@octave.org>
parents:
16867
diff
changeset
|
43 |
| 5903 | 44 for (i = 0; i < nz; i++) |
| 6580 | 45 { |
| 46 pr2[i] = 2 * pr[i]; | |
| 47 pi2[i] = 2 * pi[i]; | |
| 48 ir2[i] = ir[i]; | |
| 49 } | |
| 5903 | 50 for (i = 0; i < n + 1; i++) |
| 6580 | 51 jc2[i] = jc[i]; |
| 5903 | 52 |
| 53 if (nlhs > 0) | |
| 6580 | 54 plhs[0] = v; |
| 5903 | 55 } |
| 56 else if (mxIsLogical (prhs[0])) | |
| 57 { | |
|
15352
9a8dbd6b6b20
* mysparse.c (mexFunction): Declare pbr and pbr2 as mxLogical*, not bool*.
John W. Eaton <jwe@octave.org>
parents:
9932
diff
changeset
|
58 mxLogical *pbr, *pbr2; |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
59 mexPrintf ("Matrix is %d-by-%d logical sparse matrix", m, n); |
| 6580 | 60 mexPrintf (" with %d elements\n", nz); |
| 5903 | 61 |
| 62 pbr = mxGetLogicals (prhs[0]); | |
| 63 ir = mxGetIr (prhs[0]); | |
| 64 jc = mxGetJc (prhs[0]); | |
| 65 | |
| 66 i = n; | |
| 67 while (jc[i] == jc[i-1] && i != 0) i--; | |
|
19007
9ac2357f19bc
doc: Replace "non-zero" with "nonzero" to match existing usage.
Rik <rik@octave.org>
parents:
17289
diff
changeset
|
68 mexPrintf ("last nonzero element (%d, %d) = %d\n", |
| 7081 | 69 ir[nz-1]+ 1, i, pbr[nz-1]); |
| 5903 | 70 |
| 71 v = mxCreateSparseLogicalMatrix (m, n, nz); | |
| 72 pbr2 = mxGetLogicals (v); | |
| 73 ir2 = mxGetIr (v); | |
| 74 jc2 = mxGetJc (v); | |
|
17289
4d7f95eb8bfe
doc: Miscellaneous small tweaks to documentation for consistency.
Rik <rik@octave.org>
parents:
16867
diff
changeset
|
75 |
| 5903 | 76 for (i = 0; i < nz; i++) |
| 6580 | 77 { |
| 78 pbr2[i] = pbr[i]; | |
| 79 ir2[i] = ir[i]; | |
| 80 } | |
| 5903 | 81 for (i = 0; i < n + 1; i++) |
| 6580 | 82 jc2[i] = jc[i]; |
| 5903 | 83 |
| 84 if (nlhs > 0) | |
| 6580 | 85 plhs[0] = v; |
| 5903 | 86 } |
| 87 else | |
| 88 { | |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
89 mexPrintf ("Matrix is %d-by-%d real sparse matrix", m, n); |
| 6580 | 90 mexPrintf (" with %d elements\n", nz); |
| 5903 | 91 |
| 92 pr = mxGetPr (prhs[0]); | |
| 93 ir = mxGetIr (prhs[0]); | |
| 94 jc = mxGetJc (prhs[0]); | |
| 95 | |
| 96 i = n; | |
| 97 while (jc[i] == jc[i-1] && i != 0) i--; | |
|
19007
9ac2357f19bc
doc: Replace "non-zero" with "nonzero" to match existing usage.
Rik <rik@octave.org>
parents:
17289
diff
changeset
|
98 mexPrintf ("last nonzero element (%d, %d) = %g\n", |
|
16867
be41c30bcb44
Re-write documentation and all examples of dynamically linked functions.
Rik <rik@octave.org>
parents:
15352
diff
changeset
|
99 ir[nz-1]+ 1, i, pr[nz-1]); |
| 5903 | 100 |
| 101 v = mxCreateSparse (m, n, nz, mxREAL); | |
| 102 pr2 = mxGetPr (v); | |
| 103 ir2 = mxGetIr (v); | |
| 104 jc2 = mxGetJc (v); | |
|
17289
4d7f95eb8bfe
doc: Miscellaneous small tweaks to documentation for consistency.
Rik <rik@octave.org>
parents:
16867
diff
changeset
|
105 |
| 5903 | 106 for (i = 0; i < nz; i++) |
| 6580 | 107 { |
| 108 pr2[i] = 2 * pr[i]; | |
| 109 ir2[i] = ir[i]; | |
| 110 } | |
| 5903 | 111 for (i = 0; i < n + 1; i++) |
| 6580 | 112 jc2[i] = jc[i]; |
| 5903 | 113 |
| 114 if (nlhs > 0) | |
| 6580 | 115 plhs[0] = v; |
| 5903 | 116 } |
| 117 } |