lbio_03082005.c

基于格子Boltzmann方法开源可视化软件的源代码 具有很高的实用价值。对学习LBM方法及其软件开发非常游泳

  int    *node_ptr;
  int    n;
  double *macro_vars_ptr;
  int    frame;
  double max_u[2],  ave_u[2];
  double max_rho,   ave_rho;
  double rho_ratio, u_x_ratio, u_y_ratio;
  int    subs;

 for( subs=0; subs<NUM_FLUID_COMPONENTS; subs++)
 {

  frame = (int)((double)time/(double)lattice->param.FrameRate);
#if VERBOSITY_LEVEL > 0
  printf("read_macro_vars() -- frame = %d/%d = %d\n", 
      time,
      lattice->param.FrameRate,
      frame);
#endif /* VERBOSITY_LEVEL > 0 */

  // R E A D   R H O   A N D   U 
  //
  //  - Read the density and velocity values at the active nodes to
  //    the rho and u dat files.
  //
  sprintf( filename, "./out/rho_frame%04d_subs%02d.dat", frame, subs);
  if( !( rho_in = fopen( filename, "r+")))
  {
    printf("ERROR: fopen( \"%s\", \"r+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  sprintf( filename, "./out/u_frame%04d_subs%02d.dat", frame, subs);
  if( !( u_in = fopen( filename, "r+")))
  {
    printf("ERROR: fopen( \"%s\", \"r+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  macro_vars_ptr = &( lattice->macro_vars[subs][0].rho);

  for( n=0; n<lattice->NumNodes; n++)
  {

    fscanf( rho_in, "%lf\n", macro_vars_ptr++);

    fscanf( u_in,   "%lf ",  macro_vars_ptr++);

    fscanf( u_in,   "%lf\n", macro_vars_ptr++);

  }

  fclose(u_in);
#if VERBOSITY_LEVEL > 0
  sprintf( filename, "./out/u_frame%04d_subs%02d.dat", frame, subs);
#endif /* VERBOSITY_LEVEL > 0 */
  printf("read_macro_vars() -- Read file \"%s\"\n", filename);
  fclose(rho_in);
#if VERBOSITY_LEVEL > 0
  sprintf( filename, "./out/rho_frame%04d_subs%02d.dat", frame, subs);
#endif /* VERBOSITY_LEVEL > 0 */
  printf("read_macro_vars() -- Read file \"%s\"\n", filename);
  
 } /* for( subs=0; subs<NUM_FLUID_COMPONENTS; subs++) */

} /* void read_macro_vars( struct lattice_struct *lattice, int time) */

#endif
// void dump_pdf( struct lattice_struct *lattice, int time)
//##############################################################################
//
// D U M P   P D F 
//
//  - Output the particle distribution functions to a text file.
//
//  - This is useful mainly for debugging with small problems.
//
void dump_pdf( struct lattice_struct *lattice, int time)
{
  char   filename[1024];
  FILE   *o_feq, 
         *o_f, 
         *o_ftemp;
  double *fptr, 
         *end_ptr;
  bc_ptr bc;
  int    frame;
  int    subs;
#if WRITE_PDF_TO_TXT
  int    i, j, n;
#endif /* WRITE_PDF_TO_TXT */

 for( subs=0; subs<NUM_FLUID_COMPONENTS; subs++)
 {

  frame = time/lattice->param.FrameRate;

  sprintf( filename, "./out/feq_frame%04d_subs%02d.dat", frame, subs);
  if( !( o_feq = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  sprintf( filename, "./out/f_frame%04d_subs%02d.dat", frame, subs);
  if( !( o_f = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  sprintf( filename, "./out/ftemp_frame%04d_subs%02d.dat", frame, subs);
  if( !( o_ftemp = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  bc      =   lattice->bc[subs];
  fptr    =   lattice->pdf[subs][0].feq;
  end_ptr = &(lattice->pdf[subs][ lattice->NumNodes-1].ftemp[8]) + 1;
  while( fptr!=end_ptr)
  {
    if( 1 || !( bc++->bc_type & BC_SOLID_NODE))
    {
      fprintf( o_feq  , "%10.7f ", *fptr++);
      fprintf( o_feq  , "%10.7f ", *fptr++);
      fprintf( o_feq  , "%10.7f ", *fptr++);
      fprintf( o_feq  , "%10.7f ", *fptr++);
      fprintf( o_feq  , "%10.7f ", *fptr++);
      fprintf( o_feq  , "%10.7f ", *fptr++);
      fprintf( o_feq  , "%10.7f ", *fptr++);
      fprintf( o_feq  , "%10.7f ", *fptr++);
      fprintf( o_feq  , "%10.7f ", *fptr++);
    }
    else
    {
      fprintf( o_feq  , "%10.7f ", 0.);
      fprintf( o_feq  , "%10.7f ", 0.);
      fprintf( o_feq  , "%10.7f ", 0.);
      fprintf( o_feq  , "%10.7f ", 0.);
      fprintf( o_feq  , "%10.7f ", 0.);
      fprintf( o_feq  , "%10.7f ", 0.);
      fprintf( o_feq  , "%10.7f ", 0.);
      fprintf( o_feq  , "%10.7f ", 0.);
      fprintf( o_feq  , "%10.7f ", 0.);
      fptr+=9;
    }
    fprintf( o_feq  , "\n");

    fprintf( o_f    , "%10.7f ", *fptr++);
    fprintf( o_f    , "%10.7f ", *fptr++);
    fprintf( o_f    , "%10.7f ", *fptr++);
    fprintf( o_f    , "%10.7f ", *fptr++);
    fprintf( o_f    , "%10.7f ", *fptr++);
    fprintf( o_f    , "%10.7f ", *fptr++);
    fprintf( o_f    , "%10.7f ", *fptr++);
    fprintf( o_f    , "%10.7f ", *fptr++);
    fprintf( o_f    , "%10.7f ", *fptr++);
    fprintf( o_f    , "\n");

    fprintf( o_ftemp, "%10.7f ", *fptr++);
    fprintf( o_ftemp, "%10.7f ", *fptr++);
    fprintf( o_ftemp, "%10.7f ", *fptr++);
    fprintf( o_ftemp, "%10.7f ", *fptr++);
    fprintf( o_ftemp, "%10.7f ", *fptr++);
    fprintf( o_ftemp, "%10.7f ", *fptr++);
    fprintf( o_ftemp, "%10.7f ", *fptr++);
    fprintf( o_ftemp, "%10.7f ", *fptr++);
    fprintf( o_ftemp, "%10.7f ", *fptr++);
    fprintf( o_ftemp, "\n");

  } /* while( fptr!=end_ptr) */

  fclose( o_feq);
  fclose( o_f);
  fclose( o_ftemp);

 } /* for( subs=0; subs<NUM_FLUID_COMPONENTS; subs++) */

#if WRITE_PDF_TO_TXT
 for( subs=0; subs<NUM_FLUID_COMPONENTS; subs++)
 {

  frame = time/lattice->param.FrameRate;

  sprintf( filename, "./out/feq_frame%04d_subs%02d.txt", frame, subs);
  if( !( o_feq = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  sprintf( filename, "./out/f_frame%04d_subs%02d.txt", frame, subs);
  if( !( o_f = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  sprintf( filename, "./out/ftemp_frame%04d_subs%02d.txt", frame, subs);
  if( !( o_ftemp = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  for( i=0; i<lattice->param.LX; i++)
  {
    fprintf( o_feq  , "-----------");
    fprintf( o_feq  , "-----------");
    fprintf( o_feq  , "----------|");
    fprintf( o_f    , "-----------");
    fprintf( o_f    , "-----------");
    fprintf( o_f    , "----------|");
    fprintf( o_ftemp, "-----------");
    fprintf( o_ftemp, "-----------");
    fprintf( o_ftemp, "----------|");
  }
  fprintf( o_feq  , "\n");
  fprintf( o_f    , "\n");
  fprintf( o_ftemp, "\n");

  for( j=lattice->param.LY-1; j>=0; j--)
  {
    n = j*lattice->param.LX;
    for( i=0; i<lattice->param.LX; i++, n++)
    {
      fprintf( o_feq,  "%10.6f ", lattice->pdf[subs][n].feq[6]  );
      fprintf( o_feq,  "%10.6f ", lattice->pdf[subs][n].feq[2]  );
      fprintf( o_feq,  "%10.6f|", lattice->pdf[subs][n].feq[5]  );
      fprintf( o_f,    "%10.6f ", lattice->pdf[subs][n].f[6]    );
      fprintf( o_f,    "%10.6f ", lattice->pdf[subs][n].f[2]    );
      fprintf( o_f,    "%10.6f|", lattice->pdf[subs][n].f[5]    );
      fprintf( o_ftemp,"%10.6f ", lattice->pdf[subs][n].ftemp[6]);
      fprintf( o_ftemp,"%10.6f ", lattice->pdf[subs][n].ftemp[2]);
      fprintf( o_ftemp,"%10.6f|", lattice->pdf[subs][n].ftemp[5]);
    } /* for( i=0; i<lattice->param.LX; i++, n++) */

    fprintf( o_feq,  "\n");
    fprintf( o_f,    "\n");
    fprintf( o_ftemp,"\n");

    n = j*lattice->param.LX;
    for( i=0; i<lattice->param.LX; i++, n++)
    {
      fprintf( o_feq,  "%10.6f ", lattice->pdf[subs][n].feq[3]  );
      fprintf( o_feq,  "%10.6f ", lattice->pdf[subs][n].feq[0]  );
      fprintf( o_feq,  "%10.6f|", lattice->pdf[subs][n].feq[1]  );
      fprintf( o_f,    "%10.6f ", lattice->pdf[subs][n].f[3]    );
      fprintf( o_f,    "%10.6f ", lattice->pdf[subs][n].f[0]    );
      fprintf( o_f,    "%10.6f|", lattice->pdf[subs][n].f[1]    );
      fprintf( o_ftemp,"%10.6f ", lattice->pdf[subs][n].ftemp[3]);
      fprintf( o_ftemp,"%10.6f ", lattice->pdf[subs][n].ftemp[0]);
      fprintf( o_ftemp,"%10.6f|", lattice->pdf[subs][n].ftemp[1]);
    } /* for( i=0; i<lattice->param.LX; i++, n++) */

    fprintf( o_feq,  "\n");
    fprintf( o_f,    "\n");
    fprintf( o_ftemp,"\n");

    n = j*lattice->param.LX;
    for( i=0; i<lattice->param.LX; i++, n++)
    {
      fprintf( o_feq,  "%10.6f ", lattice->pdf[subs][n].feq[7]  );
      fprintf( o_feq,  "%10.6f ", lattice->pdf[subs][n].feq[4]  );
      fprintf( o_feq,  "%10.6f|", lattice->pdf[subs][n].feq[8]  );
      fprintf( o_f,    "%10.6f ", lattice->pdf[subs][n].f[7]    );
      fprintf( o_f,    "%10.6f ", lattice->pdf[subs][n].f[4]    );
      fprintf( o_f,    "%10.6f|", lattice->pdf[subs][n].f[8]    );
      fprintf( o_ftemp,"%10.6f ", lattice->pdf[subs][n].ftemp[7]);
      fprintf( o_ftemp,"%10.6f ", lattice->pdf[subs][n].ftemp[4]);
      fprintf( o_ftemp,"%10.6f|", lattice->pdf[subs][n].ftemp[8]);
    } /* for( i=0; i<lattice->param.LX; i++, n++) */

    fprintf( o_feq,  "\n");
    fprintf( o_f,    "\n");
    fprintf( o_ftemp,"\n");

    for( i=0; i<lattice->param.LX; i++)
    {
      fprintf( o_feq  , "-----------");
      fprintf( o_feq  , "-----------");
      fprintf( o_feq  , "----------|");
      fprintf( o_f    , "-----------");
      fprintf( o_f    , "-----------");
      fprintf( o_f    , "----------|");
      fprintf( o_ftemp, "-----------");
      fprintf( o_ftemp, "-----------");
      fprintf( o_ftemp, "----------|");
    }
    fprintf( o_feq  , "\n");
    fprintf( o_f    , "\n");
    fprintf( o_ftemp, "\n");


  } /* for( j=lattice->param.LY-1; j>=0; j--) */

  fclose( o_feq);
  fclose( o_f);
  fclose( o_ftemp);

 } /* for( subs=0; subs<NUM_FLUID_COMPONENTS; subs++) */


#endif /* WRITE_PDF_TO_TXT */

} /* void dump_pdf( struct lattice_struct *lattice, int time) */

#if NON_LOCAL_FORCES
// void dump_forces( struct lattice_struct *lattice)
//##############################################################################
// 
// D U M P   F O R C E S
//
//  - Output the interactive force values to file.
//
void dump_forces( struct lattice_struct *lattice)
{
  char   filename[1024];
  FILE   *ox, *oy;
  int    n;
  double *force;
  int    frame;
#if WRITE_RHO_AND_U_TO_TXT
  int    i, j;
#endif /* WRITE_RHO_AND_U_TO_TXT */
  int    subs;

 for( subs = 0; subs < NUM_FLUID_COMPONENTS; subs++)
 {
  frame = (int)((double)lattice->time/(double)lattice->param.FrameRate);
#if VERBOSITY_LEVEL > 0
  printf("dump_forces() -- frame = %d/%d = %d\n", 
      lattice->time,
      lattice->param.FrameRate,
      frame);
#endif /* VERBOSITY_LEVEL > 0 */

  sprintf( filename, "./out/force_x_frame%04d_subs%02d.dat", frame, subs);
  if( !( ox = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }
  sprintf( filename, "./out/force_y_frame%04d_subs%02d.dat", frame, subs);
  if( !( oy = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  force = lattice->force[subs][0].force;

  for( n=0; n<lattice->NumNodes; n++)
  {
    fprintf( ox, "%20.17f\n", *force++);
    fprintf( oy, "%20.17f\n", *force++);
    force += ( sizeof( struct force_struct)/8 - 2);
  }

  fclose(ox);
#if VERBOSITY_LEVEL > 0
  sprintf( filename, "./out/force_x_frame%04d_subs%02d.dat", frame, subs);
#endif /* VERBOSITY_LEVEL > 0 */
  printf("dump_forces() -- Wrote file \"%s\"\n", filename);
  fclose(oy);
#if VERBOSITY_LEVEL > 0
  sprintf( filename, "./out/force_y_frame%04d_subs%02d.dat", frame, subs);
#endif /* VERBOSITY_LEVEL > 0 */
  printf("dump_forces() -- Wrote file \"%s\"\n", filename);

#if WRITE_RHO_AND_U_TO_TXT
  // NOTE: This is very inefficient.  But it's only intended
  // for debugging purposes on small problems.
  sprintf( filename, "./out/force_x_frame%04d_subs%02d.txt", frame, subs);
  if( !( ox = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }
  sprintf( filename, "./out/force_y_frame%04d_subs%02d.txt", frame, subs);
  if( !( oy = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }
  for( j=0; j<lattice->param.LY; j++)
  {
    n = j*lattice->param.LX;
    for( i=0; i<lattice->param.LX; i++, n++)
    {
      fprintf( ox, "%10.7f ", lattice->force[subs][n].force[0]);
      fprintf( oy, "%10.7f ", lattice->force[subs][n].force[1]);
      if( n==lattice->NumNodes)
      {
        fprintf( ox, "%10.7f ", 0.);
        fprintf( oy, "%10.7f ", 0.);
      }
    }
    fprintf( ox, "\n");
    fprintf( oy, "\n");
  }

  fclose(ox);
#if VERBOSITY_LEVEL > 0
  sprintf( filename, "./out/force_x_frame%04d_subs%02d.txt", frame, subs);
#endif /* VERBOSITY_LEVEL > 0 */
  printf("dump_forces() -- Wrote file \"%s\"\n", filename);
  fclose(oy);
#if VERBOSITY_LEVEL > 0
  sprintf( filename, "./out/force_y_frame%04d_subs%02d.txt", frame, subs);
#endif /* VERBOSITY_LEVEL > 0 */
  printf("dump_forces() -- Wrote file \"%s\"\n", filename);

#endif /* WRITE_RHO_AND_U_TO_TXT */

 } /* for( subs = 0; subs < NUM_FLUID_COMPONENTS; subs++) */

} /* void dump_forces( struct lattice_struct *lattice) */
#endif /* NON_LOCAL_FORCES */

// void dump_checkpoint( struct lattice_struct *lattice, int time, char *fn)
//##############################################################################
//
// D U M P   C H E C K P O I N T 
//
//  - Write lattice to a checkpoint file.
//
//  - Should be binary and store all information necessary to
//    restart the current run at this point.
//
void dump_checkpoint( struct lattice_struct *lattice, int time, char *fn)
{

} /* void dump_checkpoint( struct lattice_struct *lattice, ...) */

// void read_checkpoint( struct lattice_struct *lattice)
//##############################################################################
//
// R E A D   C H E C K P O I N T 
//
//  - Read lattice from a checkpoint file (as written by dump_checkpoint).
//
//  - With this information, should be able to restart where
//    the previous run stopped.
//
void read_checkpoint( struct lattice_struct *lattice)
{

} /* void read_checkpoint( struct lattice_struct *lattice) */

// void spy_bmp( char *filename, int ***spy)
//##############################################################################