Actual source code: svdsetup.c

slepc-main 2024-11-09
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_SIZ,"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_SIZ,"Local size of signature (%" PetscInt_FMT ") does not match the local row size of A (%" PetscInt_FMT ")",n,ma);
141:   }

143:   PetscCall(VecDestroy(&svd->omega));
144:   if (omega) {
145:     PetscCall(VecDuplicate(omega,&svd->omega));
146:     PetscCall(VecCopy(omega,svd->omega));
147:   }
148:   svd->state = SVD_STATE_INITIAL;
149:   PetscFunctionReturn(PETSC_SUCCESS);
150: }

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

155:    Not Collective

157:    Input Parameter:
158: .  svd - the singular value solver context

160:    Output Parameter:
161: .  omega - a vector containing the diagonal elements of the signature matrix

163:    Notes:
164:    The signature matrix is relevant only for hyperbolic problems (HSVD).
165:    If no signature has been set, this function will return a vector of all ones.

167:    The user should pass a previously created Vec with the appropriate dimension.

169:    Level: intermediate

171: .seealso: SVDSetSignature()
172: @*/
173: PetscErrorCode SVDGetSignature(SVD svd,Vec omega)
174: {
175:   PetscFunctionBegin;
178:   if (svd->omega) PetscCall(VecCopy(svd->omega,omega));
179:   else PetscCall(VecSet(omega,1.0));
180:   PetscFunctionReturn(PETSC_SUCCESS);
181: }

183: /*@
184:    SVDSetDSType - Sets the type of the internal DS object based on the current
185:    settings of the singular value solver.

187:    Collective

189:    Input Parameter:
190: .  svd - singular value solver context

192:    Note:
193:    This function need not be called explicitly, since it will be called at
194:    both SVDSetFromOptions() and SVDSetUp().

196:    Level: developer

198: .seealso: SVDSetFromOptions(), SVDSetUp()
199: @*/
200: PetscErrorCode SVDSetDSType(SVD svd)
201: {
202:   PetscFunctionBegin;
204:   PetscTryTypeMethod(svd,setdstype);
205:   PetscFunctionReturn(PETSC_SUCCESS);
206: }

