Actual source code: svdsetup.c

slepc-3.19.2 2023-09-05
Report Typos and Errors
  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:    SVD routines for setting up the solver
 12: */

 14: #include <slepc/private/svdimpl.h>

 16: /*@
 17:    SVDSetOperators - Set the matrices associated with the singular value problem.

 19:    Collective

 21:    Input Parameters:
 22: +  svd - the singular value solver context
 23: .  A   - the matrix associated with the singular value problem
 24: -  B   - the second matrix in the case of GSVD

 26:    Level: beginner

 28: .seealso: SVDSolve(), SVDGetOperators()
 29: @*/
 30: PetscErrorCode SVDSetOperators(SVD svd,Mat A,Mat B)
 31: {
 32:   PetscInt       Ma,Na,Mb,Nb,ma,na,mb,nb,M0,N0,m0,n0;
 33:   PetscBool      samesize=PETSC_TRUE;

 35:   PetscFunctionBegin;
 39:   PetscCheckSameComm(svd,1,A,2);
 40:   if (B) PetscCheckSameComm(svd,1,B,3);

 42:   /* Check matrix sizes */
 43:   PetscCall(MatGetSize(A,&Ma,&Na));
 44:   PetscCall(MatGetLocalSize(A,&ma,&na));
 45:   if (svd->OP) {
 46:     PetscCall(MatGetSize(svd->OP,&M0,&N0));
 47:     PetscCall(MatGetLocalSize(svd->OP,&m0,&n0));
 48:     if (M0!=Ma || N0!=Na || m0!=ma || n0!=na) samesize = PETSC_FALSE;
 49:   }
 50:   if (B) {
 51:     PetscCall(MatGetSize(B,&Mb,&Nb));
 52:     PetscCall(MatGetLocalSize(B,&mb,&nb));
 53:     PetscCheck(Na==Nb,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_WRONG,"Different number of columns in A (%" PetscInt_FMT ") and B (%" PetscInt_FMT ")",Na,Nb);
 54:     PetscCheck(na==nb,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_WRONG,"Different local column size in A (%" PetscInt_FMT ") and B (%" PetscInt_FMT ")",na,nb);
 55:     if (svd->OPb) {
 56:       PetscCall(MatGetSize(svd->OPb,&M0,&N0));
 57:       PetscCall(MatGetLocalSize(svd->OPb,&m0,&n0));
 58:       if (M0!=Mb || N0!=Nb || m0!=mb || n0!=nb) samesize = PETSC_FALSE;
 59:     }
 60:   }

 62:   PetscCall(PetscObjectReference((PetscObject)A));
 63:   if (B) PetscCall(PetscObjectReference((PetscObject)B));
 64:   if (svd->state && !samesize) PetscCall(SVDReset(svd));
 65:   else {
 66:     PetscCall(MatDestroy(&svd->OP));
 67:     PetscCall(MatDestroy(&svd->OPb));
 68:     PetscCall(MatDestroy(&svd->A));
 69:     PetscCall(MatDestroy(&svd->B));
 70:     PetscCall(MatDestroy(&svd->AT));
 71:     PetscCall(MatDestroy(&svd->BT));
 72:   }
 73:   svd->nrma = 0.0;
 74:   svd->nrmb = 0.0;
 75:   svd->OP   = A;
 76:   svd->OPb  = B;
 77:   svd->state = SVD_STATE_INITIAL;
 78:   PetscFunctionReturn(PETSC_SUCCESS);
 79: }

 81: /*@
 82:    SVDGetOperators - Get the matrices associated with the singular value problem.

 84:    Not Collective

 86:    Input Parameter:
 87: .  svd - the singular value solver context

 89:    Output Parameters:
 90: +  A  - the matrix associated with the singular value problem
 91: -  B  - the second matrix in the case of GSVD

 93:    Level: intermediate

 95: .seealso: SVDSolve(), SVDSetOperators()
 96: @*/
 97: PetscErrorCode SVDGetOperators(SVD svd,Mat *A,Mat *B)
 98: {
 99:   PetscFunctionBegin;
101:   if (A) *A = svd->OP;
102:   if (B) *B = svd->OPb;
103:   PetscFunctionReturn(PETSC_SUCCESS);
104: }

