1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
|
// 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>
//
// This Source Code Form is subject to the terms of the Mozilla
// Public License v. 2.0. If a copy of the MPL was not distributed
// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
/*
* NOTE: This file is the modified version of [s,d,c,z]pruneL.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_PRUNEL_H
#define SPARSELU_PRUNEL_H
/**
* \brief Prunes the L-structure.
*
* It prunes the L-structure of supernodes whose L-structure contains the current pivot row "pivrow"
*
*
* \param jcol The current column of L
* \param [in]perm_r Row permutation
* \param [out]pivrow The pivot row
* \param nseg Number of segments
* \param segrep
* \param repfnz
* \param [out]xprune
* \param glu Global LU data
*
*/
template <typename Scalar, typename Index>
void SparseLUBase<Scalar,Index>::LU_pruneL(const int jcol, const IndexVector& perm_r, const int pivrow, const int nseg, const IndexVector& segrep, BlockIndexVector& repfnz, IndexVector& xprune, GlobalLU_t& glu)
{
// For each supernode-rep irep in U(*,j]
int jsupno = glu.supno(jcol);
int i,irep,irep1;
bool movnum, do_prune = false;
Index kmin, kmax, minloc, maxloc,krow;
for (i = 0; i < nseg; i++)
{
irep = segrep(i);
irep1 = irep + 1;
do_prune = false;
// Don't prune with a zero U-segment
if (repfnz(irep) == IND_EMPTY) continue;
// If a snode overlaps with the next panel, then the U-segment
// is fragmented into two parts -- irep and irep1. We should let
// pruning occur at the rep-column in irep1s snode.
if (glu.supno(irep) == glu.supno(irep1) ) continue; // don't prune
// If it has not been pruned & it has a nonz in row L(pivrow,i)
if (glu.supno(irep) != jsupno )
{
if ( xprune (irep) >= glu.xlsub(irep1) )
{
kmin = glu.xlsub(irep);
kmax = glu.xlsub(irep1) - 1;
for (krow = kmin; krow <= kmax; krow++)
{
if (glu.lsub(krow) == pivrow)
{
do_prune = true;
break;
}
}
}
if (do_prune)
{
// do a quicksort-type partition
// movnum=true means that the num values have to be exchanged
movnum = false;
if (irep == glu.xsup(glu.supno(irep)) ) // Snode of size 1
movnum = true;
while (kmin <= kmax)
{
if (perm_r(glu.lsub(kmax)) == IND_EMPTY)
kmax--;
else if ( perm_r(glu.lsub(kmin)) != IND_EMPTY)
kmin++;
else
{
// kmin below pivrow (not yet pivoted), and kmax
// above pivrow: interchange the two suscripts
std::swap(glu.lsub(kmin), glu.lsub(kmax));
// If the supernode has only one column, then we
// only keep one set of subscripts. For any subscript
// intercnahge performed, similar interchange must be
// done on the numerical values.
if (movnum)
{
minloc = glu.xlusup(irep) + ( kmin - glu.xlsub(irep) );
maxloc = glu.xlusup(irep) + ( kmax - glu.xlsub(irep) );
std::swap(glu.lusup(minloc), glu.lusup(maxloc));
}
kmin++;
kmax--;
}
} // end while
xprune(irep) = kmin; //Pruning
} // end if do_prune
} // end pruning
} // End for each U-segment
}
#endif
|
