248 {
249 using DESCRIPTOR = typename CELLS::template value_t<names::NavierStokes>::descriptor_t;
250 using DESCRIPTOR_ADE = typename CELLS::template value_t<names::Concentration0>::descriptor_t;
251
252 auto u = cells.template get<names::Concentration0>().template getField<descriptors::VELOCITY>();
254 cells.template get<names::NavierStokes>().computeAllMomenta(rho, u.data(), pi);
255 cells.template get<names::Concentration0>().template setField<descriptors::VELOCITY>(u);
256 cells.template get<names::Concentration1>().template setField<descriptors::VELOCITY>(u);
257 cells.template get<names::Concentration2>().template setField<descriptors::VELOCITY>(u);
258
259 T PiNeqNormSqr = pi[0]*pi[0] + 2.0*pi[1]*pi[1] + pi[2]*pi[2];
261 PiNeqNormSqr += pi[2]*pi[2] + pi[3]*pi[3] + 2*pi[4]*pi[4] +pi[5]*pi[5];
262 }
265 T tau_mol_NS = T{1} /
parameters.template get<OMEGA_NSE>();
266 auto omegasAD =
parameters.template get<OMEGAS_ADE>();
267 T tau_mol_AD0 = T{1} / omegasAD[0];
268 T tau_mol_AD1 = T{1} / omegasAD[1];
269 T tau_mol_AD2 = T{1} / omegasAD[2];
271 T smagoPrefactor =
parameters.template get<SMAGORINSKY_PREFACTOR>();
272 T tau_turb_NS = T{0.5}*(
util::sqrt(tau_mol_NS*tau_mol_NS + smagoPrefactor/rho*PiNeqNorm) - tau_mol_NS);
273
275 T tauTurbADPrefactor0 = descriptors::invCs2<T,DESCRIPTOR_ADE>() / descriptors::invCs2<T,DESCRIPTOR>() / Sc[0];
276 T tauTurbADPrefactor1 = descriptors::invCs2<T,DESCRIPTOR_ADE>() / descriptors::invCs2<T,DESCRIPTOR>() / Sc[1];
277 T tauTurbADPrefactor2 = descriptors::invCs2<T,DESCRIPTOR_ADE>() / descriptors::invCs2<T,DESCRIPTOR>() / Sc[2];
278 T tau_turb_AD0 = tau_turb_NS * tauTurbADPrefactor0;
279 T tau_turb_AD1 = tau_turb_NS * tauTurbADPrefactor1;
280 T tau_turb_AD2 = tau_turb_NS * tauTurbADPrefactor2;
281 cells.template get<names::Concentration0>().template setField<descriptors::OMEGA>(T{1} / (tau_mol_AD0 + tau_turb_AD0));
282 cells.template get<names::Concentration1>().template setField<descriptors::OMEGA>(T{1} / (tau_mol_AD1 + tau_turb_AD1));
283 cells.template get<names::Concentration2>().template setField<descriptors::OMEGA>(T{1} / (tau_mol_AD2 + tau_turb_AD2));
284
285 auto stochCoeff =
parameters.template get<STOCH_COEFF>();
286 auto reactionOrder =
parameters.template get<REACTION_ORDER>();
287 T source =
parameters.template get<LATTICE_REACTION_COEFF>();
288 {
289 T conc = cells.template get<names::Concentration0>().computeRho();
290 source *=
util::pow(conc, reactionOrder[0]);
291 }
292 {
293 T conc = cells.template get<names::Concentration1>().computeRho();
294 source *=
util::pow(conc, reactionOrder[1]);
295 }
296 {
297 T conc = cells.template get<names::Concentration2>().computeRho();
298 source *=
util::pow(conc, reactionOrder[2]);
299 }
300 {
301 cells.template get<names::Concentration0>().template setField<descriptors::SOURCE>(stochCoeff[0]*source);
302 cells.template get<names::Concentration1>().template setField<descriptors::SOURCE>(stochCoeff[1]*source);
303 cells.template get<names::Concentration2>().template setField<descriptors::SOURCE>(stochCoeff[2]*source);
304 }
305 }
cpu::simd::Pack< T > sqrt(cpu::simd::Pack< T > value)
cpu::simd::Pack< T > pow(cpu::simd::Pack< T > base, cpu::simd::Pack< T > exp)
meta::list< LATTICE_REACTION_COEFF, STOCH_COEFF, REACTION_ORDER, SMAGORINSKY_PREFACTOR, SCHMIDT, OMEGA_NSE, OMEGAS_ADE > parameters
static constexpr int n
result stored in n