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; b1[1] = b; \
10: MPIU_Allreduce(b1,b2,2,MPIU_INT,MPI_MAX,a); \
12: } while (0)
15: do { \
16: PetscMPIInt b1[2],b2[2]; \
17: b1[0] = -(PetscMPIInt)b; b1[1] = (PetscMPIInt)b; \
18: MPIU_Allreduce(b1,b2,2,MPI_INT,MPI_MAX,a); \
20: } while (0)
23: do { \
24: PetscReal b1[3],b2[3]; \
25: if (PetscIsNanReal(b)) {b1[2] = 1;} else {b1[2] = 0;}; \
26: b1[0] = -b; b1[1] = b; \
27: MPI_Allreduce(b1,b2,3,MPIU_REAL,MPIU_MAX,a); \
29: } while (0)
31: #else
37: #endif
39: struct _n_TSHistory {
40: MPI_Comm comm; /* used for runtime collective checks */
41: PetscReal *hist; /* time history */
42: PetscInt *hist_id; /* stores the stepid in time history */
43: size_t n; /* current number of steps registered */
44: PetscBool sorted; /* if the history is sorted in ascending order */
45: size_t c; /* current capacity of history */
46: size_t s; /* reallocation size */
47: };
49: PetscErrorCode TSHistoryGetNumSteps(TSHistory tsh, PetscInt *n)
50: {
52: *n = tsh->n;
53: return 0;
54: }
56: PetscErrorCode TSHistoryUpdate(TSHistory tsh, PetscInt id, PetscReal time)
57: {
58: if (tsh->n == tsh->c) { /* reallocation */
59: tsh->c += tsh->s;
60: PetscRealloc(tsh->c*sizeof(*tsh->hist),&tsh->hist);
61: PetscRealloc(tsh->c*sizeof(*tsh->hist_id),&tsh->hist_id);
62: }
63: tsh->sorted = (PetscBool)(tsh->sorted && (tsh->n ? time >= tsh->hist[tsh->n-1] : PETSC_TRUE));
64: #if defined(PETSC_USE_DEBUG)
65: if (tsh->n) { /* id should be unique */
66: PetscInt loc,*ids;
68: PetscMalloc1(tsh->n,&ids);
69: PetscArraycpy(ids,tsh->hist_id,tsh->n);
70: PetscSortInt(tsh->n,ids);
71: PetscFindInt(id,tsh->n,ids,&loc);
72: PetscFree(ids);
74: }
75: #endif
76: tsh->hist[tsh->n] = time;
77: tsh->hist_id[tsh->n] = id;
78: tsh->n += 1;
79: return 0;
80: }
82: PetscErrorCode TSHistoryGetTime(TSHistory tsh, PetscBool backward, PetscInt step, PetscReal *t)
83: {
84: if (!t) return 0;
86: if (!tsh->sorted) {
88: PetscSortRealWithArrayInt(tsh->n,tsh->hist,tsh->hist_id);
89: tsh->sorted = PETSC_TRUE;
90: }
92: if (!backward) *t = tsh->hist[step];
93: else *t = tsh->hist[tsh->n-step-1];
94: return 0;
95: }
97: PetscErrorCode TSHistoryGetTimeStep(TSHistory tsh, PetscBool backward, PetscInt step, PetscReal *dt)
98: {
99: if (!dt) return 0;
101: if (!tsh->sorted) {
103: PetscSortRealWithArrayInt(tsh->n,tsh->hist,tsh->hist_id);
104: tsh->sorted = PETSC_TRUE;
105: }
107: if (!backward) *dt = tsh->hist[PetscMin(step+1,(PetscInt)tsh->n-1)] - tsh->hist[PetscMin(step,(PetscInt)tsh->n-1)];
108: else *dt = tsh->hist[PetscMax((PetscInt)tsh->n-step-1,0)] - tsh->hist[PetscMax((PetscInt)tsh->n-step-2,0)];
109: return 0;
110: }
112: PetscErrorCode TSHistoryGetLocFromTime(TSHistory tsh, PetscReal time, PetscInt *loc)
113: {
115: if (!tsh->sorted) {
116: PetscSortRealWithArrayInt(tsh->n,tsh->hist,tsh->hist_id);
117: tsh->sorted = PETSC_TRUE;
118: }
119: PetscFindReal(time,tsh->n,tsh->hist,PETSC_SMALL,loc);
120: return 0;
121: }
123: PetscErrorCode TSHistorySetHistory(TSHistory tsh, PetscInt n, PetscReal hist[], PetscInt hist_id[], PetscBool sorted)
124: {
128: PetscFree(tsh->hist);
129: PetscFree(tsh->hist_id);
130: tsh->n = (size_t) n;
131: tsh->c = (size_t) n;
132: PetscMalloc1(tsh->n,&tsh->hist);
133: PetscMalloc1(tsh->n,&tsh->hist_id);
134: for (PetscInt i = 0; i < (PetscInt)tsh->n; i++) {
135: tsh->hist[i] = hist[i];
136: tsh->hist_id[i] = hist_id ? hist_id[i] : i;
137: }
138: if (!sorted) PetscSortRealWithArrayInt(tsh->n,tsh->hist,tsh->hist_id);
139: tsh->sorted = PETSC_TRUE;
140: return 0;
141: }
143: PetscErrorCode TSHistoryGetHistory(TSHistory tsh, PetscInt *n, const PetscReal* hist[], const PetscInt* hist_id[], PetscBool *sorted)
144: {
145: if (n) *n = tsh->n;
146: if (hist) *hist = tsh->hist;
147: if (hist_id) *hist_id = tsh->hist_id;
148: if (sorted) *sorted = tsh->sorted;
149: return 0;
150: }
152: PetscErrorCode TSHistoryDestroy(TSHistory *tsh)
153: {
154: if (!*tsh) return 0;
155: PetscFree((*tsh)->hist);
156: PetscFree((*tsh)->hist_id);
157: PetscCommDestroy(&((*tsh)->comm));
158: PetscFree((*tsh));
159: *tsh = NULL;
160: return 0;
161: }
163: PetscErrorCode TSHistoryCreate(MPI_Comm comm, TSHistory *hst)
164: {
165: TSHistory tsh;
168: *hst = NULL;
169: PetscNew(&tsh);
170: PetscCommDuplicate(comm,&tsh->comm,NULL);
172: tsh->c = 1024; /* capacity */
173: tsh->s = 1024; /* reallocation size */
174: tsh->sorted = PETSC_TRUE;
176: PetscMalloc1(tsh->c,&tsh->hist);
177: PetscMalloc1(tsh->c,&tsh->hist_id);
178: *hst = tsh;
179: return 0;
180: }