From b9d69d32bd7b93792bc9296c8558001b5f0a64d2 Mon Sep 17 00:00:00 2001 From: Holger Vogt Date: Fri, 28 Jun 2024 22:26:55 +0200 Subject: [PATCH] Re-formatting --- src/maths/KLU/klusmp.c | 1928 +++++++++++++++++++++------------------- 1 file changed, 992 insertions(+), 936 deletions(-) diff --git a/src/maths/KLU/klusmp.c b/src/maths/KLU/klusmp.c index d22c36f94..67cf106fa 100644 --- a/src/maths/KLU/klusmp.c +++ b/src/maths/KLU/klusmp.c @@ -26,435 +26,441 @@ extern double logb(double); #include "ngspice/evt.h" static int -CircuitIsDigital (void) +CircuitIsDigital(void) { #ifdef XSPICE - return g_mif_info.ckt && g_mif_info.ckt->evt && g_mif_info.ckt->evt->counts.num_nodes != 0 ; + return g_mif_info.ckt && g_mif_info.ckt->evt && g_mif_info.ckt->evt->counts.num_nodes != 0; #else - return 0 ; + return 0; #endif } -static void LoadGmin_CSC (double **diag, unsigned int n, double Gmin) ; -static void LoadGmin (SMPmatrix *eMatrix, double Gmin) ; +static void LoadGmin_CSC(double** diag, unsigned int n, double Gmin); +static void LoadGmin(SMPmatrix* eMatrix, double Gmin); typedef struct sElement { - unsigned int row ; - unsigned int col ; - double *pointer ; - unsigned int group ; -} Element ; + unsigned int row; + unsigned int col; + double* pointer; + unsigned int group; +} Element; static int -CompareRow (const void *a, const void *b) +CompareRow(const void* a, const void* b) { - Element *A = (Element *) a ; - Element *B = (Element *) b ; + Element* A = (Element*)a; + Element* B = (Element*)b; return - (A->row > B->row) ? 1 : + (A->row > B->row) ? 1 : (A->row < B->row) ? -1 : - 0 ; + 0; } static int -CompareColumn (const void *a, const void *b) +CompareColumn(const void* a, const void* b) { - Element *A = (Element *) a ; - Element *B = (Element *) b ; + Element* A = (Element*)a; + Element* B = (Element*)b; return - (A->col > B->col) ? 1 : + (A->col > B->col) ? 1 : (A->col < B->col) ? -1 : - 0 ; + 0; } static void -Compress (unsigned int *Ai, unsigned int *Bp, unsigned int n, unsigned int nz) +Compress(unsigned int* Ai, unsigned int* Bp, unsigned int n, unsigned int nz) { - unsigned int i, j ; + unsigned int i, j; - for (i = 0 ; i <= Ai [0] ; i++) - Bp [i] = 0 ; + for (i = 0; i <= Ai[0]; i++) + Bp[i] = 0; - j = Ai [0] + 1 ; - for (i = 1 ; i < nz ; i++) + j = Ai[0] + 1; + for (i = 1; i < nz; i++) { - if (Ai [i] == Ai [i - 1] + 1) + if (Ai[i] == Ai[i - 1] + 1) { - Bp [j] = i ; - j++ ; + Bp[j] = i; + j++; } - else if (Ai [i] > Ai [i - 1] + 1) + else if (Ai[i] > Ai[i - 1] + 1) { - for ( ; j <= Ai [i] ; j++) - Bp [j] = i ; + for (; j <= Ai[i]; j++) + Bp[j] = i; } } - for ( ; j <= n ; j++) - Bp [j] = i ; + for (; j <= n; j++) + Bp[j] = i; } int -BindCompare (const void *a, const void *b) +BindCompare(const void* a, const void* b) { - BindElement *A = (BindElement *) a ; - BindElement *B = (BindElement *) b ; + BindElement* A = (BindElement*)a; + BindElement* B = (BindElement*)b; return - (A->COO > B->COO) ? 1 : + (A->COO > B->COO) ? 1 : (A->COO < B->COO) ? -1 : - 0 ; + 0; } #ifdef CIDER int -BindCompareKLUforCIDER (const void *a, const void *b) +BindCompareKLUforCIDER(const void* a, const void* b) { - BindElementKLUforCIDER *A = (BindElementKLUforCIDER *) a ; - BindElementKLUforCIDER *B = (BindElementKLUforCIDER *) b ; + BindElementKLUforCIDER* A = (BindElementKLUforCIDER*)a; + BindElementKLUforCIDER* B = (BindElementKLUforCIDER*)b; return - (A->COO > B->COO) ? 1 : + (A->COO > B->COO) ? 1 : (A->COO < B->COO) ? -1 : - 0 ; + 0; } int -BindKluCompareCSCKLUforCIDER (const void *a, const void *b) +BindKluCompareCSCKLUforCIDER(const void* a, const void* b) { - BindElementKLUforCIDER *A = (BindElementKLUforCIDER *) a ; - BindElementKLUforCIDER *B = (BindElementKLUforCIDER *) b ; + BindElementKLUforCIDER* A = (BindElementKLUforCIDER*)a; + BindElementKLUforCIDER* B = (BindElementKLUforCIDER*)b; return - (A->CSC_Complex > B->CSC_Complex) ? 1 : + (A->CSC_Complex > B->CSC_Complex) ? 1 : (A->CSC_Complex < B->CSC_Complex) ? -1 : - 0 ; + 0; } #endif -void SMPconvertCOOtoCSC (SMPmatrix *Matrix) +void SMPconvertCOOtoCSC(SMPmatrix* Matrix) { - Element *MatrixCOO ; - KluLinkedListCOO *current, *temp ; - unsigned int *Ap_COO, current_group, i, j ; + Element* MatrixCOO; + KluLinkedListCOO* current, * temp; + unsigned int* Ap_COO, current_group, i, j; if (Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ == 0) { /* Assign N and NZ */ - Matrix->SMPkluMatrix->KLUmatrixN = 0 ; - Matrix->SMPkluMatrix->KLUmatrixNrhs = 0 ; - Matrix->SMPkluMatrix->KLUmatrixNZ = 0 ; + Matrix->SMPkluMatrix->KLUmatrixN = 0; + Matrix->SMPkluMatrix->KLUmatrixNrhs = 0; + Matrix->SMPkluMatrix->KLUmatrixNZ = 0; /* Allocate Diag Gmin CSC Vector */ - Matrix->SMPkluMatrix->KLUmatrixDiag = NULL ; + Matrix->SMPkluMatrix->KLUmatrixDiag = NULL; /* Allocate the temporary COO Column Index */ - Ap_COO = NULL ; + Ap_COO = NULL; /* Allocate the needed KLU data structures */ - Matrix->SMPkluMatrix->KLUmatrixAp = (int *) malloc (sizeof (int)) ; - Matrix->SMPkluMatrix->KLUmatrixAi = NULL ; - Matrix->SMPkluMatrix->KLUmatrixBindStructCOO = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAx = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL ; - Matrix->SMPkluMatrix->KLUmatrixIntermediate = NULL ; - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL ; + Matrix->SMPkluMatrix->KLUmatrixAp = (int*)malloc(sizeof(int)); + Matrix->SMPkluMatrix->KLUmatrixAi = NULL; + Matrix->SMPkluMatrix->KLUmatrixBindStructCOO = NULL; + Matrix->SMPkluMatrix->KLUmatrixAx = NULL; + Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL; + Matrix->SMPkluMatrix->KLUmatrixIntermediate = NULL; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL; /* Set the Matrix as Real */ - Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUmatrixReal ; + Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUmatrixReal; - return ; + return; } /* Allocate the compressed COO elements */ - MatrixCOO = (Element *) malloc (Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ * sizeof (Element)) ; + MatrixCOO = (Element*)malloc(Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ * sizeof(Element)); /* Populate the compressed COO elements and COO value of Binding Table */ /* Delete the Linked List in the meantime */ - i = 0 ; - temp = Matrix->SMPkluMatrix->KLUmatrixLinkedListCOO ; + i = 0; + temp = Matrix->SMPkluMatrix->KLUmatrixLinkedListCOO; while (temp != NULL) { - MatrixCOO [i].row = temp->row ; - MatrixCOO [i].col = temp->col ; - MatrixCOO [i].pointer = temp->pointer ; - MatrixCOO [i].group = 0 ; - current = temp ; - temp = temp->next ; - free (current->pointer) ; // We need only the memory address, we don't need to access it - free (current) ; - current = NULL ; - i++ ; + MatrixCOO[i].row = temp->row; + MatrixCOO[i].col = temp->col; + MatrixCOO[i].pointer = temp->pointer; + MatrixCOO[i].group = 0; + current = temp; + temp = temp->next; + free(current->pointer); // We need only the memory address, we don't need to access it + free(current); + current = NULL; + i++; } /* Order the MatrixCOO along the columns */ - qsort (MatrixCOO, Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ, sizeof (Element), CompareColumn) ; + qsort(MatrixCOO, Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ, sizeof(Element), CompareColumn); /* Order the MatrixCOO along the rows */ - i = 0 ; + i = 0; while (i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ) { /* Look for the next column */ - for (j = i + 1 ; j < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ ; j++) + for (j = i + 1; j < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ; j++) { - if (MatrixCOO [j].col != MatrixCOO [i].col) + if (MatrixCOO[j].col != MatrixCOO[i].col) { - break ; + break; } } - qsort (MatrixCOO + i, j - i, sizeof (Element), CompareRow) ; + qsort(MatrixCOO + i, j - i, sizeof(Element), CompareRow); - i = j ; + i = j; } /* Look for columns with all structural zeroes, not numerical zeroes - Circuits Matrix are NxN, so checking for columns is sufficient */ - unsigned int n = MatrixCOO [Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1].col + 1 ; -// Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew = (unsigned int *) malloc ((n + 1) * sizeof (unsigned int)) ; - Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld = (unsigned int *) malloc ((n + 1) * sizeof (unsigned int)) ; - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ ; i++) { + unsigned int n = MatrixCOO[Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1].col + 1; + // Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew = (unsigned int *) malloc ((n + 1) * sizeof (unsigned int)) ; + Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld = (unsigned int*)malloc((n + 1) * sizeof(unsigned int)); + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ; i++) { // Initialization // Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew [MatrixCOO [i].col] = MatrixCOO [i].col ; - Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld [MatrixCOO [i].col] = MatrixCOO [i].col ; + Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld[MatrixCOO[i].col] = MatrixCOO[i].col; } -// Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew [n] = n ; - Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld [n] = n ; + // Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew [n] = n ; + Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld[n] = n; - unsigned int search_index = 0 ; + unsigned int search_index = 0; while (search_index < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1) { - int col_index = -1 ; - unsigned int col_diff = 0 ; - for (i = search_index, j = search_index + 1 ; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1 ; i++, j++) + int col_index = -1; + unsigned int col_diff = 0; + for (i = search_index, j = search_index + 1; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1; i++, j++) { - col_diff = MatrixCOO [j].col - MatrixCOO [i].col ; + col_diff = MatrixCOO[j].col - MatrixCOO[i].col; if (col_diff > 1) { - col_index = (int)(MatrixCOO [i].col) ; - break ; + col_index = (int)(MatrixCOO[i].col); + break; } } - search_index = i ; + search_index = i; /* If col_index != -1 --> Row/Col elimination of all nodes greater than this one - compact by col_diff */ if (col_index != -1) { - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ ; i++) + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ; i++) { - if (MatrixCOO [i].col > (unsigned int)col_index) { -// Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew [MatrixCOO [i].col] = MatrixCOO [i].col - col_diff + 1 ; - Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld [MatrixCOO [i].col - col_diff + 1] = MatrixCOO [i].col ; - MatrixCOO [i].col = MatrixCOO [i].col - col_diff + 1 ; + if (MatrixCOO[i].col > (unsigned int)col_index) { + // Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew [MatrixCOO [i].col] = MatrixCOO [i].col - col_diff + 1 ; + Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld[MatrixCOO[i].col - col_diff + 1] = MatrixCOO[i].col; + MatrixCOO[i].col = MatrixCOO[i].col - col_diff + 1; } - if (MatrixCOO [i].row > (unsigned int)col_index) { - MatrixCOO [i].row = MatrixCOO [i].row - col_diff + 1 ; + if (MatrixCOO[i].row > (unsigned int)col_index) { + MatrixCOO[i].row = MatrixCOO[i].row - col_diff + 1; } } } } /* Assign labels to avoid duplicates */ - for (i = 0, j = 1 ; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1 ; i++, j++) { - if ((MatrixCOO [i].col == MatrixCOO [j].col) && (MatrixCOO [i].row == MatrixCOO [j].row)) { + for (i = 0, j = 1; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1; i++, j++) { + if ((MatrixCOO[i].col == MatrixCOO[j].col) && (MatrixCOO[i].row == MatrixCOO[j].row)) { // If col and row are the same - MatrixCOO [j].group = MatrixCOO [i].group ; - } else if ((MatrixCOO [i].col != MatrixCOO [j].col) || (MatrixCOO [i].row != MatrixCOO [j].row)) { + MatrixCOO[j].group = MatrixCOO[i].group; + } + else if ((MatrixCOO[i].col != MatrixCOO[j].col) || (MatrixCOO[i].row != MatrixCOO[j].row)) { // If or col either row are different, it isn't a duplicate, so assign the next label and store it in 'nz' - MatrixCOO [j].group = MatrixCOO [i].group + 1 ; - } else { - printf ("Error: Strange behavior during label assignment\n") ; + MatrixCOO[j].group = MatrixCOO[i].group + 1; + } + else { + printf("Error: Strange behavior during label assignment\n"); } } /* Assign N and NZ */ - Matrix->SMPkluMatrix->KLUmatrixN = (unsigned int)MatrixCOO [Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1].col + 1 ; - Matrix->SMPkluMatrix->KLUmatrixNrhs = Matrix->SMPkluMatrix->KLUmatrixN + 1 ; - Matrix->SMPkluMatrix->KLUmatrixNZ = MatrixCOO [Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1].group + 1 ; + Matrix->SMPkluMatrix->KLUmatrixN = (unsigned int)MatrixCOO[Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1].col + 1; + Matrix->SMPkluMatrix->KLUmatrixNrhs = Matrix->SMPkluMatrix->KLUmatrixN + 1; + Matrix->SMPkluMatrix->KLUmatrixNZ = MatrixCOO[Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ - 1].group + 1; /* Allocate Diag Gmin CSC Vector */ - Matrix->SMPkluMatrix->KLUmatrixDiag = (double **) malloc (Matrix->SMPkluMatrix->KLUmatrixN * sizeof (double *)) ; - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) { - Matrix->SMPkluMatrix->KLUmatrixDiag [i] = NULL ; + Matrix->SMPkluMatrix->KLUmatrixDiag = (double**)malloc(Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double*)); + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { + Matrix->SMPkluMatrix->KLUmatrixDiag[i] = NULL; } /* Allocate the temporary COO Column Index */ - Ap_COO = (unsigned int *) malloc (Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof (unsigned int)) ; + Ap_COO = (unsigned int*)malloc(Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(unsigned int)); /* Allocate the needed KLU data structures */ - Matrix->SMPkluMatrix->KLUmatrixAp = (int *) malloc ((Matrix->SMPkluMatrix->KLUmatrixN + 1) * sizeof (int)) ; - Matrix->SMPkluMatrix->KLUmatrixAi = (int *) malloc (Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof (int)) ; - Matrix->SMPkluMatrix->KLUmatrixBindStructCOO = (BindElement *) malloc (Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ * sizeof (BindElement)) ; - Matrix->SMPkluMatrix->KLUmatrixAx = (double *) malloc (Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof (double)) ; - Matrix->SMPkluMatrix->KLUmatrixAxComplex = (double *) malloc (2 * Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof (double)) ; - Matrix->SMPkluMatrix->KLUmatrixIntermediate = (double *) malloc (Matrix->SMPkluMatrix->KLUmatrixN * sizeof (double)) ; - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = (double *) malloc (2 * Matrix->SMPkluMatrix->KLUmatrixN * sizeof (double)) ; + Matrix->SMPkluMatrix->KLUmatrixAp = (int*)malloc((Matrix->SMPkluMatrix->KLUmatrixN + 1) * sizeof(int)); + Matrix->SMPkluMatrix->KLUmatrixAi = (int*)malloc(Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(int)); + Matrix->SMPkluMatrix->KLUmatrixBindStructCOO = (BindElement*)malloc(Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ * sizeof(BindElement)); + Matrix->SMPkluMatrix->KLUmatrixAx = (double*)malloc(Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(double)); + Matrix->SMPkluMatrix->KLUmatrixAxComplex = (double*)malloc(2 * Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(double)); + Matrix->SMPkluMatrix->KLUmatrixIntermediate = (double*)malloc(Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)); + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = (double*)malloc(2 * Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)); /* Copy back the Matrix in partial CSC */ - for (i = 0, current_group = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ ; i++) + for (i = 0, current_group = 0; i < Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ; i++) { - if (MatrixCOO [i].group > current_group) { - current_group = MatrixCOO [i].group ; + if (MatrixCOO[i].group > current_group) { + current_group = MatrixCOO[i].group; } - Ap_COO [current_group] = MatrixCOO [i].col ; - Matrix->SMPkluMatrix->KLUmatrixAi [current_group] = (int)MatrixCOO [i].row ; - Matrix->SMPkluMatrix->KLUmatrixBindStructCOO [i].COO = MatrixCOO [i].pointer ; - Matrix->SMPkluMatrix->KLUmatrixBindStructCOO [i].CSC = &(Matrix->SMPkluMatrix->KLUmatrixAx [current_group]) ; - Matrix->SMPkluMatrix->KLUmatrixBindStructCOO [i].CSC_Complex = &(Matrix->SMPkluMatrix->KLUmatrixAxComplex [2 * current_group]) ; - if (MatrixCOO [i].col == MatrixCOO [i].row) { - Matrix->SMPkluMatrix->KLUmatrixDiag [MatrixCOO [i].col] = Matrix->SMPkluMatrix->KLUmatrixBindStructCOO [i].CSC ; + Ap_COO[current_group] = MatrixCOO[i].col; + Matrix->SMPkluMatrix->KLUmatrixAi[current_group] = (int)MatrixCOO[i].row; + Matrix->SMPkluMatrix->KLUmatrixBindStructCOO[i].COO = MatrixCOO[i].pointer; + Matrix->SMPkluMatrix->KLUmatrixBindStructCOO[i].CSC = &(Matrix->SMPkluMatrix->KLUmatrixAx[current_group]); + Matrix->SMPkluMatrix->KLUmatrixBindStructCOO[i].CSC_Complex = &(Matrix->SMPkluMatrix->KLUmatrixAxComplex[2 * current_group]); + if (MatrixCOO[i].col == MatrixCOO[i].row) { + Matrix->SMPkluMatrix->KLUmatrixDiag[MatrixCOO[i].col] = Matrix->SMPkluMatrix->KLUmatrixBindStructCOO[i].CSC; } } /* Compress the COO Column Index to CSC Column Index */ - Compress (Ap_COO, (unsigned int *)Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixNZ) ; + Compress(Ap_COO, (unsigned int*)Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixNZ); /* Free the temporary stuff */ - free (Ap_COO) ; - free (MatrixCOO) ; + free(Ap_COO); + free(MatrixCOO); /* Sort the Binding Table */ - qsort (Matrix->SMPkluMatrix->KLUmatrixBindStructCOO, Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ, sizeof (BindElement), BindCompare) ; + qsort(Matrix->SMPkluMatrix->KLUmatrixBindStructCOO, Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ, sizeof(BindElement), BindCompare); /* Set the Matrix as Real */ - Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUmatrixReal ; + Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUmatrixReal; - return ; + return; } #ifdef CIDER typedef struct sElementKLUforCIDER { - unsigned int row ; - unsigned int col ; - double *pointer ; -} ElementKLUforCIDER ; + unsigned int row; + unsigned int col; + double* pointer; +} ElementKLUforCIDER; -void SMPconvertCOOtoCSCKLUforCIDER (SMPmatrix *Matrix) +void SMPconvertCOOtoCSCKLUforCIDER(SMPmatrix* Matrix) { - ElementKLUforCIDER *MatrixCOO ; - unsigned int *Ap_COO, i, j, nz ; + ElementKLUforCIDER* MatrixCOO; + unsigned int* Ap_COO, i, j, nz; /* Count the non-zero elements and store it */ - nz = 0 ; - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN ; i++) { - if ((Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER [i] != -1) && (Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER [i] != -1)) { - nz++ ; + nz = 0; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN; i++) { + if ((Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER[i] != -1) && (Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER[i] != -1)) { + nz++; } } - Matrix->SMPkluMatrix->KLUmatrixNZ = nz ; + Matrix->SMPkluMatrix->KLUmatrixNZ = nz; /* Allocate the compressed COO elements */ - MatrixCOO = (ElementKLUforCIDER *) malloc (nz * sizeof (ElementKLUforCIDER)) ; + MatrixCOO = (ElementKLUforCIDER*)malloc(nz * sizeof(ElementKLUforCIDER)); /* Allocate the temporary COO Column Index */ - Ap_COO = (unsigned int *) malloc (nz * sizeof (unsigned int)) ; + Ap_COO = (unsigned int*)malloc(nz * sizeof(unsigned int)); /* Allocate the needed KLU data structures */ - Matrix->SMPkluMatrix->KLUmatrixAp = (int *) malloc ((Matrix->SMPkluMatrix->KLUmatrixN + 1) * sizeof (int)) ; - Matrix->SMPkluMatrix->KLUmatrixAi = (int *) malloc (nz * sizeof (int)) ; - Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER = (BindElementKLUforCIDER *) malloc (nz * sizeof (BindElementKLUforCIDER)) ; - Matrix->SMPkluMatrix->KLUmatrixAxComplex = (double *) malloc (2 * nz * sizeof (double)) ; - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = (double *) malloc (2 * Matrix->SMPkluMatrix->KLUmatrixN * sizeof (double)) ; + Matrix->SMPkluMatrix->KLUmatrixAp = (int*)malloc((Matrix->SMPkluMatrix->KLUmatrixN + 1) * sizeof(int)); + Matrix->SMPkluMatrix->KLUmatrixAi = (int*)malloc(nz * sizeof(int)); + Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER = (BindElementKLUforCIDER*)malloc(nz * sizeof(BindElementKLUforCIDER)); + Matrix->SMPkluMatrix->KLUmatrixAxComplex = (double*)malloc(2 * nz * sizeof(double)); + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = (double*)malloc(2 * Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)); /* Populate the compressed COO elements and COO value of Binding Table */ - j = 0 ; - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN ; i++) { - if ((Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER [i] != -1) && (Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER [i] != -1)) { - MatrixCOO [j].row = (unsigned int)Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER [i] ; - MatrixCOO [j].col = (unsigned int)Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER [i] ; - MatrixCOO [j].pointer = &(Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER [2 * i]) ; - j++ ; + j = 0; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN; i++) { + if ((Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER[i] != -1) && (Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER[i] != -1)) { + MatrixCOO[j].row = (unsigned int)Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER[i]; + MatrixCOO[j].col = (unsigned int)Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER[i]; + MatrixCOO[j].pointer = &(Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER[2 * i]); + j++; } } /* Order the MatrixCOO along the columns */ - qsort (MatrixCOO, nz, sizeof (ElementKLUforCIDER), CompareColumn) ; + qsort(MatrixCOO, nz, sizeof(ElementKLUforCIDER), CompareColumn); /* Order the MatrixCOO along the rows */ - i = 0 ; + i = 0; while (i < nz) { /* Look for the next column */ - for (j = i + 1 ; j < nz ; j++) + for (j = i + 1; j < nz; j++) { - if (MatrixCOO [j].col != MatrixCOO [i].col) + if (MatrixCOO[j].col != MatrixCOO[i].col) { - break ; + break; } } - qsort (MatrixCOO + i, j - i, sizeof (ElementKLUforCIDER), CompareRow) ; + qsort(MatrixCOO + i, j - i, sizeof(ElementKLUforCIDER), CompareRow); - i = j ; + i = j; } /* Copy back the Matrix in partial CSC */ - for (i = 0 ; i < nz ; i++) + for (i = 0; i < nz; i++) { - Ap_COO [i] = MatrixCOO [i].col ; - Matrix->SMPkluMatrix->KLUmatrixAi [i] = (int)MatrixCOO [i].row ; - Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER [i].COO = MatrixCOO [i].pointer ; - Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER [i].CSC_Complex = &(Matrix->SMPkluMatrix->KLUmatrixAxComplex [2 * i]) ; + Ap_COO[i] = MatrixCOO[i].col; + Matrix->SMPkluMatrix->KLUmatrixAi[i] = (int)MatrixCOO[i].row; + Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER[i].COO = MatrixCOO[i].pointer; + Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER[i].CSC_Complex = &(Matrix->SMPkluMatrix->KLUmatrixAxComplex[2 * i]); } /* Compress the COO Column Index to CSC Column Index */ - Compress (Ap_COO, (unsigned int *)Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixN, nz) ; + Compress(Ap_COO, (unsigned int*)Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixN, nz); /* Free the temporary stuff */ - free (Ap_COO) ; - free (MatrixCOO) ; + free(Ap_COO); + free(MatrixCOO); /* Sort the Binding Table */ - qsort (Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER, nz, sizeof (BindElementKLUforCIDER), BindCompareKLUforCIDER) ; + qsort(Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER, nz, sizeof(BindElementKLUforCIDER), BindCompareKLUforCIDER); - return ; + return; } #endif /* * SMPmakeElt() */ -double * -SMPmakeElt (SMPmatrix *Matrix, int Row, int Col) +double* +SMPmakeElt(SMPmatrix* Matrix, int Row, int Col) { - KluLinkedListCOO *temp ; + KluLinkedListCOO* temp; if (Matrix->CKTkluMODE) { if ((Row > 0) && (Col > 0)) { - Row = Row - 1 ; - Col = Col - 1 ; - temp = (KluLinkedListCOO *) malloc (sizeof (KluLinkedListCOO)) ; - temp->row = (unsigned int)Row ; - temp->col = (unsigned int)Col ; - temp->pointer = (double *) malloc (sizeof (double)) ; - temp->next = Matrix->SMPkluMatrix->KLUmatrixLinkedListCOO ; - Matrix->SMPkluMatrix->KLUmatrixLinkedListCOO = temp ; - Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ++ ; - return temp->pointer ; - } else { - return Matrix->SMPkluMatrix->KLUmatrixTrashCOO ; + Row = Row - 1; + Col = Col - 1; + temp = (KluLinkedListCOO*)malloc(sizeof(KluLinkedListCOO)); + temp->row = (unsigned int)Row; + temp->col = (unsigned int)Col; + temp->pointer = (double*)malloc(sizeof(double)); + temp->next = Matrix->SMPkluMatrix->KLUmatrixLinkedListCOO; + Matrix->SMPkluMatrix->KLUmatrixLinkedListCOO = temp; + Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ++; + return temp->pointer; } - } else { - return spGetElement (Matrix->SPmatrix, Row, Col) ; + else { + return Matrix->SMPkluMatrix->KLUmatrixTrashCOO; + } + } + else { + return spGetElement(Matrix->SPmatrix, Row, Col); } } #ifdef CIDER -double * -SMPmakeEltKLUforCIDER (SMPmatrix *Matrix, int Row, int Col) +double* +SMPmakeEltKLUforCIDER(SMPmatrix* Matrix, int Row, int Col) { if (Matrix->CKTkluMODE) { if ((Row > 0) && (Col > 0)) { - Row = Row - 1 ; - Col = Col - 1 ; - Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER [Row * (int)Matrix->SMPkluMatrix->KLUmatrixN + Col] = Row ; - Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER [Row * (int)Matrix->SMPkluMatrix->KLUmatrixN + Col] = Col ; - return &(Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER [2 * (Row * (int)Matrix->SMPkluMatrix->KLUmatrixN + Col)]) ; - } else { - return Matrix->SMPkluMatrix->KLUmatrixTrashCOO ; + Row = Row - 1; + Col = Col - 1; + Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER[Row * (int)Matrix->SMPkluMatrix->KLUmatrixN + Col] = Row; + Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER[Row * (int)Matrix->SMPkluMatrix->KLUmatrixN + Col] = Col; + return &(Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER[2 * (Row * (int)Matrix->SMPkluMatrix->KLUmatrixN + Col)]); } - } else { - return spGetElement (Matrix->SPmatrix, Row, Col) ; + else { + return Matrix->SMPkluMatrix->KLUmatrixTrashCOO; + } + } + else { + return spGetElement(Matrix->SPmatrix, Row, Col); } } #endif @@ -464,17 +470,18 @@ SMPmakeEltKLUforCIDER (SMPmatrix *Matrix, int Row, int Col) */ void -SMPcClear (SMPmatrix *Matrix) +SMPcClear(SMPmatrix* Matrix) { - unsigned int i ; + unsigned int i; if (Matrix->CKTkluMODE) { - for (i = 0 ; i < 2 * Matrix->SMPkluMatrix->KLUmatrixNZ ; i++) { - Matrix->SMPkluMatrix->KLUmatrixAxComplex [i] = 0 ; + for (i = 0; i < 2 * Matrix->SMPkluMatrix->KLUmatrixNZ; i++) { + Matrix->SMPkluMatrix->KLUmatrixAxComplex[i] = 0; } - } else { - spClear (Matrix->SPmatrix) ; + } + else { + spClear(Matrix->SPmatrix); } } @@ -483,28 +490,29 @@ SMPcClear (SMPmatrix *Matrix) */ void -SMPclear (SMPmatrix *Matrix) +SMPclear(SMPmatrix* Matrix) { - unsigned int i ; + unsigned int i; if (Matrix->CKTkluMODE) { - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixNZ ; i++) { - Matrix->SMPkluMatrix->KLUmatrixAx [i] = 0 ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixNZ; i++) { + Matrix->SMPkluMatrix->KLUmatrixAx[i] = 0; } - } else { - spClear (Matrix->SPmatrix) ; + } + else { + spClear(Matrix->SPmatrix); } } #ifdef CIDER void -SMPclearKLUforCIDER (SMPmatrix *Matrix) +SMPclearKLUforCIDER(SMPmatrix* Matrix) { - unsigned int i ; + unsigned int i; - for (i = 0 ; i < 2 * Matrix->SMPkluMatrix->KLUmatrixNZ ; i++) { - Matrix->SMPkluMatrix->KLUmatrixAxComplex [i] = 0 ; + for (i = 0; i < 2 * Matrix->SMPkluMatrix->KLUmatrixNZ; i++) { + Matrix->SMPkluMatrix->KLUmatrixAxComplex[i] = 0; } } #endif @@ -514,30 +522,30 @@ SMPclearKLUforCIDER (SMPmatrix *Matrix) /* * SMPcLUfac() */ -/*ARGSUSED*/ + /*ARGSUSED*/ int -SMPcLUfac (SMPmatrix *Matrix, double PivTol) +SMPcLUfac(SMPmatrix* Matrix, double PivTol) { - int ret ; + int ret; - NG_IGNORE (PivTol) ; + NG_IGNORE(PivTol); if (Matrix->CKTkluMODE) { if (CircuitIsDigital() && Matrix->SMPkluMatrix->KLUmatrixN == 0) { - // XSPICE pure digital circuits produce empty KLU matrix - return 0 ; + // XSPICE pure digital circuits produce empty KLU matrix + return 0; } - ret = klu_z_refactor (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, - Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + ret = klu_z_refactor(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, + Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixCommon); if (ret == 0) { if (Matrix->SMPkluMatrix->KLUmatrixCommon == NULL) { fprintf(stderr, "Error (ReFactor Complex): KLUcommon object is NULL. A problem occurred\n"); - return 0 ; + return 0; } if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_SINGULAR) { if (ft_ngdebug) { @@ -545,51 +553,53 @@ SMPcLUfac (SMPmatrix *Matrix, double PivTol) fprintf(stderr, " Numerical Rank: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->numerical_rank); fprintf(stderr, " Singular Node: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1); } - return E_SINGULAR ; + return E_SINGULAR; } if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_EMPTY_MATRIX) { - fprintf (stderr, "Error (ReFactor Complex): KLU Matrix is empty\n") ; - return 0 ; + fprintf(stderr, "Error (ReFactor Complex): KLU Matrix is empty\n"); + return 0; } if (Matrix->SMPkluMatrix->KLUmatrixNumeric == NULL) { - fprintf (stderr, "Error (ReFactor Complex): KLUnumeric object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (ReFactor Complex): KLUnumeric object is NULL. A problem occurred\n"); } - return 1 ; - } else { - return 0 ; + return 1; } - } else { - spSetComplex (Matrix->SPmatrix) ; - return spFactor (Matrix->SPmatrix) ; + else { + return 0; + } + } + else { + spSetComplex(Matrix->SPmatrix); + return spFactor(Matrix->SPmatrix); } } /* * SMPluFac() */ -/*ARGSUSED*/ + /*ARGSUSED*/ int -SMPluFac (SMPmatrix *Matrix, double PivTol, double Gmin) +SMPluFac(SMPmatrix* Matrix, double PivTol, double Gmin) { - int ret ; + int ret; - NG_IGNORE (PivTol) ; + NG_IGNORE(PivTol); if (Matrix->CKTkluMODE) { if (CircuitIsDigital() && Matrix->SMPkluMatrix->KLUmatrixN == 0) { - // XSPICE pure digital circuits produce empty KLU matrix - return 0 ; + // XSPICE pure digital circuits produce empty KLU matrix + return 0; } if (Matrix->SMPkluMatrix->KLUloadDiagGmin) { - LoadGmin_CSC (Matrix->SMPkluMatrix->KLUmatrixDiag, Matrix->SMPkluMatrix->KLUmatrixN, Gmin) ; + LoadGmin_CSC(Matrix->SMPkluMatrix->KLUmatrixDiag, Matrix->SMPkluMatrix->KLUmatrixN, Gmin); } - ret = klu_refactor (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAx, - Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + ret = klu_refactor(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAx, + Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixCommon); if (ret == 0) { @@ -599,63 +609,66 @@ SMPluFac (SMPmatrix *Matrix, double PivTol, double Gmin) fprintf(stderr, " Numerical Rank: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->numerical_rank); fprintf(stderr, " Singular Node: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1); } - return E_SINGULAR ; + return E_SINGULAR; } if (Matrix->SMPkluMatrix->KLUmatrixCommon == NULL) { - fprintf (stderr, "Error (ReFactor): KLUcommon object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (ReFactor): KLUcommon object is NULL. A problem occurred\n"); } if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_EMPTY_MATRIX) { - fprintf (stderr, "Error (ReFactor): KLU Matrix is empty\n") ; - return 0 ; + fprintf(stderr, "Error (ReFactor): KLU Matrix is empty\n"); + return 0; } if (Matrix->SMPkluMatrix->KLUmatrixNumeric == NULL) { - fprintf (stderr, "Error (ReFactor): KLUnumeric object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (ReFactor): KLUnumeric object is NULL. A problem occurred\n"); } - return 1 ; - } else { - return 0 ; + return 1; } - } else { - spSetReal (Matrix->SPmatrix) ; - LoadGmin (Matrix, Gmin) ; - return spFactor (Matrix->SPmatrix) ; + else { + return 0; + } + } + else { + spSetReal(Matrix->SPmatrix); + LoadGmin(Matrix, Gmin); + return spFactor(Matrix->SPmatrix); } } #ifdef CIDER int -SMPluFacKLUforCIDER (SMPmatrix *Matrix) +SMPluFacKLUforCIDER(SMPmatrix* Matrix) { - unsigned int i ; - double *KLUmatrixAx ; - int ret ; + unsigned int i; + double* KLUmatrixAx; + int ret; if (Matrix->CKTkluMODE) { if (CircuitIsDigital() && Matrix->SMPkluMatrix->KLUmatrixN == 0) { - // XSPICE pure digital circuits produce empty KLU matrix - return 0 ; + // XSPICE pure digital circuits produce empty KLU matrix + return 0; } if (Matrix->SMPkluMatrix->KLUmatrixIsComplex) { - ret = klu_z_refactor (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, - Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixCommon) ; - } else { + ret = klu_z_refactor(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, + Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixCommon); + } + else { /* Allocate the Real Matrix */ - KLUmatrixAx = (double *) malloc (Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(double)) ; + KLUmatrixAx = (double*)malloc(Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(double)); /* Copy the Complex Matrix into the Real Matrix */ - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixNZ ; i++) { - KLUmatrixAx [i] = Matrix->SMPkluMatrix->KLUmatrixAxComplex [2 * i] ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixNZ; i++) { + KLUmatrixAx[i] = Matrix->SMPkluMatrix->KLUmatrixAxComplex[2 * i]; } /* Re-Factor the Real Matrix */ - ret = klu_refactor (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, KLUmatrixAx, - Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + ret = klu_refactor(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, KLUmatrixAx, + Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixCommon); /* Free the Real Matrix Storage */ - free (KLUmatrixAx) ; + free(KLUmatrixAx); } if (ret == 0) @@ -666,25 +679,27 @@ SMPluFacKLUforCIDER (SMPmatrix *Matrix) fprintf(stderr, " Numerical Rank: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->numerical_rank); fprintf(stderr, " Singular Node: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1); } - return E_SINGULAR ; + return E_SINGULAR; } if (Matrix->SMPkluMatrix->KLUmatrixCommon == NULL) { - fprintf (stderr, "Error (ReFactor for CIDER): KLUcommon object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (ReFactor for CIDER): KLUcommon object is NULL. A problem occurred\n"); } if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_EMPTY_MATRIX) { - fprintf (stderr, "Error (ReFactor for CIDER): KLU Matrix is empty\n") ; - return 0 ; + fprintf(stderr, "Error (ReFactor for CIDER): KLU Matrix is empty\n"); + return 0; } if (Matrix->SMPkluMatrix->KLUmatrixNumeric == NULL) { - fprintf (stderr, "Error (ReFactor for CIDER): KLUnumeric object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (ReFactor for CIDER): KLUnumeric object is NULL. A problem occurred\n"); } - return 1 ; - } else { - return 0 ; + return 1; } - } else { - return spFactor (Matrix->SPmatrix) ; + else { + return 0; + } + } + else { + return spFactor(Matrix->SPmatrix); } } #endif @@ -694,23 +709,23 @@ SMPluFacKLUforCIDER (SMPmatrix *Matrix) */ int -SMPcReorder (SMPmatrix *Matrix, double PivTol, double PivRel, int *NumSwaps) +SMPcReorder(SMPmatrix* Matrix, double PivTol, double PivRel, int* NumSwaps) { if (Matrix->CKTkluMODE) { if (CircuitIsDigital() && Matrix->SMPkluMatrix->KLUmatrixN == 0) { - // XSPICE pure digital circuits produce empty KLU matrix - return 0 ; + // XSPICE pure digital circuits produce empty KLU matrix + return 0; } - Matrix->SMPkluMatrix->KLUmatrixCommon->tol = PivRel ; + Matrix->SMPkluMatrix->KLUmatrixCommon->tol = PivRel; if (Matrix->SMPkluMatrix->KLUmatrixNumeric != NULL) { - klu_free_numeric (&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon) ; + klu_free_numeric(&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon); } - Matrix->SMPkluMatrix->KLUmatrixNumeric = klu_z_factor (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, - Matrix->SMPkluMatrix->KLUmatrixAxComplex, Matrix->SMPkluMatrix->KLUmatrixSymbolic, - Matrix->SMPkluMatrix->KLUmatrixCommon) ; + Matrix->SMPkluMatrix->KLUmatrixNumeric = klu_z_factor(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, + Matrix->SMPkluMatrix->KLUmatrixAxComplex, Matrix->SMPkluMatrix->KLUmatrixSymbolic, + Matrix->SMPkluMatrix->KLUmatrixCommon); if (Matrix->SMPkluMatrix->KLUmatrixNumeric == NULL) { @@ -720,28 +735,30 @@ SMPcReorder (SMPmatrix *Matrix, double PivTol, double PivRel, int *NumSwaps) fprintf(stderr, " Numerical Rank: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->numerical_rank); fprintf(stderr, " Singular Node: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1); } - return E_SINGULAR ; + return E_SINGULAR; } if (Matrix->SMPkluMatrix->KLUmatrixCommon == NULL) { - fprintf (stderr, "Error (Factor Complex): KLUnumeric object is NULL. A problem occurred\n") ; - fprintf (stderr, "Error (Factor Complex): KLUcommon object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (Factor Complex): KLUnumeric object is NULL. A problem occurred\n"); + fprintf(stderr, "Error (Factor Complex): KLUcommon object is NULL. A problem occurred\n"); } if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_EMPTY_MATRIX) { - fprintf (stderr, "Error (Factor Complex): KLU Matrix is empty\n") ; - return 0 ; + fprintf(stderr, "Error (Factor Complex): KLU Matrix is empty\n"); + return 0; } if (Matrix->SMPkluMatrix->KLUmatrixSymbolic == NULL) { - fprintf (stderr, "Error (Factor Complex): KLUnumeric object is NULL. A problem occurred\n") ; - fprintf (stderr, "Error (Factor Complex): KLUsymbolic object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (Factor Complex): KLUnumeric object is NULL. A problem occurred\n"); + fprintf(stderr, "Error (Factor Complex): KLUsymbolic object is NULL. A problem occurred\n"); } - return 1 ; - } else { - return 0 ; + return 1; } - } else { - *NumSwaps = 1 ; - spSetComplex (Matrix->SPmatrix) ; - return spOrderAndFactor (Matrix->SPmatrix, NULL, (spREAL)PivRel, (spREAL)PivTol, YES) ; + else { + return 0; + } + } + else { + *NumSwaps = 1; + spSetComplex(Matrix->SPmatrix); + return spOrderAndFactor(Matrix->SPmatrix, NULL, (spREAL)PivRel, (spREAL)PivTol, YES); } } @@ -750,27 +767,27 @@ SMPcReorder (SMPmatrix *Matrix, double PivTol, double PivRel, int *NumSwaps) */ int -SMPreorder (SMPmatrix *Matrix, double PivTol, double PivRel, double Gmin) +SMPreorder(SMPmatrix* Matrix, double PivTol, double PivRel, double Gmin) { if (Matrix->CKTkluMODE) { if (CircuitIsDigital() && Matrix->SMPkluMatrix->KLUmatrixN == 0) { - // XSPICE pure digital circuits produce empty KLU matrix - return 0 ; + // XSPICE pure digital circuits produce empty KLU matrix + return 0; } if (Matrix->SMPkluMatrix->KLUloadDiagGmin) { - LoadGmin_CSC (Matrix->SMPkluMatrix->KLUmatrixDiag, Matrix->SMPkluMatrix->KLUmatrixN, Gmin) ; + LoadGmin_CSC(Matrix->SMPkluMatrix->KLUmatrixDiag, Matrix->SMPkluMatrix->KLUmatrixN, Gmin); } - Matrix->SMPkluMatrix->KLUmatrixCommon->tol = PivRel ; + Matrix->SMPkluMatrix->KLUmatrixCommon->tol = PivRel; if (Matrix->SMPkluMatrix->KLUmatrixNumeric != NULL) { - klu_free_numeric (&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon) ; + klu_free_numeric(&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon); } - Matrix->SMPkluMatrix->KLUmatrixNumeric = klu_factor (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, - Matrix->SMPkluMatrix->KLUmatrixAx, Matrix->SMPkluMatrix->KLUmatrixSymbolic, - Matrix->SMPkluMatrix->KLUmatrixCommon) ; + Matrix->SMPkluMatrix->KLUmatrixNumeric = klu_factor(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, + Matrix->SMPkluMatrix->KLUmatrixAx, Matrix->SMPkluMatrix->KLUmatrixSymbolic, + Matrix->SMPkluMatrix->KLUmatrixCommon); if (Matrix->SMPkluMatrix->KLUmatrixNumeric == NULL) { @@ -780,68 +797,71 @@ SMPreorder (SMPmatrix *Matrix, double PivTol, double PivRel, double Gmin) fprintf(stderr, " Numerical Rank: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->numerical_rank); fprintf(stderr, " Singular Node: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1); } - return E_SINGULAR ; + return E_SINGULAR; } if (Matrix->SMPkluMatrix->KLUmatrixCommon == NULL) { - fprintf (stderr, "Error (Factor): KLUnumeric object is NULL. A problem occurred\n") ; - fprintf (stderr, "Error (Factor): KLUcommon object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (Factor): KLUnumeric object is NULL. A problem occurred\n"); + fprintf(stderr, "Error (Factor): KLUcommon object is NULL. A problem occurred\n"); } if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_EMPTY_MATRIX) { - fprintf (stderr, "Error (Factor): KLU Matrix is empty\n") ; - return 0 ; + fprintf(stderr, "Error (Factor): KLU Matrix is empty\n"); + return 0; } if (Matrix->SMPkluMatrix->KLUmatrixSymbolic == NULL) { - fprintf (stderr, "Error (Factor): KLUnumeric object is NULL. A problem occurred\n") ; - fprintf (stderr, "Error (Factor): KLUsymbolic object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (Factor): KLUnumeric object is NULL. A problem occurred\n"); + fprintf(stderr, "Error (Factor): KLUsymbolic object is NULL. A problem occurred\n"); } - return 1 ; - } else { - return 0 ; + return 1; } - } else { - spSetReal (Matrix->SPmatrix) ; - LoadGmin (Matrix, Gmin) ; - return spOrderAndFactor (Matrix->SPmatrix, NULL, (spREAL)PivRel, (spREAL)PivTol, YES) ; + else { + return 0; + } + } + else { + spSetReal(Matrix->SPmatrix); + LoadGmin(Matrix, Gmin); + return spOrderAndFactor(Matrix->SPmatrix, NULL, (spREAL)PivRel, (spREAL)PivTol, YES); } } #ifdef CIDER int -SMPreorderKLUforCIDER (SMPmatrix *Matrix) +SMPreorderKLUforCIDER(SMPmatrix* Matrix) { - unsigned int i ; - double *KLUmatrixAx ; + unsigned int i; + double* KLUmatrixAx; if (Matrix->CKTkluMODE) { if (CircuitIsDigital() && Matrix->SMPkluMatrix->KLUmatrixN == 0) { - // XSPICE pure digital circuits produce empty KLU matrix - return 0 ; + // XSPICE pure digital circuits produce empty KLU matrix + return 0; } if (Matrix->SMPkluMatrix->KLUmatrixNumeric != NULL) { - klu_free_numeric (&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon) ; + klu_free_numeric(&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon); } if (Matrix->SMPkluMatrix->KLUmatrixIsComplex) { - Matrix->SMPkluMatrix->KLUmatrixNumeric = klu_z_factor (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, - Matrix->SMPkluMatrix->KLUmatrixAxComplex, Matrix->SMPkluMatrix->KLUmatrixSymbolic, - Matrix->SMPkluMatrix->KLUmatrixCommon) ; - } else { + Matrix->SMPkluMatrix->KLUmatrixNumeric = klu_z_factor(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, + Matrix->SMPkluMatrix->KLUmatrixAxComplex, Matrix->SMPkluMatrix->KLUmatrixSymbolic, + Matrix->SMPkluMatrix->KLUmatrixCommon); + } + else { /* Allocate the Real Matrix */ - KLUmatrixAx = (double *) malloc (Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(double)) ; + KLUmatrixAx = (double*)malloc(Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(double)); /* Copy the Complex Matrix into the Real Matrix */ - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixNZ ; i++) { - KLUmatrixAx [i] = Matrix->SMPkluMatrix->KLUmatrixAxComplex [2 * i] ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixNZ; i++) { + KLUmatrixAx[i] = Matrix->SMPkluMatrix->KLUmatrixAxComplex[2 * i]; } /* Factor the Real Matrix */ - Matrix->SMPkluMatrix->KLUmatrixNumeric = klu_factor (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, - KLUmatrixAx, Matrix->SMPkluMatrix->KLUmatrixSymbolic, - Matrix->SMPkluMatrix->KLUmatrixCommon) ; + Matrix->SMPkluMatrix->KLUmatrixNumeric = klu_factor(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, + KLUmatrixAx, Matrix->SMPkluMatrix->KLUmatrixSymbolic, + Matrix->SMPkluMatrix->KLUmatrixCommon); /* Free the Real Matrix Storage */ - free (KLUmatrixAx) ; + free(KLUmatrixAx); } if (Matrix->SMPkluMatrix->KLUmatrixNumeric == NULL) @@ -852,26 +872,28 @@ SMPreorderKLUforCIDER (SMPmatrix *Matrix) fprintf(stderr, " Numerical Rank: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->numerical_rank); fprintf(stderr, " Singular Node: %d\n", Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1); } - return E_SINGULAR ; + return E_SINGULAR; } if (Matrix->SMPkluMatrix->KLUmatrixCommon == NULL) { - fprintf (stderr, "Error (Factor for CIDER): KLUnumeric object is NULL. A problem occurred\n") ; - fprintf (stderr, "Error (Factor for CIDER): KLUcommon object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (Factor for CIDER): KLUnumeric object is NULL. A problem occurred\n"); + fprintf(stderr, "Error (Factor for CIDER): KLUcommon object is NULL. A problem occurred\n"); } if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_EMPTY_MATRIX) { - fprintf (stderr, "Error (Factor for CIDER): KLU Matrix is empty\n") ; - return 0 ; + fprintf(stderr, "Error (Factor for CIDER): KLU Matrix is empty\n"); + return 0; } if (Matrix->SMPkluMatrix->KLUmatrixSymbolic == NULL) { - fprintf (stderr, "Error (Factor for CIDER): KLUnumeric object is NULL. A problem occurred\n") ; - fprintf (stderr, "Error (Factor for CIDER): KLUsymbolic object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (Factor for CIDER): KLUnumeric object is NULL. A problem occurred\n"); + fprintf(stderr, "Error (Factor for CIDER): KLUsymbolic object is NULL. A problem occurred\n"); } - return 1 ; - } else { - return 0 ; + return 1; } - } else { - return spFactor (Matrix->SPmatrix) ; + else { + return 0; + } + } + else { + return spFactor(Matrix->SPmatrix); } } #endif @@ -880,32 +902,33 @@ SMPreorderKLUforCIDER (SMPmatrix *Matrix) * SMPcaSolve() */ void -SMPcaSolve (SMPmatrix *Matrix, double RHS[], double iRHS[], double Spare[], double iSpare[]) +SMPcaSolve(SMPmatrix* Matrix, double RHS[], double iRHS[], double Spare[], double iSpare[]) { - int ret ; - unsigned int i ; + int ret; + unsigned int i; - NG_IGNORE (iSpare) ; - NG_IGNORE (Spare) ; + NG_IGNORE(iSpare); + NG_IGNORE(Spare); if (Matrix->CKTkluMODE) { - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i] = RHS [i + 1] ; - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i + 1] = iRHS [i + 1] ; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i] = RHS[i + 1]; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i + 1] = iRHS[i + 1]; } - ret = klu_z_solve (Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + ret = klu_z_solve(Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex, Matrix->SMPkluMatrix->KLUmatrixCommon); - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { - RHS [i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i] ; - iRHS [i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i + 1] ; + RHS[i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i]; + iRHS[i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i + 1]; } - } else { - spSolveTransposed (Matrix->SPmatrix, RHS, RHS, iRHS, iRHS) ; + } + else { + spSolveTransposed(Matrix->SPmatrix, RHS, RHS, iRHS, iRHS); } } @@ -914,32 +937,33 @@ SMPcaSolve (SMPmatrix *Matrix, double RHS[], double iRHS[], double Spare[], doub */ void -SMPcSolve (SMPmatrix *Matrix, double RHS[], double iRHS[], double Spare[], double iSpare[]) +SMPcSolve(SMPmatrix* Matrix, double RHS[], double iRHS[], double Spare[], double iSpare[]) { - int ret ; - unsigned int i ; + int ret; + unsigned int i; - NG_IGNORE (iSpare) ; - NG_IGNORE (Spare) ; + NG_IGNORE(iSpare); + NG_IGNORE(Spare); if (Matrix->CKTkluMODE) { - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i] = RHS [i + 1] ; - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i + 1] = iRHS [i + 1] ; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i] = RHS[i + 1]; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i + 1] = iRHS[i + 1]; } - ret = klu_z_solve (Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + ret = klu_z_solve(Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex, Matrix->SMPkluMatrix->KLUmatrixCommon); - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { - RHS [i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i] ; - iRHS [i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i + 1] ; + RHS[i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i]; + iRHS[i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i + 1]; } - } else { - spSolve (Matrix->SPmatrix, RHS, RHS, iRHS, iRHS) ; + } + else { + spSolve(Matrix->SPmatrix, RHS, RHS, iRHS, iRHS); } } @@ -948,28 +972,28 @@ SMPcSolve (SMPmatrix *Matrix, double RHS[], double iRHS[], double Spare[], doubl */ void -SMPsolve (SMPmatrix *Matrix, double RHS[], double Spare[]) +SMPsolve(SMPmatrix* Matrix, double RHS[], double Spare[]) { - int ret ; - unsigned int i ; + int ret; + unsigned int i; - NG_IGNORE (Spare) ; + NG_IGNORE(Spare); if (Matrix->CKTkluMODE) { if (CircuitIsDigital() && Matrix->SMPkluMatrix->KLUmatrixN == 0) { - // XSPICE pure digital circuits produce empty KLU matrix - return ; + // XSPICE pure digital circuits produce empty KLU matrix + return; } - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) { - if (Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld [i + 1] != 0) { - Matrix->SMPkluMatrix->KLUmatrixIntermediate [i] = RHS [Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld [i + 1]] ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { + if (Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld[i + 1] != 0) { + Matrix->SMPkluMatrix->KLUmatrixIntermediate[i] = RHS[Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld[i + 1]]; } } - ret = klu_solve (Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, - Matrix->SMPkluMatrix->KLUmatrixIntermediate, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + ret = klu_solve(Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, + Matrix->SMPkluMatrix->KLUmatrixIntermediate, Matrix->SMPkluMatrix->KLUmatrixCommon); if (ret == 0) { @@ -982,81 +1006,84 @@ SMPsolve (SMPmatrix *Matrix, double RHS[], double Spare[]) /* FIXME: Do we need a 'return E_SINGULAR' here? */ } if (Matrix->SMPkluMatrix->KLUmatrixCommon == NULL) { - fprintf (stderr, "Error (Solve): KLUcommon object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (Solve): KLUcommon object is NULL. A problem occurred\n"); } if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_EMPTY_MATRIX) { - fprintf (stderr, "Error (Solve): KLU Matrix is empty\n") ; + fprintf(stderr, "Error (Solve): KLU Matrix is empty\n"); } if (Matrix->SMPkluMatrix->KLUmatrixNumeric == NULL) { - fprintf (stderr, "Error (Solve): KLUnumeric object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (Solve): KLUnumeric object is NULL. A problem occurred\n"); } if (Matrix->SMPkluMatrix->KLUmatrixSymbolic == NULL) { - fprintf (stderr, "Error (Solve): KLUsymbolic object is NULL. A problem occurred\n") ; + fprintf(stderr, "Error (Solve): KLUsymbolic object is NULL. A problem occurred\n"); } } - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixNrhs ; i++) { - RHS [i] = 0 ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixNrhs; i++) { + RHS[i] = 0; } - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) { - if (Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld [i + 1] != 0) { - RHS [Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld [i + 1]] = Matrix->SMPkluMatrix->KLUmatrixIntermediate [i] ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { + if (Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld[i + 1] != 0) { + RHS[Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld[i + 1]] = Matrix->SMPkluMatrix->KLUmatrixIntermediate[i]; } } - } else { - spSolve (Matrix->SPmatrix, RHS, RHS, NULL, NULL) ; + } + else { + spSolve(Matrix->SPmatrix, RHS, RHS, NULL, NULL); } } #ifdef CIDER void -SMPsolveKLUforCIDER (SMPmatrix *Matrix, double RHS[], double RHSsolution[], double iRHS[], double iRHSsolution[]) +SMPsolveKLUforCIDER(SMPmatrix* Matrix, double RHS[], double RHSsolution[], double iRHS[], double iRHSsolution[]) { - int ret ; - unsigned int i ; - double *KLUmatrixIntermediate ; + int ret; + unsigned int i; + double* KLUmatrixIntermediate; if (Matrix->CKTkluMODE) { if (Matrix->SMPkluMatrix->KLUmatrixIsComplex) { - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i] = RHS [i + 1] ; - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i + 1] = iRHS [i + 1] ; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i] = RHS[i + 1]; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i + 1] = iRHS[i + 1]; } - ret = klu_z_solve (Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + ret = klu_z_solve(Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex, Matrix->SMPkluMatrix->KLUmatrixCommon); - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { - RHSsolution [i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i] ; - iRHSsolution [i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex [2 * i + 1] ; + RHSsolution[i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i]; + iRHSsolution[i + 1] = Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex[2 * i + 1]; } - } else { + } + else { /* Allocate the Intermediate Vector */ - KLUmatrixIntermediate = (double *) malloc (Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)) ; + KLUmatrixIntermediate = (double*)malloc(Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)); - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) { - KLUmatrixIntermediate [i] = RHS [i + 1] ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { + KLUmatrixIntermediate[i] = RHS[i + 1]; } - ret = klu_solve (Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, - KLUmatrixIntermediate, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + ret = klu_solve(Matrix->SMPkluMatrix->KLUmatrixSymbolic, Matrix->SMPkluMatrix->KLUmatrixNumeric, (int)Matrix->SMPkluMatrix->KLUmatrixN, 1, + KLUmatrixIntermediate, Matrix->SMPkluMatrix->KLUmatrixCommon); - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN ; i++) { - RHSsolution [i + 1] = KLUmatrixIntermediate [i] ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN; i++) { + RHSsolution[i + 1] = KLUmatrixIntermediate[i]; } /* Free the Intermediate Vector */ - free (KLUmatrixIntermediate) ; + free(KLUmatrixIntermediate); } - } else { + } + else { - spSolve (Matrix->SPmatrix, RHS, RHSsolution, iRHS, iRHSsolution) ; + spSolve(Matrix->SPmatrix, RHS, RHSsolution, iRHS, iRHSsolution); } } #endif @@ -1065,12 +1092,13 @@ SMPsolveKLUforCIDER (SMPmatrix *Matrix, double RHS[], double RHSsolution[], doub * SMPmatSize() */ int -SMPmatSize (SMPmatrix *Matrix) +SMPmatSize(SMPmatrix* Matrix) { if (Matrix->CKTkluMODE) { - return (int)Matrix->SMPkluMatrix->KLUmatrixN ; - } else { - return spGetSize (Matrix->SPmatrix, 1) ; + return (int)Matrix->SMPkluMatrix->KLUmatrixN; + } + else { + return spGetSize(Matrix->SPmatrix, 1); } } @@ -1078,98 +1106,101 @@ SMPmatSize (SMPmatrix *Matrix) * SMPnewMatrix() */ int -SMPnewMatrix (SMPmatrix *Matrix, int size) +SMPnewMatrix(SMPmatrix* Matrix, int size) { - int Error ; + int Error; if (Matrix->CKTkluMODE) { /* Allocate the KLU Matrix Data Structure */ - Matrix->SMPkluMatrix = (KLUmatrix *) malloc (sizeof (KLUmatrix)) ; - Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ = 0 ; - Matrix->SMPkluMatrix->KLUmatrixLinkedListCOO = NULL ; + Matrix->SMPkluMatrix = (KLUmatrix*)malloc(sizeof(KLUmatrix)); + Matrix->SMPkluMatrix->KLUmatrixLinkedListNZ = 0; + Matrix->SMPkluMatrix->KLUmatrixLinkedListCOO = NULL; /* Initialize the KLU Matrix Internal Pointers */ - Matrix->SMPkluMatrix->KLUmatrixCommon = (klu_common *) malloc (sizeof (klu_common)) ; ; - Matrix->SMPkluMatrix->KLUmatrixSymbolic = NULL ; - Matrix->SMPkluMatrix->KLUmatrixNumeric = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAp = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAi = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAx = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL ; - Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUmatrixReal ; - Matrix->SMPkluMatrix->KLUmatrixIntermediate = NULL ; - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL ; - Matrix->SMPkluMatrix->KLUmatrixNZ = 0 ; - Matrix->SMPkluMatrix->KLUmatrixBindStructCOO = NULL ; -// Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew = NULL ; - Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld = NULL ; - Matrix->SMPkluMatrix->KLUmatrixDiag = NULL ; + Matrix->SMPkluMatrix->KLUmatrixCommon = (klu_common*)malloc(sizeof(klu_common)); ; + Matrix->SMPkluMatrix->KLUmatrixSymbolic = NULL; + Matrix->SMPkluMatrix->KLUmatrixNumeric = NULL; + Matrix->SMPkluMatrix->KLUmatrixAp = NULL; + Matrix->SMPkluMatrix->KLUmatrixAi = NULL; + Matrix->SMPkluMatrix->KLUmatrixAx = NULL; + Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL; + Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUmatrixReal; + Matrix->SMPkluMatrix->KLUmatrixIntermediate = NULL; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL; + Matrix->SMPkluMatrix->KLUmatrixNZ = 0; + Matrix->SMPkluMatrix->KLUmatrixBindStructCOO = NULL; + // Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew = NULL ; + Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld = NULL; + Matrix->SMPkluMatrix->KLUmatrixDiag = NULL; /* Initialize the KLU Common Data Structure */ - klu_defaults (Matrix->SMPkluMatrix->KLUmatrixCommon) ; - Matrix->SMPkluMatrix->KLUmatrixCommon->memgrow = Matrix->CKTkluMemGrowFactor ; + klu_defaults(Matrix->SMPkluMatrix->KLUmatrixCommon); + Matrix->SMPkluMatrix->KLUmatrixCommon->memgrow = Matrix->CKTkluMemGrowFactor; /* Allocate KLU data structures */ - Matrix->SMPkluMatrix->KLUmatrixN = (unsigned int)size ; - Matrix->SMPkluMatrix->KLUmatrixTrashCOO = (double *) malloc (2 * sizeof (double)) ; + Matrix->SMPkluMatrix->KLUmatrixN = (unsigned int)size; + Matrix->SMPkluMatrix->KLUmatrixTrashCOO = (double*)malloc(2 * sizeof(double)); - return spOKAY ; - } else { - Matrix->SPmatrix = spCreate (size, 1, &Error) ; - return Error ; + return spOKAY; + } + else { + Matrix->SPmatrix = spCreate(size, 1, &Error); + return Error; } } #ifdef CIDER int -SMPnewMatrixKLUforCIDER (SMPmatrix *Matrix, int size, unsigned int KLUmatrixIsComplex) +SMPnewMatrixKLUforCIDER(SMPmatrix* Matrix, int size, unsigned int KLUmatrixIsComplex) { - int Error ; - unsigned int i ; + int Error; + unsigned int i; if (Matrix->CKTkluMODE) { /* Allocate the KLU Matrix Data Structure */ - Matrix->SMPkluMatrix = (KLUmatrix *) malloc (sizeof (KLUmatrix)) ; + Matrix->SMPkluMatrix = (KLUmatrix*)malloc(sizeof(KLUmatrix)); /* Initialize the KLU Matrix Internal Pointers */ - Matrix->SMPkluMatrix->KLUmatrixCommon = (klu_common *) malloc (sizeof (klu_common)) ; ; - Matrix->SMPkluMatrix->KLUmatrixSymbolic = NULL ; - Matrix->SMPkluMatrix->KLUmatrixNumeric = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAp = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAi = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL ; + Matrix->SMPkluMatrix->KLUmatrixCommon = (klu_common*)malloc(sizeof(klu_common)); ; + Matrix->SMPkluMatrix->KLUmatrixSymbolic = NULL; + Matrix->SMPkluMatrix->KLUmatrixNumeric = NULL; + Matrix->SMPkluMatrix->KLUmatrixAp = NULL; + Matrix->SMPkluMatrix->KLUmatrixAi = NULL; + Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL; if (KLUmatrixIsComplex) { - Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUMatrixComplex ; - } else { - Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUmatrixReal ; + Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUMatrixComplex; } - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL ; - Matrix->SMPkluMatrix->KLUmatrixNZ = 0 ; - Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER = NULL ; - Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER = NULL ; + else { + Matrix->SMPkluMatrix->KLUmatrixIsComplex = KLUmatrixReal; + } + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL; + Matrix->SMPkluMatrix->KLUmatrixNZ = 0; + Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER = NULL; + Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER = NULL; - Matrix->SMPkluMatrix->KLUmatrixDiag = NULL ; + Matrix->SMPkluMatrix->KLUmatrixDiag = NULL; /* Initialize the KLU Common Data Structure */ - klu_defaults (Matrix->SMPkluMatrix->KLUmatrixCommon) ; + klu_defaults(Matrix->SMPkluMatrix->KLUmatrixCommon); /* Allocate KLU data structures */ - Matrix->SMPkluMatrix->KLUmatrixN = (unsigned int)size ; - Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER = (int *) malloc (Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN * sizeof(int)) ; - Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER = (int *) malloc (Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN * sizeof(int)) ; - Matrix->SMPkluMatrix->KLUmatrixTrashCOO = (double *) malloc (2 * sizeof(double)) ; - Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER = (double *) malloc (2 * Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)) ; + Matrix->SMPkluMatrix->KLUmatrixN = (unsigned int)size; + Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER = (int*)malloc(Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN * sizeof(int)); + Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER = (int*)malloc(Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN * sizeof(int)); + Matrix->SMPkluMatrix->KLUmatrixTrashCOO = (double*)malloc(2 * sizeof(double)); + Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER = (double*)malloc(2 * Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)); /* Pre-set the values of Row and Col */ - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN ; i++) { - Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER [i] = -1 ; - Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER [i] = -1 ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixN * Matrix->SMPkluMatrix->KLUmatrixN; i++) { + Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER[i] = -1; + Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER[i] = -1; } - return spOKAY ; - } else { - Matrix->SPmatrix = spCreate (size, (int)KLUmatrixIsComplex, &Error) ; - return Error ; + return spOKAY; + } + else { + Matrix->SPmatrix = spCreate(size, (int)KLUmatrixIsComplex, &Error); + return Error; } } #endif @@ -1179,71 +1210,73 @@ SMPnewMatrixKLUforCIDER (SMPmatrix *Matrix, int size, unsigned int KLUmatrixIsCo */ void -SMPdestroy (SMPmatrix *Matrix) +SMPdestroy(SMPmatrix* Matrix) { if (Matrix->CKTkluMODE) { - klu_free_numeric (&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon) ; - klu_free_symbolic (&(Matrix->SMPkluMatrix->KLUmatrixSymbolic), Matrix->SMPkluMatrix->KLUmatrixCommon) ; - free (Matrix->SMPkluMatrix->KLUmatrixAp) ; - free (Matrix->SMPkluMatrix->KLUmatrixAi) ; - free (Matrix->SMPkluMatrix->KLUmatrixAx) ; - free (Matrix->SMPkluMatrix->KLUmatrixAxComplex) ; - free (Matrix->SMPkluMatrix->KLUmatrixIntermediate) ; - free (Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex) ; - free (Matrix->SMPkluMatrix->KLUmatrixBindStructCOO) ; -// free (Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew) ; - free (Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld) ; - free (Matrix->SMPkluMatrix->KLUmatrixTrashCOO) ; - Matrix->SMPkluMatrix->KLUmatrixAp = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAi = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAx = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL ; - Matrix->SMPkluMatrix->KLUmatrixIntermediate = NULL ; - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL ; - Matrix->SMPkluMatrix->KLUmatrixBindStructCOO = NULL ; -// Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew = NULL ; - Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld = NULL ; - Matrix->SMPkluMatrix->KLUmatrixTrashCOO = NULL ; - free (Matrix->SMPkluMatrix->KLUmatrixDiag) ; - free (Matrix->SMPkluMatrix->KLUmatrixCommon) ; - free (Matrix->SMPkluMatrix) ; - } else { - spDestroy (Matrix->SPmatrix) ; + klu_free_numeric(&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon); + klu_free_symbolic(&(Matrix->SMPkluMatrix->KLUmatrixSymbolic), Matrix->SMPkluMatrix->KLUmatrixCommon); + free(Matrix->SMPkluMatrix->KLUmatrixAp); + free(Matrix->SMPkluMatrix->KLUmatrixAi); + free(Matrix->SMPkluMatrix->KLUmatrixAx); + free(Matrix->SMPkluMatrix->KLUmatrixAxComplex); + free(Matrix->SMPkluMatrix->KLUmatrixIntermediate); + free(Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex); + free(Matrix->SMPkluMatrix->KLUmatrixBindStructCOO); + // free (Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew) ; + free(Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld); + free(Matrix->SMPkluMatrix->KLUmatrixTrashCOO); + Matrix->SMPkluMatrix->KLUmatrixAp = NULL; + Matrix->SMPkluMatrix->KLUmatrixAi = NULL; + Matrix->SMPkluMatrix->KLUmatrixAx = NULL; + Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL; + Matrix->SMPkluMatrix->KLUmatrixIntermediate = NULL; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL; + Matrix->SMPkluMatrix->KLUmatrixBindStructCOO = NULL; + // Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingOldToNew = NULL ; + Matrix->SMPkluMatrix->KLUmatrixNodeCollapsingNewToOld = NULL; + Matrix->SMPkluMatrix->KLUmatrixTrashCOO = NULL; + free(Matrix->SMPkluMatrix->KLUmatrixDiag); + free(Matrix->SMPkluMatrix->KLUmatrixCommon); + free(Matrix->SMPkluMatrix); + } + else { + spDestroy(Matrix->SPmatrix); } } #ifdef CIDER void -SMPdestroyKLUforCIDER (SMPmatrix *Matrix) +SMPdestroyKLUforCIDER(SMPmatrix* Matrix) { if (Matrix->CKTkluMODE) { - klu_free_numeric (&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon) ; - klu_free_symbolic (&(Matrix->SMPkluMatrix->KLUmatrixSymbolic), Matrix->SMPkluMatrix->KLUmatrixCommon) ; - free (Matrix->SMPkluMatrix->KLUmatrixAp) ; - free (Matrix->SMPkluMatrix->KLUmatrixAi) ; - free (Matrix->SMPkluMatrix->KLUmatrixAxComplex) ; - free (Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex) ; - free (Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER) ; - free (Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER) ; - free (Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER) ; - free (Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER) ; - free (Matrix->SMPkluMatrix->KLUmatrixTrashCOO) ; - Matrix->SMPkluMatrix->KLUmatrixAp = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAi = NULL ; - Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL ; - Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL ; - Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER = NULL ; - Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER = NULL ; - Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER = NULL ; - Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER = NULL ; - Matrix->SMPkluMatrix->KLUmatrixTrashCOO = NULL ; - free (Matrix->SMPkluMatrix->KLUmatrixDiag) ; - free (Matrix->SMPkluMatrix->KLUmatrixCommon) ; - free (Matrix->SMPkluMatrix) ; - } else { - spDestroy (Matrix->SPmatrix) ; + klu_free_numeric(&(Matrix->SMPkluMatrix->KLUmatrixNumeric), Matrix->SMPkluMatrix->KLUmatrixCommon); + klu_free_symbolic(&(Matrix->SMPkluMatrix->KLUmatrixSymbolic), Matrix->SMPkluMatrix->KLUmatrixCommon); + free(Matrix->SMPkluMatrix->KLUmatrixAp); + free(Matrix->SMPkluMatrix->KLUmatrixAi); + free(Matrix->SMPkluMatrix->KLUmatrixAxComplex); + free(Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex); + free(Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER); + free(Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER); + free(Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER); + free(Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER); + free(Matrix->SMPkluMatrix->KLUmatrixTrashCOO); + Matrix->SMPkluMatrix->KLUmatrixAp = NULL; + Matrix->SMPkluMatrix->KLUmatrixAi = NULL; + Matrix->SMPkluMatrix->KLUmatrixAxComplex = NULL; + Matrix->SMPkluMatrix->KLUmatrixIntermediateComplex = NULL; + Matrix->SMPkluMatrix->KLUmatrixBindStructForCIDER = NULL; + Matrix->SMPkluMatrix->KLUmatrixColCOOforCIDER = NULL; + Matrix->SMPkluMatrix->KLUmatrixRowCOOforCIDER = NULL; + Matrix->SMPkluMatrix->KLUmatrixValueComplexCOOforCIDER = NULL; + Matrix->SMPkluMatrix->KLUmatrixTrashCOO = NULL; + free(Matrix->SMPkluMatrix->KLUmatrixDiag); + free(Matrix->SMPkluMatrix->KLUmatrixCommon); + free(Matrix->SMPkluMatrix); + } + else { + spDestroy(Matrix->SPmatrix); } } #endif @@ -1253,34 +1286,37 @@ SMPdestroyKLUforCIDER (SMPmatrix *Matrix) */ int -SMPpreOrder (SMPmatrix *Matrix) +SMPpreOrder(SMPmatrix* Matrix) { if (Matrix->CKTkluMODE) { if (CircuitIsDigital() && Matrix->SMPkluMatrix->KLUmatrixN == 0) { - // XSPICE pure digital circuits produce empty KLU matrix - return 0 ; + // XSPICE pure digital circuits produce empty KLU matrix + return 0; } - Matrix->SMPkluMatrix->KLUmatrixSymbolic = klu_analyze ((int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixAp, - Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + Matrix->SMPkluMatrix->KLUmatrixSymbolic = klu_analyze((int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixAp, + Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixCommon); if (Matrix->SMPkluMatrix->KLUmatrixSymbolic == NULL) { if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_EMPTY_MATRIX) { - fprintf (stderr, "Error (PreOrder): KLU Matrix is empty\n") ; - return 0 ; - } else { - fprintf (stderr, "Error (PreOrder): KLUsymbolic object is NULL. A problem occurred\n") ; - return 1 ; + fprintf(stderr, "Error (PreOrder): KLU Matrix is empty\n"); + return 0; + } + else { + fprintf(stderr, "Error (PreOrder): KLUsymbolic object is NULL. A problem occurred\n"); + return 1; } - } else { - return 0 ; } - } else { - spMNA_Preorder (Matrix->SPmatrix) ; - return spError (Matrix->SPmatrix) ; + else { + return 0; + } + } + else { + spMNA_Preorder(Matrix->SPmatrix); + return spError(Matrix->SPmatrix); } } @@ -1289,10 +1325,10 @@ SMPpreOrder (SMPmatrix *Matrix) */ void -SMPprintRHS (SMPmatrix *Matrix, char *Filename, RealVector RHS, RealVector iRHS) +SMPprintRHS(SMPmatrix* Matrix, char* Filename, RealVector RHS, RealVector iRHS) { if (!Matrix->CKTkluMODE) - spFileVector (Matrix->SPmatrix, Filename, RHS, iRHS) ; + spFileVector(Matrix->SPmatrix, Filename, RHS, iRHS); } /* @@ -1300,59 +1336,63 @@ SMPprintRHS (SMPmatrix *Matrix, char *Filename, RealVector RHS, RealVector iRHS) */ void -SMPprint (SMPmatrix *Matrix, char *Filename) +SMPprint(SMPmatrix* Matrix, char* Filename) { if (Matrix->CKTkluMODE) { if (Matrix->SMPkluMatrix->KLUmatrixIsComplex) { - klu_z_print (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, - (int)Matrix->SMPkluMatrix->KLUmatrixN, NULL, NULL) ; - } else { - klu_print (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAx, - (int)Matrix->SMPkluMatrix->KLUmatrixN, NULL, NULL) ; + klu_z_print(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, + (int)Matrix->SMPkluMatrix->KLUmatrixN, NULL, NULL); } - } else { + else { + klu_print(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAx, + (int)Matrix->SMPkluMatrix->KLUmatrixN, NULL, NULL); + } + } + else { if (Filename) - spFileMatrix (Matrix->SPmatrix, Filename, "Circuit Matrix", 0, 1, 1) ; + spFileMatrix(Matrix->SPmatrix, Filename, "Circuit Matrix", 0, 1, 1); else - spPrint (Matrix->SPmatrix, 0, 1, 1) ; + spPrint(Matrix->SPmatrix, 0, 1, 1); } } #ifdef CIDER void -SMPprintKLUforCIDER (SMPmatrix *Matrix, char *Filename) +SMPprintKLUforCIDER(SMPmatrix* Matrix, char* Filename) { - unsigned int i ; - double *KLUmatrixAx ; + unsigned int i; + double* KLUmatrixAx; if (Matrix->CKTkluMODE) { if (Matrix->SMPkluMatrix->KLUmatrixIsComplex) { - klu_z_print (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, - (int)Matrix->SMPkluMatrix->KLUmatrixN, NULL, NULL) ; - } else { + klu_z_print(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, + (int)Matrix->SMPkluMatrix->KLUmatrixN, NULL, NULL); + } + else { /* Allocate the Real Matrix */ - KLUmatrixAx = (double *) malloc (Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(double)) ; + KLUmatrixAx = (double*)malloc(Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(double)); /* Copy the Complex Matrix into the Real Matrix */ - for (i = 0 ; i < Matrix->SMPkluMatrix->KLUmatrixNZ ; i++) { - KLUmatrixAx [i] = Matrix->SMPkluMatrix->KLUmatrixAxComplex [2 * i] ; + for (i = 0; i < Matrix->SMPkluMatrix->KLUmatrixNZ; i++) { + KLUmatrixAx[i] = Matrix->SMPkluMatrix->KLUmatrixAxComplex[2 * i]; } /* Print the Real Matrix */ - klu_print (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, KLUmatrixAx, (int)Matrix->SMPkluMatrix->KLUmatrixN, NULL, NULL) ; + klu_print(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, KLUmatrixAx, (int)Matrix->SMPkluMatrix->KLUmatrixN, NULL, NULL); /* Free the Real Matrix Storage */ - free (KLUmatrixAx) ; + free(KLUmatrixAx); } - } else { + } + else { if (Filename) - spFileMatrix (Matrix->SPmatrix, Filename, "Circuit Matrix", 0, 1, 1) ; + spFileMatrix(Matrix->SPmatrix, Filename, "Circuit Matrix", 0, 1, 1); else - spPrint (Matrix->SPmatrix, 0, 1, 1) ; + spPrint(Matrix->SPmatrix, 0, 1, 1); } } #endif @@ -1361,277 +1401,279 @@ SMPprintKLUforCIDER (SMPmatrix *Matrix, char *Filename) * SMPgetError() */ void -SMPgetError (SMPmatrix *Matrix, int *Col, int *Row) +SMPgetError(SMPmatrix* Matrix, int* Col, int* Row) { if (Matrix->CKTkluMODE) { if (Matrix->SMPkluMatrix->KLUmatrixNZ == 0) { - *Row = 0 ; - *Col = 0 ; - } else { - *Row = Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1 ; - *Col = Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1 ; + *Row = 0; + *Col = 0; } - } else { - spWhereSingular (Matrix->SPmatrix, Row, Col) ; + else { + *Row = Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1; + *Col = Matrix->SMPkluMatrix->KLUmatrixCommon->singular_col + 1; + } + } + else { + spWhereSingular(Matrix->SPmatrix, Row, Col); } } void -spDeterminant_KLU (SMPmatrix *Matrix, int *pExponent, RealNumber *pDeterminant, RealNumber *piDeterminant) +spDeterminant_KLU(SMPmatrix* Matrix, int* pExponent, RealNumber* pDeterminant, RealNumber* piDeterminant) { - int I, Size ; - RealNumber Norm, nr, ni ; - ComplexNumber Pivot, cDeterminant, Udiag ; + int I, Size; + RealNumber Norm, nr, ni; + ComplexNumber Pivot, cDeterminant, Udiag; - int *P, *Q ; - double *Rs, *Ux, *Uz ; - unsigned int nSwap, nSwapP, nSwapQ ; + int* P, * Q; + double* Rs, * Ux, * Uz; + unsigned int nSwap, nSwapP, nSwapQ; #define NORM(a) (nr = ABS((a).Real), ni = ABS((a).Imag), MAX (nr,ni)) - *pExponent = 0 ; + * pExponent = 0; if (Matrix->SMPkluMatrix->KLUmatrixCommon->status == KLU_SINGULAR) { - *pDeterminant = 0.0 ; + *pDeterminant = 0.0; if (Matrix->SMPkluMatrix->KLUmatrixIsComplex == KLUMatrixComplex) { - *piDeterminant = 0.0 ; + *piDeterminant = 0.0; } - return ; + return; } - Size = (int)Matrix->SMPkluMatrix->KLUmatrixN ; - I = 0 ; + Size = (int)Matrix->SMPkluMatrix->KLUmatrixN; + I = 0; - P = (int *) malloc ((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof (int)) ; - Q = (int *) malloc ((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof (int)) ; + P = (int*)malloc((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof(int)); + Q = (int*)malloc((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof(int)); - Ux = (double *) malloc ((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof (double)) ; + Ux = (double*)malloc((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)); - Rs = (double *) malloc ((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof (double)) ; + Rs = (double*)malloc((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)); if (Matrix->SMPkluMatrix->KLUmatrixIsComplex == KLUMatrixComplex) /* Complex Case. */ { - cDeterminant.Real = 1.0 ; - cDeterminant.Imag = 0.0 ; + cDeterminant.Real = 1.0; + cDeterminant.Imag = 0.0; - Uz = (double *) malloc ((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof (double)) ; -/* - int *Lp, *Li, *Up, *Ui, *Fp, *Fi, *P, *Q ; - double *Lx, *Lz, *Ux, *Uz, *Fx, *Fz, *Rs ; - Lp = (int *) malloc (((size_t)Matrix->CKTkluN + 1) * sizeof (int)) ; - Li = (int *) malloc ((size_t)Matrix->CKTkluNumeric->lnz * sizeof (int)) ; - Lx = (double *) malloc ((size_t)Matrix->CKTkluNumeric->lnz * sizeof (double)) ; - Lz = (double *) malloc ((size_t)Matrix->CKTkluNumeric->lnz * sizeof (double)) ; - Up = (int *) malloc (((size_t)Matrix->CKTkluN + 1) * sizeof (int)) ; - Ui = (int *) malloc ((size_t)Matrix->CKTkluNumeric->unz * sizeof (int)) ; - Ux = (double *) malloc ((size_t)Matrix->CKTkluNumeric->unz * sizeof (double)) ; - Uz = (double *) malloc ((size_t)Matrix->CKTkluNumeric->unz * sizeof (double)) ; - Fp = (int *) malloc (((size_t)Matrix->CKTkluN + 1) * sizeof (int)) ; - Fi = (int *) malloc ((size_t)Matrix->CKTkluNumeric->Offp [Matrix->CKTkluN] * sizeof (int)) ; - Fx = (double *) malloc ((size_t)Matrix->CKTkluNumeric->Offp [Matrix->CKTkluN] * sizeof (double)) ; - Fz = (double *) malloc ((size_t)Matrix->CKTkluNumeric->Offp [Matrix->CKTkluN] * sizeof (double)) ; - klu_z_extract (Matrix->CKTkluNumeric, Matrix->CKTkluSymbolic, - Lp, Li, Lx, Lz, - Up, Ui, Ux, Uz, - Fp, Fi, Fx, Fz, - P, Q, Rs, NULL, - Matrix->CKTkluCommon) ; -*/ - klu_z_extract_Udiag (Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixSymbolic, Ux, Uz, P, Q, Rs, Matrix->SMPkluMatrix->KLUmatrixCommon) ; -/* - for (I = 0 ; I < Matrix->CKTkluNumeric->lnz ; I++) - { - printf ("L - Value: %-.9g\t%-.9g\n", Lx [I], Lz [I]) ; - } - for (I = 0 ; I < Matrix->CKTkluNumeric->unz ; I++) - { - printf ("U - Value: %-.9g\t%-.9g\n", Ux [I], Uz [I]) ; - } - for (I = 0 ; I < Matrix->CKTkluNumeric->Offp [Matrix->CKTkluN] ; I++) - { - printf ("F - Value: %-.9g\t%-.9g\n", Fx [I], Fz [I]) ; - } + Uz = (double*)malloc((size_t)Matrix->SMPkluMatrix->KLUmatrixN * sizeof(double)); + /* + int *Lp, *Li, *Up, *Ui, *Fp, *Fi, *P, *Q ; + double *Lx, *Lz, *Ux, *Uz, *Fx, *Fz, *Rs ; + Lp = (int *) malloc (((size_t)Matrix->CKTkluN + 1) * sizeof (int)) ; + Li = (int *) malloc ((size_t)Matrix->CKTkluNumeric->lnz * sizeof (int)) ; + Lx = (double *) malloc ((size_t)Matrix->CKTkluNumeric->lnz * sizeof (double)) ; + Lz = (double *) malloc ((size_t)Matrix->CKTkluNumeric->lnz * sizeof (double)) ; + Up = (int *) malloc (((size_t)Matrix->CKTkluN + 1) * sizeof (int)) ; + Ui = (int *) malloc ((size_t)Matrix->CKTkluNumeric->unz * sizeof (int)) ; + Ux = (double *) malloc ((size_t)Matrix->CKTkluNumeric->unz * sizeof (double)) ; + Uz = (double *) malloc ((size_t)Matrix->CKTkluNumeric->unz * sizeof (double)) ; + Fp = (int *) malloc (((size_t)Matrix->CKTkluN + 1) * sizeof (int)) ; + Fi = (int *) malloc ((size_t)Matrix->CKTkluNumeric->Offp [Matrix->CKTkluN] * sizeof (int)) ; + Fx = (double *) malloc ((size_t)Matrix->CKTkluNumeric->Offp [Matrix->CKTkluN] * sizeof (double)) ; + Fz = (double *) malloc ((size_t)Matrix->CKTkluNumeric->Offp [Matrix->CKTkluN] * sizeof (double)) ; + klu_z_extract (Matrix->CKTkluNumeric, Matrix->CKTkluSymbolic, + Lp, Li, Lx, Lz, + Up, Ui, Ux, Uz, + Fp, Fi, Fx, Fz, + P, Q, Rs, NULL, + Matrix->CKTkluCommon) ; + */ + klu_z_extract_Udiag(Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixSymbolic, Ux, Uz, P, Q, Rs, Matrix->SMPkluMatrix->KLUmatrixCommon); + /* + for (I = 0 ; I < Matrix->CKTkluNumeric->lnz ; I++) + { + printf ("L - Value: %-.9g\t%-.9g\n", Lx [I], Lz [I]) ; + } + for (I = 0 ; I < Matrix->CKTkluNumeric->unz ; I++) + { + printf ("U - Value: %-.9g\t%-.9g\n", Ux [I], Uz [I]) ; + } + for (I = 0 ; I < Matrix->CKTkluNumeric->Offp [Matrix->CKTkluN] ; I++) + { + printf ("F - Value: %-.9g\t%-.9g\n", Fx [I], Fz [I]) ; + } - for (I = 0 ; I < Matrix->CKTkluN ; I++) + for (I = 0 ; I < Matrix->CKTkluN ; I++) + { + printf ("U - Value: %-.9g\t%-.9g\n", Ux [I], Uz [I]) ; + } + */ + nSwapP = 0; + for (I = 0; I < (int)Matrix->SMPkluMatrix->KLUmatrixN; I++) { - printf ("U - Value: %-.9g\t%-.9g\n", Ux [I], Uz [I]) ; - } -*/ - nSwapP = 0 ; - for (I = 0 ; I < (int)Matrix->SMPkluMatrix->KLUmatrixN ; I++) - { - if (P [I] != I) + if (P[I] != I) { - nSwapP++ ; + nSwapP++; } } - nSwapP /= 2 ; + nSwapP /= 2; - nSwapQ = 0 ; - for (I = 0 ; I < (int)Matrix->SMPkluMatrix->KLUmatrixN ; I++) + nSwapQ = 0; + for (I = 0; I < (int)Matrix->SMPkluMatrix->KLUmatrixN; I++) { - if (Q [I] != I) + if (Q[I] != I) { - nSwapQ++ ; + nSwapQ++; } } - nSwapQ /= 2 ; + nSwapQ /= 2; - nSwap = nSwapP + nSwapQ ; -/* - free (Lp) ; - free (Li) ; - free (Lx) ; - free (Lz) ; - free (Up) ; - free (Ui) ; - free (Fp) ; - free (Fi) ; - free (Fx) ; - free (Fz) ; -*/ - I = 0 ; + nSwap = nSwapP + nSwapQ; + /* + free (Lp) ; + free (Li) ; + free (Lx) ; + free (Lz) ; + free (Up) ; + free (Ui) ; + free (Fp) ; + free (Fi) ; + free (Fx) ; + free (Fz) ; + */ + I = 0; while (I < Size) { - Udiag.Real = 1 / (Ux [I] * Rs [I]) ; - Udiag.Imag = Uz [I] * Rs [I] ; + Udiag.Real = 1 / (Ux[I] * Rs[I]); + Udiag.Imag = Uz[I] * Rs[I]; -// printf ("Udiag.Real: %-.9g\tUdiag.Imag %-.9g\n", Udiag.Real, Udiag.Imag) ; + // printf ("Udiag.Real: %-.9g\tUdiag.Imag %-.9g\n", Udiag.Real, Udiag.Imag) ; - CMPLX_RECIPROCAL (Pivot, Udiag) ; - CMPLX_MULT_ASSIGN (cDeterminant, Pivot) ; + CMPLX_RECIPROCAL(Pivot, Udiag); + CMPLX_MULT_ASSIGN(cDeterminant, Pivot); -// printf ("cDeterminant.Real: %-.9g\tcDeterminant.Imag %-.9g\n", cDeterminant.Real, cDeterminant.Imag) ; + // printf ("cDeterminant.Real: %-.9g\tcDeterminant.Imag %-.9g\n", cDeterminant.Real, cDeterminant.Imag) ; - /* Scale Determinant. */ - Norm = NORM (cDeterminant) ; + /* Scale Determinant. */ + Norm = NORM(cDeterminant); if (Norm != 0.0) { - while (Norm >= 1.0e12) + while (Norm >= 1.0e12) { - cDeterminant.Real *= 1.0e-12 ; - cDeterminant.Imag *= 1.0e-12 ; - *pExponent += 12 ; - Norm = NORM (cDeterminant) ; + cDeterminant.Real *= 1.0e-12; + cDeterminant.Imag *= 1.0e-12; + *pExponent += 12; + Norm = NORM(cDeterminant); } while (Norm < 1.0e-12) { - cDeterminant.Real *= 1.0e12 ; - cDeterminant.Imag *= 1.0e12 ; - *pExponent -= 12 ; - Norm = NORM (cDeterminant) ; + cDeterminant.Real *= 1.0e12; + cDeterminant.Imag *= 1.0e12; + *pExponent -= 12; + Norm = NORM(cDeterminant); } } - I++ ; + I++; } - /* Scale Determinant again, this time to be between 1.0 <= x < 10.0. */ - Norm = NORM (cDeterminant) ; + /* Scale Determinant again, this time to be between 1.0 <= x < 10.0. */ + Norm = NORM(cDeterminant); if (Norm != 0.0) { - while (Norm >= 10.0) + while (Norm >= 10.0) { - cDeterminant.Real *= 0.1 ; - cDeterminant.Imag *= 0.1 ; - (*pExponent)++ ; - Norm = NORM (cDeterminant) ; + cDeterminant.Real *= 0.1; + cDeterminant.Imag *= 0.1; + (*pExponent)++; + Norm = NORM(cDeterminant); } while (Norm < 1.0) { - cDeterminant.Real *= 10.0 ; - cDeterminant.Imag *= 10.0 ; - (*pExponent)-- ; - Norm = NORM (cDeterminant) ; + cDeterminant.Real *= 10.0; + cDeterminant.Imag *= 10.0; + (*pExponent)--; + Norm = NORM(cDeterminant); } } if (nSwap % 2 != 0) { - CMPLX_NEGATE (cDeterminant) ; + CMPLX_NEGATE(cDeterminant); } - *pDeterminant = cDeterminant.Real ; - *piDeterminant = cDeterminant.Imag ; + *pDeterminant = cDeterminant.Real; + *piDeterminant = cDeterminant.Imag; - free (Uz) ; + free(Uz); } else { - /* Real Case. */ - *pDeterminant = 1.0 ; + /* Real Case. */ + *pDeterminant = 1.0; - klu_extract_Udiag (Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixSymbolic, Ux, P, Q, Rs, Matrix->SMPkluMatrix->KLUmatrixCommon) ; + klu_extract_Udiag(Matrix->SMPkluMatrix->KLUmatrixNumeric, Matrix->SMPkluMatrix->KLUmatrixSymbolic, Ux, P, Q, Rs, Matrix->SMPkluMatrix->KLUmatrixCommon); - nSwapP = 0 ; - for (I = 0 ; I < (int)Matrix->SMPkluMatrix->KLUmatrixN ; I++) + nSwapP = 0; + for (I = 0; I < (int)Matrix->SMPkluMatrix->KLUmatrixN; I++) { - if (P [I] != I) + if (P[I] != I) { - nSwapP++ ; + nSwapP++; } } - nSwapP /= 2 ; + nSwapP /= 2; - nSwapQ = 0 ; - for (I = 0 ; I < (int)Matrix->SMPkluMatrix->KLUmatrixN ; I++) + nSwapQ = 0; + for (I = 0; I < (int)Matrix->SMPkluMatrix->KLUmatrixN; I++) { - if (Q [I] != I) + if (Q[I] != I) { - nSwapQ++ ; + nSwapQ++; } } - nSwapQ /= 2 ; + nSwapQ /= 2; - nSwap = nSwapP + nSwapQ ; + nSwap = nSwapP + nSwapQ; while (I < Size) { - *pDeterminant /= (Ux [I] * Rs [I]) ; + *pDeterminant /= (Ux[I] * Rs[I]); - /* Scale Determinant. */ + /* Scale Determinant. */ if (*pDeterminant != 0.0) { - while (ABS(*pDeterminant) >= 1.0e12) + while (ABS(*pDeterminant) >= 1.0e12) { - *pDeterminant *= 1.0e-12 ; - *pExponent += 12 ; + *pDeterminant *= 1.0e-12; + *pExponent += 12; } while (ABS(*pDeterminant) < 1.0e-12) { - *pDeterminant *= 1.0e12 ; - *pExponent -= 12 ; + *pDeterminant *= 1.0e12; + *pExponent -= 12; } } - I++ ; + I++; } - /* Scale Determinant again, this time to be between 1.0 <= x < - 10.0. */ + /* Scale Determinant again, this time to be between 1.0 <= x < + 10.0. */ if (*pDeterminant != 0.0) { - while (ABS(*pDeterminant) >= 10.0) + while (ABS(*pDeterminant) >= 10.0) { - *pDeterminant *= 0.1 ; - (*pExponent)++ ; + *pDeterminant *= 0.1; + (*pExponent)++; } while (ABS(*pDeterminant) < 1.0) { - *pDeterminant *= 10.0 ; - (*pExponent)-- ; + *pDeterminant *= 10.0; + (*pExponent)--; } } if (nSwap % 2 != 0) { - *pDeterminant = -*pDeterminant ; + *pDeterminant = -*pDeterminant; } } - free (P) ; - free (Q) ; - free (Ux) ; - free (Rs) ; + free(P); + free(Q); + free(Ux); + free(Rs); } /* @@ -1639,31 +1681,33 @@ spDeterminant_KLU (SMPmatrix *Matrix, int *pExponent, RealNumber *pDeterminant, * note: obsolete for Spice3d2 and later */ int -SMPcProdDiag (SMPmatrix *Matrix, SPcomplex *pMantissa, int *pExponent) +SMPcProdDiag(SMPmatrix* Matrix, SPcomplex* pMantissa, int* pExponent) { if (Matrix->CKTkluMODE) { - spDeterminant_KLU (Matrix, pExponent, &(pMantissa->real), &(pMantissa->imag)) ; - } else { - spDeterminant (Matrix->SPmatrix, pExponent, &(pMantissa->real), &(pMantissa->imag)) ; + spDeterminant_KLU(Matrix, pExponent, &(pMantissa->real), &(pMantissa->imag)); } - return spError (Matrix->SPmatrix) ; + else { + spDeterminant(Matrix->SPmatrix, pExponent, &(pMantissa->real), &(pMantissa->imag)); + } + return spError(Matrix->SPmatrix); } /* * SMPcDProd() */ int -SMPcDProd (SMPmatrix *Matrix, SPcomplex *pMantissa, int *pExponent) +SMPcDProd(SMPmatrix* Matrix, SPcomplex* pMantissa, int* pExponent) { double re, im, x, y, z; int p; if (Matrix->CKTkluMODE) { - spDeterminant_KLU (Matrix, &p, &re, &im) ; - } else { - spDeterminant (Matrix->SPmatrix, &p, &re, &im) ; + spDeterminant_KLU(Matrix, &p, &re, &im); + } + else { + spDeterminant(Matrix->SPmatrix, &p, &re, &im); } #ifndef M_LN2 @@ -1674,71 +1718,74 @@ SMPcDProd (SMPmatrix *Matrix, SPcomplex *pMantissa, int *pExponent) #endif #ifdef debug_print - printf ("Determinant 10: (%20g,%20g)^%d\n", re, im, p) ; + printf("Determinant 10: (%20g,%20g)^%d\n", re, im, p); #endif /* Convert base 10 numbers to base 2 numbers, for comparison */ y = p * M_LN10 / M_LN2; - x = (int) y; + x = (int)y; y -= x; /* ASSERT * x = integral part of exponent, y = fraction part of exponent */ - /* Fold in the fractional part */ + /* Fold in the fractional part */ #ifdef debug_print - printf (" ** base10 -> base2 int = %g, frac = %20g\n", x, y) ; + printf(" ** base10 -> base2 int = %g, frac = %20g\n", x, y); #endif - z = pow (2.0, y) ; - re *= z ; - im *= z ; + z = pow(2.0, y); + re *= z; + im *= z; #ifdef debug_print - printf (" ** multiplier = %20g\n", z) ; + printf(" ** multiplier = %20g\n", z); #endif /* Re-normalize (re or im may be > 2.0 or both < 1.0 */ if (re != 0.0) { - y = logb (re) ; - if (im != 0.0) - z = logb (im) ; - else - z = 0 ; - } else if (im != 0.0) { - z = logb (im) ; - y = 0 ; - } else { - /* Singular */ - /*printf("10 -> singular\n");*/ - y = 0 ; - z = 0 ; + y = logb(re); + if (im != 0.0) + z = logb(im); + else + z = 0; + } + else if (im != 0.0) { + z = logb(im); + y = 0; + } + else { + /* Singular */ + /*printf("10 -> singular\n");*/ + y = 0; + z = 0; } #ifdef debug_print - printf (" ** renormalize changes = %g,%g\n", y, z) ; + printf(" ** renormalize changes = %g,%g\n", y, z); #endif if (y < z) - y = z ; + y = z; - *pExponent = (int)(x + y) ; - x = scalbn (re, (int) -y) ; - z = scalbn (im, (int) -y) ; + *pExponent = (int)(x + y); + x = scalbn(re, (int)-y); + z = scalbn(im, (int)-y); #ifdef debug_print - printf (" ** values are: re %g, im %g, y %g, re' %g, im' %g\n", re, im, y, x, z) ; + printf(" ** values are: re %g, im %g, y %g, re' %g, im' %g\n", re, im, y, x, z); #endif - pMantissa->real = scalbn (re, (int) -y) ; - pMantissa->imag = scalbn (im, (int) -y) ; + pMantissa->real = scalbn(re, (int)-y); + pMantissa->imag = scalbn(im, (int)-y); #ifdef debug_print - printf ("Determinant 10->2: (%20g,%20g)^%d\n", pMantissa->real, pMantissa->imag, *pExponent) ; + printf("Determinant 10->2: (%20g,%20g)^%d\n", pMantissa->real, pMantissa->imag, *pExponent); #endif if (Matrix->CKTkluMODE) { - return 0 ; - } else { - return spError (Matrix->SPmatrix) ; + return 0; + } + else { + return spError(Matrix->SPmatrix); } } @@ -1749,53 +1796,53 @@ SMPcDProd (SMPmatrix *Matrix, SPcomplex *pMantissa, int *pExponent) * structures. */ -/* - * LOAD GMIN - * - * This routine adds Gmin to each diagonal element. Because Gmin is - * added to the current diagonal, which may bear little relation to - * what the outside world thinks is a diagonal, and because the - * elements that are diagonals may change after calling spOrderAndFactor, - * use of this routine is not recommended. It is included here simply - * for compatibility with Spice3. - */ + /* + * LOAD GMIN + * + * This routine adds Gmin to each diagonal element. Because Gmin is + * added to the current diagonal, which may bear little relation to + * what the outside world thinks is a diagonal, and because the + * elements that are diagonals may change after calling spOrderAndFactor, + * use of this routine is not recommended. It is included here simply + * for compatibility with Spice3. + */ static void -LoadGmin_CSC (double **diag, unsigned int n, double Gmin) +LoadGmin_CSC(double** diag, unsigned int n, double Gmin) { - unsigned int i ; + unsigned int i; if (Gmin != 0.0) { - for (i = 0 ; i < n ; i++) { - if (diag [i] != NULL) { + for (i = 0; i < n; i++) { + if (diag[i] != NULL) { // Not all the elements on the diagonal are present, when the circuit is parsed - *(diag [i]) += Gmin ; + *(diag[i]) += Gmin; } } } } static void -LoadGmin (SMPmatrix *eMatrix, double Gmin) +LoadGmin(SMPmatrix* eMatrix, double Gmin) { - MatrixPtr Matrix = eMatrix->SPmatrix ; - int I ; - ArrayOfElementPtrs Diag ; - ElementPtr diag ; + MatrixPtr Matrix = eMatrix->SPmatrix; + int I; + ArrayOfElementPtrs Diag; + ElementPtr diag; /* Begin `LoadGmin'. */ - assert (IS_SPARSE (Matrix)) ; + assert(IS_SPARSE(Matrix)); if (Gmin != 0.0) { - Diag = Matrix->Diag ; - for (I = Matrix->Size ; I > 0 ; I--) + Diag = Matrix->Diag; + for (I = Matrix->Size; I > 0; I--) { - if ((diag = Diag [I]) != NULL) - diag->Real += Gmin ; - } + if ((diag = Diag[I]) != NULL) + diag->Real += Gmin; + } } - return ; + return; } @@ -1810,44 +1857,47 @@ LoadGmin (SMPmatrix *eMatrix, double Gmin) * pointer to the new element is returned. */ -SMPelement * -SMPfindElt (SMPmatrix *eMatrix, int Row, int Col, int CreateIfMissing) +SMPelement* +SMPfindElt(SMPmatrix* eMatrix, int Row, int Col, int CreateIfMissing) { - MatrixPtr Matrix = eMatrix->SPmatrix ; + MatrixPtr Matrix = eMatrix->SPmatrix; if (eMatrix->CKTkluMODE) { - int i ; + int i; - Row = Row - 1 ; - Col = Col - 1 ; + Row = Row - 1; + Col = Col - 1; if (Col < 0) { -// printf ("Information: Cannot find an element with row '%d' and column '%d' in the KLU matrix\n", Row, Col) ; - return NULL ; + // printf ("Information: Cannot find an element with row '%d' and column '%d' in the KLU matrix\n", Row, Col) ; + return NULL; } - for (i = eMatrix->SMPkluMatrix->KLUmatrixAp [Col] ; i < eMatrix->SMPkluMatrix->KLUmatrixAp [Col + 1] ; i++) { - if (eMatrix->SMPkluMatrix->KLUmatrixAi [i] == Row) { + for (i = eMatrix->SMPkluMatrix->KLUmatrixAp[Col]; i < eMatrix->SMPkluMatrix->KLUmatrixAp[Col + 1]; i++) { + if (eMatrix->SMPkluMatrix->KLUmatrixAi[i] == Row) { if (eMatrix->SMPkluMatrix->KLUmatrixIsComplex == KLUmatrixReal) { - return (SMPelement *) &(eMatrix->SMPkluMatrix->KLUmatrixAx [i]) ; - } else if (eMatrix->SMPkluMatrix->KLUmatrixIsComplex == KLUMatrixComplex) { - return (SMPelement *) &(eMatrix->SMPkluMatrix->KLUmatrixAxComplex [2 * i]) ; - } else { - printf ("Information: Cannot find an element with row '%d' and column '%d' in the KLU matrix\n", Row, Col) ; - return NULL ; + return (SMPelement*)&(eMatrix->SMPkluMatrix->KLUmatrixAx[i]); + } + else if (eMatrix->SMPkluMatrix->KLUmatrixIsComplex == KLUMatrixComplex) { + return (SMPelement*)&(eMatrix->SMPkluMatrix->KLUmatrixAxComplex[2 * i]); + } + else { + printf("Information: Cannot find an element with row '%d' and column '%d' in the KLU matrix\n", Row, Col); + return NULL; } } } - return NULL ; - } else { - ElementPtr Element ; + return NULL; + } + else { + ElementPtr Element; /* Begin `SMPfindElt'. */ - assert (IS_SPARSE (Matrix)) ; - Row = Matrix->ExtToIntRowMap [Row] ; - Col = Matrix->ExtToIntColMap [Col] ; - Element = Matrix->FirstInCol [Col] ; - Element = spcFindElementInCol (Matrix, &Element, Row, Col, CreateIfMissing) ; - return (SMPelement *)Element ; + assert(IS_SPARSE(Matrix)); + Row = Matrix->ExtToIntRowMap[Row]; + Col = Matrix->ExtToIntColMap[Col]; + Element = Matrix->FirstInCol[Col]; + Element = spcFindElementInCol(Matrix, &Element, Row, Col, CreateIfMissing); + return (SMPelement*)Element; } } @@ -1857,28 +1907,29 @@ SMPfindElt (SMPmatrix *eMatrix, int Row, int Col, int CreateIfMissing) * SMPcZeroCol() */ int -SMPcZeroCol (SMPmatrix *eMatrix, int Col) +SMPcZeroCol(SMPmatrix* eMatrix, int Col) { - MatrixPtr Matrix = eMatrix->SPmatrix ; - ElementPtr Element ; + MatrixPtr Matrix = eMatrix->SPmatrix; + ElementPtr Element; if (eMatrix->CKTkluMODE) { - int i ; - for (i = eMatrix->SMPkluMatrix->KLUmatrixAp [Col - 1] ; i < eMatrix->SMPkluMatrix->KLUmatrixAp [Col] ; i++) + int i; + for (i = eMatrix->SMPkluMatrix->KLUmatrixAp[Col - 1]; i < eMatrix->SMPkluMatrix->KLUmatrixAp[Col]; i++) { - eMatrix->SMPkluMatrix->KLUmatrixAxComplex [2 * i] = 0 ; - eMatrix->SMPkluMatrix->KLUmatrixAxComplex [2 * i + 1] = 0 ; + eMatrix->SMPkluMatrix->KLUmatrixAxComplex[2 * i] = 0; + eMatrix->SMPkluMatrix->KLUmatrixAxComplex[2 * i + 1] = 0; } - return 0 ; - } else { - Col = Matrix->ExtToIntColMap [Col] ; - for (Element = Matrix->FirstInCol [Col] ; Element != NULL ; Element = Element->NextInCol) + return 0; + } + else { + Col = Matrix->ExtToIntColMap[Col]; + for (Element = Matrix->FirstInCol[Col]; Element != NULL; Element = Element->NextInCol) { - Element->Real = 0.0 ; - Element->Imag = 0.0 ; + Element->Real = 0.0; + Element->Imag = 0.0; } - return spError (Matrix) ; + return spError(Matrix); } } @@ -1886,86 +1937,89 @@ SMPcZeroCol (SMPmatrix *eMatrix, int Col) * SMPcAddCol() */ int -SMPcAddCol (SMPmatrix *eMatrix, int Accum_Col, int Addend_Col) +SMPcAddCol(SMPmatrix* eMatrix, int Accum_Col, int Addend_Col) { - MatrixPtr Matrix = eMatrix->SPmatrix ; - ElementPtr Accum, Addend, *Prev ; + MatrixPtr Matrix = eMatrix->SPmatrix; + ElementPtr Accum, Addend, * Prev; - Accum_Col = Matrix->ExtToIntColMap [Accum_Col] ; - Addend_Col = Matrix->ExtToIntColMap [Addend_Col] ; + Accum_Col = Matrix->ExtToIntColMap[Accum_Col]; + Addend_Col = Matrix->ExtToIntColMap[Addend_Col]; - Addend = Matrix->FirstInCol [Addend_Col] ; - Prev = &Matrix->FirstInCol [Accum_Col] ; + Addend = Matrix->FirstInCol[Addend_Col]; + Prev = &Matrix->FirstInCol[Accum_Col]; Accum = *Prev; while (Addend != NULL) { - while (Accum && Accum->Row < Addend->Row) + while (Accum && Accum->Row < Addend->Row) { - Prev = &Accum->NextInCol ; - Accum = *Prev ; - } - if (!Accum || Accum->Row > Addend->Row) + Prev = &Accum->NextInCol; + Accum = *Prev; + } + if (!Accum || Accum->Row > Addend->Row) { - Accum = spcCreateElement (Matrix, Addend->Row, Accum_Col, Prev, 0) ; - } - Accum->Real += Addend->Real ; - Accum->Imag += Addend->Imag ; - Addend = Addend->NextInCol ; + Accum = spcCreateElement(Matrix, Addend->Row, Accum_Col, Prev, 0); + } + Accum->Real += Addend->Real; + Accum->Imag += Addend->Imag; + Addend = Addend->NextInCol; } - return spError (Matrix) ; + return spError(Matrix); } /* * SMPzeroRow() */ int -SMPzeroRow (SMPmatrix *eMatrix, int Row) +SMPzeroRow(SMPmatrix* eMatrix, int Row) { - MatrixPtr Matrix = eMatrix->SPmatrix ; - ElementPtr Element ; + MatrixPtr Matrix = eMatrix->SPmatrix; + ElementPtr Element; - Row = Matrix->ExtToIntColMap [Row] ; + Row = Matrix->ExtToIntColMap[Row]; if (Matrix->RowsLinked == NO) - spcLinkRows (Matrix) ; + spcLinkRows(Matrix); if (Matrix->PreviousMatrixWasComplex || Matrix->Complex) { - for (Element = Matrix->FirstInRow[Row] ; Element != NULL; Element = Element->NextInRow) - { - Element->Real = 0.0 ; - Element->Imag = 0.0 ; - } - } else { - for (Element = Matrix->FirstInRow [Row] ; Element != NULL ; Element = Element->NextInRow) - { - Element->Real = 0.0 ; - } + for (Element = Matrix->FirstInRow[Row]; Element != NULL; Element = Element->NextInRow) + { + Element->Real = 0.0; + Element->Imag = 0.0; + } + } + else { + for (Element = Matrix->FirstInRow[Row]; Element != NULL; Element = Element->NextInRow) + { + Element->Real = 0.0; + } } - return spError (Matrix) ; + return spError(Matrix); } /* * SMPconstMult() */ void -SMPconstMult (SMPmatrix *Matrix, double constant) +SMPconstMult(SMPmatrix* Matrix, double constant) { if (Matrix->CKTkluMODE) { if (Matrix->SMPkluMatrix->KLUmatrixIsComplex) { - klu_z_constant_multiply (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAxComplex, - (int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixCommon, constant) ; - } else { - klu_constant_multiply (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAx, - (int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixCommon, constant) ; + klu_z_constant_multiply(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAxComplex, + (int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixCommon, constant); } - } else { - spConstMult (Matrix->SPmatrix, constant) ; + else { + klu_constant_multiply(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAx, + (int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixCommon, constant); + } + } + else { + spConstMult(Matrix->SPmatrix, constant); } } @@ -1973,37 +2027,39 @@ SMPconstMult (SMPmatrix *Matrix, double constant) * SMPmultiply() */ void -SMPmultiply (SMPmatrix *Matrix, double *RHS, double *Solution, double *iRHS, double *iSolution) +SMPmultiply(SMPmatrix* Matrix, double* RHS, double* Solution, double* iRHS, double* iSolution) { if (Matrix->CKTkluMODE) { - int *Ap_CSR, *Ai_CSR ; - double *Ax_CSR ; + int* Ap_CSR, * Ai_CSR; + double* Ax_CSR; - Ap_CSR = (int *) malloc ((size_t)(Matrix->SMPkluMatrix->KLUmatrixN + 1) * sizeof (int)) ; - Ai_CSR = (int *) malloc ((size_t)Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof (int)) ; + Ap_CSR = (int*)malloc((size_t)(Matrix->SMPkluMatrix->KLUmatrixN + 1) * sizeof(int)); + Ai_CSR = (int*)malloc((size_t)Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(int)); if (Matrix->SMPkluMatrix->KLUmatrixIsComplex) { - Ax_CSR = (double *) malloc ((size_t)(2 * Matrix->SMPkluMatrix->KLUmatrixNZ) * sizeof (double)) ; - klu_z_convert_matrix_in_CSR (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, Ap_CSR, - Ai_CSR, Ax_CSR, (int)Matrix->SMPkluMatrix->KLUmatrixN, (int)Matrix->SMPkluMatrix->KLUmatrixNZ, Matrix->SMPkluMatrix->KLUmatrixCommon) ; - klu_z_matrix_vector_multiply (Ap_CSR, Ai_CSR, Ax_CSR, RHS, Solution, iRHS, iSolution, NULL, NULL, - (int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixCommon) ; - } else { - Ax_CSR = (double *) malloc ((size_t)Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof (double)) ; - klu_convert_matrix_in_CSR (Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAx, Ap_CSR, Ai_CSR, - Ax_CSR, (int)Matrix->SMPkluMatrix->KLUmatrixN, (int)Matrix->SMPkluMatrix->KLUmatrixNZ, Matrix->SMPkluMatrix->KLUmatrixCommon) ; - klu_matrix_vector_multiply (Ap_CSR, Ai_CSR, Ax_CSR, RHS, Solution, NULL, NULL, - (int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixCommon) ; - iSolution = iRHS ; + Ax_CSR = (double*)malloc((size_t)(2 * Matrix->SMPkluMatrix->KLUmatrixNZ) * sizeof(double)); + klu_z_convert_matrix_in_CSR(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAxComplex, Ap_CSR, + Ai_CSR, Ax_CSR, (int)Matrix->SMPkluMatrix->KLUmatrixN, (int)Matrix->SMPkluMatrix->KLUmatrixNZ, Matrix->SMPkluMatrix->KLUmatrixCommon); + klu_z_matrix_vector_multiply(Ap_CSR, Ai_CSR, Ax_CSR, RHS, Solution, iRHS, iSolution, NULL, NULL, + (int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixCommon); + } + else { + Ax_CSR = (double*)malloc((size_t)Matrix->SMPkluMatrix->KLUmatrixNZ * sizeof(double)); + klu_convert_matrix_in_CSR(Matrix->SMPkluMatrix->KLUmatrixAp, Matrix->SMPkluMatrix->KLUmatrixAi, Matrix->SMPkluMatrix->KLUmatrixAx, Ap_CSR, Ai_CSR, + Ax_CSR, (int)Matrix->SMPkluMatrix->KLUmatrixN, (int)Matrix->SMPkluMatrix->KLUmatrixNZ, Matrix->SMPkluMatrix->KLUmatrixCommon); + klu_matrix_vector_multiply(Ap_CSR, Ai_CSR, Ax_CSR, RHS, Solution, NULL, NULL, + (int)Matrix->SMPkluMatrix->KLUmatrixN, Matrix->SMPkluMatrix->KLUmatrixCommon); + iSolution = iRHS; } - free (Ap_CSR) ; - free (Ai_CSR) ; - free (Ax_CSR) ; - } else { - spMultiply (Matrix->SPmatrix, RHS, Solution, iRHS, iSolution) ; + free(Ap_CSR); + free(Ai_CSR); + free(Ax_CSR); + } + else { + spMultiply(Matrix->SPmatrix, RHS, Solution, iRHS, iSolution); } }