/* nag_sparse_sym_chol_fac (f11jac) Example Program. * * Copyright 1998 Numerical Algorithms Group. * * Mark 5, 1998. * */ #include #include #include #include #include int main(void) { double dtol; double *a; double dscale; Integer *irow, *icol; Integer *ipiv, nnzc, *istr; Integer exit_status = 0, i, n, lfill, npivm; Integer nnz; Integer num; char nag_enum_arg[40]; Nag_SparseSym_Piv pstrat; Nag_SparseSym_Fact mic; Nag_Sparse_Comm comm; NagError fail; INIT_FAIL(fail); printf( "nag_sparse_sym_chol_fac (f11jac) Example Program Results\n\n"); /* Skip heading in data file */ scanf("%*[^\n]"); /* Read algorithmic parameters */ scanf("%ld", &n); scanf("%*[^\n]"); scanf("%ld%*[^\n]", &nnz); scanf("%ld%lf%*[^\n]", &lfill, &dtol); scanf("%s%lf%*[^\n]", nag_enum_arg, &dscale); /* nag_enum_name_to_value(x04nac). * Converts NAG enum member name to value */ mic = (Nag_SparseSym_Fact) nag_enum_name_to_value(nag_enum_arg); scanf("%s%*[^\n]", nag_enum_arg); pstrat = (Nag_SparseNsym_Piv) nag_enum_name_to_value(nag_enum_arg); /* Allocate memory */ num = 2 * nnz; ipiv = NAG_ALLOC(n, Integer); istr = NAG_ALLOC(n+1, Integer); irow = NAG_ALLOC(num, Integer); icol = NAG_ALLOC(num, Integer); a = NAG_ALLOC(num, double); if (!ipiv || !istr || !irow || !icol || !a) { printf("Allocation failure\n"); exit_status = -1; goto END; } /* Read the matrix a */ for (i = 1; i <= nnz; ++i) scanf("%lf%ld%ld%*[^\n]", &a[i-1], &irow[i-1], &icol[i-1]); /* Calculate incomplete Cholesky factorization */ /* nag_sparse_sym_chol_fac (f11jac). * Incomplete Cholesky factorization (symmetric) */ nag_sparse_sym_chol_fac(n, nnz, &a, &num, &irow, &icol, lfill, dtol, mic, dscale, pstrat, ipiv, istr, &nnzc, &npivm, &comm, &fail); if (fail.code != NE_NOERROR) { printf("Error from nag_sparse_sym_chol_fac (f11jac).\n%s\n", fail.message); exit_status = 1; goto END; } /* Output original matrix */ printf(" Original Matrix \n"); printf(" n = %6ld\n", n); printf(" nnz = %6ld\n\n", nnz); for (i = 1; i <= nnz; ++i) printf(" %8ld%16.4e%8ld%8ld\n", i, a[i-1], irow[i-1], icol[i-1]); printf("\n"); /* Output details of the factorization */ printf(" Factorization\n n = %6ld \n nnz = %6ld\n", n, nnzc); printf(" npivm = %6ld\n\n", npivm); for (i = nnz + 1; i <= nnz + nnzc; ++i) printf(" %8ld%16.4e%8ld%8ld\n", i, a[i-1], irow[i-1], icol[i-1]); printf("\n i ipiv(i) \n"); for (i = 1; i <= n; ++i) printf(" %8ld%8ld\n", i, ipiv[i-1]); END: if (irow) NAG_FREE(irow); if (icol) NAG_FREE(icol); if (a) NAG_FREE(a); if (istr) NAG_FREE(istr); if (ipiv) NAG_FREE(ipiv); return exit_status; }