Background
I'm writing code that searches for the row number of the maximum value in a partial column of a matrix. The partial column corresponds with the non-zero values of an equivalent lower triangular matrix. That is, starting from the input row index idx
through the final row. (This is my own version of subroutine for Gaussian Elimination with Partial Pivoting, for those who are interested).
Basically, the behavior of MaxLocColumnWise
is the search of columnwise maximum location
in the submatrices of A. Passing by reference, MaxLocColumnWise
changes the value of k
.
For example, for the following matrix $$ A = \begin{bmatrix}1 & 2 & 3 \\ 4 & 5 & 6 \\ 7 & 8 & 9 \end{bmatrix} $$
- when
k=0
,MaxLocColumnWise(A,k)
convertsk
to the index of maximum element in[1,4,7]
, which is2
, which is the index of7
i.e. we get2
as a result of
int k = 0;
`MaxLocColumnWise(A,k)`;
std::cout << k <<"\n";
when
k=1
,MaxLocColumnWise(A,k)
convertsk
to the index of maximum element in[5,8]
, which is2
, which is the index of8
. (This case is highlighted in red here)when
k=2
,MaxLocColumnWise(A,k)
convertsk
to the index of maximum element in[9]
, which is2
, which is the index of9
Code
Interestingly it turns out that the following code takes quite long time.
void MaxLocColumnWise(Matrix A, int &idx){
int n = A.size();
int col = idx;
double currentAbsMax = abs(A[col][col]);
for (int i=(col+1); i<n; i++){
double currentVal = abs(A[i][col]);
if (currentVal > currentAbsMax){
currentAbsMax = currentVal;
idx = i;
}
}
}
I have implemented this subroutine to the following Gaussian elimination with partial pivoting
routine:
void GaussElimPartialPivot(Matrix& A, Vector& b){
int n = A.size();
for (int j=0; j<(n-1); j++){
int index = j;
MaxLocColumnWise(A, index);
SwapMatrixRows(A, j, index);
SwapVector(b, j, index);
// main loop
for (int i=(j+1); i<n; i++){
double m = A[i][j]/A[j][j];
b[i] -= m*b[j];
for(int k=j; k<n; k++){
A[i][k] -= m*A[j][k];
}
}
}
}
But when A gets large, the program gets slower exactly due to the subroutine MaxLocColumnWise
, which was verified from excluding each subroutine in the main code.
But I'm not sure exactly where in MaxLocColumnWise()
to blame. Any help will be appreciated.
(An excuse for the exclusion of the code : Matrix
is just from typedef std::vector<Vector> Matrix;
, and the Vector
is typedef std::vector<double> Vector;
)