External force on specific cells

[luka]

Dear community,

I’m new to LBM method as well ass cpp but i have some experience in other programming and CFD (mostly FEM) methods.

I have to create simulation that will be very simplified version of ionic wind (electrically driven flow). My approach for now is to define area between my “electrodes” – two cylinders and apply external force to them. I’m trying for now only to avoid transport of ions, as well as transients at startup of simulation. For now i have box in which force should be applied. I started setup of my case from laminar cylinder2d.

Other approach I was thinking of was to add material and than apply force to it but definition of material that is correctly representing parts of domain in which force is applied doesn’t seem feasible for now.

Full code is available here: https://drive.google.com/file/d/18q9-TVMC7OFMiknYUoT6MrG5EYQfborg/view?usp=sharing

My function of interest for now is this:
void applyExternalForce(olb::SuperLattice<double, olb::descriptors::D2Q9<olb::descriptors::FORCE>>& sLattice,
const olb::UnitConverter<double, olb::descriptors::D2Q9<olb::descriptors::FORCE>>& converter,
olb::SuperGeometry<double, 2>& superGeometry) {
// Parameters for external force
const double Fx = 0.0001; // Constant force magnitude in x-direction
const double Fy = 0.0000; // Constant force magnitude in y-direction

// Define region of interest
Vector<double, 2> minPosition(centerCylinderX1 + radiusCylinder, centerCylinderY1 – radiusCylinder);
Vector<double, 2> maxPosition(centerCylinderX2 – radiusCylinder, centerCylinderY2 + radiusCylinder);

// Create a force vector
olb::Vector<double, 2> forceVector({0.0000, 0.0000});
olb::AnalyticalConst2D<double, double> force(forceVector);

// Iterate over all lattice cells
//for (std::size_t iC = 0; iC < sLattice.getNblock(); ++iC) {
for (std::size_t iC = 0; iC < 1; ++iC) {
auto& block = sLattice.getBlock(iC);
for (int iX = 0; iX < block.getNx(); ++iX) {
for (int iY = 0; iY < block.getNy(); ++iY) {
// Retrieve the cell
auto cell = sLattice.get(static_cast<int>(iX), static_cast<int>(iY), static_cast<int>(iC));
auto force = cell.template getFieldPointer<descriptors::FORCE>();

olb::Vector<double, 2> position(iX, iY);

std::cout << “cell: ” << iC << std::endl;

std::cout << “Position: ” << position << std::endl;
std::cout << “ForceVector before loop: ” << forceVector << std::endl;

if (position[0] > minPosition[0] && position[0] < maxPosition[0] &&
position[1] > minPosition[1] && position[1] < maxPosition[1]) {

for (int iPop = 0; iPop < olb::descriptors::D2Q9<olb::descriptors::FORCE>::q; ++iPop) {
forceVector[0] += force[iPop] + Fx;
forceVector[1] += force[iPop] + Fy;
}

std::cout << “ForceVector after loop: ” << forceVector << std::endl;

cell.template setField<descriptors::FORCE>(forceVector);
}

// Add external force to the cell
//olb::lbm<olb::descriptors::D2Q9<olb::descriptors::FORCE>>::defineField(cell, rho, u, omega, forceVector);
}
}
}
}

[FBukreev]

Hello Luka,

you can try to use gpu appropriate postprocessor:
struct PostProcessor {
static constexpr OperatorScope scope = OperatorScope::PerCellWithParameters;

using parameters = meta::list;

int getPriority() const {
return 0;
}

template
void apply(CELL& cell, PARAMETERS& params) any_platform {
using V = typename CELL::value_t;
std::size_t iT = params.template get();
….
V force = …;
cell.template setField(force);
}
};

then in main class:
sLattice.addPostProcessor(meta::id{});

and in the simulation for loop:
sLattice.executePostProcessors(stage::PreCoupling());

Please do not forget to convert force from SI units to lattice units, whereby in lattice units in openlb force must be also divided by density.

Best wishes
Fedor

[Author]

Posts