1c4762a1bSJed Brown static char help[] = "Plots the various potentials used in the examples.\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscdmda.h> 4c4762a1bSJed Brown #include <petscts.h> 5c4762a1bSJed Brown #include <petscdraw.h> 6c4762a1bSJed Brown 7d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv) 8d71ae5a4SJacob Faibussowitsch { 9c4762a1bSJed Brown PetscDrawLG lg; 10c4762a1bSJed Brown PetscInt Mx = 100, i; 11c4762a1bSJed Brown PetscReal x, hx = .1 / Mx, pause, xx[3], yy[3]; 12c4762a1bSJed Brown PetscDraw draw; 13c4762a1bSJed Brown const char *const legend[] = {"(1 - u^2)^2", "1 - u^2", "-(1 - u)log(1 - u)"}; 14c4762a1bSJed Brown PetscDrawAxis axis; 15c4762a1bSJed Brown PetscDrawViewPorts *ports; 16c4762a1bSJed Brown 17c4762a1bSJed Brown PetscFunctionBegin; 18327415f7SBarry Smith PetscFunctionBeginUser; 199566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &argv, 0, help)); 209566063dSJacob Faibussowitsch PetscCall(PetscViewerDrawResize(PETSC_VIEWER_DRAW_(PETSC_COMM_WORLD), 1200, 800)); 219566063dSJacob Faibussowitsch PetscCall(PetscViewerDrawGetDrawLG(PETSC_VIEWER_DRAW_(PETSC_COMM_WORLD), 0, &lg)); 229566063dSJacob Faibussowitsch PetscCall(PetscDrawLGGetDraw(lg, &draw)); 239566063dSJacob Faibussowitsch PetscCall(PetscDrawCheckResizedWindow(draw)); 249566063dSJacob Faibussowitsch PetscCall(PetscDrawViewPortsCreateRect(draw, 1, 2, &ports)); 259566063dSJacob Faibussowitsch PetscCall(PetscDrawLGGetAxis(lg, &axis)); 269566063dSJacob Faibussowitsch PetscCall(PetscDrawLGReset(lg)); 27c4762a1bSJed Brown 28c4762a1bSJed Brown /* 29c4762a1bSJed Brown Plot the energies 30c4762a1bSJed Brown */ 319566063dSJacob Faibussowitsch PetscCall(PetscDrawLGSetDimension(lg, 3)); 329566063dSJacob Faibussowitsch PetscCall(PetscDrawViewPortsSet(ports, 1)); 33c4762a1bSJed Brown x = .9; 34c4762a1bSJed Brown for (i = 0; i < Mx; i++) { 35c4762a1bSJed Brown xx[0] = xx[1] = xx[2] = x; 36c4762a1bSJed Brown yy[0] = (1. - x * x) * (1. - x * x); 37c4762a1bSJed Brown yy[1] = (1. - x * x); 38c4762a1bSJed Brown yy[2] = -(1. - x) * PetscLogReal(1. - x); 399566063dSJacob Faibussowitsch PetscCall(PetscDrawLGAddPoint(lg, xx, yy)); 40c4762a1bSJed Brown x += hx; 41c4762a1bSJed Brown } 429566063dSJacob Faibussowitsch PetscCall(PetscDrawGetPause(draw, &pause)); 439566063dSJacob Faibussowitsch PetscCall(PetscDrawSetPause(draw, 0.0)); 449566063dSJacob Faibussowitsch PetscCall(PetscDrawAxisSetLabels(axis, "Energy", "", "")); 459566063dSJacob Faibussowitsch PetscCall(PetscDrawLGSetLegend(lg, legend)); 469566063dSJacob Faibussowitsch PetscCall(PetscDrawLGDraw(lg)); 47c4762a1bSJed Brown 48c4762a1bSJed Brown /* 49c4762a1bSJed Brown Plot the forces 50c4762a1bSJed Brown */ 519566063dSJacob Faibussowitsch PetscCall(PetscDrawViewPortsSet(ports, 0)); 529566063dSJacob Faibussowitsch PetscCall(PetscDrawLGReset(lg)); 53c4762a1bSJed Brown x = .9; 54c4762a1bSJed Brown for (i = 0; i < Mx; i++) { 55c4762a1bSJed Brown xx[0] = xx[1] = xx[2] = x; 56c4762a1bSJed Brown yy[0] = x * x * x - x; 57c4762a1bSJed Brown yy[1] = -x; 58c4762a1bSJed Brown yy[2] = 1.0 + PetscLogReal(1. - x); 599566063dSJacob Faibussowitsch PetscCall(PetscDrawLGAddPoint(lg, xx, yy)); 60c4762a1bSJed Brown x += hx; 61c4762a1bSJed Brown } 629566063dSJacob Faibussowitsch PetscCall(PetscDrawAxisSetLabels(axis, "Derivative", "", "")); 639566063dSJacob Faibussowitsch PetscCall(PetscDrawLGSetLegend(lg, NULL)); 649566063dSJacob Faibussowitsch PetscCall(PetscDrawLGDraw(lg)); 65c4762a1bSJed Brown 669566063dSJacob Faibussowitsch PetscCall(PetscDrawSetPause(draw, pause)); 679566063dSJacob Faibussowitsch PetscCall(PetscDrawPause(draw)); 689566063dSJacob Faibussowitsch PetscCall(PetscDrawViewPortsDestroy(ports)); 699566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 70b122ec5aSJacob Faibussowitsch return 0; 71c4762a1bSJed Brown } 72c4762a1bSJed Brown 73c4762a1bSJed Brown /*TEST 74c4762a1bSJed Brown 75c4762a1bSJed Brown test: 76c4762a1bSJed Brown requires: x 77*3886731fSPierre Jolivet output_file: output/empty.out 78c4762a1bSJed Brown 79c4762a1bSJed Brown TEST*/ 80