106: /*@
107:    SVDSetSignature - Set the signature matrix defining a hyperbolic singular value problem.

109:    Collective

111:    Input Parameters:
112: +  svd   - the singular value solver context
113: -  omega - a vector containing the diagonal elements of the signature matrix (or NULL)

115:    Notes:
116:    The signature matrix is relevant only for hyperbolic problems (HSVD).
117:    Use NULL to reset a previously set signature.

119:    Level: intermediate

121: .seealso: SVDSetProblemType(), SVDSetOperators(), SVDGetSignature()
122: @*/
123: PetscErrorCode SVDSetSignature(SVD svd,Vec omega)
124: {
125:   PetscInt N,Ma,n,ma;

127:   PetscFunctionBegin;
129:   if (omega) {
131:     PetscCheckSameComm(svd,1,omega,2);
132:   }

134:   if (omega && svd->OP) {  /* Check sizes */
135:     PetscCall(VecGetSize(omega,&N));
136:     PetscCall(VecGetLocalSize(omega,&n));
137:     PetscCall(MatGetSize(svd->OP,&Ma,NULL));
138:     PetscCall(MatGetLocalSize(svd->OP,&ma,NULL));
139:     PetscCheck(N==Ma,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_WRONG,"Global size of signature (%" PetscInt_FMT ") does not match the row size of A (%" PetscInt_FMT ")",N,Ma);
140:     PetscCheck(n==ma,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_WRONG,"Local size of signature (%" PetscInt_FMT ") does not match the local row size of A (%" PetscInt_FMT ")",n,ma);
141:   }

143:   if (omega) PetscCall(PetscObjectReference((PetscObject)omega));
144:   PetscCall(VecDestroy(&svd->omega));
145:   svd->omega = omega;
146:   svd->state = SVD_STATE_INITIAL;
147:   PetscFunctionReturn(PETSC_SUCCESS);
148: }

150: /*@
151:    SVDGetSignature - Get the signature matrix defining a hyperbolic singular value problem.

153:    Not Collective

155:    Input Parameter:
156: .  svd - the singular value solver context

158:    Output Parameter:
159: .  omega - a vector containing the diagonal elements of the signature matrix

161:    Level: intermediate

163: .seealso: SVDSetSignature()
164: @*/
165: PetscErrorCode SVDGetSignature(SVD svd,Vec *omega)
166: {
167:   PetscFunctionBegin;
170:   *omega = svd->omega;
171:   PetscFunctionReturn(PETSC_SUCCESS);
172: }

174: /*@
175:    SVDSetDSType - Sets the type of the internal DS object based on the current
176:    settings of the singular value solver.

178:    Collective

180:    Input Parameter:
181: .  svd - singular value solver context

183:    Note:
184:    This function need not be called explicitly, since it will be called at
185:    both SVDSetFromOptions() and SVDSetUp().

187:    Level: developer

189: .seealso: SVDSetFromOptions(), SVDSetUp()
190: @*/
191: PetscErrorCode SVDSetDSType(SVD svd)
192: {
193:   PetscFunctionBegin;
195:   PetscTryTypeMethod(svd,setdstype);
196:   PetscFunctionReturn(PETSC_SUCCESS);
197: }

