![]() |
Eigen
3.3.3
|
00001 // This file is part of Eigen, a lightweight C++ template library 00002 // for linear algebra. 00003 // 00004 // Copyright (C) 2012 Désiré Nuentsa-Wakam <desire.nuentsa_wakam@inria.fr> 00005 // 00006 // This Source Code Form is subject to the terms of the Mozilla 00007 // Public License v. 2.0. If a copy of the MPL was not distributed 00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00009 00010 /* This file is a modified version of heap_relax_snode.c file in SuperLU 00011 * -- SuperLU routine (version 3.0) -- 00012 * Univ. of California Berkeley, Xerox Palo Alto Research Center, 00013 * and Lawrence Berkeley National Lab. 00014 * October 15, 2003 00015 * 00016 * Copyright (c) 1994 by Xerox Corporation. All rights reserved. 00017 * 00018 * THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY 00019 * EXPRESSED OR IMPLIED. ANY USE IS AT YOUR OWN RISK. 00020 * 00021 * Permission is hereby granted to use or copy this program for any 00022 * purpose, provided the above notices are retained on all copies. 00023 * Permission to modify the code and to distribute modified code is 00024 * granted, provided the above notices are retained, and a notice that 00025 * the code was modified is included with the above copyright notice. 00026 */ 00027 00028 #ifndef SPARSELU_HEAP_RELAX_SNODE_H 00029 #define SPARSELU_HEAP_RELAX_SNODE_H 00030 00031 namespace Eigen { 00032 namespace internal { 00033 00045 template <typename Scalar, typename StorageIndex> 00046 void SparseLUImpl<Scalar,StorageIndex>::heap_relax_snode (const Index n, IndexVector& et, const Index relax_columns, IndexVector& descendants, IndexVector& relax_end) 00047 { 00048 00049 // The etree may not be postordered, but its heap ordered 00050 IndexVector post; 00051 internal::treePostorder(StorageIndex(n), et, post); // Post order etree 00052 IndexVector inv_post(n+1); 00053 for (StorageIndex i = 0; i < n+1; ++i) inv_post(post(i)) = i; // inv_post = post.inverse()??? 00054 00055 // Renumber etree in postorder 00056 IndexVector iwork(n); 00057 IndexVector et_save(n+1); 00058 for (Index i = 0; i < n; ++i) 00059 { 00060 iwork(post(i)) = post(et(i)); 00061 } 00062 et_save = et; // Save the original etree 00063 et = iwork; 00064 00065 // compute the number of descendants of each node in the etree 00066 relax_end.setConstant(emptyIdxLU); 00067 Index j, parent; 00068 descendants.setZero(); 00069 for (j = 0; j < n; j++) 00070 { 00071 parent = et(j); 00072 if (parent != n) // not the dummy root 00073 descendants(parent) += descendants(j) + 1; 00074 } 00075 // Identify the relaxed supernodes by postorder traversal of the etree 00076 Index snode_start; // beginning of a snode 00077 StorageIndex k; 00078 Index nsuper_et_post = 0; // Number of relaxed snodes in postordered etree 00079 Index nsuper_et = 0; // Number of relaxed snodes in the original etree 00080 StorageIndex l; 00081 for (j = 0; j < n; ) 00082 { 00083 parent = et(j); 00084 snode_start = j; 00085 while ( parent != n && descendants(parent) < relax_columns ) 00086 { 00087 j = parent; 00088 parent = et(j); 00089 } 00090 // Found a supernode in postordered etree, j is the last column 00091 ++nsuper_et_post; 00092 k = StorageIndex(n); 00093 for (Index i = snode_start; i <= j; ++i) 00094 k = (std::min)(k, inv_post(i)); 00095 l = inv_post(j); 00096 if ( (l - k) == (j - snode_start) ) // Same number of columns in the snode 00097 { 00098 // This is also a supernode in the original etree 00099 relax_end(k) = l; // Record last column 00100 ++nsuper_et; 00101 } 00102 else 00103 { 00104 for (Index i = snode_start; i <= j; ++i) 00105 { 00106 l = inv_post(i); 00107 if (descendants(i) == 0) 00108 { 00109 relax_end(l) = l; 00110 ++nsuper_et; 00111 } 00112 } 00113 } 00114 j++; 00115 // Search for a new leaf 00116 while (descendants(j) != 0 && j < n) j++; 00117 } // End postorder traversal of the etree 00118 00119 // Recover the original etree 00120 et = et_save; 00121 } 00122 00123 } // end namespace internal 00124 00125 } // end namespace Eigen 00126 #endif // SPARSELU_HEAP_RELAX_SNODE_H