25#ifndef FREE_ENERGY_DYNAMICS_H
26#define FREE_ENERGY_DYNAMICS_H
48 template <
typename DESCRIPTOR,
typename MOMENTA>
52 template <
typename RHO,
typename U,
typename V=RHO>
56 eq += (V{1} - descriptors::t<V,DESCRIPTOR>(0)) * rho;
58 eq -= descriptors::t<V,DESCRIPTOR>(iPop) * rho;
63 template <
typename CELL,
typename PARAMETERS,
typename FEQ,
typename V=
typename CELL::value_t>
65 auto u = cell.template getField<descriptors::FORCE>();
67 const V uSqr = util::normSqr<V,DESCRIPTOR::d>(u);
68 for (
int iPop=0; iPop < DESCRIPTOR::q; ++iPop) {
89 template <
typename DESCRIPTOR,
typename MOMENTA,
typename EQUILIBRIUM>
93 template <
typename CELL,
typename PARAMETERS,
typename V=
typename CELL::value_t>
95 V fEq[DESCRIPTOR::q] { };
97 const V omega =
parameters.template get<descriptors::OMEGA>();
98 const V gamma =
parameters.template get<GAMMA>();
99 for (
int iPop=0; iPop < DESCRIPTOR::q; ++iPop) {
100 cell[iPop] *= V{1} - omega;
101 cell[iPop] += omega * fEq[iPop];
103 const V tmp = gamma * descriptors::invCs2<V,DESCRIPTOR>()
104 * cell.template getField<descriptors::CHEM_POTENTIAL>();
105 for (
int iPop=1; iPop < DESCRIPTOR::q; ++iPop) {
106 cell[iPop] -= omega * descriptors::t<V,DESCRIPTOR>(iPop) * (statistic.rho - tmp);
108 cell[0] += omega * (V{1} - descriptors::t<V,DESCRIPTOR>(0)) * (statistic.rho - tmp);
114template <
int direction,
int orientation>
119 return "FreeEnergyInletOutlet<"
120 + std::to_string(direction) +
"," + std::to_string(orientation) +
124 template <
typename DESCRIPTOR,
typename MOMENTA,
typename EQUILIBRIUM>
128 template <
typename CELL,
typename PARAMETERS,
typename V=
typename CELL::value_t>
131 V rho, u[DESCRIPTOR::d];
132 MomentaF().computeRhoU(cell, rho, u);
133 const V omega =
parameters.template get<descriptors::OMEGA>();
135 for (
int iPop=1; iPop < DESCRIPTOR::q; ++iPop) {
136 cell[iPop] -= omega * descriptors::t<V,DESCRIPTOR>(iPop) * rho;
138 cell[0] += omega * (V{1} - descriptors::t<V,DESCRIPTOR>(0)) * rho;
140 constexpr auto missingIndices = util::subIndexOutgoing<DESCRIPTOR,direction,orientation>();
141 V missingRho = rho - V{1};
142 V missingWeightSum = 0;
143 for (
int iPop=0; iPop < DESCRIPTOR::q; ++iPop) {
152 bool contains =
false;
153 for(
long unsigned int i=0; i < missingIndices.size(); i++){
154 if(missingIndices[i] == (
long unsigned int)iPop){
160 missingWeightSum += descriptors::t<V,DESCRIPTOR>(iPop);
162 missingRho -= cell[iPop];
167 for (
unsigned iPop=0; iPop < missingIndices.size(); ++iPop) {
168 cell[missingIndices[iPop]] = missingRho * descriptors::t<V,DESCRIPTOR>(missingIndices[iPop]) / missingWeightSum;
177template <
typename T,
typename DESCRIPTOR,
typename MOMENTA=momenta::FreeEnergyBulkTuple>
185template <
typename T,
typename DESCRIPTOR>
198template <
typename T,
typename DESCRIPTOR,
int direction,
int orientation>
Instantiation of Momenta tuples which define the computation and definition of momenta in a specific ...
Top level namespace for all of OpenLB.
Return value of any collision.
typename MOMENTA::template type< DESCRIPTOR > MomentaF
CellStatistic< V > apply(CELL &cell, PARAMETERS ¶meters) any_platform
static std::string getName()
typename meta::list< descriptors::OMEGA > parameters
CellStatistic< V > apply(CELL &cell, PARAMETERS ¶meters) any_platform
typename EQUILIBRIUM::template type< DESCRIPTOR, MOMENTA > EquilibriumF
static std::string getName()
typename meta::list< descriptors::OMEGA, GAMMA > parameters
Base of a field whose size is defined by [C,U_1,...,U_N]^T * [1,V_1,...V_N].
Dynamics constructed as a tuple of momenta, equilibrium and collision.
auto compute(int iPop, const RHO &rho, const U &u) any_platform
CellStatistic< V > compute(CELL &cell, PARAMETERS ¶meters, FEQ &fEq) any_platform
typename MOMENTA::template type< DESCRIPTOR > MomentaF
static std::string getName()
static V secondOrder(int iPop, const RHO &rho, const U &u, const USQR &uSqr) any_platform
Computation of equilibrium distribution, second order in u.
static V bgkCollision(CELL &cell, const RHO &rho, const VELOCITY &u, const OMEGA &omega) any_platform
BGK collision step.
Standard computation for density in the bulk as zeroth moment of the population.
The momenta are defined one after the other.
The density is stored in descriptors::FORCE[0] (TODO: absurd, to be changed)
Computation of the stress tensor for regularized boundary nodes.
Momentum is zero at solid material.
The stress is always zero.