CppAD: A C++ Algorithmic Differentiation Package
20130918
|
00001 /* $Id$ */ 00002 /* -------------------------------------------------------------------------- 00003 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-14 Bradley M. Bell 00004 00005 CppAD is distributed under multiple licenses. This distribution is under 00006 the terms of the 00007 Eclipse Public License Version 1.0. 00008 00009 A copy of this license is included in the COPYING file of this distribution. 00010 Please visit http://www.coin-or.org/CppAD/ for information on other licenses. 00011 -------------------------------------------------------------------------- */ 00012 00013 # include <vector> 00014 # include <cppad/vector.hpp> 00015 # include <ColPack/ColPackHeaders.h> 00016 00017 namespace CppAD { 00018 void cppad_colpack_general( 00019 CppAD::vector<size_t>& color , 00020 size_t m , 00021 size_t n , 00022 const CppAD::vector<unsigned int*>& adolc_pattern ) 00023 { size_t i, k; 00024 CPPAD_ASSERT_UNKNOWN( adolc_pattern.size() == m ); 00025 00026 // Use adolc sparsity pattern to create corresponding bipartite graph 00027 ColPack::BipartiteGraphPartialColoringInterface graph( 00028 SRC_MEM_ADOLC, 00029 adolc_pattern.data(), 00030 m, 00031 n 00032 ); 00033 00034 // row ordered Partial-Distance-Two-Coloring of the bipartite graph 00035 graph.PartialDistanceTwoColoring( 00036 "SMALLEST_LAST", "ROW_PARTIAL_DISTANCE_TWO" 00037 ); 00038 00039 // Use coloring information to create seed matrix 00040 int n_seed_row; 00041 int n_seed_col; 00042 double** seed_matrix = graph.GetSeedMatrix(&n_seed_row, &n_seed_col); 00043 CPPAD_ASSERT_UNKNOWN( size_t(n_seed_col) == m ); 00044 00045 // now return coloring in format required by CppAD 00046 for(i = 0; i < m; i++) 00047 color[i] = m; 00048 for(k = 0; k < size_t(n_seed_row); k++) 00049 { for(i = 0; i < m; i++) 00050 { if( seed_matrix[k][i] != 0.0 ) 00051 { // check that no row appears twice in the coloring 00052 CPPAD_ASSERT_UNKNOWN( color[i] == m ); 00053 color[i] = k; 00054 } 00055 } 00056 } 00057 # ifndef NDEBUG 00058 // check that all non-zero rows appear in the coloring 00059 for(i = 0; i < m; i++) 00060 CPPAD_ASSERT_UNKNOWN(color[i] < m || adolc_pattern[i][0] == 0); 00061 00062 // check that no rows with the same color have overlapping entries 00063 CppAD::vector<bool> found(n); 00064 for(k = 0; k < size_t(n_seed_row); k++) 00065 { size_t j, ell; 00066 for(j = 0; j < n; j++) 00067 found[j] = false; 00068 for(i = 0; i < m; i++) if( color[i] == k ) 00069 { for(ell = 0; ell < adolc_pattern[i][0]; ell++) 00070 { j = adolc_pattern[i][1 + ell]; 00071 CPPAD_ASSERT_UNKNOWN( ! found[j] ); 00072 found[j] = true; 00073 } 00074 } 00075 } 00076 # endif 00077 return; 00078 } 00079 // The following routine is not yet used or tested. 00080 void cppad_colpack_symmetric( 00081 CppAD::vector<size_t>& color , 00082 size_t n , 00083 const CppAD::vector<unsigned int*>& adolc_pattern ) 00084 { size_t i, k; 00085 CPPAD_ASSERT_UNKNOWN( adolc_pattern.size() == n ); 00086 00087 // Use adolc sparsity pattern to create corresponding bipartite graph 00088 ColPack::GraphColoringInterface graph( 00089 SRC_MEM_ADOLC, 00090 adolc_pattern.data(), 00091 n 00092 ); 00093 00094 // Color the graph with the speciied ordering 00095 // graph.Coloring("SMALLEST_LAST", "STAR") is slower in adolc testing 00096 graph.Coloring("SMALLEST_LAST", "ACYCLIC_FOR_INDIRECT_RECOVERY"); 00097 00098 // Use coloring information to create seed matrix 00099 int n_seed_row; 00100 int n_seed_col; 00101 double** seed_matrix = graph.GetSeedMatrix(&n_seed_row, &n_seed_col); 00102 CPPAD_ASSERT_UNKNOWN( size_t(n_seed_col) == n ); 00103 00104 // now return coloring in format required by CppAD 00105 for(i = 0; i < n; i++) 00106 color[i] = n; 00107 for(k = 0; k < size_t(n_seed_row); k++) 00108 { for(i = 0; i < n; i++) 00109 { if( seed_matrix[k][i] != 0.0 ) 00110 { CPPAD_ASSERT_UNKNOWN( color[i] == n ); 00111 color[i] = k; 00112 } 00113 } 00114 } 00115 # ifndef NDEBUG 00116 for(i = 0; i < n; i++) 00117 CPPAD_ASSERT_UNKNOWN(color[i] < n || adolc_pattern[i][0] == 0); 00118 00119 // The coloring above will probably fail this test. 00120 // Check that no rows with the same color have overlapping entries: 00121 CppAD::vector<bool> found(n); 00122 for(k = 0; k < size_t(n_seed_row); k++) 00123 { size_t j, ell; 00124 for(j = 0; j < n; j++) 00125 found[j] = false; 00126 for(i = 0; i < n; i++) if( color[i] == k ) 00127 { for(ell = 0; ell < adolc_pattern[i][0]; ell++) 00128 { j = adolc_pattern[i][1 + ell]; 00129 CPPAD_ASSERT_UNKNOWN( ! found[j] ); 00130 found[j] = true; 00131 } 00132 } 00133 } 00134 # endif 00135 return; 00136 } 00137 }