Actual source code: dviewp.c

  1: #include <petsc/private/drawimpl.h>

  3: /*@
  4:   PetscDrawSetViewPort - Sets the portion of the window (page) to which draw
  5:   routines will write.

  7:   Collective

  9:   Input Parameters:
 10: + xl   - the horizontal coordinate of the lower left corner of the subwindow.
 11: . yl   - the vertical coordinate of the lower left corner of the subwindow.
 12: . xr   - the horizontal coordinate of the upper right corner of the subwindow.
 13: . yr   - the vertical coordinate of the upper right corner of the subwindow.
 14: - draw - the drawing context

 16:   Level: advanced

 18:   Notes:
 19:   These numbers must always be between 0.0 and 1.0.

 21:   Lower left corner is (0,0).

 23: .seealso: `PetscDrawGetViewPort()`, `PetscDraw`, `PetscDrawSplitViewPort()`, `PetscDrawViewPortsCreate()`
 24: @*/
 25: PetscErrorCode PetscDrawSetViewPort(PetscDraw draw, PetscReal xl, PetscReal yl, PetscReal xr, PetscReal yr)
 26: {
 27:   PetscFunctionBegin;
 29:   PetscCheck(xl >= 0.0 && xr <= 1.0 && yl >= 0.0 && yr <= 1.0 && xr > xl && yr > yl, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "ViewPort values must be >= 0 and <= 1: Instead %g %g %g %g", (double)xl, (double)yl, (double)xr, (double)yr);
 30:   draw->port_xl = xl;
 31:   draw->port_yl = yl;
 32:   draw->port_xr = xr;
 33:   draw->port_yr = yr;
 34:   PetscTryTypeMethod(draw, setviewport, xl, yl, xr, yr);
 35:   PetscFunctionReturn(PETSC_SUCCESS);
 36: }

 38: /*@
 39:   PetscDrawGetViewPort - Gets the portion of the window (page) to which draw
 40:   routines will write.

 42:   Collective

 44:   Input Parameter:
 45: . draw - the drawing context

 47:   Output Parameters:
 48: + xl - the horizontal coordinate of the lower left corner of the subwindow.
 49: . yl - the vertical coordinate of the lower left corner of the subwindow.
 50: . xr - the horizontal coordinate of the upper right corner of the subwindow.
 51: - yr - the vertical coordinate of the upper right corner of the subwindow.

 53:   Level: advanced

 55:   Notes:
 56:   These numbers must always be between 0.0 and 1.0.

 58:   Lower left corner is (0,0).

 60: .seealso: `PetscDraw`, `PetscDrawSplitViewPort()`, `PetscDrawSetViewPort()`
 61: @*/
 62: PetscErrorCode PetscDrawGetViewPort(PetscDraw draw, PetscReal *xl, PetscReal *yl, PetscReal *xr, PetscReal *yr)
 63: {
 64:   PetscFunctionBegin;
 66:   PetscAssertPointer(xl, 2);
 67:   PetscAssertPointer(yl, 3);
 68:   PetscAssertPointer(xr, 4);
 69:   PetscAssertPointer(yr, 5);
 70:   *xl = draw->port_xl;
 71:   *yl = draw->port_yl;
 72:   *xr = draw->port_xr;
 73:   *yr = draw->port_yr;
 74:   PetscFunctionReturn(PETSC_SUCCESS);
 75: }

 77: /*@
 78:   PetscDrawSplitViewPort - Splits a window shared by several processes into smaller
 79:   view ports. One for each process.

 81:   Collective

 83:   Input Parameter:
 84: . draw - the drawing context

 86:   Level: advanced

 88: .seealso: `PetscDrawDivideViewPort()`, `PetscDrawSetViewPort()`
 89: @*/
 90: PetscErrorCode PetscDrawSplitViewPort(PetscDraw draw)
 91: {
 92:   PetscMPIInt rank, size;
 93:   PetscInt    n;
 94:   PetscBool   isnull;
 95:   PetscReal   xl, xr, yl, yr, h;

 97:   PetscFunctionBegin;
 99:   PetscCall(PetscDrawIsNull(draw, &isnull));
100:   if (isnull) PetscFunctionReturn(PETSC_SUCCESS);
101:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)draw), &rank));
102:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)draw), &size));

