Line data Source code
1 : /*
2 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3 : SLEPc - Scalable Library for Eigenvalue Problem Computations
4 : Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5 :
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 with a dense Mat
12 : */
13 :
14 : #include <slepc/private/bvimpl.h>
15 : #include "bvmat.h"
16 :
17 6562 : static PetscErrorCode BVMult_Mat(BV Y,PetscScalar alpha,PetscScalar beta,BV X,Mat Q)
18 : {
19 6562 : BV_MAT *y = (BV_MAT*)Y->data,*x = (BV_MAT*)X->data;
20 6562 : PetscScalar *py;
21 6562 : const PetscScalar *px,*q;
22 6562 : PetscInt ldq;
23 :
24 6562 : PetscFunctionBegin;
25 6562 : PetscCall(MatDenseGetArrayRead(x->A,&px));
26 6562 : PetscCall(MatDenseGetArray(y->A,&py));
27 6562 : if (Q) {
28 5144 : PetscCall(MatDenseGetLDA(Q,&ldq));
29 5144 : PetscCall(MatDenseGetArrayRead(Q,&q));
30 5144 : 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 5144 : PetscCall(MatDenseRestoreArrayRead(Q,&q));
32 1418 : } 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 6562 : PetscCall(MatDenseRestoreArrayRead(x->A,&px));
34 6562 : PetscCall(MatDenseRestoreArray(y->A,&py));
35 6562 : PetscFunctionReturn(PETSC_SUCCESS);
36 : }
37 :
38 722969 : static PetscErrorCode BVMultVec_Mat(BV X,PetscScalar alpha,PetscScalar beta,Vec y,PetscScalar *q)
39 : {
40 722969 : BV_MAT *x = (BV_MAT*)X->data;
41 722969 : PetscScalar *py,*qq=q;
42 722969 : const PetscScalar *px;
43 :
44 722969 : PetscFunctionBegin;
45 722969 : PetscCall(MatDenseGetArrayRead(x->A,&px));
46 722969 : PetscCall(VecGetArray(y,&py));
47 722969 : if (!q) PetscCall(VecGetArray(X->buffer,&qq));
48 722969 : PetscCall(BVMultVec_BLAS_Private(X,X->n,X->k-X->l,alpha,px+(X->nc+X->l)*X->ld,X->ld,qq,beta,py));
49 722969 : if (!q) PetscCall(VecRestoreArray(X->buffer,&qq));
50 722969 : PetscCall(MatDenseRestoreArrayRead(x->A,&px));
51 722969 : PetscCall(VecRestoreArray(y,&py));
52 722969 : PetscFunctionReturn(PETSC_SUCCESS);
53 : }
54 :
55 23738 : static PetscErrorCode BVMultInPlace_Mat(BV V,Mat Q,PetscInt s,PetscInt e)
56 : {
57 23738 : BV_MAT *ctx = (BV_MAT*)V->data;
58 23738 : PetscScalar *pv;
59 23738 : const PetscScalar *q;
60 23738 : PetscInt ldq;
61 :
62 23738 : PetscFunctionBegin;
63 23738 : if (s>=e || !V->n) PetscFunctionReturn(PETSC_SUCCESS);
64 21937 : PetscCall(MatDenseGetLDA(Q,&ldq));
65 21937 : PetscCall(MatDenseGetArray(ctx->A,&pv));
66 21937 : PetscCall(MatDenseGetArrayRead(Q,&q));
67 21937 : 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));
68 21937 : PetscCall(MatDenseRestoreArrayRead(Q,&q));
69 21937 : PetscCall(MatDenseRestoreArray(ctx->A,&pv));
70 21937 : PetscFunctionReturn(PETSC_SUCCESS);
71 : }
72 :
73 1 : static PetscErrorCode BVMultInPlaceHermitianTranspose_Mat(BV V,Mat Q,PetscInt s,PetscInt e)
74 : {
75 1 : BV_MAT *ctx = (BV_MAT*)V->data;
76 1 : PetscScalar *pv;
77 1 : const PetscScalar *q;
78 1 : PetscInt ldq;
79 :
80 1 : PetscFunctionBegin;
81 1 : if (s>=e || !V->n) PetscFunctionReturn(PETSC_SUCCESS);
82 1 : PetscCall(MatDenseGetLDA(Q,&ldq));
83 1 : PetscCall(MatDenseGetArray(ctx->A,&pv));
84 1 : PetscCall(MatDenseGetArrayRead(Q,&q));
85 1 : 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));
86 1 : PetscCall(MatDenseRestoreArrayRead(Q,&q));
87 1 : PetscCall(MatDenseRestoreArray(ctx->A,&pv));
88 1 : PetscFunctionReturn(PETSC_SUCCESS);
89 : }
90 :
91 29616 : static PetscErrorCode BVDot_Mat(BV X,BV Y,Mat M)
92 : {
93 29616 : BV_MAT *x = (BV_MAT*)X->data,*y = (BV_MAT*)Y->data;
94 29616 : PetscScalar *m;
95 29616 : const PetscScalar *px,*py;
96 29616 : PetscInt ldm;
97 :
98 29616 : PetscFunctionBegin;
99 29616 : PetscCall(MatDenseGetLDA(M,&ldm));
100 29616 : PetscCall(MatDenseGetArrayRead(x->A,&px));
101 29616 : PetscCall(MatDenseGetArrayRead(y->A,&py));
102 29616 : PetscCall(MatDenseGetArray(M,&m));
103 29616 : 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));
104 29616 : PetscCall(MatDenseRestoreArray(M,&m));
105 29616 : PetscCall(MatDenseRestoreArrayRead(x->A,&px));
106 29616 : PetscCall(MatDenseRestoreArrayRead(y->A,&py));
107 29616 : PetscFunctionReturn(PETSC_SUCCESS);
108 : }
109 :
110 424752 : static PetscErrorCode BVDotVec_Mat(BV X,Vec y,PetscScalar *q)
111 : {
112 424752 : BV_MAT *x = (BV_MAT*)X->data;
113 424752 : PetscScalar *qq=q;
114 424752 : const PetscScalar *px,*py;
115 424752 : Vec z = y;
116 :
117 424752 : PetscFunctionBegin;
118 424752 : if (PetscUnlikely(X->matrix)) {
119 47381 : PetscCall(BV_IPMatMult(X,y));
120 47381 : z = X->Bx;
121 : }
122 424752 : PetscCall(MatDenseGetArrayRead(x->A,&px));
123 424752 : PetscCall(VecGetArrayRead(z,&py));
124 424752 : if (!q) PetscCall(VecGetArray(X->buffer,&qq));
125 424752 : PetscCall(BVDotVec_BLAS_Private(X,X->n,X->k-X->l,px+(X->nc+X->l)*X->ld,X->ld,py,qq,x->mpi));
126 424752 : if (!q) PetscCall(VecRestoreArray(X->buffer,&qq));
127 424752 : PetscCall(VecRestoreArrayRead(z,&py));
128 424752 : PetscCall(MatDenseRestoreArrayRead(x->A,&px));
129 424752 : PetscFunctionReturn(PETSC_SUCCESS);
130 : }
131 :
132 2351 : static PetscErrorCode BVDotVec_Local_Mat(BV X,Vec y,PetscScalar *m)
133 : {
134 2351 : BV_MAT *x = (BV_MAT*)X->data;
135 2351 : const PetscScalar *px,*py;
136 2351 : Vec z = y;
137 :
138 2351 : PetscFunctionBegin;
139 2351 : if (PetscUnlikely(X->matrix)) {
140 0 : PetscCall(BV_IPMatMult(X,y));
141 0 : z = X->Bx;
142 : }
143 2351 : PetscCall(MatDenseGetArrayRead(x->A,&px));
144 2351 : PetscCall(VecGetArrayRead(z,&py));
145 2351 : PetscCall(BVDotVec_BLAS_Private(X,X->n,X->k-X->l,px+(X->nc+X->l)*X->ld,X->ld,py,m,PETSC_FALSE));
146 2351 : PetscCall(VecRestoreArrayRead(z,&py));
147 2351 : PetscCall(MatDenseRestoreArrayRead(x->A,&px));
148 2351 : PetscFunctionReturn(PETSC_SUCCESS);
149 : }
150 :
151 179251 : static PetscErrorCode BVScale_Mat(BV bv,PetscInt j,PetscScalar alpha)
152 : {
153 179251 : BV_MAT *ctx = (BV_MAT*)bv->data;
154 179251 : PetscScalar *array;
155 :
156 179251 : PetscFunctionBegin;
157 179251 : if (!bv->n) PetscFunctionReturn(PETSC_SUCCESS);
158 179251 : PetscCall(MatDenseGetArray(ctx->A,&array));
159 179251 : if (PetscUnlikely(j<0)) PetscCall(BVScale_BLAS_Private(bv,(bv->k-bv->l)*bv->ld,array+(bv->nc+bv->l)*bv->ld,alpha));
160 178897 : else PetscCall(BVScale_BLAS_Private(bv,bv->n,array+(bv->nc+j)*bv->ld,alpha));
161 179251 : PetscCall(MatDenseRestoreArray(ctx->A,&array));
162 179251 : PetscFunctionReturn(PETSC_SUCCESS);
163 : }
164 :
165 8379 : static PetscErrorCode BVNorm_Mat(BV bv,PetscInt j,NormType type,PetscReal *val)
166 : {
167 8379 : BV_MAT *ctx = (BV_MAT*)bv->data;
168 8379 : const PetscScalar *array;
169 :
170 8379 : PetscFunctionBegin;
171 8379 : PetscCall(MatDenseGetArrayRead(ctx->A,&array));
172 8379 : 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));
173 8307 : else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->ld,bv->ld,type,val,ctx->mpi));
174 8379 : PetscCall(MatDenseRestoreArrayRead(ctx->A,&array));
175 8379 : PetscFunctionReturn(PETSC_SUCCESS);
176 : }
177 :
178 2545 : static PetscErrorCode BVNorm_Local_Mat(BV bv,PetscInt j,NormType type,PetscReal *val)
179 : {
180 2545 : BV_MAT *ctx = (BV_MAT*)bv->data;
181 2545 : const PetscScalar *array;
182 :
183 2545 : PetscFunctionBegin;
184 2545 : PetscCall(MatDenseGetArrayRead(ctx->A,&array));
185 2545 : 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));
186 2545 : else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->ld,bv->ld,type,val,PETSC_FALSE));
187 2545 : PetscCall(MatDenseRestoreArrayRead(ctx->A,&array));
188 2545 : PetscFunctionReturn(PETSC_SUCCESS);
189 : }
190 :
191 3868 : static PetscErrorCode BVNormalize_Mat(BV bv,PetscScalar *eigi)
192 : {
193 3868 : BV_MAT *ctx = (BV_MAT*)bv->data;
194 3868 : PetscScalar *array,*wi=NULL;
195 :
196 3868 : PetscFunctionBegin;
197 3868 : PetscCall(MatDenseGetArray(ctx->A,&array));
198 3868 : if (eigi) wi = eigi+bv->l;
199 3868 : PetscCall(BVNormalize_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->ld,bv->ld,wi,ctx->mpi));
200 3868 : PetscCall(MatDenseRestoreArray(ctx->A,&array));
201 3868 : PetscFunctionReturn(PETSC_SUCCESS);
202 : }
203 :
204 11024 : static PetscErrorCode BVMatMult_Mat(BV V,Mat A,BV W)
205 : {
206 11024 : PetscInt j;
207 11024 : Mat Vmat,Wmat;
208 11024 : Vec vv,ww;
209 :
210 11024 : PetscFunctionBegin;
211 11024 : if (V->vmm) {
212 11021 : PetscCall(BVGetMat(V,&Vmat));
213 11021 : PetscCall(BVGetMat(W,&Wmat));
214 11021 : PetscCall(MatProductCreateWithMat(A,Vmat,NULL,Wmat));
215 11021 : PetscCall(MatProductSetType(Wmat,MATPRODUCT_AB));
216 11021 : PetscCall(MatProductSetFromOptions(Wmat));
217 11021 : PetscCall(MatProductSymbolic(Wmat));
218 11021 : PetscCall(MatProductNumeric(Wmat));
219 11021 : PetscCall(MatProductClear(Wmat));
220 11021 : PetscCall(BVRestoreMat(V,&Vmat));
221 11021 : PetscCall(BVRestoreMat(W,&Wmat));
222 : } else {
223 26 : for (j=0;j<V->k-V->l;j++) {
224 23 : PetscCall(BVGetColumn(V,V->l+j,&vv));
225 23 : PetscCall(BVGetColumn(W,W->l+j,&ww));
226 23 : PetscCall(MatMult(A,vv,ww));
227 23 : PetscCall(BVRestoreColumn(V,V->l+j,&vv));
228 23 : PetscCall(BVRestoreColumn(W,W->l+j,&ww));
229 : }
230 : }
231 11024 : PetscFunctionReturn(PETSC_SUCCESS);
232 : }
233 :
234 13206 : static PetscErrorCode BVCopy_Mat(BV V,BV W)
235 : {
236 13206 : BV_MAT *v = (BV_MAT*)V->data,*w = (BV_MAT*)W->data;
237 13206 : const PetscScalar *pv;
238 13206 : PetscScalar *pw;
239 13206 : PetscInt j;
240 :
241 13206 : PetscFunctionBegin;
242 13206 : PetscCall(MatDenseGetArrayRead(v->A,&pv));
243 13206 : PetscCall(MatDenseGetArray(w->A,&pw));
244 84495 : 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));
245 13206 : PetscCall(MatDenseRestoreArrayRead(v->A,&pv));
246 13206 : PetscCall(MatDenseRestoreArray(w->A,&pw));
247 13206 : PetscFunctionReturn(PETSC_SUCCESS);
248 : }
249 :
250 5909 : static PetscErrorCode BVCopyColumn_Mat(BV V,PetscInt j,PetscInt i)
251 : {
252 5909 : BV_MAT *v = (BV_MAT*)V->data;
253 5909 : PetscScalar *pv;
254 :
255 5909 : PetscFunctionBegin;
256 5909 : PetscCall(MatDenseGetArray(v->A,&pv));
257 5909 : PetscCall(PetscArraycpy(pv+(V->nc+i)*V->ld,pv+(V->nc+j)*V->ld,V->n));
258 5909 : PetscCall(MatDenseRestoreArray(v->A,&pv));
259 5909 : PetscFunctionReturn(PETSC_SUCCESS);
260 : }
261 :
262 162 : static PetscErrorCode BVResize_Mat(BV bv,PetscInt m,PetscBool copy)
263 : {
264 162 : BV_MAT *ctx = (BV_MAT*)bv->data;
265 162 : Mat A,Msrc,Mdst;
266 162 : char str[50];
267 :
268 162 : PetscFunctionBegin;
269 162 : PetscCall(MatCreateDenseFromVecType(PetscObjectComm((PetscObject)bv),bv->vtype,bv->n,PETSC_DECIDE,bv->N,m,bv->ld,NULL,&A));
270 162 : if (((PetscObject)bv)->name) {
271 15 : PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
272 15 : PetscCall(PetscObjectSetName((PetscObject)A,str));
273 : }
274 162 : if (copy) {
275 20 : PetscCall(MatDenseGetSubMatrix(ctx->A,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PetscMin(m,bv->m),&Msrc));
276 20 : PetscCall(MatDenseGetSubMatrix(A,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PetscMin(m,bv->m),&Mdst));
277 20 : PetscCall(MatCopy(Msrc,Mdst,SAME_NONZERO_PATTERN));
278 20 : PetscCall(MatDenseRestoreSubMatrix(ctx->A,&Msrc));
279 20 : PetscCall(MatDenseRestoreSubMatrix(A,&Mdst));
280 : }
281 162 : PetscCall(MatDestroy(&ctx->A));
282 162 : ctx->A = A;
283 162 : PetscFunctionReturn(PETSC_SUCCESS);
284 : }
285 :
286 1198086 : static PetscErrorCode BVGetColumn_Mat(BV bv,PetscInt j,Vec *v)
287 : {
288 1198086 : BV_MAT *ctx = (BV_MAT*)bv->data;
289 1198086 : PetscScalar *pA;
290 1198086 : PetscInt l;
291 :
292 1198086 : PetscFunctionBegin;
293 1198086 : l = BVAvailableVec;
294 1198086 : PetscCall(MatDenseGetArray(ctx->A,&pA));
295 1198086 : PetscCall(VecPlaceArray(bv->cv[l],pA+(bv->nc+j)*bv->ld));
296 1198086 : PetscFunctionReturn(PETSC_SUCCESS);
297 : }
298 :
299 1198086 : static PetscErrorCode BVRestoreColumn_Mat(BV bv,PetscInt j,Vec *v)
300 : {
301 1198086 : BV_MAT *ctx = (BV_MAT*)bv->data;
302 1198086 : PetscScalar *pA;
303 1198086 : PetscInt l;
304 :
305 1198086 : PetscFunctionBegin;
306 1198086 : l = (j==bv->ci[0])? 0: 1;
307 1198086 : PetscCall(VecResetArray(bv->cv[l]));
308 1198086 : PetscCall(MatDenseRestoreArray(ctx->A,&pA));
309 1198086 : PetscFunctionReturn(PETSC_SUCCESS);
310 : }
311 :
312 26187 : static PetscErrorCode BVGetArray_Mat(BV bv,PetscScalar **a)
313 : {
314 26187 : BV_MAT *ctx = (BV_MAT*)bv->data;
315 :
316 26187 : PetscFunctionBegin;
317 26187 : PetscCall(MatDenseGetArray(ctx->A,a));
318 26187 : PetscFunctionReturn(PETSC_SUCCESS);
319 : }
320 :
321 26187 : static PetscErrorCode BVRestoreArray_Mat(BV bv,PetscScalar **a)
322 : {
323 26187 : BV_MAT *ctx = (BV_MAT*)bv->data;
324 :
325 26187 : PetscFunctionBegin;
326 26187 : if (a) PetscCall(MatDenseRestoreArray(ctx->A,a));
327 26187 : PetscFunctionReturn(PETSC_SUCCESS);
328 : }
329 :
330 166 : static PetscErrorCode BVGetArrayRead_Mat(BV bv,const PetscScalar **a)
331 : {
332 166 : BV_MAT *ctx = (BV_MAT*)bv->data;
333 :
334 166 : PetscFunctionBegin;
335 166 : PetscCall(MatDenseGetArrayRead(ctx->A,a));
336 166 : PetscFunctionReturn(PETSC_SUCCESS);
337 : }
338 :
339 166 : static PetscErrorCode BVRestoreArrayRead_Mat(BV bv,const PetscScalar **a)
340 : {
341 166 : BV_MAT *ctx = (BV_MAT*)bv->data;
342 :
343 166 : PetscFunctionBegin;
344 166 : if (a) PetscCall(MatDenseRestoreArrayRead(ctx->A,a));
345 166 : PetscFunctionReturn(PETSC_SUCCESS);
346 : }
347 :
348 29 : static PetscErrorCode BVView_Mat(BV bv,PetscViewer viewer)
349 : {
350 29 : BV_MAT *ctx = (BV_MAT*)bv->data;
351 29 : PetscViewerFormat format;
352 29 : PetscBool isascii;
353 29 : const char *bvname,*name;
354 :
355 29 : PetscFunctionBegin;
356 29 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
357 29 : if (isascii) {
358 29 : PetscCall(PetscViewerGetFormat(viewer,&format));
359 29 : if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) PetscFunctionReturn(PETSC_SUCCESS);
360 10 : PetscCall(MatView(ctx->A,viewer));
361 10 : if (format == PETSC_VIEWER_ASCII_MATLAB) {
362 0 : PetscCall(PetscObjectGetName((PetscObject)bv,&bvname));
363 0 : PetscCall(PetscObjectGetName((PetscObject)ctx->A,&name));
364 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"%s=%s;clear %s\n",bvname,name,name));
365 0 : if (bv->nc) PetscCall(PetscViewerASCIIPrintf(viewer,"%s=%s(:,%" PetscInt_FMT ":end);\n",bvname,bvname,bv->nc+1));
366 : }
367 0 : } else PetscCall(MatView(ctx->A,viewer));
368 10 : PetscFunctionReturn(PETSC_SUCCESS);
369 : }
370 :
371 14014 : static PetscErrorCode BVDestroy_Mat(BV bv)
372 : {
373 14014 : BV_MAT *ctx = (BV_MAT*)bv->data;
374 :
375 14014 : PetscFunctionBegin;
376 14014 : PetscCall(MatDestroy(&ctx->A));
377 14014 : PetscCall(VecDestroy(&bv->cv[0]));
378 14014 : PetscCall(VecDestroy(&bv->cv[1]));
379 14014 : PetscCall(PetscFree(bv->data));
380 14014 : PetscFunctionReturn(PETSC_SUCCESS);
381 : }
382 :
383 14014 : SLEPC_EXTERN PetscErrorCode BVCreate_Mat(BV bv)
384 : {
385 14014 : BV_MAT *ctx;
386 14014 : PetscInt nloc,lsplit;
387 14014 : PetscBool seq;
388 14014 : char str[50];
389 14014 : PetscScalar *array,*ptr=NULL;
390 14014 : BV parent;
391 14014 : Mat Apar;
392 :
393 14014 : PetscFunctionBegin;
394 14014 : PetscCall(PetscNew(&ctx));
395 14014 : bv->data = (void*)ctx;
396 :
397 14014 : PetscCall(PetscStrcmpAny(bv->vtype,&bv->cuda,VECSEQCUDA,VECMPICUDA,""));
398 14014 : PetscCall(PetscStrcmpAny(bv->vtype,&ctx->mpi,VECMPI,VECMPICUDA,""));
399 :
400 14014 : PetscCall(PetscStrcmp(bv->vtype,VECSEQ,&seq));
401 14014 : PetscCheck(seq || ctx->mpi || bv->cuda,PetscObjectComm((PetscObject)bv),PETSC_ERR_SUP,"BVMAT does not support the requested vector type: %s",bv->vtype);
402 :
403 14014 : PetscCall(PetscLayoutGetLocalSize(bv->map,&nloc));
404 14014 : PetscCall(BV_SetDefaultLD(bv,nloc));
405 :
406 14014 : if (PetscUnlikely(bv->issplit)) {
407 : /* split BV: share the memory of the parent BV */
408 80 : parent = bv->splitparent;
409 80 : lsplit = parent->lsplit;
410 80 : Apar = ((BV_MAT*)parent->data)->A;
411 80 : if (bv->cuda) {
412 : #if defined(PETSC_HAVE_CUDA)
413 : PetscCall(MatDenseCUDAGetArray(Apar,&array));
414 : ptr = (bv->issplit==1)? array: array+lsplit*bv->ld;
415 : PetscCall(MatDenseCUDARestoreArray(Apar,&array));
416 : #endif
417 : } else {
418 80 : PetscCall(MatDenseGetArray(Apar,&array));
419 80 : ptr = (bv->issplit==1)? array: array+lsplit*bv->ld;
420 80 : PetscCall(MatDenseRestoreArray(Apar,&array));
421 : }
422 : }
423 :
424 14014 : PetscCall(MatCreateDenseFromVecType(PetscObjectComm((PetscObject)bv),bv->vtype,nloc,PETSC_DECIDE,bv->N,bv->m,bv->ld,ptr,&ctx->A));
425 14014 : if (((PetscObject)bv)->name) {
426 114 : PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
427 114 : PetscCall(PetscObjectSetName((PetscObject)ctx->A,str));
428 : }
429 :
430 14014 : if (PetscUnlikely(bv->Acreate)) {
431 150 : PetscCall(MatConvert(bv->Acreate,bv->cuda?MATDENSECUDA:MATDENSE,MAT_REUSE_MATRIX,&ctx->A));
432 75 : PetscCall(MatDestroy(&bv->Acreate));
433 : }
434 :
435 14014 : PetscCall(BVCreateVecEmpty(bv,&bv->cv[0]));
436 14014 : PetscCall(BVCreateVecEmpty(bv,&bv->cv[1]));
437 :
438 14014 : if (bv->cuda) {
439 : #if defined(PETSC_HAVE_CUDA)
440 : bv->ops->mult = BVMult_Mat_CUDA;
441 : bv->ops->multvec = BVMultVec_Mat_CUDA;
442 : bv->ops->multinplace = BVMultInPlace_Mat_CUDA;
443 : bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Mat_CUDA;
444 : bv->ops->dot = BVDot_Mat_CUDA;
445 : bv->ops->dotvec = BVDotVec_Mat_CUDA;
446 : bv->ops->dotvec_local = BVDotVec_Local_Mat_CUDA;
447 : bv->ops->scale = BVScale_Mat_CUDA;
448 : bv->ops->norm = BVNorm_Mat_CUDA;
449 : bv->ops->norm_local = BVNorm_Local_Mat_CUDA;
450 : bv->ops->normalize = BVNormalize_Mat_CUDA;
451 : bv->ops->matmult = BVMatMult_Mat_CUDA;
452 : bv->ops->copy = BVCopy_Mat_CUDA;
453 : bv->ops->copycolumn = BVCopyColumn_Mat_CUDA;
454 : bv->ops->getcolumn = BVGetColumn_Mat_CUDA;
455 : bv->ops->restorecolumn = BVRestoreColumn_Mat_CUDA;
456 : bv->ops->restoresplit = BVRestoreSplit_Mat_CUDA;
457 : bv->ops->getmat = BVGetMat_Mat_CUDA;
458 : bv->ops->restoremat = BVRestoreMat_Mat_CUDA;
459 : #endif
460 : } else {
461 14014 : bv->ops->mult = BVMult_Mat;
462 14014 : bv->ops->multvec = BVMultVec_Mat;
463 14014 : bv->ops->multinplace = BVMultInPlace_Mat;
464 14014 : bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Mat;
465 14014 : bv->ops->dot = BVDot_Mat;
466 14014 : bv->ops->dotvec = BVDotVec_Mat;
467 14014 : bv->ops->dotvec_local = BVDotVec_Local_Mat;
468 14014 : bv->ops->scale = BVScale_Mat;
469 14014 : bv->ops->norm = BVNorm_Mat;
470 14014 : bv->ops->norm_local = BVNorm_Local_Mat;
471 14014 : bv->ops->normalize = BVNormalize_Mat;
472 14014 : bv->ops->matmult = BVMatMult_Mat;
473 14014 : bv->ops->copy = BVCopy_Mat;
474 14014 : bv->ops->copycolumn = BVCopyColumn_Mat;
475 14014 : bv->ops->getcolumn = BVGetColumn_Mat;
476 14014 : bv->ops->restorecolumn = BVRestoreColumn_Mat;
477 14014 : bv->ops->getmat = BVGetMat_Default;
478 14014 : bv->ops->restoremat = BVRestoreMat_Default;
479 : }
480 14014 : bv->ops->resize = BVResize_Mat;
481 14014 : bv->ops->getarray = BVGetArray_Mat;
482 14014 : bv->ops->restorearray = BVRestoreArray_Mat;
483 14014 : bv->ops->getarrayread = BVGetArrayRead_Mat;
484 14014 : bv->ops->restorearrayread = BVRestoreArrayRead_Mat;
485 14014 : bv->ops->destroy = BVDestroy_Mat;
486 14014 : if (!ctx->mpi) bv->ops->view = BVView_Mat;
487 14014 : PetscFunctionReturn(PETSC_SUCCESS);
488 : }
|