|
|
|
@ -18,33 +18,33 @@ |
|
|
|
*/ |
|
|
|
|
|
|
|
|
|
|
|
/* |
|
|
|
* Revision and copyright information. |
|
|
|
* |
|
|
|
* Copyright (c) 1985,86,87,88,89,90 |
|
|
|
* by Kenneth S. Kundert and the University of California. |
|
|
|
* |
|
|
|
* Permission to use, copy, modify, and distribute this software and |
|
|
|
* its documentation for any purpose and without fee is hereby granted, |
|
|
|
* provided that the copyright notices appear in all copies and |
|
|
|
* supporting documentation and that the authors and the University of |
|
|
|
* California are properly credited. The authors and the University of |
|
|
|
* California make no representations as to the suitability of this |
|
|
|
* software for any purpose. It is provided `as is', without express |
|
|
|
* or implied warranty. |
|
|
|
*/ |
|
|
|
|
|
|
|
/* |
|
|
|
* IMPORTS |
|
|
|
* |
|
|
|
* >>> Import descriptions: |
|
|
|
* spConfig.h |
|
|
|
* Macros that customize the sparse matrix routines. |
|
|
|
* spMatrix.h |
|
|
|
* Macros and declarations to be imported by the user. |
|
|
|
* spDefs.h |
|
|
|
* Matrix type and macro definitions for the sparse matrix routines. |
|
|
|
*/ |
|
|
|
/* |
|
|
|
* Revision and copyright information. |
|
|
|
* |
|
|
|
* Copyright (c) 1985,86,87,88,89,90 |
|
|
|
* by Kenneth S. Kundert and the University of California. |
|
|
|
* |
|
|
|
* Permission to use, copy, modify, and distribute this software and |
|
|
|
* its documentation for any purpose and without fee is hereby granted, |
|
|
|
* provided that the copyright notices appear in all copies and |
|
|
|
* supporting documentation and that the authors and the University of |
|
|
|
* California are properly credited. The authors and the University of |
|
|
|
* California make no representations as to the suitability of this |
|
|
|
* software for any purpose. It is provided `as is', without express |
|
|
|
* or implied warranty. |
|
|
|
*/ |
|
|
|
|
|
|
|
/* |
|
|
|
* IMPORTS |
|
|
|
* |
|
|
|
* >>> Import descriptions: |
|
|
|
* spConfig.h |
|
|
|
* Macros that customize the sparse matrix routines. |
|
|
|
* spMatrix.h |
|
|
|
* Macros and declarations to be imported by the user. |
|
|
|
* spDefs.h |
|
|
|
* Matrix type and macro definitions for the sparse matrix routines. |
|
|
|
*/ |
|
|
|
#include <assert.h> |
|
|
|
#include <stdlib.h> |
|
|
|
|
|
|
|
@ -63,7 +63,7 @@ int Printer_Width = PRINTER_WIDTH; |
|
|
|
|
|
|
|
|
|
|
|
#if DOCUMENTATION |
|
|
|
|
|
|
|
|
|
|
|
/* |
|
|
|
* PRINT MATRIX |
|
|
|
* |
|
|
|
@ -140,22 +140,22 @@ spPrint(MatrixPtr Matrix, int PrintReordered, int Data, int Header) |
|
|
|
int J = 0; |
|
|
|
int I, Row, Col, Size, Top; |
|
|
|
int StartCol = 1, StopCol, Columns, ElementCount = 0; |
|
|
|
double Magnitude; |
|
|
|
double SmallestDiag = 0; |
|
|
|
double Magnitude; |
|
|
|
double SmallestDiag = 0; |
|
|
|
double SmallestElement = 0; |
|
|
|
double LargestElement = 0.0, LargestDiag = 0.0; |
|
|
|
ElementPtr pElement, *pImagElements; |
|
|
|
int *PrintOrdToIntRowMap, *PrintOrdToIntColMap; |
|
|
|
ElementPtr pElement, * pImagElements; |
|
|
|
int* PrintOrdToIntRowMap, * PrintOrdToIntColMap; |
|
|
|
|
|
|
|
/* Begin `spPrint'. */ |
|
|
|
assert( IS_SPARSE( Matrix ) ); |
|
|
|
assert(IS_SPARSE(Matrix)); |
|
|
|
Size = Matrix->Size; |
|
|
|
SP_CALLOC(pImagElements, ElementPtr, Printer_Width / 10 + 1); |
|
|
|
if ( pImagElements == NULL) |
|
|
|
if (pImagElements == NULL) |
|
|
|
{ |
|
|
|
Matrix->Error = spNO_MEMORY; |
|
|
|
SP_FREE(pImagElements); |
|
|
|
return; |
|
|
|
Matrix->Error = spNO_MEMORY; |
|
|
|
SP_FREE(pImagElements); |
|
|
|
return; |
|
|
|
} |
|
|
|
|
|
|
|
/* Create a packed external to internal row and column translation |
|
|
|
@ -165,49 +165,49 @@ spPrint(MatrixPtr Matrix, int PrintReordered, int Data, int Header) |
|
|
|
#else |
|
|
|
Top = Matrix->AllocatedSize; |
|
|
|
#endif |
|
|
|
SP_CALLOC( PrintOrdToIntRowMap, int, Top + 1 ); |
|
|
|
if ( PrintOrdToIntRowMap == NULL) |
|
|
|
SP_CALLOC(PrintOrdToIntRowMap, int, Top + 1); |
|
|
|
if (PrintOrdToIntRowMap == NULL) |
|
|
|
{ |
|
|
|
Matrix->Error = spNO_MEMORY; |
|
|
|
SP_FREE(pImagElements); |
|
|
|
Matrix->Error = spNO_MEMORY; |
|
|
|
SP_FREE(pImagElements); |
|
|
|
return; |
|
|
|
} |
|
|
|
SP_CALLOC( PrintOrdToIntColMap, int, Top + 1 ); |
|
|
|
SP_CALLOC(PrintOrdToIntColMap, int, Top + 1); |
|
|
|
if (PrintOrdToIntColMap == NULL) |
|
|
|
{ |
|
|
|
Matrix->Error = spNO_MEMORY; |
|
|
|
SP_FREE(pImagElements); |
|
|
|
Matrix->Error = spNO_MEMORY; |
|
|
|
SP_FREE(pImagElements); |
|
|
|
SP_FREE(PrintOrdToIntRowMap); |
|
|
|
return; |
|
|
|
} |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
{ |
|
|
|
PrintOrdToIntRowMap[ Matrix->IntToExtRowMap[I] ] = I; |
|
|
|
PrintOrdToIntColMap[ Matrix->IntToExtColMap[I] ] = I; |
|
|
|
PrintOrdToIntRowMap[Matrix->IntToExtRowMap[I]] = I; |
|
|
|
PrintOrdToIntColMap[Matrix->IntToExtColMap[I]] = I; |
|
|
|
} |
|
|
|
|
|
|
|
/* Pack the arrays. */ |
|
|
|
for (J = 1, I = 1; I <= Top; I++) |
|
|
|
{ |
|
|
|
if (PrintOrdToIntRowMap[I] != 0) |
|
|
|
PrintOrdToIntRowMap[ J++ ] = PrintOrdToIntRowMap[ I ]; |
|
|
|
if (PrintOrdToIntRowMap[I] != 0) |
|
|
|
PrintOrdToIntRowMap[J++] = PrintOrdToIntRowMap[I]; |
|
|
|
} |
|
|
|
for (J = 1, I = 1; I <= Top; I++) |
|
|
|
{ |
|
|
|
if (PrintOrdToIntColMap[I] != 0) |
|
|
|
PrintOrdToIntColMap[ J++ ] = PrintOrdToIntColMap[ I ]; |
|
|
|
if (PrintOrdToIntColMap[I] != 0) |
|
|
|
PrintOrdToIntColMap[J++] = PrintOrdToIntColMap[I]; |
|
|
|
} |
|
|
|
|
|
|
|
/* Print header. */ |
|
|
|
if (Header) |
|
|
|
{ |
|
|
|
printf("MATRIX SUMMARY\n\n"); |
|
|
|
printf("MATRIX SUMMARY\n\n"); |
|
|
|
printf("Size of matrix = %1d x %1d.\n", Size, Size); |
|
|
|
if ( Matrix->Reordered && PrintReordered ) |
|
|
|
if (Matrix->Reordered && PrintReordered) |
|
|
|
printf("Matrix has been reordered.\n"); |
|
|
|
putchar('\n'); |
|
|
|
|
|
|
|
if ( Matrix->Factored ) |
|
|
|
if (Matrix->Factored) |
|
|
|
printf("Matrix after factorization:\n"); |
|
|
|
else |
|
|
|
printf("Matrix before factorization:\n"); |
|
|
|
@ -219,74 +219,74 @@ spPrint(MatrixPtr Matrix, int PrintReordered, int Data, int Header) |
|
|
|
/* Determine how many columns to use. */ |
|
|
|
Columns = Printer_Width; |
|
|
|
if (Header) Columns -= 5; |
|
|
|
if (Data) Columns = (Columns+1) / 10; |
|
|
|
if (Data) Columns = (Columns + 1) / 10; |
|
|
|
|
|
|
|
/* Print matrix by printing groups of complete columns until all |
|
|
|
* the columns are printed. */ |
|
|
|
J = 0; |
|
|
|
while ( J <= Size ) |
|
|
|
while (J <= Size) |
|
|
|
{ |
|
|
|
/* Calculatestrchr of last column to printed in this group. */ |
|
|
|
StopCol = StartCol + Columns - 1; |
|
|
|
/* Calculatestrchr of last column to printed in this group. */ |
|
|
|
StopCol = StartCol + Columns - 1; |
|
|
|
if (StopCol > Size) |
|
|
|
StopCol = Size; |
|
|
|
|
|
|
|
/* Label the columns. */ |
|
|
|
/* Label the columns. */ |
|
|
|
if (Header) |
|
|
|
{ |
|
|
|
if (Data) |
|
|
|
if (Data) |
|
|
|
{ |
|
|
|
printf(" "); |
|
|
|
printf(" "); |
|
|
|
for (I = StartCol; I <= StopCol; I++) |
|
|
|
{ |
|
|
|
if (PrintReordered) |
|
|
|
if (PrintReordered) |
|
|
|
Col = I; |
|
|
|
else |
|
|
|
Col = PrintOrdToIntColMap[I]; |
|
|
|
printf(" %9d", Matrix->IntToExtColMap[ Col ]); |
|
|
|
printf(" %9d", Matrix->IntToExtColMap[Col]); |
|
|
|
} |
|
|
|
printf("\n\n"); |
|
|
|
} |
|
|
|
else |
|
|
|
{ |
|
|
|
if (PrintReordered) |
|
|
|
printf("Columns %1d to %1d.\n",StartCol,StopCol); |
|
|
|
if (PrintReordered) |
|
|
|
printf("Columns %1d to %1d.\n", StartCol, StopCol); |
|
|
|
else |
|
|
|
{ |
|
|
|
printf("Columns %1d to %1d.\n", |
|
|
|
Matrix->IntToExtColMap[ PrintOrdToIntColMap[StartCol] ], |
|
|
|
Matrix->IntToExtColMap[ PrintOrdToIntColMap[StopCol] ]); |
|
|
|
printf("Columns %1d to %1d.\n", |
|
|
|
Matrix->IntToExtColMap[PrintOrdToIntColMap[StartCol]], |
|
|
|
Matrix->IntToExtColMap[PrintOrdToIntColMap[StopCol]]); |
|
|
|
} |
|
|
|
} |
|
|
|
} |
|
|
|
|
|
|
|
/* Print every row ... */ |
|
|
|
/* Print every row ... */ |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
{ |
|
|
|
if (PrintReordered) |
|
|
|
if (PrintReordered) |
|
|
|
Row = I; |
|
|
|
else |
|
|
|
Row = PrintOrdToIntRowMap[I]; |
|
|
|
|
|
|
|
if (Header) |
|
|
|
{ |
|
|
|
if (PrintReordered && !Data) |
|
|
|
if (PrintReordered && !Data) |
|
|
|
printf("%4d", I); |
|
|
|
else |
|
|
|
printf("%4d", Matrix->IntToExtRowMap[ Row ]); |
|
|
|
printf("%4d", Matrix->IntToExtRowMap[Row]); |
|
|
|
if (!Data) putchar(' '); |
|
|
|
} |
|
|
|
|
|
|
|
/* ... in each column of the group. */ |
|
|
|
/* ... in each column of the group. */ |
|
|
|
for (J = StartCol; J <= StopCol; J++) |
|
|
|
{ |
|
|
|
if (PrintReordered) |
|
|
|
if (PrintReordered) |
|
|
|
Col = J; |
|
|
|
else |
|
|
|
Col = PrintOrdToIntColMap[J]; |
|
|
|
|
|
|
|
pElement = Matrix->FirstInCol[Col]; |
|
|
|
while(pElement != NULL && pElement->Row != Row) |
|
|
|
while (pElement != NULL && pElement->Row != Row) |
|
|
|
pElement = pElement->NextInCol; |
|
|
|
|
|
|
|
if (Data) |
|
|
|
@ -294,24 +294,24 @@ spPrint(MatrixPtr Matrix, int PrintReordered, int Data, int Header) |
|
|
|
|
|
|
|
if (pElement != NULL) |
|
|
|
{ |
|
|
|
/* Case where element exists */ |
|
|
|
if (Data) |
|
|
|
/* Case where element exists */ |
|
|
|
if (Data) |
|
|
|
printf(" %9.3g", pElement->Real); |
|
|
|
else |
|
|
|
putchar('x'); |
|
|
|
|
|
|
|
/* Update status variables */ |
|
|
|
if ( (Magnitude = ELEMENT_MAG(pElement)) > LargestElement ) |
|
|
|
/* Update status variables */ |
|
|
|
if ((Magnitude = ELEMENT_MAG(pElement)) > LargestElement) |
|
|
|
LargestElement = Magnitude; |
|
|
|
if ((Magnitude < SmallestElement) && (Magnitude != 0.0)) |
|
|
|
SmallestElement = Magnitude; |
|
|
|
ElementCount++; |
|
|
|
} |
|
|
|
|
|
|
|
/* Case where element is structurally zero */ |
|
|
|
/* Case where element is structurally zero */ |
|
|
|
else |
|
|
|
{ |
|
|
|
if (Data) |
|
|
|
if (Data) |
|
|
|
printf(" ..."); |
|
|
|
else |
|
|
|
putchar('.'); |
|
|
|
@ -321,13 +321,13 @@ spPrint(MatrixPtr Matrix, int PrintReordered, int Data, int Header) |
|
|
|
|
|
|
|
if (Matrix->Complex && Data) |
|
|
|
{ |
|
|
|
printf(" "); |
|
|
|
printf(" "); |
|
|
|
for (J = StartCol; J <= StopCol; J++) |
|
|
|
{ |
|
|
|
if (pImagElements[J - StartCol] != NULL) |
|
|
|
if (pImagElements[J - StartCol] != NULL) |
|
|
|
{ |
|
|
|
printf(" %8.2gj", |
|
|
|
pImagElements[J-StartCol]->Imag); |
|
|
|
printf(" %8.2gj", |
|
|
|
pImagElements[J - StartCol]->Imag); |
|
|
|
} |
|
|
|
else printf(" "); |
|
|
|
} |
|
|
|
@ -335,44 +335,44 @@ spPrint(MatrixPtr Matrix, int PrintReordered, int Data, int Header) |
|
|
|
} |
|
|
|
} |
|
|
|
|
|
|
|
/* Calculatestrchr of first column in next group. */ |
|
|
|
/* Calculatestrchr of first column in next group. */ |
|
|
|
StartCol = StopCol; |
|
|
|
StartCol++; |
|
|
|
putchar('\n'); |
|
|
|
} |
|
|
|
if (Header) |
|
|
|
{ |
|
|
|
printf("\nLargest element in matrix = %-1.4g.\n", LargestElement); |
|
|
|
printf("\nLargest element in matrix = %-1.4g.\n", LargestElement); |
|
|
|
printf("Smallest element in matrix = %-1.4g.\n", SmallestElement); |
|
|
|
|
|
|
|
/* Search for largest and smallest diagonal values */ |
|
|
|
/* Search for largest and smallest diagonal values */ |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
{ |
|
|
|
if (Matrix->Diag[I] != NULL) |
|
|
|
if (Matrix->Diag[I] != NULL) |
|
|
|
{ |
|
|
|
Magnitude = ELEMENT_MAG( Matrix->Diag[I] ); |
|
|
|
if ( Magnitude > LargestDiag ) LargestDiag = Magnitude; |
|
|
|
if ( Magnitude < SmallestDiag ) SmallestDiag = Magnitude; |
|
|
|
Magnitude = ELEMENT_MAG(Matrix->Diag[I]); |
|
|
|
if (Magnitude > LargestDiag) LargestDiag = Magnitude; |
|
|
|
if (Magnitude < SmallestDiag) SmallestDiag = Magnitude; |
|
|
|
} |
|
|
|
} |
|
|
|
|
|
|
|
/* Print the largest and smallest diagonal values */ |
|
|
|
if ( Matrix->Factored ) |
|
|
|
/* Print the largest and smallest diagonal values */ |
|
|
|
if (Matrix->Factored) |
|
|
|
{ |
|
|
|
printf("\nLargest diagonal element = %-1.4g.\n", LargestDiag); |
|
|
|
printf("\nLargest diagonal element = %-1.4g.\n", LargestDiag); |
|
|
|
printf("Smallest diagonal element = %-1.4g.\n", SmallestDiag); |
|
|
|
} |
|
|
|
else |
|
|
|
{ |
|
|
|
printf("\nLargest pivot element = %-1.4g.\n", LargestDiag); |
|
|
|
printf("\nLargest pivot element = %-1.4g.\n", LargestDiag); |
|
|
|
printf("Smallest pivot element = %-1.4g.\n", SmallestDiag); |
|
|
|
} |
|
|
|
|
|
|
|
/* Calculate and print sparsity and number of fill-ins created. */ |
|
|
|
/* Calculate and print sparsity and number of fill-ins created. */ |
|
|
|
printf("\nDensity = %2.2f%%.\n", ((double)(ElementCount * 100)) / |
|
|
|
((double)(Size * Size))); |
|
|
|
((double)(Size * Size))); |
|
|
|
|
|
|
|
printf("Number of originals = %1d.\n", Matrix->Originals); |
|
|
|
printf("Number of originals = %1d.\n", Matrix->Originals); |
|
|
|
if (!Matrix->NeedsOrdering) |
|
|
|
printf("Number of fill-ins = %1d.\n", Matrix->Fillins); |
|
|
|
} |
|
|
|
@ -393,7 +393,7 @@ spPrint(MatrixPtr Matrix, int PrintReordered, int Data, int Header) |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* |
|
|
|
* OUTPUT MATRIX TO FILE |
|
|
|
* |
|
|
|
@ -437,16 +437,16 @@ spPrint(MatrixPtr Matrix, int PrintReordered, int Data, int Header) |
|
|
|
*/ |
|
|
|
|
|
|
|
int |
|
|
|
spFileMatrix(MatrixPtr Matrix, char *File, char *Label, int Reordered, |
|
|
|
int Data, int Header) |
|
|
|
spFileMatrix(MatrixPtr Matrix, char* File, char* Label, int Reordered, |
|
|
|
int Data, int Header) |
|
|
|
{ |
|
|
|
int I, Size; |
|
|
|
ElementPtr pElement; |
|
|
|
int Row, Col, Err; |
|
|
|
FILE *pMatrixFile; |
|
|
|
FILE* pMatrixFile; |
|
|
|
|
|
|
|
/* Begin `spFileMatrix'. */ |
|
|
|
assert( IS_SPARSE( Matrix ) ); |
|
|
|
assert(IS_SPARSE(Matrix)); |
|
|
|
|
|
|
|
/* Open file matrix file in write mode. */ |
|
|
|
if ((pMatrixFile = fopen(File, "w")) == NULL) |
|
|
|
@ -456,100 +456,100 @@ spFileMatrix(MatrixPtr Matrix, char *File, char *Label, int Reordered, |
|
|
|
Size = Matrix->Size; |
|
|
|
if (Header) |
|
|
|
{ |
|
|
|
if (Matrix->Factored && Data) |
|
|
|
if (Matrix->Factored && Data) |
|
|
|
{ |
|
|
|
Err = fprintf(pMatrixFile, |
|
|
|
"Warning : The following matrix is " |
|
|
|
"factored in to LU form.\n"); |
|
|
|
if (Err < 0) |
|
|
|
return 0; |
|
|
|
Err = fprintf(pMatrixFile, |
|
|
|
"Warning : The following matrix is " |
|
|
|
"factored in to LU form.\n"); |
|
|
|
if (Err < 0) |
|
|
|
return 0; |
|
|
|
} |
|
|
|
if (fprintf(pMatrixFile, "%s\n", Label) < 0) |
|
|
|
return 0; |
|
|
|
Err = fprintf( pMatrixFile, "%d\t%s\n", Size, |
|
|
|
(Matrix->Complex ? "complex" : "real")); |
|
|
|
return 0; |
|
|
|
Err = fprintf(pMatrixFile, "%d\t%s\n", Size, |
|
|
|
(Matrix->Complex ? "complex" : "real")); |
|
|
|
if (Err < 0) |
|
|
|
return 0; |
|
|
|
return 0; |
|
|
|
} |
|
|
|
|
|
|
|
/* Output matrix. */ |
|
|
|
if (!Data) |
|
|
|
{ |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
{ |
|
|
|
pElement = Matrix->FirstInCol[I]; |
|
|
|
pElement = Matrix->FirstInCol[I]; |
|
|
|
while (pElement != NULL) |
|
|
|
{ |
|
|
|
if (Reordered) |
|
|
|
if (Reordered) |
|
|
|
{ |
|
|
|
Row = pElement->Row; |
|
|
|
Row = pElement->Row; |
|
|
|
Col = I; |
|
|
|
} |
|
|
|
else |
|
|
|
{ |
|
|
|
Row = Matrix->IntToExtRowMap[pElement->Row]; |
|
|
|
Row = Matrix->IntToExtRowMap[pElement->Row]; |
|
|
|
Col = Matrix->IntToExtColMap[I]; |
|
|
|
} |
|
|
|
pElement = pElement->NextInCol; |
|
|
|
if (fprintf(pMatrixFile, "%d\t%d\n", Row, Col) < 0) return 0; |
|
|
|
} |
|
|
|
} |
|
|
|
/* Output terminator, a line of zeros. */ |
|
|
|
/* Output terminator, a line of zeros. */ |
|
|
|
if (Header) |
|
|
|
if (fprintf(pMatrixFile, "0\t0\n") < 0) return 0; |
|
|
|
} |
|
|
|
|
|
|
|
if (Data && Matrix->Complex) |
|
|
|
{ |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
{ |
|
|
|
pElement = Matrix->FirstInCol[I]; |
|
|
|
pElement = Matrix->FirstInCol[I]; |
|
|
|
while (pElement != NULL) |
|
|
|
{ |
|
|
|
if (Reordered) |
|
|
|
if (Reordered) |
|
|
|
{ |
|
|
|
Row = pElement->Row; |
|
|
|
Row = pElement->Row; |
|
|
|
Col = I; |
|
|
|
} |
|
|
|
else |
|
|
|
{ |
|
|
|
Row = Matrix->IntToExtRowMap[pElement->Row]; |
|
|
|
Row = Matrix->IntToExtRowMap[pElement->Row]; |
|
|
|
Col = Matrix->IntToExtColMap[I]; |
|
|
|
} |
|
|
|
Err = fprintf |
|
|
|
( pMatrixFile,"%d\t%d\t%-.15g\t%-.15g\n", |
|
|
|
Row, Col, pElement->Real, pElement->Imag |
|
|
|
); |
|
|
|
(pMatrixFile, "%d\t%d\t%-.15g\t%-.15g\n", |
|
|
|
Row, Col, pElement->Real, pElement->Imag |
|
|
|
); |
|
|
|
if (Err < 0) return 0; |
|
|
|
pElement = pElement->NextInCol; |
|
|
|
} |
|
|
|
} |
|
|
|
/* Output terminator, a line of zeros. */ |
|
|
|
/* Output terminator, a line of zeros. */ |
|
|
|
if (Header) |
|
|
|
if (fprintf(pMatrixFile,"0\t0\t0.0\t0.0\n") < 0) return 0; |
|
|
|
if (fprintf(pMatrixFile, "0\t0\t0.0\t0.0\n") < 0) return 0; |
|
|
|
|
|
|
|
} |
|
|
|
|
|
|
|
if (Data && !Matrix->Complex) |
|
|
|
{ |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
{ |
|
|
|
pElement = Matrix->FirstInCol[I]; |
|
|
|
pElement = Matrix->FirstInCol[I]; |
|
|
|
while (pElement != NULL) |
|
|
|
{ |
|
|
|
Row = Matrix->IntToExtRowMap[pElement->Row]; |
|
|
|
Row = Matrix->IntToExtRowMap[pElement->Row]; |
|
|
|
Col = Matrix->IntToExtColMap[I]; |
|
|
|
Err = fprintf |
|
|
|
( pMatrixFile,"%d\t%d\t%-.15g\n", |
|
|
|
Row, Col, pElement->Real |
|
|
|
); |
|
|
|
(pMatrixFile, "%d\t%d\t%-.15g\n", |
|
|
|
Row, Col, pElement->Real |
|
|
|
); |
|
|
|
if (Err < 0) return 0; |
|
|
|
pElement = pElement->NextInCol; |
|
|
|
} |
|
|
|
} |
|
|
|
/* Output terminator, a line of zeros. */ |
|
|
|
/* Output terminator, a line of zeros. */ |
|
|
|
if (Header) |
|
|
|
if (fprintf(pMatrixFile,"0\t0\t0.0\n") < 0) return 0; |
|
|
|
if (fprintf(pMatrixFile, "0\t0\t0.0\n") < 0) return 0; |
|
|
|
|
|
|
|
} |
|
|
|
|
|
|
|
@ -563,7 +563,7 @@ spFileMatrix(MatrixPtr Matrix, char *File, char *Label, int Reordered, |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* |
|
|
|
* OUTPUT SOURCE VECTOR TO FILE |
|
|
|
* |
|
|
|
@ -595,22 +595,22 @@ spFileMatrix(MatrixPtr Matrix, char *File, char *Label, int Reordered, |
|
|
|
*/ |
|
|
|
|
|
|
|
int |
|
|
|
spFileVector(MatrixPtr Matrix, char *File, RealVector RHS, RealVector iRHS) |
|
|
|
spFileVector(MatrixPtr Matrix, char* File, RealVector RHS, RealVector iRHS) |
|
|
|
{ |
|
|
|
int I, Size, Err; |
|
|
|
FILE *pMatrixFile; |
|
|
|
FILE* pMatrixFile; |
|
|
|
|
|
|
|
/* Begin `spFileVector'. */ |
|
|
|
assert( IS_SPARSE( Matrix ) && RHS != NULL); |
|
|
|
assert(IS_SPARSE(Matrix) && RHS != NULL); |
|
|
|
|
|
|
|
if (File) { |
|
|
|
/* Open File in write mode. */ |
|
|
|
pMatrixFile = fopen(File,"w"); |
|
|
|
pMatrixFile = fopen(File, "w"); |
|
|
|
if (pMatrixFile == NULL) |
|
|
|
return 0; |
|
|
|
return 0; |
|
|
|
} |
|
|
|
else |
|
|
|
pMatrixFile=stdout; |
|
|
|
else |
|
|
|
pMatrixFile = stdout; |
|
|
|
|
|
|
|
/* Output vector. */ |
|
|
|
Size = Matrix->Size; |
|
|
|
@ -618,18 +618,18 @@ spFileVector(MatrixPtr Matrix, char *File, RealVector RHS, RealVector iRHS) |
|
|
|
{ |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
{ |
|
|
|
Err = fprintf |
|
|
|
( pMatrixFile, "%-.15g\t%-.15g\n", |
|
|
|
RHS[I], iRHS[I] |
|
|
|
); |
|
|
|
Err = fprintf |
|
|
|
(pMatrixFile, "%-.15g\t%-.15g\n", |
|
|
|
RHS[I], iRHS[I] |
|
|
|
); |
|
|
|
if (Err < 0) return 0; |
|
|
|
} |
|
|
|
} |
|
|
|
else |
|
|
|
{ |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
{ |
|
|
|
if (fprintf(pMatrixFile, "%-.15g\n", RHS[I]) < 0) |
|
|
|
if (fprintf(pMatrixFile, "%-.15g\n", RHS[I]) < 0) |
|
|
|
return 0; |
|
|
|
} |
|
|
|
} |
|
|
|
@ -647,13 +647,13 @@ spFileVector(MatrixPtr Matrix, char *File, RealVector RHS, RealVector iRHS) |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* |
|
|
|
* OUTPUT STATISTICS TO FILE |
|
|
|
* |
|
|
|
* Writes useful information concerning the matrix to a file. Should be |
|
|
|
* executed after the matrix is factored. |
|
|
|
* |
|
|
|
* |
|
|
|
* >>> Returns: |
|
|
|
* One is returned if routine was successful, otherwise zero is returned. |
|
|
|
* The calling function can query errno (the system global error variable) |
|
|
|
@ -685,16 +685,16 @@ spFileVector(MatrixPtr Matrix, char *File, RealVector RHS, RealVector iRHS) |
|
|
|
*/ |
|
|
|
|
|
|
|
int |
|
|
|
spFileStats(MatrixPtr Matrix, char *File, char *Label) |
|
|
|
spFileStats(MatrixPtr Matrix, char* File, char* Label) |
|
|
|
{ |
|
|
|
int Size, I; |
|
|
|
ElementPtr pElement; |
|
|
|
int NumberOfElements; |
|
|
|
RealNumber Data, LargestElement, SmallestElement; |
|
|
|
FILE *pStatsFile; |
|
|
|
FILE* pStatsFile; |
|
|
|
|
|
|
|
/* Begin `spFileStats'. */ |
|
|
|
assert( IS_SPARSE( Matrix ) ); |
|
|
|
assert(IS_SPARSE(Matrix)); |
|
|
|
|
|
|
|
/* Open File in append mode. */ |
|
|
|
if ((pStatsFile = fopen(File, "a")) == NULL) |
|
|
|
@ -710,7 +710,7 @@ spFileStats(MatrixPtr Matrix, char *File, char *Label) |
|
|
|
fprintf(pStatsFile, "Matrix is complex.\n"); |
|
|
|
else |
|
|
|
fprintf(pStatsFile, "Matrix is real.\n"); |
|
|
|
fprintf(pStatsFile," Size = %d\n",Size); |
|
|
|
fprintf(pStatsFile, " Size = %d\n", Size); |
|
|
|
|
|
|
|
/* Search matrix. */ |
|
|
|
NumberOfElements = 0; |
|
|
|
@ -719,10 +719,10 @@ spFileStats(MatrixPtr Matrix, char *File, char *Label) |
|
|
|
|
|
|
|
for (I = 1; I <= Size; I++) |
|
|
|
{ |
|
|
|
pElement = Matrix->FirstInCol[I]; |
|
|
|
pElement = Matrix->FirstInCol[I]; |
|
|
|
while (pElement != NULL) |
|
|
|
{ |
|
|
|
NumberOfElements++; |
|
|
|
NumberOfElements++; |
|
|
|
Data = ELEMENT_MAG(pElement); |
|
|
|
if (Data > LargestElement) |
|
|
|
LargestElement = Data; |
|
|
|
@ -732,27 +732,27 @@ spFileStats(MatrixPtr Matrix, char *File, char *Label) |
|
|
|
} |
|
|
|
} |
|
|
|
|
|
|
|
SmallestElement = MIN( SmallestElement, LargestElement ); |
|
|
|
SmallestElement = MIN(SmallestElement, LargestElement); |
|
|
|
|
|
|
|
/* Output remaining statistics. */ |
|
|
|
fprintf(pStatsFile, " Initial number of elements = %d\n", |
|
|
|
NumberOfElements - Matrix->Fillins); |
|
|
|
NumberOfElements - Matrix->Fillins); |
|
|
|
fprintf(pStatsFile, |
|
|
|
" Initial average number of elements per row = %f\n", |
|
|
|
(double)(NumberOfElements - Matrix->Fillins) / (double)Size); |
|
|
|
fprintf(pStatsFile, " Fill-ins = %d\n",Matrix->Fillins); |
|
|
|
" Initial average number of elements per row = %f\n", |
|
|
|
(double)(NumberOfElements - Matrix->Fillins) / (double)Size); |
|
|
|
fprintf(pStatsFile, " Fill-ins = %d\n", Matrix->Fillins); |
|
|
|
fprintf(pStatsFile, " Average number of fill-ins per row = %f%%\n", |
|
|
|
(double)Matrix->Fillins / (double)Size); |
|
|
|
(double)Matrix->Fillins / (double)Size); |
|
|
|
fprintf(pStatsFile, " Total number of elements = %d\n", |
|
|
|
NumberOfElements); |
|
|
|
NumberOfElements); |
|
|
|
fprintf(pStatsFile, " Average number of elements per row = %f\n", |
|
|
|
(double)NumberOfElements / (double)Size); |
|
|
|
fprintf(pStatsFile," Density = %f%%\n", |
|
|
|
(double)(100.0*NumberOfElements)/(double)(Size*Size)); |
|
|
|
fprintf(pStatsFile," Relative Threshold = %e\n", Matrix->RelThreshold); |
|
|
|
fprintf(pStatsFile," Absolute Threshold = %e\n", Matrix->AbsThreshold); |
|
|
|
fprintf(pStatsFile," Largest Element = %e\n", LargestElement); |
|
|
|
fprintf(pStatsFile," Smallest Element = %e\n\n\n", SmallestElement); |
|
|
|
(double)NumberOfElements / (double)Size); |
|
|
|
fprintf(pStatsFile, " Density = %f%%\n", |
|
|
|
(double)(100.0 * NumberOfElements) / (double)(Size * Size)); |
|
|
|
fprintf(pStatsFile, " Relative Threshold = %e\n", Matrix->RelThreshold); |
|
|
|
fprintf(pStatsFile, " Absolute Threshold = %e\n", Matrix->AbsThreshold); |
|
|
|
fprintf(pStatsFile, " Largest Element = %e\n", LargestElement); |
|
|
|
fprintf(pStatsFile, " Smallest Element = %e\n\n\n", SmallestElement); |
|
|
|
|
|
|
|
/* Close file. */ |
|
|
|
(void)fclose(pStatsFile); |
|
|
|
|