lbio.c

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

    {
      k = fread( &rgb, sizeof(struct rgb_quad), 1, in );
      if( k!=1)
      {
        printf("Error reading palette entry %d.  Exiting!\n", i);
        exit(1);
      }
    }
  }

  // Bytes per row of the bitmap.
  bytes_per_row = 
    ((int)ceil(( (((double)(ENDIAN4(*width_ptr)))*((double)(ENDIAN2(*bitcount_ptr))))/8.)));

  // Bitmaps pad rows to preserve 4-byte boundaries.
  // The length of a row in the file will be bytes_per_row + pad .
  pad = ((4) - bytes_per_row%4)%4;

  n = 0;
  m = 0;
  n+=( k = fread( &b, 1, 1, in ));
  i = 0;
  //j = *height_ptr-1;
  j = 0;//*height_ptr-1;
  while( !feof(in))
  {
    switch(ENDIAN2(*bitcount_ptr))
    {
      case 1: // Monochrome.
        printf("read_bcs() -- "
            "Support for Monochrome BMPs is pending.  "
            "Exiting!\n");
        exit(1);

        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b & 0x80) == 0); }
        i++;     
        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b & 0x40) == 0); }
        i++;     
        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b & 0x20) == 0); }
        i++;     
        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b & 0x10) == 0); }
        i++;     
        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b & 0x08) == 0); }
        i++;     
        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b & 0x04) == 0); }
        i++;     
        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b & 0x02) == 0); }
        i++;     
        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b & 0x01) == 0); }
        i++;
        break;

      case 4: // 16 colors.
        printf("read_bcs() -- "
            "Support for 16 color BMPs is pending.  "
            "Exiting!\n");
        exit(1);

        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b&0xf0)>>4 != 15); }
        i++;
        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b&0x0f) != 15); }
        i++;
        break;

      case 8: // 256 colors.
        printf("read_bcs() -- "
            "Support for 256 color BMPs is pending.  "
            "Exiting!\n");
        exit(1);

        if( i < ENDIAN4(*width_ptr)) { (bcs)[j][i] = ( (b&0xff) != 255); }
        i++;
        break;

      case 24: // 24-bit colors.
        if( i < 3*(ENDIAN4(*width_ptr))) 
        { 
          i++; n+=( k = fread( &g, 1, 1, in ));
          i++; n+=( k = fread( &r, 1, 1, in ));
          if( ( (b&0xff) == 0) &&( (g&0xff) == 0) &&( (r&0xff) == 255) )
          { // R E D ==> Inflow, Pressure boundaries.
            bcs[j][(int)floor((double)i/3.)] = 1;
          }
          if( ( (b&0xff) == 0) &&( (g&0xff) == 255) &&( (r&0xff) == 0) )
          { // G R E E N ==> Outflow, Pressure boundaries.
            bcs[j][(int)floor((double)i/3.)] = 2;
          }

        }
        i++;
        break;

      default: // 32-bit colors?
        printf("ERROR: Unhandled color depth, "
            "BitCount = %d. Exiting!\n", ENDIAN2(*bitcount_ptr));
        exit(1);
        break;

    } /* switch(*(bmih.biBitCount)) */

    if( !(n%(bytes_per_row+pad))) { m++; i=0; j++;}
    n+=( k = fread( &b, 1, 1, in ));

  } /* while( !feof(in)) */

  if( (bytes_per_row+pad)*m!=n)
  {
    printf("WARNING: Num bytes read = %d versus num bytes predicted = %d .\n",
        n, (bytes_per_row+pad)*m);
  }

  if( m != ENDIAN4(*height_ptr))
  {
    printf("WARNING: m (%d) != bmih.biHeight (%d).\n", m, ENDIAN4(*height_ptr));
  }

  fclose(in);

  ei = lattice->param.LX-1;
  ej = lattice->param.LY-1;

  for( n=0; n<lattice->NumNodes; n++)
  {
    i = n%lattice->param.LX;
    j = n/lattice->param.LX;

    if( bcs[ j][ i] != 0)
    {
//printf("read_bcs() -- n = %d, ( %d, %d) of ( %d, %d).\n", n, i, j, ei, ej);
#if 0
      if(    ( i==0  && j==0 )
          || ( i==ei && j==0 )
          || ( i==ei && j==ej)
          || ( i==0  && j==ej)  )
      {
        // Skip corners for now.
        printf("read_bcs() -- WARNING: Skipping corner ( %d, %d).", i, j);
      }
      else
      {
#endif
#if 0
        if( i==0)
        {
          // West
          if( bcs[ j][ i] == 1)
          {
            // Inflow
            lattice->bc[n].bc_type |= BC_PRESSURE_W_IN;
            //lattice->periodic_x = 0;
          }
          else if( bcs[ j][ i] == 2)
          {
            // Outflow
            lattice->bc[n].bc_type |= BC_PRESSURE_W_OUT;
            //lattice->periodic_x = 0;
          }
          else
          {
            // Unhandled case.
            printf("read_bcs() -- Unhandled case: "
                "bcs[ %d][ %d] = %d .  Exiting!\n", 
                i, j, bcs[j][i]);
            exit(1);
          }

        }
        else if( i==ei)
        {
          // East
          if( bcs[ j][ i] == 1)
          {
            // Inflow
            lattice->bc[n].bc_type |= BC_PRESSURE_E_IN;
            //lattice->periodic_x = 0;
          }
          else if( bcs[ j][ i] == 2)
          {
            // Outflow
            lattice->bc[n].bc_type |= BC_PRESSURE_E_OUT;
            //lattice->periodic_x = 0;
          }
          else
          {
            // Unhandled case.
            printf("read_bcs() -- Unhandled case: "
                "bcs[ %d][ %d] = %d .  Exiting!\n", 
                i, j, bcs[j][i]);
            exit(1);
          }

        }
        else 
#endif
        if( j==0)
        {
//printf("read_bcs() -- South at i=%d\n", i);
          // South
          if( bcs[ j][ i] == 1)
          {
            // Inflow
            lattice->bc[subs][n].bc_type |= BC_PRESSURE_S_IN;
            //lattice->periodic_y = 0;
          }
          else if( bcs[ j][ i] == 2)
          {
            // Outflow
            lattice->bc[subs][n].bc_type |= BC_PRESSURE_S_OUT;
            //lattice->periodic_y = 0;
          }
          else
          {
            // Unhandled case.
            printf("read_bcs() -- Unhandled case: "
                "bcs[ %d][ %d] = %d .  Exiting!\n", 
                i, j, bcs[j][i]);
            exit(1);
          }

        }
        else if( j==ej)
        {
//printf("read_bcs() -- North at i=%d\n", i);
          // North
          if( bcs[ j][ i] == 1)
          {
            // Inflow
            lattice->bc[subs][n].bc_type |= BC_PRESSURE_N_IN;
            //lattice->periodic_y = 0;
          }
          else if( bcs[ j][ i] == 2)
          {
            // Outflow
            lattice->bc[subs][n].bc_type |= BC_PRESSURE_N_OUT;
            //lattice->periodic_y = 0;
          }
          else
          {
            // Unhandled case.
            printf("read_bcs() -- Unhandled case: "
                "bcs[ %d][ %d] = %d .  Exiting!\n", 
                i, j, bcs[j][i]);
            exit(1);
          }

        }
        else
        {
          // Unhandled case.
          printf("read_bcs() -- WARNING: "
              "Support for interior flow bcs is pending! "
              "Skipping ( i, j) = ( %d, %d).\n", i, j);
        }
#if 0
      }
#endif

    }

  }

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

  printf("read_bcs() -- Bye!\n");
  printf("\n");

} /* read_bcs( char *filename, int ***bcs, int *height, int *width) */

