source: trunk/cppad_lib/cppad_colpack.cpp @ 3845

Last change on this file since 3845 was 3845, checked in by bradbell, 3 years ago

merge to branch: trunk
from repository: https://github.com/coin-or/CppAD
start hash code: 6b0d372ce8350b57b0c4fb3ad6e815316aa28b31
end hash code: 60923b98407ffc46249875d39ac667a2a5bdb8c9

commit 60923b98407ffc46249875d39ac667a2a5bdb8c9
Author: Brad Bell <bradbell@…>
Date: Fri Nov 18 18:33:53 2016 -0700

Advance to cppad-20161118.

commit 9dbc363f75937f64686b315320e1884bd39e89f3
Author: Brad Bell <bradbell@…>
Date: Fri Nov 18 09:02:02 2016 -0700

local branch:
declare_ad.hpp: remove some unecessary declarations.

commit bb55f01c0a486b9577a4a5c71c35a5a16b096c4c
Author: Brad Bell <bradbell@…>
Date: Fri Nov 18 08:05:45 2016 -0700

local branch:
Split hash_code and independent into core and local parts.

commit cade0ae629a5035b8a435a0d8f1d8534ae3c1293
Author: Brad Bell <bradbell@…>
Date: Fri Nov 18 07:01:47 2016 -0700

local branch:
Move op_code.hpp identifiers into local namespace.

commit 67ddec52a43933fb19e3fb85907056fb8493d5e0
Author: Brad Bell <bradbell@…>
Date: Fri Nov 18 04:45:24 2016 -0700

local branch:

  1. Move tape_link from local to core directory.
  2. Put following routines in local namespace: one_element_std_set, two_element_std_set, set_get_in_parallel.


check_makefile.sh: use ls_files.sh to restrict to source code makefile.am files.

commit 8da16c2145003e0fb8f133f64764d597296e486c
Author: Brad Bell <bradbell@…>
Date: Fri Nov 18 03:44:26 2016 -0700

local branch:
Move pod_vector, sparse_list, sparse_pack into local namespace.

commit fecacda7359cd23fda184a64b6afccb0785bc5fa
Author: Brad Bell <bradbell@…>
Date: Thu Nov 17 18:00:25 2016 -0700

local branch:
Move some more files from local to core.

commit 12866a5cf24f1b0a08c2be6d4afdf24b913c0796
Author: Brad Bell <bradbell@…>
Date: Thu Nov 17 07:02:54 2016 -0700

local branch:
Move add.hpp and add_eq.hpp out of local namespace directory.

commit ea43f73fcaa657a1ab6fe55f6ee473ffb6480837
Author: Brad Bell <bradbell@…>
Date: Thu Nov 17 06:52:07 2016 -0700

local branch:
Move remaining files that do not contain local namespace routines from local to core.

commit 9793adb3baf44b975a79d00d0e1265a9ed9a3393
Author: Brad Bell <bradbell@…>
Date: Thu Nov 17 05:43:14 2016 -0700

local branch:
Move more routines that are not in local namespace from local to core.

commit 324694ff288e38861cead02f785c7e736aa0d322
Author: Brad Bell <bradbell@…>
Date: Thu Nov 17 04:45:02 2016 -0700

local branch:
Move more routines not in local namespace from local to core.

commit c25e66d13a4f49456e1ed65fd7ad0b6b1d88c825
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 21:13:09 2016 -0700

local branch:
Move some more functions that are not in local namespace out of local directory.

commit 9245dd1e8621c31290f13bb01f7f52c1dab1bc7f
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 13:49:01 2016 -0700

local branch:
Move local/ad*.hpp files that declare functions in CppAD namespace to core.

commit 9dd0f3a4f4d584e3df6ae0c7b7a0352a4eabb6eb
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 13:06:05 2016 -0700

local branch:
Move all the tape operators into the local namespace.

commit 93858e16ae715f8d870192ed94dc46948c992837
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 12:53:36 2016 -0700

local branch:
Move all the sweep functions into local namespace.

commit 3e3025bd7a747252354fd49e0768e27ec59c081b
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 11:34:42 2016 -0700

local branch:
Move coloring routines into local namespace.

commit ca06cfe07cfed44b66f7a5abd47236fbb2db2e72
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 09:51:09 2016 -0700

local branch:
Move CppAD::ADTape -> CppAD::local::ADTape.

commit d3c474afcd6f87cd7130d2487a6cb5171f0a7d31
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 07:26:09 2016 -0700

check_replace.sh: remove trailing white space.
optimize.hpp: move to CppAD::local::optimize namespace.

commit 27fdf442bd0bd76d705bcbb73240da7c909eb9f0
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 07:02:32 2016 -0700

