49 using V =
typename CELL::value_t;
52 auto j = cell.template getField<descriptors::VELOCITY>();
54 cell.template setField<descriptors::VELOCITY>(j);
57 const V rho0 =
parameters.template get<RHO0>();
58 const V omega =
parameters.template get<OMEGA>();
60 const V rho = cell.template getFieldComponent<descriptors::STATISTIC>(0) * rho0;
62 const V rhoTot = rho * omega;
65 auto u = cell.template getField<descriptors::VELOCITY>();
66 uTot = (u*rho0*omega) / rhoTot;
70 const V psi = POTENTIAL().compute(rho,
parameters);
72 for (
int iPop=0; iPop < DESCRIPTOR::q; ++iPop) {
73 const V rho_ = cell.neighbor(descriptors::c<DESCRIPTOR>(iPop))
74 .template getFieldComponent<descriptors::STATISTIC>(0)
76 const V psi_ = POTENTIAL().compute(rho_,
parameters);
77 rhoContribution += psi * psi_
78 * descriptors::c<DESCRIPTOR>(iPop)
79 * descriptors::t<V,DESCRIPTOR>(iPop);
85 auto force = cell.template getField<descriptors::EXTERNAL_FORCE>();
86 cell.template setField<descriptors::VELOCITY>(uTot);
87 cell.template setField<descriptors::FORCE>(force - g*rhoContribution/rho);