// void rho2bmp( char *filename, int time)
//##############################################################################
//
// R H O 2 B M P 
//
#if 1
void rho2bmp( lattice_ptr lattice, int time)
{
  FILE   *in, 
         *o;
  int    i, j, 
         n, m;
  int    pad, 
         bytes_per_row;
  int    frame;
  char   k;
  char   b;
  struct bitmap_file_header bmfh;
  struct bitmap_info_header bmih;
  struct rgb_quad rgb;
  int    *int_ptr;
  short  int *short_int_ptr;
  int    *width_ptr;
  int    *height_ptr;
  short  int *bitcount_ptr;
  char   filename[1024];
  char   red_val, 
         green_val, 
         blue_val, 
         val;
  double fval;
  double min_rho, max_rho;
  int    subs;
  double **colormap;
  int    num_colors;

#if SAY_HI
  printf("rho2bmp() -- Hi!\n");
#endif /* SAY_HI */

  if( lattice->param.use_colormap)
  {
    count_colormap(      &num_colors);
    allocate_colormap(   &colormap, num_colors);
    read_colormap(       colormap, num_colors);
  }

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

  frame = time/lattice->param.FrameRate;

  sprintf( filename, "./in/%dx%d.bmp", lattice->param.LX, lattice->param.LY);
  if( !( in = fopen( filename, "r")))
  {
    printf("rho2bmp() -- Error opening file \"%s\".\n", filename);
    exit(1);
  }

  // n = fread( void *BUF, size_t SIZE, size_t COUNT, FILE *FP);

  n = fread( &bmfh, sizeof(struct bitmap_file_header), 1, in );
  if( strncmp(bmfh.bfType,"BM",2))
  {
    printf("ERROR: Can't process this file type.  Exiting!\n");
    printf("\n");
    exit(1);
  }
  n = fread( &bmih, sizeof(struct bitmap_info_header), 1, in );
  int_ptr = (int*)bmih.biCompression;
  if( *int_ptr != 0)
  {
    printf("ERROR: Can't handle compression.  Exiting!\n");
    printf("\n");
    exit(1);
  }

  width_ptr = (int*)bmih.biWidth;
  height_ptr = (int*)bmih.biHeight;
  bitcount_ptr = (short int*)bmih.biBitCount;

  // Read palette entries, if applicable.
  if( ENDIAN2(*bitcount_ptr) < 24)
  {
    n = (int)pow(2.,(double)ENDIAN2(*bitcount_ptr)); // Num palette entries.
    for( i=0; i<n; i++)
    {
      k = fread( &rgb, sizeof(struct rgb_quad), 1, in );
      if( k!=1)
      {
        printf("Error reading palette entry %d.  Exiting!\n", i);
        exit(1);
      }
    }
  }

  fclose(in);

  // Bytes per row of the bitmap.
  bytes_per_row = 
    ((int)ceil(( (((double)(ENDIAN4(*width_ptr)))*((double)(ENDIAN2(*bitcount_ptr))))/8.)));

  // Bitmaps pad rows to preserve 4-byte boundaries.
  // The length of a row in the file will be bytes_per_row + pad .
  pad = ((4) - bytes_per_row%4)%4;

  compute_min_rho( lattice, &min_rho, subs);
  compute_max_rho( lattice, &max_rho, subs);

  sprintf( filename, "./out/rho%dx%d_frame%04d_subs%02d.bmp", 
      lattice->param.LX, 
      lattice->param.LY, 
      frame, subs);
  if( !( o = fopen( filename, "w+")))
  {
    printf("ERROR: fopen( \"%s\", \"w+\") = NULL.  Bye, bye!\n", filename);
    exit(1);
  }

  fwrite( &bmfh, sizeof(struct bitmap_file_header), 1, o );
  fwrite( &bmih, sizeof(struct bitmap_info_header), 1, o );

  for( j=0; j<lattice->param.LY; j++)
  {
    n = j*lattice->param.LX;

    for( i=0; i<lattice->param.LX; i++, n++)
    {
      if( lattice->bc[subs][ n].bc_type == /*FLUID_NODE*/0)
      {
        if( lattice->param.use_colormap)
        {
          if( lattice->param.plot_scale_dynamic)
          {
            if( max_rho!=min_rho)
            {
              get_color(
                colormap,
                num_colors,
                (lattice->macro_vars[subs][ n].rho - min_rho)/(max_rho-min_rho),
                &red_val,
                &green_val,
                &blue_val   );
            }
            else
            {
              get_color(
                colormap,
                num_colors,
                1.,
                &red_val,
                &green_val,
                &blue_val   );
            }
          }
          else
          {
            get_color(
              colormap,
              num_colors,
              (lattice->macro_vars[subs][ n].rho
             /( (lattice->param.rho_A[subs]>lattice->param.rho_B[subs])
               ?(lattice->param.rho_A[subs])
               :(lattice->param.rho_B[subs]) )),
              &red_val,
              &green_val,
              &blue_val   );
          }

        }
        else
        {
          if( subs==0)
          {
            if( lattice->param.plot_scale_dynamic)
            {
              if( max_rho!=min_rho)