Actual source code: dtri.c
2: /*
3: Provides the calling sequences for all the basic PetscDraw routines.
4: */
5: #include <petsc/private/drawimpl.h>
7: /*@
8: PetscDrawTriangle - draws a triangle onto a drawable.
10: Not Collective
12: Input Parameters:
13: + draw - the drawing context
14: . x1,y1,x2,y2,x3,y3 - the coordinates of the vertices
15: - c1,c2,c3 - the colors of the three vertices in the same order as the xi,yi
17: Level: beginner
19: .seealso: `PetscDraw`, `PetscDrawLine()`, `PetscDrawRectangle()`, `PetscDrawEllipse()`, `PetscDrawMarker()`, `PetscDrawPoint()`, `PetscDrawArrow()`
20: @*/
21: PetscErrorCode PetscDrawTriangle(PetscDraw draw, PetscReal x1, PetscReal y_1, PetscReal x2, PetscReal y2, PetscReal x3, PetscReal y3, int c1, int c2, int c3)
22: {
23: PetscFunctionBegin;
25: PetscUseTypeMethod(draw, triangle, x1, y_1, x2, y2, x3, y3, c1, c2, c3);
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: /*@
30: PetscDrawScalePopup - draws a contour scale window.
32: Collective
34: Input Parameters:
35: + popup - the window (often a window obtained via `PetscDrawGetPopup()`
36: . min - minimum value being plotted
37: - max - maximum value being plotted
39: Level: intermediate
41: Note:
42: All processors that share the draw MUST call this routine
44: .seealso: `PetscDraw`, `PetscDrawGetPopup()`, `PetscDrawTensorContour()`
45: @*/
46: PetscErrorCode PetscDrawScalePopup(PetscDraw popup, PetscReal min, PetscReal max)
47: {
48: PetscBool isnull;
49: PetscReal xl = 0.0, yl = 0.0, xr = 1.0, yr = 1.0;
50: PetscMPIInt rank;
51: int i;
52: char string[32];
54: PetscFunctionBegin;
55: if (!popup) PetscFunctionReturn(PETSC_SUCCESS);
57: PetscCall(PetscDrawIsNull(popup, &isnull));
58: if (isnull) PetscFunctionReturn(PETSC_SUCCESS);
59: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)popup), &rank));
61: PetscCall(PetscDrawCheckResizedWindow(popup));
62: PetscCall(PetscDrawClear(popup));
63: PetscCall(PetscDrawSetTitle(popup, "Contour Scale"));
64: PetscCall(PetscDrawSetCoordinates(popup, xl, yl, xr, yr));
65: PetscDrawCollectiveBegin(popup);
66: if (rank == 0) {
67: for (i = 0; i < 10; i++) {
68: int c = PetscDrawRealToColor((PetscReal)i / 9, 0, 1);
69: PetscCall(PetscDrawRectangle(popup, xl, yl, xr, yr, c, c, c, c));
70: yl += 0.1;
71: }
72: for (i = 0; i < 10; i++) {
73: PetscReal value = min + i * (max - min) / 9;
74: /* look for a value that should be zero, but is not due to round-off */
75: if (PetscAbsReal(value) < 1.e-10 && max - min > 1.e-6) value = 0.0;
76: PetscCall(PetscSNPrintf(string, sizeof(string), "%18.16e", (double)value));
77: PetscCall(PetscDrawString(popup, 0.2, 0.02 + i / 10.0, PETSC_DRAW_BLACK, string));
78: }
79: }
80: PetscDrawCollectiveEnd(popup);
81: PetscCall(PetscDrawFlush(popup));
82: PetscCall(PetscDrawSave(popup));
83: PetscFunctionReturn(PETSC_SUCCESS);
84: }
86: typedef struct {
87: int m, n;
88: PetscReal *x, *y, min, max, *v;
89: PetscBool showgrid;
90: } ZoomCtx;
92: static PetscErrorCode PetscDrawTensorContour_Zoom(PetscDraw win, void *dctx)
93: {
94: int i;
95: ZoomCtx *ctx = (ZoomCtx *)dctx;
97: PetscFunctionBegin;
98: PetscCall(PetscDrawTensorContourPatch(win, ctx->m, ctx->n, ctx->x, ctx->y, ctx->min, ctx->max, ctx->v));
99: if (ctx->showgrid) {
100: for (i = 0; i < ctx->m; i++) PetscCall(PetscDrawLine(win, ctx->x[i], ctx->y[0], ctx->x[i], ctx->y[ctx->n - 1], PETSC_DRAW_BLACK));
101: for (i = 0; i < ctx->n; i++) PetscCall(PetscDrawLine(win, ctx->x[0], ctx->y[i], ctx->x[ctx->m - 1], ctx->y[i], PETSC_DRAW_BLACK));
102: }
103: PetscFunctionReturn(PETSC_SUCCESS);
104: }
106: /*@C
107: PetscDrawTensorContour - draws a contour plot for a two-dimensional array
109: Collective on draw, but draw must be sequential
111: Input Parameters:
112: + draw - the draw context
113: . m,n - the global number of mesh points in the x and y directions
114: . xi,yi - the locations of the global mesh points (optional, use NULL
115: to indicate uniform spacing on [0,1])
116: - V - the values
118: Options Database Keys:
119: + -draw_x_shared_colormap - Indicates use of private colormap
120: - -draw_contour_grid - draws grid contour
122: Level: intermediate
124: .seealso: `PetscDraw`, `PetscDrawTensorContourPatch()`, `PetscDrawScalePopup()`
125: @*/
126: PetscErrorCode PetscDrawTensorContour(PetscDraw draw, int m, int n, const PetscReal xi[], const PetscReal yi[], PetscReal *v)
127: {
128: int N = m * n;
129: PetscBool isnull;
130: PetscDraw popup;
131: int xin = 1, yin = 1, i;
132: PetscMPIInt size;
133: PetscReal h;
134: ZoomCtx ctx;
136: PetscFunctionBegin;
138: PetscCall(PetscDrawIsNull(draw, &isnull));
139: if (isnull) PetscFunctionReturn(PETSC_SUCCESS);
140: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)draw), &size));
141: PetscCheck(size <= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "May only be used with single processor PetscDraw");
142: PetscCheck(N > 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "n %d and m %d must be positive", m, n);
144: ctx.v = v;
145: ctx.m = m;
146: ctx.n = n;
147: ctx.max = ctx.min = v[0];
148: for (i = 0; i < N; i++) {
149: if (ctx.max < ctx.v[i]) ctx.max = ctx.v[i];
150: if (ctx.min > ctx.v[i]) ctx.min = ctx.v[i];
151: }
152: if (ctx.max - ctx.min < 1.e-7) {
153: ctx.min -= 5.e-8;
154: ctx.max += 5.e-8;
155: }
157: /* PetscDraw the scale window */
158: PetscCall(PetscDrawGetPopup(draw, &popup));
159: PetscCall(PetscDrawScalePopup(popup, ctx.min, ctx.max));
161: ctx.showgrid = PETSC_FALSE;
162: PetscCall(PetscOptionsGetBool(((PetscObject)draw)->options, NULL, "-draw_contour_grid", &ctx.showgrid, NULL));
164: /* fill up x and y coordinates */
165: if (!xi) {
166: xin = 0;
167: PetscCall(PetscMalloc1(ctx.m, &ctx.x));
168: h = 1.0 / (ctx.m - 1);
169: ctx.x[0] = 0.0;
170: for (i = 1; i < ctx.m; i++) ctx.x[i] = ctx.x[i - 1] + h;
171: } else ctx.x = (PetscReal *)xi;
173: if (!yi) {
174: yin = 0;
175: PetscCall(PetscMalloc1(ctx.n, &ctx.y));
176: h = 1.0 / (ctx.n - 1);
177: ctx.y[0] = 0.0;
178: for (i = 1; i < ctx.n; i++) ctx.y[i] = ctx.y[i - 1] + h;
179: } else ctx.y = (PetscReal *)yi;
181: PetscCall(PetscDrawZoom(draw, PetscDrawTensorContour_Zoom, &ctx));
183: if (!xin) PetscCall(PetscFree(ctx.x));
184: if (!yin) PetscCall(PetscFree(ctx.y));
185: PetscFunctionReturn(PETSC_SUCCESS);
186: }
188: /*@
189: PetscDrawTensorContourPatch - draws a rectangular patch of a contour plot
190: for a two-dimensional array.
192: Not Collective
194: Input Parameters:
195: + draw - the draw context
196: . m,n - the number of local mesh points in the x and y direction
197: . x,y - the locations of the local mesh points
198: . min,max - the minimum and maximum value in the entire contour
199: - v - the data
201: Options Database Keys:
202: . -draw_x_shared_colormap - Activates private colormap
204: Level: advanced
206: Note:
207: This is a lower level support routine, usually the user will call
208: `PetscDrawTensorContour()`.
210: .seealso: `PetscDraw`, `PetscDrawTensorContour()`
211: @*/
212: PetscErrorCode PetscDrawTensorContourPatch(PetscDraw draw, int m, int n, PetscReal *x, PetscReal *y, PetscReal min, PetscReal max, PetscReal *v)
213: {
214: int c1, c2, c3, c4, i, j;
215: PetscReal x1, x2, x3, x4, y_1, y2, y3, y4;
217: PetscFunctionBegin;
219: /* PetscDraw the contour plot patch */
220: for (j = 0; j < n - 1; j++) {
221: for (i = 0; i < m - 1; i++) {
222: x1 = x[i];
223: y_1 = y[j];
224: c1 = PetscDrawRealToColor(v[i + j * m], min, max);
225: x2 = x[i + 1];
226: y2 = y_1;
227: c2 = PetscDrawRealToColor(v[i + j * m + 1], min, max);
228: x3 = x2;
229: y3 = y[j + 1];
230: c3 = PetscDrawRealToColor(v[i + j * m + 1 + m], min, max);
231: x4 = x1;
232: y4 = y3;
233: c4 = PetscDrawRealToColor(v[i + j * m + m], min, max);
235: PetscCall(PetscDrawTriangle(draw, x1, y_1, x2, y2, x3, y3, c1, c2, c3));
236: PetscCall(PetscDrawTriangle(draw, x1, y_1, x3, y3, x4, y4, c1, c3, c4));
237: }
238: }
239: PetscFunctionReturn(PETSC_SUCCESS);
240: }