IFPACK Development
Loading...
Searching...
No Matches
Ifpack_IKLU_Utils.h
1/*@HEADER
2// ***********************************************************************
3//
4// Ifpack: Object-Oriented Algebraic Preconditioner Package
5// Copyright (2002) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ***********************************************************************
40//@HEADER
41*/
42
43#ifndef IFPACK_IKLU_UTILS_H
44#define IFPACK_IKLU_UTILS_H
45
46#include <stdlib.h>
47#include <limits.h>
48#include <math.h>
49#include <stdio.h>
50
51/* The code found in this file is adapted from CSparse Version 2.0.0
52 written by Tim Davis, UFL
53*/
54
55/* --- primary CSparse routines and data structures ------------------------- */
56
57typedef struct row_matrix /* matrix in compressed-row or triplet form */
58{
59 int nzmax ; /* maximum number of entries */
60 int m ; /* number of rows */
61 int n ; /* number of columns */
62 int *p ; /* row pointers (size m+1) or col indices (size nzmax) */
63 int *j ; /* col indices, size nzmax */
64 double *x ; /* numerical values, size nzmax */
65 int nz ; /* # of entries in triplet matrix, -1 for compressed-row */
66} csr ;
67
68csr *csr_add (const csr *A, const csr *B, double alpha, double beta) ;
69csr *csr_multiply (const csr *A, const csr *B) ;
70double csr_norm (const csr *A) ;
71int csr_print (const csr *A, int brief) ;
72csr *csr_transpose (const csr *A, int values) ;
73
74/* utilities */
75void *csr_realloc (void *p, int n, size_t size, int *ok) ;
76
77/* csr utilities */
78csr *csr_spalloc (int m, int n, int nzmax, int values, int triplet) ;
79csr *csr_spfree (csr *A) ;
80int csr_sprealloc (csr *A, int nzmax) ;
81
82
83
84/* --- secondary CSparse routines and data structures ----------------------- */
85typedef struct cs_symbolic /* symbolic Cholesky, LU, or QR analysis */
86{
87 int *pinv ; /* inverse row perm. for QR, fill red. perm for Chol */
88 int *q ; /* fill-reducing column permutation for LU and QR */
89 int *parent ; /* elimination tree for Cholesky and QR */
90 int *cp ; /* column pointers for Cholesky, row counts for QR */
91 int *leftmost ; /* leftmost[i] = min(find(A(i,:))), for QR */
92 int m2 ; /* # of rows for QR, after adding fictitious rows */
93 double lnz ; /* # entries in L for LU or Cholesky; in V for QR */
94 double unz ; /* # entries in U for LU; in R for QR */
95} css ;
96
97typedef struct csr_numeric /* numeric Cholesky, LU, or QR factorization */
98{
99 csr *L ; /* L for LU and Cholesky, V for QR */
100 csr *U ; /* U for LU, R for QR, not used for Cholesky */
101 int *pinv ; /* partial pivoting for LU */
102 int *perm ; /* partial pivoting for LU */
103 double *B ; /* beta [0..n-1] for QR */
104} csrn ;
105
106typedef struct csr_dmperm_results /* csr_dmperm or csr_scc output */
107{
108 int *p ; /* size m, row permutation */ /* may be back wards */
109 int *q ; /* size n, column permutation */
110 int *r ; /* size nb+1, block k is rows r[k] to r[k+1]-1 in A(p,q) */
111 int *s ; /* size nb+1, block k is cols s[k] to s[k+1]-1 in A(p,q) */
112 int nb ; /* # of blocks in fine dmperm decomposition */
113 int rr [5] ; /* coarse row decomposition */
114 int cc [5] ; /* coarse column decomposition */
115} csrd ;
116
117int *csr_amd (int order, const csr *A) ;
118
119int csr_droptol (csr *A, double tol);
120int csr_dropzeros (csr *A);
121int csr_lsolve (const csr *L, double *x);
122csrn *csr_lu (const csr *A, const css *S, double tol);
123csr *csr_permute (const csr *A, const int *pinv, const int *q, int values);
124css *csr_sqr (int order, const csr *A);
125int csr_usolve (const csr *U, double *x);
126
127/* utilities */
128css *csr_sfree (css *S) ;
129csrd *csr_dfree (csrd *D);
130csrn *csr_nfree (csrn *N);
131
132
133
134/* --- tertiary CSparse routines -------------------------------------------- */
135double csr_cumsum (int *p, int *c, int n) ;
136int csr_dfs (int j, csr *G, int top, int *xi, int *pstack, const int *pinv);
137int csr_reach (csr *G, const csr *B, int k, int *xi, const int *pinv);
138int csr_scatter (const csr *A, int j, double beta, int *w, double *x, int mark,
139 csr *C, int nz) ;
140csrd *csr_scc (csr *A);
141int csr_spsolve (csr *G, const csr *B, int k, int *xi,
142 double *x, const int *pinv, int up);
143int csr_tdfs (int j, int k, int *head, const int *next, int *post,
144 int *stack) ;
145/* utilities */
146csrd *csr_dalloc (int m, int n);
147/* csr utilities */
148csrd *csr_ddone (csrd *D, csr *C, void *w, int ok) ;
149csr *csr_done (csr *C, void *w, void *x, int ok) ;
150int *csr_idone (int *p, csr *C, void *w, int ok) ;
151csrn *csr_ndone (csrn *N, csr *C, void *w, void *x, int ok) ;
152
153int csr_fkeep (csr *A, int (*fkeep) (int, int, double, void *), void *other);
154
155#define CS_MAX(a,b) (((a) > (b)) ? (a) : (b))
156#define CS_MIN(a,b) (((a) < (b)) ? (a) : (b))
157#define CS_FLIP(i) (-(i)-2)
158#define CS_UNFLIP(i) (((i) < 0) ? CS_FLIP(i) : (i))
159#define CS_MARKED(w,j) (w [j] < 0)
160#define CS_MARK(w,j) { w [j] = CS_FLIP (w [j]) ; }
161#define CS_CSC(A) (A && (A->nz == -1))
162#define CS_TRIPLET(A) (A && (A->nz >= 0))
163
164#endif // IFPACK_IKLU_UTILS_H