Actual source code: svdview.c

slepc-main 2024-11-15
Report Typos and Errors
  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:    The SVD routines related to various viewers
 12: */

 14: #include <slepc/private/svdimpl.h>
 15: #include <petscdraw.h>

 17: /*@
 18:    SVDView - Prints the SVD data structure.

 20:    Collective

 22:    Input Parameters:
 23: +  svd - the singular value solver context
 24: -  viewer - optional visualization context

 26:    Options Database Key:
 27: .  -svd_view -  Calls SVDView() at end of SVDSolve()

 29:    Note:
 30:    The available visualization contexts include
 31: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
 32: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
 33:          output where only the first processor opens
 34:          the file.  All other processors send their
 35:          data to the first processor to print.

 37:    The user can open an alternative visualization context with
 38:    PetscViewerASCIIOpen() - output to a specified file.

 40:    Level: beginner

 42: .seealso: EPSView()
 43: @*/
 44: PetscErrorCode SVDView(SVD svd,PetscViewer viewer)
 45: {
 46:   const char     *type=NULL;
 47:   PetscBool      isascii,isshell,isexternal;

 49:   PetscFunctionBegin;
 51:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
 53:   PetscCheckSameComm(svd,1,viewer,2);

 55:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
 56:   if (isascii) {
 57:     PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)svd,viewer));
 58:     PetscCall(PetscViewerASCIIPushTab(viewer));
 59:     PetscTryTypeMethod(svd,view,viewer);
 60:     PetscCall(PetscViewerASCIIPopTab(viewer));
 61:     if (svd->problem_type) {
 62:       switch (svd->problem_type) {
 63:         case SVD_STANDARD:    type = "(standard) singular value problem"; break;
 64:         case SVD_GENERALIZED: type = "generalized singular value problem"; break;
 65:         case SVD_HYPERBOLIC:  type = "hyperbolic singular value problem"; break;
 66:       }
 67:     } else type = "not yet set";
 68:     PetscCall(PetscViewerASCIIPrintf(viewer,"  problem type: %s\n",type));
 69:     PetscCall(PetscViewerASCIIPrintf(viewer,"  transpose mode: %s\n",svd->impltrans?"implicit":"explicit"));
 70:     if (svd->which == SVD_LARGEST) PetscCall(PetscViewerASCIIPrintf(viewer,"  selected portion of the spectrum: largest\n"));
 71:     else PetscCall(PetscViewerASCIIPrintf(viewer,"  selected portion of the spectrum: smallest\n"));
 72:     PetscCall(PetscViewerASCIIPrintf(viewer,"  number of singular values (nsv): %" PetscInt_FMT "\n",svd->nsv));
 73:     PetscCall(PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %" PetscInt_FMT "\n",svd->ncv));
 74:     PetscCall(PetscViewerASCIIPrintf(viewer,"  maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",svd->mpd));
 75:     PetscCall(PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %" PetscInt_FMT "\n",svd->max_it));
 76:     PetscCall(PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)svd->tol));
 77:     PetscCall(PetscViewerASCIIPrintf(viewer,"  convergence test: "));
 78:     PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
 79:     switch (svd->conv) {
 80:     case SVD_CONV_ABS:
 81:       PetscCall(PetscViewerASCIIPrintf(viewer,"absolute\n"));break;
 82:     case SVD_CONV_REL:
 83:       PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the singular value\n"));break;
 84:     case SVD_CONV_NORM:
 85:       PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n"));
 86:       PetscCall(PetscViewerASCIIPrintf(viewer,"  computed matrix norms: norm(A)=%g",(double)svd->nrma));
 87:       if (svd->isgeneralized) PetscCall(PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)svd->nrmb));
 88:       PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
 89:       break;
 90:     case SVD_CONV_MAXIT:
 91:       PetscCall(PetscViewerASCIIPrintf(viewer,"maximum number of iterations\n"));break;
 92:     case SVD_CONV_USER:
 93:       PetscCall(PetscViewerASCIIPrintf(viewer,"user-defined\n"));break;
 94:     }
 95:     PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
 96:     if (svd->nini) PetscCall(PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(svd->nini)));
 97:     if (svd->ninil) PetscCall(PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial left space: %" PetscInt_FMT "\n",PetscAbs(svd->ninil)));
 98:   } else PetscTryTypeMethod(svd,view,viewer);
 99:   PetscCall(PetscObjectTypeCompareAny((PetscObject)svd,&isshell,SVDCROSS,SVDCYCLIC,""));
