Actual source code: potentials.c

petsc-3.8.4 2018-03-24
Report Typos and Errors

  2: static char help[] = "Plots the various potentials used in the examples.\n";


  5:  #include <petscdmda.h>
  6:  #include <petscts.h>
  7:  #include <petscdraw.h>


 10: int main(int argc,char **argv)
 11: {
 12:   PetscDrawLG               lg;
 13:   PetscErrorCode            ierr;
 14:   PetscInt                  Mx = 100,i;
 15:   PetscReal                 x,hx = .1/Mx,pause,xx[3],yy[3];
 16:   PetscDraw                 draw;
 17:   const char *const         legend[] = {"(1 - u^2)^2","1 - u^2","-(1 - u)log(1 - u)"};
 18:   PetscDrawAxis             axis;
 19:   static PetscDrawViewPorts *ports = 0;


 23:   PetscInitialize(&argc,&argv,0,help);if (ierr) return ierr;
 24:   PetscViewerDrawResize(PETSC_VIEWER_DRAW_(PETSC_COMM_WORLD),1200,800);
 25:   PetscViewerDrawGetDrawLG(PETSC_VIEWER_DRAW_(PETSC_COMM_WORLD),0,&lg);
 26:   PetscDrawLGGetDraw(lg,&draw);
 27:   PetscDrawCheckResizedWindow(draw);
 28:   if (!ports) {
 29:     PetscDrawViewPortsCreateRect(draw,1,2,&ports);
 30:   }
 31:   PetscDrawLGGetAxis(lg,&axis);
 32:   PetscDrawLGReset(lg);

 34:   /*
 35:       Plot the  energies
 36:   */
 37:   PetscDrawLGSetDimension(lg,3);
 38:   PetscDrawViewPortsSet(ports,1);
 39:   x    = .9;
 40:   for (i=0; i<Mx; i++) {
 41:     xx[0] = xx[1] = xx[2] = x;
 42:     yy[0] = (1.-x*x)*(1. - x*x);
 43:     yy[1] = (1. - x*x);
 44:     yy[2] = -(1.-x)*PetscLogReal(1.-x);
 45:     PetscDrawLGAddPoint(lg,xx,yy);
 46:     x    += hx;
 47:   }
 48:   PetscDrawGetPause(draw,&pause);
 49:   PetscDrawSetPause(draw,0.0);
 50:   PetscDrawAxisSetLabels(axis,"Energy","","");
 51:   PetscDrawLGSetLegend(lg,legend);
 52:   PetscDrawLGDraw(lg);

 54:   /*
 55:       Plot the  forces
 56:   */
 57:   PetscDrawViewPortsSet(ports,0);
 58:   PetscDrawLGReset(lg);
 59:   x    = .9;
 60:   for (i=0; i<Mx; i++) {
 61:     xx[0] = xx[1] = xx[2] = x;
 62:     yy[0] = x*x*x - x;
 63:     yy[1] = -x;
 64:     yy[2] = 1.0 + PetscLogReal(1. - x);
 65:     PetscDrawLGAddPoint(lg,xx,yy);
 66:     x    += hx;
 67:   }
 68:   PetscDrawAxisSetLabels(axis,"Derivative","","");
 69:   PetscDrawLGSetLegend(lg,NULL);
 70:   PetscDrawLGDraw(lg);

 72:   PetscDrawSetPause(draw,pause);
 73:   PetscDrawPause(draw);
 74:   return(0);
 75: }