Actual source code: tshistory.c
1: #include <petsc/private/tshistoryimpl.h>
3: /* These macros can be moved to petscimpl.h eventually */
4: #if defined(PETSC_USE_DEBUG)
7: do { \
8: PetscInt b1[2], b2[2]; \
9: b1[0] = -b; \
10: b1[1] = b; \
11: MPIU_Allreduce(b1, b2, 2, MPIU_INT, MPI_MAX, a); \
13: } while (0)
16: do { \
17: PetscMPIInt b1[2], b2[2]; \
18: b1[0] = -(PetscMPIInt)b; \
19: b1[1] = (PetscMPIInt)b; \
20: MPIU_Allreduce(b1, b2, 2, MPI_INT, MPI_MAX, a); \
22: } while (0)
25: do { \
26: PetscReal b1[3], b2[3]; \
27: if (PetscIsNanReal(b)) { \
28: b1[2] = 1; \
29: } else { \
30: b1[2] = 0; \
31: }; \
32: b1[0] = -b; \
33: b1[1] = b; \
34: MPI_Allreduce(b1, b2, 3, MPIU_REAL, MPIU_MAX, a); \
36: } while (0)
38: #else
41: do { \
42: } while (0)
44: do { \
45: } while (0)
47: do { \
48: } while (0)
50: #endif
52: struct _n_TSHistory {
53: MPI_Comm comm; /* used for runtime collective checks */
54: PetscReal *hist; /* time history */
55: PetscInt *hist_id; /* stores the stepid in time history */
56: size_t n; /* current number of steps registered */
57: PetscBool sorted; /* if the history is sorted in ascending order */
58: size_t c; /* current capacity of history */
59: size_t s; /* reallocation size */
60: };
62: PetscErrorCode TSHistoryGetNumSteps(TSHistory tsh, PetscInt *n)
63: {
65: *n = tsh->n;
66: return 0;
67: }
69: PetscErrorCode TSHistoryUpdate(TSHistory tsh, PetscInt id, PetscReal time)
70: {
71: if (tsh->n == tsh->c) { /* reallocation */
72: tsh->c += tsh->s;
73: PetscRealloc(tsh->c * sizeof(*tsh->hist), &tsh->hist);
74: PetscRealloc(tsh->c * sizeof(*tsh->hist_id), &tsh->hist_id);
75: }
76: tsh->sorted = (PetscBool)(tsh->sorted && (tsh->n ? time >= tsh->hist[tsh->n - 1] : PETSC_TRUE));
77: #if defined(PETSC_USE_DEBUG)
78: if (tsh->n) { /* id should be unique */
79: PetscInt loc, *ids;
81: PetscMalloc1(tsh->n, &ids);
82: PetscArraycpy(ids, tsh->hist_id, tsh->n);
83: PetscSortInt(tsh->n, ids);
84: PetscFindInt(id, tsh->n, ids, &loc);
85: PetscFree(ids);
87: }
88: #endif
89: tsh->hist[tsh->n] = time;
90: tsh->hist_id[tsh->n] = id;
91: tsh->n += 1;
92: return 0;
93: }
95: PetscErrorCode TSHistoryGetTime(TSHistory tsh, PetscBool backward, PetscInt step, PetscReal *t)
96: {
97: if (!t) return 0;
99: if (!tsh->sorted) {
100: PetscSortRealWithArrayInt(tsh->n, tsh->hist, tsh->hist_id);
101: tsh->sorted = PETSC_TRUE;
102: }
104: if (!backward) *t = tsh->hist[step];
105: else *t = tsh->hist[tsh->n - step - 1];
106: return 0;
107: }
109: PetscErrorCode TSHistoryGetTimeStep(TSHistory tsh, PetscBool backward, PetscInt step, PetscReal *dt)
110: {
111: if (!dt) return 0;
113: if (!tsh->sorted) {
114: PetscSortRealWithArrayInt(tsh->n, tsh->hist, tsh->hist_id);
115: tsh->sorted = PETSC_TRUE;
116: }
118: if (!backward) *dt = tsh->hist[PetscMin(step + 1, (PetscInt)tsh->n - 1)] - tsh->hist[PetscMin(step, (PetscInt)tsh->n - 1)];
119: else *dt = tsh->hist[PetscMax((PetscInt)tsh->n - step - 1, 0)] - tsh->hist[PetscMax((PetscInt)tsh->n - step - 2, 0)];
120: return 0;
121: }
123: PetscErrorCode TSHistoryGetLocFromTime(TSHistory tsh, PetscReal time, PetscInt *loc)
124: {
126: if (!tsh->sorted) {
127: PetscSortRealWithArrayInt(tsh->n, tsh->hist, tsh->hist_id);
128: tsh->sorted = PETSC_TRUE;
129: }
130: PetscFindReal(time, tsh->n, tsh->hist, PETSC_SMALL, loc);
131: return 0;
132: }
134: PetscErrorCode TSHistorySetHistory(TSHistory tsh, PetscInt n, PetscReal hist[], PetscInt hist_id[], PetscBool sorted)
135: {
139: PetscFree(tsh->hist);
140: PetscFree(tsh->hist_id);
141: tsh->n = (size_t)n;
142: tsh->c = (size_t)n;
143: PetscMalloc1(tsh->n, &tsh->hist);
144: PetscMalloc1(tsh->n, &tsh->hist_id);
145: for (PetscInt i = 0; i < (PetscInt)tsh->n; i++) {
146: tsh->hist[i] = hist[i];
147: tsh->hist_id[i] = hist_id ? hist_id[i] : i;
148: }
149: if (!sorted) PetscSortRealWithArrayInt(tsh->n, tsh->hist, tsh->hist_id);
150: tsh->sorted = PETSC_TRUE;
151: return 0;
152: }
154: PetscErrorCode TSHistoryGetHistory(TSHistory tsh, PetscInt *n, const PetscReal *hist[], const PetscInt *hist_id[], PetscBool *sorted)
155: {
156: if (n) *n = tsh->n;
157: if (hist) *hist = tsh->hist;
158: if (hist_id) *hist_id = tsh->hist_id;
159: if (sorted) *sorted = tsh->sorted;
160: return 0;
161: }
163: PetscErrorCode TSHistoryDestroy(TSHistory *tsh)
164: {
165: if (!*tsh) return 0;
166: PetscFree((*tsh)->hist);
167: PetscFree((*tsh)->hist_id);
168: PetscCommDestroy(&((*tsh)->comm));
169: PetscFree((*tsh));
170: *tsh = NULL;
171: return 0;
172: }
174: PetscErrorCode TSHistoryCreate(MPI_Comm comm, TSHistory *hst)
175: {
176: TSHistory tsh;
179: *hst = NULL;
180: PetscNew(&tsh);
181: PetscCommDuplicate(comm, &tsh->comm, NULL);
183: tsh->c = 1024; /* capacity */
184: tsh->s = 1024; /* reallocation size */
185: tsh->sorted = PETSC_TRUE;
187: PetscMalloc1(tsh->c, &tsh->hist);
188: PetscMalloc1(tsh->c, &tsh->hist_id);
189: *hst = tsh;
190: return 0;
191: }