Actual source code: cross.c

slepc-3.17.1 2022-04-11
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:    SLEPc singular value solver: "cross"

 13:    Method: Uses a Hermitian eigensolver for A^T*A
 14: */

 16: #include <slepc/private/svdimpl.h>
 17: #include <slepc/private/epsimpl.h>

 19: typedef struct {
 20:   PetscBool explicitmatrix;
 21:   EPS       eps;
 22:   PetscBool usereps;
 23:   Mat       C,D;
 24: } SVD_CROSS;

 26: typedef struct {
 27:   Mat       A,AT;
 28:   Vec       w,diag;
 29:   PetscBool swapped;
 30: } SVD_CROSS_SHELL;

 32: static PetscErrorCode MatMult_Cross(Mat B,Vec x,Vec y)
 33: {
 34:   SVD_CROSS_SHELL *ctx;

 36:   MatShellGetContext(B,&ctx);
 37:   MatMult(ctx->A,x,ctx->w);
 38:   MatMult(ctx->AT,ctx->w,y);
 39:   PetscFunctionReturn(0);
 40: }

 42: static PetscErrorCode MatGetDiagonal_Cross(Mat B,Vec d)
 43: {
 44:   SVD_CROSS_SHELL   *ctx;
 45:   PetscMPIInt       len;
 46:   PetscInt          N,n,i,j,start,end,ncols;
 47:   PetscScalar       *work1,*work2,*diag;
 48:   const PetscInt    *cols;
 49:   const PetscScalar *vals;

 51:   MatShellGetContext(B,&ctx);
 52:   if (!ctx->diag) {
 53:     /* compute diagonal from rows and store in ctx->diag */
 54:     VecDuplicate(d,&ctx->diag);
 55:     MatGetSize(ctx->A,NULL,&N);
 56:     MatGetLocalSize(ctx->A,NULL,&n);
 57:     PetscCalloc2(N,&work1,N,&work2);
 58:     if (ctx->swapped) {
 59:       MatGetOwnershipRange(ctx->AT,&start,&end);
 60:       for (i=start;i<end;i++) {
 61:         MatGetRow(ctx->AT,i,&ncols,NULL,&vals);
 62:         for (j=0;j<ncols;j++) work1[i] += vals[j]*vals[j];
 63:         MatRestoreRow(ctx->AT,i,&ncols,NULL,&vals);
 64:       }
 65:     } else {
 66:       MatGetOwnershipRange(ctx->A,&start,&end);
 67:       for (i=start;i<end;i++) {
 68:         MatGetRow(ctx->A,i,&ncols,&cols,&vals);
 69:         for (j=0;j<ncols;j++) work1[cols[j]] += vals[j]*vals[j];
 70:         MatRestoreRow(ctx->A,i,&ncols,&cols,&vals);
 71:       }
 72:     }
 73:     PetscMPIIntCast(N,&len);
 74:     MPIU_Allreduce(work1,work2,len,MPIU_SCALAR,MPIU_SUM,PetscObjectComm((PetscObject)B));
 75:     VecGetOwnershipRange(ctx->diag,&start,&end);
 76:     VecGetArrayWrite(ctx->diag,&diag);
 77:     for (i=start;i<end;i++) diag[i-start] = work2[i];
 78:     VecRestoreArrayWrite(ctx->diag,&diag);
 79:     PetscFree2(work1,work2);
 80:   }
 81:   VecCopy(ctx->diag,d);
 82:   PetscFunctionReturn(0);
 83: }

 85: static PetscErrorCode MatDestroy_Cross(Mat B)
 86: {
 87:   SVD_CROSS_SHELL *ctx;

 89:   MatShellGetContext(B,&ctx);
 90:   VecDestroy(&ctx->w);
 91:   VecDestroy(&ctx->diag);
 92:   PetscFree(ctx);
 93:   PetscFunctionReturn(0);
 94: }

 96: static PetscErrorCode SVDCrossGetProductMat(SVD svd,Mat A,Mat AT,Mat *C)
 97: {
 98:   SVD_CROSS       *cross = (SVD_CROSS*)svd->data;
 99:   SVD_CROSS_SHELL *ctx;
100:   PetscInt        n;
101:   VecType         vtype;

103:   if (cross->explicitmatrix) {
104:     if (svd->expltrans) {  /* explicit transpose */
105:       MatProductCreate(AT,A,NULL,C);
106:       MatProductSetType(*C,MATPRODUCT_AB);
107:     } else {  /* implicit transpose */
108: #if defined(PETSC_USE_COMPLEX)
109:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"Must use explicit transpose with complex scalars");
110: #else
111:       if (!svd->swapped) {
112:         MatProductCreate(A,A,NULL,C);
113:         MatProductSetType(*C,MATPRODUCT_AtB);
114:       } else {
115:         MatProductCreate(AT,AT,NULL,C);
116:         MatProductSetType(*C,MATPRODUCT_ABt);
117:       }
118: #endif
119:     }
120:     MatProductSetFromOptions(*C);
121:     MatProductSymbolic(*C);
122:     MatProductNumeric(*C);
123:   } else {
124:     PetscNew(&ctx);
125:     ctx->A       = A;
126:     ctx->AT      = AT;
127:     ctx->swapped = svd->swapped;
128:     MatCreateVecs(A,NULL,&ctx->w);
129:     PetscLogObjectParent((PetscObject)svd,(PetscObject)ctx->w);
130:     MatGetLocalSize(A,NULL,&n);
131:     MatCreateShell(PetscObjectComm((PetscObject)svd),n,n,PETSC_DETERMINE,PETSC_DETERMINE,(void*)ctx,C);
132:     MatShellSetOperation(*C,MATOP_MULT,(void(*)(void))MatMult_Cross);
133:     MatShellSetOperation(*C,MATOP_GET_DIAGONAL,(void(*)(void))MatGetDiagonal_Cross);
134:     MatShellSetOperation(*C,MATOP_DESTROY,(void(*)(void))MatDestroy_Cross);
135:     MatGetVecType(A,&vtype);
136:     MatSetVecType(*C,vtype);
137:   }
138:   PetscLogObjectParent((PetscObject)svd,(PetscObject)*C);
139:   PetscFunctionReturn(0);
140: }

