Actual source code: tshistory.c

petsc-3.11.4 2019-09-28
Report Typos and Errors
  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:     PetscErrorCode _7_ierr;                                             \
  9:     PetscInt b1[2],b2[2];                                               \
 10:     b1[0] = -b; b1[1] = b;                                              \
 11:     _7_MPIU_Allreduce(b1,b2,2,MPIU_INT,MPI_MAX,a);CHKERRQ(_7_ierr); \
 12:     if (-b2[0] != b2[1]) SETERRQ1(a,PETSC_ERR_ARG_WRONG,"Int value must be same on all processes, argument # %d",c); \
 13:   } while (0)

 16:   do {                                                                  \
 17:     PetscErrorCode _7_ierr;                                             \
 18:     PetscMPIInt b1[2],b2[2];                                            \
 19:     b1[0] = -(PetscMPIInt)b; b1[1] = (PetscMPIInt)b;                    \
 20:     _7_MPIU_Allreduce(b1,b2,2,MPI_INT,MPI_MAX,a);CHKERRQ(_7_ierr); \
 21:     if (-b2[0] != b2[1]) SETERRQ1(a,PETSC_ERR_ARG_WRONG,"Bool value must be same on all processes, argument # %d",c); \
 22:   } while (0)

 25:   do {                                                                  \
 26:     PetscErrorCode _7_ierr;                                             \
 27:     PetscReal b1[3],b2[3];                                              \
 28:     if (PetscIsNanReal(b)) {b1[2] = 1;} else {b1[2] = 0;};              \
 29:     b1[0] = -b; b1[1] = b;                                              \
 30:     _7_MPI_Allreduce(b1,b2,3,MPIU_REAL,MPIU_MAX,a);CHKERRQ(_7_ierr); \
 31:     if (!(b2[2] > 0) && !PetscEqualReal(-b2[0],b2[1])) SETERRQ1(a,PETSC_ERR_ARG_WRONG,"Real value must be same on all processes, argument # %d",c); \
 32:   } while (0)

 34: #else


 40: #endif

 42: struct _n_TSHistory {
 43:   MPI_Comm  comm;     /* used for runtime collective checks */
 44:   PetscReal *hist;    /* time history */
 45:   PetscInt  *hist_id; /* stores the stepid in time history */
 46:   PetscInt  n;        /* current number of steps registered */
 47:   PetscBool sorted;   /* if the history is sorted in ascending order */
 48:   PetscInt  c;        /* current capacity of hist */
 49:   PetscInt  s;        /* reallocation size */
 50: };

 52: PetscErrorCode TSHistoryGetNumSteps(TSHistory tsh, PetscInt *n)
 53: {
 56:   *n = tsh->n;
 57:   return(0);
 58: }

 60: PetscErrorCode TSHistoryUpdate(TSHistory tsh, PetscInt id, PetscReal time)
 61: {

 67:   if (tsh->n == tsh->c) { /* reallocation */
 68:     tsh->c += tsh->s;
 69:     PetscRealloc(tsh->c*sizeof(*tsh->hist),&tsh->hist);
 70:     PetscRealloc(tsh->c*sizeof(*tsh->hist_id),&tsh->hist_id);
 71:   }
 72:   tsh->sorted = (PetscBool)(tsh->sorted && (tsh->n ? time >= tsh->hist[tsh->n-1] : PETSC_TRUE));
 73: #if defined(PETSC_USE_DEBUG)
 74:   if (tsh->n) { /* id should be unique */
 75:     PetscInt loc,*ids;

 77:     PetscMalloc1(tsh->n,&ids);
 78:     PetscMemcpy(ids,tsh->hist_id,tsh->n*sizeof(PetscInt));
 79:     PetscSortInt(tsh->n,ids);
 80:     PetscFindInt(id,tsh->n,ids,&loc);
 81:     PetscFree(ids);
 82:     if (loc >=0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"History id should be unique");
 83:   }
 84: #endif
 85:   tsh->hist[tsh->n]    = time;
 86:   tsh->hist_id[tsh->n] = id;
 87:   tsh->n += 1;
 88:   return(0);
 89: }

 91: PetscErrorCode TSHistoryGetTime(TSHistory tsh, PetscBool backward, PetscInt step, PetscReal *t)
 92: {
 96:   if (!t) return(0);
 98:   if (!tsh->sorted) {

101:     PetscSortRealWithArrayInt(tsh->n,tsh->hist,tsh->hist_id);
102:     tsh->sorted = PETSC_TRUE;
103:   }
104:   if (step < 0 || step >= tsh->n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Given time step %D does not match any in history [0,%D]",step,tsh->n);
105:   if (!backward) *t = tsh->hist[step];
106:   else           *t = tsh->hist[tsh->n-step-1];
107:   return(0);
108: }

110: PetscErrorCode TSHistoryGetTimeStep(TSHistory tsh, PetscBool backward, PetscInt step, PetscReal *dt)
111: {
115:   if (!dt) return(0);
117:   if (!tsh->sorted) {

120:     PetscSortRealWithArrayInt(tsh->n,tsh->hist,tsh->hist_id);
121:     tsh->sorted = PETSC_TRUE;
122:   }
123:   if (step < 0 || step > tsh->n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Given time step %D does not match any in history [0,%D]",step,tsh->n);
124:   if (!backward) *dt = tsh->hist[PetscMin(step+1,tsh->n-1)] - tsh->hist[PetscMin(step,tsh->n-1)];
125:   else           *dt = tsh->hist[PetscMax(tsh->n-step-1,0)] - tsh->hist[PetscMax(tsh->n-step-2,0)];
126:   return(0);
127: }

129: PetscErrorCode TSHistoryGetLocFromTime(TSHistory tsh, PetscReal time, PetscInt *loc)
130: {

136:   if (!tsh->sorted) {
137:     PetscSortRealWithArrayInt(tsh->n,tsh->hist,tsh->hist_id);
138:     tsh->sorted = PETSC_TRUE;
139:   }
140:   PetscFindReal(time,tsh->n,tsh->hist,PETSC_SMALL,loc);
141:   return(0);
142: }

144: PetscErrorCode TSHistorySetHistory(TSHistory tsh, PetscInt n, PetscReal hist[], PetscInt hist_id[], PetscBool sorted)
145: {
146:   PetscInt       i;

152:   PetscFree(tsh->hist);
153:   PetscFree(tsh->hist_id);
154:   tsh->n = n;
155:   tsh->c = n;
156:   PetscMalloc1(tsh->n,&tsh->hist);
157:   PetscMalloc1(tsh->n,&tsh->hist_id);
158:   for (i = 0; i < tsh->n; i++) {
159:     tsh->hist[i]    = hist[i];
160:     tsh->hist_id[i] = hist_id ? hist_id[i] : i;
161:   }
162:   if (!sorted) {
163:     PetscSortRealWithArrayInt(tsh->n,tsh->hist,tsh->hist_id);
164:   }
165:   tsh->sorted = PETSC_TRUE;
166:   return(0);
167: }

169: PetscErrorCode TSHistoryGetHistory(TSHistory tsh, PetscInt *n, const PetscReal* hist[], const PetscInt* hist_id[], PetscBool *sorted)
170: {
172:   if (n)             *n = tsh->n;
173:   if (hist)       *hist = tsh->hist;
174:   if (hist_id) *hist_id = tsh->hist_id;
175:   if (sorted)   *sorted = tsh->sorted;
176:   return(0);
177: }

179: PetscErrorCode TSHistoryDestroy(TSHistory *tsh)
180: {

184:   if (!*tsh) return(0);
185:   PetscFree((*tsh)->hist);
186:   PetscFree((*tsh)->hist_id);
187:   PetscCommDestroy(&((*tsh)->comm));
188:   PetscFree((*tsh));
189:   *tsh = NULL;
190:   return(0);
191: }

193: PetscErrorCode TSHistoryCreate(MPI_Comm comm, TSHistory *hst)
194: {
195:   TSHistory      tsh;

200:   *hst = NULL;
201:   PetscNew(&tsh);
202:   PetscCommDuplicate(comm,&tsh->comm,NULL);

204:   tsh->c      = 1024; /* capacity */
205:   tsh->s      = 1024; /* reallocation size */
206:   tsh->sorted = PETSC_TRUE;

208:   PetscMalloc1(tsh->c,&tsh->hist);
209:   PetscMalloc1(tsh->c,&tsh->hist_id);
210:   *hst = tsh;
211:   return(0);
212: }