HARM
harm and utilities
 All Data Structures Files Functions Variables Typedefs Macros Pages
restart.checks.c
Go to the documentation of this file.
1 
7 #include "decs.h"
8 
9 
10 static int restart_init_point_check_pglobal(int which, int i, int j, int k);
11 static int restart_init_point_check_unewglobal(int which, int i, int j, int k);
12 static int restart_init_point_check_pstagglobal(int which, int i, int j, int k);
13 
14 
16 int restart_init_simple_checks(int which)
17 {
18  int gotnan;
19  int i,j,k;
20 
22  //
23  // make sure all zones are not nan just as read-in from file
24  //
26  // make sure all zones are not nan
27  gotnan=0;
28  // OPENMPOPTMARK: Don't optimize since critical region
29 
30  if(which<=2){ // see initbase.c
31  LOOP{ gotnan+=restart_init_point_check_pglobal(which,i,j,k); }
32  LOOP{ gotnan+=restart_init_point_check_unewglobal(which,i,j,k); }
33  }
34  else if(which<=3){
35  FULLLOOP{ gotnan+=restart_init_point_check_pglobal(which,i,j,k); }
36  LOOP{ gotnan+=restart_init_point_check_unewglobal(which,i,j,k); }
37  }
38  else{
39  FULLLOOP{ gotnan+=restart_init_point_check_pglobal(which,i,j,k); }
40  LOOP{ gotnan+=restart_init_point_check_unewglobal(which,i,j,k); }
41  FULLLOOP{ gotnan+=restart_init_point_check_pstagglobal(which,i,j,k); }
42  }
43 
44 
45 
46  if(gotnan) myexit(39476346);
47 
48  return(0);
49 
50 }
51 
52 
54 static int restart_init_point_check_pglobal(int which, int i, int j, int k)
55 {
56  int pliter,pl;
57  int gotnan;
58 
59  gotnan=0;
60 
61  PDUMPLOOP(pliter,pl){
62  if(!finite(GLOBALMACP0A1(pglobal,i,j,k,pl)) ){
63  dualfprintf(fail_file,"restart_init(%d): restart data has NaN at i=%d j=%d k=%d ti=%d tj=%d tk=%d :: pl=%d : pglobal=%21.15g\n",which,i,j,k,startpos[1]+i,startpos[2]+j,startpos[3]+k,pl,GLOBALMACP0A1(pglobal,i,j,k,pl));
64  if(SCALARPL(pl)){
65  dualfprintf(fail_file,"scalar went nan, reset to floor: pl=%d\n",pl);
66  GLOBALMACP0A1(pglobal,i,j,k,pl)=NUMEPSILON;
67  }
68  else{
69  // myexit(24968346);
70  gotnan++;
71  }
72  }
73  }
74 
75  return(gotnan);
76 
77 }
78 
79 
81 static int restart_init_point_check_unewglobal(int which, int i, int j, int k)
82 {
83  int pliter,pl;
84  int gotnan;
85 
86  gotnan=0;
87 
88  PDUMPLOOP(pliter,pl){
89  if(DOENOFLUX != NOENOFLUX || (FLUXB==FLUXCTSTAG &&(pl==B1 && i>=-N1BND+SHIFT1 || pl==B2 && j>=-N2BND+SHIFT2 || pl==B3 && k>=-N3BND+SHIFT3)) ){
90  if(!finite(GLOBALMACP0A1(unewglobal,i,j,k,pl)) ){
91  dualfprintf(fail_file,"restart_init(%d): restart data has NaN at i=%d j=%d k=%d ti=%d tj=%d tk=%d :: pl=%d : unewglobal=%21.15g\n",which,i,j,k,startpos[1]+i,startpos[2]+j,startpos[3]+k,pl,GLOBALMACP0A1(unewglobal,i,j,k,pl));
92  // myexit(24968346);
93  gotnan++;
94  }
95  }
96  }
97 
98  return(gotnan);
99 
100 }
101 
102 
104 static int restart_init_point_check_pstagglobal(int which, int i, int j, int k)
105 {
106  int pliter,pl;
107  int gotnan;
108 
109  gotnan=0;
110 
111 
112  if(FLUXB==FLUXCTSTAG){
113  PDUMPLOOP(pliter,pl){
114  if(pl==B1 && i>=-N1BND+SHIFT1 || pl==B2 && j>=-N2BND+SHIFT2 || pl==B3 && k>=-N3BND+SHIFT3){
115  if(!finite(GLOBALMACP0A1(pstagglobal,i,j,k,pl)) ){
116  dualfprintf(fail_file,"restart_init(%d): restart data has NaN at i=%d j=%d k=%d ti=%d tj=%d tk=%d :: pl=%d : pstagglobal=%21.15g\n",which,i,j,k,startpos[1]+i,startpos[2]+j,startpos[3]+k,pl,GLOBALMACP0A1(pstagglobal,i,j,k,pl));
117  // myexit(24968346);
118  gotnan++;
119  }
120  }
121  }
122  }
123 
124 
125  return(gotnan);
126 
127 }
128 
129 
130 
131 
141 int restart_init_checks(int which, FTYPE (*prim)[NSTORE2][NSTORE3][NPR], FTYPE (*pstag)[NSTORE2][NSTORE3][NPR], FTYPE (*ucons)[NSTORE2][NSTORE3][NPR])
142 {
143  char ans[100];
144  struct of_geom geomdontuse;
145  struct of_geom *ptrgeom=&geomdontuse;
146  int failreturn;
147  FTYPE ucon[NDIM];
148  FTYPE utmax=0.0;
149  int i,j,k,pl,pliter;
150  int failflag=0;
151  extern int set_dt(FTYPE (*prim)[NSTORE2][NSTORE3][NPR], SFTYPE *dt);
152  int gotnan;
153 
154 
155 
157  //
158  // NOW have all parameters and data.
159  //
161 
162 
163 
165  //
166  // make sure all zones are not nan before fixup and bound
167  //
169  // make sure all zones are not nan
170  gotnan=0;
171  LOOP{// OPENMPOPTMARK: Don't optimize since critical region
172  PDUMPLOOP(pliter,pl){
173  if(!finite(MACP0A1(prim,i,j,k,pl))){
174  dualfprintf(fail_file,"before fixup & bound: restart data has NaN at i=%d j=%d k=%d ti=%d tj=%d tk=%d :: pl=%d\n",i,j,k,startpos[1]+i,startpos[2]+j,startpos[3]+k,pl);
175  // myexit(24968346);
176  gotnan=1;
177  }
178  }
179  }
180  if(gotnan) myexit(24968341);
181 
182 
183 
185  // NOW CHECK THE DATA and apply "fixups" if necessary
186  //
188 
189 #if(CHECKRHONEGZERORESTART)
190 
191  // report if data has rho<0 or u<0
192  ZLOOP{
193  if(STEPOVERNEGRHO!=NEGDENSITY_NEVERFIXUP){
194  if(MACP0A1(prim,i,j,k,RHO)<=0.0){
195  dualfprintf(fail_file,"restart data has negative mass density at i=%d j=%d k=%d\n",startpos[1]+i,startpos[2]+j,startpos[3]+k);
196  failflag++;
197  }
198  }
199  if(STEPOVERNEGU!=NEGDENSITY_NEVERFIXUP){
200  if(MACP0A1(prim,i,j,k,UU)<=0.0){
201  dualfprintf(fail_file,"restart data has negative ie density at i=%d j=%d k=%d\n",startpos[1]+i,startpos[2]+j,startpos[3]+k);
202  failflag++;
203  }
204  }
205 
206  // check 3-velocity
207 #if(WHICHVEL==VEL3)
208  if(jonchecks){
209  get_geometry(i,j,k,CENT,ptrgeom);
210  failreturn=check_pr(MAC(prim,i,j,k),MAC(prim,i,j,k),MAC(ucons,i,j,k), ptrgeom,-2,-1);
211  if(failreturn){
212  dualfprintf(fail_file,"restart data has large or imaginary u^t=%21.15g at i=%d j=%d k=%d. I will attempt to correct.\n",1.0/sqrt(uttdiscr),startpos[1]+i,startpos[2]+j,startpos[3]+k);
213  }
214  if(1.0/sqrt(uttdiscr)>utmax) utmax=1.0/sqrt(uttdiscr);
215  // need to settle over limit u^t's
216  failreturn=check_pr(MAC(prim,i,j,k),MAC(prim,i,j,k),MAC(ucons,i,j,k), ptrgeom,-1,-1);
217  if(failreturn){
218  dualfprintf(fail_file,"restart data has imaginary u^t at i=%d j=%d k=%d. Unable to correct.\n",startpos[1]+i,startpos[2]+j,startpos[3]+k);
219  return(1);
220  }
221  }
222 #endif
223  }
224  if(WHICHVEL==VEL3){
225  if(jonchecks){
226  dualfprintf(fail_file,"max u^t of restart data=%21.15g\n",utmax);
227  }
228  }
229 #endif
230 
231  if(failflag>0){
232  dualfprintf(fail_file,"Restart data has at least %d failures, please correct or accept that fixup will fix them.\n",failflag);
233  //myexit(1); // aborts MPI
234  // return(1); // currently doesn't abort properly in MPI
235  }
236 
237 
238 
239 
241  //
242  // fixup() during restart
243  //
245 #if(FIXUPAFTERRESTART)
246  if(fixup(STAGEM1,prim,ucons,-1)>=1)
247  FAILSTATEMENT("restart.c:restart_init()", "fixup()", 1);
248 
249  trifprintf( "proc: %d fixup restart completed: failed=%d\n", myid,failed);
250 #endif
251 
252 
253 
255  //
256  // make sure all zones are not nan after fixup
257  //
259  // make sure all zones are not nan
260  gotnan=0;
261  LOOP{// OPENMPOPTMARK: Don't optimize since critical region
262  PDUMPLOOP(pliter,pl){
263  if(!finite(MACP0A1(prim,i,j,k,pl))){
264  dualfprintf(fail_file,"after fixup & before bound: restart data has NaN at i=%d j=%d k=%d ti=%d tj=%d tk=%d :: pl=%d\n",i,j,k,startpos[1]+i,startpos[2]+j,startpos[3]+k,pl);
265  // myexit(24968346);
266  gotnan=1;
267  }
268  }
269  }
270  if(gotnan) myexit(24968341);
271 
272 
273 
275  //
276  // BOUND during restart
277  //
279  int finalstep=1; // user would want to know about changes to conserved quants during restart
280  if (bound_allprim(STAGEM1,finalstep,t,prim,pstag,ucons, USEMPI) >= 1) {
281  dualfprintf(fail_file, "restart_init:bound_allprim: failure\n");
282  fflush(fail_file);
283  return (1);
284  }
285 
286  trifprintf( "proc: %d bound restart completed: failed=%d\n", myid,failed);
287 
288 
290  //
291  // make sure all zones are not nan after bound
292  //
294  // make sure all zones are not nan
295  gotnan=0;
296  FULLLOOP{// OPENMPOPTMARK: Don't optimize since critical region
297  PDUMPLOOP(pliter,pl){
298  if(!finite(MACP0A1(prim,i,j,k,pl))){
299  dualfprintf(fail_file,"after fixup & bound: restart data has NaN at i=%d j=%d k=%d ti=%d tj=%d tk=%d :: pl=%d\n",i,j,k,startpos[1]+i,startpos[2]+j,startpos[3]+k,pl);
300  // myexit(24968346);
301  gotnan=1;
302  }
303  }
304  }
305  if(gotnan) myexit(24968346);
306 
308  //
309  // pre_fixup() during restart
310  //
312 
313  if(pre_fixup(STAGEM1,prim)>=1)
314  FAILSTATEMENT("init.c:init()", "postbc_fixup()", 1);
315 
316 
317 
319  //
320  // make sure all zones are good now
321  //
323 #if(CHECKRHONEGZERORESTART)
324  failflag=0;
325 
326  if(STEPOVERNEGRHO!=NEGDENSITY_NEVERFIXUP){
327  FULLLOOP{
328  if(MACP0A1(prim,i,j,k,RHO)<=0.0){
329  dualfprintf(fail_file,"restart data has negative mass density at i=%d j=%d k=%d : %21.15g\n",startpos[1]+i,startpos[2]+j,startpos[3]+k,MACP0A1(prim,i,j,k,RHO));
330  // return(1);
331  failflag++;
332  }
333  }// end fullloop
334  }// end negrho check
335 
336  if(STEPOVERNEGU!=NEGDENSITY_NEVERFIXUP){
337  FULLLOOP{
338  if(MACP0A1(prim,i,j,k,UU)<=0.0){
339  dualfprintf(fail_file,"restart data has negative ie density at i=%d j=%d k=%d : %21.15g\n",startpos[1]+i,startpos[2]+j,startpos[3]+k,MACP0A1(prim,i,j,k,UU));
340  // return(1);
341  failflag++;
342  }
343  }// end fullloop
344  }
345 
346  if(failflag>0){
347  dualfprintf(fail_file,"Restart data has at least %d failures -- even after fixup should have been applied!\n",failflag);
348  myexit(1); // aborts MPI (abort, not return, so no stall on other cores)
349  }
350 
351 #endif
352 
353 
354 
356  //
357  // test read by looking at images
358  //
360  if(image_dump(-3)>=1) return(1);
361 
362 
363 
365  //
366  // Now adjust timestep since may be dt=0 or whatever if last run ended on tf
367  // NOW done whether restarting or not at end of init()
368  //
370  // set_dt(prim,&dt);
371 
372 
373  trifprintf("end restart_init_checks\n");
374 
375  /* done! */
376  return (0);
377 
378 }
Generated on Fri May 20 2016 15:52:35 for HARM by   doxygen 1.8.3.1