Actual source code: svec.c
slepc-3.21.1 2024-04-26
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: BV implemented as a single Vec
12: */
14: #include <slepc/private/bvimpl.h>
15: #include "svec.h"
17: static PetscErrorCode BVMult_Svec(BV Y,PetscScalar alpha,PetscScalar beta,BV X,Mat Q)
18: {
19: BV_SVEC *y = (BV_SVEC*)Y->data,*x = (BV_SVEC*)X->data;
20: const PetscScalar *px,*q;
21: PetscScalar *py;
22: PetscInt ldq;
24: PetscFunctionBegin;
25: PetscCall(VecGetArrayRead(x->v,&px));
26: PetscCall(VecGetArray(y->v,&py));
27: if (Q) {
28: PetscCall(MatDenseGetLDA(Q,&ldq));
29: PetscCall(MatDenseGetArrayRead(Q,&q));
30: PetscCall(BVMult_BLAS_Private(Y,Y->n,Y->k-Y->l,X->k-X->l,alpha,px+(X->nc+X->l)*X->ld,X->ld,q+Y->l*ldq+X->l,ldq,beta,py+(Y->nc+Y->l)*Y->ld,Y->ld));
31: PetscCall(MatDenseRestoreArrayRead(Q,&q));
32: } else PetscCall(BVAXPY_BLAS_Private(Y,Y->n,Y->k-Y->l,alpha,px+(X->nc+X->l)*X->ld,X->ld,beta,py+(Y->nc+Y->l)*Y->ld,Y->ld));
33: PetscCall(VecRestoreArrayRead(x->v,&px));
34: PetscCall(VecRestoreArray(y->v,&py));
35: PetscFunctionReturn(PETSC_SUCCESS);
36: }
38: static PetscErrorCode BVMultVec_Svec(BV X,PetscScalar alpha,PetscScalar beta,Vec y,PetscScalar *q)
39: {
40: BV_SVEC *x = (BV_SVEC*)X->data;
41: PetscScalar *px,*py,*qq=q;
43: PetscFunctionBegin;
44: PetscCall(VecGetArray(x->v,&px));
45: PetscCall(VecGetArray(y,&py));
46: if (!q) PetscCall(VecGetArray(X->buffer,&qq));
47: PetscCall(BVMultVec_BLAS_Private(X,X->n,X->k-X->l,alpha,px+(X->nc+X->l)*X->ld,X->ld,qq,beta,py));
48: if (!q) PetscCall(VecRestoreArray(X->buffer,&qq));
49: PetscCall(VecRestoreArray(x->v,&px));
50: PetscCall(VecRestoreArray(y,&py));
51: PetscFunctionReturn(PETSC_SUCCESS);
52: }
54: static PetscErrorCode BVMultInPlace_Svec(BV V,Mat Q,PetscInt s,PetscInt e)
55: {
56: BV_SVEC *ctx = (BV_SVEC*)V->data;
57: PetscScalar *pv;
58: const PetscScalar *q;
59: PetscInt ldq;
61: PetscFunctionBegin;
62: if (s>=e || !V->n) PetscFunctionReturn(PETSC_SUCCESS);
63: PetscCall(MatDenseGetLDA(Q,&ldq));
64: PetscCall(VecGetArray(ctx->v,&pv));
65: PetscCall(MatDenseGetArrayRead(Q,&q));
66: PetscCall(BVMultInPlace_BLAS_Private(V,V->n,V->k-V->l,s-V->l,e-V->l,pv+(V->nc+V->l)*V->ld,V->ld,q+V->l*ldq+V->l,ldq,PETSC_FALSE));
67: PetscCall(MatDenseRestoreArrayRead(Q,&q));
68: PetscCall(VecRestoreArray(ctx->v,&pv));
69: PetscFunctionReturn(PETSC_SUCCESS);
70: }
72: static PetscErrorCode BVMultInPlaceHermitianTranspose_Svec(BV V,Mat Q,PetscInt s,PetscInt e)
73: {
74: BV_SVEC *ctx = (BV_SVEC*)V->data;
75: PetscScalar *pv;
76: const PetscScalar *q;
77: PetscInt ldq;
79: PetscFunctionBegin;
80: if (s>=e || !V->n) PetscFunctionReturn(PETSC_SUCCESS);
81: PetscCall(MatDenseGetLDA(Q,&ldq));
82: PetscCall(VecGetArray(ctx->v,&pv));
83: PetscCall(MatDenseGetArrayRead(Q,&q));
84: PetscCall(BVMultInPlace_BLAS_Private(V,V->n,V->k-V->l,s-V->l,e-V->l,pv+(V->nc+V->l)*V->ld,V->ld,q+V->l*ldq+V->l,ldq,PETSC_TRUE));
85: PetscCall(MatDenseRestoreArrayRead(Q,&q));
86: PetscCall(VecRestoreArray(ctx->v,&pv));
87: PetscFunctionReturn(PETSC_SUCCESS);
88: }
90: static PetscErrorCode BVDot_Svec(BV X,BV Y,Mat M)
91: {
92: BV_SVEC *x = (BV_SVEC*)X->data,*y = (BV_SVEC*)Y->data;
93: const PetscScalar *px,*py;
94: PetscScalar *m;
95: PetscInt ldm;
97: PetscFunctionBegin;
98: PetscCall(MatDenseGetLDA(M,&ldm));
99: PetscCall(VecGetArrayRead(x->v,&px));
100: PetscCall(VecGetArrayRead(y->v,&py));
101: PetscCall(MatDenseGetArray(M,&m));
102: PetscCall(BVDot_BLAS_Private(X,Y->k-Y->l,X->k-X->l,X->n,py+(Y->nc+Y->l)*Y->ld,Y->ld,px+(X->nc+X->l)*X->ld,X->ld,m+X->l*ldm+Y->l,ldm,x->mpi));
103: PetscCall(MatDenseRestoreArray(M,&m));
104: PetscCall(VecRestoreArrayRead(x->v,&px));
105: PetscCall(VecRestoreArrayRead(y->v,&py));
106: PetscFunctionReturn(PETSC_SUCCESS);
107: }
109: static PetscErrorCode BVDotVec_Svec(BV X,Vec y,PetscScalar *q)
110: {
111: BV_SVEC *x = (BV_SVEC*)X->data;
112: const PetscScalar *px,*py;
113: PetscScalar *qq=q;
114: Vec z = y;
116: PetscFunctionBegin;
117: if (PetscUnlikely(X->matrix)) {
118: PetscCall(BV_IPMatMult(X,y));
119: z = X->Bx;
120: }
121: PetscCall(VecGetArrayRead(x->v,&px));
122: PetscCall(VecGetArrayRead(z,&py));
123: if (!q) PetscCall(VecGetArray(X->buffer,&qq));
124: PetscCall(BVDotVec_BLAS_Private(X,X->n,X->k-X->l,px+(X->nc+X->l)*X->ld,X->ld,py,qq,x->mpi));
125: if (!q) PetscCall(VecRestoreArray(X->buffer,&qq));
126: PetscCall(VecRestoreArrayRead(z,&py));
127: PetscCall(VecRestoreArrayRead(x->v,&px));
128: PetscFunctionReturn(PETSC_SUCCESS);
129: }
131: static PetscErrorCode BVDotVec_Local_Svec(BV X,Vec y,PetscScalar *m)
132: {
133: BV_SVEC *x = (BV_SVEC*)X->data;
134: PetscScalar *px,*py;
135: Vec z = y;
137: PetscFunctionBegin;
138: if (PetscUnlikely(X->matrix)) {
139: PetscCall(BV_IPMatMult(X,y));
140: z = X->Bx;
141: }
142: PetscCall(VecGetArray(x->v,&px));
143: PetscCall(VecGetArray(z,&py));
144: PetscCall(BVDotVec_BLAS_Private(X,X->n,X->k-X->l,px+(X->nc+X->l)*X->ld,X->ld,py,m,PETSC_FALSE));
145: PetscCall(VecRestoreArray(z,&py));
146: PetscCall(VecRestoreArray(x->v,&px));
147: PetscFunctionReturn(PETSC_SUCCESS);
148: }
150: static PetscErrorCode BVScale_Svec(BV bv,PetscInt j,PetscScalar alpha)
151: {
152: BV_SVEC *ctx = (BV_SVEC*)bv->data;
153: PetscScalar *array;
155: PetscFunctionBegin;
156: if (!bv->n) PetscFunctionReturn(PETSC_SUCCESS);
157: PetscCall(VecGetArray(ctx->v,&array));
158: if (PetscUnlikely(j<0)) PetscCall(BVScale_BLAS_Private(bv,(bv->k-bv->l)*bv->ld,array+(bv->nc+bv->l)*bv->ld,alpha));
159: else PetscCall(BVScale_BLAS_Private(bv,bv->n,array+(bv->nc+j)*bv->ld,alpha));
160: PetscCall(VecRestoreArray(ctx->v,&array));
161: PetscFunctionReturn(PETSC_SUCCESS);
162: }
164: static PetscErrorCode BVNorm_Svec(BV bv,PetscInt j,NormType type,PetscReal *val)
165: {
166: BV_SVEC *ctx = (BV_SVEC*)bv->data;
167: const PetscScalar *array;
169: PetscFunctionBegin;
170: PetscCall(VecGetArrayRead(ctx->v,&array));
171: if (PetscUnlikely(j<0)) PetscCall(BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->ld,bv->ld,type,val,ctx->mpi));
172: else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->ld,bv->ld,type,val,ctx->mpi));
173: PetscCall(VecRestoreArrayRead(ctx->v,&array));
174: PetscFunctionReturn(PETSC_SUCCESS);
175: }
177: static PetscErrorCode BVNorm_Local_Svec(BV bv,PetscInt j,NormType type,PetscReal *val)
178: {
179: BV_SVEC *ctx = (BV_SVEC*)bv->data;
180: const PetscScalar *array;
182: PetscFunctionBegin;
183: PetscCall(VecGetArrayRead(ctx->v,&array));
184: if (PetscUnlikely(j<0)) PetscCall(BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->ld,bv->ld,type,val,PETSC_FALSE));
185: else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->ld,bv->ld,type,val,PETSC_FALSE));
186: PetscCall(VecRestoreArrayRead(ctx->v,&array));
187: PetscFunctionReturn(PETSC_SUCCESS);
188: }
190: static PetscErrorCode BVNormalize_Svec(BV bv,PetscScalar *eigi)
191: {
192: BV_SVEC *ctx = (BV_SVEC*)bv->data;
193: PetscScalar *array,*wi=NULL;
195: PetscFunctionBegin;
196: PetscCall(VecGetArray(ctx->v,&array));
197: if (eigi) wi = eigi+bv->l;
198: PetscCall(BVNormalize_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->ld,bv->ld,wi,ctx->mpi));
199: PetscCall(VecRestoreArray(ctx->v,&array));
200: PetscFunctionReturn(PETSC_SUCCESS);
201: }
203: static PetscErrorCode BVMatMult_Svec(BV V,Mat A,BV W)
204: {
205: PetscInt j;
206: Mat Vmat,Wmat;
207: Vec vv,ww;
209: PetscFunctionBegin;
210: if (V->vmm) {
211: PetscCall(BVGetMat(V,&Vmat));
212: PetscCall(BVGetMat(W,&Wmat));
213: PetscCall(MatProductCreateWithMat(A,Vmat,NULL,Wmat));
214: PetscCall(MatProductSetType(Wmat,MATPRODUCT_AB));
215: PetscCall(MatProductSetFromOptions(Wmat));
216: PetscCall(MatProductSymbolic(Wmat));
217: PetscCall(MatProductNumeric(Wmat));
218: PetscCall(MatProductClear(Wmat));
219: PetscCall(BVRestoreMat(V,&Vmat));
220: PetscCall(BVRestoreMat(W,&Wmat));
221: } else {
222: for (j=0;j<V->k-V->l;j++) {
223: PetscCall(BVGetColumn(V,V->l+j,&vv));
224: PetscCall(BVGetColumn(W,W->l+j,&ww));
225: PetscCall(MatMult(A,vv,ww));
226: PetscCall(BVRestoreColumn(V,V->l+j,&vv));
227: PetscCall(BVRestoreColumn(W,W->l+j,&ww));
228: }
229: }
230: PetscFunctionReturn(PETSC_SUCCESS);
231: }
233: static PetscErrorCode BVCopy_Svec(BV V,BV W)
234: {
235: BV_SVEC *v = (BV_SVEC*)V->data,*w = (BV_SVEC*)W->data;
236: const PetscScalar *pv;
237: PetscScalar *pw;
238: PetscInt j;
240: PetscFunctionBegin;
241: PetscCall(VecGetArrayRead(v->v,&pv));
242: PetscCall(VecGetArray(w->v,&pw));
243: for (j=0;j<V->k-V->l;j++) PetscCall(PetscArraycpy(pw+(W->nc+W->l+j)*W->ld,pv+(V->nc+V->l+j)*V->ld,V->n));
244: PetscCall(VecRestoreArrayRead(v->v,&pv));
245: PetscCall(VecRestoreArray(w->v,&pw));
246: PetscFunctionReturn(PETSC_SUCCESS);
247: }
249: static PetscErrorCode BVCopyColumn_Svec(BV V,PetscInt j,PetscInt i)
250: {
251: BV_SVEC *v = (BV_SVEC*)V->data;
252: PetscScalar *pv;
254: PetscFunctionBegin;
255: PetscCall(VecGetArray(v->v,&pv));
256: PetscCall(PetscArraycpy(pv+(V->nc+i)*V->ld,pv+(V->nc+j)*V->ld,V->n));
257: PetscCall(VecRestoreArray(v->v,&pv));
258: PetscFunctionReturn(PETSC_SUCCESS);
259: }
261: static PetscErrorCode BVResize_Svec(BV bv,PetscInt m,PetscBool copy)
262: {
263: BV_SVEC *ctx = (BV_SVEC*)bv->data;
264: PetscScalar *pnew;
265: const PetscScalar *pv;
266: PetscInt bs;
267: Vec vnew;
268: char str[50];
270: PetscFunctionBegin;
271: PetscCall(PetscLayoutGetBlockSize(bv->map,&bs));
272: PetscCall(VecCreate(PetscObjectComm((PetscObject)bv),&vnew));
273: PetscCall(VecSetType(vnew,bv->vtype));
274: PetscCall(VecSetSizes(vnew,m*bv->ld,PETSC_DECIDE));
275: PetscCall(VecSetBlockSize(vnew,bs));
276: if (((PetscObject)bv)->name) {
277: PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
278: PetscCall(PetscObjectSetName((PetscObject)vnew,str));
279: }
280: if (copy) {
281: PetscCall(VecGetArrayRead(ctx->v,&pv));
282: PetscCall(VecGetArray(vnew,&pnew));
283: PetscCall(PetscArraycpy(pnew,pv,PetscMin(m,bv->m)*bv->ld));
284: PetscCall(VecRestoreArrayRead(ctx->v,&pv));
285: PetscCall(VecRestoreArray(vnew,&pnew));
286: }
287: PetscCall(VecDestroy(&ctx->v));
288: ctx->v = vnew;
289: PetscFunctionReturn(PETSC_SUCCESS);
290: }
292: static PetscErrorCode BVGetColumn_Svec(BV bv,PetscInt j,Vec *v)
293: {
294: BV_SVEC *ctx = (BV_SVEC*)bv->data;
295: PetscScalar *pv;
296: PetscInt l;
298: PetscFunctionBegin;
299: l = BVAvailableVec;
300: PetscCall(VecGetArray(ctx->v,&pv));
301: PetscCall(VecPlaceArray(bv->cv[l],pv+(bv->nc+j)*bv->ld));
302: PetscFunctionReturn(PETSC_SUCCESS);
303: }
305: static PetscErrorCode BVRestoreColumn_Svec(BV bv,PetscInt j,Vec *v)
306: {
307: BV_SVEC *ctx = (BV_SVEC*)bv->data;
308: PetscInt l;
310: PetscFunctionBegin;
311: l = (j==bv->ci[0])? 0: 1;
312: PetscCall(VecResetArray(bv->cv[l]));
313: PetscCall(VecRestoreArray(ctx->v,NULL));
314: PetscFunctionReturn(PETSC_SUCCESS);
315: }
317: static PetscErrorCode BVGetArray_Svec(BV bv,PetscScalar **a)
318: {
319: BV_SVEC *ctx = (BV_SVEC*)bv->data;
321: PetscFunctionBegin;
322: PetscCall(VecGetArray(ctx->v,a));
323: PetscFunctionReturn(PETSC_SUCCESS);
324: }
326: static PetscErrorCode BVRestoreArray_Svec(BV bv,PetscScalar **a)
327: {
328: BV_SVEC *ctx = (BV_SVEC*)bv->data;
330: PetscFunctionBegin;
331: PetscCall(VecRestoreArray(ctx->v,a));
332: PetscFunctionReturn(PETSC_SUCCESS);
333: }
335: static PetscErrorCode BVGetArrayRead_Svec(BV bv,const PetscScalar **a)
336: {
337: BV_SVEC *ctx = (BV_SVEC*)bv->data;
339: PetscFunctionBegin;
340: PetscCall(VecGetArrayRead(ctx->v,a));
341: PetscFunctionReturn(PETSC_SUCCESS);
342: }
344: static PetscErrorCode BVRestoreArrayRead_Svec(BV bv,const PetscScalar **a)
345: {
346: BV_SVEC *ctx = (BV_SVEC*)bv->data;
348: PetscFunctionBegin;
349: PetscCall(VecRestoreArrayRead(ctx->v,a));
350: PetscFunctionReturn(PETSC_SUCCESS);
351: }
353: static PetscErrorCode BVView_Svec(BV bv,PetscViewer viewer)
354: {
355: PetscInt j;
356: Vec v;
357: PetscViewerFormat format;
358: PetscBool isascii,ismatlab=PETSC_FALSE;
359: const char *bvname,*name;
361: PetscFunctionBegin;
362: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
363: if (isascii) {
364: PetscCall(PetscViewerGetFormat(viewer,&format));
365: if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) PetscFunctionReturn(PETSC_SUCCESS);
366: if (format == PETSC_VIEWER_ASCII_MATLAB) ismatlab = PETSC_TRUE;
367: }
368: if (ismatlab) {
369: PetscCall(PetscObjectGetName((PetscObject)bv,&bvname));
370: PetscCall(PetscViewerASCIIPrintf(viewer,"%s=[];\n",bvname));
371: }
372: for (j=0;j<bv->m;j++) {
373: PetscCall(BVGetColumn(bv,j,&v));
374: PetscCall(VecView(v,viewer));
375: if (ismatlab) {
376: PetscCall(PetscObjectGetName((PetscObject)v,&name));
377: PetscCall(PetscViewerASCIIPrintf(viewer,"%s=[%s,%s];clear %s\n",bvname,bvname,name,name));
378: }
379: PetscCall(BVRestoreColumn(bv,j,&v));
380: }
381: PetscFunctionReturn(PETSC_SUCCESS);
382: }
384: static PetscErrorCode BVDestroy_Svec(BV bv)
385: {
386: BV_SVEC *ctx = (BV_SVEC*)bv->data;
388: PetscFunctionBegin;
389: PetscCall(VecDestroy(&ctx->v));
390: PetscCall(VecDestroy(&bv->cv[0]));
391: PetscCall(VecDestroy(&bv->cv[1]));
392: PetscCall(PetscFree(bv->data));
393: bv->cuda = PETSC_FALSE;
394: PetscFunctionReturn(PETSC_SUCCESS);
395: }
397: SLEPC_EXTERN PetscErrorCode BVCreate_Svec(BV bv)
398: {
399: BV_SVEC *ctx;
400: PetscInt nloc,N,bs,tglobal=0,tlocal,lsplit,j,lda;
401: PetscBool seq,isdense;
402: PetscScalar *vv;
403: const PetscScalar *aa,*array,*ptr;
404: char str[50];
405: BV parent;
406: Vec vpar;
407: MatType mtype;
409: PetscFunctionBegin;
410: PetscCall(PetscNew(&ctx));
411: bv->data = (void*)ctx;
413: PetscCall(PetscStrcmpAny(bv->vtype,&bv->cuda,VECSEQCUDA,VECMPICUDA,""));
414: PetscCall(PetscStrcmpAny(bv->vtype,&ctx->mpi,VECMPI,VECMPICUDA,""));
416: PetscCall(PetscStrcmp(bv->vtype,VECSEQ,&seq));
417: PetscCheck(seq || ctx->mpi || bv->cuda,PetscObjectComm((PetscObject)bv),PETSC_ERR_SUP,"BVSVEC does not support the requested vector type: %s",bv->vtype);
419: PetscCall(PetscLayoutGetLocalSize(bv->map,&nloc));
420: PetscCall(PetscLayoutGetSize(bv->map,&N));
421: PetscCall(PetscLayoutGetBlockSize(bv->map,&bs));
422: PetscCall(BV_SetDefaultLD(bv,nloc));
423: tlocal = bv->m*bv->ld;
424: PetscCall(PetscIntMultError(bv->m,N,&tglobal)); /* just to check integer overflow */
426: if (PetscUnlikely(bv->issplit)) {
427: /* split BV: create Vec sharing the memory of the parent BV */
428: parent = bv->splitparent;
429: lsplit = parent->lsplit;
430: vpar = ((BV_SVEC*)parent->data)->v;
431: if (bv->cuda) {
432: #if defined(PETSC_HAVE_CUDA)
433: PetscCall(VecCUDAGetArrayRead(vpar,&array));
434: ptr = (bv->issplit==1)? array: array+lsplit*bv->ld;
435: PetscCall(VecCUDARestoreArrayRead(vpar,&array));
436: if (ctx->mpi) PetscCall(VecCreateMPICUDAWithArray(PetscObjectComm((PetscObject)bv),bs,tlocal,PETSC_DECIDE,NULL,&ctx->v));
437: else PetscCall(VecCreateSeqCUDAWithArray(PetscObjectComm((PetscObject)bv),bs,tlocal,NULL,&ctx->v));
438: PetscCall(VecCUDAPlaceArray(ctx->v,ptr));
439: #endif
440: } else {
441: PetscCall(VecGetArrayRead(vpar,&array));
442: ptr = (bv->issplit==1)? array: array+lsplit*bv->ld;
443: PetscCall(VecRestoreArrayRead(vpar,&array));
444: if (ctx->mpi) PetscCall(VecCreateMPIWithArray(PetscObjectComm((PetscObject)bv),bs,tlocal,PETSC_DECIDE,NULL,&ctx->v));
445: else PetscCall(VecCreateSeqWithArray(PetscObjectComm((PetscObject)bv),bs,tlocal,NULL,&ctx->v));
446: PetscCall(VecPlaceArray(ctx->v,ptr));
447: }
448: } else {
449: /* regular BV: create Vec to store the BV entries */
450: PetscCall(VecCreate(PetscObjectComm((PetscObject)bv),&ctx->v));
451: PetscCall(VecSetType(ctx->v,bv->vtype));
452: PetscCall(VecSetSizes(ctx->v,tlocal,PETSC_DECIDE));
453: PetscCall(VecSetBlockSize(ctx->v,bs));
454: }
455: if (((PetscObject)bv)->name) {
456: PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
457: PetscCall(PetscObjectSetName((PetscObject)ctx->v,str));
458: }
460: if (PetscUnlikely(bv->Acreate)) {
461: PetscCall(MatGetType(bv->Acreate,&mtype));
462: PetscCall(PetscStrcmpAny(mtype,&isdense,MATSEQDENSE,MATMPIDENSE,""));
463: PetscCheck(isdense,PetscObjectComm((PetscObject)bv->Acreate),PETSC_ERR_SUP,"BVSVEC requires a dense matrix in BVCreateFromMat()");
464: PetscCall(MatDenseGetArrayRead(bv->Acreate,&aa));
465: PetscCall(MatDenseGetLDA(bv->Acreate,&lda));
466: PetscCall(VecGetArray(ctx->v,&vv));
467: for (j=0;j<bv->m;j++) PetscCall(PetscArraycpy(vv+j*bv->ld,aa+j*lda,bv->n));
468: PetscCall(VecRestoreArray(ctx->v,&vv));
469: PetscCall(MatDenseRestoreArrayRead(bv->Acreate,&aa));
470: PetscCall(MatDestroy(&bv->Acreate));
471: }
473: PetscCall(BVCreateVecEmpty(bv,&bv->cv[0]));
474: PetscCall(BVCreateVecEmpty(bv,&bv->cv[1]));
476: if (bv->cuda) {
477: #if defined(PETSC_HAVE_CUDA)
478: bv->ops->mult = BVMult_Svec_CUDA;
479: bv->ops->multvec = BVMultVec_Svec_CUDA;
480: bv->ops->multinplace = BVMultInPlace_Svec_CUDA;
481: bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Svec_CUDA;
482: bv->ops->dot = BVDot_Svec_CUDA;
483: bv->ops->dotvec = BVDotVec_Svec_CUDA;
484: bv->ops->dotvec_local = BVDotVec_Local_Svec_CUDA;
485: bv->ops->scale = BVScale_Svec_CUDA;
486: bv->ops->norm = BVNorm_Svec_CUDA;
487: bv->ops->norm_local = BVNorm_Local_Svec_CUDA;
488: bv->ops->normalize = BVNormalize_Svec_CUDA;
489: bv->ops->matmult = BVMatMult_Svec_CUDA;
490: bv->ops->copy = BVCopy_Svec_CUDA;
491: bv->ops->copycolumn = BVCopyColumn_Svec_CUDA;
492: bv->ops->resize = BVResize_Svec_CUDA;
493: bv->ops->getcolumn = BVGetColumn_Svec_CUDA;
494: bv->ops->restorecolumn = BVRestoreColumn_Svec_CUDA;
495: bv->ops->restoresplit = BVRestoreSplit_Svec_CUDA;
496: bv->ops->getmat = BVGetMat_Svec_CUDA;
497: bv->ops->restoremat = BVRestoreMat_Svec_CUDA;
498: #endif
499: } else {
500: bv->ops->mult = BVMult_Svec;
501: bv->ops->multvec = BVMultVec_Svec;
502: bv->ops->multinplace = BVMultInPlace_Svec;
503: bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Svec;
504: bv->ops->dot = BVDot_Svec;
505: bv->ops->dotvec = BVDotVec_Svec;
506: bv->ops->dotvec_local = BVDotVec_Local_Svec;
507: bv->ops->scale = BVScale_Svec;
508: bv->ops->norm = BVNorm_Svec;
509: bv->ops->norm_local = BVNorm_Local_Svec;
510: bv->ops->normalize = BVNormalize_Svec;
511: bv->ops->matmult = BVMatMult_Svec;
512: bv->ops->copy = BVCopy_Svec;
513: bv->ops->copycolumn = BVCopyColumn_Svec;
514: bv->ops->resize = BVResize_Svec;
515: bv->ops->getcolumn = BVGetColumn_Svec;
516: bv->ops->restorecolumn = BVRestoreColumn_Svec;
517: bv->ops->getmat = BVGetMat_Default;
518: bv->ops->restoremat = BVRestoreMat_Default;
519: }
520: bv->ops->getarray = BVGetArray_Svec;
521: bv->ops->restorearray = BVRestoreArray_Svec;
522: bv->ops->getarrayread = BVGetArrayRead_Svec;
523: bv->ops->restorearrayread = BVRestoreArrayRead_Svec;
524: bv->ops->destroy = BVDestroy_Svec;
525: if (!ctx->mpi) bv->ops->view = BVView_Svec;
526: PetscFunctionReturn(PETSC_SUCCESS);
527: }