142: /* Convergence test relative to the norm of R (used in GSVD only) */
143: static PetscErrorCode EPSConv_Cross(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
144: {
145:   SVD svd = (SVD)ctx;

147:   *errest = res/PetscMax(svd->nrma,svd->nrmb);
148:   PetscFunctionReturn(0);
149: }

151: PetscErrorCode SVDSetUp_Cross(SVD svd)
152: {
153:   SVD_CROSS      *cross = (SVD_CROSS*)svd->data;
154:   ST             st;
155:   PetscBool      trackall,issinv;

157:   if (!cross->eps) SVDCrossGetEPS(svd,&cross->eps);
158:   MatDestroy(&cross->C);
159:   MatDestroy(&cross->D);
160:   SVDCrossGetProductMat(svd,svd->A,svd->AT,&cross->C);
161:   if (svd->isgeneralized) {
162:     SVDCrossGetProductMat(svd,svd->B,svd->BT,&cross->D);
163:     EPSSetOperators(cross->eps,cross->C,cross->D);
164:     EPSSetProblemType(cross->eps,EPS_GHEP);
165:   } else {
166:     EPSSetOperators(cross->eps,cross->C,NULL);
167:     EPSSetProblemType(cross->eps,EPS_HEP);
168:   }
169:   if (!cross->usereps) {
170:     EPSGetST(cross->eps,&st);
171:     if (svd->isgeneralized && svd->which==SVD_SMALLEST) {
172:       STSetType(st,STSINVERT);
173:       EPSSetTarget(cross->eps,0.0);
174:       EPSSetWhichEigenpairs(cross->eps,EPS_TARGET_REAL);
175:     } else {
176:       PetscObjectTypeCompare((PetscObject)st,STSINVERT,&issinv);
177:       if (issinv) EPSSetWhichEigenpairs(cross->eps,EPS_TARGET_MAGNITUDE);
178:       else EPSSetWhichEigenpairs(cross->eps,svd->which==SVD_LARGEST?EPS_LARGEST_REAL:EPS_SMALLEST_REAL);
179:     }
180:     EPSSetDimensions(cross->eps,svd->nsv,svd->ncv,svd->mpd);
181:     EPSSetTolerances(cross->eps,svd->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL/10.0:svd->tol,svd->max_it);
182:     switch (svd->conv) {
183:     case SVD_CONV_ABS:
184:       EPSSetConvergenceTest(cross->eps,EPS_CONV_ABS);break;
185:     case SVD_CONV_REL:
186:       EPSSetConvergenceTest(cross->eps,EPS_CONV_REL);break;
187:     case SVD_CONV_NORM:
188:       if (svd->isgeneralized) {
189:         if (!svd->nrma) MatNorm(svd->OP,NORM_INFINITY,&svd->nrma);
190:         if (!svd->nrmb) MatNorm(svd->OPb,NORM_INFINITY,&svd->nrmb);
191:         EPSSetConvergenceTestFunction(cross->eps,EPSConv_Cross,svd,NULL);
192:       } else {
193:         EPSSetConvergenceTest(cross->eps,EPS_CONV_NORM);break;
194:       }
195:       break;
196:     case SVD_CONV_MAXIT:
197:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"Maxit convergence test not supported in this solver");
198:     case SVD_CONV_USER:
199:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"User-defined convergence test not supported in this solver");
200:     }
201:   }
202:   SVDCheckUnsupported(svd,SVD_FEATURE_STOPPING);
203:   /* Transfer the trackall option from svd to eps */
204:   SVDGetTrackAll(svd,&trackall);
205:   EPSSetTrackAll(cross->eps,trackall);
206:   /* Transfer the initial space from svd to eps */
207:   if (svd->nini<0) {
208:     EPSSetInitialSpace(cross->eps,-svd->nini,svd->IS);
209:     SlepcBasisDestroy_Private(&svd->nini,&svd->IS);
210:   }
211:   EPSSetUp(cross->eps);
212:   EPSGetDimensions(cross->eps,NULL,&svd->ncv,&svd->mpd);
213:   EPSGetTolerances(cross->eps,NULL,&svd->max_it);
214:   if (svd->tol==PETSC_DEFAULT) svd->tol = SLEPC_DEFAULT_TOL;

