Actual source code: linearp.h
1: /*
3: Private header for QEPLINEAR.
5: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
6: SLEPc - Scalable Library for Eigenvalue Problem Computations
7: Copyright (c) 2002-2012, Universitat Politecnica de Valencia, Spain
9: This file is part of SLEPc.
10:
11: SLEPc is free software: you can redistribute it and/or modify it under the
12: terms of version 3 of the GNU Lesser General Public License as published by
13: the Free Software Foundation.
15: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
16: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
18: more details.
20: You should have received a copy of the GNU Lesser General Public License
21: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
22: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
23: */
28: typedef struct {
29: PetscBool explicitmatrix;
30: PetscInt cform; /* companion form */
31: PetscReal sfactor; /* scaling factor */
32: Mat A,B; /* matrices of generalized eigenproblem */
33: EPS eps; /* linear eigensolver for Az=lBz */
34: Mat M,C,K; /* copy of QEP coefficient matrices */
35: Vec x1,x2,y1,y2; /* work vectors */
36: PetscBool setfromoptionscalled;
37: } QEP_LINEAR;
39: /* N1 */
40: extern PetscErrorCode MatMult_Linear_N1A(Mat,Vec,Vec);
41: extern PetscErrorCode MatMult_Linear_N1B(Mat,Vec,Vec);
42: extern PetscErrorCode MatGetDiagonal_Linear_N1A(Mat,Vec);
43: extern PetscErrorCode MatGetDiagonal_Linear_N1B(Mat,Vec);
44: extern PetscErrorCode MatCreateExplicit_Linear_N1A(MPI_Comm,QEP_LINEAR*,Mat*);
45: extern PetscErrorCode MatCreateExplicit_Linear_N1B(MPI_Comm,QEP_LINEAR*,Mat*);
47: /* N2 */
48: extern PetscErrorCode MatMult_Linear_N2A(Mat,Vec,Vec);
49: extern PetscErrorCode MatMult_Linear_N2B(Mat,Vec,Vec);
50: extern PetscErrorCode MatGetDiagonal_Linear_N2A(Mat,Vec);
51: extern PetscErrorCode MatGetDiagonal_Linear_N2B(Mat,Vec);
52: extern PetscErrorCode MatCreateExplicit_Linear_N2A(MPI_Comm,QEP_LINEAR*,Mat*);
53: extern PetscErrorCode MatCreateExplicit_Linear_N2B(MPI_Comm,QEP_LINEAR*,Mat*);
55: /* S1 */
56: extern PetscErrorCode MatMult_Linear_S1A(Mat,Vec,Vec);
57: extern PetscErrorCode MatMult_Linear_S1B(Mat,Vec,Vec);
58: extern PetscErrorCode MatGetDiagonal_Linear_S1A(Mat,Vec);
59: extern PetscErrorCode MatGetDiagonal_Linear_S1B(Mat,Vec);
60: extern PetscErrorCode MatCreateExplicit_Linear_S1A(MPI_Comm,QEP_LINEAR*,Mat*);
61: extern PetscErrorCode MatCreateExplicit_Linear_S1B(MPI_Comm,QEP_LINEAR*,Mat*);
63: /* S2 */
64: extern PetscErrorCode MatMult_Linear_S2A(Mat,Vec,Vec);
65: extern PetscErrorCode MatMult_Linear_S2B(Mat,Vec,Vec);
66: extern PetscErrorCode MatGetDiagonal_Linear_S2A(Mat,Vec);
67: extern PetscErrorCode MatGetDiagonal_Linear_S2B(Mat,Vec);
68: extern PetscErrorCode MatCreateExplicit_Linear_S2A(MPI_Comm,QEP_LINEAR*,Mat*);
69: extern PetscErrorCode MatCreateExplicit_Linear_S2B(MPI_Comm,QEP_LINEAR*,Mat*);
71: /* H1 */
72: extern PetscErrorCode MatMult_Linear_H1A(Mat,Vec,Vec);
73: extern PetscErrorCode MatMult_Linear_H1B(Mat,Vec,Vec);
74: extern PetscErrorCode MatGetDiagonal_Linear_H1A(Mat,Vec);
75: extern PetscErrorCode MatGetDiagonal_Linear_H1B(Mat,Vec);
76: extern PetscErrorCode MatCreateExplicit_Linear_H1A(MPI_Comm,QEP_LINEAR*,Mat*);
77: extern PetscErrorCode MatCreateExplicit_Linear_H1B(MPI_Comm,QEP_LINEAR*,Mat*);
79: /* H2 */
80: extern PetscErrorCode MatMult_Linear_H2A(Mat,Vec,Vec);
81: extern PetscErrorCode MatMult_Linear_H2B(Mat,Vec,Vec);
82: extern PetscErrorCode MatGetDiagonal_Linear_H2A(Mat,Vec);
83: extern PetscErrorCode MatGetDiagonal_Linear_H2B(Mat,Vec);
84: extern PetscErrorCode MatCreateExplicit_Linear_H2A(MPI_Comm,QEP_LINEAR*,Mat*);
85: extern PetscErrorCode MatCreateExplicit_Linear_H2B(MPI_Comm,QEP_LINEAR*,Mat*);
87: #endif