100:   PetscCall(PetscObjectTypeCompareAny((PetscObject)svd,&isexternal,SVDSCALAPACK,SVDKSVD,SVDELEMENTAL,SVDPRIMME,""));
101:   if (!isshell && !isexternal) {
102:     PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
103:     if (!svd->V) PetscCall(SVDGetBV(svd,&svd->V,NULL));
104:     PetscCall(BVView(svd->V,viewer));
105:     if (!svd->ds) PetscCall(SVDGetDS(svd,&svd->ds));
106:     PetscCall(DSView(svd->ds,viewer));
107:     PetscCall(PetscViewerPopFormat(viewer));
108:   }
109:   PetscFunctionReturn(PETSC_SUCCESS);
110: }

112: /*@
113:    SVDViewFromOptions - View from options

115:    Collective

117:    Input Parameters:
118: +  svd  - the singular value solver context
119: .  obj  - optional object
120: -  name - command line option

122:    Level: intermediate

124: .seealso: SVDView(), SVDCreate()
125: @*/
126: PetscErrorCode SVDViewFromOptions(SVD svd,PetscObject obj,const char name[])
127: {
128:   PetscFunctionBegin;
130:   PetscCall(PetscObjectViewFromOptions((PetscObject)svd,obj,name));
131:   PetscFunctionReturn(PETSC_SUCCESS);
132: }