216:   svd->leftbasis = PETSC_FALSE;
217:   SVDAllocateSolution(svd,0);
218:   PetscFunctionReturn(0);
219: }

221: PetscErrorCode SVDSolve_Cross(SVD svd)
222: {
223:   SVD_CROSS      *cross = (SVD_CROSS*)svd->data;
224:   PetscInt       i;
225:   PetscScalar    lambda;
226:   PetscReal      sigma;

228:   EPSSolve(cross->eps);
229:   EPSGetConverged(cross->eps,&svd->nconv);
230:   EPSGetIterationNumber(cross->eps,&svd->its);
231:   EPSGetConvergedReason(cross->eps,(EPSConvergedReason*)&svd->reason);
232:   for (i=0;i<svd->nconv;i++) {
233:     EPSGetEigenvalue(cross->eps,i,&lambda,NULL);
234:     sigma = PetscRealPart(lambda);
236:     if (sigma<0.0) {
237:       PetscInfo(svd,"Negative eigenvalue computed by EPS: %g, resetting to 0\n",(double)sigma);
238:       sigma = 0.0;
239:     }
240:     svd->sigma[i] = PetscSqrtReal(sigma);
241:   }
242:   PetscFunctionReturn(0);
243: }

245: PetscErrorCode SVDComputeVectors_Cross(SVD svd)
246: {
247:   SVD_CROSS         *cross = (SVD_CROSS*)svd->data;
248:   PetscInt          i,mloc,ploc;
249:   Vec               u,v,x,uv;
250:   PetscScalar       *dst,alpha,lambda;
251:   const PetscScalar *src;
252:   PetscReal         nrm;

254:   if (svd->isgeneralized) {
255:     MatCreateVecs(svd->A,NULL,&u);
256:     VecGetLocalSize(u,&mloc);
257:     MatCreateVecs(svd->B,NULL,&v);
258:     VecGetLocalSize(v,&ploc);
259:     for (i=0;i<svd->nconv;i++) {
260:       BVGetColumn(svd->V,i,&x);
261:       EPSGetEigenpair(cross->eps,i,&lambda,NULL,x,NULL);
262:       MatMult(svd->A,x,u);     /* u_i*c_i/alpha = A*x_i */
263:       VecNormalize(u,NULL);
264:       MatMult(svd->B,x,v);     /* v_i*s_i/alpha = B*x_i */
265:       VecNormalize(v,&nrm);    /* ||v||_2 = s_i/alpha   */
266:       alpha = 1.0/(PetscSqrtReal(1.0+PetscRealPart(lambda))*nrm);    /* alpha=s_i/||v||_2 */
267:       VecScale(x,alpha);
268:       BVRestoreColumn(svd->V,i,&x);
269:       /* copy [u;v] to U[i] */
270:       BVGetColumn(svd->U,i,&uv);
271:       VecGetArrayWrite(uv,&dst);
272:       VecGetArrayRead(u,&src);
273:       PetscArraycpy(dst,src,mloc);
274:       VecRestoreArrayRead(u,&src);
275:       VecGetArrayRead(v,&src);
276:       PetscArraycpy(dst+mloc,src,ploc);
277:       VecRestoreArrayRead(v,&src);
278:       VecRestoreArrayWrite(uv,&dst);
279:       BVRestoreColumn(svd->U,i,&uv);
280:     }
281:     VecDestroy(&v);
282:     VecDestroy(&u);
283:   } else {
284:     for (i=0;i<svd->nconv;i++) {
285:       BVGetColumn(svd->V,i,&v);
286:       EPSGetEigenvector(cross->eps,i,v,NULL);
287:       BVRestoreColumn(svd->V,i,&v);
288:     }
289:     SVDComputeVectors_Left(svd);
290:   }
291:   PetscFunctionReturn(0);
292: }

