LCOV - code coverage report
Current view: top level - sys/classes/bv/interface - bvfunc.c (source / functions) Hit Total Coverage
Test: SLEPc Lines: 274 309 88.7 %
Date: 2024-03-28 00:28:38 Functions: 13 16 81.2 %
Legend: Lines: hit not hit

          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 (basis vectors) interface routines, callable by users
      12             : */
      13             : 
      14             : #include <slepc/private/bvimpl.h>            /*I "slepcbv.h" I*/
      15             : 
      16             : PetscClassId     BV_CLASSID = 0;
      17             : PetscLogEvent    BV_Create = 0,BV_Copy = 0,BV_Mult = 0,BV_MultVec = 0,BV_MultInPlace = 0,BV_Dot = 0,BV_DotVec = 0,BV_Orthogonalize = 0,BV_OrthogonalizeVec = 0,BV_Scale = 0,BV_Norm = 0,BV_NormVec = 0,BV_Normalize = 0,BV_SetRandom = 0,BV_MatMult = 0,BV_MatMultVec = 0,BV_MatProject = 0,BV_SVDAndRank = 0;
      18             : static PetscBool BVPackageInitialized = PETSC_FALSE;
      19             : MPI_Op MPIU_TSQR = 0,MPIU_LAPY2;
      20             : 
      21             : const char *BVOrthogTypes[] = {"CGS","MGS","BVOrthogType","BV_ORTHOG_",NULL};
      22             : const char *BVOrthogRefineTypes[] = {"IFNEEDED","NEVER","ALWAYS","BVOrthogRefineType","BV_ORTHOG_REFINE_",NULL};
      23             : const char *BVOrthogBlockTypes[] = {"GS","CHOL","TSQR","TSQRCHOL","SVQB","BVOrthogBlockType","BV_ORTHOG_BLOCK_",NULL};
      24             : const char *BVMatMultTypes[] = {"VECS","MAT","MAT_SAVE","BVMatMultType","BV_MATMULT_",NULL};
      25             : const char *BVSVDMethods[] = {"REFINE","QR","QR_CAA","BVSVDMethod","BV_SVD_METHOD_",NULL};
      26             : 
      27             : /*@C
      28             :    BVFinalizePackage - This function destroys everything in the Slepc interface
      29             :    to the BV package. It is called from SlepcFinalize().
      30             : 
      31             :    Level: developer
      32             : 
      33             : .seealso: SlepcFinalize()
      34             : @*/
      35        1321 : PetscErrorCode BVFinalizePackage(void)
      36             : {
      37        1321 :   PetscFunctionBegin;
      38        1321 :   PetscCall(PetscFunctionListDestroy(&BVList));
      39        1321 :   PetscCallMPI(MPI_Op_free(&MPIU_TSQR));
      40        1321 :   PetscCallMPI(MPI_Op_free(&MPIU_LAPY2));
      41        1321 :   BVPackageInitialized = PETSC_FALSE;
      42        1321 :   BVRegisterAllCalled  = PETSC_FALSE;
      43        1321 :   PetscFunctionReturn(PETSC_SUCCESS);
      44             : }
      45             : 
      46             : /*@C
      47             :    BVInitializePackage - This function initializes everything in the BV package.
      48             :    It is called from PetscDLLibraryRegister() when using dynamic libraries, and
      49             :    on the first call to BVCreate() when using static libraries.
      50             : 
      51             :    Level: developer
      52             : 
      53             : .seealso: SlepcInitialize()
      54             : @*/
      55       21777 : PetscErrorCode BVInitializePackage(void)
      56             : {
      57       21777 :   char           logList[256];
      58       21777 :   PetscBool      opt,pkg;
      59       21777 :   PetscClassId   classids[1];
      60             : 
      61       21777 :   PetscFunctionBegin;
      62       21777 :   if (BVPackageInitialized) PetscFunctionReturn(PETSC_SUCCESS);
      63        1321 :   BVPackageInitialized = PETSC_TRUE;
      64             :   /* Register Classes */
      65        1321 :   PetscCall(PetscClassIdRegister("Basis Vectors",&BV_CLASSID));
      66             :   /* Register Constructors */
      67        1321 :   PetscCall(BVRegisterAll());
      68             :   /* Register Events */
      69        1321 :   PetscCall(PetscLogEventRegister("BVCreate",BV_CLASSID,&BV_Create));
      70        1321 :   PetscCall(PetscLogEventRegister("BVCopy",BV_CLASSID,&BV_Copy));
      71        1321 :   PetscCall(PetscLogEventRegister("BVMult",BV_CLASSID,&BV_Mult));
      72        1321 :   PetscCall(PetscLogEventRegister("BVMultVec",BV_CLASSID,&BV_MultVec));
      73        1321 :   PetscCall(PetscLogEventRegister("BVMultInPlace",BV_CLASSID,&BV_MultInPlace));
      74        1321 :   PetscCall(PetscLogEventRegister("BVDot",BV_CLASSID,&BV_Dot));
      75        1321 :   PetscCall(PetscLogEventRegister("BVDotVec",BV_CLASSID,&BV_DotVec));
      76        1321 :   PetscCall(PetscLogEventRegister("BVOrthogonalize",BV_CLASSID,&BV_Orthogonalize));
      77        1321 :   PetscCall(PetscLogEventRegister("BVOrthogonalizeV",BV_CLASSID,&BV_OrthogonalizeVec));
      78        1321 :   PetscCall(PetscLogEventRegister("BVScale",BV_CLASSID,&BV_Scale));
      79        1321 :   PetscCall(PetscLogEventRegister("BVNorm",BV_CLASSID,&BV_Norm));
      80        1321 :   PetscCall(PetscLogEventRegister("BVNormVec",BV_CLASSID,&BV_NormVec));
      81        1321 :   PetscCall(PetscLogEventRegister("BVNormalize",BV_CLASSID,&BV_Normalize));
      82        1321 :   PetscCall(PetscLogEventRegister("BVSetRandom",BV_CLASSID,&BV_SetRandom));
      83        1321 :   PetscCall(PetscLogEventRegister("BVMatMult",BV_CLASSID,&BV_MatMult));
      84        1321 :   PetscCall(PetscLogEventRegister("BVMatMultVec",BV_CLASSID,&BV_MatMultVec));
      85        1321 :   PetscCall(PetscLogEventRegister("BVMatProject",BV_CLASSID,&BV_MatProject));
      86        1321 :   PetscCall(PetscLogEventRegister("BVSVDAndRank",BV_CLASSID,&BV_SVDAndRank));
      87             :   /* MPI reduction operation used in BVOrthogonalize */
      88        1321 :   PetscCallMPI(MPI_Op_create(SlepcGivensPacked,PETSC_FALSE,&MPIU_TSQR));
      89        1321 :   PetscCallMPI(MPI_Op_create(SlepcPythag,PETSC_TRUE,&MPIU_LAPY2));
      90             :   /* Process Info */
      91        1321 :   classids[0] = BV_CLASSID;
      92        1321 :   PetscCall(PetscInfoProcessClass("bv",1,&classids[0]));
      93             :   /* Process summary exclusions */
      94        1321 :   PetscCall(PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt));
      95        1321 :   if (opt) {
      96          11 :     PetscCall(PetscStrInList("bv",logList,',',&pkg));
      97          11 :     if (pkg) PetscCall(PetscLogEventDeactivateClass(BV_CLASSID));
      98             :   }
      99             :   /* Register package finalizer */
     100        1321 :   PetscCall(PetscRegisterFinalize(BVFinalizePackage));
     101        1321 :   PetscFunctionReturn(PETSC_SUCCESS);
     102             : }
     103             : 
     104             : /*@C
     105             :    BVDestroy - Destroys BV context that was created with BVCreate().
     106             : 
     107             :    Collective
     108             : 
     109             :    Input Parameter:
     110             : .  bv - the basis vectors context
     111             : 
     112             :    Level: beginner
     113             : 
     114             : .seealso: BVCreate()
     115             : @*/
     116       62522 : PetscErrorCode BVDestroy(BV *bv)
     117             : {
     118       62522 :   PetscFunctionBegin;
     119       62522 :   if (!*bv) PetscFunctionReturn(PETSC_SUCCESS);
     120       15416 :   PetscValidHeaderSpecific(*bv,BV_CLASSID,1);
     121       15416 :   PetscCheck(!(*bv)->lsplit,PetscObjectComm((PetscObject)*bv),PETSC_ERR_ARG_WRONGSTATE,"Must call BVRestoreSplit before destroying the BV");
     122       15416 :   if (--((PetscObject)*bv)->refct > 0) { *bv = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
     123       15171 :   PetscTryTypeMethod(*bv,destroy);
     124       15171 :   PetscCall(PetscLayoutDestroy(&(*bv)->map));
     125       15171 :   PetscCall(PetscFree((*bv)->vtype));
     126       15171 :   PetscCall(MatDestroy(&(*bv)->matrix));
     127       15171 :   PetscCall(VecDestroy(&(*bv)->Bx));
     128       15171 :   PetscCall(VecDestroy(&(*bv)->buffer));
     129       15171 :   PetscCall(BVDestroy(&(*bv)->cached));
     130       15171 :   PetscCall(BVDestroy(&(*bv)->L));
     131       15171 :   PetscCall(BVDestroy(&(*bv)->R));
     132       15171 :   PetscCall(PetscFree((*bv)->work));
     133       15171 :   PetscCall(PetscFree2((*bv)->h,(*bv)->c));
     134       15171 :   PetscCall(VecDestroy(&(*bv)->omega));
     135       15171 :   PetscCall(MatDestroy(&(*bv)->Acreate));
     136       15171 :   PetscCall(MatDestroy(&(*bv)->Aget));
     137       15171 :   PetscCall(MatDestroy(&(*bv)->Abuffer));
     138       15171 :   PetscCall(PetscRandomDestroy(&(*bv)->rand));
     139       15171 :   PetscCall(PetscHeaderDestroy(bv));
     140       15171 :   PetscFunctionReturn(PETSC_SUCCESS);
     141             : }
     142             : 
     143             : /*@
     144             :    BVCreate - Creates a basis vectors context.
     145             : 
     146             :    Collective
     147             : 
     148             :    Input Parameter:
     149             : .  comm - MPI communicator
     150             : 
     151             :    Output Parameter:
     152             : .  newbv - location to put the basis vectors context
     153             : 
     154             :    Level: beginner
     155             : 
     156             : .seealso: BVSetUp(), BVDestroy(), BV
     157             : @*/
     158       15171 : PetscErrorCode BVCreate(MPI_Comm comm,BV *newbv)
     159             : {
     160       15171 :   BV             bv;
     161             : 
     162       15171 :   PetscFunctionBegin;
     163       15171 :   PetscAssertPointer(newbv,2);
     164       15171 :   *newbv = NULL;
     165       15171 :   PetscCall(BVInitializePackage());
     166       15171 :   PetscCall(SlepcHeaderCreate(bv,BV_CLASSID,"BV","Basis Vectors","BV",comm,BVDestroy,BVView));
     167             : 
     168       15171 :   bv->map          = NULL;
     169       15171 :   bv->vtype        = NULL;
     170       15171 :   bv->n            = -1;
     171       15171 :   bv->N            = -1;
     172       15171 :   bv->m            = 0;
     173       15171 :   bv->l            = 0;
     174       15171 :   bv->k            = 0;
     175       15171 :   bv->nc           = 0;
     176       15171 :   bv->ld           = 0;
     177       15171 :   bv->orthog_type  = BV_ORTHOG_CGS;
     178       15171 :   bv->orthog_ref   = BV_ORTHOG_REFINE_IFNEEDED;
     179       15171 :   bv->orthog_eta   = 0.7071;
     180       15171 :   bv->orthog_block = BV_ORTHOG_BLOCK_GS;
     181       15171 :   bv->matrix       = NULL;
     182       15171 :   bv->indef        = PETSC_FALSE;
     183       15171 :   bv->vmm          = BV_MATMULT_MAT;
     184       15171 :   bv->rrandom      = PETSC_FALSE;
     185       15171 :   bv->deftol       = 10*PETSC_MACHINE_EPSILON;
     186             : 
     187       15171 :   bv->buffer       = NULL;
     188       15171 :   bv->Abuffer      = NULL;
     189       15171 :   bv->Bx           = NULL;
     190       15171 :   bv->xid          = 0;
     191       15171 :   bv->xstate       = 0;
     192       15171 :   bv->cv[0]        = NULL;
     193       15171 :   bv->cv[1]        = NULL;
     194       15171 :   bv->ci[0]        = -1;
     195       15171 :   bv->ci[1]        = -1;
     196       15171 :   bv->st[0]        = -1;
     197       15171 :   bv->st[1]        = -1;
     198       15171 :   bv->id[0]        = 0;
     199       15171 :   bv->id[1]        = 0;
     200       15171 :   bv->h            = NULL;
     201       15171 :   bv->c            = NULL;
     202       15171 :   bv->omega        = NULL;
     203       15171 :   bv->defersfo     = PETSC_FALSE;
     204       15171 :   bv->cached       = NULL;
     205       15171 :   bv->bvstate      = 0;
     206       15171 :   bv->L            = NULL;
     207       15171 :   bv->R            = NULL;
     208       15171 :   bv->lstate       = 0;
     209       15171 :   bv->rstate       = 0;
     210       15171 :   bv->lsplit       = 0;
     211       15171 :   bv->issplit      = 0;
     212       15171 :   bv->splitparent  = NULL;
     213       15171 :   bv->rand         = NULL;
     214       15171 :   bv->Acreate      = NULL;
     215       15171 :   bv->Aget         = NULL;
     216       15171 :   bv->cuda         = PETSC_FALSE;
     217       15171 :   bv->sfocalled    = PETSC_FALSE;
     218       15171 :   bv->work         = NULL;
     219       15171 :   bv->lwork        = 0;
     220       15171 :   bv->data         = NULL;
     221             : 
     222       15171 :   *newbv = bv;
     223       15171 :   PetscFunctionReturn(PETSC_SUCCESS);
     224             : }
     225             : 
     226             : /*@
     227             :    BVCreateFromMat - Creates a basis vectors object from a dense Mat object.
     228             : 
     229             :    Collective
     230             : 
     231             :    Input Parameter:
     232             : .  A - a dense tall-skinny matrix
     233             : 
     234             :    Output Parameter:
     235             : .  bv - the new basis vectors context
     236             : 
     237             :    Notes:
     238             :    The matrix values are copied to the BV data storage, memory is not shared.
     239             : 
     240             :    The communicator of the BV object will be the same as A, and so will be
     241             :    the dimensions.
     242             : 
     243             :    Level: intermediate
     244             : 
     245             : .seealso: BVCreate(), BVDestroy(), BVCreateMat()
     246             : @*/
     247          81 : PetscErrorCode BVCreateFromMat(Mat A,BV *bv)
     248             : {
     249          81 :   PetscInt  n,N,k;
     250          81 :   VecType   vtype;
     251             : 
     252          81 :   PetscFunctionBegin;
     253          81 :   PetscValidHeaderSpecific(A,MAT_CLASSID,1);
     254             : 
     255          81 :   PetscCall(MatGetSize(A,&N,&k));
     256          81 :   PetscCall(MatGetLocalSize(A,&n,NULL));
     257          81 :   PetscCall(MatGetVecType(A,&vtype));
     258          81 :   PetscCall(BVCreate(PetscObjectComm((PetscObject)A),bv));
     259          81 :   PetscCall(BVSetSizes(*bv,n,N,k));
     260          81 :   PetscCall(BVSetVecType(*bv,vtype));
     261             : 
     262          81 :   (*bv)->Acreate = A;
     263          81 :   PetscCall(PetscObjectReference((PetscObject)A));
     264          81 :   PetscFunctionReturn(PETSC_SUCCESS);
     265             : }
     266             : 
     267             : /*@
     268             :    BVInsertVec - Insert a vector into the specified column.
     269             : 
     270             :    Logically Collective
     271             : 
     272             :    Input Parameters:
     273             : +  V - basis vectors
     274             : .  j - the column of V to be overwritten
     275             : -  w - the vector to be copied
     276             : 
     277             :    Level: intermediate
     278             : 
     279             : .seealso: BVInsertVecs()
     280             : @*/
     281       55452 : PetscErrorCode BVInsertVec(BV V,PetscInt j,Vec w)
     282             : {
     283       55452 :   PetscInt       n,N;
     284       55452 :   Vec            v;
     285             : 
     286       55452 :   PetscFunctionBegin;
     287       55452 :   PetscValidHeaderSpecific(V,BV_CLASSID,1);
     288      221808 :   PetscValidLogicalCollectiveInt(V,j,2);
     289       55452 :   PetscValidHeaderSpecific(w,VEC_CLASSID,3);
     290       55452 :   PetscValidType(V,1);
     291       55452 :   BVCheckSizes(V,1);
     292       55452 :   PetscCheckSameComm(V,1,w,3);
     293             : 
     294       55452 :   PetscCall(VecGetSize(w,&N));
     295       55452 :   PetscCall(VecGetLocalSize(w,&n));
     296       55452 :   PetscCheck(N==V->N && n==V->n,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_INCOMP,"Vec sizes (global %" PetscInt_FMT ", local %" PetscInt_FMT ") do not match BV sizes (global %" PetscInt_FMT ", local %" PetscInt_FMT ")",N,n,V->N,V->n);
     297       55452 :   PetscCheck(j>=-V->nc && j<V->m,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Argument j has wrong value %" PetscInt_FMT ", should be between %" PetscInt_FMT " and %" PetscInt_FMT,j,-V->nc,V->m-1);
     298             : 
     299       55452 :   PetscCall(BVGetColumn(V,j,&v));
     300       55452 :   PetscCall(VecCopy(w,v));
     301       55452 :   PetscCall(BVRestoreColumn(V,j,&v));
     302       55452 :   PetscCall(PetscObjectStateIncrease((PetscObject)V));
     303       55452 :   PetscFunctionReturn(PETSC_SUCCESS);
     304             : }
     305             : 
     306             : /*@
     307             :    BVInsertVecs - Insert a set of vectors into the specified columns.
     308             : 
     309             :    Collective
     310             : 
     311             :    Input Parameters:
     312             : +  V - basis vectors
     313             : .  s - first column of V to be overwritten
     314             : .  W - set of vectors to be copied
     315             : -  orth - flag indicating if the vectors must be orthogonalized
     316             : 
     317             :    Input/Output Parameter:
     318             : .  m - number of input vectors, on output the number of linearly independent
     319             :        vectors
     320             : 
     321             :    Notes:
     322             :    Copies the contents of vectors W to V(:,s:s+n). If the orthogonalization
     323             :    flag is set, then the vectors are copied one by one and then orthogonalized
     324             :    against the previous ones. If any of them is linearly dependent then it
     325             :    is discarded and the value of m is decreased.
     326             : 
     327             :    Level: intermediate
     328             : 
     329             : .seealso: BVInsertVec(), BVOrthogonalizeColumn()
     330             : @*/
     331         331 : PetscErrorCode BVInsertVecs(BV V,PetscInt s,PetscInt *m,Vec *W,PetscBool orth)
     332             : {
     333         331 :   PetscInt       n,N,i,ndep;
     334         331 :   PetscBool      lindep;
     335         331 :   PetscReal      norm;
     336         331 :   Vec            v;
     337             : 
     338         331 :   PetscFunctionBegin;
     339         331 :   PetscValidHeaderSpecific(V,BV_CLASSID,1);
     340        1324 :   PetscValidLogicalCollectiveInt(V,s,2);
     341         331 :   PetscAssertPointer(m,3);
     342        1324 :   PetscValidLogicalCollectiveInt(V,*m,3);
     343         331 :   if (!*m) PetscFunctionReturn(PETSC_SUCCESS);
     344         331 :   PetscCheck(*m>0,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %" PetscInt_FMT ") cannot be negative",*m);
     345         331 :   PetscAssertPointer(W,4);
     346         331 :   PetscValidHeaderSpecific(*W,VEC_CLASSID,4);
     347        1324 :   PetscValidLogicalCollectiveBool(V,orth,5);
     348         331 :   PetscValidType(V,1);
     349         331 :   BVCheckSizes(V,1);
     350         331 :   PetscCheckSameComm(V,1,*W,4);
     351             : 
     352         331 :   PetscCall(VecGetSize(*W,&N));
     353         331 :   PetscCall(VecGetLocalSize(*W,&n));
     354         331 :   PetscCheck(N==V->N && n==V->n,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_INCOMP,"Vec sizes (global %" PetscInt_FMT ", local %" PetscInt_FMT ") do not match BV sizes (global %" PetscInt_FMT ", local %" PetscInt_FMT ")",N,n,V->N,V->n);
     355         331 :   PetscCheck(s>=0 && s<V->m,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Argument s has wrong value %" PetscInt_FMT ", should be between 0 and %" PetscInt_FMT,s,V->m-1);
     356         331 :   PetscCheck(s+(*m)<=V->m,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Too many vectors provided, there is only room for %" PetscInt_FMT,V->m);
     357             : 
     358             :   ndep = 0;
     359         710 :   for (i=0;i<*m;i++) {
     360         379 :     PetscCall(BVGetColumn(V,s+i-ndep,&v));
     361         379 :     PetscCall(VecCopy(W[i],v));
     362         379 :     PetscCall(BVRestoreColumn(V,s+i-ndep,&v));
     363         379 :     if (orth) {
     364         379 :       PetscCall(BVOrthogonalizeColumn(V,s+i-ndep,NULL,&norm,&lindep));
     365         379 :       if (norm==0.0 || lindep) {
     366           0 :         PetscCall(PetscInfo(V,"Removing linearly dependent vector %" PetscInt_FMT "\n",i));
     367           0 :         ndep++;
     368         379 :       } else PetscCall(BVScaleColumn(V,s+i-ndep,1.0/norm));
     369             :     }
     370             :   }
     371         331 :   *m -= ndep;
     372         331 :   PetscCall(PetscObjectStateIncrease((PetscObject)V));
     373         331 :   PetscFunctionReturn(PETSC_SUCCESS);
     374             : }
     375             : 
     376             : /*@
     377             :    BVInsertConstraints - Insert a set of vectors as constraints.
     378             : 
     379             :    Collective
     380             : 
     381             :    Input Parameters:
     382             : +  V - basis vectors
     383             : -  C - set of vectors to be inserted as constraints
     384             : 
     385             :    Input/Output Parameter:
     386             : .  nc - number of input vectors, on output the number of linearly independent
     387             :        vectors
     388             : 
     389             :    Notes:
     390             :    The constraints are relevant only during orthogonalization. Constraint
     391             :    vectors span a subspace that is deflated in every orthogonalization
     392             :    operation, so they are intended for removing those directions from the
     393             :    orthogonal basis computed in regular BV columns.
     394             : 
     395             :    Constraints are not stored in regular BV columns, but in a special part of
     396             :    the storage. They can be accessed with negative indices in BVGetColumn().
     397             : 
     398             :    This operation is DESTRUCTIVE, meaning that all data contained in the
     399             :    columns of V is lost. This is typically invoked just after creating the BV.
     400             :    Once a set of constraints has been set, it is not allowed to call this
     401             :    function again.
     402             : 
     403             :    The vectors are copied one by one and then orthogonalized against the
     404             :    previous ones. If any of them is linearly dependent then it is discarded
     405             :    and the value of nc is decreased. The behaviour is similar to BVInsertVecs().
     406             : 
     407             :    Level: advanced
     408             : 
     409             : .seealso: BVInsertVecs(), BVOrthogonalizeColumn(), BVGetColumn(), BVGetNumConstraints()
     410             : @*/
     411          38 : PetscErrorCode BVInsertConstraints(BV V,PetscInt *nc,Vec *C)
     412             : {
     413          38 :   PetscInt       msave;
     414             : 
     415          38 :   PetscFunctionBegin;
     416          38 :   PetscValidHeaderSpecific(V,BV_CLASSID,1);
     417          38 :   PetscAssertPointer(nc,2);
     418         152 :   PetscValidLogicalCollectiveInt(V,*nc,2);
     419          38 :   if (!*nc) PetscFunctionReturn(PETSC_SUCCESS);
     420          38 :   PetscCheck(*nc>0,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of constraints (given %" PetscInt_FMT ") cannot be negative",*nc);
     421          38 :   PetscAssertPointer(C,3);
     422          38 :   PetscValidHeaderSpecific(*C,VEC_CLASSID,3);
     423          38 :   PetscValidType(V,1);
     424          38 :   BVCheckSizes(V,1);
     425          38 :   PetscCheckSameComm(V,1,*C,3);
     426          38 :   PetscCheck(!V->issplit,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Operation not permitted for a BV obtained from BVGetSplit");
     427          38 :   PetscCheck(!V->nc,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Constraints already present in this BV object");
     428          38 :   PetscCheck(V->ci[0]==-1 && V->ci[1]==-1,PetscObjectComm((PetscObject)V),PETSC_ERR_SUP,"Cannot call BVInsertConstraints after BVGetColumn");
     429             : 
     430          38 :   msave = V->m;
     431          38 :   PetscCall(BVResize(V,*nc+V->m,PETSC_FALSE));
     432          38 :   PetscCall(BVInsertVecs(V,0,nc,C,PETSC_TRUE));
     433          38 :   V->nc = *nc;
     434          38 :   V->m  = msave;
     435          38 :   V->ci[0] = -V->nc-1;
     436          38 :   V->ci[1] = -V->nc-1;
     437          38 :   PetscCall(PetscObjectStateIncrease((PetscObject)V));
     438          38 :   PetscFunctionReturn(PETSC_SUCCESS);
     439             : }
     440             : 
     441             : /*@C
     442             :    BVSetOptionsPrefix - Sets the prefix used for searching for all
     443             :    BV options in the database.
     444             : 
     445             :    Logically Collective
     446             : 
     447             :    Input Parameters:
     448             : +  bv     - the basis vectors context
     449             : -  prefix - the prefix string to prepend to all BV option requests
     450             : 
     451             :    Notes:
     452             :    A hyphen (-) must NOT be given at the beginning of the prefix name.
     453             :    The first character of all runtime options is AUTOMATICALLY the
     454             :    hyphen.
     455             : 
     456             :    Level: advanced
     457             : 
     458             : .seealso: BVAppendOptionsPrefix(), BVGetOptionsPrefix()
     459             : @*/
     460         222 : PetscErrorCode BVSetOptionsPrefix(BV bv,const char *prefix)
     461             : {
     462         222 :   PetscFunctionBegin;
     463         222 :   PetscValidHeaderSpecific(bv,BV_CLASSID,1);
     464         222 :   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)bv,prefix));
     465         222 :   PetscFunctionReturn(PETSC_SUCCESS);
     466             : }
     467             : 
     468             : /*@C
     469             :    BVAppendOptionsPrefix - Appends to the prefix used for searching for all
     470             :    BV options in the database.
     471             : 
     472             :    Logically Collective
     473             : 
     474             :    Input Parameters:
     475             : +  bv     - the basis vectors context
     476             : -  prefix - the prefix string to prepend to all BV option requests
     477             : 
     478             :    Notes:
     479             :    A hyphen (-) must NOT be given at the beginning of the prefix name.
     480             :    The first character of all runtime options is AUTOMATICALLY the
     481             :    hyphen.
     482             : 
     483             :    Level: advanced
     484             : 
     485             : .seealso: BVSetOptionsPrefix(), BVGetOptionsPrefix()
     486             : @*/
     487         186 : PetscErrorCode BVAppendOptionsPrefix(BV bv,const char *prefix)
     488             : {
     489         186 :   PetscFunctionBegin;
     490         186 :   PetscValidHeaderSpecific(bv,BV_CLASSID,1);
     491         186 :   PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)bv,prefix));
     492         186 :   PetscFunctionReturn(PETSC_SUCCESS);
     493             : }
     494             : 
     495             : /*@C
     496             :    BVGetOptionsPrefix - Gets the prefix used for searching for all
     497             :    BV options in the database.
     498             : 
     499             :    Not Collective
     500             : 
     501             :    Input Parameters:
     502             : .  bv - the basis vectors context
     503             : 
     504             :    Output Parameters:
     505             : .  prefix - pointer to the prefix string used, is returned
     506             : 
     507             :    Note:
     508             :    On the Fortran side, the user should pass in a string 'prefix' of
     509             :    sufficient length to hold the prefix.
     510             : 
     511             :    Level: advanced
     512             : 
     513             : .seealso: BVSetOptionsPrefix(), BVAppendOptionsPrefix()
     514             : @*/
     515           0 : PetscErrorCode BVGetOptionsPrefix(BV bv,const char *prefix[])
     516             : {
     517           0 :   PetscFunctionBegin;
     518           0 :   PetscValidHeaderSpecific(bv,BV_CLASSID,1);
     519           0 :   PetscAssertPointer(prefix,2);
     520           0 :   PetscCall(PetscObjectGetOptionsPrefix((PetscObject)bv,prefix));
     521           0 :   PetscFunctionReturn(PETSC_SUCCESS);
     522             : }
     523             : 
     524             : /*@C
     525             :    BVView - Prints the BV data structure.
     526             : 
     527             :    Collective
     528             : 
     529             :    Input Parameters:
     530             : +  bv     - the BV context
     531             : -  viewer - optional visualization context
     532             : 
     533             :    Note:
     534             :    The available visualization contexts include
     535             : +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
     536             : -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
     537             :          output where only the first processor opens
     538             :          the file.  All other processors send their
     539             :          data to the first processor to print.
     540             : 
     541             :    The user can open an alternative visualization contexts with
     542             :    PetscViewerASCIIOpen() (output to a specified file).
     543             : 
     544             :    Level: beginner
     545             : 
     546             : .seealso: BVCreate()
     547             : @*/
     548          75 : PetscErrorCode BVView(BV bv,PetscViewer viewer)
     549             : {
     550          75 :   PetscBool         isascii;
     551          75 :   PetscViewerFormat format;
     552          75 :   const char        *orthname[2] = {"classical","modified"};
     553          75 :   const char        *refname[3] = {"if needed","never","always"};
     554             : 
     555          75 :   PetscFunctionBegin;
     556          75 :   PetscValidHeaderSpecific(bv,BV_CLASSID,1);
     557          75 :   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)bv),&viewer));
     558          75 :   PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2);
     559             : 
     560          75 :   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
     561          75 :   if (isascii) {
     562          75 :     PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)bv,viewer));
     563          75 :     PetscCall(PetscViewerGetFormat(viewer,&format));
     564          75 :     if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
     565          70 :       PetscCall(PetscViewerASCIIPrintf(viewer,"  %" PetscInt_FMT " columns of global length %" PetscInt_FMT "%s\n",bv->m,bv->N,bv->cuda?" (CUDA)":""));
     566          35 :       if (bv->nc>0) PetscCall(PetscViewerASCIIPrintf(viewer,"  number of constraints: %" PetscInt_FMT "\n",bv->nc));
     567          35 :       PetscCall(PetscViewerASCIIPrintf(viewer,"  vector orthogonalization method: %s Gram-Schmidt\n",orthname[bv->orthog_type]));
     568          35 :       switch (bv->orthog_ref) {
     569          35 :         case BV_ORTHOG_REFINE_IFNEEDED:
     570          35 :           PetscCall(PetscViewerASCIIPrintf(viewer,"  orthogonalization refinement: %s (eta: %g)\n",refname[bv->orthog_ref],(double)bv->orthog_eta));
     571             :           break;
     572           0 :         case BV_ORTHOG_REFINE_NEVER:
     573             :         case BV_ORTHOG_REFINE_ALWAYS:
     574           0 :           PetscCall(PetscViewerASCIIPrintf(viewer,"  orthogonalization refinement: %s\n",refname[bv->orthog_ref]));
     575             :           break;
     576             :       }
     577          35 :       PetscCall(PetscViewerASCIIPrintf(viewer,"  block orthogonalization method: %s\n",BVOrthogBlockTypes[bv->orthog_block]));
     578          35 :       if (bv->matrix) {
     579           0 :         if (bv->indef) PetscCall(PetscViewerASCIIPrintf(viewer,"  indefinite inner product\n"));
     580           0 :         else PetscCall(PetscViewerASCIIPrintf(viewer,"  non-standard inner product\n"));
     581           0 :         PetscCall(PetscViewerASCIIPrintf(viewer,"  tolerance for definite inner product: %g\n",(double)bv->deftol));
     582           0 :         PetscCall(PetscViewerASCIIPrintf(viewer,"  inner product matrix:\n"));
     583           0 :         PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
     584           0 :         PetscCall(PetscViewerASCIIPushTab(viewer));
     585           0 :         PetscCall(MatView(bv->matrix,viewer));
     586           0 :         PetscCall(PetscViewerASCIIPopTab(viewer));
     587           0 :         PetscCall(PetscViewerPopFormat(viewer));
     588             :       }
     589          35 :       switch (bv->vmm) {
     590           4 :         case BV_MATMULT_VECS:
     591           4 :           PetscCall(PetscViewerASCIIPrintf(viewer,"  doing matmult as matrix-vector products\n"));
     592             :           break;
     593          31 :         case BV_MATMULT_MAT:
     594          31 :           PetscCall(PetscViewerASCIIPrintf(viewer,"  doing matmult as a single matrix-matrix product\n"));
     595             :           break;
     596           0 :         case BV_MATMULT_MAT_SAVE:
     597           0 :           PetscCall(PetscViewerASCIIPrintf(viewer,"  mat_save is deprecated, use mat\n"));
     598             :           break;
     599             :       }
     600          35 :       if (bv->rrandom) PetscCall(PetscViewerASCIIPrintf(viewer,"  generating random vectors independent of the number of processes\n"));
     601             :     }
     602             :   }
     603          75 :   PetscTryTypeMethod(bv,view,viewer);
     604          75 :   PetscFunctionReturn(PETSC_SUCCESS);
     605             : }
     606             : 
     607             : /*@C
     608             :    BVViewFromOptions - View from options
     609             : 
     610             :    Collective
     611             : 
     612             :    Input Parameters:
     613             : +  bv   - the basis vectors context
     614             : .  obj  - optional object
     615             : -  name - command line option
     616             : 
     617             :    Level: intermediate
     618             : 
     619             : .seealso: BVView(), BVCreate()
     620             : @*/
     621           0 : PetscErrorCode BVViewFromOptions(BV bv,PetscObject obj,const char name[])
     622             : {
     623           0 :   PetscFunctionBegin;
     624           0 :   PetscValidHeaderSpecific(bv,BV_CLASSID,1);
     625           0 :   PetscCall(PetscObjectViewFromOptions((PetscObject)bv,obj,name));
     626           0 :   PetscFunctionReturn(PETSC_SUCCESS);
     627             : }
     628             : 
     629             : /*@C
     630             :    BVRegister - Adds a new storage format to the BV package.
     631             : 
     632             :    Not Collective
     633             : 
     634             :    Input Parameters:
     635             : +  name     - name of a new user-defined BV
     636             : -  function - routine to create context
     637             : 
     638             :    Notes:
     639             :    BVRegister() may be called multiple times to add several user-defined
     640             :    basis vectors.
     641             : 
     642             :    Level: advanced
     643             : 
     644             : .seealso: BVRegisterAll()
     645             : @*/
     646        6605 : PetscErrorCode BVRegister(const char *name,PetscErrorCode (*function)(BV))
     647             : {
     648        6605 :   PetscFunctionBegin;
     649        6605 :   PetscCall(BVInitializePackage());
     650        6605 :   PetscCall(PetscFunctionListAdd(&BVList,name,function));
     651        6605 :   PetscFunctionReturn(PETSC_SUCCESS);
     652             : }
     653             : 
     654       70844 : PetscErrorCode BVAllocateWork_Private(BV bv,PetscInt s)
     655             : {
     656       70844 :   PetscFunctionBegin;
     657       70844 :   if (s>bv->lwork) {
     658       17530 :     PetscCall(PetscFree(bv->work));
     659       17530 :     PetscCall(PetscMalloc1(s,&bv->work));
     660       17530 :     bv->lwork = s;
     661             :   }
     662       70844 :   PetscFunctionReturn(PETSC_SUCCESS);
     663             : }
     664             : 
     665             : #if defined(PETSC_USE_DEBUG) && !defined(PETSC_CLANG_STATIC_ANALYZER)
     666             : /*
     667             :    SlepcDebugBVView - partially view a BV object, to be used from within a debugger.
     668             : 
     669             :      ini, end: columns to be viewed
     670             :      s: name of Matlab variable
     671             :      filename: optionally write output to a file
     672             :  */
     673           0 : PETSC_UNUSED PetscErrorCode SlepcDebugBVView(BV bv,PetscInt ini,PetscInt end,const char *s,const char *filename)
     674             : {
     675           0 :   PetscInt       N,m;
     676           0 :   PetscScalar    *array;
     677             : 
     678           0 :   PetscFunctionBegin;
     679           0 :   PetscCall(BVGetArray(bv,&array));
     680           0 :   PetscCall(BVGetSizes(bv,NULL,&N,&m));
     681           0 :   PetscCall(SlepcDebugViewMatrix(N,end-ini+1,array+ini*N,NULL,bv->ld,s,filename));
     682           0 :   PetscCall(BVRestoreArray(bv,&array));
     683           0 :   PetscFunctionReturn(PETSC_SUCCESS);
     684             : }
     685             : #endif

Generated by: LCOV version 1.14