Actual source code: potentials.c
petsc-3.12.2 2019-11-22
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: PetscDrawViewPorts *ports;
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: PetscDrawViewPortsCreateRect(draw,1,2,&ports);
29: PetscDrawLGGetAxis(lg,&axis);
30: PetscDrawLGReset(lg);
32: /*
33: Plot the energies
34: */
35: PetscDrawLGSetDimension(lg,3);
36: PetscDrawViewPortsSet(ports,1);
37: x = .9;
38: for (i=0; i<Mx; i++) {
39: xx[0] = xx[1] = xx[2] = x;
40: yy[0] = (1.-x*x)*(1. - x*x);
41: yy[1] = (1. - x*x);
42: yy[2] = -(1.-x)*PetscLogReal(1.-x);
43: PetscDrawLGAddPoint(lg,xx,yy);
44: x += hx;
45: }
46: PetscDrawGetPause(draw,&pause);
47: PetscDrawSetPause(draw,0.0);
48: PetscDrawAxisSetLabels(axis,"Energy","","");
49: PetscDrawLGSetLegend(lg,legend);
50: PetscDrawLGDraw(lg);
52: /*
53: Plot the forces
54: */
55: PetscDrawViewPortsSet(ports,0);
56: PetscDrawLGReset(lg);
57: x = .9;
58: for (i=0; i<Mx; i++) {
59: xx[0] = xx[1] = xx[2] = x;
60: yy[0] = x*x*x - x;
61: yy[1] = -x;
62: yy[2] = 1.0 + PetscLogReal(1. - x);
63: PetscDrawLGAddPoint(lg,xx,yy);
64: x += hx;
65: }
66: PetscDrawAxisSetLabels(axis,"Derivative","","");
67: PetscDrawLGSetLegend(lg,NULL);
68: PetscDrawLGDraw(lg);
70: PetscDrawSetPause(draw,pause);
71: PetscDrawPause(draw);
72: PetscDrawViewPortsDestroy(ports);
73: PetscFinalize();
74: return ierr;
75: }
77: /*TEST
79: test:
80: requires: x
82: TEST*/