104:   n = (PetscInt)(.1 + PetscSqrtReal((PetscReal)size));
105:   while (n * n < size) n++;

107:   h  = 1.0 / n;
108:   xl = (rank % n) * h;
109:   xr = xl + h;
110:   yl = (rank / n) * h;
111:   yr = yl + h;

113:   PetscDrawCollectiveBegin(draw);
114:   PetscCall(PetscDrawLine(draw, xl, yl, xl, yr, PETSC_DRAW_BLACK));
115:   PetscCall(PetscDrawLine(draw, xl, yr, xr, yr, PETSC_DRAW_BLACK));
116:   PetscCall(PetscDrawLine(draw, xr, yr, xr, yl, PETSC_DRAW_BLACK));
117:   PetscCall(PetscDrawLine(draw, xr, yl, xl, yl, PETSC_DRAW_BLACK));
118:   PetscDrawCollectiveEnd(draw);
119:   PetscCall(PetscDrawFlush(draw));

121:   draw->port_xl = xl + .05 * h;
122:   draw->port_xr = xr - .05 * h;
123:   draw->port_yl = yl + .05 * h;
124:   draw->port_yr = yr - .05 * h;

126:   PetscTryTypeMethod(draw, setviewport, xl, yl, xr, yr);
127:   PetscFunctionReturn(PETSC_SUCCESS);
128: }

130: /*@C
131:   PetscDrawViewPortsCreate - Splits a window into smaller view ports. Each processor shares all the viewports.

133:   Collective

135:   Input Parameters:
136: + draw   - the drawing context
137: - nports - the number of ports

139:   Output Parameter:
140: . newports - a `PetscDrawViewPorts` context (C structure)

142:   Options Database Key:
143: . -draw_ports - display multiple fields in the same window with PetscDrawPorts() instead of in separate windows

145:   Level: advanced

147: .seealso: `PetscDrawSplitViewPort()`, `PetscDrawSetViewPort()`, `PetscDrawViewPortsSet()`, `PetscDrawViewPortsDestroy()`
148: @*/
149: PetscErrorCode PetscDrawViewPortsCreate(PetscDraw draw, PetscInt nports, PetscDrawViewPorts **newports)
150: {
151:   PetscDrawViewPorts *ports;
152:   PetscInt            i, n;
153:   PetscBool           isnull;
154:   PetscMPIInt         rank;
155:   PetscReal          *xl, *xr, *yl, *yr, h;

157:   PetscFunctionBegin;
159:   PetscCheck(nports >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Number of divisions must be positive: %" PetscInt_FMT, nports);
160:   PetscAssertPointer(newports, 3);
161:   PetscCall(PetscDrawIsNull(draw, &isnull));
162:   if (isnull) {
163:     *newports = NULL;
164:     PetscFunctionReturn(PETSC_SUCCESS);
165:   }
166:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)draw), &rank));

168:   PetscCall(PetscNew(&ports));
169:   *newports     = ports;
170:   ports->draw   = draw;
171:   ports->nports = nports;
172:   PetscCall(PetscObjectReference((PetscObject)draw));
173:   /* save previous drawport of window */
174:   PetscCall(PetscDrawGetViewPort(draw, &ports->port_xl, &ports->port_yl, &ports->port_xr, &ports->port_yr));

176:   n = (PetscInt)(.1 + PetscSqrtReal((PetscReal)nports));
177:   while (n * n < nports) n++;
178:   h = 1.0 / n;

180:   PetscCall(PetscMalloc4(n * n, &xl, n * n, &xr, n * n, &yl, n * n, &yr));
181:   ports->xl = xl;
182:   ports->xr = xr;
183:   ports->yl = yl;
184:   ports->yr = yr;

