ammarhakim / gkylzero

Lowest, compiled layer of Gkeyll infrastructure.
MIT License
22 stars 5 forks source link

Package BGK operator assembly in GK app by using zero/gkyl_bgk_collisions #287

Closed manauref closed 7 months ago

manauref commented 8 months ago

currently the BGK operator

(nu*f_M)_sum - nu_sum * f

is put together, and its cflrate computed, with these 4 lines in app/gk_species_bgk (L174-181)

  // Accumulate RHS of BGK collisions
  gkyl_dg_mul_conf_phase_op_range(&app->confBasis, &app->basis, bgk->nu_sum_f, 
    bgk->nu_sum, fin, &app->local, &species->local);  
  gkyl_array_accumulate(rhs, 1.0, bgk->nu_fmax);
  gkyl_array_accumulate(rhs, -1.0, bgk->nu_sum_f);

  // Determine the maximum collision frequency in each cell for finding the stable time-step
  gkyl_dg_calc_average_range(app->confBasis, 0, bgk->max_nu, 0, bgk->nu_sum, app->local);

Rather than having 4 lines, 4 kernel launches, 4 places for potential mistakes and possible doubts about whether cflrate is computed correctly (although it seems to me that it is for the nu=const case), it'd be better to simply use the gkyl_bgk_collisions object in zero/. That way we have 1 line, 1 kernel launch, only 1 place to make a mistake (in the app) and the cflfreq is also computed within it like it is for other objects. Further more, if we do that, these lines in apps/gk_species (L399-410)

  // If we are doing BGK collisions, find the maximum collision frequency and accumulate it
  // onto the CFL. Collision frequency is only a configuration space quantity, so do the reduce
  // here instead of utilizing species->cflrate, which is a phase space array. 
  if (species->collision_id == GKYL_BGK_COLLISIONS) {
    gkyl_array_reduce_range(species->omega_cfl, species->bgk.max_nu, GKYL_MAX, &app->local);
    double omega_cfl_bgk_ho[1];
    if (app->use_gpu)
      gkyl_cu_memcpy(omega_cfl_bgk_ho, species->omega_cfl, sizeof(double), GKYL_CU_MEMCPY_D2H);
    else
      omega_cfl_bgk_ho[0] = species->omega_cfl[0]; 
    omega_cfl_ho[0] += omega_cfl_bgk_ho[0];
  }

can be removed, and we'd be treating the BGK cflrate with the same code/interface we already treat the LBO cflrate.

JunoRavin commented 7 months ago

This issue was addressed by commit https://github.com/ammarhakim/gkylzero/commit/9b22e4c8d7a0fc65f324ffc6fd0258a5193882e0

to more accurately compute the stable time step in BGK collisions. The use of this operator is now uniform across GK and Vlasov BGK collisions within the app.