Actual source code: extchemfield.c

  1: static const char help[] = "Integrate chemistry using TChem.\n";

  3: #include <petscts.h>
  4: #include <petscdmda.h>

  6: #if defined(PETSC_HAVE_TCHEM)
  7: #if defined(MAX)
  8: #undef MAX
  9: #endif
 10: #if defined(MIN)
 11: #undef MIN
 12: #endif
 13: #  include <TC_params.h>
 14: #  include <TC_interface.h>
 15: #else
 16: #  error TChem is required for this example.  Reconfigure PETSc using --download-tchem.
 17: #endif
 18: /*

 20:     This is an extension of extchem.c to solve the reaction equations independently in each cell of a one dimensional field

 22:     Obtain the three files into this directory

 24:        curl http://combustion.berkeley.edu/gri_mech/version30/files30/grimech30.dat > chem.inp
 25:        curl http://combustion.berkeley.edu/gri_mech/version30/files30/thermo30.dat > therm.dat
 26:        cp $PETSC_DIR/$PETSC_ARCH/externalpackages/tchem/data/periodictable.dat .

 28:     Run with
 29:      ./extchemfield  -ts_arkimex_fully_implicit -ts_max_snes_failures -1 -ts_adapt_monitor -ts_adapt_dt_max 1e-4 -ts_arkimex_type 4 -ts_max_time .005

 31:      Options for visualizing the solution:
 32:         Watch certain variables in each cell evolve with time
 33:         -draw_solution 1 -ts_monitor_lg_solution_variables Temp,H2,O2,H2O,CH4,CO,CO2,C2H2,N2 -lg_use_markers false  -draw_pause -2

 35:         Watch certain variables in all cells evolve with time
 36:         -da_refine 4 -ts_monitor_draw_solution -draw_fields_by_name Temp,H2 -draw_vec_mark_points  -draw_pause -2

 38:         Keep the initial temperature distribution as one monitors the current temperature distribution
 39:         -ts_monitor_draw_solution_initial -draw_bounds .9,1.7 -draw_fields_by_name Temp

 41:         Save the images in a .gif (movie) file
 42:         -draw_save -draw_save_single_file

 44:         Compute the sensitivies of the solution of the first tempature on the initial conditions
 45:         -ts_adjoint_solve  -ts_dt 1.e-5 -ts_type cn -ts_adjoint_view_solution draw

 47:         Turn off diffusion
 48:         -diffusion no

 50:         Turn off reactions
 51:         -reactions no


 54:     The solution for component i = 0 is the temperature.

 56:     The solution, i > 0, is the mass fraction, massf[i], of species i, i.e. mass of species i/ total mass of all species

 58:     The mole fraction molef[i], i > 0, is the number of moles of a species/ total number of moles of all species
 59:         Define M[i] = mass per mole of species i then
 60:         molef[i] = massf[i]/(M[i]*(sum_j massf[j]/M[j]))

 62:     FormMoleFraction(User,massf,molef) converts the mass fraction solution of each species to the mole fraction of each species.

 64: */
 65: typedef struct _User *User;
 66: struct _User {
 67:   PetscReal pressure;
 68:   int       Nspec;
 69:   int       Nreac;
 70:   PetscReal Tini,dx;
 71:   PetscReal diffus;
 72:   DM        dm;
 73:   PetscBool diffusion,reactions;
 74:   double    *tchemwork;
 75:   double    *Jdense;        /* Dense array workspace where Tchem computes the Jacobian */
 76:   PetscInt  *rows;
 77: };

 79: static PetscErrorCode MonitorCell(TS,User,PetscInt);
 80: static PetscErrorCode FormRHSFunction(TS,PetscReal,Vec,Vec,void*);
 81: static PetscErrorCode FormRHSJacobian(TS,PetscReal,Vec,Mat,Mat,void*);
 82: static PetscErrorCode FormInitialSolution(TS,Vec,void*);

 84: #define TCCHKERRQ(ierr) do {if (ierr) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in TChem library, return code %d",ierr);} while (0)

 86: int main(int argc,char **argv)
 87: {
 88:   TS                ts;         /* time integrator */
 89:   TSAdapt           adapt;
 90:   Vec               X;          /* solution vector */
 91:   Mat               J;          /* Jacobian matrix */
 92:   PetscInt          steps,ncells,xs,xm,i;
 93:   PetscErrorCode    ierr;
 94:   PetscReal         ftime,dt;
 95:   char              chemfile[PETSC_MAX_PATH_LEN] = "chem.inp",thermofile[PETSC_MAX_PATH_LEN] = "therm.dat";
 96:   struct _User      user;
 97:   TSConvergedReason reason;
 98:   PetscBool         showsolutions = PETSC_FALSE;
 99:   char              **snames,*names;
100:   Vec               lambda;     /* used with TSAdjoint for sensitivities */

102:   PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
103:   PetscOptionsBegin(PETSC_COMM_WORLD,NULL,"Chemistry solver options","");
104:   PetscOptionsString("-chem","CHEMKIN input file","",chemfile,chemfile,sizeof(chemfile),NULL);
105:   PetscOptionsString("-thermo","NASA thermo input file","",thermofile,thermofile,sizeof(thermofile),NULL);
106:   user.pressure = 1.01325e5;    /* Pascal */
107:   PetscOptionsReal("-pressure","Pressure of reaction [Pa]","",user.pressure,&user.pressure,NULL);
108:   user.Tini   = 1550;
109:   PetscOptionsReal("-Tini","Initial temperature [K]","",user.Tini,&user.Tini,NULL);
110:   user.diffus = 100;
111:   PetscOptionsReal("-diffus","Diffusion constant","",user.diffus,&user.diffus,NULL);
112:   PetscOptionsBool("-draw_solution","Plot the solution for each cell","",showsolutions,&showsolutions,NULL);
113:   user.diffusion = PETSC_TRUE;
114:   PetscOptionsBool("-diffusion","Have diffusion","",user.diffusion,&user.diffusion,NULL);
115:   user.reactions = PETSC_TRUE;
116:   PetscOptionsBool("-reactions","Have reactions","",user.reactions,&user.reactions,NULL);
117:   PetscOptionsEnd();

119:   TC_initChem(chemfile, thermofile, 0, 1.0);TC
120:   user.Nspec = TC_getNspec();
121:   user.Nreac = TC_getNreac();

123:   DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_PERIODIC,10,user.Nspec+1,1,NULL,&user.dm);
124:   DMSetFromOptions(user.dm);
125:   DMSetUp(user.dm);
126:   DMDAGetInfo(user.dm,NULL,&ncells,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
127:   user.dx = 1.0/ncells;  /* Set the coordinates of the cell centers; note final ghost cell is at x coordinate 1.0 */
128:   DMDASetUniformCoordinates(user.dm,0.0,1.0,0.0,1.0,0.0,1.0);

130:   /* set the names of each field in the DMDA based on the species name */
131:   PetscMalloc1((user.Nspec+1)*LENGTHOFSPECNAME,&names);
132:   PetscStrcpy(names,"Temp");
133:   TC_getSnames(user.Nspec,names+LENGTHOFSPECNAME);
134:   PetscMalloc1((user.Nspec+2),&snames);
135:   for (i=0; i<user.Nspec+1; i++) snames[i] = names+i*LENGTHOFSPECNAME;
136:   snames[user.Nspec+1] = NULL;
137:   DMDASetFieldNames(user.dm,(const char * const *)snames);
138:   PetscFree(snames);
139:   PetscFree(names);


142:   DMCreateMatrix(user.dm,&J);
143:   DMCreateGlobalVector(user.dm,&X);

145:   PetscMalloc3(user.Nspec+1,&user.tchemwork,PetscSqr(user.Nspec+1),&user.Jdense,user.Nspec+1,&user.rows);

147:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
148:      Create timestepping solver context
149:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
150:   TSCreate(PETSC_COMM_WORLD,&ts);
151:   TSSetDM(ts,user.dm);
152:   TSSetType(ts,TSARKIMEX);
153:   TSARKIMEXSetFullyImplicit(ts,PETSC_TRUE);
154:   TSARKIMEXSetType(ts,TSARKIMEX4);
155:   TSSetRHSFunction(ts,NULL,FormRHSFunction,&user);
156:   TSSetRHSJacobian(ts,J,J,FormRHSJacobian,&user);

158:   ftime = 1.0;
159:   TSSetMaxTime(ts,ftime);
160:   TSSetExactFinalTime(ts,TS_EXACTFINALTIME_STEPOVER);

162:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
163:      Set initial conditions
164:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
165:   FormInitialSolution(ts,X,&user);
166:   TSSetSolution(ts,X);
167:   dt   = 1e-10;                 /* Initial time step */
168:   TSSetTimeStep(ts,dt);
169:   TSGetAdapt(ts,&adapt);
170:   TSAdaptSetStepLimits(adapt,1e-12,1e-4); /* Also available with -ts_adapt_dt_min/-ts_adapt_dt_max */
171:   TSSetMaxSNESFailures(ts,-1);            /* Retry step an unlimited number of times */


174:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
175:      Pass information to graphical monitoring routine
176:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
177:   if (showsolutions) {
178:     DMDAGetCorners(user.dm,&xs,NULL,NULL,&xm,NULL,NULL);
179:     for (i=xs;i<xs+xm;i++) {
180:       MonitorCell(ts,&user,i);
181:     }
182:   }

184:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
185:      Set runtime options
186:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
187:   TSSetFromOptions(ts);

189:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
190:      Set final conditions for sensitivities
191:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
192:   DMCreateGlobalVector(user.dm,&lambda);
193:   TSSetCostGradients(ts,1,&lambda,NULL);
194:   VecSetValue(lambda,0,1.0,INSERT_VALUES);
195:   VecAssemblyBegin(lambda);
196:   VecAssemblyEnd(lambda);

198:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
199:      Solve ODE
200:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
201:   TSSolve(ts,X);
202:   TSGetSolveTime(ts,&ftime);
203:   TSGetStepNumber(ts,&steps);
204:   TSGetConvergedReason(ts,&reason);
205:   PetscPrintf(PETSC_COMM_WORLD,"%s at time %g after %D steps\n",TSConvergedReasons[reason],(double)ftime,steps);

207:   {
208:     Vec                max;
209:     const char * const *names;
210:     PetscInt           i;
211:     const PetscReal    *bmax;

213:     TSMonitorEnvelopeGetBounds(ts,&max,NULL);
214:     if (max) {
215:       TSMonitorLGGetVariableNames(ts,&names);
216:       if (names) {
217:         VecGetArrayRead(max,&bmax);
218:         PetscPrintf(PETSC_COMM_SELF,"Species - maximum mass fraction\n");
219:         for (i=1; i<user.Nspec; i++) {
220:           if (bmax[i] > .01) {PetscPrintf(PETSC_COMM_SELF,"%s %g\n",names[i],bmax[i]);}
221:         }
222:         VecRestoreArrayRead(max,&bmax);
223:       }
224:     }
225:   }

227:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
228:      Free work space.
229:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
230:   TC_reset();
231:   DMDestroy(&user.dm);
232:   MatDestroy(&J);
233:   VecDestroy(&X);
234:   VecDestroy(&lambda);
235:   TSDestroy(&ts);
236:   PetscFree3(user.tchemwork,user.Jdense,user.rows);
237:   PetscFinalize();
238:   return ierr;
239: }

