Actual source code: arpackp.h
slepc-3.6.1 2015-09-03
1: /*
2: Private data structure used by the ARPACK interface
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
27: typedef struct {
28: PetscBool *select;
29: PetscScalar *workev;
30: PetscScalar *workd;
31: PetscScalar *workl;
32: PetscBLASInt lworkl;
33: PetscReal *rwork;
34: } EPS_ARPACK;
36: /*
37: Definition of routines from the ARPACK package
38: */
40: #if defined(PETSC_HAVE_MPIUNI)
42: #if defined(PETSC_USE_COMPLEX)
44: #if defined(PETSC_USE_REAL_SINGLE)
46: #if defined(SLEPC_ARPACK_HAVE_UNDERSCORE)
47: #define SLEPC_ARPACK(lcase,ucase) c##lcase##_
48: #elif defined(SLEPC_ARPACK_HAVE_CAPS)
49: #define SLEPC_ARPACK(lcase,ucase) C##ucase
50: #else
51: #define SLEPC_ARPACK(lcase,ucase) c##lcase
52: #endif
54: #else
56: #if defined(SLEPC_ARPACK_HAVE_UNDERSCORE)
57: #define SLEPC_ARPACK(lcase,ucase) z##lcase##_
58: #elif defined(SLEPC_ARPACK_HAVE_CAPS)
59: #define SLEPC_ARPACK(lcase,ucase) Z##ucase
60: #else
61: #define SLEPC_ARPACK(lcase,ucase) z##lcase
62: #endif
64: #endif
66: #else
68: #if defined(PETSC_USE_REAL_SINGLE)
70: #if defined(SLEPC_ARPACK_HAVE_UNDERSCORE)
71: #define SLEPC_ARPACK(lcase,ucase) s##lcase##_
72: #elif defined(SLEPC_ARPACK_HAVE_CAPS)
73: #define SLEPC_ARPACK(lcase,ucase) S##ucase
74: #else
75: #define SLEPC_ARPACK(lcase,ucase) s##lcase
76: #endif
78: #else
80: #if defined(SLEPC_ARPACK_HAVE_UNDERSCORE)
81: #define SLEPC_ARPACK(lcase,ucase) d##lcase##_
82: #elif defined(SLEPC_ARPACK_HAVE_CAPS)
83: #define SLEPC_ARPACK(lcase,ucase) D##ucase
84: #else
85: #define SLEPC_ARPACK(lcase,ucase) d##lcase
86: #endif
88: #endif
90: #endif
92: #else /* not MPIUNI */
94: #if defined(PETSC_USE_COMPLEX)
96: #if defined(PETSC_USE_REAL_SINGLE)
98: #if defined(SLEPC_ARPACK_HAVE_UNDERSCORE)
99: #define SLEPC_ARPACK(lcase,ucase) pc##lcase##_
100: #elif defined(SLEPC_ARPACK_HAVE_CAPS)
101: #define SLEPC_ARPACK(lcase,ucase) PC##ucase
102: #else
103: #define SLEPC_ARPACK(lcase,ucase) pc##lcase
104: #endif
106: #else
108: #if defined(SLEPC_ARPACK_HAVE_UNDERSCORE)
109: #define SLEPC_ARPACK(lcase,ucase) pz##lcase##_
110: #elif defined(SLEPC_ARPACK_HAVE_CAPS)
111: #define SLEPC_ARPACK(lcase,ucase) PZ##ucase
112: #else
113: #define SLEPC_ARPACK(lcase,ucase) pz##lcase
114: #endif
116: #endif
118: #else
120: #if defined(PETSC_USE_REAL_SINGLE)
122: #if defined(SLEPC_ARPACK_HAVE_UNDERSCORE)
123: #define SLEPC_ARPACK(lcase,ucase) ps##lcase##_
124: #elif defined(SLEPC_ARPACK_HAVE_CAPS)
125: #define SLEPC_ARPACK(lcase,ucase) PS##ucase
126: #else
127: #define SLEPC_ARPACK(lcase,ucase) ps##lcase
128: #endif
130: #else
132: #if defined(SLEPC_ARPACK_HAVE_UNDERSCORE)
133: #define SLEPC_ARPACK(lcase,ucase) pd##lcase##_
134: #elif defined(SLEPC_ARPACK_HAVE_CAPS)
135: #define SLEPC_ARPACK(lcase,ucase) PD##ucase
136: #else
137: #define SLEPC_ARPACK(lcase,ucase) pd##lcase
138: #endif
140: #endif
142: #endif
144: #endif
146: #if defined(PETSC_HAVE_MPIUNI)
148: #define COMM_ARG
150: #if !defined(PETSC_USE_COMPLEX)
152: #define ARPACKnaupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p) SLEPC_ARPACK(naupd,NAUPD) ((a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),1,2)
153: #define ARPACKneupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q,r,s,t,u,v,w,x,y) SLEPC_ARPACK(neupd,NEUPD) ((a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),(q),(r),(s),(t),(u),(v),(w),(x),(y),1,1,2)
154: #define ARPACKsaupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p) SLEPC_ARPACK(saupd,SAUPD) ((a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),1,2)
155: #define ARPACKseupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,o,p,q,r,s,t,u,v,w) SLEPC_ARPACK(seupd,SEUPD) ((a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(o),(p),(q),(r),(s),(t),(u),(v),(w),1,1,2)
157: #else
159: #define ARPACKnaupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q) SLEPC_ARPACK(naupd,NAUPD) ((a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),(q),1,2)
160: #define ARPACKneupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q,r,s,t,u,v,w,x) SLEPC_ARPACK(neupd,NEUPD) ((a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),(q),(r),(s),(t),(u),(v),(w),(x),1,1,2)
162: #endif
164: #else /* not MPIUNI */
166: #define COMM_ARG MPI_Fint*,
168: #if !defined(PETSC_USE_COMPLEX)
170: #define ARPACKnaupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p) SLEPC_ARPACK(naupd,NAUPD) ((comm),(a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),1,2)
171: #define ARPACKneupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q,r,s,t,u,v,w,x,y) SLEPC_ARPACK(neupd,NEUPD) ((comm),(a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),(q),(r),(s),(t),(u),(v),(w),(x),(y),1,1,2)
172: #define ARPACKsaupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p) SLEPC_ARPACK(saupd,SAUPD) ((comm),(a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),1,2)
173: #define ARPACKseupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,o,p,q,r,s,t,u,v,w) SLEPC_ARPACK(seupd,SEUPD) ((comm),(a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(o),(p),(q),(r),(s),(t),(u),(v),(w),1,1,2)
175: #else
177: #define ARPACKnaupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q) SLEPC_ARPACK(naupd,NAUPD) ((comm),(a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),(q),1,2)
178: #define ARPACKneupd_(comm,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q,r,s,t,u,v,w,x) SLEPC_ARPACK(neupd,NEUPD) ((comm),(a),(b),(c),(d),(e),(f),(g),(h),(i),(j),(k),(l),(m),(n),(o),(p),(q),(r),(s),(t),(u),(v),(w),(x),1,1,2)
180: #endif
182: #endif
184: PETSC_EXTERN void SLEPC_ARPACK(saupd,SAUPD)(COMM_ARG PetscBLASInt*,char*,PetscBLASInt*,const char*,PetscBLASInt*,PetscReal*,PetscScalar*,PetscBLASInt*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,PetscBLASInt*,PetscScalar*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,int,int);
185: PETSC_EXTERN void SLEPC_ARPACK(seupd,SEUPD)(COMM_ARG PetscBool*,char*,PetscBool*,PetscReal*,PetscReal*,PetscBLASInt*,PetscReal*,char*,PetscBLASInt*,const char*,PetscBLASInt*,PetscReal*,PetscScalar*,PetscBLASInt*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,PetscBLASInt*,PetscScalar*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,int,int,int);
187: #if !defined(PETSC_USE_COMPLEX)
188: PETSC_EXTERN void SLEPC_ARPACK(naupd,NAUPD)(COMM_ARG PetscBLASInt*,char*,PetscBLASInt*,const char*,PetscBLASInt*,PetscReal*,PetscScalar*,PetscBLASInt*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,PetscBLASInt*,PetscScalar*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,int,int);
189: PETSC_EXTERN void SLEPC_ARPACK(neupd,NEUPD)(COMM_ARG PetscBool*,char*,PetscBool*,PetscReal*,PetscReal*,PetscReal*,PetscBLASInt*,PetscReal*,PetscReal*,PetscReal*,char*,PetscBLASInt*,const char*,PetscBLASInt*,PetscReal*,PetscScalar*,PetscBLASInt*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,PetscBLASInt*,PetscScalar*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,int,int,int);
190: #else
191: PETSC_EXTERN void SLEPC_ARPACK(naupd,NAUPD)(COMM_ARG PetscBLASInt*,char*,PetscBLASInt*,const char*,PetscBLASInt*,PetscReal*,PetscScalar*,PetscBLASInt*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,PetscBLASInt*,PetscScalar*,PetscScalar*,PetscBLASInt*,PetscReal*,PetscBLASInt*,int,int);
192: PETSC_EXTERN void SLEPC_ARPACK(neupd,NEUPD)(COMM_ARG PetscBool*,char*,PetscBool*,PetscScalar*,PetscScalar*,PetscBLASInt*,PetscScalar*,PetscScalar*,char*,PetscBLASInt*,const char*,PetscBLASInt*,PetscReal*,PetscScalar*,PetscBLASInt*,PetscScalar*,PetscBLASInt*,PetscBLASInt*,PetscBLASInt*,PetscScalar*,PetscScalar*,PetscBLASInt*,PetscReal*,PetscBLASInt*,int,int,int);
193: #endif
195: #endif