288 {
289 V rho, u[DESCRIPTOR::d],eta[DESCRIPTOR::q],uPlus[DESCRIPTOR::d],tmp_cell[(DESCRIPTOR::q+1)/2];
290 MomentaF().computeRhoU(cell, rho, u);
291 const V omega =
parameters.template get<descriptors::OMEGA>();
292 V tmpMomentumLoss[DESCRIPTOR::d]{ };
293 auto velNumerator = cell.template getFieldPointer<descriptors::VELOCITY_NUMERATOR>();
294 auto zeta = cell.template getFieldPointer<descriptors::ZETA>();
295 auto velDenominator = cell.template getField<descriptors::VELOCITY_DENOMINATOR>();
296
297 if (velDenominator > 1) {
298 rho /= velDenominator;
299 }
300 for (int tmp_iPop=1; 2*tmp_iPop<DESCRIPTOR::q; tmp_iPop++) {
301 eta[tmp_iPop]=6.*descriptors::t<V,DESCRIPTOR>(tmp_iPop)*rho*(descriptors::c<DESCRIPTOR>(tmp_iPop,0)*(velNumerator[0])+descriptors::c<DESCRIPTOR>(tmp_iPop,1)*(velNumerator[1]));
302 tmp_cell[tmp_iPop]=(zeta[tmp_iPop]*(-cell[tmp_iPop]+cell[descriptors::opposite<DESCRIPTOR>(tmp_iPop)]+eta[tmp_iPop]));
303 cell[tmp_iPop]+=tmp_cell[tmp_iPop]/(1.+2.*(zeta[tmp_iPop]));
304 cell[descriptors::opposite<DESCRIPTOR>(tmp_iPop)]-=tmp_cell[tmp_iPop]/(1.+2.*(zeta[tmp_iPop]));
305 zeta[tmp_iPop] = 0.;
306 zeta[descriptors::opposite<DESCRIPTOR>(tmp_iPop)] = 0.;
307 }
308
309 cell.template setField<descriptors::POROSITY>(1.);
310 cell.template setField<descriptors::VELOCITY_DENOMINATOR>(0);
311
312 MomentaF().computeRhoU(cell, rho, uPlus);
314
315 V diff[DESCRIPTOR::q] = {};
316 for (int tmp_iPop=0; tmp_iPop<DESCRIPTOR::q; tmp_iPop++) {
317 diff[tmp_iPop] +=
EquilibriumF().compute(tmp_iPop, rho, uPlus)
319
320 for (int iDim=0; iDim<DESCRIPTOR::d; iDim++) {
321 tmpMomentumLoss[iDim] -= descriptors::c<DESCRIPTOR>(tmp_iPop,iDim) * diff[tmp_iPop];
322 }
323 }
324
325 for (int i_dim=0; i_dim<DESCRIPTOR::d; i_dim++) {
326 velNumerator[i_dim] = tmpMomentumLoss[i_dim];
327 }
328
329 return {rho, util::normSqr<V,DESCRIPTOR::d>(u)};
330 }
typename EQUILIBRIUM::template type< DESCRIPTOR, MOMENTA > EquilibriumF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
typename meta::list< descriptors::OMEGA > parameters
static V bgkCollision(CELL &cell, const RHO &rho, const VELOCITY &u, const OMEGA &omega) any_platform
BGK collision step.