Mercurial > hg > octave-nkf
annotate liboctave/numeric/sparse-base-lu.cc @ 20830:b65888ec820e draft default tip gccjit
dmalcom gcc jit import
author | Stefan Mahr <dac922@gmx.de> |
---|---|
date | Fri, 27 Feb 2015 16:59:36 +0100 |
parents | 4197fc428c7d |
children |
rev | line source |
---|---|
5164 | 1 /* |
2 | |
19898
4197fc428c7d
maint: Update copyright notices for 2015.
John W. Eaton <jwe@octave.org>
parents:
19461
diff
changeset
|
3 Copyright (C) 2004-2015 David Bateman |
11523 | 4 Copyright (C) 1998-2004 Andy Adler |
7016 | 5 |
6 This file is part of Octave. | |
5164 | 7 |
8 Octave is free software; you can redistribute it and/or modify it | |
9 under the terms of the GNU General Public License as published by the | |
7016 | 10 Free Software Foundation; either version 3 of the License, or (at your |
11 option) any later version. | |
5164 | 12 |
13 Octave is distributed in the hope that it will be useful, but WITHOUT | |
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |
16 for more details. | |
17 | |
18 You should have received a copy of the GNU General Public License | |
7016 | 19 along with Octave; see the file COPYING. If not, see |
20 <http://www.gnu.org/licenses/>. | |
5164 | 21 |
22 */ | |
23 | |
24 #ifdef HAVE_CONFIG_H | |
25 #include <config.h> | |
26 #endif | |
27 | |
28 #include "sparse-base-lu.h" | |
29 | |
19461
65554f5847ac
don't include oct-locbuf.h in header files unnecessarily
John W. Eaton <jwe@octave.org>
parents:
18493
diff
changeset
|
30 #include "PermMatrix.h" |
65554f5847ac
don't include oct-locbuf.h in header files unnecessarily
John W. Eaton <jwe@octave.org>
parents:
18493
diff
changeset
|
31 |
5164 | 32 template <class lu_type, class lu_elt_type, class p_type, class p_elt_type> |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
33 lu_type |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
34 sparse_base_lu <lu_type, lu_elt_type, p_type, p_elt_type> :: Y (void) const |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
35 { |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
36 octave_idx_type nr = Lfact.rows (); |
18493
2a45b6b87bee
correct numerical errors in sparse LU factorization (bug #41116).
Michael C. Grant <mcg@cvxr.com>
parents:
17744
diff
changeset
|
37 octave_idx_type nz = Lfact.cols (); |
2a45b6b87bee
correct numerical errors in sparse LU factorization (bug #41116).
Michael C. Grant <mcg@cvxr.com>
parents:
17744
diff
changeset
|
38 octave_idx_type nc = Ufact.cols (); |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
39 |
18493
2a45b6b87bee
correct numerical errors in sparse LU factorization (bug #41116).
Michael C. Grant <mcg@cvxr.com>
parents:
17744
diff
changeset
|
40 lu_type Yout (nr, nc, Lfact.nnz () + Ufact.nnz () - (nr<nz?nr:nz)); |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
41 octave_idx_type ii = 0; |
15018
3d8ace26c5b4
maint: Use Octave coding conventions for cuddled parentheses in liboctave/.
Rik <rik@octave.org>
parents:
14846
diff
changeset
|
42 Yout.xcidx (0) = 0; |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
43 |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
44 for (octave_idx_type j = 0; j < nc; j++) |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
45 { |
15018
3d8ace26c5b4
maint: Use Octave coding conventions for cuddled parentheses in liboctave/.
Rik <rik@octave.org>
parents:
14846
diff
changeset
|
46 for (octave_idx_type i = Ufact.cidx (j); i < Ufact.cidx (j + 1); i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
47 { |
15018
3d8ace26c5b4
maint: Use Octave coding conventions for cuddled parentheses in liboctave/.
Rik <rik@octave.org>
parents:
14846
diff
changeset
|
48 Yout.xridx (ii) = Ufact.ridx (i); |
3d8ace26c5b4
maint: Use Octave coding conventions for cuddled parentheses in liboctave/.
Rik <rik@octave.org>
parents:
14846
diff
changeset
|
49 Yout.xdata (ii++) = Ufact.data (i); |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
50 } |
18493
2a45b6b87bee
correct numerical errors in sparse LU factorization (bug #41116).
Michael C. Grant <mcg@cvxr.com>
parents:
17744
diff
changeset
|
51 if (j < nz) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
52 { |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
53 // Note the +1 skips the 1.0 on the diagonal |
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
54 for (octave_idx_type i = Lfact.cidx (j) + 1; |
15018
3d8ace26c5b4
maint: Use Octave coding conventions for cuddled parentheses in liboctave/.
Rik <rik@octave.org>
parents:
14846
diff
changeset
|
55 i < Lfact.cidx (j +1); i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
56 { |
15018
3d8ace26c5b4
maint: Use Octave coding conventions for cuddled parentheses in liboctave/.
Rik <rik@octave.org>
parents:
14846
diff
changeset
|
57 Yout.xridx (ii) = Lfact.ridx (i); |
3d8ace26c5b4
maint: Use Octave coding conventions for cuddled parentheses in liboctave/.
Rik <rik@octave.org>
parents:
14846
diff
changeset
|
58 Yout.xdata (ii++) = Lfact.data (i); |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
59 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
60 } |
15018
3d8ace26c5b4
maint: Use Octave coding conventions for cuddled parentheses in liboctave/.
Rik <rik@octave.org>
parents:
14846
diff
changeset
|
61 Yout.xcidx (j + 1) = ii; |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
62 } |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
63 |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
64 return Yout; |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
65 } |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
66 |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
67 template <class lu_type, class lu_elt_type, class p_type, class p_elt_type> |
5164 | 68 p_type |
69 sparse_base_lu <lu_type, lu_elt_type, p_type, p_elt_type> :: Pr (void) const | |
70 { | |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
71 |
5275 | 72 octave_idx_type nr = Lfact.rows (); |
5164 | 73 |
74 p_type Pout (nr, nr, nr); | |
75 | |
5275 | 76 for (octave_idx_type i = 0; i < nr; i++) |
5164 | 77 { |
78 Pout.cidx (i) = i; | |
79 Pout.ridx (P (i)) = i; | |
80 Pout.data (i) = 1; | |
81 } | |
82 Pout.cidx (nr) = nr; | |
83 | |
84 return Pout; | |
85 } | |
86 | |
87 template <class lu_type, class lu_elt_type, class p_type, class p_elt_type> | |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
88 ColumnVector |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
89 sparse_base_lu <lu_type, lu_elt_type, p_type, p_elt_type> :: Pr_vec (void) const |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
90 { |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
91 |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
92 octave_idx_type nr = Lfact.rows (); |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
93 |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
94 ColumnVector Pout (nr); |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
95 |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
96 for (octave_idx_type i = 0; i < nr; i++) |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
97 Pout.xelem (i) = static_cast<double> (P(i) + 1); |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
98 |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
99 return Pout; |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
100 } |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
101 |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
102 template <class lu_type, class lu_elt_type, class p_type, class p_elt_type> |
8969
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
103 PermMatrix |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
104 sparse_base_lu <lu_type, lu_elt_type, p_type, p_elt_type> :: Pr_mat (void) const |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
105 { |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
106 return PermMatrix (P, false); |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
107 } |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
108 |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
109 template <class lu_type, class lu_elt_type, class p_type, class p_elt_type> |
5164 | 110 p_type |
111 sparse_base_lu <lu_type, lu_elt_type, p_type, p_elt_type> :: Pc (void) const | |
112 { | |
5275 | 113 octave_idx_type nc = Ufact.cols (); |
5164 | 114 |
115 p_type Pout (nc, nc, nc); | |
116 | |
5275 | 117 for (octave_idx_type i = 0; i < nc; i++) |
5164 | 118 { |
119 Pout.cidx (i) = i; | |
120 Pout.ridx (i) = Q (i); | |
121 Pout.data (i) = 1; | |
122 } | |
123 Pout.cidx (nc) = nc; | |
124 | |
125 return Pout; | |
126 } | |
127 | |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
128 template <class lu_type, class lu_elt_type, class p_type, class p_elt_type> |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
129 ColumnVector |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
130 sparse_base_lu <lu_type, lu_elt_type, p_type, p_elt_type> :: Pc_vec (void) const |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
131 { |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
132 |
8954
97c84c4c2247
Make the column permutation vector in sparse LU cols()-long.
Jason Riedy <jason@acm.org>
parents:
8920
diff
changeset
|
133 octave_idx_type nc = Ufact.cols (); |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
134 |
8954
97c84c4c2247
Make the column permutation vector in sparse LU cols()-long.
Jason Riedy <jason@acm.org>
parents:
8920
diff
changeset
|
135 ColumnVector Pout (nc); |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
136 |
8954
97c84c4c2247
Make the column permutation vector in sparse LU cols()-long.
Jason Riedy <jason@acm.org>
parents:
8920
diff
changeset
|
137 for (octave_idx_type i = 0; i < nc; i++) |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
138 Pout.xelem (i) = static_cast<double> (Q(i) + 1); |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
139 |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
140 return Pout; |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
141 } |
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
142 |
8969
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
143 template <class lu_type, class lu_elt_type, class p_type, class p_elt_type> |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
144 PermMatrix |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
145 sparse_base_lu <lu_type, lu_elt_type, p_type, p_elt_type> :: Pc_mat (void) const |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
146 { |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
147 return PermMatrix (Q, true); |
3ecbc236e2e0
Have sparse LU return permutation matrices rather than sparse matrices.
Jason Riedy <jason@acm.org>
parents:
8954
diff
changeset
|
148 } |