294: static PetscErrorCode EPSMonitor_Cross(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)
295: {
296:   PetscInt       i;
297:   SVD            svd = (SVD)ctx;
298:   PetscScalar    er,ei;

300:   for (i=0;i<PetscMin(nest,svd->ncv);i++) {
301:     er = eigr[i]; ei = eigi[i];
302:     STBackTransform(eps->st,1,&er,&ei);
303:     svd->sigma[i] = PetscSqrtReal(PetscAbsReal(PetscRealPart(er)));
304:     svd->errest[i] = errest[i];
305:   }
306:   SVDMonitor(svd,its,nconv,svd->sigma,svd->errest,nest);
307:   PetscFunctionReturn(0);
308: }

310: PetscErrorCode SVDSetFromOptions_Cross(PetscOptionItems *PetscOptionsObject,SVD svd)
311: {
312:   PetscBool      set,val;
313:   SVD_CROSS      *cross = (SVD_CROSS*)svd->data;
314:   ST             st;

316:   PetscOptionsHead(PetscOptionsObject,"SVD Cross Options");

318:     PetscOptionsBool("-svd_cross_explicitmatrix","Use cross explicit matrix","SVDCrossSetExplicitMatrix",cross->explicitmatrix,&val,&set);
319:     if (set) SVDCrossSetExplicitMatrix(svd,val);

321:   PetscOptionsTail();

323:   if (!cross->eps) SVDCrossGetEPS(svd,&cross->eps);
324:   if (!cross->explicitmatrix && !cross->usereps) {
325:     /* use as default an ST with shell matrix and Jacobi */
326:     EPSGetST(cross->eps,&st);
327:     STSetMatMode(st,ST_MATMODE_SHELL);
328:   }
329:   EPSSetFromOptions(cross->eps);
330:   PetscFunctionReturn(0);
331: }