186:   PetscCall(PetscDrawSetCoordinates(draw, 0.0, 0.0, 1.0, 1.0));
187:   PetscDrawCollectiveBegin(draw);
188:   for (i = 0; i < n * n; i++) {
189:     xl[i] = (i % n) * h;
190:     xr[i] = xl[i] + h;
191:     yl[i] = (i / n) * h;
192:     yr[i] = yl[i] + h;

194:     if (rank == 0) {
195:       PetscCall(PetscDrawLine(draw, xl[i], yl[i], xl[i], yr[i], PETSC_DRAW_BLACK));
196:       PetscCall(PetscDrawLine(draw, xl[i], yr[i], xr[i], yr[i], PETSC_DRAW_BLACK));
197:       PetscCall(PetscDrawLine(draw, xr[i], yr[i], xr[i], yl[i], PETSC_DRAW_BLACK));
198:       PetscCall(PetscDrawLine(draw, xr[i], yl[i], xl[i], yl[i], PETSC_DRAW_BLACK));
199:     }

201:     xl[i] += .05 * h;
202:     xr[i] -= .05 * h;
203:     yl[i] += .05 * h;
204:     yr[i] -= .05 * h;
205:   }
206:   PetscDrawCollectiveEnd(draw);
207:   PetscCall(PetscDrawFlush(draw));
208:   PetscFunctionReturn(PETSC_SUCCESS);
209: }

211: /*@C
212:   PetscDrawViewPortsCreateRect - Splits a window into smaller
213:   view ports. Each processor shares all the viewports. The number
214:   of views in the x- and y-directions is specified.

216:   Collective

218:   Input Parameters:
219: + draw - the drawing context
220: . nx   - the number of x divisions
221: - ny   - the number of y divisions

223:   Output Parameter:
224: . newports - a `PetscDrawViewPorts` context (C structure)

226:   Level: advanced

228: .seealso: `PetscDrawSplitViewPort()`, `PetscDrawSetViewPort()`, `PetscDrawViewPortsSet()`, `PetscDrawViewPortsDestroy()`, `PetscDrawViewPorts`
229: @*/
230: PetscErrorCode PetscDrawViewPortsCreateRect(PetscDraw draw, PetscInt nx, PetscInt ny, PetscDrawViewPorts **newports)
231: {
232:   PetscDrawViewPorts *ports;
233:   PetscReal          *xl, *xr, *yl, *yr, hx, hy;
234:   PetscInt            i, j, k, n;
235:   PetscBool           isnull;
236:   PetscMPIInt         rank;

238:   PetscFunctionBegin;
240:   PetscCheck(nx >= 1 && ny >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Number of divisions must be positive: %" PetscInt_FMT " x %" PetscInt_FMT, nx, ny);
241:   PetscAssertPointer(newports, 4);
242:   PetscCall(PetscDrawIsNull(draw, &isnull));
243:   if (isnull) {
244:     *newports = NULL;
245:     PetscFunctionReturn(PETSC_SUCCESS);
246:   }
247:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)draw), &rank));

249:   n  = nx * ny;
250:   hx = 1.0 / nx;
251:   hy = 1.0 / ny;
252:   PetscCall(PetscNew(&ports));
253:   *newports     = ports;
254:   ports->draw   = draw;
255:   ports->nports = n;
256:   PetscCall(PetscObjectReference((PetscObject)draw));
257:   /* save previous drawport of window */
258:   PetscCall(PetscDrawGetViewPort(draw, &ports->port_xl, &ports->port_yl, &ports->port_xr, &ports->port_yr));

260:   PetscCall(PetscMalloc4(n, &xl, n, &xr, n, &yl, n, &yr));
261:   ports->xr = xr;
262:   ports->xl = xl;
263:   ports->yl = yl;
264:   ports->yr = yr;

