Actual source code: slepcbv.h

  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    User interface for the basis vectors object in SLEPc
 12: */

 14: #pragma once

 16: #include <slepcsys.h>

 18: /* SUBMANSEC = BV */

 20: SLEPC_EXTERN PetscErrorCode BVInitializePackage(void);
 21: SLEPC_EXTERN PetscErrorCode BVFinalizePackage(void);

 23: /*S
 24:     BV - Basis vectors, SLEPc object representing a collection of vectors
 25:     that typically constitute a basis of a subspace.

 27:     Level: beginner

 29: .seealso: [](sec:bv), `BVCreate()`
 30: S*/
 31: typedef struct _p_BV* BV;

 33: /*J
 34:     BVType - String with the name of the type of `BV`. Each type differs in
 35:     the way data is stored internally.

 37:     Level: beginner

 39: .seealso: [](sec:bv), `BVSetType()`, `BV`
 40: J*/
 41: typedef const char *BVType;
 42: #define BVMAT        "mat"
 43: #define BVSVEC       "svec"
 44: #define BVVECS       "vecs"
 45: #define BVCONTIGUOUS "contiguous"
 46: #define BVTENSOR     "tensor"

 48: /* Logging support */
 49: SLEPC_EXTERN PetscClassId BV_CLASSID;

 51: /*E
 52:    BVOrthogType - Determines the method used in the orthogonalization of vectors.

 54:    Values:
 55: +  `BV_ORTHOG_CGS` - Classical Gram-Schmidt
 56: -  `BV_ORTHOG_MGS` - Modified Gram-Schmidt

 58:    Note:
 59:    The default is to use CGS. Both CGS and MGS provide similar numerical accuracy when
 60:    combined with iterative refinement (the default). MGS results in poorer performance,
 61:    both sequential and in parallel, so it is not recommended unless iterative
 62:    refinement is disabled.

 64:    Level: advanced

 66: .seealso: [](sec:bv), `BVSetOrthogonalization()`, `BVGetOrthogonalization()`, `BVOrthogonalizeColumn()`, `BVOrthogRefineType`
 67: E*/
 68: typedef enum { BV_ORTHOG_CGS,
 69:                BV_ORTHOG_MGS } BVOrthogType;
 70: SLEPC_EXTERN const char *BVOrthogTypes[];

 72: /*MC
 73:    BV_ORTHOG_CGS - Use Classical Gram-Schmidt for orthogonalization of vectors.

 75:    Level: advanced

 77:    Note:
 78:    Although CGS is not numerically stable, when combined with iterative refinement
 79:    it yields a sufficiently accurate result.

 81: .seealso: [](sec:bv), `BVOrthogType`, `BVSetOrthogonalization()`, `BV_ORTHOG_MGS`
 82: M*/

 84: /*MC
 85:    BV_ORTHOG_MGS - Use Modified Gram-Schmidt for orthogonalization of vectors.

 87:    Level: advanced

 89:    Note:
 90:    MGS results in poorer performance, both sequential and in parallel, so it is
 91:    not recommended unless iterative refinement is disabled.

 93: .seealso: [](sec:bv), `BVOrthogType`, `BVSetOrthogonalization()`, `BV_ORTHOG_CGS`
 94: M*/

 96: /*E
 97:    BVOrthogRefineType - Determines what type of iterative refinement to use
 98:    during orthogonalization of vectors.

100:    Values:
101: +  `BV_ORTHOG_REFINE_IFNEEDED` - Refine only if a certain criterion is satisfied
102: .  `BV_ORTHOG_REFINE_NEVER`    - Never refine, do plain CGS or MGS
103: -  `BV_ORTHOG_REFINE_ALWAYS`   - Always refine, i.e., run CGS2 or MGS2

105:    Notes:
106:    The default is to perform one step of iterative refinement if a certain
107:    numerical criterion holds. In ill-conditioned cases, a third orthogonalization
108:    can also be done.

110:    Never refining is not recommended because it will generally make the eigensolver
111:    numerically unstable.

113:    Always refining is numerically stable, but it often performs unnecessary
114:    computation.

116:    Level: advanced

118: .seealso: [](sec:bv), `BVSetOrthogonalization()`, `BVGetOrthogonalization()`, `BVOrthogonalizeColumn()`
119: E*/
120: typedef enum { BV_ORTHOG_REFINE_IFNEEDED,
121:                BV_ORTHOG_REFINE_NEVER,
122:                BV_ORTHOG_REFINE_ALWAYS } BVOrthogRefineType;
123: SLEPC_EXTERN const char *BVOrthogRefineTypes[];