333: static PetscErrorCode SVDCrossSetExplicitMatrix_Cross(SVD svd,PetscBool explicitmatrix)
334: {
335:   SVD_CROSS *cross = (SVD_CROSS*)svd->data;

337:   if (cross->explicitmatrix != explicitmatrix) {
338:     cross->explicitmatrix = explicitmatrix;
339:     svd->state = SVD_STATE_INITIAL;
340:   }
341:   PetscFunctionReturn(0);
342: }

344: /*@
345:    SVDCrossSetExplicitMatrix - Indicate if the eigensolver operator A^T*A must
346:    be computed explicitly.

348:    Logically Collective on svd

350:    Input Parameters:
351: +  svd         - singular value solver
352: -  explicitmat - boolean flag indicating if A^T*A is built explicitly

354:    Options Database Key:
355: .  -svd_cross_explicitmatrix <boolean> - Indicates the boolean flag

357:    Level: advanced

359: .seealso: SVDCrossGetExplicitMatrix()
360: @*/
361: PetscErrorCode SVDCrossSetExplicitMatrix(SVD svd,PetscBool explicitmat)
362: {
365:   PetscTryMethod(svd,"SVDCrossSetExplicitMatrix_C",(SVD,PetscBool),(svd,explicitmat));
366:   PetscFunctionReturn(0);
367: }

369: static PetscErrorCode SVDCrossGetExplicitMatrix_Cross(SVD svd,PetscBool *explicitmat)
370: {
371:   SVD_CROSS *cross = (SVD_CROSS*)svd->data;

373:   *explicitmat = cross->explicitmatrix;
374:   PetscFunctionReturn(0);
375: }

377: /*@
378:    SVDCrossGetExplicitMatrix - Returns the flag indicating if A^T*A is built explicitly.

380:    Not Collective

382:    Input Parameter:
383: .  svd  - singular value solver

385:    Output Parameter:
386: .  explicitmat - the mode flag

388:    Level: advanced

390: .seealso: SVDCrossSetExplicitMatrix()
391: @*/
392: PetscErrorCode SVDCrossGetExplicitMatrix(SVD svd,PetscBool *explicitmat)
393: {
396:   PetscUseMethod(svd,"SVDCrossGetExplicitMatrix_C",(SVD,PetscBool*),(svd,explicitmat));
397:   PetscFunctionReturn(0);
398: }

400: static PetscErrorCode SVDCrossSetEPS_Cross(SVD svd,EPS eps)
401: {
402:   SVD_CROSS      *cross = (SVD_CROSS*)svd->data;

404:   PetscObjectReference((PetscObject)eps);
405:   EPSDestroy(&cross->eps);
406:   cross->eps = eps;
407:   cross->usereps = PETSC_TRUE;
408:   PetscLogObjectParent((PetscObject)svd,(PetscObject)cross->eps);
409:   svd->state = SVD_STATE_INITIAL;
410:   PetscFunctionReturn(0);
411: }

413: /*@
414:    SVDCrossSetEPS - Associate an eigensolver object (EPS) to the
415:    singular value solver.

417:    Collective on svd

419:    Input Parameters:
420: +  svd - singular value solver
421: -  eps - the eigensolver object

423:    Level: advanced

425: .seealso: SVDCrossGetEPS()
426: @*/
427: PetscErrorCode SVDCrossSetEPS(SVD svd,EPS eps)
428: {
432:   PetscTryMethod(svd,"SVDCrossSetEPS_C",(SVD,EPS),(svd,eps));
433:   PetscFunctionReturn(0);
434: }