241: /*
242:    Applies the second order centered difference diffusion operator on a one dimensional periodic domain
243: */
244: static PetscErrorCode FormDiffusionFunction(TS ts,PetscReal t,Vec X,Vec F,void *ptr)
245: {
246:   User              user = (User)ptr;
247:   PetscErrorCode    ierr;
248:   PetscScalar       **f;
249:   const PetscScalar **x;
250:   DM                dm;
251:   PetscInt          i,xs,xm,j,dof;
252:   Vec               Xlocal;
253:   PetscReal         idx;

256:   TSGetDM(ts,&dm);
257:   DMDAGetInfo(dm,NULL,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
258:   DMGetLocalVector(dm,&Xlocal);
259:   DMGlobalToLocalBegin(dm,X,INSERT_VALUES,Xlocal);
260:   DMGlobalToLocalEnd(dm,X,INSERT_VALUES,Xlocal);
261:   DMDAVecGetArrayDOFRead(dm,Xlocal,&x);
262:   DMDAVecGetArrayDOF(dm,F,&f);
263:   DMDAGetCorners(dm,&xs,NULL,NULL,&xm,NULL,NULL);

265:   idx = 1.0*user->diffus/user->dx;
266:   for (i=xs; i<xs+xm; i++) {
267:     for (j=0; j<dof; j++) {
268:       f[i][j] += idx*(x[i+1][j] - 2.0*x[i][j] + x[i-1][j]);
269:     }
270:   }
271:   DMDAVecRestoreArrayDOFRead(dm,Xlocal,&x);
272:   DMDAVecRestoreArrayDOF(dm,F,&f);
273:   DMRestoreLocalVector(dm,&Xlocal);
274:   return(0);
275: }

277: /*
278:    Produces the second order centered difference diffusion operator on a one dimensional periodic domain
279: */
280: static PetscErrorCode FormDiffusionJacobian(TS ts,PetscReal t,Vec X,Mat Amat,Mat Pmat,void *ptr)
281: {
282:   User              user = (User)ptr;
283:   PetscErrorCode    ierr;
284:   DM                dm;
285:   PetscInt          i,xs,xm,j,dof;
286:   PetscReal         idx,values[3];
287:   MatStencil        row,col[3];

290:   TSGetDM(ts,&dm);
291:   DMDAGetInfo(dm,NULL,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
292:   DMDAGetCorners(dm,&xs,NULL,NULL,&xm,NULL,NULL);

294:   idx = 1.0*user->diffus/user->dx;
295:   values[0] = idx;
296:   values[1] = -2.0*idx;
297:   values[2] = idx;
298:   for (i=xs; i<xs+xm; i++) {
299:     for (j=0; j<dof; j++) {
300:       row.i = i;      row.c = j;
301:       col[0].i = i-1; col[0].c = j;
302:       col[1].i = i;   col[1].c = j;
303:       col[2].i = i+1; col[2].c = j;
304:       MatSetValuesStencil(Pmat,1,&row,3,col,values,ADD_VALUES);
305:     }
306:   }
307:   MatAssemblyBegin(Pmat,MAT_FINAL_ASSEMBLY);
308:   MatAssemblyEnd(Pmat,MAT_FINAL_ASSEMBLY);
309:   return(0);
310: }

312: static PetscErrorCode FormRHSFunction(TS ts,PetscReal t,Vec X,Vec F,void *ptr)
313: {
314:   User              user = (User)ptr;
315:   PetscErrorCode    ierr;
316:   PetscScalar       **f;
317:   const PetscScalar **x;
318:   DM                dm;
319:   PetscInt          i,xs,xm;

322:   if (user->reactions) {
323:     TSGetDM(ts,&dm);
324:     DMDAVecGetArrayDOFRead(dm,X,&x);
325:     DMDAVecGetArrayDOF(dm,F,&f);
326:     DMDAGetCorners(dm,&xs,NULL,NULL,&xm,NULL,NULL);

328:     for (i=xs; i<xs+xm; i++) {
329:       PetscArraycpy(user->tchemwork,x[i],user->Nspec+1);
330:       user->tchemwork[0] *= user->Tini; /* Dimensionalize */
331:       TC_getSrc(user->tchemwork,user->Nspec+1,f[i]);TC
332:       f[i][0] /= user->Tini;           /* Non-dimensionalize */
333:     }

335:     DMDAVecRestoreArrayDOFRead(dm,X,&x);
336:     DMDAVecRestoreArrayDOF(dm,F,&f);
337:   } else {
338:     VecZeroEntries(F);
339:   }
340:   if (user->diffusion) {
341:     FormDiffusionFunction(ts,t,X,F,ptr);
342:   }
343:   return(0);
344: }

346: static PetscErrorCode FormRHSJacobian(TS ts,PetscReal t,Vec X,Mat Amat,Mat Pmat,void *ptr)
347: {
348:   User              user = (User)ptr;
349:   PetscErrorCode    ierr;
350:   const PetscScalar **x;
351:   PetscInt          M = user->Nspec+1,i,j,xs,xm;
352:   DM                dm;

355:   if (user->reactions) {
356:     TSGetDM(ts,&dm);
357:     MatZeroEntries(Pmat);
358:     MatSetOption(Pmat,MAT_ROW_ORIENTED,PETSC_FALSE);
359:     MatSetOption(Pmat,MAT_IGNORE_ZERO_ENTRIES,PETSC_TRUE);
360:     DMDAVecGetArrayDOFRead(dm,X,&x);
361:     DMDAGetCorners(dm,&xs,NULL,NULL,&xm,NULL,NULL);

363:     for (i=xs; i<xs+xm; i++) {
364:       PetscArraycpy(user->tchemwork,x[i],user->Nspec+1);
365:       user->tchemwork[0] *= user->Tini;  /* Dimensionalize temperature (first row) because that is what Tchem wants */
366:       TC_getJacTYN(user->tchemwork,user->Nspec,user->Jdense,1);

368:       for (j=0; j<M; j++) user->Jdense[j + 0*M] /= user->Tini; /* Non-dimensionalize first column */
369:       for (j=0; j<M; j++) user->Jdense[0 + j*M] /= user->Tini; /* Non-dimensionalize first row */
370:       for (j=0; j<M; j++) user->rows[j] = i*M+j;
371:       MatSetValues(Pmat,M,user->rows,M,user->rows,user->Jdense,INSERT_VALUES);
372:     }
373:     DMDAVecRestoreArrayDOFRead(dm,X,&x);
374:     MatAssemblyBegin(Pmat,MAT_FINAL_ASSEMBLY);
375:     MatAssemblyEnd(Pmat,MAT_FINAL_ASSEMBLY);
376:   } else {
377:     MatZeroEntries(Pmat);
378:   }
379:   if (user->diffusion) {
380:     FormDiffusionJacobian(ts,t,X,Amat,Pmat,ptr);
381:   }
382:   if (Amat != Pmat) {
383:     MatAssemblyBegin(Amat,MAT_FINAL_ASSEMBLY);
384:     MatAssemblyEnd(Amat,MAT_FINAL_ASSEMBLY);
385:   }
386:   return(0);
387: }

389: PetscErrorCode FormInitialSolution(TS ts,Vec X,void *ctx)
390: {
391:   PetscScalar    **x,*xc;
393:   struct {const char *name; PetscReal massfrac;} initial[] = {
394:     {"CH4", 0.0948178320887},
395:     {"O2", 0.189635664177},
396:     {"N2", 0.706766236705},
397:     {"AR", 0.00878026702874}
398:   };
399:   PetscInt       i,j,xs,xm;
400:   DM             dm;

403:   VecZeroEntries(X);
404:   TSGetDM(ts,&dm);
405:   DMDAGetCorners(dm,&xs,NULL,NULL,&xm,NULL,NULL);

407:   DMDAGetCoordinateArray(dm,&xc);
408:   DMDAVecGetArrayDOF(dm,X,&x);
409:   for (i=xs; i<xs+xm; i++) {
410:     x[i][0] = 1.0 + .05*PetscSinScalar(2.*PETSC_PI*xc[i]);  /* Non-dimensionalized by user->Tini */
411:     for (j=0; j<sizeof(initial)/sizeof(initial[0]); j++) {
412:       int ispec = TC_getSpos(initial[j].name, strlen(initial[j].name));
413:       if (ispec < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Could not find species %s",initial[j].name);
414:       PetscPrintf(PETSC_COMM_SELF,"Species %d: %s %g\n",j,initial[j].name,initial[j].massfrac);
415:       x[i][1+ispec] = initial[j].massfrac;
416:     }
417:   }
418:   DMDAVecRestoreArrayDOF(dm,X,&x);
419:   DMDARestoreCoordinateArray(dm,&xc);
420:   return(0);
421: }

423: /*
424:     Routines for displaying the solutions
425: */
426: typedef struct {
427:   PetscInt cell;
428:   User     user;
429: } UserLGCtx;

431: static PetscErrorCode FormMoleFraction(UserLGCtx *ctx,Vec massf,Vec *molef)
432: {
433:   User              user = ctx->user;
434:   PetscErrorCode    ierr;
435:   PetscReal         *M,tM=0;
436:   PetscInt          i,n = user->Nspec+1;
437:   PetscScalar       *mof;
438:   const PetscScalar **maf;

441:   VecCreateSeq(PETSC_COMM_SELF,n,molef);
442:   PetscMalloc1(user->Nspec,&M);
443:   TC_getSmass(user->Nspec, M);
444:   DMDAVecGetArrayDOFRead(user->dm,massf,&maf);
445:   VecGetArray(*molef,&mof);
446:   mof[0] = maf[ctx->cell][0]; /* copy over temperature */
447:   for (i=1; i<n; i++) tM += maf[ctx->cell][i]/M[i-1];
448:   for (i=1; i<n; i++) {
449:     mof[i] = maf[ctx->cell][i]/(M[i-1]*tM);
450:   }
451:   DMDAVecRestoreArrayDOFRead(user->dm,massf,&maf);
452:   VecRestoreArray(*molef,&mof);
453:   PetscFree(M);
454:   return(0);
455: }

457: static PetscErrorCode MonitorCellDestroy(UserLGCtx *uctx)
458: {

462:   PetscFree(uctx);
463:   return(0);
464: }

466: /*
467:    Use TSMonitorLG to monitor the reactions in a particular cell
468: */
469: static PetscErrorCode MonitorCell(TS ts,User user,PetscInt cell)
470: {
472:   TSMonitorLGCtx ctx;
473:   char           **snames;
474:   UserLGCtx      *uctx;
475:   char           label[128];
476:   PetscReal      temp,*xc;
477:   PetscMPIInt    rank;

480:   DMDAGetCoordinateArray(user->dm,&xc);
481:   temp = 1.0 + .05*PetscSinScalar(2.*PETSC_PI*xc[cell]);  /* Non-dimensionalized by user->Tini */
482:   DMDARestoreCoordinateArray(user->dm,&xc);
483:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
484:   PetscSNPrintf(label,sizeof(label),"Initial Temperature %g Cell %d Rank %d",(double)user->Tini*temp,(int)cell,rank);
485:   TSMonitorLGCtxCreate(PETSC_COMM_SELF,NULL,label,PETSC_DECIDE,PETSC_DECIDE,600,400,1,&ctx);
486:   DMDAGetFieldNames(user->dm,(const char * const **)&snames);
487:   TSMonitorLGCtxSetVariableNames(ctx,(const char * const *)snames);
488:   PetscNew(&uctx);
489:   uctx->cell = cell;
490:   uctx->user = user;
491:   TSMonitorLGCtxSetTransform(ctx,(PetscErrorCode (*)(void*,Vec,Vec*))FormMoleFraction,(PetscErrorCode (*)(void*))MonitorCellDestroy,uctx);
492:   TSMonitorSet(ts,TSMonitorLGSolution,ctx,(PetscErrorCode (*)(void**))TSMonitorLGCtxDestroy);
493:   return(0);
494: }