134: /*@
135:    SVDConvergedReasonView - Displays the reason an SVD solve converged or diverged.

137:    Collective

139:    Input Parameters:
140: +  svd - the singular value solver context
141: -  viewer - the viewer to display the reason

143:    Options Database Keys:
144: .  -svd_converged_reason - print reason for convergence, and number of iterations

146:    Note:
147:    To change the format of the output call PetscViewerPushFormat(viewer,format) before
148:    this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
149:    display a reason if it fails. The latter can be set in the command line with
150:    -svd_converged_reason ::failed

152:    Level: intermediate

154: .seealso: SVDSetTolerances(), SVDGetIterationNumber(), SVDConvergedReasonViewFromOptions()
155: @*/
156: PetscErrorCode SVDConvergedReasonView(SVD svd,PetscViewer viewer)
157: {
158:   PetscBool         isAscii;
159:   PetscViewerFormat format;

161:   PetscFunctionBegin;
162:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)svd));
163:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
164:   if (isAscii) {
165:     PetscCall(PetscViewerGetFormat(viewer,&format));
166:     PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
167:     if (svd->reason > 0 && format != PETSC_VIEWER_FAILED) PetscCall(PetscViewerASCIIPrintf(viewer,"%s SVD solve converged (%" PetscInt_FMT " singular triplet%s) due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)svd)->prefix?((PetscObject)svd)->prefix:"",svd->nconv,(svd->nconv>1)?"s":"",SVDConvergedReasons[svd->reason],svd->its));
168:     else if (svd->reason <= 0) PetscCall(PetscViewerASCIIPrintf(viewer,"%s SVD solve did not converge due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)svd)->prefix?((PetscObject)svd)->prefix:"",SVDConvergedReasons[svd->reason],svd->its));
169:     PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
170:   }
171:   PetscFunctionReturn(PETSC_SUCCESS);
172: }

174: /*@
175:    SVDConvergedReasonViewFromOptions - Processes command line options to determine if/how
176:    the SVD converged reason is to be viewed.

178:    Collective

180:    Input Parameter:
181: .  svd - the singular value solver context

183:    Level: developer

185: .seealso: SVDConvergedReasonView()
186: @*/
187: PetscErrorCode SVDConvergedReasonViewFromOptions(SVD svd)
188: {
189:   PetscViewer       viewer;
190:   PetscBool         flg;
191:   static PetscBool  incall = PETSC_FALSE;
192:   PetscViewerFormat format;

194:   PetscFunctionBegin;
195:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
196:   incall = PETSC_TRUE;
197:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_converged_reason",&viewer,&format,&flg));
198:   if (flg) {
199:     PetscCall(PetscViewerPushFormat(viewer,format));
200:     PetscCall(SVDConvergedReasonView(svd,viewer));
201:     PetscCall(PetscViewerPopFormat(viewer));
202:     PetscCall(PetscViewerDestroy(&viewer));
203:   }
204:   incall = PETSC_FALSE;
205:   PetscFunctionReturn(PETSC_SUCCESS);
206: }

208: static PetscErrorCode SVDErrorView_ASCII(SVD svd,SVDErrorType etype,PetscViewer viewer)
209: {
210:   PetscReal      error,sigma;
211:   PetscInt       i,j;

213:   PetscFunctionBegin;
214:   if (svd->nconv<svd->nsv) {
215:     PetscCall(PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " singular values converged\n\n",svd->nsv));
216:     PetscFunctionReturn(PETSC_SUCCESS);
217:   }
218:   for (i=0;i<svd->nsv;i++) {
219:     PetscCall(SVDComputeError(svd,i,etype,&error));
220:     if (error>=5.0*svd->tol) {
221:       PetscCall(PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",svd->nsv));
222:       PetscFunctionReturn(PETSC_SUCCESS);
223:     }
224:   }
225:   PetscCall(PetscViewerASCIIPrintf(viewer," All requested %ssingular values computed up to the required tolerance:",svd->isgeneralized?"generalized ":""));
226:   for (i=0;i<=(svd->nsv-1)/8;i++) {
227:     PetscCall(PetscViewerASCIIPrintf(viewer,"\n     "));
228:     for (j=0;j<PetscMin(8,svd->nsv-8*i);j++) {
229:       PetscCall(SVDGetSingularTriplet(svd,8*i+j,&sigma,NULL,NULL));
230:       PetscCall(PetscViewerASCIIPrintf(viewer,"%.5f",(double)sigma));
231:       if (8*i+j+1<svd->nsv) PetscCall(PetscViewerASCIIPrintf(viewer,", "));
232:     }
233:   }
234:   PetscCall(PetscViewerASCIIPrintf(viewer,"\n\n"));
235:   PetscFunctionReturn(PETSC_SUCCESS);
236: }

238: static PetscErrorCode SVDErrorView_DETAIL(SVD svd,SVDErrorType etype,PetscViewer viewer)
239: {
240:   PetscReal      error,sigma;
241:   PetscInt       i;
242:   char           ex[30],sep[]=" ---------------------- --------------------\n";

244:   PetscFunctionBegin;
245:   if (!svd->nconv) PetscFunctionReturn(PETSC_SUCCESS);
246:   switch (etype) {
247:     case SVD_ERROR_ABSOLUTE:
248:       PetscCall(PetscSNPrintf(ex,sizeof(ex)," absolute error"));
249:       break;
250:     case SVD_ERROR_RELATIVE:
251:       PetscCall(PetscSNPrintf(ex,sizeof(ex)," relative error"));
252:       break;
253:     case SVD_ERROR_NORM:
254:       if (svd->isgeneralized) PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||r||/||[A;B]||"));
255:       else PetscCall(PetscSNPrintf(ex,sizeof(ex),"  ||r||/||A||"));
256:       break;
257:   }
258:   PetscCall(PetscViewerASCIIPrintf(viewer,"%s          sigma           %s\n%s",sep,ex,sep));
259:   for (i=0;i<svd->nconv;i++) {
260:     PetscCall(SVDGetSingularTriplet(svd,i,&sigma,NULL,NULL));
261:     PetscCall(SVDComputeError(svd,i,etype,&error));
262:     PetscCall(PetscViewerASCIIPrintf(viewer,"       % 6f          %12g\n",(double)sigma,(double)error));
263:   }
264:   PetscCall(PetscViewerASCIIPrintf(viewer,"%s",sep));
265:   PetscFunctionReturn(PETSC_SUCCESS);
266: }

268: static PetscErrorCode SVDErrorView_MATLAB(SVD svd,SVDErrorType etype,PetscViewer viewer)
269: {
270:   PetscReal      error;
271:   PetscInt       i;
272:   const char     *name;

274:   PetscFunctionBegin;
275:   PetscCall(PetscObjectGetName((PetscObject)svd,&name));
276:   PetscCall(PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name));
277:   for (i=0;i<svd->nconv;i++) {
278:     PetscCall(SVDComputeError(svd,i,etype,&error));
279:     PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error));
280:   }
281:   PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
282:   PetscFunctionReturn(PETSC_SUCCESS);
283: }

