LCOV - code coverage report
Current view: top level - sys/classes/bv/impls/contiguous - contig.c (source / functions) Hit Total Coverage
Test: SLEPc Lines: 266 269 98.9 %
Date: 2021-08-02 00:35:43 Functions: 20 20 100.0 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /*
       2             :    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
       3             :    SLEPc - Scalable Library for Eigenvalue Problem Computations
       4             :    Copyright (c) 2002-2021, 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 as an array of Vecs sharing a contiguous array for elements
      12             : */
      13             : 
      14             : #include <slepc/private/bvimpl.h>
      15             : 
      16             : typedef struct {
      17             :   Vec         *V;
      18             :   PetscScalar *array;
      19             :   PetscBool   mpi;
      20             : } BV_CONTIGUOUS;
      21             : 
      22          99 : PetscErrorCode BVMult_Contiguous(BV Y,PetscScalar alpha,PetscScalar beta,BV X,Mat Q)
      23             : {
      24          99 :   PetscErrorCode    ierr;
      25          99 :   BV_CONTIGUOUS     *y = (BV_CONTIGUOUS*)Y->data,*x = (BV_CONTIGUOUS*)X->data;
      26          99 :   const PetscScalar *q;
      27          99 :   PetscInt          ldq;
      28             : 
      29          99 :   PetscFunctionBegin;
      30          99 :   if (Q) {
      31          79 :     ierr = MatGetSize(Q,&ldq,NULL);CHKERRQ(ierr);
      32          79 :     ierr = MatDenseGetArrayRead(Q,&q);CHKERRQ(ierr);
      33          79 :     ierr = BVMult_BLAS_Private(Y,Y->n,Y->k-Y->l,X->k-X->l,ldq,alpha,x->array+(X->nc+X->l)*X->n,q+Y->l*ldq+X->l,beta,y->array+(Y->nc+Y->l)*Y->n);CHKERRQ(ierr);
      34          79 :     ierr = MatDenseRestoreArrayRead(Q,&q);CHKERRQ(ierr);
      35             :   } else {
      36          20 :     ierr = BVAXPY_BLAS_Private(Y,Y->n,Y->k-Y->l,alpha,x->array+(X->nc+X->l)*X->n,beta,y->array+(Y->nc+Y->l)*Y->n);CHKERRQ(ierr);
      37             :   }
      38          99 :   PetscFunctionReturn(0);
      39             : }
      40             : 
      41         345 : PetscErrorCode BVMultVec_Contiguous(BV X,PetscScalar alpha,PetscScalar beta,Vec y,PetscScalar *q)
      42             : {
      43         345 :   PetscErrorCode ierr;
      44         345 :   BV_CONTIGUOUS  *x = (BV_CONTIGUOUS*)X->data;
      45         345 :   PetscScalar    *py,*qq=q;
      46             : 
      47         345 :   PetscFunctionBegin;
      48         345 :   ierr = VecGetArray(y,&py);CHKERRQ(ierr);
      49         345 :   if (!q) { ierr = VecGetArray(X->buffer,&qq);CHKERRQ(ierr); }
      50         345 :   ierr = BVMultVec_BLAS_Private(X,X->n,X->k-X->l,alpha,x->array+(X->nc+X->l)*X->n,qq,beta,py);CHKERRQ(ierr);
      51         345 :   if (!q) { ierr = VecRestoreArray(X->buffer,&qq);CHKERRQ(ierr); }
      52         345 :   ierr = VecRestoreArray(y,&py);CHKERRQ(ierr);
      53         345 :   PetscFunctionReturn(0);
      54             : }
      55             : 
      56          45 : PetscErrorCode BVMultInPlace_Contiguous(BV V,Mat Q,PetscInt s,PetscInt e)
      57             : {
      58          45 :   PetscErrorCode    ierr;
      59          45 :   BV_CONTIGUOUS     *ctx = (BV_CONTIGUOUS*)V->data;
      60          45 :   const PetscScalar *q;
      61          45 :   PetscInt          ldq;
      62             : 
      63          45 :   PetscFunctionBegin;
      64          45 :   ierr = MatGetSize(Q,&ldq,NULL);CHKERRQ(ierr);
      65          45 :   ierr = MatDenseGetArrayRead(Q,&q);CHKERRQ(ierr);
      66          45 :   ierr = BVMultInPlace_BLAS_Private(V,V->n,V->k-V->l,ldq,s-V->l,e-V->l,ctx->array+(V->nc+V->l)*V->n,q+V->l*ldq+V->l,PETSC_FALSE);CHKERRQ(ierr);
      67          45 :   ierr = MatDenseRestoreArrayRead(Q,&q);CHKERRQ(ierr);
      68          45 :   PetscFunctionReturn(0);
      69             : }
      70             : 
      71           1 : PetscErrorCode BVMultInPlaceTranspose_Contiguous(BV V,Mat Q,PetscInt s,PetscInt e)
      72             : {
      73           1 :   PetscErrorCode    ierr;
      74           1 :   BV_CONTIGUOUS     *ctx = (BV_CONTIGUOUS*)V->data;
      75           1 :   const PetscScalar *q;
      76           1 :   PetscInt          ldq;
      77             : 
      78           1 :   PetscFunctionBegin;
      79           1 :   ierr = MatGetSize(Q,&ldq,NULL);CHKERRQ(ierr);
      80           1 :   ierr = MatDenseGetArrayRead(Q,&q);CHKERRQ(ierr);
      81           1 :   ierr = BVMultInPlace_BLAS_Private(V,V->n,V->k-V->l,ldq,s-V->l,e-V->l,ctx->array+(V->nc+V->l)*V->n,q+V->l*ldq+V->l,PETSC_TRUE);CHKERRQ(ierr);
      82           1 :   ierr = MatDenseRestoreArrayRead(Q,&q);CHKERRQ(ierr);
      83           1 :   PetscFunctionReturn(0);
      84             : }
      85             : 
      86         236 : PetscErrorCode BVDot_Contiguous(BV X,BV Y,Mat M)
      87             : {
      88         236 :   PetscErrorCode ierr;
      89         236 :   BV_CONTIGUOUS  *x = (BV_CONTIGUOUS*)X->data,*y = (BV_CONTIGUOUS*)Y->data;
      90         236 :   PetscScalar    *m;
      91         236 :   PetscInt       ldm;
      92             : 
      93         236 :   PetscFunctionBegin;
      94         236 :   ierr = MatGetSize(M,&ldm,NULL);CHKERRQ(ierr);
      95         236 :   ierr = MatDenseGetArray(M,&m);CHKERRQ(ierr);
      96         236 :   ierr = BVDot_BLAS_Private(X,Y->k-Y->l,X->k-X->l,X->n,ldm,y->array+(Y->nc+Y->l)*Y->n,x->array+(X->nc+X->l)*X->n,m+X->l*ldm+Y->l,x->mpi);CHKERRQ(ierr);
      97         236 :   ierr = MatDenseRestoreArray(M,&m);CHKERRQ(ierr);
      98         236 :   PetscFunctionReturn(0);
      99             : }
     100             : 
     101         353 : PetscErrorCode BVDotVec_Contiguous(BV X,Vec y,PetscScalar *q)
     102             : {
     103         353 :   PetscErrorCode    ierr;
     104         353 :   BV_CONTIGUOUS     *x = (BV_CONTIGUOUS*)X->data;
     105         353 :   const PetscScalar *py;
     106         353 :   PetscScalar       *qq=q;
     107         353 :   Vec               z = y;
     108             : 
     109         353 :   PetscFunctionBegin;
     110         353 :   if (X->matrix) {
     111          57 :     ierr = BV_IPMatMult(X,y);CHKERRQ(ierr);
     112          57 :     z = X->Bx;
     113             :   }
     114         353 :   ierr = VecGetArrayRead(z,&py);CHKERRQ(ierr);
     115         353 :   if (!q) { ierr = VecGetArray(X->buffer,&qq);CHKERRQ(ierr); }
     116         353 :   ierr = BVDotVec_BLAS_Private(X,X->n,X->k-X->l,x->array+(X->nc+X->l)*X->n,py,qq,x->mpi);CHKERRQ(ierr);
     117         353 :   if (!q) { ierr = VecRestoreArray(X->buffer,&qq);CHKERRQ(ierr); }
     118         353 :   ierr = VecRestoreArrayRead(z,&py);CHKERRQ(ierr);
     119         353 :   PetscFunctionReturn(0);
     120             : }
     121             : 
     122           4 : PetscErrorCode BVDotVec_Local_Contiguous(BV X,Vec y,PetscScalar *m)
     123             : {
     124           4 :   PetscErrorCode ierr;
     125           4 :   BV_CONTIGUOUS  *x = (BV_CONTIGUOUS*)X->data;
     126           4 :   PetscScalar    *py;
     127           4 :   Vec            z = y;
     128             : 
     129           4 :   PetscFunctionBegin;
     130           4 :   if (X->matrix) {
     131           0 :     ierr = BV_IPMatMult(X,y);CHKERRQ(ierr);
     132           0 :     z = X->Bx;
     133             :   }
     134           4 :   ierr = VecGetArray(z,&py);CHKERRQ(ierr);
     135           4 :   ierr = BVDotVec_BLAS_Private(X,X->n,X->k-X->l,x->array+(X->nc+X->l)*X->n,py,m,PETSC_FALSE);CHKERRQ(ierr);
     136           4 :   ierr = VecRestoreArray(z,&py);CHKERRQ(ierr);
     137           4 :   PetscFunctionReturn(0);
     138             : }
     139             : 
     140         323 : PetscErrorCode BVScale_Contiguous(BV bv,PetscInt j,PetscScalar alpha)
     141             : {
     142         323 :   PetscErrorCode ierr;
     143         323 :   BV_CONTIGUOUS  *ctx = (BV_CONTIGUOUS*)bv->data;
     144             : 
     145         323 :   PetscFunctionBegin;
     146         323 :   if (j<0) {
     147           3 :     ierr = BVScale_BLAS_Private(bv,(bv->k-bv->l)*bv->n,ctx->array+(bv->nc+bv->l)*bv->n,alpha);CHKERRQ(ierr);
     148             :   } else {
     149         320 :     ierr = BVScale_BLAS_Private(bv,bv->n,ctx->array+(bv->nc+j)*bv->n,alpha);CHKERRQ(ierr);
     150             :   }
     151         323 :   PetscFunctionReturn(0);
     152             : }
     153             : 
     154         108 : PetscErrorCode BVNorm_Contiguous(BV bv,PetscInt j,NormType type,PetscReal *val)
     155             : {
     156         108 :   PetscErrorCode ierr;
     157         108 :   BV_CONTIGUOUS  *ctx = (BV_CONTIGUOUS*)bv->data;
     158             : 
     159         108 :   PetscFunctionBegin;
     160         108 :   if (j<0) {
     161          69 :     ierr = BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,ctx->array+(bv->nc+bv->l)*bv->n,type,val,ctx->mpi);CHKERRQ(ierr);
     162             :   } else {
     163          39 :     ierr = BVNorm_LAPACK_Private(bv,bv->n,1,ctx->array+(bv->nc+j)*bv->n,type,val,ctx->mpi);CHKERRQ(ierr);
     164             :   }
     165         108 :   PetscFunctionReturn(0);
     166             : }
     167             : 
     168           2 : PetscErrorCode BVNorm_Local_Contiguous(BV bv,PetscInt j,NormType type,PetscReal *val)
     169             : {
     170           2 :   PetscErrorCode ierr;
     171           2 :   BV_CONTIGUOUS  *ctx = (BV_CONTIGUOUS*)bv->data;
     172             : 
     173           2 :   PetscFunctionBegin;
     174           2 :   if (j<0) {
     175           0 :     ierr = BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,ctx->array+(bv->nc+bv->l)*bv->n,type,val,PETSC_FALSE);CHKERRQ(ierr);
     176             :   } else {
     177           2 :     ierr = BVNorm_LAPACK_Private(bv,bv->n,1,ctx->array+(bv->nc+j)*bv->n,type,val,PETSC_FALSE);CHKERRQ(ierr);
     178             :   }
     179           2 :   PetscFunctionReturn(0);
     180             : }
     181             : 
     182           3 : PetscErrorCode BVNormalize_Contiguous(BV bv,PetscScalar *eigi)
     183             : {
     184           3 :   PetscErrorCode ierr;
     185           3 :   BV_CONTIGUOUS  *ctx = (BV_CONTIGUOUS*)bv->data;
     186           3 :   PetscScalar    *wi=NULL;
     187             : 
     188           3 :   PetscFunctionBegin;
     189           3 :   if (eigi) wi = eigi+bv->l;
     190           3 :   ierr = BVNormalize_LAPACK_Private(bv,bv->n,bv->k-bv->l,ctx->array+(bv->nc+bv->l)*bv->n,wi,ctx->mpi);CHKERRQ(ierr);
     191           3 :   PetscFunctionReturn(0);
     192             : }
     193             : 
     194         105 : PetscErrorCode BVMatMult_Contiguous(BV V,Mat A,BV W)
     195             : {
     196         105 :   PetscErrorCode ierr;
     197         105 :   BV_CONTIGUOUS  *v = (BV_CONTIGUOUS*)V->data,*w = (BV_CONTIGUOUS*)W->data;
     198         105 :   PetscInt       j;
     199         105 :   Mat            Vmat,Wmat;
     200             : 
     201         105 :   PetscFunctionBegin;
     202         105 :   if (V->vmm) {
     203         102 :     ierr = BVGetMat(V,&Vmat);CHKERRQ(ierr);
     204         102 :     ierr = BVGetMat(W,&Wmat);CHKERRQ(ierr);
     205         102 :     ierr = MatProductCreateWithMat(A,Vmat,NULL,Wmat);CHKERRQ(ierr);
     206         102 :     ierr = MatProductSetType(Wmat,MATPRODUCT_AB);CHKERRQ(ierr);
     207         102 :     ierr = MatProductSetFromOptions(Wmat);CHKERRQ(ierr);
     208         102 :     ierr = MatProductSymbolic(Wmat);CHKERRQ(ierr);
     209         102 :     ierr = MatProductNumeric(Wmat);CHKERRQ(ierr);
     210         102 :     ierr = MatProductClear(Wmat);CHKERRQ(ierr);
     211         102 :     ierr = BVRestoreMat(V,&Vmat);CHKERRQ(ierr);
     212         102 :     ierr = BVRestoreMat(W,&Wmat);CHKERRQ(ierr);
     213             :   } else {
     214          26 :     for (j=0;j<V->k-V->l;j++) {
     215          23 :       ierr = MatMult(A,v->V[V->nc+V->l+j],w->V[W->nc+W->l+j]);CHKERRQ(ierr);
     216             :     }
     217             :   }
     218         105 :   PetscFunctionReturn(0);
     219             : }
     220             : 
     221          99 : PetscErrorCode BVCopy_Contiguous(BV V,BV W)
     222             : {
     223          99 :   PetscErrorCode ierr;
     224          99 :   BV_CONTIGUOUS  *v = (BV_CONTIGUOUS*)V->data,*w = (BV_CONTIGUOUS*)W->data;
     225          99 :   PetscScalar    *pvc,*pwc;
     226             : 
     227          99 :   PetscFunctionBegin;
     228          99 :   pvc = v->array+(V->nc+V->l)*V->n;
     229          99 :   pwc = w->array+(W->nc+W->l)*W->n;
     230          99 :   ierr = PetscArraycpy(pwc,pvc,(V->k-V->l)*V->n);CHKERRQ(ierr);
     231          99 :   PetscFunctionReturn(0);
     232             : }
     233             : 
     234           2 : PetscErrorCode BVCopyColumn_Contiguous(BV V,PetscInt j,PetscInt i)
     235             : {
     236           2 :   PetscErrorCode ierr;
     237           2 :   BV_CONTIGUOUS  *v = (BV_CONTIGUOUS*)V->data;
     238             : 
     239           2 :   PetscFunctionBegin;
     240           2 :   ierr = PetscArraycpy(v->array+(V->nc+i)*V->n,v->array+(V->nc+j)*V->n,V->n);CHKERRQ(ierr);
     241           2 :   PetscFunctionReturn(0);
     242             : }
     243             : 
     244           8 : PetscErrorCode BVResize_Contiguous(BV bv,PetscInt m,PetscBool copy)
     245             : {
     246           8 :   PetscErrorCode ierr;
     247           8 :   BV_CONTIGUOUS  *ctx = (BV_CONTIGUOUS*)bv->data;
     248           8 :   PetscInt       j,bs;
     249           8 :   PetscScalar    *newarray;
     250           8 :   Vec            *newV;
     251           8 :   char           str[50];
     252             : 
     253           8 :   PetscFunctionBegin;
     254           8 :   ierr = VecGetBlockSize(bv->t,&bs);CHKERRQ(ierr);
     255           8 :   ierr = PetscMalloc1(m*bv->n,&newarray);CHKERRQ(ierr);
     256           8 :   ierr = PetscArrayzero(newarray,m*bv->n);CHKERRQ(ierr);
     257           8 :   ierr = PetscMalloc1(m,&newV);CHKERRQ(ierr);
     258          94 :   for (j=0;j<m;j++) {
     259          86 :     if (ctx->mpi) {
     260          34 :       ierr = VecCreateMPIWithArray(PetscObjectComm((PetscObject)bv->t),bs,bv->n,PETSC_DECIDE,newarray+j*bv->n,newV+j);CHKERRQ(ierr);
     261             :     } else {
     262          86 :       ierr = VecCreateSeqWithArray(PetscObjectComm((PetscObject)bv->t),bs,bv->n,newarray+j*bv->n,newV+j);CHKERRQ(ierr);
     263             :     }
     264             :   }
     265          94 :   ierr = PetscLogObjectParents(bv,m,newV);CHKERRQ(ierr);
     266           8 :   if (((PetscObject)bv)->name) {
     267          94 :     for (j=0;j<m;j++) {
     268          86 :       ierr = PetscSNPrintf(str,sizeof(str),"%s_%d",((PetscObject)bv)->name,(int)j);CHKERRQ(ierr);
     269          86 :       ierr = PetscObjectSetName((PetscObject)newV[j],str);CHKERRQ(ierr);
     270             :     }
     271             :   }
     272           8 :   if (copy) {
     273           2 :     ierr = PetscArraycpy(newarray,ctx->array,PetscMin(m,bv->m)*bv->n);CHKERRQ(ierr);
     274             :   }
     275           8 :   ierr = VecDestroyVecs(bv->m,&ctx->V);CHKERRQ(ierr);
     276           8 :   ctx->V = newV;
     277           8 :   ierr = PetscFree(ctx->array);CHKERRQ(ierr);
     278           8 :   ctx->array = newarray;
     279           8 :   PetscFunctionReturn(0);
     280             : }
     281             : 
     282        2503 : PetscErrorCode BVGetColumn_Contiguous(BV bv,PetscInt j,Vec *v)
     283             : {
     284        2503 :   BV_CONTIGUOUS *ctx = (BV_CONTIGUOUS*)bv->data;
     285        2503 :   PetscInt      l;
     286             : 
     287        2503 :   PetscFunctionBegin;
     288        2503 :   l = BVAvailableVec;
     289        2503 :   bv->cv[l] = ctx->V[bv->nc+j];
     290        2503 :   PetscFunctionReturn(0);
     291             : }
     292             : 
     293         251 : PetscErrorCode BVGetArray_Contiguous(BV bv,PetscScalar **a)
     294             : {
     295         251 :   BV_CONTIGUOUS *ctx = (BV_CONTIGUOUS*)bv->data;
     296             : 
     297         251 :   PetscFunctionBegin;
     298         251 :   *a = ctx->array;
     299         251 :   PetscFunctionReturn(0);
     300             : }
     301             : 
     302           9 : PetscErrorCode BVGetArrayRead_Contiguous(BV bv,const PetscScalar **a)
     303             : {
     304           9 :   BV_CONTIGUOUS *ctx = (BV_CONTIGUOUS*)bv->data;
     305             : 
     306           9 :   PetscFunctionBegin;
     307           9 :   *a = ctx->array;
     308           9 :   PetscFunctionReturn(0);
     309             : }
     310             : 
     311         329 : PetscErrorCode BVDestroy_Contiguous(BV bv)
     312             : {
     313         329 :   PetscErrorCode ierr;
     314         329 :   BV_CONTIGUOUS  *ctx = (BV_CONTIGUOUS*)bv->data;
     315             : 
     316         329 :   PetscFunctionBegin;
     317         329 :   if (!bv->issplit) {
     318         249 :     ierr = VecDestroyVecs(bv->nc+bv->m,&ctx->V);CHKERRQ(ierr);
     319         249 :     ierr = PetscFree(ctx->array);CHKERRQ(ierr);
     320             :   }
     321         329 :   ierr = PetscFree(bv->data);CHKERRQ(ierr);
     322         329 :   PetscFunctionReturn(0);
     323             : }
     324             : 
     325         329 : SLEPC_EXTERN PetscErrorCode BVCreate_Contiguous(BV bv)
     326             : {
     327         329 :   PetscErrorCode ierr;
     328         329 :   BV_CONTIGUOUS  *ctx;
     329         329 :   PetscInt       j,nloc,bs,lsplit;
     330         329 :   PetscBool      seq;
     331         329 :   PetscScalar    *aa;
     332         329 :   char           str[50];
     333         329 :   PetscScalar    *array;
     334         329 :   BV             parent;
     335         329 :   Vec            *Vpar;
     336             : 
     337         329 :   PetscFunctionBegin;
     338         329 :   ierr = PetscNewLog(bv,&ctx);CHKERRQ(ierr);
     339         329 :   bv->data = (void*)ctx;
     340             : 
     341         329 :   ierr = PetscObjectTypeCompare((PetscObject)bv->t,VECMPI,&ctx->mpi);CHKERRQ(ierr);
     342         329 :   if (!ctx->mpi) {
     343          59 :     ierr = PetscObjectTypeCompare((PetscObject)bv->t,VECSEQ,&seq);CHKERRQ(ierr);
     344          59 :     if (!seq) SETERRQ(PetscObjectComm((PetscObject)bv),PETSC_ERR_SUP,"Cannot create a contiguous BV from a non-standard template vector");
     345             :   }
     346             : 
     347         329 :   ierr = VecGetLocalSize(bv->t,&nloc);CHKERRQ(ierr);
     348         329 :   ierr = VecGetBlockSize(bv->t,&bs);CHKERRQ(ierr);
     349             : 
     350         329 :   if (bv->issplit) {
     351             :     /* split BV: share memory and Vecs of the parent BV */
     352          80 :     parent = bv->splitparent;
     353          80 :     lsplit = parent->lsplit;
     354          80 :     Vpar   = ((BV_CONTIGUOUS*)parent->data)->V;
     355          80 :     ctx->V = (bv->issplit==1)? Vpar: Vpar+lsplit;
     356          80 :     array  = ((BV_CONTIGUOUS*)parent->data)->array;
     357          80 :     ctx->array = (bv->issplit==1)? array: array+lsplit*nloc;
     358             :   } else {
     359             :     /* regular BV: allocate memory and Vecs for the BV entries */
     360         249 :     ierr = PetscCalloc1(bv->m*nloc,&ctx->array);CHKERRQ(ierr);
     361         249 :     ierr = PetscMalloc1(bv->m,&ctx->V);CHKERRQ(ierr);
     362        2167 :     for (j=0;j<bv->m;j++) {
     363        1918 :       if (ctx->mpi) {
     364        1500 :         ierr = VecCreateMPIWithArray(PetscObjectComm((PetscObject)bv->t),bs,nloc,PETSC_DECIDE,ctx->array+j*nloc,ctx->V+j);CHKERRQ(ierr);
     365             :       } else {
     366        1918 :         ierr = VecCreateSeqWithArray(PetscObjectComm((PetscObject)bv->t),bs,nloc,ctx->array+j*nloc,ctx->V+j);CHKERRQ(ierr);
     367             :       }
     368             :     }
     369        2167 :     ierr = PetscLogObjectParents(bv,bv->m,ctx->V);CHKERRQ(ierr);
     370             :   }
     371         329 :   if (((PetscObject)bv)->name) {
     372         933 :     for (j=0;j<bv->m;j++) {
     373         828 :       ierr = PetscSNPrintf(str,sizeof(str),"%s_%d",((PetscObject)bv)->name,(int)j);CHKERRQ(ierr);
     374         828 :       ierr = PetscObjectSetName((PetscObject)ctx->V[j],str);CHKERRQ(ierr);
     375             :     }
     376             :   }
     377             : 
     378         329 :   if (bv->Acreate) {
     379           2 :     ierr = MatDenseGetArray(bv->Acreate,&aa);CHKERRQ(ierr);
     380           2 :     ierr = PetscArraycpy(ctx->array,aa,bv->m*nloc);CHKERRQ(ierr);
     381           2 :     ierr = MatDenseRestoreArray(bv->Acreate,&aa);CHKERRQ(ierr);
     382           2 :     ierr = MatDestroy(&bv->Acreate);CHKERRQ(ierr);
     383             :   }
     384             : 
     385         329 :   bv->ops->mult             = BVMult_Contiguous;
     386         329 :   bv->ops->multvec          = BVMultVec_Contiguous;
     387         329 :   bv->ops->multinplace      = BVMultInPlace_Contiguous;
     388         329 :   bv->ops->multinplacetrans = BVMultInPlaceTranspose_Contiguous;
     389         329 :   bv->ops->dot              = BVDot_Contiguous;
     390         329 :   bv->ops->dotvec           = BVDotVec_Contiguous;
     391         329 :   bv->ops->dotvec_local     = BVDotVec_Local_Contiguous;
     392         329 :   bv->ops->scale            = BVScale_Contiguous;
     393         329 :   bv->ops->norm             = BVNorm_Contiguous;
     394         329 :   bv->ops->norm_local       = BVNorm_Local_Contiguous;
     395         329 :   bv->ops->normalize        = BVNormalize_Contiguous;
     396         329 :   bv->ops->matmult          = BVMatMult_Contiguous;
     397         329 :   bv->ops->copy             = BVCopy_Contiguous;
     398         329 :   bv->ops->copycolumn       = BVCopyColumn_Contiguous;
     399         329 :   bv->ops->resize           = BVResize_Contiguous;
     400         329 :   bv->ops->getcolumn        = BVGetColumn_Contiguous;
     401         329 :   bv->ops->getarray         = BVGetArray_Contiguous;
     402         329 :   bv->ops->getarrayread     = BVGetArrayRead_Contiguous;
     403         329 :   bv->ops->destroy          = BVDestroy_Contiguous;
     404         329 :   PetscFunctionReturn(0);
     405             : }
     406             : 

Generated by: LCOV version 1.14