diff options
author | Desire NUENTSA <desire.nuentsa_wakam@inria.fr> | 2012-05-31 17:10:29 +0200 |
---|---|---|
committer | Desire NUENTSA <desire.nuentsa_wakam@inria.fr> | 2012-05-31 17:10:29 +0200 |
commit | b26d6b02de24f2c96f4bdfb6bf1c42afc80693c6 (patch) | |
tree | b899e6f43d2b7dfd6926eda6c514d7c50d317f0d /Eigen/src/SparseLU/SparseLU_copy_to_ucol.h | |
parent | 8608d08d658b09bfd92057d752eb80d59462cdc8 (diff) |
Eliminate and prune columns in a panel
Diffstat (limited to 'Eigen/src/SparseLU/SparseLU_copy_to_ucol.h')
-rw-r--r-- | Eigen/src/SparseLU/SparseLU_copy_to_ucol.h | 123 |
1 files changed, 123 insertions, 0 deletions
diff --git a/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h b/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h new file mode 100644 index 000000000..3f8d8abe2 --- /dev/null +++ b/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h @@ -0,0 +1,123 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2012 Désiré Nuentsa-Wakam <desire.nuentsa_wakam@inria.fr> +// +// Eigen is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 3 of the License, or (at your option) any later version. +// +// Alternatively, you can redistribute it and/or +// modify it under the terms of the GNU General Public License as +// published by the Free Software Foundation; either version 2 of +// the License, or (at your option) any later version. +// +// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY +// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS +// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the +// GNU General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License and a copy of the GNU General Public License along with +// Eigen. If not, see <http://www.gnu.org/licenses/>. + +/* + + * NOTE: This file is the modified version of xcopy_to_ucol.c file in SuperLU + + * -- SuperLU routine (version 2.0) -- + * Univ. of California Berkeley, Xerox Palo Alto Research Center, + * and Lawrence Berkeley National Lab. + * November 15, 1997 + * + * Copyright (c) 1994 by Xerox Corporation. All rights reserved. + * + * THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY + * EXPRESSED OR IMPLIED. ANY USE IS AT YOUR OWN RISK. + * + * Permission is hereby granted to use or copy this program for any + * purpose, provided the above notices are retained on all copies. + * Permission to modify the code and to distribute modified code is + * granted, provided the above notices are retained, and a notice that + * the code was modified is included with the above copyright notice. + */ +#ifndef SPARSELU_COPY_TO_UCOL_H +#define SPARSELU_COPY_TO_UCOL_H +/** + * \brief Performs numeric block updates (sup-col) in topological order + * + * \param jcol current column to update + * \param nseg Number of segments in the U part + * \param segrep segment representative ... + * \param repfnz ??? First nonzero column in each row ??? ... + * \param perm_r Row permutation + * \param dense Store the full representation of the column + * \param Glu Global LU data. + * \return 0 - successful return + * > 0 - number of bytes allocated when run out of space + * + */ +template <typename VectorType> +int SparseLU::LU_copy_to_ucol(const int jcol, const int nseg, VectorXi& segrep, VectorXi& repfnz, VectorXi& perm_r, VectorType& dense, LU_GlobalLu_t& Glu) +{ + int ksupno, k, ksub, krep, ksupno; + + VectorXi& xsup = Glu.xsup; + VectorXi& supno = Glu.supno; + VectorXi& lsub = Glu.lsub; + VectorXi& xlsub = Glu.xlsub; + VectorType& ucol = GLu.ucol; + VectorXi& usub = Glu.usub; + VectorXi& xusub = Glu.xusub; + int nzumax = GLu.nzumax; + int jsupno = supno(jcol); + + // For each nonzero supernode segment of U[*,j] in topological order + k = nseg - 1; + int nextu = xusub(jcol); + int kfnz, isub, segsize; + int new_next,irow; + for (ksub = 0; ksub < nseg; ksub++) + { + krep = segrep(k); k--; + ksupno = supno(krep); + if (jsupno != ksupno ) // should go into ucol(); + { + kfnz = repfnz(krep); + if (kfnz != IND_EMPTY) + { // Nonzero U-segment + fsupc = xsup(ksupno); + isub = xlsub(fsupc) + kfnz - fsupc; + segsize = krep - kfnz + 1; + new_next = nextu + segsize; + while (new_next > nzumax) + { + Glu.ucol = LU_MemXpand<Scalar>(jcol, nextu, UCOL, nzumax); //FIXME try and catch errors + ucol = Glu.ucol; + Glu.nzumax = nzumax; + Glu.usub = LU_MemXpand<Index>(jcol, nextu, USUB, nzumax); //FIXME try and catch errors + Glu.nzumax = nzumax; + usub = Glu.usub; + lsub = Glu.lsub; + } + + for (i = 0; i < segsize; i++) + { + irow = lsub(isub); + usub(nextu) = perm_r(irow); // Unlike teh L part, the U part is stored in its final order + ucol(nextu) = dense(irow); + dense(irow) = Scalar(0.0); + nextu++; + isub++; + } + + } // end nonzero U-segment + + } // end if jsupno + + } // end for each segment + xusub(jcol + 1) = nextu; // close U(*,jcol) + return 0; +} +#endif
\ No newline at end of file |