125: /*MC
126:    BV_ORTHOG_REFINE_IFNEEDED - In the orthogonalization of vectors, do iterative
127:    refinement only if a certain criterion is satisfied.

129:    Note:
130:    This is the default. One step of iterative refinement will be done if a certain
131:    numerical criterion holds. In ill-conditioned cases, a third orthogonalization
132:    can also be done. The criterion is similar to the one proposed in {cite:p}`Dan76`.
133:    The parameter $\eta$ for the criterion can be set in `BVSetOrthogonalization()`.

135:    Level: advanced

137: .seealso: [](sec:bv), `BVOrthogRefineType`, `BVSetOrthogonalization()`, `BV_ORTHOG_REFINE_NEVER`, `BV_ORTHOG_REFINE_ALWAYS`
138: M*/

140: /*MC
141:    BV_ORTHOG_REFINE_NEVER - In the orthogonalization of vectors, never do an iterative
142:    refinement step, i.e., do plain CGS or MGS.

144:    Note:
145:    Never refining is not recommended because it will generally make the eigensolver
146:    numerically unstable.

148:    Level: advanced

150: .seealso: [](sec:bv), `BVOrthogRefineType`, `BVSetOrthogonalization()`, `BV_ORTHOG_REFINE_IFNEEDED`, `BV_ORTHOG_REFINE_ALWAYS`
151: M*/

153: /*MC
154:    BV_ORTHOG_REFINE_ALWAYS - In the orthogonalization of vectors, always do an iterative
155:    refinement step, i.e., run CGS2 or MGS2.

157:    Note:
158:    Always refining is numerically stable, but it often performs unnecessary
159:    computation.

161:    Level: advanced

163: .seealso: [](sec:bv), `BVOrthogRefineType`, `BVSetOrthogonalization()`, `BV_ORTHOG_REFINE_IFNEEDED`, `BV_ORTHOG_REFINE_NEVER`
164: M*/

166: /*E
167:    BVOrthogBlockType - Determines the method used in block
168:    orthogonalization (simultaneous orthogonalization of a set of vectors).

170:    Values:
171: +  `BV_ORTHOG_BLOCK_GS`       - Gram-Schmidt, column by column
172: .  `BV_ORTHOG_BLOCK_CHOL`     - Cholesky QR method
173: .  `BV_ORTHOG_BLOCK_TSQR`     - Tall-Skinny QR method
174: .  `BV_ORTHOG_BLOCK_TSQRCHOL` - TSQR, but computing the triangular factor only
175: -  `BV_ORTHOG_BLOCK_SVQB`     - SVQB method

177:    Level: advanced

179: .seealso: [](sec:bv), `BVSetOrthogonalization()`, `BVGetOrthogonalization()`, `BVOrthogonalize()`
180: E*/
181: typedef enum { BV_ORTHOG_BLOCK_GS,
182:                BV_ORTHOG_BLOCK_CHOL,
183:                BV_ORTHOG_BLOCK_TSQR,
184:                BV_ORTHOG_BLOCK_TSQRCHOL,
185:                BV_ORTHOG_BLOCK_SVQB     } BVOrthogBlockType;
186: SLEPC_EXTERN const char *BVOrthogBlockTypes[];

188: /*E
189:    BVMatMultType - Different ways of performing the `BVMatMult()` operation.

191:    Values:
192:    Allowed values are
193: +  `BV_MATMULT_VECS` - perform a matrix-vector multiply per each column
194: .  `BV_MATMULT_MAT` - carry out a Mat-Mat product with a dense matrix
195: -  `BV_MATMULT_MAT_SAVE` - this case is deprecated

197:    Note:
198:    The default is `BV_MATMULT_MAT` except in the case of `BVVECS`.

200:    Level: advanced

202: .seealso: [](sec:bv), `BVSetMatMultMethod()`, `BVMatMult()`
203: E*/
204: typedef enum { BV_MATMULT_VECS,
205:                BV_MATMULT_MAT,
206:                BV_MATMULT_MAT_SAVE } BVMatMultType;
207: SLEPC_EXTERN const char *BVMatMultTypes[];

