1: #include <petsc-private/kspimpl.h>  /*I "petscksp.h" I*/
  2: #include <petscviewerams.h>

  4: typedef struct {
  5:   PetscViewer viewer;
  6:   PetscInt    neigs;
  7:   PetscReal   *eigi;
  8:   PetscReal   *eigr;
  9:   AMS_Memory amem;
 10: } KSPMonitor_AMS;

 14: /*@C
 15:    KSPMonitorAMSCreate - create an AMS monitor context

 17:    Collective

 19:    Input Arguments:
 20: +  ksp - KSP to monitor
 21: -  amscommname - name of AMS communicator to use, if NULL uses default "PETSc" communicator

 23:    Output Arguments:
 24: .  ctx - context for monitor

 26:    Level: developer

 28: .seealso: KSPMonitorAMS(), KSPMonitorAMSDestroy()
 29: @*/
 30: PetscErrorCode KSPMonitorAMSCreate(KSP ksp,const char *amscommname,void **ctx)
 31: {
 33:   KSPMonitor_AMS *mon;

 36:   PetscNewLog(ksp,KSPMonitor_AMS,&mon);
 37:   if (!amscommname) {
 38:     mon->viewer = PETSC_VIEWER_AMS_(PetscObjectComm((PetscObject)ksp));
 39:     if (!mon->viewer) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_PLIB,"Cannot create AMS default communicator");
 40:   } else {
 41:     PetscViewerAMSOpen(PetscObjectComm((PetscObject)ksp),amscommname,&mon->viewer);
 42:   }
 43:   mon->amem = -1;
 44:   *ctx = (void*)mon;
 45:   return(0);
 46: }

 50: /*@C
 51:    KSPMonitorAMSDestroy - destroy a monitor context created with KSPMonitorAMSCreate()

 53:    Collective

 55:    Input Arguments:
 56: .  ctx - monitor context

 58:    Level: developer

 60: .seealso: KSPMonitorAMSCreate()
 61: @*/
 62: PetscErrorCode KSPMonitorAMSDestroy(void **ctx)
 63: {
 64:   KSPMonitor_AMS *mon = (KSPMonitor_AMS*)*ctx;

 68:   if (mon->amem != -1) {
 69:     PetscStackCallAMS(AMS_Memory_destroy,(mon->amem));
 70:     mon->amem = -1;
 71:   }
 72:   PetscViewerDestroy(&mon->viewer);
 73:   PetscFree(mon->eigr);
 74:   mon->eigi = NULL;
 75:   PetscFree(*ctx);
 76:   return(0);
 77: }

 81: /*@C
 82:    KSPMonitorAMS - monitor solution using AMS

 84:    Logically Collective on KSP

 86:    Input Parameters:
 87: +  ksp   - iterative context
 88: .  n     - iteration number
 89: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
 90: -  ctx -  PetscViewer of type AMS

 92:    Level: advanced

 94: .keywords: KSP, CG, monitor, AMS, singular values

 96: .seealso: KSPMonitorSingularValue(), KSPComputeExtremeSingularValues(), PetscViewerAMSOpen()
 97: @*/
 98: PetscErrorCode KSPMonitorAMS(KSP ksp,PetscInt n,PetscReal rnorm,void *ctx)
 99: {
101:   KSPMonitor_AMS *mon   = (KSPMonitor_AMS*)ctx;
102:   PetscViewer    viewer = mon->viewer;
103:   PetscReal      emax,emin;;
104:   AMS_Comm       acomm;

109:   KSPComputeExtremeSingularValues(ksp,&emax,&emin);

111:   /* UnPublish  */
112:   if (mon->amem != -1) PetscStackCallAMS(AMS_Memory_destroy,(mon->amem));
113:   mon->amem = -1;

115:   PetscFree(mon->eigr);
116:   PetscMalloc(2*n*sizeof(PetscReal),&mon->eigr);
117:   mon->eigi = mon->eigr + n;
118:   if (n) {KSPComputeEigenvalues(ksp,n,mon->eigr,mon->eigi,&mon->neigs);}

120:   PetscViewerAMSGetAMSComm(viewer,&acomm);
121:   PetscStackCallAMS(AMS_Memory_create,(acomm,"ksp_monitor_ams",&mon->amem));
122:   PetscStackCallAMS(AMS_Memory_take_access,(mon->amem));

124:   PetscStackCallAMS(AMS_Memory_add_field,(mon->amem,"rnorm",&ksp->rnorm,1,AMS_DOUBLE,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF));
125:   PetscStackCallAMS(AMS_Memory_add_field,(mon->amem,"neigs",&mon->neigs,1,AMS_INT,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF));
126:   if (mon->neigs > 0) {
127:     PetscStackCallAMS(AMS_Memory_add_field,(mon->amem,"eigr",&mon->eigr,mon->neigs,AMS_DOUBLE,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF));
128:     PetscStackCallAMS(AMS_Memory_add_field,(mon->amem,"eigi",&mon->eigr,mon->neigs,AMS_DOUBLE,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF));
129:   }
130:   PetscStackCallAMS(AMS_Memory_publish,(mon->amem));
131:   PetscStackCallAMS(AMS_Memory_grant_access,(mon->amem));

133:   PetscInfo2(ksp,"KSP extreme singular values min=%G max=%G\n",emin,emax);
134:   return(0);
135: }