208: /*@
209:    SVDSetUp - Sets up all the internal data structures necessary for the
210:    execution of the singular value solver.

212:    Collective

214:    Input Parameter:
215: .  svd   - singular value solver context

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

222:    Level: developer

224: .seealso: SVDCreate(), SVDSolve(), SVDDestroy()
225: @*/
226: PetscErrorCode SVDSetUp(SVD svd)
227: {
228:   PetscBool      flg;
229:   PetscInt       M,N,P=0,k,maxnsol,m,Nom,nom;
230:   SlepcSC        sc;
231:   Vec            *T;
232:   BV             bv;

234:   PetscFunctionBegin;
236:   if (svd->state) PetscFunctionReturn(PETSC_SUCCESS);
237:   PetscCall(PetscLogEventBegin(SVD_SetUp,svd,0,0,0));

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

242:   /* set default solver type (SVDSetFromOptions was not called) */
243:   if (!((PetscObject)svd)->type_name) PetscCall(SVDSetType(svd,SVDCROSS));
244:   if (!svd->ds) PetscCall(SVDGetDS(svd,&svd->ds));

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

249:   if (!svd->problem_type) {  /* set default problem type */
250:     if (svd->OPb) {
251:       PetscCheck(!svd->omega,PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"There is no support yet for generalized hyperbolic problems");
252:       PetscCall(SVDSetProblemType(svd,SVD_GENERALIZED));
253:     } else {
254:       if (svd->omega) PetscCall(SVDSetProblemType(svd,SVD_HYPERBOLIC));
255:       else PetscCall(SVDSetProblemType(svd,SVD_STANDARD));
256:     }
257:   } else {  /* check consistency of problem type set by user */
258:     if (svd->OPb) {
259:       PetscCheck(svd->isgeneralized,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_INCOMP,"Inconsistent SVD state: the problem type does not match the number of matrices");
260:       PetscCheck(!svd->omega,PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"There is no support yet for generalized hyperbolic problems");
261:     } else {
262:       PetscCheck(!svd->isgeneralized,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_INCOMP,"Inconsistent SVD state: the problem type does not match the number of matrices");
263:       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()");
264:       else PetscCheck(!svd->ishyperbolic,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_INCOMP,"Inconsistent SVD state: a hyperbolic problem requires passing a signature with SVDSetSignature()");
265:     }
266:   }

268:   /* set DS type once the default problem type has been determined */
269:   PetscCall(SVDSetDSType(svd));

271:   /* determine how to handle the transpose */
272:   svd->expltrans = PETSC_TRUE;
273:   if (svd->impltrans) svd->expltrans = PETSC_FALSE;
274:   else {
275:     PetscCall(MatHasOperation(svd->OP,MATOP_TRANSPOSE,&flg));
276:     if (!flg) svd->expltrans = PETSC_FALSE;
277:     else {
278:       PetscCall(PetscObjectTypeCompareAny((PetscObject)svd,&flg,SVDLAPACK,SVDSCALAPACK,SVDKSVD,SVDELEMENTAL,""));
279:       if (flg) svd->expltrans = PETSC_FALSE;
280:     }
281:   }

283:   /* get matrix dimensions */
284:   PetscCall(MatGetSize(svd->OP,&M,&N));
285:   if (svd->isgeneralized) {
286:     PetscCall(MatGetSize(svd->OPb,&P,NULL));
287:     PetscCheck(M+P>=N,PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"The case when [A;B] has less rows than columns is not supported");
288:   } else if (svd->ishyperbolic) {
289:     PetscCall(VecGetSize(svd->omega,&Nom));
290:     PetscCall(VecGetLocalSize(svd->omega,&nom));
291:     PetscCall(MatGetLocalSize(svd->OP,&m,NULL));
292:     PetscCheck(Nom==M,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_SIZ,"Global size of signature (%" PetscInt_FMT ") does not match the row size of A (%" PetscInt_FMT ")",Nom,M);
293:     PetscCheck(nom==m,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_SIZ,"Local size of signature (%" PetscInt_FMT ") does not match the local row size of A (%" PetscInt_FMT ")",nom,m);
294:   }

296:   /* build transpose matrix */
297:   PetscCall(MatDestroy(&svd->A));
298:   PetscCall(MatDestroy(&svd->AT));
299:   PetscCall(PetscObjectReference((PetscObject)svd->OP));
300:   if (svd->expltrans) {
301:     if (svd->isgeneralized || M>=N) {
302:       svd->A = svd->OP;
303:       PetscCall(MatHermitianTranspose(svd->OP,MAT_INITIAL_MATRIX,&svd->AT));
304:     } else {
305:       PetscCall(MatHermitianTranspose(svd->OP,MAT_INITIAL_MATRIX,&svd->A));
306:       svd->AT = svd->OP;
307:     }
308:   } else {
309:     if (svd->isgeneralized || M>=N) {
310:       svd->A = svd->OP;
311:       PetscCall(MatCreateHermitianTranspose(svd->OP,&svd->AT));
312:     } else {
313:       PetscCall(MatCreateHermitianTranspose(svd->OP,&svd->A));
314:       svd->AT = svd->OP;
315:     }
316:   }

318:   /* build transpose matrix B for GSVD */
319:   if (svd->isgeneralized) {
320:     PetscCall(MatDestroy(&svd->B));
321:     PetscCall(MatDestroy(&svd->BT));
322:     PetscCall(PetscObjectReference((PetscObject)svd->OPb));
323:     if (svd->expltrans) {
324:       svd->B = svd->OPb;
325:       PetscCall(MatHermitianTranspose(svd->OPb,MAT_INITIAL_MATRIX,&svd->BT));
326:     } else {
327:       svd->B = svd->OPb;
328:       PetscCall(MatCreateHermitianTranspose(svd->OPb,&svd->BT));
329:     }
330:   }

332:   if (!svd->isgeneralized && M<N) {
333:     /* swap initial vectors */
334:     if (svd->nini || svd->ninil) {
335:       T=svd->ISL; svd->ISL=svd->IS; svd->IS=T;
336:       k=svd->ninil; svd->ninil=svd->nini; svd->nini=k;
337:     }
338:     /* swap basis vectors */
339:     if (!svd->swapped) {  /* only the first time in case of multiple calls */
340:       bv=svd->V; svd->V=svd->U; svd->U=bv;
341:       svd->swapped = PETSC_TRUE;
342:     }
343:   }

345:   maxnsol = svd->isgeneralized? PetscMin(PetscMin(M,N),P): PetscMin(M,N);
346:   svd->ncv = PetscMin(svd->ncv,maxnsol);
347:   svd->nsv = PetscMin(svd->nsv,maxnsol);
348:   PetscCheck(svd->ncv==PETSC_DETERMINE || svd->nsv<=svd->ncv,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nsv bigger than ncv");

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

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

357:   /* call specific solver setup */
358:   PetscUseTypeMethod(svd,setup);

360:   /* set tolerance if not yet set */
361:   if (svd->tol==(PetscReal)PETSC_DETERMINE) svd->tol = SLEPC_DEFAULT_TOL;

363:   /* fill sorting criterion context */
364:   PetscCall(DSGetSlepcSC(svd->ds,&sc));
365:   sc->comparison    = (svd->which==SVD_LARGEST)? SlepcCompareLargestReal: SlepcCompareSmallestReal;
366:   sc->comparisonctx = NULL;
367:   sc->map           = NULL;
368:   sc->mapobj        = NULL;

370:   /* process initial vectors */
371:   if (svd->nini<0) {
372:     k = -svd->nini;
373:     PetscCheck(k<=svd->ncv,PetscObjectComm((PetscObject)svd),PETSC_ERR_USER_INPUT,"The number of initial vectors is larger than ncv");
374:     PetscCall(BVInsertVecs(svd->V,0,&k,svd->IS,PETSC_TRUE));
375:     PetscCall(SlepcBasisDestroy_Private(&svd->nini,&svd->IS));
376:     svd->nini = k;
377:   }
378:   if (svd->ninil<0) {
379:     k = 0;
380:     if (svd->leftbasis) {
381:       k = -svd->ninil;
382:       PetscCheck(k<=svd->ncv,PetscObjectComm((PetscObject)svd),PETSC_ERR_USER_INPUT,"The number of left initial vectors is larger than ncv");
383:       PetscCall(BVInsertVecs(svd->U,0,&k,svd->ISL,PETSC_TRUE));
384:     } else PetscCall(PetscInfo(svd,"Ignoring initial left vectors\n"));
385:     PetscCall(SlepcBasisDestroy_Private(&svd->ninil,&svd->ISL));
386:     svd->ninil = k;
387:   }

389:   PetscCall(PetscLogEventEnd(SVD_SetUp,svd,0,0,0));
390:   svd->state = SVD_STATE_SETUP;
391:   PetscFunctionReturn(PETSC_SUCCESS);
392: }

394: /*@
395:    SVDSetInitialSpaces - Specify two basis of vectors that constitute the initial
396:    right and/or left spaces.

398:    Collective

400:    Input Parameters:
401: +  svd   - the singular value solver context
402: .  nr    - number of right vectors
403: .  isr   - set of basis vectors of the right initial space
404: .  nl    - number of left vectors
405: -  isl   - set of basis vectors of the left initial space

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

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

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

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

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

424:    Level: intermediate

426: .seealso: SVDSetUp()
427: @*/
428: PetscErrorCode SVDSetInitialSpaces(SVD svd,PetscInt nr,Vec isr[],PetscInt nl,Vec isl[])
429: {
430:   PetscFunctionBegin;
434:   PetscCheck(nr>=0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Argument nr cannot be negative");
435:   if (nr>0) {
436:     PetscAssertPointer(isr,3);
438:   }
439:   PetscCheck(nl>=0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Argument nl cannot be negative");
440:   if (nl>0) {
441:     PetscAssertPointer(isl,5);
443:   }
444:   PetscCall(SlepcBasisReference_Private(nr,isr,&svd->nini,&svd->IS));
445:   PetscCall(SlepcBasisReference_Private(nl,isl,&svd->ninil,&svd->ISL));
446:   if (nr>0 || nl>0) svd->state = SVD_STATE_INITIAL;
447:   PetscFunctionReturn(PETSC_SUCCESS);
448: }

450: /*
451:   SVDSetDimensions_Default - Set reasonable values for ncv, mpd if not set
452:   by the user. This is called at setup.
453:  */
454: PetscErrorCode SVDSetDimensions_Default(SVD svd)
455: {
456:   PetscInt       N,M,P,maxnsol;

458:   PetscFunctionBegin;
459:   PetscCall(MatGetSize(svd->OP,&M,&N));
460:   maxnsol = PetscMin(M,N);
461:   if (svd->isgeneralized) {
462:     PetscCall(MatGetSize(svd->OPb,&P,NULL));
463:     maxnsol = PetscMin(maxnsol,P);
464:   }
465:   if (svd->ncv!=PETSC_DETERMINE) { /* ncv set */
466:     PetscCheck(svd->ncv>=svd->nsv,PetscObjectComm((PetscObject)svd),PETSC_ERR_USER_INPUT,"The value of ncv must be at least nsv");
467:   } else if (svd->mpd!=PETSC_DETERMINE) { /* mpd set */
468:     svd->ncv = PetscMin(maxnsol,svd->nsv+svd->mpd);
469:   } else { /* neither set: defaults depend on nsv being small or large */
470:     if (svd->nsv<500) svd->ncv = PetscMin(maxnsol,PetscMax(2*svd->nsv,10));
471:     else {
472:       svd->mpd = 500;
473:       svd->ncv = PetscMin(maxnsol,svd->nsv+svd->mpd);
474:     }
475:   }
476:   if (svd->mpd==PETSC_DETERMINE) svd->mpd = svd->ncv;
477:   PetscFunctionReturn(PETSC_SUCCESS);
478: }

480: /*@
481:    SVDAllocateSolution - Allocate memory storage for common variables such
482:    as the singular values and the basis vectors.

484:    Collective

486:    Input Parameters:
487: +  svd   - eigensolver context
488: -  extra - number of additional positions, used for methods that require a
489:            working basis slightly larger than ncv

491:    Developer Notes:
492:    This is SLEPC_EXTERN because it may be required by user plugin SVD
493:    implementations.

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

497:    Level: developer

499: .seealso: SVDSetUp()
500: @*/
501: PetscErrorCode SVDAllocateSolution(SVD svd,PetscInt extra)
502: {
503:   PetscInt       oldsize,requested;
504:   Vec            tr,tl;

506:   PetscFunctionBegin;
507:   requested = svd->ncv + extra;

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

512:   /* allocate sigma */
513:   if (requested != oldsize || !svd->sigma) {
514:     PetscCall(PetscFree3(svd->sigma,svd->perm,svd->errest));
515:     if (svd->sign) PetscCall(PetscFree(svd->sign));
516:     PetscCall(PetscMalloc3(requested,&svd->sigma,requested,&svd->perm,requested,&svd->errest));
517:     if (svd->ishyperbolic) PetscCall(PetscMalloc1(requested,&svd->sign));
518:   }
519:   /* allocate V */
520:   if (!svd->V) PetscCall(SVDGetBV(svd,&svd->V,NULL));
521:   if (!oldsize) {
522:     if (!((PetscObject)svd->V)->type_name) PetscCall(BVSetType(svd->V,BVMAT));
523:     PetscCall(MatCreateVecsEmpty(svd->A,&tr,NULL));
524:     PetscCall(BVSetSizesFromVec(svd->V,tr,requested));
525:     PetscCall(VecDestroy(&tr));
526:   } else PetscCall(BVResize(svd->V,requested,PETSC_FALSE));
527:   /* allocate U */
528:   if (svd->leftbasis && !svd->isgeneralized) {
529:     if (!svd->U) PetscCall(SVDGetBV(svd,NULL,&svd->U));
530:     if (!oldsize) {
531:       if (!((PetscObject)svd->U)->type_name) PetscCall(BVSetType(svd->U,((PetscObject)svd->V)->type_name));
532:       PetscCall(MatCreateVecsEmpty(svd->A,NULL,&tl));
533:       PetscCall(BVSetSizesFromVec(svd->U,tl,requested));
534:       PetscCall(VecDestroy(&tl));
535:     } else PetscCall(BVResize(svd->U,requested,PETSC_FALSE));
536:   } else if (svd->isgeneralized) {  /* left basis for the GSVD */
537:     if (!svd->U) PetscCall(SVDGetBV(svd,NULL,&svd->U));
538:     if (!oldsize) {
539:       if (!((PetscObject)svd->U)->type_name) PetscCall(BVSetType(svd->U,((PetscObject)svd->V)->type_name));
540:       PetscCall(SVDCreateLeftTemplate(svd,&tl));
541:       PetscCall(BVSetSizesFromVec(svd->U,tl,requested));
542:       PetscCall(VecDestroy(&tl));
543:     } else PetscCall(BVResize(svd->U,requested,PETSC_FALSE));
544:   }
545:   PetscFunctionReturn(PETSC_SUCCESS);
546: }