1 /*! \file
2 Copyright (c) 2003, The Regents of the University of California, through
3 Lawrence Berkeley National Laboratory (subject to receipt of any required
4 approvals from U.S. Dept. of Energy)
5
6 All rights reserved.
7
8 The source code is distributed under BSD license, see the file License.txt
9 at the top-level directory.
10 */
11
12 /*
13 * -- SuperLU routine (version 3.0) --
14 * Univ. of California Berkeley, Xerox Palo Alto Research Center,
15 * and Lawrence Berkeley National Lab.
16 * October 15, 2003
17 *
18 */
19 #include "slu_cdefs.h"
20
main(int argc,char * argv[])21 int main(int argc, char *argv[])
22 {
23 SuperMatrix A;
24 NCformat *Astore;
25 complex *a;
26 int *asub, *xa;
27 int *perm_c; /* column permutation vector */
28 int *perm_r; /* row permutations from partial pivoting */
29 SuperMatrix L; /* factor L */
30 SCformat *Lstore;
31 SuperMatrix U; /* factor U */
32 NCformat *Ustore;
33 SuperMatrix B;
34 int nrhs, ldx, info, m, n, nnz;
35 complex *xact, *rhs;
36 mem_usage_t mem_usage;
37 superlu_options_t options;
38 SuperLUStat_t stat;
39 FILE *fp = stdin;
40
41 #if ( DEBUGlevel>=1 )
42 CHECK_MALLOC("Enter main()");
43 #endif
44
45 /* Set the default input options:
46 options.Fact = DOFACT;
47 options.Equil = YES;
48 options.ColPerm = COLAMD;
49 options.DiagPivotThresh = 1.0;
50 options.Trans = NOTRANS;
51 options.IterRefine = NOREFINE;
52 options.SymmetricMode = NO;
53 options.PivotGrowth = NO;
54 options.ConditionNumber = NO;
55 options.PrintStat = YES;
56 */
57 set_default_options(&options);
58
59 /* Read the matrix in Harwell-Boeing format. */
60 creadhb(fp, &m, &n, &nnz, &a, &asub, &xa);
61
62 cCreate_CompCol_Matrix(&A, m, n, nnz, a, asub, xa, SLU_NC, SLU_C, SLU_GE);
63 Astore = A.Store;
64 printf("Dimension %dx%d; # nonzeros %d\n", A.nrow, A.ncol, Astore->nnz);
65
66 nrhs = 1;
67 if ( !(rhs = complexMalloc(m * nrhs)) ) ABORT("Malloc fails for rhs[].");
68 cCreate_Dense_Matrix(&B, m, nrhs, rhs, m, SLU_DN, SLU_C, SLU_GE);
69 xact = complexMalloc(n * nrhs);
70 ldx = n;
71 cGenXtrue(n, nrhs, xact, ldx);
72 cFillRHS(options.Trans, nrhs, xact, ldx, &A, &B);
73
74 if ( !(perm_c = intMalloc(n)) ) ABORT("Malloc fails for perm_c[].");
75 if ( !(perm_r = intMalloc(m)) ) ABORT("Malloc fails for perm_r[].");
76
77 /* Initialize the statistics variables. */
78 StatInit(&stat);
79
80 cgssv(&options, &A, perm_c, perm_r, &L, &U, &B, &stat, &info);
81
82 if ( info == 0 ) {
83
84 /* This is how you could access the solution matrix. */
85 complex *sol = (complex*) ((DNformat*) B.Store)->nzval;
86
87 /* Compute the infinity norm of the error. */
88 cinf_norm_error(nrhs, &B, xact);
89
90 Lstore = (SCformat *) L.Store;
91 Ustore = (NCformat *) U.Store;
92 printf("No of nonzeros in factor L = %d\n", Lstore->nnz);
93 printf("No of nonzeros in factor U = %d\n", Ustore->nnz);
94 printf("No of nonzeros in L+U = %d\n", Lstore->nnz + Ustore->nnz - n);
95 printf("FILL ratio = %.1f\n", (float)(Lstore->nnz + Ustore->nnz - n)/nnz);
96
97 cQuerySpace(&L, &U, &mem_usage);
98 printf("L\\U MB %.3f\ttotal MB needed %.3f\n",
99 mem_usage.for_lu/1e6, mem_usage.total_needed/1e6);
100
101 } else {
102 printf("cgssv() error returns INFO= %d\n", info);
103 if ( info <= n ) { /* factorization completes */
104 cQuerySpace(&L, &U, &mem_usage);
105 printf("L\\U MB %.3f\ttotal MB needed %.3f\n",
106 mem_usage.for_lu/1e6, mem_usage.total_needed/1e6);
107 }
108 }
109
110 if ( options.PrintStat ) StatPrint(&stat);
111 StatFree(&stat);
112
113 SUPERLU_FREE (rhs);
114 SUPERLU_FREE (xact);
115 SUPERLU_FREE (perm_r);
116 SUPERLU_FREE (perm_c);
117 Destroy_CompCol_Matrix(&A);
118 Destroy_SuperMatrix_Store(&B);
119 Destroy_SuperNode_Matrix(&L);
120 Destroy_CompCol_Matrix(&U);
121
122 #if ( DEBUGlevel>=1 )
123 CHECK_MALLOC("Exit main()");
124 #endif
125 }
126
127