Actual source code: cross.c
slepc-3.18.0 2022-10-01
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,omega;
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: if (ctx->omega && !ctx->swapped) VecPointwiseMult(ctx->w,ctx->w,ctx->omega);
39: MatMult(ctx->AT,ctx->w,y);
40: return 0;
41: }
43: static PetscErrorCode MatGetDiagonal_Cross(Mat B,Vec d)
44: {
45: SVD_CROSS_SHELL *ctx;
46: PetscMPIInt len;
47: PetscInt N,n,i,j,start,end,ncols;
48: PetscScalar *work1,*work2,*diag;
49: const PetscInt *cols;
50: const PetscScalar *vals;
52: MatShellGetContext(B,&ctx);
53: if (!ctx->diag) {
54: /* compute diagonal from rows and store in ctx->diag */
55: VecDuplicate(d,&ctx->diag);
56: MatGetSize(ctx->A,NULL,&N);
57: MatGetLocalSize(ctx->A,NULL,&n);
58: PetscCalloc2(N,&work1,N,&work2);
59: if (ctx->swapped) {
60: MatGetOwnershipRange(ctx->AT,&start,&end);
61: for (i=start;i<end;i++) {
62: MatGetRow(ctx->AT,i,&ncols,NULL,&vals);
63: for (j=0;j<ncols;j++) work1[i] += vals[j]*vals[j];
64: MatRestoreRow(ctx->AT,i,&ncols,NULL,&vals);
65: }
66: } else {
67: MatGetOwnershipRange(ctx->A,&start,&end);
68: for (i=start;i<end;i++) {
69: MatGetRow(ctx->A,i,&ncols,&cols,&vals);
70: for (j=0;j<ncols;j++) work1[cols[j]] += vals[j]*vals[j];
71: MatRestoreRow(ctx->A,i,&ncols,&cols,&vals);
72: }
73: }
74: PetscMPIIntCast(N,&len);
75: MPIU_Allreduce(work1,work2,len,MPIU_SCALAR,MPIU_SUM,PetscObjectComm((PetscObject)B));
76: VecGetOwnershipRange(ctx->diag,&start,&end);
77: VecGetArrayWrite(ctx->diag,&diag);
78: for (i=start;i<end;i++) diag[i-start] = work2[i];
79: VecRestoreArrayWrite(ctx->diag,&diag);
80: PetscFree2(work1,work2);
81: }
82: VecCopy(ctx->diag,d);
83: return 0;
84: }
86: static PetscErrorCode MatDestroy_Cross(Mat B)
87: {
88: SVD_CROSS_SHELL *ctx;
90: MatShellGetContext(B,&ctx);
91: VecDestroy(&ctx->w);
92: VecDestroy(&ctx->diag);
93: PetscFree(ctx);
94: return 0;
95: }
97: static PetscErrorCode SVDCrossGetProductMat(SVD svd,Mat A,Mat AT,Mat *C)
98: {
99: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
100: SVD_CROSS_SHELL *ctx;
101: PetscInt n;
102: VecType vtype;
103: Mat B;
105: if (cross->explicitmatrix) {
106: if (!svd->ishyperbolic || svd->swapped) B = (!svd->expltrans && svd->swapped)? AT: A;
107: else { /* duplicate A and scale by signature */
108: MatDuplicate(A,MAT_COPY_VALUES,&B);
109: MatDiagonalScale(B,svd->omega,NULL);
110: }
111: if (svd->expltrans) { /* explicit transpose */
112: MatProductCreate(AT,B,NULL,C);
113: MatProductSetType(*C,MATPRODUCT_AB);
114: } else { /* implicit transpose */
115: #if defined(PETSC_USE_COMPLEX)
116: SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"Must use explicit transpose with complex scalars");
117: #else
118: if (!svd->swapped) {
119: MatProductCreate(A,B,NULL,C);
120: MatProductSetType(*C,MATPRODUCT_AtB);
121: } else {
122: MatProductCreate(B,AT,NULL,C);
123: MatProductSetType(*C,MATPRODUCT_ABt);
124: }
125: #endif
126: }
127: MatProductSetFromOptions(*C);
128: MatProductSymbolic(*C);
129: MatProductNumeric(*C);
130: if (svd->ishyperbolic && !svd->swapped) MatDestroy(&B);
131: } else {
132: PetscNew(&ctx);
133: ctx->A = A;
134: ctx->AT = AT;
135: ctx->omega = svd->omega;
136: ctx->swapped = svd->swapped;
137: MatCreateVecs(A,NULL,&ctx->w);
138: MatGetLocalSize(A,NULL,&n);
139: MatCreateShell(PetscObjectComm((PetscObject)svd),n,n,PETSC_DETERMINE,PETSC_DETERMINE,(void*)ctx,C);
140: MatShellSetOperation(*C,MATOP_MULT,(void(*)(void))MatMult_Cross);
141: if (!svd->ishyperbolic || svd->swapped) MatShellSetOperation(*C,MATOP_GET_DIAGONAL,(void(*)(void))MatGetDiagonal_Cross);
142: MatShellSetOperation(*C,MATOP_DESTROY,(void(*)(void))MatDestroy_Cross);
143: MatGetVecType(A,&vtype);
144: MatSetVecType(*C,vtype);
145: }
146: return 0;
147: }
149: /* Convergence test relative to the norm of R (used in GSVD only) */
150: static PetscErrorCode EPSConv_Cross(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
151: {
152: SVD svd = (SVD)ctx;
154: *errest = res/PetscMax(svd->nrma,svd->nrmb);
155: return 0;
156: }
158: PetscErrorCode SVDSetUp_Cross(SVD svd)
159: {
160: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
161: ST st;
162: PetscBool trackall,issinv,isks;
163: EPSProblemType ptype;
164: EPSWhich which;
165: Mat Omega;
166: MatType Atype;
167: PetscInt n,N;
169: if (!cross->eps) SVDCrossGetEPS(svd,&cross->eps);
170: MatDestroy(&cross->C);
171: MatDestroy(&cross->D);
172: SVDCrossGetProductMat(svd,svd->A,svd->AT,&cross->C);
173: if (svd->isgeneralized) {
174: SVDCrossGetProductMat(svd,svd->B,svd->BT,&cross->D);
175: EPSSetOperators(cross->eps,cross->C,cross->D);
176: EPSGetProblemType(cross->eps,&ptype);
177: if (!ptype) EPSSetProblemType(cross->eps,EPS_GHEP);
178: } else if (svd->ishyperbolic && svd->swapped) {
179: MatGetType(svd->OP,&Atype);
180: MatGetSize(svd->A,NULL,&N);
181: MatGetLocalSize(svd->A,NULL,&n);
182: MatCreate(PetscObjectComm((PetscObject)svd),&Omega);
183: MatSetSizes(Omega,n,n,N,N);
184: MatSetType(Omega,Atype);
185: MatSetUp(Omega);
186: MatDiagonalSet(Omega,svd->omega,INSERT_VALUES);
187: EPSSetOperators(cross->eps,cross->C,Omega);
188: EPSSetProblemType(cross->eps,EPS_GHIEP);
189: MatDestroy(&Omega);
190: } else {
191: EPSSetOperators(cross->eps,cross->C,NULL);
192: EPSSetProblemType(cross->eps,EPS_HEP);
193: }
194: if (!cross->usereps) {
195: EPSGetST(cross->eps,&st);
196: PetscObjectTypeCompare((PetscObject)st,STSINVERT,&issinv);
197: PetscObjectTypeCompare((PetscObject)cross->eps,EPSKRYLOVSCHUR,&isks);
198: if (svd->isgeneralized && svd->which==SVD_SMALLEST) {
199: if (cross->explicitmatrix && isks && !issinv) { /* default to shift-and-invert */
200: STSetType(st,STSINVERT);
201: EPSSetTarget(cross->eps,0.0);
202: which = EPS_TARGET_REAL;
203: } else which = issinv?EPS_TARGET_REAL:EPS_SMALLEST_REAL;
204: } else {
205: if (issinv) which = EPS_TARGET_MAGNITUDE;
206: else if (svd->ishyperbolic) which = svd->which==SVD_LARGEST?EPS_LARGEST_MAGNITUDE:EPS_SMALLEST_MAGNITUDE;
207: else which = svd->which==SVD_LARGEST?EPS_LARGEST_REAL:EPS_SMALLEST_REAL;
208: }
209: EPSSetWhichEigenpairs(cross->eps,which);
210: EPSSetDimensions(cross->eps,svd->nsv,svd->ncv,svd->mpd);
211: EPSSetTolerances(cross->eps,svd->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL/10.0:svd->tol,svd->max_it);
212: switch (svd->conv) {
213: case SVD_CONV_ABS:
214: EPSSetConvergenceTest(cross->eps,EPS_CONV_ABS);break;
215: case SVD_CONV_REL:
216: EPSSetConvergenceTest(cross->eps,EPS_CONV_REL);break;
217: case SVD_CONV_NORM:
218: if (svd->isgeneralized) {
219: if (!svd->nrma) MatNorm(svd->OP,NORM_INFINITY,&svd->nrma);
220: if (!svd->nrmb) MatNorm(svd->OPb,NORM_INFINITY,&svd->nrmb);
221: EPSSetConvergenceTestFunction(cross->eps,EPSConv_Cross,svd,NULL);
222: } else {
223: EPSSetConvergenceTest(cross->eps,EPS_CONV_NORM);break;
224: }
225: break;
226: case SVD_CONV_MAXIT:
227: SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"Maxit convergence test not supported in this solver");
228: case SVD_CONV_USER:
229: SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"User-defined convergence test not supported in this solver");
230: }
231: }
232: SVDCheckUnsupported(svd,SVD_FEATURE_STOPPING);
233: /* Transfer the trackall option from svd to eps */
234: SVDGetTrackAll(svd,&trackall);
235: EPSSetTrackAll(cross->eps,trackall);
236: /* Transfer the initial space from svd to eps */
237: if (svd->nini<0) {
238: EPSSetInitialSpace(cross->eps,-svd->nini,svd->IS);
239: SlepcBasisDestroy_Private(&svd->nini,&svd->IS);
240: }
241: EPSSetUp(cross->eps);
242: EPSGetDimensions(cross->eps,NULL,&svd->ncv,&svd->mpd);
243: EPSGetTolerances(cross->eps,NULL,&svd->max_it);
244: if (svd->tol==PETSC_DEFAULT) svd->tol = SLEPC_DEFAULT_TOL;
246: svd->leftbasis = PETSC_FALSE;
247: SVDAllocateSolution(svd,0);
248: return 0;
249: }
251: PetscErrorCode SVDSolve_Cross(SVD svd)
252: {
253: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
254: PetscInt i;
255: PetscScalar lambda;
256: PetscReal sigma;
258: EPSSolve(cross->eps);
259: EPSGetConverged(cross->eps,&svd->nconv);
260: EPSGetIterationNumber(cross->eps,&svd->its);
261: EPSGetConvergedReason(cross->eps,(EPSConvergedReason*)&svd->reason);
262: for (i=0;i<svd->nconv;i++) {
263: EPSGetEigenvalue(cross->eps,i,&lambda,NULL);
264: sigma = PetscRealPart(lambda);
265: if (svd->ishyperbolic) svd->sigma[i] = PetscSqrtReal(PetscAbsReal(sigma));
266: else {
268: if (sigma<0.0) {
269: PetscInfo(svd,"Negative eigenvalue computed by EPS: %g, resetting to 0\n",(double)sigma);
270: sigma = 0.0;
271: }
272: svd->sigma[i] = PetscSqrtReal(sigma);
273: }
274: }
275: return 0;
276: }
278: PetscErrorCode SVDComputeVectors_Cross(SVD svd)
279: {
280: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
281: PetscInt i,mloc,ploc;
282: Vec u,v,x,uv,w,omega2=NULL;
283: Mat Omega;
284: PetscScalar *dst,alpha,lambda,*varray;
285: const PetscScalar *src;
286: PetscReal nrm;
288: if (svd->isgeneralized) {
289: MatCreateVecs(svd->A,NULL,&u);
290: VecGetLocalSize(u,&mloc);
291: MatCreateVecs(svd->B,NULL,&v);
292: VecGetLocalSize(v,&ploc);
293: for (i=0;i<svd->nconv;i++) {
294: BVGetColumn(svd->V,i,&x);
295: EPSGetEigenpair(cross->eps,i,&lambda,NULL,x,NULL);
296: MatMult(svd->A,x,u); /* u_i*c_i/alpha = A*x_i */
297: VecNormalize(u,NULL);
298: MatMult(svd->B,x,v); /* v_i*s_i/alpha = B*x_i */
299: VecNormalize(v,&nrm); /* ||v||_2 = s_i/alpha */
300: alpha = 1.0/(PetscSqrtReal(1.0+PetscRealPart(lambda))*nrm); /* alpha=s_i/||v||_2 */
301: VecScale(x,alpha);
302: BVRestoreColumn(svd->V,i,&x);
303: /* copy [u;v] to U[i] */
304: BVGetColumn(svd->U,i,&uv);
305: VecGetArrayWrite(uv,&dst);
306: VecGetArrayRead(u,&src);
307: PetscArraycpy(dst,src,mloc);
308: VecRestoreArrayRead(u,&src);
309: VecGetArrayRead(v,&src);
310: PetscArraycpy(dst+mloc,src,ploc);
311: VecRestoreArrayRead(v,&src);
312: VecRestoreArrayWrite(uv,&dst);
313: BVRestoreColumn(svd->U,i,&uv);
314: }
315: VecDestroy(&v);
316: VecDestroy(&u);
317: } else if (svd->ishyperbolic && svd->swapped) { /* was solved as GHIEP, set u=Omega*u and normalize */
318: EPSGetOperators(cross->eps,NULL,&Omega);
319: MatCreateVecs(Omega,&w,NULL);
320: VecCreateSeq(PETSC_COMM_SELF,svd->ncv,&omega2);
321: VecGetArrayWrite(omega2,&varray);
322: for (i=0;i<svd->nconv;i++) {
323: BVGetColumn(svd->V,i,&v);
324: EPSGetEigenvector(cross->eps,i,v,NULL);
325: MatMult(Omega,v,w);
326: VecDot(v,w,&alpha);
327: svd->sign[i] = PetscSign(PetscRealPart(alpha));
328: varray[i] = svd->sign[i];
329: alpha = 1.0/PetscSqrtScalar(PetscAbsScalar(alpha));
330: VecScale(w,alpha);
331: VecCopy(w,v);
332: BVRestoreColumn(svd->V,i,&v);
333: }
334: BVSetSignature(svd->V,omega2);
335: VecRestoreArrayWrite(omega2,&varray);
336: VecDestroy(&omega2);
337: VecDestroy(&w);
338: SVDComputeVectors_Left(svd);
339: } else {
340: for (i=0;i<svd->nconv;i++) {
341: BVGetColumn(svd->V,i,&v);
342: EPSGetEigenvector(cross->eps,i,v,NULL);
343: BVRestoreColumn(svd->V,i,&v);
344: }
345: SVDComputeVectors_Left(svd);
346: }
347: return 0;
348: }
350: static PetscErrorCode EPSMonitor_Cross(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)
351: {
352: PetscInt i;
353: SVD svd = (SVD)ctx;
354: PetscScalar er,ei;
356: for (i=0;i<PetscMin(nest,svd->ncv);i++) {
357: er = eigr[i]; ei = eigi[i];
358: STBackTransform(eps->st,1,&er,&ei);
359: svd->sigma[i] = PetscSqrtReal(PetscAbsReal(PetscRealPart(er)));
360: svd->errest[i] = errest[i];
361: }
362: SVDMonitor(svd,its,nconv,svd->sigma,svd->errest,nest);
363: return 0;
364: }
366: PetscErrorCode SVDSetFromOptions_Cross(SVD svd,PetscOptionItems *PetscOptionsObject)
367: {
368: PetscBool set,val;
369: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
370: ST st;
372: PetscOptionsHeadBegin(PetscOptionsObject,"SVD Cross Options");
374: PetscOptionsBool("-svd_cross_explicitmatrix","Use cross explicit matrix","SVDCrossSetExplicitMatrix",cross->explicitmatrix,&val,&set);
375: if (set) SVDCrossSetExplicitMatrix(svd,val);
377: PetscOptionsHeadEnd();
379: if (!cross->eps) SVDCrossGetEPS(svd,&cross->eps);
380: if (!cross->explicitmatrix && !cross->usereps) {
381: /* use as default an ST with shell matrix and Jacobi */
382: EPSGetST(cross->eps,&st);
383: STSetMatMode(st,ST_MATMODE_SHELL);
384: }
385: EPSSetFromOptions(cross->eps);
386: return 0;
387: }
389: static PetscErrorCode SVDCrossSetExplicitMatrix_Cross(SVD svd,PetscBool explicitmatrix)
390: {
391: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
393: if (cross->explicitmatrix != explicitmatrix) {
394: cross->explicitmatrix = explicitmatrix;
395: svd->state = SVD_STATE_INITIAL;
396: }
397: return 0;
398: }
400: /*@
401: SVDCrossSetExplicitMatrix - Indicate if the eigensolver operator A^T*A must
402: be computed explicitly.
404: Logically Collective on svd
406: Input Parameters:
407: + svd - singular value solver
408: - explicitmat - boolean flag indicating if A^T*A is built explicitly
410: Options Database Key:
411: . -svd_cross_explicitmatrix <boolean> - Indicates the boolean flag
413: Level: advanced
415: .seealso: SVDCrossGetExplicitMatrix()
416: @*/
417: PetscErrorCode SVDCrossSetExplicitMatrix(SVD svd,PetscBool explicitmat)
418: {
421: PetscTryMethod(svd,"SVDCrossSetExplicitMatrix_C",(SVD,PetscBool),(svd,explicitmat));
422: return 0;
423: }
425: static PetscErrorCode SVDCrossGetExplicitMatrix_Cross(SVD svd,PetscBool *explicitmat)
426: {
427: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
429: *explicitmat = cross->explicitmatrix;
430: return 0;
431: }
433: /*@
434: SVDCrossGetExplicitMatrix - Returns the flag indicating if A^T*A is built explicitly.
436: Not Collective
438: Input Parameter:
439: . svd - singular value solver
441: Output Parameter:
442: . explicitmat - the mode flag
444: Level: advanced
446: .seealso: SVDCrossSetExplicitMatrix()
447: @*/
448: PetscErrorCode SVDCrossGetExplicitMatrix(SVD svd,PetscBool *explicitmat)
449: {
452: PetscUseMethod(svd,"SVDCrossGetExplicitMatrix_C",(SVD,PetscBool*),(svd,explicitmat));
453: return 0;
454: }
456: static PetscErrorCode SVDCrossSetEPS_Cross(SVD svd,EPS eps)
457: {
458: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
460: PetscObjectReference((PetscObject)eps);
461: EPSDestroy(&cross->eps);
462: cross->eps = eps;
463: cross->usereps = PETSC_TRUE;
464: svd->state = SVD_STATE_INITIAL;
465: return 0;
466: }
468: /*@
469: SVDCrossSetEPS - Associate an eigensolver object (EPS) to the
470: singular value solver.
472: Collective on svd
474: Input Parameters:
475: + svd - singular value solver
476: - eps - the eigensolver object
478: Level: advanced
480: .seealso: SVDCrossGetEPS()
481: @*/
482: PetscErrorCode SVDCrossSetEPS(SVD svd,EPS eps)
483: {
487: PetscTryMethod(svd,"SVDCrossSetEPS_C",(SVD,EPS),(svd,eps));
488: return 0;
489: }
491: static PetscErrorCode SVDCrossGetEPS_Cross(SVD svd,EPS *eps)
492: {
493: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
495: if (!cross->eps) {
496: EPSCreate(PetscObjectComm((PetscObject)svd),&cross->eps);
497: PetscObjectIncrementTabLevel((PetscObject)cross->eps,(PetscObject)svd,1);
498: EPSSetOptionsPrefix(cross->eps,((PetscObject)svd)->prefix);
499: EPSAppendOptionsPrefix(cross->eps,"svd_cross_");
500: PetscObjectSetOptions((PetscObject)cross->eps,((PetscObject)svd)->options);
501: EPSSetWhichEigenpairs(cross->eps,EPS_LARGEST_REAL);
502: EPSMonitorSet(cross->eps,EPSMonitor_Cross,svd,NULL);
503: }
504: *eps = cross->eps;
505: return 0;
506: }
508: /*@
509: SVDCrossGetEPS - Retrieve the eigensolver object (EPS) associated
510: to the singular value solver.
512: Not Collective
514: Input Parameter:
515: . svd - singular value solver
517: Output Parameter:
518: . eps - the eigensolver object
520: Level: advanced
522: .seealso: SVDCrossSetEPS()
523: @*/
524: PetscErrorCode SVDCrossGetEPS(SVD svd,EPS *eps)
525: {
528: PetscUseMethod(svd,"SVDCrossGetEPS_C",(SVD,EPS*),(svd,eps));
529: return 0;
530: }
532: PetscErrorCode SVDView_Cross(SVD svd,PetscViewer viewer)
533: {
534: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
535: PetscBool isascii;
537: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
538: if (isascii) {
539: if (!cross->eps) SVDCrossGetEPS(svd,&cross->eps);
540: PetscViewerASCIIPrintf(viewer," %s matrix\n",cross->explicitmatrix?"explicit":"implicit");
541: PetscViewerASCIIPushTab(viewer);
542: EPSView(cross->eps,viewer);
543: PetscViewerASCIIPopTab(viewer);
544: }
545: return 0;
546: }
548: PetscErrorCode SVDReset_Cross(SVD svd)
549: {
550: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
552: EPSReset(cross->eps);
553: MatDestroy(&cross->C);
554: MatDestroy(&cross->D);
555: return 0;
556: }
558: PetscErrorCode SVDDestroy_Cross(SVD svd)
559: {
560: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
562: EPSDestroy(&cross->eps);
563: PetscFree(svd->data);
564: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetEPS_C",NULL);
565: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetEPS_C",NULL);
566: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetExplicitMatrix_C",NULL);
567: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetExplicitMatrix_C",NULL);
568: return 0;
569: }
571: SLEPC_EXTERN PetscErrorCode SVDCreate_Cross(SVD svd)
572: {
573: SVD_CROSS *cross;
575: PetscNew(&cross);
576: svd->data = (void*)cross;
578: svd->ops->solve = SVDSolve_Cross;
579: svd->ops->solveg = SVDSolve_Cross;
580: svd->ops->solveh = SVDSolve_Cross;
581: svd->ops->setup = SVDSetUp_Cross;
582: svd->ops->setfromoptions = SVDSetFromOptions_Cross;
583: svd->ops->destroy = SVDDestroy_Cross;
584: svd->ops->reset = SVDReset_Cross;
585: svd->ops->view = SVDView_Cross;
586: svd->ops->computevectors = SVDComputeVectors_Cross;
587: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetEPS_C",SVDCrossSetEPS_Cross);
588: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetEPS_C",SVDCrossGetEPS_Cross);
589: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetExplicitMatrix_C",SVDCrossSetExplicitMatrix_Cross);
590: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetExplicitMatrix_C",SVDCrossGetExplicitMatrix_Cross);
591: return 0;
592: }