199: /*@
200:    SVDSetUp - Sets up all the internal data structures necessary for the
201:    execution of the singular value solver.

203:    Collective

205:    Input Parameter:
206: .  svd   - singular value solver context

208:    Notes:
209:    This function need not be called explicitly in most cases, since SVDSolve()
210:    calls it. It can be useful when one wants to measure the set-up time
211:    separately from the solve time.

213:    Level: developer

215: .seealso: SVDCreate(), SVDSolve(), SVDDestroy()
216: @*/
217: PetscErrorCode SVDSetUp(SVD svd)
218: {
219:   PetscBool      flg;
220:   PetscInt       M,N,P=0,k,maxnsol;
221:   SlepcSC        sc;
222:   Vec            *T;
223:   BV             bv;

225:   PetscFunctionBegin;
227:   if (svd->state) PetscFunctionReturn(PETSC_SUCCESS);
228:   PetscCall(PetscLogEventBegin(SVD_SetUp,svd,0,0,0));

230:   /* reset the convergence flag from the previous solves */
231:   svd->reason = SVD_CONVERGED_ITERATING;

233:   /* set default solver type (SVDSetFromOptions was not called) */
234:   if (!((PetscObject)svd)->type_name) PetscCall(SVDSetType(svd,SVDCROSS));
235:   if (!svd->ds) PetscCall(SVDGetDS(svd,&svd->ds));
236:   PetscCall(SVDSetDSType(svd));

238:   /* check matrices */
239:   PetscCheck(svd->OP,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_WRONGSTATE,"SVDSetOperators() must be called first");

241:   if (!svd->problem_type) {  /* set default problem type */
242:     if (svd->OPb) {
243:       PetscCheck(!svd->omega,PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"There is no support yet for generalized hyperbolic problems");
244:       PetscCall(SVDSetProblemType(svd,SVD_GENERALIZED));
245:     } else {
246:       if (svd->omega) PetscCall(SVDSetProblemType(svd,SVD_HYPERBOLIC));
247:       else PetscCall(SVDSetProblemType(svd,SVD_STANDARD));
248:     }
249:   } else {  /* check consistency of problem type set by user */
250:     if (svd->OPb) {
251:       PetscCheck(svd->isgeneralized,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_INCOMP,"Inconsistent SVD state: the problem type does not match the number of matrices");
252:       PetscCheck(!svd->omega,PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"There is no support yet for generalized hyperbolic problems");
253:     } else {
254:       PetscCheck(!svd->isgeneralized,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_INCOMP,"Inconsistent SVD state: the problem type does not match the number of matrices");
255:       if (svd->omega) PetscCheck(svd->ishyperbolic,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_INCOMP,"Inconsistent SVD state: the problem type must be set to hyperbolic when passing a signature with SVDSetSignature()");
256:       else PetscCheck(!svd->ishyperbolic,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_INCOMP,"Inconsistent SVD state: a hyperbolic problem requires passing a signature with SVDSetSignature()");
257:     }
258:   }

260:   /* determine how to handle the transpose */
261:   svd->expltrans = PETSC_TRUE;
262:   if (svd->impltrans) svd->expltrans = PETSC_FALSE;
263:   else {
264:     PetscCall(MatHasOperation(svd->OP,MATOP_TRANSPOSE,&flg));
265:     if (!flg) svd->expltrans = PETSC_FALSE;
266:     else {
267:       PetscCall(PetscObjectTypeCompareAny((PetscObject)svd,&flg,SVDLAPACK,SVDSCALAPACK,SVDKSVD,SVDELEMENTAL,""));
268:       if (flg) svd->expltrans = PETSC_FALSE;
269:     }
270:   }

272:   /* get matrix dimensions */
273:   PetscCall(MatGetSize(svd->OP,&M,&N));
274:   if (svd->isgeneralized) {
275:     PetscCall(MatGetSize(svd->OPb,&P,NULL));
276:     PetscCheck(M+P>=N,PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"The case when [A;B] has less rows than columns is not supported");
277:   }

279:   /* build transpose matrix */
280:   PetscCall(MatDestroy(&svd->A));
281:   PetscCall(MatDestroy(&svd->AT));
282:   PetscCall(PetscObjectReference((PetscObject)svd->OP));
283:   if (svd->expltrans) {
284:     if (svd->isgeneralized || M>=N) {
285:       svd->A = svd->OP;
286:       PetscCall(MatHermitianTranspose(svd->OP,MAT_INITIAL_MATRIX,&svd->AT));
287:     } else {
288:       PetscCall(MatHermitianTranspose(svd->OP,MAT_INITIAL_MATRIX,&svd->A));
289:       svd->AT = svd->OP;
290:     }
291:   } else {
292:     if (svd->isgeneralized || M>=N) {
293:       svd->A = svd->OP;
294:       PetscCall(MatCreateHermitianTranspose(svd->OP,&svd->AT));
295:     } else {
296:       PetscCall(MatCreateHermitianTranspose(svd->OP,&svd->A));
297:       svd->AT = svd->OP;
298:     }
299:   }

301:   /* build transpose matrix B for GSVD */
302:   if (svd->isgeneralized) {
303:     PetscCall(MatDestroy(&svd->B));
304:     PetscCall(MatDestroy(&svd->BT));
305:     PetscCall(PetscObjectReference((PetscObject)svd->OPb));
306:     if (svd->expltrans) {
307:       svd->B = svd->OPb;
308:       PetscCall(MatHermitianTranspose(svd->OPb,MAT_INITIAL_MATRIX,&svd->BT));
309:     } else {
310:       svd->B = svd->OPb;
311:       PetscCall(MatCreateHermitianTranspose(svd->OPb,&svd->BT));
312:     }
313:   }

315:   if (!svd->isgeneralized && M<N) {
316:     /* swap initial vectors */
317:     if (svd->nini || svd->ninil) {
318:       T=svd->ISL; svd->ISL=svd->IS; svd->IS=T;
319:       k=svd->ninil; svd->ninil=svd->nini; svd->nini=k;
320:     }
321:     /* swap basis vectors */
322:     if (!svd->swapped) {  /* only the first time in case of multiple calls */
323:       bv=svd->V; svd->V=svd->U; svd->U=bv;
324:       svd->swapped = PETSC_TRUE;
325:     }
326:   }

328:   maxnsol = svd->isgeneralized? PetscMin(PetscMin(M,N),P): PetscMin(M,N);
329:   svd->ncv = PetscMin(svd->ncv,maxnsol);
330:   svd->nsv = PetscMin(svd->nsv,maxnsol);
331:   PetscCheck(svd->ncv==PETSC_DEFAULT || svd->nsv<=svd->ncv,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nsv bigger than ncv");

333:   /* relative convergence criterion is not allowed in GSVD */
334:   if (svd->conv==(SVDConv)-1) PetscCall(SVDSetConvergenceTest(svd,svd->isgeneralized?SVD_CONV_NORM:SVD_CONV_REL));
335:   PetscCheck(!svd->isgeneralized || svd->conv!=SVD_CONV_REL,PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"Relative convergence criterion is not allowed in GSVD");

337:   /* initialization of matrix norm (stardard case only, for GSVD it is done inside setup()) */
338:   if (!svd->isgeneralized && svd->conv==SVD_CONV_NORM && !svd->nrma) PetscCall(MatNorm(svd->OP,NORM_INFINITY,&svd->nrma));

340:   /* call specific solver setup */
341:   PetscUseTypeMethod(svd,setup);

343:   /* set tolerance if not yet set */
344:   if (svd->tol==(PetscReal)PETSC_DEFAULT) svd->tol = SLEPC_DEFAULT_TOL;

346:   /* fill sorting criterion context */
347:   PetscCall(DSGetSlepcSC(svd->ds,&sc));
348:   sc->comparison    = (svd->which==SVD_LARGEST)? SlepcCompareLargestReal: SlepcCompareSmallestReal;
349:   sc->comparisonctx = NULL;
350:   sc->map           = NULL;
351:   sc->mapobj        = NULL;

353:   /* process initial vectors */
354:   if (svd->nini<0) {
355:     k = -svd->nini;
356:     PetscCheck(k<=svd->ncv,PetscObjectComm((PetscObject)svd),PETSC_ERR_USER_INPUT,"The number of initial vectors is larger than ncv");
357:     PetscCall(BVInsertVecs(svd->V,0,&k,svd->IS,PETSC_TRUE));
358:     PetscCall(SlepcBasisDestroy_Private(&svd->nini,&svd->IS));
359:     svd->nini = k;
360:   }
361:   if (svd->ninil<0) {
362:     k = 0;
363:     if (svd->leftbasis) {
364:       k = -svd->ninil;
365:       PetscCheck(k<=svd->ncv,PetscObjectComm((PetscObject)svd),PETSC_ERR_USER_INPUT,"The number of left initial vectors is larger than ncv");
366:       PetscCall(BVInsertVecs(svd->U,0,&k,svd->ISL,PETSC_TRUE));
367:     } else PetscCall(PetscInfo(svd,"Ignoring initial left vectors\n"));
368:     PetscCall(SlepcBasisDestroy_Private(&svd->ninil,&svd->ISL));
369:     svd->ninil = k;
370:   }

372:   PetscCall(PetscLogEventEnd(SVD_SetUp,svd,0,0,0));
373:   svd->state = SVD_STATE_SETUP;
374:   PetscFunctionReturn(PETSC_SUCCESS);
375: }

377: /*@C
378:    SVDSetInitialSpaces - Specify two basis of vectors that constitute the initial
379:    right and/or left spaces.

381:    Collective

383:    Input Parameters:
384: +  svd   - the singular value solver context
385: .  nr    - number of right vectors
386: .  isr   - set of basis vectors of the right initial space
387: .  nl    - number of left vectors
388: -  isl   - set of basis vectors of the left initial space

390:    Notes:
391:    The initial right and left spaces are rough approximations to the right and/or
392:    left singular subspaces from which the solver starts to iterate.
393:    It is not necessary to provide both sets of vectors.

395:    Some solvers start to iterate on a single vector (initial vector). In that case,
396:    the other vectors are ignored.

398:    These vectors do not persist from one SVDSolve() call to the other, so the
399:    initial space should be set every time.

401:    The vectors do not need to be mutually orthonormal, since they are explicitly
402:    orthonormalized internally.

404:    Common usage of this function is when the user can provide a rough approximation
405:    of the wanted singular space. Then, convergence may be faster.

407:    Level: intermediate

409: .seealso: SVDSetUp()
410: @*/
411: PetscErrorCode SVDSetInitialSpaces(SVD svd,PetscInt nr,Vec isr[],PetscInt nl,Vec isl[])
412: {
413:   PetscFunctionBegin;
417:   PetscCheck(nr>=0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Argument nr cannot be negative");
418:   if (nr>0) {
421:   }
422:   PetscCheck(nl>=0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Argument nl cannot be negative");
423:   if (nl>0) {
426:   }
427:   PetscCall(SlepcBasisReference_Private(nr,isr,&svd->nini,&svd->IS));
428:   PetscCall(SlepcBasisReference_Private(nl,isl,&svd->ninil,&svd->ISL));
429:   if (nr>0 || nl>0) svd->state = SVD_STATE_INITIAL;
430:   PetscFunctionReturn(PETSC_SUCCESS);
431: }

433: /*
434:   SVDSetDimensions_Default - Set reasonable values for ncv, mpd if not set
435:   by the user. This is called at setup.
436:  */
437: PetscErrorCode SVDSetDimensions_Default(SVD svd)
438: {
439:   PetscInt       N,M,P,maxnsol;

441:   PetscFunctionBegin;
442:   PetscCall(MatGetSize(svd->OP,&M,&N));
443:   maxnsol = PetscMin(M,N);
444:   if (svd->isgeneralized) {
445:     PetscCall(MatGetSize(svd->OPb,&P,NULL));
446:     maxnsol = PetscMin(maxnsol,P);
447:   }
448:   if (svd->ncv!=PETSC_DEFAULT) { /* ncv set */
449:     PetscCheck(svd->ncv>=svd->nsv,PetscObjectComm((PetscObject)svd),PETSC_ERR_USER_INPUT,"The value of ncv must be at least nsv");
450:   } else if (svd->mpd!=PETSC_DEFAULT) { /* mpd set */
451:     svd->ncv = PetscMin(maxnsol,svd->nsv+svd->mpd);
452:   } else { /* neither set: defaults depend on nsv being small or large */
453:     if (svd->nsv<500) svd->ncv = PetscMin(maxnsol,PetscMax(2*svd->nsv,10));
454:     else {
455:       svd->mpd = 500;
456:       svd->ncv = PetscMin(maxnsol,svd->nsv+svd->mpd);
457:     }
458:   }
459:   if (svd->mpd==PETSC_DEFAULT) svd->mpd = svd->ncv;
460:   PetscFunctionReturn(PETSC_SUCCESS);
461: }

463: /*@
464:    SVDAllocateSolution - Allocate memory storage for common variables such
465:    as the singular values and the basis vectors.

467:    Collective

469:    Input Parameters:
470: +  svd   - eigensolver context
471: -  extra - number of additional positions, used for methods that require a
472:            working basis slightly larger than ncv

474:    Developer Notes:
475:    This is SLEPC_EXTERN because it may be required by user plugin SVD
476:    implementations.

478:    This is called at setup after setting the value of ncv and the flag leftbasis.

480:    Level: developer

482: .seealso: SVDSetUp()
483: @*/
484: PetscErrorCode SVDAllocateSolution(SVD svd,PetscInt extra)
485: {
486:   PetscInt       oldsize,requested;
487:   Vec            tr,tl;

489:   PetscFunctionBegin;
490:   requested = svd->ncv + extra;

492:   /* oldsize is zero if this is the first time setup is called */
493:   PetscCall(BVGetSizes(svd->V,NULL,NULL,&oldsize));

495:   /* allocate sigma */
496:   if (requested != oldsize || !svd->sigma) {
497:     PetscCall(PetscFree3(svd->sigma,svd->perm,svd->errest));
498:     if (svd->sign) PetscCall(PetscFree(svd->sign));
499:     PetscCall(PetscMalloc3(requested,&svd->sigma,requested,&svd->perm,requested,&svd->errest));
500:     if (svd->ishyperbolic) PetscCall(PetscMalloc1(requested,&svd->sign));
501:   }
502:   /* allocate V */
503:   if (!svd->V) PetscCall(SVDGetBV(svd,&svd->V,NULL));
504:   if (!oldsize) {
505:     if (!((PetscObject)(svd->V))->type_name) PetscCall(BVSetType(svd->V,BVSVEC));
506:     PetscCall(MatCreateVecsEmpty(svd->A,&tr,NULL));
507:     PetscCall(BVSetSizesFromVec(svd->V,tr,requested));
508:     PetscCall(VecDestroy(&tr));
509:   } else PetscCall(BVResize(svd->V,requested,PETSC_FALSE));
510:   /* allocate U */
511:   if (svd->leftbasis && !svd->isgeneralized) {
512:     if (!svd->U) PetscCall(SVDGetBV(svd,NULL,&svd->U));
513:     if (!oldsize) {
514:       if (!((PetscObject)(svd->U))->type_name) PetscCall(BVSetType(svd->U,((PetscObject)(svd->V))->type_name));
515:       PetscCall(MatCreateVecsEmpty(svd->A,NULL,&tl));
516:       PetscCall(BVSetSizesFromVec(svd->U,tl,requested));
517:       PetscCall(VecDestroy(&tl));
518:     } else PetscCall(BVResize(svd->U,requested,PETSC_FALSE));
519:   } else if (svd->isgeneralized) {  /* left basis for the GSVD */
520:     if (!svd->U) PetscCall(SVDGetBV(svd,NULL,&svd->U));
521:     if (!oldsize) {
522:       if (!((PetscObject)(svd->U))->type_name) PetscCall(BVSetType(svd->U,((PetscObject)(svd->V))->type_name));
523:       PetscCall(SVDCreateLeftTemplate(svd,&tl));
524:       PetscCall(BVSetSizesFromVec(svd->U,tl,requested));
525:       PetscCall(VecDestroy(&tl));
526:     } else PetscCall(BVResize(svd->U,requested,PETSC_FALSE));
527:   }
528:   PetscFunctionReturn(PETSC_SUCCESS);
529: }