209: /*MC
210:    BV_MATMULT_VECS - Perform a matrix-vector multiply per each column.

212:    Level: advanced

214: .seealso: [](sec:bv), `BVMatMultType`, `BVSetMatMultMethod()`, `BVMatMult()`, `BV_MATMULT_MAT`, `BV_MATMULT_MAT_SAVE`
215: M*/

217: /*MC
218:    BV_MATMULT_MAT - Perform a `Mat`-`Mat` product with a dense matrix.

220:    Level: advanced

222: .seealso: [](sec:bv), `BVMatMultType`, `BVSetMatMultMethod()`, `BVMatMult()`, `BV_MATMULT_VECS`, `BV_MATMULT_MAT_SAVE`
223: M*/

225: /*MC
226:    BV_MATMULT_MAT_SAVE - This case is deprecated, do not use.

228:    Level: advanced

230: .seealso: [](sec:bv), `BVMatMultType`, `BVSetMatMultMethod()`, `BVMatMult()`, `BV_MATMULT_VECS`, `BV_MATMULT_MAT`
231: M*/

233: SLEPC_EXTERN PetscErrorCode BVCreate(MPI_Comm,BV*);
234: SLEPC_EXTERN PetscErrorCode BVDestroy(BV*);
235: SLEPC_EXTERN PetscErrorCode BVSetType(BV,BVType);
236: SLEPC_EXTERN PetscErrorCode BVGetType(BV,BVType*);
237: SLEPC_EXTERN PetscErrorCode BVSetSizes(BV,PetscInt,PetscInt,PetscInt);
238: SLEPC_EXTERN PetscErrorCode BVSetSizesFromVec(BV,Vec,PetscInt);
239: SLEPC_EXTERN PetscErrorCode BVSetLeadingDimension(BV,PetscInt);
240: SLEPC_EXTERN PetscErrorCode BVGetLeadingDimension(BV,PetscInt*);
241: SLEPC_EXTERN PetscErrorCode BVGetSizes(BV,PetscInt*,PetscInt*,PetscInt*);
242: SLEPC_EXTERN PetscErrorCode BVResize(BV,PetscInt,PetscBool);
243: SLEPC_EXTERN PetscErrorCode BVSetFromOptions(BV);
244: SLEPC_EXTERN PetscErrorCode BVView(BV,PetscViewer);
245: SLEPC_EXTERN PetscErrorCode BVViewFromOptions(BV,PetscObject,const char[]);

247: SLEPC_EXTERN PetscErrorCode BVGetColumn(BV,PetscInt,Vec*);
248: SLEPC_EXTERN PetscErrorCode BVRestoreColumn(BV,PetscInt,Vec*);
249: SLEPC_EXTERN PetscErrorCode BVGetSplit(BV,BV*,BV*);
250: SLEPC_EXTERN PetscErrorCode BVRestoreSplit(BV,BV*,BV*);
251: SLEPC_EXTERN PetscErrorCode BVGetSplitRows(BV,IS,IS,BV*,BV*);
252: SLEPC_EXTERN PetscErrorCode BVRestoreSplitRows(BV,IS,IS,BV*,BV*);
253: SLEPC_EXTERN PetscErrorCode BVGetArray(BV,PetscScalar*[]);
254: SLEPC_EXTERN PetscErrorCode BVRestoreArray(BV,PetscScalar*[]);
255: SLEPC_EXTERN PetscErrorCode BVGetArrayRead(BV,const PetscScalar*[]);
256: SLEPC_EXTERN PetscErrorCode BVRestoreArrayRead(BV,const PetscScalar*[]);
257: SLEPC_EXTERN PetscErrorCode BVCreateVec(BV,Vec*);
258: SLEPC_EXTERN PetscErrorCode BVCreateVecEmpty(BV,Vec*);
259: SLEPC_EXTERN PetscErrorCode BVSetVecType(BV,VecType);
260: SLEPC_EXTERN PetscErrorCode BVGetVecType(BV,VecType*);
261: SLEPC_EXTERN PetscErrorCode BVSetActiveColumns(BV,PetscInt,PetscInt);
262: SLEPC_EXTERN PetscErrorCode BVGetActiveColumns(BV,PetscInt*,PetscInt*);
263: SLEPC_EXTERN PetscErrorCode BVInsertVec(BV,PetscInt,Vec);
264: SLEPC_EXTERN PetscErrorCode BVInsertVecs(BV,PetscInt,PetscInt*,Vec[],PetscBool);
265: SLEPC_EXTERN PetscErrorCode BVInsertConstraints(BV,PetscInt*,Vec[]);
266: SLEPC_EXTERN PetscErrorCode BVSetNumConstraints(BV,PetscInt);
267: SLEPC_EXTERN PetscErrorCode BVGetNumConstraints(BV,PetscInt*);
268: SLEPC_EXTERN PetscErrorCode BVSetDefiniteTolerance(BV,PetscReal);
269: SLEPC_EXTERN PetscErrorCode BVGetDefiniteTolerance(BV,PetscReal*);
270: SLEPC_EXTERN PetscErrorCode BVDuplicate(BV,BV*);
271: SLEPC_EXTERN PetscErrorCode BVDuplicateResize(BV,PetscInt,BV*);
272: SLEPC_EXTERN PetscErrorCode BVCopy(BV,BV);
273: SLEPC_EXTERN PetscErrorCode BVCopyVec(BV,PetscInt,Vec);
274: SLEPC_EXTERN PetscErrorCode BVCopyColumn(BV,PetscInt,PetscInt);
275: SLEPC_EXTERN PetscErrorCode BVSetMatrix(BV,Mat,PetscBool);
276: SLEPC_EXTERN PetscErrorCode BVGetMatrix(BV,Mat*,PetscBool*);
277: SLEPC_EXTERN PetscErrorCode BVApplyMatrix(BV,Vec,Vec);
278: SLEPC_EXTERN PetscErrorCode BVApplyMatrixBV(BV,BV);
279: SLEPC_EXTERN PetscErrorCode BVGetCachedBV(BV,BV*);
280: SLEPC_EXTERN PetscErrorCode BVSetSignature(BV,Vec);
281: SLEPC_EXTERN PetscErrorCode BVGetSignature(BV,Vec);
282: SLEPC_EXTERN PetscErrorCode BVSetBufferVec(BV,Vec);
283: SLEPC_EXTERN PetscErrorCode BVGetBufferVec(BV,Vec*);