local branch:
Move cppad/local/optimize.hpp -> cppad/local/optimize/optimize.hpp

commit bca9607442ea768a064853dd524e63b4d695b621
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 06:12:17 2016 -0700

local branch:
Move recorder to local namespace.

commit 66360e9a4d832296ceb610f354e389b1d954ea62
Author: Brad Bell <bradbell@…>
Date: Wed Nov 16 05:54:40 2016 -0700

local branch:
Move player to local namespace.

File size: 6.6 KB
Line 
1// $Id$
2/* --------------------------------------------------------------------------
3CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-16 Bradley M. Bell
4
5CppAD is distributed under multiple licenses. This distribution is under
6the terms of the
7                    Eclipse Public License Version 1.0.
8
9A copy of this license is included in the COPYING file of this distribution.
10Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
11-------------------------------------------------------------------------- */
12# include <vector>
13# include <cppad/utility/vector.hpp>
14# include <cppad/configure.hpp>
15# include <cppad/core/define.hpp>
16# include <cppad/core/cppad_assert.hpp>
17
18# if CPPAD_HAS_COLPACK == 0
19namespace CppAD { namespace local {
20        CPPAD_LIB_EXPORT void this_routine_should_never_get_called(void)
21        {       CPPAD_ASSERT_UNKNOWN(false); }
22}
23# else // CPPAD_HAS_COLPACK
24# include <ColPack/ColPackHeaders.h>
25
26namespace CppAD { namespace local { // BEGIN_CPPAD_LOCAL_NAMESPACE
27
28/*!
29\file cppad_colpack.cpp
30The CppAD interface to the Colpack coloring algorithms.
31*/
32
33/*!
34Determine which rows of a general sparse matrix can be computed together.
35
36\param color
37is a vector with color.size() == m.
38For i = 0 , ... , m-1, color[i] is the color for the corresponding row
39of the matrix. If color[i1]==color[i2], (i1, j1) is in sparsity pattern,
40and (i2, j2) is in sparsity pattern, then j1 is not equal to j2.
41
42\param m
43is the number of rows in the matrix.
44
45\param n
46is the number of columns in the matrix.
47
48\param adolc_pattern
49is a vector with adolc_pattern.size() == m.
50For i = 0 , ... , m-1, and for k = 1, ... ,adolc_pattern[i][0],
51the entry with index (i, adolc_pattern[i][k]) is a non-zero
52in the sparsity pattern for the matrix.
53*/
54// ----------------------------------------------------------------------
55CPPAD_LIB_EXPORT void cppad_colpack_general(
56        CppAD::vector<size_t>&               color         ,
57        size_t                               m             ,
58        size_t                               n             ,
59        const CppAD::vector<unsigned int*>&  adolc_pattern )
60{       size_t i, k;
61        CPPAD_ASSERT_UNKNOWN( adolc_pattern.size() == m );
62        CPPAD_ASSERT_UNKNOWN( color.size() == m );
63
64        // Use adolc sparsity pattern to create corresponding bipartite graph
65        ColPack::BipartiteGraphPartialColoringInterface graph(
66                        SRC_MEM_ADOLC,
67                        adolc_pattern.data(),
68                        m,
69                        n
70        );
71
72        // row ordered Partial-Distance-Two-Coloring of the bipartite graph
73        graph.PartialDistanceTwoColoring(
74                "SMALLEST_LAST", "ROW_PARTIAL_DISTANCE_TWO"
75        );
76
77        // Use coloring information to create seed matrix
78        int n_seed_row;
79        int n_seed_col;
80        double** seed_matrix = graph.GetSeedMatrix(&n_seed_row, &n_seed_col);
81        CPPAD_ASSERT_UNKNOWN( size_t(n_seed_col) == m );
82
83        // now return coloring in format required by CppAD
84        for(i = 0; i < m; i++)
85                color[i] = m;
86        for(k = 0; k < size_t(n_seed_row); k++)
87        {       for(i = 0; i < m; i++)
88                {       if( seed_matrix[k][i] != 0.0 )
89                        {       // check that no row appears twice in the coloring
90                                CPPAD_ASSERT_UNKNOWN( color[i] == m );
91                                color[i] = k;
92                        }
93                }
94        }
95# ifndef NDEBUG
96        // check that all non-zero rows appear in the coloring
97        for(i = 0; i < m; i++)
98                CPPAD_ASSERT_UNKNOWN(color[i] < m || adolc_pattern[i][0] == 0);
99
100        // check that no rows with the same color have overlapping entries
101        CppAD::vector<bool> found(n);
102        for(k = 0; k < size_t(n_seed_row); k++)
103        {       size_t j, ell;
104                for(j = 0; j < n; j++)
105                        found[j] = false;
106                for(i = 0; i < m; i++) if( color[i] == k )
107                {       for(ell = 0; ell < adolc_pattern[i][0]; ell++)
108                        {       j = adolc_pattern[i][1 + ell];
109                                CPPAD_ASSERT_UNKNOWN( ! found[j] );
110                                found[j] = true;
111                        }
112                }
113        }
114# endif
115        return;
116}
117// ----------------------------------------------------------------------
118/*!
119Determine which rows of a symmetrix sparse matrix can be computed together.
120
121\param color
122is a vector with color.size() == m.
123For i = 0 , ... , m-1, color[i] is the color for the corresponding row
124of the matrix. We say that a sparsity pattern entry (i, j) is valid if
125for all i1, such that i1 != i and color[i1]==color[i],
126and all j1, such that (i1, j1) is in sparsity pattern, j1 != j.
127The coloring is chosen so that for all (i, j) in the sparsity pattern;
128either (i, j) or (j, i) is valid (possibly both).
129
130\param m
131is the number of rows (and columns) in the matrix.
132
133\param adolc_pattern
134is a vector with adolc_pattern.size() == m.
135For i = 0 , ... , m-1, and for k = 1, ... ,adolc_pattern[i][0],
136the entry with index (i, adolc_pattern[i][k]) is
137in the sparsity pattern for the symmetric matrix.
138*/
139CPPAD_LIB_EXPORT void cppad_colpack_symmetric(
140        CppAD::vector<size_t>&               color         ,
141        size_t                               m             ,
142        const CppAD::vector<unsigned int*>&  adolc_pattern )
143{       size_t i;
144        CPPAD_ASSERT_UNKNOWN( adolc_pattern.size() == m );
145        CPPAD_ASSERT_UNKNOWN( color.size() == m );
146
147        // Use adolc sparsity pattern to create corresponding bipartite graph
148        ColPack::GraphColoringInterface graph(
149                        SRC_MEM_ADOLC,
150                        adolc_pattern.data(),
151                        m
152        );
153
154        // Use STAR coloring because it has a direct recovery scheme; i.e.,
155        // not necessary to solve equations to extract values.
156        graph.Coloring("SMALLEST_LAST", "STAR");
157
158        // Use coloring information to create seed matrix
159        int n_seed_row;
160        int n_seed_col;
161        double** seed_matrix = graph.GetSeedMatrix(&n_seed_row, &n_seed_col);
162        CPPAD_ASSERT_UNKNOWN( size_t(n_seed_row) == m );
163
164        // now return coloring for each row in format required by CppAD
165        for(i = 0; i < m; i++)
166                color[i] = m;
167        for(i = 0; i < m; i++)
168        {       for(size_t k = 0; k < size_t(n_seed_col); k++)
169                {       if( seed_matrix[i][k] != 0.0 )
170                        {       CPPAD_ASSERT_UNKNOWN( color[i] == m );
171                                color[i] = k;
172                        }
173                }
174        }
175
176# ifndef NDEBUG
177        // check that every entry in the symetric matrix can be direclty recovered
178        size_t i1, i2, j1, j2, k1, k2, nz1, nz2;
179        for(i1 = 0; i1 < m; i1++)
180        {       nz1 = size_t(adolc_pattern[i1][0]);
181                for(k1 = 1; k1 <= nz1; k1++)
182                {       j1 = adolc_pattern[i1][k1];
183
184                        // check of a forward on color[i1] followed by a reverse
185                        // can recover entry (i1, j1)
186                        bool color_i1_ok = true;
187                        for(i2 = 0; i2 < m; i2++) if( i1 != i2 && color[i1] == color[i2] )
188                        {       nz2 = adolc_pattern[i2][0];
189                                for(k2 = 1; k2 <= nz2; k2++)
190                                {       j2 = adolc_pattern[i2][k2];
191                                        color_i1_ok &= (j1 != j2);
192                                }
193                        }
194
195                        // check of a forward on color[j1] followed by a reverse
196                        // can recover entry (j1, i1)
197                        bool color_j1_ok = true;
198                        for(j2 = 0; j2 < m; j2++) if( j1 != j2 && color[j1] == color[j2] )
199                        {       nz2 = adolc_pattern[j2][0];
200                                for(k2 = 1; k2 <= nz2; k2++)
201                                {       i2 = adolc_pattern[j2][k2];
202                                        color_j1_ok &= (i1 != i2);
203                                }
204                        }
205
206                        CPPAD_ASSERT_UNKNOWN( color_i1_ok || color_j1_ok );
207                }
208        }
209# endif
210        return;
211}
212
213} } // END_CPPAD_LOCAL_NAMESPACE
214
215# endif // CPPAD_HAS_COLPACK
Note: See TracBrowser for help on using the repository browser.