285: /*@
286:    SVDErrorView - Displays the errors associated with the computed solution
287:    (as well as the singular values).

289:    Collective

291:    Input Parameters:
292: +  svd    - the singular value solver context
293: .  etype  - error type
294: -  viewer - optional visualization context

296:    Options Database Keys:
297: +  -svd_error_absolute - print absolute errors of each singular triplet
298: .  -svd_error_relative - print relative errors of each singular triplet
299: -  -svd_error_norm     - print errors relative to the matrix norms of each singular triplet

301:    Notes:
302:    By default, this function checks the error of all singular triplets and prints
303:    the singular values if all of them are below the requested tolerance.
304:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
305:    singular values and corresponding errors is printed.

307:    Level: intermediate

309: .seealso: SVDSolve(), SVDValuesView(), SVDVectorsView()
310: @*/
311: PetscErrorCode SVDErrorView(SVD svd,SVDErrorType etype,PetscViewer viewer)
312: {
313:   PetscBool         isascii;
314:   PetscViewerFormat format;

316:   PetscFunctionBegin;
318:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
320:   PetscCheckSameComm(svd,1,viewer,3);
321:   SVDCheckSolved(svd,1);
322:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
323:   if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);

325:   PetscCall(PetscViewerGetFormat(viewer,&format));
326:   switch (format) {
327:     case PETSC_VIEWER_DEFAULT:
328:     case PETSC_VIEWER_ASCII_INFO:
329:       PetscCall(SVDErrorView_ASCII(svd,etype,viewer));
330:       break;
331:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
332:       PetscCall(SVDErrorView_DETAIL(svd,etype,viewer));
333:       break;
334:     case PETSC_VIEWER_ASCII_MATLAB:
335:       PetscCall(SVDErrorView_MATLAB(svd,etype,viewer));
336:       break;
337:     default:
338:       PetscCall(PetscInfo(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
339:   }
340:   PetscFunctionReturn(PETSC_SUCCESS);
341: }

343: /*@
344:    SVDErrorViewFromOptions - Processes command line options to determine if/how
345:    the errors of the computed solution are to be viewed.

347:    Collective

349:    Input Parameter:
350: .  svd - the singular value solver context

352:    Level: developer

354: .seealso: SVDErrorView()
355: @*/
356: PetscErrorCode SVDErrorViewFromOptions(SVD svd)
357: {
358:   PetscViewer       viewer;
359:   PetscBool         flg;
360:   static PetscBool  incall = PETSC_FALSE;
361:   PetscViewerFormat format;

363:   PetscFunctionBegin;
364:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
365:   incall = PETSC_TRUE;
366:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_absolute",&viewer,&format,&flg));
367:   if (flg) {
368:     PetscCall(PetscViewerPushFormat(viewer,format));
369:     PetscCall(SVDErrorView(svd,SVD_ERROR_ABSOLUTE,viewer));
370:     PetscCall(PetscViewerPopFormat(viewer));
371:     PetscCall(PetscViewerDestroy(&viewer));
372:   }
373:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_relative",&viewer,&format,&flg));
374:   if (flg) {
375:     PetscCall(PetscViewerPushFormat(viewer,format));
376:     PetscCall(SVDErrorView(svd,SVD_ERROR_RELATIVE,viewer));
377:     PetscCall(PetscViewerPopFormat(viewer));
378:     PetscCall(PetscViewerDestroy(&viewer));
379:   }
380:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_norm",&viewer,&format,&flg));
381:   if (flg) {
382:     PetscCall(PetscViewerPushFormat(viewer,format));
383:     PetscCall(SVDErrorView(svd,SVD_ERROR_NORM,viewer));
384:     PetscCall(PetscViewerPopFormat(viewer));
385:     PetscCall(PetscViewerDestroy(&viewer));
386:   }
387:   incall = PETSC_FALSE;
388:   PetscFunctionReturn(PETSC_SUCCESS);
389: }

391: static PetscErrorCode SVDValuesView_DRAW(SVD svd,PetscViewer viewer)
392: {
393:   PetscDraw      draw;
394:   PetscDrawSP    drawsp;
395:   PetscReal      re,im=0.0;
396:   PetscInt       i;

398:   PetscFunctionBegin;
399:   if (!svd->nconv) PetscFunctionReturn(PETSC_SUCCESS);
400:   PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
401:   PetscCall(PetscDrawSetTitle(draw,"Computed singular values"));
402:   PetscCall(PetscDrawSPCreate(draw,1,&drawsp));
403:   for (i=0;i<svd->nconv;i++) {
404:     re = svd->sigma[svd->perm[i]];
405:     PetscCall(PetscDrawSPAddPoint(drawsp,&re,&im));
406:   }
407:   PetscCall(PetscDrawSPDraw(drawsp,PETSC_TRUE));
408:   PetscCall(PetscDrawSPSave(drawsp));
409:   PetscCall(PetscDrawSPDestroy(&drawsp));
410:   PetscFunctionReturn(PETSC_SUCCESS);
411: }

413: static PetscErrorCode SVDValuesView_BINARY(SVD svd,PetscViewer viewer)
414: {
415:   PetscInt       i,k;
416:   PetscReal      *sv;

418:   PetscFunctionBegin;
419:   PetscCall(PetscMalloc1(svd->nconv,&sv));
420:   for (i=0;i<svd->nconv;i++) {
421:     k = svd->perm[i];
422:     sv[i] = svd->sigma[k];
423:   }
424:   PetscCall(PetscViewerBinaryWrite(viewer,sv,svd->nconv,PETSC_REAL));
425:   PetscCall(PetscFree(sv));
426:   PetscFunctionReturn(PETSC_SUCCESS);
427: }

429: #if defined(PETSC_HAVE_HDF5)
430: static PetscErrorCode SVDValuesView_HDF5(SVD svd,PetscViewer viewer)
431: {
432:   PetscInt       i,k,n,N;
433:   PetscMPIInt    rank;
434:   Vec            v;
435:   char           vname[30];
436:   const char     *ename;

438:   PetscFunctionBegin;
439:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)svd),&rank));
440:   N = svd->nconv;
441:   n = rank? 0: N;
442:   /* create a vector containing the singular values */
443:   PetscCall(VecCreateMPI(PetscObjectComm((PetscObject)svd),n,N,&v));
444:   PetscCall(PetscObjectGetName((PetscObject)svd,&ename));
445:   PetscCall(PetscSNPrintf(vname,sizeof(vname),"sigma_%s",ename));
446:   PetscCall(PetscObjectSetName((PetscObject)v,vname));
447:   if (!rank) {
448:     for (i=0;i<svd->nconv;i++) {
449:       k = svd->perm[i];
450:       PetscCall(VecSetValue(v,i,svd->sigma[k],INSERT_VALUES));
451:     }
452:   }
453:   PetscCall(VecAssemblyBegin(v));
454:   PetscCall(VecAssemblyEnd(v));
455:   PetscCall(VecView(v,viewer));
456:   PetscCall(VecDestroy(&v));
457:   PetscFunctionReturn(PETSC_SUCCESS);
458: }
459: #endif