285: SLEPC_EXTERN PetscErrorCode BVMult(BV,PetscScalar,PetscScalar,BV,Mat);
286: SLEPC_EXTERN PetscErrorCode BVMultVec(BV,PetscScalar,PetscScalar,Vec,PetscScalar[]);
287: SLEPC_EXTERN PetscErrorCode BVMultColumn(BV,PetscScalar,PetscScalar,PetscInt,PetscScalar[]);
288: SLEPC_EXTERN PetscErrorCode BVMultInPlace(BV,Mat,PetscInt,PetscInt);
289: SLEPC_EXTERN PetscErrorCode BVMultInPlaceHermitianTranspose(BV,Mat,PetscInt,PetscInt);
290: PETSC_DEPRECATED_FUNCTION(3, 16, 0, "BVMultInPlaceHermitianTranspose()", ) static inline PetscErrorCode BVMultInPlaceTranspose(BV bv,Mat A,PetscInt s,PetscInt e) {return BVMultInPlaceHermitianTranspose(bv,A,s,e);}
291: SLEPC_EXTERN PetscErrorCode BVMatMult(BV,Mat,BV);
292: SLEPC_EXTERN PetscErrorCode BVMatMultTranspose(BV,Mat,BV);
293: SLEPC_EXTERN PetscErrorCode BVMatMultHermitianTranspose(BV,Mat,BV);
294: SLEPC_EXTERN PetscErrorCode BVMatMultColumn(BV,Mat,PetscInt);
295: SLEPC_EXTERN PetscErrorCode BVMatMultTransposeColumn(BV,Mat,PetscInt);
296: SLEPC_EXTERN PetscErrorCode BVMatMultHermitianTransposeColumn(BV,Mat,PetscInt);
297: SLEPC_EXTERN PetscErrorCode BVMatProject(BV,Mat,BV,Mat);
298: SLEPC_EXTERN PetscErrorCode BVMatArnoldi(BV,Mat,Mat,PetscInt,PetscInt*,PetscReal*,PetscBool*);
299: SLEPC_EXTERN PetscErrorCode BVMatLanczos(BV,Mat,Mat,PetscInt,PetscInt*,PetscReal*,PetscBool*);

