Line data Source code
1 : /*
2 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3 : SLEPc - Scalable Library for Eigenvalue Problem Computations
4 : Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5 :
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 : */
13 :
14 : #include <slepc/private/svdimpl.h> /*I "slepcsvd.h" I*/
15 :
16 : /* Logging support */
17 : PetscClassId SVD_CLASSID = 0;
18 : PetscLogEvent SVD_SetUp = 0,SVD_Solve = 0;
19 :
20 : /* List of registered SVD routines */
21 : PetscFunctionList SVDList = NULL;
22 : PetscBool SVDRegisterAllCalled = PETSC_FALSE;
23 :
24 : /* List of registered SVD monitors */
25 : PetscFunctionList SVDMonitorList = NULL;
26 : PetscFunctionList SVDMonitorCreateList = NULL;
27 : PetscFunctionList SVDMonitorDestroyList = NULL;
28 : PetscBool SVDMonitorRegisterAllCalled = PETSC_FALSE;
29 :
30 : /*@
31 : SVDCreate - Creates the default SVD context.
32 :
33 : Collective
34 :
35 : Input Parameter:
36 : . comm - MPI communicator
37 :
38 : Output Parameter:
39 : . outsvd - location to put the SVD context
40 :
41 : Note:
42 : The default SVD type is SVDCROSS
43 :
44 : Level: beginner
45 :
46 : .seealso: SVDSetUp(), SVDSolve(), SVDDestroy(), SVD
47 : @*/
48 239 : PetscErrorCode SVDCreate(MPI_Comm comm,SVD *outsvd)
49 : {
50 239 : SVD svd;
51 :
52 239 : PetscFunctionBegin;
53 239 : PetscAssertPointer(outsvd,2);
54 239 : PetscCall(SVDInitializePackage());
55 239 : PetscCall(SlepcHeaderCreate(svd,SVD_CLASSID,"SVD","Singular Value Decomposition","SVD",comm,SVDDestroy,SVDView));
56 :
57 239 : svd->OP = NULL;
58 239 : svd->OPb = NULL;
59 239 : svd->omega = NULL;
60 239 : svd->max_it = PETSC_DETERMINE;
61 239 : svd->nsv = 0;
62 239 : svd->ncv = PETSC_DETERMINE;
63 239 : svd->mpd = PETSC_DETERMINE;
64 239 : svd->nini = 0;
65 239 : svd->ninil = 0;
66 239 : svd->tol = PETSC_DETERMINE;
67 239 : svd->thres = 0.0;
68 239 : svd->threlative = PETSC_FALSE;
69 239 : svd->conv = (SVDConv)-1;
70 239 : svd->stop = SVD_STOP_BASIC;
71 239 : svd->which = SVD_LARGEST;
72 239 : svd->problem_type = (SVDProblemType)0;
73 239 : svd->impltrans = PETSC_FALSE;
74 239 : svd->trackall = PETSC_FALSE;
75 :
76 239 : svd->converged = NULL;
77 239 : svd->convergeduser = NULL;
78 239 : svd->convergeddestroy = NULL;
79 239 : svd->stopping = SVDStoppingBasic;
80 239 : svd->stoppinguser = NULL;
81 239 : svd->stoppingdestroy = NULL;
82 239 : svd->convergedctx = NULL;
83 239 : svd->stoppingctx = NULL;
84 239 : svd->numbermonitors = 0;
85 :
86 239 : svd->ds = NULL;
87 239 : svd->U = NULL;
88 239 : svd->V = NULL;
89 239 : svd->A = NULL;
90 239 : svd->B = NULL;
91 239 : svd->AT = NULL;
92 239 : svd->BT = NULL;
93 239 : svd->IS = NULL;
94 239 : svd->ISL = NULL;
95 239 : svd->sigma = NULL;
96 239 : svd->errest = NULL;
97 239 : svd->sign = NULL;
98 239 : svd->perm = NULL;
99 239 : svd->nworkl = 0;
100 239 : svd->nworkr = 0;
101 239 : svd->workl = NULL;
102 239 : svd->workr = NULL;
103 239 : svd->data = NULL;
104 :
105 239 : svd->state = SVD_STATE_INITIAL;
106 239 : svd->nconv = 0;
107 239 : svd->its = 0;
108 239 : svd->leftbasis = PETSC_FALSE;
109 239 : svd->swapped = PETSC_FALSE;
110 239 : svd->expltrans = PETSC_FALSE;
111 239 : svd->nrma = 0.0;
112 239 : svd->nrmb = 0.0;
113 239 : svd->isgeneralized = PETSC_FALSE;
114 239 : svd->reason = SVD_CONVERGED_ITERATING;
115 :
116 239 : PetscCall(PetscNew(&svd->sc));
117 239 : *outsvd = svd;
118 239 : PetscFunctionReturn(PETSC_SUCCESS);
119 : }
120 :
121 : /*@
122 : SVDReset - Resets the SVD context to the initial state (prior to setup)
123 : and destroys any allocated Vecs and Mats.
124 :
125 : Collective
126 :
127 : Input Parameter:
128 : . svd - singular value solver context obtained from SVDCreate()
129 :
130 : Level: advanced
131 :
132 : .seealso: SVDDestroy()
133 : @*/
134 245 : PetscErrorCode SVDReset(SVD svd)
135 : {
136 245 : PetscFunctionBegin;
137 245 : if (svd) PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
138 245 : if (!svd) PetscFunctionReturn(PETSC_SUCCESS);
139 245 : PetscTryTypeMethod(svd,reset);
140 245 : PetscCall(MatDestroy(&svd->OP));
141 245 : PetscCall(MatDestroy(&svd->OPb));
142 245 : PetscCall(VecDestroy(&svd->omega));
143 245 : PetscCall(MatDestroy(&svd->A));
144 245 : PetscCall(MatDestroy(&svd->B));
145 245 : PetscCall(MatDestroy(&svd->AT));
146 245 : PetscCall(MatDestroy(&svd->BT));
147 245 : PetscCall(BVDestroy(&svd->U));
148 245 : PetscCall(BVDestroy(&svd->V));
149 245 : PetscCall(VecDestroyVecs(svd->nworkl,&svd->workl));
150 245 : svd->nworkl = 0;
151 245 : PetscCall(VecDestroyVecs(svd->nworkr,&svd->workr));
152 245 : svd->nworkr = 0;
153 245 : svd->swapped = PETSC_FALSE;
154 245 : svd->state = SVD_STATE_INITIAL;
155 245 : PetscFunctionReturn(PETSC_SUCCESS);
156 : }
157 :
158 : /*@
159 : SVDDestroy - Destroys the SVD context.
160 :
161 : Collective
162 :
163 : Input Parameter:
164 : . svd - singular value solver context obtained from SVDCreate()
165 :
166 : Level: beginner
167 :
168 : .seealso: SVDCreate(), SVDSetUp(), SVDSolve()
169 : @*/
170 239 : PetscErrorCode SVDDestroy(SVD *svd)
171 : {
172 239 : PetscFunctionBegin;
173 239 : if (!*svd) PetscFunctionReturn(PETSC_SUCCESS);
174 239 : PetscValidHeaderSpecific(*svd,SVD_CLASSID,1);
175 239 : if (--((PetscObject)*svd)->refct > 0) { *svd = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
176 239 : PetscCall(SVDReset(*svd));
177 239 : PetscTryTypeMethod(*svd,destroy);
178 239 : if ((*svd)->sigma) PetscCall(PetscFree3((*svd)->sigma,(*svd)->perm,(*svd)->errest));
179 239 : if ((*svd)->sign) PetscCall(PetscFree((*svd)->sign));
180 239 : PetscCall(DSDestroy(&(*svd)->ds));
181 239 : PetscCall(PetscFree((*svd)->sc));
182 : /* just in case the initial vectors have not been used */
183 239 : PetscCall(SlepcBasisDestroy_Private(&(*svd)->nini,&(*svd)->IS));
184 239 : PetscCall(SlepcBasisDestroy_Private(&(*svd)->ninil,&(*svd)->ISL));
185 239 : if ((*svd)->convergeddestroy) PetscCall((*(*svd)->convergeddestroy)(&(*svd)->convergedctx));
186 239 : if ((*svd)->stoppingdestroy) PetscCall((*(*svd)->stoppingdestroy)(&(*svd)->stoppingctx));
187 239 : PetscCall(SVDMonitorCancel(*svd));
188 239 : PetscCall(PetscHeaderDestroy(svd));
189 239 : PetscFunctionReturn(PETSC_SUCCESS);
190 : }
191 :
192 : /*@
193 : SVDSetType - Selects the particular solver to be used in the SVD object.
194 :
195 : Logically Collective
196 :
197 : Input Parameters:
198 : + svd - the singular value solver context
199 : - type - a known method
200 :
201 : Options Database Key:
202 : . -svd_type <method> - Sets the method; use -help for a list
203 : of available methods
204 :
205 : Notes:
206 : See "slepc/include/slepcsvd.h" for available methods. The default
207 : is SVDCROSS.
208 :
209 : Normally, it is best to use the SVDSetFromOptions() command and
210 : then set the SVD type from the options database rather than by using
211 : this routine. Using the options database provides the user with
212 : maximum flexibility in evaluating the different available methods.
213 : The SVDSetType() routine is provided for those situations where it
214 : is necessary to set the iterative solver independently of the command
215 : line or options database.
216 :
217 : Level: intermediate
218 :
219 : .seealso: SVDType
220 : @*/
221 240 : PetscErrorCode SVDSetType(SVD svd,SVDType type)
222 : {
223 240 : PetscErrorCode (*r)(SVD);
224 240 : PetscBool match;
225 :
226 240 : PetscFunctionBegin;
227 240 : PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
228 240 : PetscAssertPointer(type,2);
229 :
230 240 : PetscCall(PetscObjectTypeCompare((PetscObject)svd,type,&match));
231 240 : if (match) PetscFunctionReturn(PETSC_SUCCESS);
232 :
233 240 : PetscCall(PetscFunctionListFind(SVDList,type,&r));
234 240 : PetscCheck(r,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown SVD type given: %s",type);
235 :
236 240 : PetscTryTypeMethod(svd,destroy);
237 240 : PetscCall(PetscMemzero(svd->ops,sizeof(struct _SVDOps)));
238 :
239 240 : svd->state = SVD_STATE_INITIAL;
240 240 : PetscCall(PetscObjectChangeTypeName((PetscObject)svd,type));
241 240 : PetscCall((*r)(svd));
242 240 : PetscFunctionReturn(PETSC_SUCCESS);
243 : }
244 :
245 : /*@
246 : SVDGetType - Gets the SVD type as a string from the SVD object.
247 :
248 : Not Collective
249 :
250 : Input Parameter:
251 : . svd - the singular value solver context
252 :
253 : Output Parameter:
254 : . type - name of SVD method
255 :
256 : Level: intermediate
257 :
258 : .seealso: SVDSetType()
259 : @*/
260 56 : PetscErrorCode SVDGetType(SVD svd,SVDType *type)
261 : {
262 56 : PetscFunctionBegin;
263 56 : PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
264 56 : PetscAssertPointer(type,2);
265 56 : *type = ((PetscObject)svd)->type_name;
266 56 : PetscFunctionReturn(PETSC_SUCCESS);
267 : }
268 :
269 : /*@C
270 : SVDRegister - Adds a method to the singular value solver package.
271 :
272 : Not Collective
273 :
274 : Input Parameters:
275 : + name - name of a new user-defined solver
276 : - function - routine to create the solver context
277 :
278 : Notes:
279 : SVDRegister() may be called multiple times to add several user-defined solvers.
280 :
281 : Example Usage:
282 : .vb
283 : SVDRegister("my_solver",MySolverCreate);
284 : .ve
285 :
286 : Then, your solver can be chosen with the procedural interface via
287 : $ SVDSetType(svd,"my_solver")
288 : or at runtime via the option
289 : $ -svd_type my_solver
290 :
291 : Level: advanced
292 :
293 : .seealso: SVDRegisterAll()
294 : @*/
295 1680 : PetscErrorCode SVDRegister(const char *name,PetscErrorCode (*function)(SVD))
296 : {
297 1680 : PetscFunctionBegin;
298 1680 : PetscCall(SVDInitializePackage());
299 1680 : PetscCall(PetscFunctionListAdd(&SVDList,name,function));
300 1680 : PetscFunctionReturn(PETSC_SUCCESS);
301 : }
302 :
303 : /*@C
304 : SVDMonitorRegister - Adds SVD monitor routine.
305 :
306 : Not Collective
307 :
308 : Input Parameters:
309 : + name - name of a new monitor routine
310 : . vtype - a PetscViewerType for the output
311 : . format - a PetscViewerFormat for the output
312 : . monitor - monitor routine
313 : . create - creation routine, or NULL
314 : - destroy - destruction routine, or NULL
315 :
316 : Notes:
317 : SVDMonitorRegister() may be called multiple times to add several user-defined monitors.
318 :
319 : Example Usage:
320 : .vb
321 : SVDMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
322 : .ve
323 :
324 : Then, your monitor can be chosen with the procedural interface via
325 : $ SVDMonitorSetFromOptions(svd,"-svd_monitor_my_monitor","my_monitor",NULL)
326 : or at runtime via the option
327 : $ -svd_monitor_my_monitor
328 :
329 : Level: advanced
330 :
331 : .seealso: SVDMonitorRegisterAll()
332 : @*/
333 1680 : PetscErrorCode SVDMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,PetscViewerAndFormat*),PetscErrorCode (*create)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**),PetscErrorCode (*destroy)(PetscViewerAndFormat**))
334 : {
335 1680 : char key[PETSC_MAX_PATH_LEN];
336 :
337 1680 : PetscFunctionBegin;
338 1680 : PetscCall(SVDInitializePackage());
339 1680 : PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
340 1680 : PetscCall(PetscFunctionListAdd(&SVDMonitorList,key,monitor));
341 1680 : if (create) PetscCall(PetscFunctionListAdd(&SVDMonitorCreateList,key,create));
342 1680 : if (destroy) PetscCall(PetscFunctionListAdd(&SVDMonitorDestroyList,key,destroy));
343 1680 : PetscFunctionReturn(PETSC_SUCCESS);
344 : }
345 :
346 : /*@
347 : SVDSetBV - Associates basis vectors objects to the singular value solver.
348 :
349 : Collective
350 :
351 : Input Parameters:
352 : + svd - singular value solver context obtained from SVDCreate()
353 : . V - the basis vectors object for right singular vectors
354 : - U - the basis vectors object for left singular vectors
355 :
356 : Note:
357 : Use SVDGetBV() to retrieve the basis vectors contexts (for example,
358 : to free them at the end of the computations).
359 :
360 : Level: advanced
361 :
362 : .seealso: SVDGetBV()
363 : @*/
364 3 : PetscErrorCode SVDSetBV(SVD svd,BV V,BV U)
365 : {
366 3 : PetscFunctionBegin;
367 3 : PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
368 3 : if (V) {
369 3 : PetscValidHeaderSpecific(V,BV_CLASSID,2);
370 3 : PetscCheckSameComm(svd,1,V,2);
371 3 : PetscCall(PetscObjectReference((PetscObject)V));
372 3 : PetscCall(BVDestroy(&svd->V));
373 3 : svd->V = V;
374 : }
375 3 : if (U) {
376 3 : PetscValidHeaderSpecific(U,BV_CLASSID,3);
377 3 : PetscCheckSameComm(svd,1,U,3);
378 3 : PetscCall(PetscObjectReference((PetscObject)U));
379 3 : PetscCall(BVDestroy(&svd->U));
380 3 : svd->U = U;
381 : }
382 3 : PetscFunctionReturn(PETSC_SUCCESS);
383 : }
384 :
385 : /*@
386 : SVDGetBV - Obtain the basis vectors objects associated to the singular
387 : value solver object.
388 :
389 : Not Collective
390 :
391 : Input Parameter:
392 : . svd - singular value solver context obtained from SVDCreate()
393 :
394 : Output Parameters:
395 : + V - basis vectors context for right singular vectors
396 : - U - basis vectors context for left singular vectors
397 :
398 : Level: advanced
399 :
400 : .seealso: SVDSetBV()
401 : @*/
402 482 : PetscErrorCode SVDGetBV(SVD svd,BV *V,BV *U)
403 : {
404 482 : PetscFunctionBegin;
405 482 : PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
406 482 : if (V) {
407 242 : if (!svd->V) {
408 242 : PetscCall(BVCreate(PetscObjectComm((PetscObject)svd),&svd->V));
409 242 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)svd->V,(PetscObject)svd,0));
410 242 : PetscCall(PetscObjectSetOptions((PetscObject)svd->V,((PetscObject)svd)->options));
411 : }
412 242 : *V = svd->V;
413 : }
414 482 : if (U) {
415 240 : if (!svd->U) {
416 240 : PetscCall(BVCreate(PetscObjectComm((PetscObject)svd),&svd->U));
417 240 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)svd->U,(PetscObject)svd,0));
418 240 : PetscCall(PetscObjectSetOptions((PetscObject)svd->U,((PetscObject)svd)->options));
419 : }
420 240 : *U = svd->U;
421 : }
422 482 : PetscFunctionReturn(PETSC_SUCCESS);
423 : }
424 :
425 : /*@
426 : SVDSetDS - Associates a direct solver object to the singular value solver.
427 :
428 : Collective
429 :
430 : Input Parameters:
431 : + svd - singular value solver context obtained from SVDCreate()
432 : - ds - the direct solver object
433 :
434 : Note:
435 : Use SVDGetDS() to retrieve the direct solver context (for example,
436 : to free it at the end of the computations).
437 :
438 : Level: advanced
439 :
440 : .seealso: SVDGetDS()
441 : @*/
442 0 : PetscErrorCode SVDSetDS(SVD svd,DS ds)
443 : {
444 0 : PetscFunctionBegin;
445 0 : PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
446 0 : PetscValidHeaderSpecific(ds,DS_CLASSID,2);
447 0 : PetscCheckSameComm(svd,1,ds,2);
448 0 : PetscCall(PetscObjectReference((PetscObject)ds));
449 0 : PetscCall(DSDestroy(&svd->ds));
450 0 : svd->ds = ds;
451 0 : PetscFunctionReturn(PETSC_SUCCESS);
452 : }
453 :
454 : /*@
455 : SVDGetDS - Obtain the direct solver object associated to the singular value
456 : solver object.
457 :
458 : Not Collective
459 :
460 : Input Parameters:
461 : . svd - singular value solver context obtained from SVDCreate()
462 :
463 : Output Parameter:
464 : . ds - direct solver context
465 :
466 : Level: advanced
467 :
468 : .seealso: SVDSetDS()
469 : @*/
470 239 : PetscErrorCode SVDGetDS(SVD svd,DS *ds)
471 : {
472 239 : PetscFunctionBegin;
473 239 : PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
474 239 : PetscAssertPointer(ds,2);
475 239 : if (!svd->ds) {
476 239 : PetscCall(DSCreate(PetscObjectComm((PetscObject)svd),&svd->ds));
477 239 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)svd->ds,(PetscObject)svd,0));
478 239 : PetscCall(PetscObjectSetOptions((PetscObject)svd->ds,((PetscObject)svd)->options));
479 : }
480 239 : *ds = svd->ds;
481 239 : PetscFunctionReturn(PETSC_SUCCESS);
482 : }
|