461: static PetscErrorCode SVDValuesView_ASCII(SVD svd,PetscViewer viewer)
462: {
463:   PetscInt       i;

465:   PetscFunctionBegin;
466:   PetscCall(PetscViewerASCIIPrintf(viewer,"Singular values = \n"));
467:   for (i=0;i<svd->nconv;i++) PetscCall(PetscViewerASCIIPrintf(viewer,"   %.5f\n",(double)svd->sigma[svd->perm[i]]));
468:   PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
469:   PetscFunctionReturn(PETSC_SUCCESS);
470: }

472: static PetscErrorCode SVDValuesView_MATLAB(SVD svd,PetscViewer viewer)
473: {
474:   PetscInt       i;
475:   const char     *name;

477:   PetscFunctionBegin;
478:   PetscCall(PetscObjectGetName((PetscObject)svd,&name));
479:   PetscCall(PetscViewerASCIIPrintf(viewer,"Sigma_%s = [\n",name));
480:   for (i=0;i<svd->nconv;i++) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)svd->sigma[svd->perm[i]]));
481:   PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
482:   PetscFunctionReturn(PETSC_SUCCESS);
483: }

485: /*@
486:    SVDValuesView - Displays the computed singular values in a viewer.

488:    Collective

490:    Input Parameters:
491: +  svd    - the singular value solver context
492: -  viewer - the viewer

494:    Options Database Key:
495: .  -svd_view_values - print computed singular values

497:    Level: intermediate

499: .seealso: SVDSolve(), SVDVectorsView(), SVDErrorView()
500: @*/
501: PetscErrorCode SVDValuesView(SVD svd,PetscViewer viewer)
502: {
503:   PetscBool         isascii,isdraw,isbinary;
504:   PetscViewerFormat format;
505: #if defined(PETSC_HAVE_HDF5)
506:   PetscBool         ishdf5;
507: #endif

509:   PetscFunctionBegin;
511:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
513:   PetscCheckSameComm(svd,1,viewer,2);
514:   SVDCheckSolved(svd,1);
515:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
516:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
517: #if defined(PETSC_HAVE_HDF5)
518:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
519: #endif
520:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
521:   if (isdraw) PetscCall(SVDValuesView_DRAW(svd,viewer));
522:   else if (isbinary) PetscCall(SVDValuesView_BINARY(svd,viewer));
523: #if defined(PETSC_HAVE_HDF5)
524:   else if (ishdf5) PetscCall(SVDValuesView_HDF5(svd,viewer));
525: #endif
526:   else if (isascii) {
527:     PetscCall(PetscViewerGetFormat(viewer,&format));
528:     switch (format) {
529:       case PETSC_VIEWER_DEFAULT:
530:       case PETSC_VIEWER_ASCII_INFO:
531:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
532:         PetscCall(SVDValuesView_ASCII(svd,viewer));
533:         break;
534:       case PETSC_VIEWER_ASCII_MATLAB:
535:         PetscCall(SVDValuesView_MATLAB(svd,viewer));
536:         break;
537:       default:
538:         PetscCall(PetscInfo(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
539:     }
540:   }
541:   PetscFunctionReturn(PETSC_SUCCESS);
542: }

544: /*@
545:    SVDValuesViewFromOptions - Processes command line options to determine if/how
546:    the computed singular values are to be viewed.

548:    Collective

550:    Input Parameter:
551: .  svd - the singular value solver context

553:    Level: developer

555: .seealso: SVDValuesView()
556: @*/
557: PetscErrorCode SVDValuesViewFromOptions(SVD svd)
558: {
559:   PetscViewer       viewer;
560:   PetscBool         flg;
561:   static PetscBool  incall = PETSC_FALSE;
562:   PetscViewerFormat format;

564:   PetscFunctionBegin;
565:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
566:   incall = PETSC_TRUE;
567:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_values",&viewer,&format,&flg));
568:   if (flg) {
569:     PetscCall(PetscViewerPushFormat(viewer,format));
570:     PetscCall(SVDValuesView(svd,viewer));
571:     PetscCall(PetscViewerPopFormat(viewer));
572:     PetscCall(PetscViewerDestroy(&viewer));
573:   }
574:   incall = PETSC_FALSE;
575:   PetscFunctionReturn(PETSC_SUCCESS);
576: }

578: /*@
579:    SVDVectorsView - Outputs computed singular vectors to a viewer.

581:    Collective

583:    Input Parameters:
584: +  svd    - the singular value solver context
585: -  viewer - the viewer

587:    Options Database Key:
588: .  -svd_view_vectors - output singular vectors

590:    Note:
591:    Right and left singular vectors are interleaved, that is, the vectors are
592:    output in the following order V0, U0, V1, U1, V2, U2, ...

594:    Level: intermediate

596: .seealso: SVDSolve(), SVDValuesView(), SVDErrorView()
597: @*/
598: PetscErrorCode SVDVectorsView(SVD svd,PetscViewer viewer)
599: {
600:   PetscInt       i,k;
601:   Vec            x;
602:   char           vname[30];
603:   const char     *ename;

605:   PetscFunctionBegin;
607:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
609:   PetscCheckSameComm(svd,1,viewer,2);
610:   SVDCheckSolved(svd,1);
611:   if (svd->nconv) {
612:     PetscCall(PetscObjectGetName((PetscObject)svd,&ename));
613:     PetscCall(SVDComputeVectors(svd));
614:     for (i=0;i<svd->nconv;i++) {
615:       k = svd->perm[i];
616:       PetscCall(PetscSNPrintf(vname,sizeof(vname),"V%" PetscInt_FMT "_%s",i,ename));
617:       PetscCall(BVGetColumn(svd->V,k,&x));
618:       PetscCall(PetscObjectSetName((PetscObject)x,vname));
619:       PetscCall(VecView(x,viewer));
620:       PetscCall(BVRestoreColumn(svd->V,k,&x));
621:       PetscCall(PetscSNPrintf(vname,sizeof(vname),"U%" PetscInt_FMT "_%s",i,ename));
622:       PetscCall(BVGetColumn(svd->U,k,&x));
623:       PetscCall(PetscObjectSetName((PetscObject)x,vname));
624:       PetscCall(VecView(x,viewer));
625:       PetscCall(BVRestoreColumn(svd->U,k,&x));
626:     }
627:   }
628:   PetscFunctionReturn(PETSC_SUCCESS);
629: }

631: /*@
632:    SVDVectorsViewFromOptions - Processes command line options to determine if/how
633:    the computed singular vectors are to be viewed.

635:    Collective

637:    Input Parameter:
638: .  svd - the singular value solver context

640:    Level: developer

642: .seealso: SVDVectorsView()
643: @*/
644: PetscErrorCode SVDVectorsViewFromOptions(SVD svd)
645: {
646:   PetscViewer       viewer;
647:   PetscBool         flg = PETSC_FALSE;
648:   static PetscBool  incall = PETSC_FALSE;
649:   PetscViewerFormat format;

651:   PetscFunctionBegin;
652:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
653:   incall = PETSC_TRUE;
654:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_vectors",&viewer,&format,&flg));
655:   if (flg) {
656:     PetscCall(PetscViewerPushFormat(viewer,format));
657:     PetscCall(SVDVectorsView(svd,viewer));
658:     PetscCall(PetscViewerPopFormat(viewer));
659:     PetscCall(PetscViewerDestroy(&viewer));
660:   }
661:   incall = PETSC_FALSE;
662:   PetscFunctionReturn(PETSC_SUCCESS);
663: }