Actual source code: potentials.c
petsc-3.7.5 2017-01-01
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>
12: int main(int argc,char **argv)
13: {
14: PetscDrawLG lg;
15: PetscErrorCode ierr;
16: PetscInt Mx = 100,i;
17: PetscReal x,hx = .1/Mx,pause,xx[3],yy[3];
18: PetscDraw draw;
19: const char *const legend[] = {"(1 - u^2)^2","1 - u^2","-(1 - u)log(1 - u)"};
20: PetscDrawAxis axis;
21: static PetscDrawViewPorts *ports = 0;
25: PetscInitialize(&argc,&argv,0,help);
26: PetscViewerDrawResize(PETSC_VIEWER_DRAW_(PETSC_COMM_WORLD),1200,800);
27: PetscViewerDrawGetDrawLG(PETSC_VIEWER_DRAW_(PETSC_COMM_WORLD),0,&lg);
28: PetscDrawLGGetDraw(lg,&draw);
29: PetscDrawCheckResizedWindow(draw);
30: if (!ports) {
31: PetscDrawViewPortsCreateRect(draw,1,2,&ports);
32: }
33: PetscDrawLGGetAxis(lg,&axis);
34: PetscDrawLGReset(lg);
36: /*
37: Plot the energies
38: */
39: PetscDrawLGSetDimension(lg,3);
40: PetscDrawViewPortsSet(ports,1);
41: x = .9;
42: for (i=0; i<Mx; i++) {
43: xx[0] = xx[1] = xx[2] = x;
44: yy[0] = (1.-x*x)*(1. - x*x);
45: yy[1] = (1. - x*x);
46: yy[2] = -(1.-x)*PetscLogReal(1.-x);
47: PetscDrawLGAddPoint(lg,xx,yy);
48: x += hx;
49: }
50: PetscDrawGetPause(draw,&pause);
51: PetscDrawSetPause(draw,0.0);
52: PetscDrawAxisSetLabels(axis,"Energy","","");
53: PetscDrawLGSetLegend(lg,legend);
54: PetscDrawLGDraw(lg);
56: /*
57: Plot the forces
58: */
59: PetscDrawViewPortsSet(ports,0);
60: PetscDrawLGReset(lg);
61: x = .9;
62: for (i=0; i<Mx; i++) {
63: xx[0] = xx[1] = xx[2] = x;
64: yy[0] = x*x*x - x;
65: yy[1] = -x;
66: yy[2] = 1.0 + PetscLogReal(1. - x);
67: PetscDrawLGAddPoint(lg,xx,yy);
68: x += hx;
69: }
70: PetscDrawAxisSetLabels(axis,"Derivative","","");
71: PetscDrawLGSetLegend(lg,NULL);
72: PetscDrawLGDraw(lg);
74: PetscDrawSetPause(draw,pause);
75: PetscDrawPause(draw);
76: return(0);
77: }