Actual source code: pointqueue.c
1: #include <petsc/private/dmpleximpl.h>
3: PetscErrorCode DMPlexPointQueueCreate(PetscInt size, DMPlexPointQueue *queue)
4: {
5: DMPlexPointQueue q;
8: PetscCalloc1(1, &q);
9: q->size = size;
10: PetscMalloc1(q->size, &q->points);
11: q->num = 0;
12: q->front = 0;
13: q->back = q->size - 1;
14: *queue = q;
15: return 0;
16: }
18: PetscErrorCode DMPlexPointQueueDestroy(DMPlexPointQueue *queue)
19: {
20: DMPlexPointQueue q = *queue;
22: PetscFree(q->points);
23: PetscFree(q);
24: *queue = NULL;
25: return 0;
26: }
28: PetscErrorCode DMPlexPointQueueEnsureSize(DMPlexPointQueue queue)
29: {
30: if (queue->num < queue->size) return 0;
31: queue->size *= 2;
32: PetscRealloc(queue->size * sizeof(PetscInt), &queue->points);
33: return 0;
34: }
36: PetscErrorCode DMPlexPointQueueEnqueue(DMPlexPointQueue queue, PetscInt p)
37: {
38: DMPlexPointQueueEnsureSize(queue);
39: queue->back = (queue->back + 1) % queue->size;
40: queue->points[queue->back] = p;
41: ++queue->num;
42: return 0;
43: }
45: PetscErrorCode DMPlexPointQueueDequeue(DMPlexPointQueue queue, PetscInt *p)
46: {
48: *p = queue->points[queue->front];
49: queue->front = (queue->front + 1) % queue->size;
50: --queue->num;
51: return 0;
52: }
54: PetscErrorCode DMPlexPointQueueFront(DMPlexPointQueue queue, PetscInt *p)
55: {
57: *p = queue->points[queue->front];
58: return 0;
59: }
61: PetscErrorCode DMPlexPointQueueBack(DMPlexPointQueue queue, PetscInt *p)
62: {
64: *p = queue->points[queue->back];
65: return 0;
66: }
68: PetscBool DMPlexPointQueueEmpty(DMPlexPointQueue queue)
69: {
70: if (!queue->num) return PETSC_TRUE;
71: return PETSC_FALSE;
72: }
74: PetscErrorCode DMPlexPointQueueEmptyCollective(PetscObject obj, DMPlexPointQueue queue, PetscBool *empty)
75: {
77: *empty = DMPlexPointQueueEmpty(queue);
78: MPI_Allreduce(MPI_IN_PLACE, empty, 1, MPIU_BOOL, MPI_LAND, PetscObjectComm(obj));
79: return 0;
80: }