Actual source code: cross.c
slepc-3.18.3 2023-03-24
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>
18: typedef struct {
19: PetscBool explicitmatrix;
20: EPS eps;
21: PetscBool usereps;
22: Mat C,D;
23: } SVD_CROSS;
25: typedef struct {
26: Mat A,AT;
27: Vec w,diag,omega;
28: PetscBool swapped;
29: } SVD_CROSS_SHELL;
31: static PetscErrorCode MatMult_Cross(Mat B,Vec x,Vec y)
32: {
33: SVD_CROSS_SHELL *ctx;
35: MatShellGetContext(B,&ctx);
36: MatMult(ctx->A,x,ctx->w);
37: if (ctx->omega && !ctx->swapped) VecPointwiseMult(ctx->w,ctx->w,ctx->omega);
38: MatMult(ctx->AT,ctx->w,y);
39: return 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: return 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: return 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;
102: Mat B;
104: if (cross->explicitmatrix) {
105: if (!svd->ishyperbolic || svd->swapped) B = (!svd->expltrans && svd->swapped)? AT: A;
106: else { /* duplicate A and scale by signature */
107: MatDuplicate(A,MAT_COPY_VALUES,&B);
108: MatDiagonalScale(B,svd->omega,NULL);
109: }
110: if (svd->expltrans) { /* explicit transpose */
111: MatProductCreate(AT,B,NULL,C);
112: MatProductSetType(*C,MATPRODUCT_AB);
113: } else { /* implicit transpose */
114: #if defined(PETSC_USE_COMPLEX)
115: SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"Must use explicit transpose with complex scalars");
116: #else
117: if (!svd->swapped) {
118: MatProductCreate(A,B,NULL,C);
119: MatProductSetType(*C,MATPRODUCT_AtB);
120: } else {
121: MatProductCreate(B,AT,NULL,C);
122: MatProductSetType(*C,MATPRODUCT_ABt);
123: }
124: #endif
125: }
126: MatProductSetFromOptions(*C);
127: MatProductSymbolic(*C);
128: MatProductNumeric(*C);
129: if (svd->ishyperbolic && !svd->swapped) MatDestroy(&B);
130: } else {
131: PetscNew(&ctx);
132: ctx->A = A;
133: ctx->AT = AT;
134: ctx->omega = svd->omega;
135: ctx->swapped = svd->swapped;
136: MatCreateVecs(A,NULL,&ctx->w);
137: MatGetLocalSize(A,NULL,&n);
138: MatCreateShell(PetscObjectComm((PetscObject)svd),n,n,PETSC_DETERMINE,PETSC_DETERMINE,(void*)ctx,C);
139: MatShellSetOperation(*C,MATOP_MULT,(void(*)(void))MatMult_Cross);
140: if (!svd->ishyperbolic || svd->swapped) MatShellSetOperation(*C,MATOP_GET_DIAGONAL,(void(*)(void))MatGetDiagonal_Cross);
141: MatShellSetOperation(*C,MATOP_DESTROY,(void(*)(void))MatDestroy_Cross);
142: MatGetVecType(A,&vtype);
143: MatSetVecType(*C,vtype);
144: }
145: return 0;
146: }
148: /* Convergence test relative to the norm of R (used in GSVD only) */
149: static PetscErrorCode EPSConv_Cross(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
150: {
151: SVD svd = (SVD)ctx;
153: *errest = res/PetscMax(svd->nrma,svd->nrmb);
154: return 0;
155: }
157: PetscErrorCode SVDSetUp_Cross(SVD svd)
158: {
159: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
160: ST st;
161: PetscBool trackall,issinv,isks;
162: EPSProblemType ptype;
163: EPSWhich which;
164: Mat Omega;
165: MatType Atype;
166: PetscInt n,N;
168: if (!cross->eps) SVDCrossGetEPS(svd,&cross->eps);
169: MatDestroy(&cross->C);
170: MatDestroy(&cross->D);
171: SVDCrossGetProductMat(svd,svd->A,svd->AT,&cross->C);
172: if (svd->isgeneralized) {
173: SVDCrossGetProductMat(svd,svd->B,svd->BT,&cross->D);
174: EPSSetOperators(cross->eps,cross->C,cross->D);
175: EPSGetProblemType(cross->eps,&ptype);
176: if (!ptype) EPSSetProblemType(cross->eps,EPS_GHEP);
177: } else if (svd->ishyperbolic && svd->swapped) {
178: MatGetType(svd->OP,&Atype);
179: MatGetSize(svd->A,NULL,&N);
180: MatGetLocalSize(svd->A,NULL,&n);
181: MatCreate(PetscObjectComm((PetscObject)svd),&Omega);
182: MatSetSizes(Omega,n,n,N,N);
183: MatSetType(Omega,Atype);
184: MatSetUp(Omega);
185: MatDiagonalSet(Omega,svd->omega,INSERT_VALUES);
186: EPSSetOperators(cross->eps,cross->C,Omega);
187: EPSSetProblemType(cross->eps,EPS_GHIEP);
188: MatDestroy(&Omega);
189: } else {
190: EPSSetOperators(cross->eps,cross->C,NULL);
191: EPSSetProblemType(cross->eps,EPS_HEP);
192: }
193: if (!cross->usereps) {
194: EPSGetST(cross->eps,&st);
195: PetscObjectTypeCompare((PetscObject)st,STSINVERT,&issinv);
196: PetscObjectTypeCompare((PetscObject)cross->eps,EPSKRYLOVSCHUR,&isks);
197: if (svd->isgeneralized && svd->which==SVD_SMALLEST) {
198: if (cross->explicitmatrix && isks && !issinv) { /* default to shift-and-invert */
199: STSetType(st,STSINVERT);
200: EPSSetTarget(cross->eps,0.0);
201: which = EPS_TARGET_REAL;
202: } else which = issinv?EPS_TARGET_REAL:EPS_SMALLEST_REAL;
203: } else {
204: if (issinv) which = EPS_TARGET_MAGNITUDE;
205: else if (svd->ishyperbolic) which = svd->which==SVD_LARGEST?EPS_LARGEST_MAGNITUDE:EPS_SMALLEST_MAGNITUDE;
206: else which = svd->which==SVD_LARGEST?EPS_LARGEST_REAL:EPS_SMALLEST_REAL;
207: }
208: EPSSetWhichEigenpairs(cross->eps,which);
209: EPSSetDimensions(cross->eps,svd->nsv,svd->ncv,svd->mpd);
210: EPSSetTolerances(cross->eps,svd->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL/10.0:svd->tol,svd->max_it);
211: switch (svd->conv) {
212: case SVD_CONV_ABS:
213: EPSSetConvergenceTest(cross->eps,EPS_CONV_ABS);break;
214: case SVD_CONV_REL:
215: EPSSetConvergenceTest(cross->eps,EPS_CONV_REL);break;
216: case SVD_CONV_NORM:
217: if (svd->isgeneralized) {
218: if (!svd->nrma) MatNorm(svd->OP,NORM_INFINITY,&svd->nrma);
219: if (!svd->nrmb) MatNorm(svd->OPb,NORM_INFINITY,&svd->nrmb);
220: EPSSetConvergenceTestFunction(cross->eps,EPSConv_Cross,svd,NULL);
221: } else {
222: EPSSetConvergenceTest(cross->eps,EPS_CONV_NORM);break;
223: }
224: break;
225: case SVD_CONV_MAXIT:
226: SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"Maxit convergence test not supported in this solver");
227: case SVD_CONV_USER:
228: SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"User-defined convergence test not supported in this solver");
229: }
230: }
231: SVDCheckUnsupported(svd,SVD_FEATURE_STOPPING);
232: /* Transfer the trackall option from svd to eps */
233: SVDGetTrackAll(svd,&trackall);
234: EPSSetTrackAll(cross->eps,trackall);
235: /* Transfer the initial space from svd to eps */
236: if (svd->nini<0) {
237: EPSSetInitialSpace(cross->eps,-svd->nini,svd->IS);
238: SlepcBasisDestroy_Private(&svd->nini,&svd->IS);
239: }
240: EPSSetUp(cross->eps);
241: EPSGetDimensions(cross->eps,NULL,&svd->ncv,&svd->mpd);
242: EPSGetTolerances(cross->eps,NULL,&svd->max_it);
243: if (svd->tol==PETSC_DEFAULT) svd->tol = SLEPC_DEFAULT_TOL;
245: svd->leftbasis = PETSC_FALSE;
246: SVDAllocateSolution(svd,0);
247: return 0;
248: }
250: PetscErrorCode SVDSolve_Cross(SVD svd)
251: {
252: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
253: PetscInt i;
254: PetscScalar lambda;
255: PetscReal sigma;
257: EPSSolve(cross->eps);
258: EPSGetConverged(cross->eps,&svd->nconv);
259: EPSGetIterationNumber(cross->eps,&svd->its);
260: EPSGetConvergedReason(cross->eps,(EPSConvergedReason*)&svd->reason);
261: for (i=0;i<svd->nconv;i++) {
262: EPSGetEigenvalue(cross->eps,i,&lambda,NULL);
263: sigma = PetscRealPart(lambda);
264: if (svd->ishyperbolic) svd->sigma[i] = PetscSqrtReal(PetscAbsReal(sigma));
265: else {
267: if (sigma<0.0) {
268: PetscInfo(svd,"Negative eigenvalue computed by EPS: %g, resetting to 0\n",(double)sigma);
269: sigma = 0.0;
270: }
271: svd->sigma[i] = PetscSqrtReal(sigma);
272: }
273: }
274: return 0;
275: }
277: PetscErrorCode SVDComputeVectors_Cross(SVD svd)
278: {
279: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
280: PetscInt i,mloc,ploc;
281: Vec u,v,x,uv,w,omega2=NULL;
282: Mat Omega;
283: PetscScalar *dst,alpha,lambda,*varray;
284: const PetscScalar *src;
285: PetscReal nrm;
287: if (svd->isgeneralized) {
288: MatCreateVecs(svd->A,NULL,&u);
289: VecGetLocalSize(u,&mloc);
290: MatCreateVecs(svd->B,NULL,&v);
291: VecGetLocalSize(v,&ploc);
292: for (i=0;i<svd->nconv;i++) {
293: BVGetColumn(svd->V,i,&x);
294: EPSGetEigenpair(cross->eps,i,&lambda,NULL,x,NULL);
295: MatMult(svd->A,x,u); /* u_i*c_i/alpha = A*x_i */
296: VecNormalize(u,NULL);
297: MatMult(svd->B,x,v); /* v_i*s_i/alpha = B*x_i */
298: VecNormalize(v,&nrm); /* ||v||_2 = s_i/alpha */
299: alpha = 1.0/(PetscSqrtReal(1.0+PetscRealPart(lambda))*nrm); /* alpha=s_i/||v||_2 */
300: VecScale(x,alpha);
301: BVRestoreColumn(svd->V,i,&x);
302: /* copy [u;v] to U[i] */
303: BVGetColumn(svd->U,i,&uv);
304: VecGetArrayWrite(uv,&dst);
305: VecGetArrayRead(u,&src);
306: PetscArraycpy(dst,src,mloc);
307: VecRestoreArrayRead(u,&src);
308: VecGetArrayRead(v,&src);
309: PetscArraycpy(dst+mloc,src,ploc);
310: VecRestoreArrayRead(v,&src);
311: VecRestoreArrayWrite(uv,&dst);
312: BVRestoreColumn(svd->U,i,&uv);
313: }
314: VecDestroy(&v);
315: VecDestroy(&u);
316: } else if (svd->ishyperbolic && svd->swapped) { /* was solved as GHIEP, set u=Omega*u and normalize */
317: EPSGetOperators(cross->eps,NULL,&Omega);
318: MatCreateVecs(Omega,&w,NULL);
319: VecCreateSeq(PETSC_COMM_SELF,svd->ncv,&omega2);
320: VecGetArrayWrite(omega2,&varray);
321: for (i=0;i<svd->nconv;i++) {
322: BVGetColumn(svd->V,i,&v);
323: EPSGetEigenvector(cross->eps,i,v,NULL);
324: MatMult(Omega,v,w);
325: VecDot(v,w,&alpha);
326: svd->sign[i] = PetscSign(PetscRealPart(alpha));
327: varray[i] = svd->sign[i];
328: alpha = 1.0/PetscSqrtScalar(PetscAbsScalar(alpha));
329: VecScale(w,alpha);
330: VecCopy(w,v);
331: BVRestoreColumn(svd->V,i,&v);
332: }
333: BVSetSignature(svd->V,omega2);
334: VecRestoreArrayWrite(omega2,&varray);
335: VecDestroy(&omega2);
336: VecDestroy(&w);
337: SVDComputeVectors_Left(svd);
338: } else {
339: for (i=0;i<svd->nconv;i++) {
340: BVGetColumn(svd->V,i,&v);
341: EPSGetEigenvector(cross->eps,i,v,NULL);
342: BVRestoreColumn(svd->V,i,&v);
343: }
344: SVDComputeVectors_Left(svd);
345: }
346: return 0;
347: }
349: static PetscErrorCode EPSMonitor_Cross(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)
350: {
351: PetscInt i;
352: SVD svd = (SVD)ctx;
353: PetscScalar er,ei;
354: ST st;
356: EPSGetST(eps,&st);
357: for (i=0;i<PetscMin(nest,svd->ncv);i++) {
358: er = eigr[i]; ei = eigi[i];
359: STBackTransform(st,1,&er,&ei);
360: svd->sigma[i] = PetscSqrtReal(PetscAbsReal(PetscRealPart(er)));
361: svd->errest[i] = errest[i];
362: }
363: SVDMonitor(svd,its,nconv,svd->sigma,svd->errest,nest);
364: return 0;
365: }
367: PetscErrorCode SVDSetFromOptions_Cross(SVD svd,PetscOptionItems *PetscOptionsObject)
368: {
369: PetscBool set,val;
370: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
371: ST st;
373: PetscOptionsHeadBegin(PetscOptionsObject,"SVD Cross Options");
375: PetscOptionsBool("-svd_cross_explicitmatrix","Use cross explicit matrix","SVDCrossSetExplicitMatrix",cross->explicitmatrix,&val,&set);
376: if (set) SVDCrossSetExplicitMatrix(svd,val);
378: PetscOptionsHeadEnd();
380: if (!cross->eps) SVDCrossGetEPS(svd,&cross->eps);
381: if (!cross->explicitmatrix && !cross->usereps) {
382: /* use as default an ST with shell matrix and Jacobi */
383: EPSGetST(cross->eps,&st);
384: STSetMatMode(st,ST_MATMODE_SHELL);
385: }
386: EPSSetFromOptions(cross->eps);
387: return 0;
388: }
390: static PetscErrorCode SVDCrossSetExplicitMatrix_Cross(SVD svd,PetscBool explicitmatrix)
391: {
392: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
394: if (cross->explicitmatrix != explicitmatrix) {
395: cross->explicitmatrix = explicitmatrix;
396: svd->state = SVD_STATE_INITIAL;
397: }
398: return 0;
399: }
401: /*@
402: SVDCrossSetExplicitMatrix - Indicate if the eigensolver operator A^T*A must
403: be computed explicitly.
405: Logically Collective on svd
407: Input Parameters:
408: + svd - singular value solver
409: - explicitmat - boolean flag indicating if A^T*A is built explicitly
411: Options Database Key:
412: . -svd_cross_explicitmatrix <boolean> - Indicates the boolean flag
414: Level: advanced
416: .seealso: SVDCrossGetExplicitMatrix()
417: @*/
418: PetscErrorCode SVDCrossSetExplicitMatrix(SVD svd,PetscBool explicitmat)
419: {
422: PetscTryMethod(svd,"SVDCrossSetExplicitMatrix_C",(SVD,PetscBool),(svd,explicitmat));
423: return 0;
424: }
426: static PetscErrorCode SVDCrossGetExplicitMatrix_Cross(SVD svd,PetscBool *explicitmat)
427: {
428: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
430: *explicitmat = cross->explicitmatrix;
431: return 0;
432: }
434: /*@
435: SVDCrossGetExplicitMatrix - Returns the flag indicating if A^T*A is built explicitly.
437: Not Collective
439: Input Parameter:
440: . svd - singular value solver
442: Output Parameter:
443: . explicitmat - the mode flag
445: Level: advanced
447: .seealso: SVDCrossSetExplicitMatrix()
448: @*/
449: PetscErrorCode SVDCrossGetExplicitMatrix(SVD svd,PetscBool *explicitmat)
450: {
453: PetscUseMethod(svd,"SVDCrossGetExplicitMatrix_C",(SVD,PetscBool*),(svd,explicitmat));
454: return 0;
455: }
457: static PetscErrorCode SVDCrossSetEPS_Cross(SVD svd,EPS eps)
458: {
459: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
461: PetscObjectReference((PetscObject)eps);
462: EPSDestroy(&cross->eps);
463: cross->eps = eps;
464: cross->usereps = PETSC_TRUE;
465: svd->state = SVD_STATE_INITIAL;
466: return 0;
467: }
469: /*@
470: SVDCrossSetEPS - Associate an eigensolver object (EPS) to the
471: singular value solver.
473: Collective on svd
475: Input Parameters:
476: + svd - singular value solver
477: - eps - the eigensolver object
479: Level: advanced
481: .seealso: SVDCrossGetEPS()
482: @*/
483: PetscErrorCode SVDCrossSetEPS(SVD svd,EPS eps)
484: {
488: PetscTryMethod(svd,"SVDCrossSetEPS_C",(SVD,EPS),(svd,eps));
489: return 0;
490: }
492: static PetscErrorCode SVDCrossGetEPS_Cross(SVD svd,EPS *eps)
493: {
494: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
496: if (!cross->eps) {
497: EPSCreate(PetscObjectComm((PetscObject)svd),&cross->eps);
498: PetscObjectIncrementTabLevel((PetscObject)cross->eps,(PetscObject)svd,1);
499: EPSSetOptionsPrefix(cross->eps,((PetscObject)svd)->prefix);
500: EPSAppendOptionsPrefix(cross->eps,"svd_cross_");
501: PetscObjectSetOptions((PetscObject)cross->eps,((PetscObject)svd)->options);
502: EPSSetWhichEigenpairs(cross->eps,EPS_LARGEST_REAL);
503: EPSMonitorSet(cross->eps,EPSMonitor_Cross,svd,NULL);
504: }
505: *eps = cross->eps;
506: return 0;
507: }
509: /*@
510: SVDCrossGetEPS - Retrieve the eigensolver object (EPS) associated
511: to the singular value solver.
513: Not Collective
515: Input Parameter:
516: . svd - singular value solver
518: Output Parameter:
519: . eps - the eigensolver object
521: Level: advanced
523: .seealso: SVDCrossSetEPS()
524: @*/
525: PetscErrorCode SVDCrossGetEPS(SVD svd,EPS *eps)
526: {
529: PetscUseMethod(svd,"SVDCrossGetEPS_C",(SVD,EPS*),(svd,eps));
530: return 0;
531: }
533: PetscErrorCode SVDView_Cross(SVD svd,PetscViewer viewer)
534: {
535: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
536: PetscBool isascii;
538: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
539: if (isascii) {
540: if (!cross->eps) SVDCrossGetEPS(svd,&cross->eps);
541: PetscViewerASCIIPrintf(viewer," %s matrix\n",cross->explicitmatrix?"explicit":"implicit");
542: PetscViewerASCIIPushTab(viewer);
543: EPSView(cross->eps,viewer);
544: PetscViewerASCIIPopTab(viewer);
545: }
546: return 0;
547: }
549: PetscErrorCode SVDReset_Cross(SVD svd)
550: {
551: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
553: EPSReset(cross->eps);
554: MatDestroy(&cross->C);
555: MatDestroy(&cross->D);
556: return 0;
557: }
559: PetscErrorCode SVDDestroy_Cross(SVD svd)
560: {
561: SVD_CROSS *cross = (SVD_CROSS*)svd->data;
563: EPSDestroy(&cross->eps);
564: PetscFree(svd->data);
565: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetEPS_C",NULL);
566: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetEPS_C",NULL);
567: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetExplicitMatrix_C",NULL);
568: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetExplicitMatrix_C",NULL);
569: return 0;
570: }
572: SLEPC_EXTERN PetscErrorCode SVDCreate_Cross(SVD svd)
573: {
574: SVD_CROSS *cross;
576: PetscNew(&cross);
577: svd->data = (void*)cross;
579: svd->ops->solve = SVDSolve_Cross;
580: svd->ops->solveg = SVDSolve_Cross;
581: svd->ops->solveh = SVDSolve_Cross;
582: svd->ops->setup = SVDSetUp_Cross;
583: svd->ops->setfromoptions = SVDSetFromOptions_Cross;
584: svd->ops->destroy = SVDDestroy_Cross;
585: svd->ops->reset = SVDReset_Cross;
586: svd->ops->view = SVDView_Cross;
587: svd->ops->computevectors = SVDComputeVectors_Cross;
588: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetEPS_C",SVDCrossSetEPS_Cross);
589: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetEPS_C",SVDCrossGetEPS_Cross);
590: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossSetExplicitMatrix_C",SVDCrossSetExplicitMatrix_Cross);
591: PetscObjectComposeFunction((PetscObject)svd,"SVDCrossGetExplicitMatrix_C",SVDCrossGetExplicitMatrix_Cross);
592: return 0;
593: }