301: SLEPC_EXTERN PetscErrorCode BVDot(BV,BV,Mat);
302: SLEPC_EXTERN PetscErrorCode BVDotVec(BV,Vec,PetscScalar[]);
303: SLEPC_EXTERN PetscErrorCode BVDotVecBegin(BV,Vec,PetscScalar[]);
304: SLEPC_EXTERN PetscErrorCode BVDotVecEnd(BV,Vec,PetscScalar[]);
305: SLEPC_EXTERN PetscErrorCode BVDotColumn(BV,PetscInt,PetscScalar[]);
306: SLEPC_EXTERN PetscErrorCode BVDotColumnBegin(BV,PetscInt,PetscScalar[]);
307: SLEPC_EXTERN PetscErrorCode BVDotColumnEnd(BV,PetscInt,PetscScalar[]);
308: SLEPC_EXTERN PetscErrorCode BVScale(BV,PetscScalar);
309: SLEPC_EXTERN PetscErrorCode BVScaleColumn(BV,PetscInt,PetscScalar);
310: SLEPC_EXTERN PetscErrorCode BVNorm(BV,NormType,PetscReal*);
311: SLEPC_EXTERN PetscErrorCode BVNormVec(BV,Vec,NormType,PetscReal*);
312: SLEPC_EXTERN PetscErrorCode BVNormVecBegin(BV,Vec,NormType,PetscReal*);
313: SLEPC_EXTERN PetscErrorCode BVNormVecEnd(BV,Vec,NormType,PetscReal*);
314: SLEPC_EXTERN PetscErrorCode BVNormColumn(BV,PetscInt,NormType,PetscReal*);
315: SLEPC_EXTERN PetscErrorCode BVNormColumnBegin(BV,PetscInt,NormType,PetscReal*);
316: SLEPC_EXTERN PetscErrorCode BVNormColumnEnd(BV,PetscInt,NormType,PetscReal*);
317: SLEPC_EXTERN PetscErrorCode BVNormalize(BV,PetscScalar[]);
318: SLEPC_EXTERN PetscErrorCode BVSetRandom(BV);
319: SLEPC_EXTERN PetscErrorCode BVSetRandomNormal(BV);
320: SLEPC_EXTERN PetscErrorCode BVSetRandomSign(BV);
321: SLEPC_EXTERN PetscErrorCode BVSetRandomColumn(BV,PetscInt);
322: SLEPC_EXTERN PetscErrorCode BVSetRandomCond(BV,PetscReal);
323: SLEPC_EXTERN PetscErrorCode BVSetRandomContext(BV,PetscRandom);
324: SLEPC_EXTERN PetscErrorCode BVGetRandomContext(BV,PetscRandom*);

326: SLEPC_EXTERN PetscErrorCode BVSetOrthogonalization(BV,BVOrthogType,BVOrthogRefineType,PetscReal,BVOrthogBlockType);
327: SLEPC_EXTERN PetscErrorCode BVGetOrthogonalization(BV,BVOrthogType*,BVOrthogRefineType*,PetscReal*,BVOrthogBlockType*);
328: SLEPC_EXTERN PetscErrorCode BVOrthogonalize(BV,Mat);
329: SLEPC_EXTERN PetscErrorCode BVOrthogonalizeVec(BV,Vec,PetscScalar[],PetscReal*,PetscBool*);
330: SLEPC_EXTERN PetscErrorCode BVOrthogonalizeColumn(BV,PetscInt,PetscScalar[],PetscReal*,PetscBool*);
331: SLEPC_EXTERN PetscErrorCode BVOrthonormalizeColumn(BV,PetscInt,PetscBool,PetscReal*,PetscBool*);
332: SLEPC_EXTERN PetscErrorCode BVOrthogonalizeSomeColumn(BV,PetscInt,PetscBool*,PetscScalar[],PetscReal*,PetscBool*);
333: SLEPC_EXTERN PetscErrorCode BVBiorthogonalizeColumn(BV,BV,PetscInt);
334: SLEPC_EXTERN PetscErrorCode BVBiorthonormalizeColumn(BV,BV,PetscInt,PetscReal*);
335: SLEPC_EXTERN PetscErrorCode BVSetMatMultMethod(BV,BVMatMultType);
336: SLEPC_EXTERN PetscErrorCode BVGetMatMultMethod(BV,BVMatMultType*);

338: SLEPC_EXTERN PetscErrorCode BVCreateFromMat(Mat,BV*);
339: SLEPC_EXTERN PetscErrorCode BVCreateMat(BV,Mat*);
340: SLEPC_EXTERN PetscErrorCode BVGetMat(BV,Mat*);
341: SLEPC_EXTERN PetscErrorCode BVRestoreMat(BV,Mat*);