266:   PetscCall(PetscDrawSetCoordinates(draw, 0.0, 0.0, 1.0, 1.0));
267:   PetscDrawCollectiveBegin(draw);
268:   for (i = 0; i < nx; i++) {
269:     for (j = 0; j < ny; j++) {
270:       k = j * nx + i;

272:       xl[k] = i * hx;
273:       xr[k] = xl[k] + hx;
274:       yl[k] = j * hy;
275:       yr[k] = yl[k] + hy;

277:       if (rank == 0) {
278:         PetscCall(PetscDrawLine(draw, xl[k], yl[k], xl[k], yr[k], PETSC_DRAW_BLACK));
279:         PetscCall(PetscDrawLine(draw, xl[k], yr[k], xr[k], yr[k], PETSC_DRAW_BLACK));
280:         PetscCall(PetscDrawLine(draw, xr[k], yr[k], xr[k], yl[k], PETSC_DRAW_BLACK));
281:         PetscCall(PetscDrawLine(draw, xr[k], yl[k], xl[k], yl[k], PETSC_DRAW_BLACK));
282:       }

284:       xl[k] += .05 * hx;
285:       xr[k] -= .05 * hx;
286:       yl[k] += .05 * hy;
287:       yr[k] -= .05 * hy;
288:     }
289:   }
290:   PetscDrawCollectiveEnd(draw);
291:   PetscCall(PetscDrawFlush(draw));
292:   PetscFunctionReturn(PETSC_SUCCESS);
293: }

295: /*@C
296:   PetscDrawViewPortsDestroy - frees a `PetscDrawViewPorts` object

298:   Collective on the `PetscDraw` inside `ports`

300:   Input Parameter:
301: . ports - the `PetscDrawViewPorts` object

303:   Level: advanced

305: .seealso: `PetscDrawViewPorts`, `PetscDrawSplitViewPort()`, `PetscDrawSetViewPort()`, `PetscDrawViewPortsSet()`, `PetscDrawViewPortsCreate()`
306: @*/
307: PetscErrorCode PetscDrawViewPortsDestroy(PetscDrawViewPorts *ports)
308: {
309:   PetscFunctionBegin;
310:   if (!ports) PetscFunctionReturn(PETSC_SUCCESS);
311:   PetscAssertPointer(ports, 1);
312:   /* reset Drawport of Window back to previous value */
313:   PetscCall(PetscDrawSetViewPort(ports->draw, ports->port_xl, ports->port_yl, ports->port_xr, ports->port_yr));
314:   PetscCall(PetscDrawDestroy(&ports->draw));
315:   PetscCall(PetscFree4(ports->xl, ports->xr, ports->yl, ports->yr));
316:   PetscCall(PetscFree(ports));
317:   PetscFunctionReturn(PETSC_SUCCESS);
318: }

320: /*@C
321:   PetscDrawViewPortsSet - sets a draw object to use a particular subport

323:   Logically Collective on the `PetscDraw` inside `ports`

325:   Input Parameters:
326: + ports - the `PetscDrawViewPorts` object
327: - port  - the port number, from 0 to nports-1

329:   Level: advanced

331: .seealso: `PetscDrawViewPorts`, `PetscDrawSplitViewPort()`, `PetscDrawSetViewPort()`, `PetscDrawViewPortsDestroy()`, `PetscDrawViewPortsCreate()`
332: @*/
333: PetscErrorCode PetscDrawViewPortsSet(PetscDrawViewPorts *ports, PetscInt port)
334: {
335:   PetscFunctionBegin;
336:   if (!ports) PetscFunctionReturn(PETSC_SUCCESS);
337:   PetscAssertPointer(ports, 1);
338:   PetscCheck(port >= 0 && (port <= ports->nports - 1), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Port is out of range requested %" PetscInt_FMT " from 0 to %" PetscInt_FMT, port, ports->nports - 1);
339:   PetscCall(PetscDrawSetViewPort(ports->draw, ports->xl[port], ports->yl[port], ports->xr[port], ports->yr[port]));
340:   PetscFunctionReturn(PETSC_SUCCESS);
341: }