Actual source code: svdbasic.c
1: /*
2: The basic SVD routines, Create, View, etc. are here.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2009, Universidad Politecnica de Valencia, Spain
8: This file is part of SLEPc.
9:
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include private/svdimpl.h
26: PetscFList SVDList = 0;
27: PetscCookie SVD_COOKIE = 0;
28: PetscLogEvent SVD_SetUp = 0, SVD_Solve = 0, SVD_Dense = 0;
32: /*@C
33: SVDInitializePackage - This function initializes everything in the SVD package. It is called
34: from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to SVDCreate()
35: when using static libraries.
37: Input Parameter:
38: path - The dynamic library path, or PETSC_NULL
40: Level: developer
42: .seealso: SlepcInitialize()
43: @*/
44: PetscErrorCode SVDInitializePackage(char *path)
45: {
46: static PetscTruth initialized = PETSC_FALSE;
47: char logList[256];
48: char *className;
49: PetscTruth opt;
50: PetscErrorCode ierr;
53: if (initialized) return(0);
54: initialized = PETSC_TRUE;
55: /* Register Classes */
56: PetscCookieRegister("Singular Value Solver",&SVD_COOKIE);
57: /* Register Constructors */
58: SVDRegisterAll(path);
59: /* Register Events */
60: PetscLogEventRegister("SVDSetUp",SVD_COOKIE,&SVD_SetUp);
61: PetscLogEventRegister("SVDSolve",SVD_COOKIE,&SVD_Solve);
62: PetscLogEventRegister("SVDDense",SVD_COOKIE,&SVD_Dense);
63: /* Process info exclusions */
64: PetscOptionsGetString(PETSC_NULL, "-log_info_exclude", logList, 256, &opt);
65: if (opt) {
66: PetscStrstr(logList, "svd", &className);
67: if (className) {
68: PetscInfoDeactivateClass(SVD_COOKIE);
69: }
70: }
71: /* Process summary exclusions */
72: PetscOptionsGetString(PETSC_NULL, "-log_summary_exclude", logList, 256, &opt);
73: if (opt) {
74: PetscStrstr(logList, "svd", &className);
75: if (className) {
76: PetscLogEventDeactivateClass(SVD_COOKIE);
77: }
78: }
79: return(0);
80: }
84: /*@C
85: SVDView - Prints the SVD data structure.
87: Collective on SVD
89: Input Parameters:
90: + svd - the singular value solver context
91: - viewer - optional visualization context
93: Options Database Key:
94: . -svd_view - Calls SVDView() at end of SVDSolve()
96: Note:
97: The available visualization contexts include
98: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
99: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
100: output where only the first processor opens
101: the file. All other processors send their
102: data to the first processor to print.
104: The user can open an alternative visualization context with
105: PetscViewerASCIIOpen() - output to a specified file.
107: Level: beginner
109: .seealso: STView(), PetscViewerASCIIOpen()
110: @*/
111: PetscErrorCode SVDView(SVD svd,PetscViewer viewer)
112: {
114: const SVDType type;
115: PetscTruth isascii;
119: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(((PetscObject)svd)->comm);
123: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
124: if (isascii) {
125: PetscViewerASCIIPrintf(viewer,"SVD Object:\n");
126: SVDGetType(svd,&type);
127: if (type) {
128: PetscViewerASCIIPrintf(viewer," method: %s\n",type);
129: } else {
130: PetscViewerASCIIPrintf(viewer," method: not yet set\n");
131: }
132: switch (svd->transmode) {
133: case SVD_TRANSPOSE_EXPLICIT:
134: PetscViewerASCIIPrintf(viewer," transpose mode: explicit\n");
135: break;
136: case SVD_TRANSPOSE_IMPLICIT:
137: PetscViewerASCIIPrintf(viewer," transpose mode: implicit\n");
138: break;
139: default:
140: PetscViewerASCIIPrintf(viewer," transpose mode: not yet set\n");
141: }
142: if (svd->which == SVD_LARGEST) {
143: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: largest\n");
144: } else {
145: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: smallest\n");
146: }
147: PetscViewerASCIIPrintf(viewer," number of singular values (nsv): %d\n",svd->nsv);
148: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %d\n",svd->ncv);
149: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %d\n",svd->mpd);
150: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %d\n",svd->max_it);
151: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",svd->tol);
152: if (svd->ops->view) {
153: PetscViewerASCIIPushTab(viewer);
154: (*svd->ops->view)(svd,viewer);
155: PetscViewerASCIIPopTab(viewer);
156: }
157: PetscViewerASCIIPushTab(viewer);
158: IPView(svd->ip,viewer);
159: PetscViewerASCIIPopTab(viewer);
160: } else {
161: if (svd->ops->view) {
162: (*svd->ops->view)(svd,viewer);
163: }
164: }
165: return(0);
166: }
170: /*@C
171: SVDCreate - Creates the default SVD context.
173: Collective on MPI_Comm
175: Input Parameter:
176: . comm - MPI communicator
178: Output Parameter:
179: . svd - location to put the SVD context
181: Note:
182: The default SVD type is SVDCROSS
184: Level: beginner
186: .seealso: SVDSetUp(), SVDSolve(), SVDDestroy(), SVD
187: @*/
188: PetscErrorCode SVDCreate(MPI_Comm comm,SVD *outsvd)
189: {
191: SVD svd;
196: PetscHeaderCreate(svd,_p_SVD,struct _SVDOps,SVD_COOKIE,-1,"SVD",comm,SVDDestroy,SVDView);
197: *outsvd = svd;
199: PetscMemzero(svd->ops,sizeof(struct _SVDOps));
201: svd->OP = PETSC_NULL;
202: svd->A = PETSC_NULL;
203: svd->AT = PETSC_NULL;
204: svd->transmode = (SVDTransposeMode)PETSC_DECIDE;
205: svd->sigma = PETSC_NULL;
206: svd->perm = PETSC_NULL;
207: svd->U = PETSC_NULL;
208: svd->V = PETSC_NULL;
209: svd->vec_initial = PETSC_NULL;
210: svd->which = SVD_LARGEST;
211: svd->n = 0;
212: svd->nconv = 0;
213: svd->nsv = 1;
214: svd->ncv = 0;
215: svd->mpd = 0;
216: svd->its = 0;
217: svd->max_it = 0;
218: svd->tol = 1e-7;
219: svd->errest = PETSC_NULL;
220: svd->data = PETSC_NULL;
221: svd->setupcalled = 0;
222: svd->reason = SVD_CONVERGED_ITERATING;
223: svd->numbermonitors = 0;
224: svd->matvecs = 0;
226: IPCreate(comm,&svd->ip);
227: IPSetOptionsPrefix(svd->ip,((PetscObject)svd)->prefix);
228: IPAppendOptionsPrefix(svd->ip,"svd_");
229: PetscLogObjectParent(svd,svd->ip);
231: PetscPublishAll(svd);
232: return(0);
233: }
234:
237: /*@
238: SVDDestroy - Destroys the SVD context.
240: Collective on SVD
242: Input Parameter:
243: . svd - singular value solver context obtained from SVDCreate()
245: Level: beginner
247: .seealso: SVDCreate(), SVDSetUp(), SVDSolve()
248: @*/
249: PetscErrorCode SVDDestroy(SVD svd)
250: {
252: PetscInt i;
253: PetscScalar *p;
254:
257: if (--((PetscObject)svd)->refct > 0) return(0);
259: /* if memory was published with AMS then destroy it */
260: PetscObjectDepublish(svd);
262: if (svd->ops->destroy) {
263: (*svd->ops->destroy)(svd);
264: }
266: if (svd->OP) { MatDestroy(svd->OP); }
267: if (svd->A) { MatDestroy(svd->A); }
268: if (svd->AT) { MatDestroy(svd->AT); }
269: if (svd->n) {
270: PetscFree(svd->sigma);
271: PetscFree(svd->perm);
272: PetscFree(svd->errest);
273: if (svd->U) {
274: VecGetArray(svd->U[0],&p);
275: for (i=0;i<svd->n;i++) {
276: VecDestroy(svd->U[i]);
277: }
278: PetscFree(p);
279: PetscFree(svd->U);
280: }
281: VecGetArray(svd->V[0],&p);
282: for (i=0;i<svd->n;i++) {
283: VecDestroy(svd->V[i]);
284: }
285: PetscFree(p);
286: PetscFree(svd->V);
287: }
288: if (svd->vec_initial) { VecDestroy(svd->vec_initial); }
289: SVDMonitorCancel(svd);
290:
291: IPDestroy(svd->ip);
292:
293: PetscHeaderDestroy(svd);
294: return(0);
295: }
299: PetscErrorCode SVDDestroy_Default(SVD svd)
300: {
305: PetscFree(svd->data);
306: return(0);
307: }
311: /*@C
312: SVDSetType - Selects the particular solver to be used in the SVD object.
314: Collective on SVD
316: Input Parameters:
317: + svd - the singular value solver context
318: - type - a known method
320: Options Database Key:
321: . -svd_type <method> - Sets the method; use -help for a list
322: of available methods
323:
324: Notes:
325: See "slepc/include/slepcsvd.h" for available methods. The default
326: is SVDCROSS.
328: Normally, it is best to use the SVDSetFromOptions() command and
329: then set the SVD type from the options database rather than by using
330: this routine. Using the options database provides the user with
331: maximum flexibility in evaluating the different available methods.
332: The SVDSetType() routine is provided for those situations where it
333: is necessary to set the iterative solver independently of the command
334: line or options database.
336: Level: intermediate
338: .seealso: SVDType
339: @*/
340: PetscErrorCode SVDSetType(SVD svd,const SVDType type)
341: {
342: PetscErrorCode ierr,(*r)(SVD);
343: PetscTruth match;
349: PetscTypeCompare((PetscObject)svd,type,&match);
350: if (match) return(0);
352: if (svd->data) {
353: /* destroy the old private SVD context */
354: (*svd->ops->destroy)(svd);
355: svd->data = 0;
356: }
358: PetscFListFind(SVDList,((PetscObject)svd)->comm,type,(void (**)(void)) &r);
360: if (!r) SETERRQ1(1,"Unknown SVD type given: %s",type);
362: svd->setupcalled = 0;
363: PetscMemzero(svd->ops,sizeof(struct _SVDOps));
364: (*r)(svd);
366: PetscObjectChangeTypeName((PetscObject)svd,type);
367: return(0);
368: }
372: /*@C
373: SVDGetType - Gets the SVD type as a string from the SVD object.
375: Not Collective
377: Input Parameter:
378: . svd - the singular value solver context
380: Output Parameter:
381: . name - name of SVD method
383: Level: intermediate
385: .seealso: SVDSetType()
386: @*/
387: PetscErrorCode SVDGetType(SVD svd,const SVDType *type)
388: {
392: *type = ((PetscObject)svd)->type_name;
393: return(0);
394: }
396: /*MC
397: SVDRegisterDynamic - Adds a method to the singular value solver package.
399: Synopsis:
400: SVDRegisterDynamic(char *name_solver,char *path,char *name_create,PetscErrorCode (*routine_create)(SVD))
402: Not Collective
404: Input Parameters:
405: + name_solver - name of a new user-defined solver
406: . path - path (either absolute or relative) the library containing this solver
407: . name_create - name of routine to create the solver context
408: - routine_create - routine to create the solver context
410: Notes:
411: SVDRegisterDynamic() may be called multiple times to add several user-defined solvers.
413: If dynamic libraries are used, then the fourth input argument (routine_create)
414: is ignored.
416: Sample usage:
417: .vb
418: SVDRegisterDynamic("my_solver",/home/username/my_lib/lib/libO/solaris/mylib.a,
419: "MySolverCreate",MySolverCreate);
420: .ve
422: Then, your solver can be chosen with the procedural interface via
423: $ SVDSetType(svd,"my_solver")
424: or at runtime via the option
425: $ -svd_type my_solver
427: Level: advanced
429: Environmental variables such as ${PETSC_ARCH}, ${SLEPC_DIR},
430: and others of the form ${any_environmental_variable} occuring in pathname will be
431: replaced with appropriate values.
433: .seealso: SVDRegisterDestroy(), SVDRegisterAll()
435: M*/
439: /*@C
440: SVDRegister - See SVDRegisterDynamic()
442: Level: advanced
443: @*/
444: PetscErrorCode SVDRegister(const char *sname,const char *path,const char *name,PetscErrorCode (*function)(SVD))
445: {
447: char fullname[256];
450: PetscFListConcat(path,name,fullname);
451: PetscFListAdd(&SVDList,sname,fullname,(void (*)(void))function);
452: return(0);
453: }
457: /*@
458: SVDRegisterDestroy - Frees the list of SVD methods that were
459: registered by SVDRegisterDynamic().
461: Not Collective
463: Level: advanced
465: .seealso: SVDRegisterDynamic(), SVDRegisterAll()
466: @*/
467: PetscErrorCode SVDRegisterDestroy(void)
468: {
472: PetscFListDestroy(&SVDList);
473: SVDRegisterAll(PETSC_NULL);
474: return(0);
475: }
479: /*@
480: SVDSetIP - Associates an inner product object to the
481: singular value solver.
483: Collective on SVD
485: Input Parameters:
486: + svd - singular value solver context obtained from SVDCreate()
487: - ip - the inner product object
489: Note:
490: Use SVDGetIP() to retrieve the inner product context (for example,
491: to free it at the end of the computations).
493: Level: advanced
495: .seealso: SVDGetIP()
496: @*/
497: PetscErrorCode SVDSetIP(SVD svd,IP ip)
498: {
505: PetscObjectReference((PetscObject)ip);
506: IPDestroy(svd->ip);
507: svd->ip = ip;
508: return(0);
509: }
513: /*@C
514: SVDGetIP - Obtain the inner product object associated
515: to the singular value solver object.
517: Not Collective
519: Input Parameters:
520: . svd - singular value solver context obtained from SVDCreate()
522: Output Parameter:
523: . ip - inner product context
525: Level: advanced
527: .seealso: SVDSetIP()
528: @*/
529: PetscErrorCode SVDGetIP(SVD svd,IP *ip)
530: {
534: *ip = svd->ip;
535: return(0);
536: }