343: SLEPC_EXTERN PetscErrorCode BVScatter(BV,BV,VecScatter,Vec);
344: SLEPC_EXTERN PetscErrorCode BVSumQuadrature(BV,BV,PetscInt,PetscInt,PetscInt,PetscScalar[],PetscScalar[],VecScatter,PetscSubcomm,PetscInt,PetscBool);
345: SLEPC_EXTERN PetscErrorCode BVDotQuadrature(BV,BV,PetscScalar[],PetscInt,PetscInt,PetscInt,PetscScalar[],PetscScalar[],PetscSubcomm,PetscInt,PetscBool);
346: SLEPC_EXTERN PetscErrorCode BVTraceQuadrature(BV,BV,PetscInt,PetscInt,PetscScalar[],VecScatter,PetscSubcomm,PetscInt,PetscBool,PetscReal*);

348: /*E
349:    BVSVDMethod - Different methods for computing the SVD of a `BV`.

351:    Values:
352: +  `BV_SVD_METHOD_REFINE` - based on the SVD of the cross product matrix $S^*S$, with refinement
353: .  `BV_SVD_METHOD_QR`     - based on the SVD of the triangular factor of qr(S)
354: -  `BV_SVD_METHOD_QR_CAA` - variant of QR intended for use in communication-avoiding Arnoldi

356:    Level: developer

358: .seealso: [](sec:bv), `BVSVDAndRank()`
359: E*/
360: typedef enum { BV_SVD_METHOD_REFINE,
361:                BV_SVD_METHOD_QR,
362:                BV_SVD_METHOD_QR_CAA } BVSVDMethod;
363: SLEPC_EXTERN const char *BVSVDMethods[];

365: /*MC
366:    BV_SVD_METHOD_REFINE - Compute the SVD of the `BV` via the cross product matrix $S^*S$,
367:    with iterative refinement.

369:    Level: advanced

371: .seealso: [](sec:bv), `BVSVDMethod`, `BVSVDAndRank()`, `BV_SVD_METHOD_QR`, `BV_SVD_METHOD_QR_CAA`
372: M*/

374: /*MC
375:    BV_SVD_METHOD_QR - Compute the SVD of the `BV` via the SVD of the triangular factor of
376:    the QR factorization of $S$.

378:    Level: advanced

380: .seealso: [](sec:bv), `BVSVDMethod`, `BVSVDAndRank()`, `BV_SVD_METHOD_REFINE`, `BV_SVD_METHOD_QR_CAA`
381: M*/

383: /*MC
384:    BV_SVD_METHOD_QR_CAA - Compute the SVD of the `BV` employing a variant of
385:    `BV_SVD_METHOD_QR` intended for use in communication-avoiding Arnoldi.

387:    Level: advanced

389: .seealso: [](sec:bv), `BVSVDMethod`, `BVSVDAndRank()`, `BV_SVD_METHOD_REFINE`, `BV_SVD_METHOD_QR`
390: M*/

392: SLEPC_EXTERN PetscErrorCode BVSVDAndRank(BV,PetscInt,PetscInt,PetscReal,BVSVDMethod,PetscScalar[],PetscReal[],PetscInt*);
393: SLEPC_EXTERN PetscErrorCode BVCISSResizeBases(BV,BV,BV,PetscInt,PetscInt,PetscInt,PetscInt);

395: SLEPC_EXTERN PetscErrorCode BVCreateTensor(BV,PetscInt,BV*);
396: SLEPC_EXTERN PetscErrorCode BVTensorBuildFirstColumn(BV,PetscInt);
397: SLEPC_EXTERN PetscErrorCode BVTensorCompress(BV,PetscInt);
398: SLEPC_EXTERN PetscErrorCode BVTensorGetDegree(BV,PetscInt*);
399: SLEPC_EXTERN PetscErrorCode BVTensorGetFactors(BV,BV*,Mat*);
400: SLEPC_EXTERN PetscErrorCode BVTensorRestoreFactors(BV,BV*,Mat*);

402: SLEPC_EXTERN PetscErrorCode BVSetOptionsPrefix(BV,const char[]);
403: SLEPC_EXTERN PetscErrorCode BVAppendOptionsPrefix(BV,const char[]);
404: SLEPC_EXTERN PetscErrorCode BVGetOptionsPrefix(BV,const char*[]);

406: SLEPC_EXTERN PetscFunctionList BVList;
407: SLEPC_EXTERN PetscErrorCode BVRegister(const char[],PetscErrorCode(*)(BV));