collide.c

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

//##############################################################################
//
// collide.c
//

void collide( lattice_ptr lattice)
{
  double *f;

  double omega;

  int    *bc_type;

  int    n, a;

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

  //dump_pdf( lattice, 998);

  f = lattice->pdf[0].f;

  bc_type = &( lattice->bc[0].bc_type);

  for( n=0; n<lattice->NumNodes; n++, f+=18, bc_type++)
  {
    if( !( *bc_type & BC_SOLID_NODE))
    {
      if( *bc_type == 0)
      {
        // C O L L I D E 

        for( a=0; a<=8; a++, f++)
        {
          // f = f - (1/tau)( ftemp - feq)
//printf("Collide (Before): n = %2d, a = %d, "
//       "*f = %10.7f, *(f+9) = %10.7f, *(f-9) = %10.7f\n",
//        n, a, *(f), *(f+9), *(f-9));


          *f = *(f+9) - ( ( *(f+9) / lattice->param.tau )
                        - ( *(f-9) / lattice->param.tau ) );

#if PUKE_NEGATIVE_DENSITIES
          if( *f < 0.)
          {
            printf("\n");
            printf(
              "collide() -- Node %d (%d,%d) has negative density %20.17f "
              "in direction %d "
              "at timestep %d. Exiting!\n", 
              n, lattice->node[n].i, 
                 lattice->node[n].j, *f, a,
                 lattice->time             );
            printf("\n");
            exit(1);
          }
#endif /* PUKE_NEGATIVE_DENSITIES */

//printf("Collide (After ): n = %2d, a = %d, "
//       "*f = %10.7f, *(f+9) = %10.7f, *(f-9) = %10.7f\n",
//        n, a, *(f-1), *(f-1+9), *(f-1-9));

        } /* for( a=0; a<=8; a++) */
      }
      else
      {
//printf("collide() -- Skipping bc %d at n = %d\n", *bc_type, n);
        *f++ = *( f + 9);
        *f++ = *( f + 9);
        *f++ = *( f + 9);
        *f++ = *( f + 9);
        *f++ = *( f + 9);
        *f++ = *( f + 9);
        *f++ = *( f + 9);
        *f++ = *( f + 9);
        *f++ = *( f + 9);
      }

    } /* if( !( *bc_type++ & BC_SOLID_NODE)) */

    else // *bc_type++ & BC_SOLID_NODE
    {
      // B O U N C E B A C K 

      f++; // Skip rest particle.

      *f++ = *( f + 9 + 2); //f++; // f[1] = ftemp[3]
      *f++ = *( f + 9 + 2); //f++; // f[2] = ftemp[4]
      *f++ = *( f + 9 - 2); //f++; // f[3] = ftemp[1]
      *f++ = *( f + 9 - 2); //f++; // f[4] = ftemp[2]
      *f++ = *( f + 9 + 2); //f++; // f[5] = ftemp[7]
      *f++ = *( f + 9 + 2); //f++; // f[6] = ftemp[8]
      *f++ = *( f + 9 - 2); //f++; // f[7] = ftemp[5]
      *f++ = *( f + 9 - 2); //f++; // f[8] = ftemp[6]

//printf("Bncback: n = %2d\n", n);

    } /* if( !( *bc_type++ & BC_SOLID_NODE)) else */

  } /* for( n=0; n<lattice_NumNodes; n++) */

  //dump_pdf( lattice, 999);

#if SAY_HI
  printf("collide() -- Bye!\n");
#endif /* SAY_HI */

} /* void collide( lattice_ptr lattice) */