CppAD: A C++ Algorithmic Differentiation Package  20130918
cppad_colpack.cpp
Go to the documentation of this file.
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 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines