I would like to modify damBreak3d.cpp in such a way that there is a vertical gate that slowly shifts upwards and releases the water gradually.
I thought that I could call
setToFunction(fields.flag …
once every 100 steps and modify flag field to match the shifted state of the gate. This however does not seem to work. Any idea how to achieve it?
Hello czb,
you could for example have a wall of bounceback nodes and every 100 or whatever timesteps change the dynamics of the lowest row of the wall to whatever dynamics object you are using via
[code=“cpp”]
defineDynamics(…)
This is just a wild guess, I am not really familiar with the free surface model. Anyway, I hope i could help.
best
PhilippeI have written the following code to do the work (although not for the DamBreak example but for lifting a cone instead):
[code=“cpp”]
#include “core/globalDefs.h”
#include “core/block3D.h”
#include “atomicBlock/dataProcessor3D.h”
#include “atomicBlock/dataProcessingFunctional3D.h”
#include “atomicBlock/blockLattice3D.h”
#include “atomicBlock/atomicContainerBlock3D.h”
namespace plb {
template< typename T,template class Descriptor>
class ConeLifter3D : public BoxProcessingFunctional3D {
public:
ConeLifter3D() { }
virtual void processGenericBlocks(Box3D domain, std::vector<AtomicBlock3D*> atomicBlocks) {
using namespace twoPhaseFlag;
FreeSurfaceProcessorParam3D<T,Descriptor> param(atomicBlocks);
typedef Descriptor D;
std::vector toGas, toWall;
for (plint iX=domain.x0; iX<=domain.x1; ++iX) {
for (plint iY=domain.y0; iY<=domain.y1; ++iY) {
for (plint iZ=domain.z0+1; iZ<=domain.z1; ++iZ) {
if (param.flag(iX,iY,iZ) == wall && (iZ < 3 || param.flag(iX,iY,iZ-1) != wall)) {
toGas.push_back(iX);
toGas.push_back(iY);
toGas.push_back(iZ);
} else if(iZ > 2 && param.flag(iX,iY,iZ) == empty && param.flag(iX,iY,iZ-1) == wall) {
toWall.push_back(iX);
toWall.push_back(iY);
toWall.push_back(iZ);
}
}
}
}
for (std::vector<plint>::size_type i = 0; i != toGas.size(); i=i+3) {
param.flag(toGas[i], toGas[i+1], toGas[i+2]) = empty;
param.attributeDynamics(toGas[i], toGas[i+1], toGas[i+2], new NoDynamics<T, Descriptor>);
}
for (std::vector<plint>::size_type i = 0; i != toWall.size(); i=i+3) {
param.flag(toWall[i], toWall[i+1], toWall[i+2]) = wall;
param.attributeDynamics(toWall[i], toWall[i+1], toWall[i+2], new BounceBack<T, Descriptor>);
}
// In the following, spot the interface cells and tag them.
for (plint iX=domain.x0; iX<=domain.x1; ++iX) {
for (plint iY=domain.y0; iY<=domain.y1; ++iY) {
for (plint iZ=domain.z0+1; iZ<=domain.z1; ++iZ) {
if (isEmpty(param.flag(iX,iY,iZ))) {
for (plint iPop=1; iPop<D::q; iPop++) {
plint nextX = iX + D::c[iPop][0];
plint nextY = iY + D::c[iPop][1];
plint nextZ = iZ + D::c[iPop][2];
if (isFullWet(param.flag(nextX,nextY,nextZ))) {
param.flag(iX,iY,iZ) = interface;
T rho = param.getDensity(nextX, nextY, nextZ);
Array<T,3> j = param.getMomentum(nextX, nextY, nextZ);
iniCellAtEquilibrium(param.cell(iX,iY,iZ), rho, j/rho);
param.setDensity(iX,iY,iZ, rho);
param.setMomentum(iX,iY,iZ, j);
param.setForce(iX,iY,iZ, param.getForce(nextX, nextY, nextZ));
param.attributeDynamics(iX,iY,iZ, param.cell(nextX, nextY, nextZ).getDynamics().clone());
param.mass(iX,iY,iZ) = 0.01 * rho;
param.volumeFraction(iX,iY,iZ) = (T)0.01;
break;
}
}
}
}
}
}
}
virtual ConeLifter3D<T,Descriptor>* clone() const {
return new ConeLifter3D(*this);
}
virtual void getTypeOfModification(std::vector<modif::ModifT>& modified) const {
std::fill(modified.begin(), modified.end(), modif::nothing);
modified[0] = modif::dataStructure; // Fluid
modified[1] = modif::staticVariables; // rhoBar.
modified[2] = modif::staticVariables; // j.
modified[3] = modif::staticVariables; // Mass
modified[4] = modif::staticVariables; // Volume-fraction
modified[5] = modif::staticVariables; // Flag-status
modified[6] = modif::nothing; // Normal.
modified[7] = modif::nothing; // Interface lists
modified[8] = modif::nothing; // Curvature.
modified[9] = modif::nothing; // Outside density.
}
};
template< typename T,template class Descriptor>
void liftCone3D(FreeSurfaceFields3D<T,Descriptor>& fields) {
executeDataProcessor(
BoxProcessorGenerator3D(
new ConeLifter3D<T,Descriptor>(), fields.lattice.getBoundingBox()
), fields.twoPhaseArgs );
}
} // namespace plb