436: static PetscErrorCode SVDCrossGetEPS_Cross(SVD svd,EPS *eps)
437: {
438:   SVD_CROSS      *cross = (SVD_CROSS*)svd->data;

440:   if (!cross->eps) {
441:     EPSCreate(PetscObjectComm((PetscObject)svd),&cross->eps);
442:     PetscObjectIncrementTabLevel((PetscObject)cross->eps,(PetscObject)svd,1);
443:     EPSSetOptionsPrefix(cross->eps,((PetscObject)svd)->prefix);
444:     EPSAppendOptionsPrefix(cross->eps,"svd_cross_");
445:     PetscLogObjectParent((PetscObject)svd,(PetscObject)cross->eps);
446:     PetscObjectSetOptions((PetscObject)cross->eps,((PetscObject)svd)->options);
447:     EPSSetWhichEigenpairs(cross->eps,EPS_LARGEST_REAL);
448:     EPSMonitorSet(cross->eps,EPSMonitor_Cross,svd,NULL);
449:   }
450:   *eps = cross->eps;
451:   PetscFunctionReturn(0);
452: }

454: /*@
455:    SVDCrossGetEPS - Retrieve the eigensolver object (EPS) associated
456:    to the singular value solver.

458:    Not Collective

460:    Input Parameter:
461: .  svd - singular value solver

463:    Output Parameter:
464: .  eps - the eigensolver object

466:    Level: advanced

468: .seealso: SVDCrossSetEPS()
469: @*/
470: PetscErrorCode SVDCrossGetEPS(SVD svd,EPS *eps)
471: {
474:   PetscUseMethod(svd,"SVDCrossGetEPS_C",(SVD,EPS*),(svd,eps));
475:   PetscFunctionReturn(0);
476: }

478: PetscErrorCode SVDView_Cross(SVD svd,PetscViewer viewer)
479: {
480:   SVD_CROSS      *cross = (SVD_CROSS*)svd->data;
481:   PetscBool      isascii;

483:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
484:   if (isascii) {
485:     if (!cross->eps) SVDCrossGetEPS(svd,&cross->eps);
486:     PetscViewerASCIIPrintf(viewer,"  %s matrix\n",cross->explicitmatrix?"explicit":"implicit");
487:     PetscViewerASCIIPushTab(viewer);
488:     EPSView(cross->eps,viewer);
489:     PetscViewerASCIIPopTab(viewer);
490:   }
491:   PetscFunctionReturn(0);
492: }

494: PetscErrorCode SVDReset_Cross(SVD svd)
495: {
496:   SVD_CROSS      *cross = (SVD_CROSS*)svd->data;

498:   EPSReset(cross->eps);
499:   MatDestroy(&cross->C);
500:   MatDestroy(&cross->D);
501:   PetscFunctionReturn(0);
502: }

504: PetscErrorCode SVDDestroy_Cross(SVD svd)
505: {
506:   SVD_CROSS      *cross = (SVD_CROSS*)svd->data;

508:   EPSDestroy(&cross->eps);
509:   PetscFree(svd->data);
510:   PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetEPS_C",NULL);
511:   PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetEPS_C",NULL);
512:   PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetExplicitMatrix_C",NULL);
513:   PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetExplicitMatrix_C",NULL);
514:   PetscFunctionReturn(0);
515: }

517: SLEPC_EXTERN PetscErrorCode SVDCreate_Cross(SVD svd)
518: {
519:   SVD_CROSS      *cross;

521:   PetscNewLog(svd,&cross);
522:   svd->data = (void*)cross;

524:   svd->ops->solve          = SVDSolve_Cross;
525:   svd->ops->solveg         = SVDSolve_Cross;
526:   svd->ops->setup          = SVDSetUp_Cross;
527:   svd->ops->setfromoptions = SVDSetFromOptions_Cross;
528:   svd->ops->destroy        = SVDDestroy_Cross;
529:   svd->ops->reset          = SVDReset_Cross;
530:   svd->ops->view           = SVDView_Cross;
531:   svd->ops->computevectors = SVDComputeVectors_Cross;
532:   PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetEPS_C",SVDCrossSetEPS_Cross);
533:   PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetEPS_C",SVDCrossGetEPS_Cross);
534:   PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetExplicitMatrix_C",SVDCrossSetExplicitMatrix_Cross);
535:   PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetExplicitMatrix_C",SVDCrossGetExplicitMatrix_Cross);
536:   PetscFunctionReturn(0);
537: }