Actual source code: svdbasic.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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: Basic SVD routines
12: */
14: #include <slepc/private/svdimpl.h> /*I "slepcsvd.h" I*/
16: PetscFunctionList SVDList = 0;
17: PetscBool SVDRegisterAllCalled = PETSC_FALSE;
18: PetscClassId SVD_CLASSID = 0;
19: PetscLogEvent SVD_SetUp = 0,SVD_Solve = 0;
21: /*@
22: SVDCreate - Creates the default SVD context.
24: Collective on MPI_Comm
26: Input Parameter:
27: . comm - MPI communicator
29: Output Parameter:
30: . svd - location to put the SVD context
32: Note:
33: The default SVD type is SVDCROSS
35: Level: beginner
37: .seealso: SVDSetUp(), SVDSolve(), SVDDestroy(), SVD
38: @*/
39: PetscErrorCode SVDCreate(MPI_Comm comm,SVD *outsvd)
40: {
42: SVD svd;
46: *outsvd = 0;
47: SVDInitializePackage();
48: SlepcHeaderCreate(svd,SVD_CLASSID,"SVD","Singular Value Decomposition","SVD",comm,SVDDestroy,SVDView);
50: svd->OP = NULL;
51: svd->max_it = 0;
52: svd->nsv = 1;
53: svd->ncv = 0;
54: svd->mpd = 0;
55: svd->nini = 0;
56: svd->ninil = 0;
57: svd->tol = PETSC_DEFAULT;
58: svd->conv = SVD_CONV_REL;
59: svd->stop = SVD_STOP_BASIC;
60: svd->which = SVD_LARGEST;
61: svd->impltrans = PETSC_FALSE;
62: svd->trackall = PETSC_FALSE;
64: svd->converged = SVDConvergedRelative;
65: svd->convergeduser = NULL;
66: svd->convergeddestroy = NULL;
67: svd->stopping = SVDStoppingBasic;
68: svd->stoppinguser = NULL;
69: svd->stoppingdestroy = NULL;
70: svd->convergedctx = NULL;
71: svd->stoppingctx = NULL;
72: svd->numbermonitors = 0;
74: svd->ds = NULL;
75: svd->U = NULL;
76: svd->V = NULL;
77: svd->A = NULL;
78: svd->AT = NULL;
79: svd->IS = NULL;
80: svd->ISL = NULL;
81: svd->sigma = NULL;
82: svd->perm = NULL;
83: svd->errest = NULL;
84: svd->data = NULL;
86: svd->state = SVD_STATE_INITIAL;
87: svd->nconv = 0;
88: svd->its = 0;
89: svd->leftbasis = PETSC_FALSE;
90: svd->reason = SVD_CONVERGED_ITERATING;
92: PetscNewLog(svd,&svd->sc);
93: *outsvd = svd;
94: return(0);
95: }
97: /*@
98: SVDReset - Resets the SVD context to the initial state (prior to setup)
99: and destroys any allocated Vecs and Mats.
101: Collective on SVD
103: Input Parameter:
104: . svd - singular value solver context obtained from SVDCreate()
106: Level: advanced
108: .seealso: SVDDestroy()
109: @*/
110: PetscErrorCode SVDReset(SVD svd)
111: {
116: if (!svd) return(0);
117: if (svd->ops->reset) { (svd->ops->reset)(svd); }
118: MatDestroy(&svd->OP);
119: MatDestroy(&svd->A);
120: MatDestroy(&svd->AT);
121: BVDestroy(&svd->U);
122: BVDestroy(&svd->V);
123: svd->state = SVD_STATE_INITIAL;
124: return(0);
125: }
127: /*@
128: SVDDestroy - Destroys the SVD context.
130: Collective on SVD
132: Input Parameter:
133: . svd - singular value solver context obtained from SVDCreate()
135: Level: beginner
137: .seealso: SVDCreate(), SVDSetUp(), SVDSolve()
138: @*/
139: PetscErrorCode SVDDestroy(SVD *svd)
140: {
144: if (!*svd) return(0);
146: if (--((PetscObject)(*svd))->refct > 0) { *svd = 0; return(0); }
147: SVDReset(*svd);
148: if ((*svd)->ops->destroy) { (*(*svd)->ops->destroy)(*svd); }
149: if ((*svd)->sigma) {
150: PetscFree3((*svd)->sigma,(*svd)->perm,(*svd)->errest);
151: }
152: DSDestroy(&(*svd)->ds);
153: PetscFree((*svd)->sc);
154: /* just in case the initial vectors have not been used */
155: SlepcBasisDestroy_Private(&(*svd)->nini,&(*svd)->IS);
156: SlepcBasisDestroy_Private(&(*svd)->ninil,&(*svd)->ISL);
157: SVDMonitorCancel(*svd);
158: PetscHeaderDestroy(svd);
159: return(0);
160: }
162: /*@C
163: SVDSetType - Selects the particular solver to be used in the SVD object.
165: Logically Collective on SVD
167: Input Parameters:
168: + svd - the singular value solver context
169: - type - a known method
171: Options Database Key:
172: . -svd_type <method> - Sets the method; use -help for a list
173: of available methods
175: Notes:
176: See "slepc/include/slepcsvd.h" for available methods. The default
177: is SVDCROSS.
179: Normally, it is best to use the SVDSetFromOptions() command and
180: then set the SVD type from the options database rather than by using
181: this routine. Using the options database provides the user with
182: maximum flexibility in evaluating the different available methods.
183: The SVDSetType() routine is provided for those situations where it
184: is necessary to set the iterative solver independently of the command
185: line or options database.
187: Level: intermediate
189: .seealso: SVDType
190: @*/
191: PetscErrorCode SVDSetType(SVD svd,SVDType type)
192: {
193: PetscErrorCode ierr,(*r)(SVD);
194: PetscBool match;
200: PetscObjectTypeCompare((PetscObject)svd,type,&match);
201: if (match) return(0);
203: PetscFunctionListFind(SVDList,type,&r);
204: if (!r) SETERRQ1(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown SVD type given: %s",type);
206: if (svd->ops->destroy) { (*svd->ops->destroy)(svd); }
207: PetscMemzero(svd->ops,sizeof(struct _SVDOps));
209: svd->state = SVD_STATE_INITIAL;
210: PetscObjectChangeTypeName((PetscObject)svd,type);
211: (*r)(svd);
212: return(0);
213: }
215: /*@C
216: SVDGetType - Gets the SVD type as a string from the SVD object.
218: Not Collective
220: Input Parameter:
221: . svd - the singular value solver context
223: Output Parameter:
224: . name - name of SVD method
226: Level: intermediate
228: .seealso: SVDSetType()
229: @*/
230: PetscErrorCode SVDGetType(SVD svd,SVDType *type)
231: {
235: *type = ((PetscObject)svd)->type_name;
236: return(0);
237: }
239: /*@C
240: SVDRegister - Adds a method to the singular value solver package.
242: Not Collective
244: Input Parameters:
245: + name - name of a new user-defined solver
246: - function - routine to create the solver context
248: Notes:
249: SVDRegister() may be called multiple times to add several user-defined solvers.
251: Sample usage:
252: .vb
253: SVDRegister("my_solver",MySolverCreate);
254: .ve
256: Then, your solver can be chosen with the procedural interface via
257: $ SVDSetType(svd,"my_solver")
258: or at runtime via the option
259: $ -svd_type my_solver
261: Level: advanced
263: .seealso: SVDRegisterAll()
264: @*/
265: PetscErrorCode SVDRegister(const char *name,PetscErrorCode (*function)(SVD))
266: {
270: SVDInitializePackage();
271: PetscFunctionListAdd(&SVDList,name,function);
272: return(0);
273: }
275: /*@
276: SVDSetBV - Associates basis vectors objects to the singular value solver.
278: Collective on SVD
280: Input Parameters:
281: + svd - singular value solver context obtained from SVDCreate()
282: . V - the basis vectors object for right singular vectors
283: - U - the basis vectors object for left singular vectors
285: Note:
286: Use SVDGetBV() to retrieve the basis vectors contexts (for example,
287: to free them at the end of the computations).
289: Level: advanced
291: .seealso: SVDGetBV()
292: @*/
293: PetscErrorCode SVDSetBV(SVD svd,BV V,BV U)
294: {
299: if (V) {
302: PetscObjectReference((PetscObject)V);
303: BVDestroy(&svd->V);
304: svd->V = V;
305: PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->V);
306: }
307: if (U) {
310: PetscObjectReference((PetscObject)U);
311: BVDestroy(&svd->U);
312: svd->U = U;
313: PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->U);
314: }
315: return(0);
316: }
318: /*@
319: SVDGetBV - Obtain the basis vectors objects associated to the singular
320: value solver object.
322: Not Collective
324: Input Parameters:
325: . svd - singular value solver context obtained from SVDCreate()
327: Output Parameter:
328: + V - basis vectors context for right singular vectors
329: - U - basis vectors context for left singular vectors
331: Level: advanced
333: .seealso: SVDSetBV()
334: @*/
335: PetscErrorCode SVDGetBV(SVD svd,BV *V,BV *U)
336: {
341: if (V) {
342: if (!svd->V) {
343: BVCreate(PetscObjectComm((PetscObject)svd),&svd->V);
344: PetscObjectIncrementTabLevel((PetscObject)svd->V,(PetscObject)svd,0);
345: PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->V);
346: PetscObjectSetOptions((PetscObject)svd->V,((PetscObject)svd)->options);
347: }
348: *V = svd->V;
349: }
350: if (U) {
351: if (!svd->U) {
352: BVCreate(PetscObjectComm((PetscObject)svd),&svd->U);
353: PetscObjectIncrementTabLevel((PetscObject)svd->U,(PetscObject)svd,0);
354: PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->U);
355: PetscObjectSetOptions((PetscObject)svd->U,((PetscObject)svd)->options);
356: }
357: *U = svd->U;
358: }
359: return(0);
360: }
362: /*@
363: SVDSetDS - Associates a direct solver object to the singular value solver.
365: Collective on SVD
367: Input Parameters:
368: + svd - singular value solver context obtained from SVDCreate()
369: - ds - the direct solver object
371: Note:
372: Use SVDGetDS() to retrieve the direct solver context (for example,
373: to free it at the end of the computations).
375: Level: advanced
377: .seealso: SVDGetDS()
378: @*/
379: PetscErrorCode SVDSetDS(SVD svd,DS ds)
380: {
387: PetscObjectReference((PetscObject)ds);
388: DSDestroy(&svd->ds);
389: svd->ds = ds;
390: PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->ds);
391: return(0);
392: }
394: /*@
395: SVDGetDS - Obtain the direct solver object associated to the singular value
396: solver object.
398: Not Collective
400: Input Parameters:
401: . svd - singular value solver context obtained from SVDCreate()
403: Output Parameter:
404: . ds - direct solver context
406: Level: advanced
408: .seealso: SVDSetDS()
409: @*/
410: PetscErrorCode SVDGetDS(SVD svd,DS *ds)
411: {
417: if (!svd->ds) {
418: DSCreate(PetscObjectComm((PetscObject)svd),&svd->ds);
419: PetscObjectIncrementTabLevel((PetscObject)svd->ds,(PetscObject)svd,0);
420: PetscLogObjectParent((PetscObject)svd,(PetscObject)svd->ds);
421: PetscObjectSetOptions((PetscObject)svd->ds,((PetscObject)svd)->options);
422: }
423: *ds = svd->ds;
424: return(0);
425: }