olb::MixedScaleBoussinesqCouplingPostProcessor3D< T, DESCRIPTOR > Class Template Reference

#include <navierStokesAdvectionDiffusionCouplingPostProcessor3D.h>

Inheritance diagram for olb::MixedScaleBoussinesqCouplingPostProcessor3D< T, DESCRIPTOR >:

Collaboration diagram for olb::MixedScaleBoussinesqCouplingPostProcessor3D< T, DESCRIPTOR >:

Public Member Functions
	MixedScaleBoussinesqCouplingPostProcessor3D (int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T gravity_, T T0_, T deltaTemp_, std::vector< T > dir_, T PrTurb_, std::vector< BlockStructureD< 3 > * > partners_)
int	extent () const override
	Extent of application area (0 for purely local operations)
int	extent (int whichDirection) const override
	Extent of application area along a direction (0 or 1)
void	process (BlockLattice< T, DESCRIPTOR > &blockLattice) override
	Execute post-processing step.
void	processSubDomain (BlockLattice< T, DESCRIPTOR > &blockLattice, int x0_, int x1_, int y0_, int y1_, int z0_, int z1_) override
	Execute post-processing step on a sublattice.
Public Member Functions inherited from olb::PostProcessor3D< T, DESCRIPTOR >
	PostProcessor3D ()
virtual	~PostProcessor3D ()
std::string &	getName ()
	read and write access to name
std::string const &	getName () const
	read only access to name
int	getPriority () const
	read only access to priority

Additional Inherited Members
Protected Attributes inherited from olb::PostProcessor3D< T, DESCRIPTOR >
int	_priority

Detailed Description

template<typename T, typename DESCRIPTOR>
class olb::MixedScaleBoussinesqCouplingPostProcessor3D< T, DESCRIPTOR >

Definition at line 304 of file navierStokesAdvectionDiffusionCouplingPostProcessor3D.h.

Constructor & Destructor Documentation

MixedScaleBoussinesqCouplingPostProcessor3D()

template<typename T , typename DESCRIPTOR >

olb::MixedScaleBoussinesqCouplingPostProcessor3D< T, DESCRIPTOR >::MixedScaleBoussinesqCouplingPostProcessor3D	(	int	x0_,
		int	x1_,
		int	y0_,
		int	y1_,
		int	z0_,
		int	z1_,
		T	gravity_,
		T	T0_,
		T	deltaTemp_,
		std::vector< T >	dir_,
		T	PrTurb_,
		std::vector< BlockStructureD< 3 > * >	partners_ )

Definition at line 692 of file navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh.

  :  x0(x0_), x1(x1_), y0(y0_), y1(y1_), z0(z0_), z1(z1_),
     gravity(gravity_), T0(T0_), deltaTemp(deltaTemp_),
     dir(dir_), PrTurb(PrTurb_), partners(partners_)
{
  // we normalize the direction of force vector
  T normDir = T();
  for (unsigned iD = 0; iD < dir.size(); ++iD) {
    normDir += dir[iD]*dir[iD];
  }
  normDir = util::sqrt(normDir);
  for (unsigned iD = 0; iD < dir.size(); ++iD) {
    dir[iD] /= normDir;
  }
 
  for (unsigned iD = 0; iD < dir.size(); ++iD) {
    forcePrefactor[iD] = gravity * dir[iD];
  }
 
  tauTurbADPrefactor = descriptors::invCs2<T,descriptors::D3Q7<>>() / descriptors::invCs2<T,DESCRIPTOR>() / PrTurb;
  tPartner = static_cast<BlockLattice<T,descriptors::D3Q7<descriptors::VELOCITY,descriptors::TAU_EFF,descriptors::CUTOFF_HEAT_FLUX>> *>(partners[0]);
}

References olb::util::sqrt().

Here is the call graph for this function:

Member Function Documentation

extent() [1/2]

template<typename T , typename DESCRIPTOR >

int olb::MixedScaleBoussinesqCouplingPostProcessor3D< T, DESCRIPTOR >::extent ( ) const

inlineoverridevirtual

Extent of application area (0 for purely local operations)

Implements olb::PostProcessor3D< T, DESCRIPTOR >.

Definition at line 309 of file navierStokesAdvectionDiffusionCouplingPostProcessor3D.h.

  {
    return 0;
  }

extent() [2/2]

template<typename T , typename DESCRIPTOR >

int olb::MixedScaleBoussinesqCouplingPostProcessor3D< T, DESCRIPTOR >::extent ( int direction ) const

inlineoverridevirtual

Extent of application area along a direction (0 or 1)

Implements olb::PostProcessor3D< T, DESCRIPTOR >.

Definition at line 313 of file navierStokesAdvectionDiffusionCouplingPostProcessor3D.h.

  {
    return 0;
  }

process()

template<typename T , typename DESCRIPTOR >

void olb::MixedScaleBoussinesqCouplingPostProcessor3D< T, DESCRIPTOR >::process ( BlockLattice< T, DESCRIPTOR > & blockLattice )

overridevirtual

Execute post-processing step.

Implements olb::PostProcessor3D< T, DESCRIPTOR >.

Definition at line 962 of file navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh.

{
  processSubDomain(blockLattice, x0, x1, y0, y1, z0, z1);
}

processSubDomain()

template<typename T , typename DESCRIPTOR >

void olb::MixedScaleBoussinesqCouplingPostProcessor3D< T, DESCRIPTOR >::processSubDomain ( BlockLattice< T, DESCRIPTOR > & blockLattice, int x0_, int x1_, int y0_, int y1_, int z0_, int z1_ )

overridevirtual

Execute post-processing step on a sublattice.

Molecular realaxation time

Effective realaxation time

Implements olb::PostProcessor3D< T, DESCRIPTOR >.

Definition at line 719 of file navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh.

{
 
  const T C_nu = 0.04;
  const T C_alpha = 0.5;
  // const T deltaT = 1.0;
 
  const T invCs2 = descriptors::invCs2<T,DESCRIPTOR>();
  const T invCs2_g = descriptors::invCs2<T,descriptors::D3Q7<>>();
 
  int newX0, newX1, newY0, newY1, newZ1, newZ0;
  if ( util::intersect (
         x0, x1, y0, y1, z0, z1,
         x0_, x1_, y0_, y1_, z0_, z1_,
         newX0, newX1, newY0, newY1, newZ1, newZ0 ) ) {
 
    auto& heatFluxCache = blockLattice.template getField<HEAT_FLUX_CACHE>()[0];
    for (int iX=newX0-1; iX<=newX1+1; ++iX) {
      for (int iY=newY0-1; iY<=newY1+1; ++iY) {
        for (int iZ=newZ0-1; iZ<=newZ1+1; ++iZ) {
          const T temperature = tPartner->get(iX,iY,iZ).computeRho();
          heatFluxCache[blockLattice.getCellId(iX, iY, iZ)] = temperature;
 
          // computation of the bousinessq force
          T temperatureDifference = temperature - T0;
          blockLattice.get(iX,iY,iZ).template setField<descriptors::FORCE>(temperatureDifference*forcePrefactor);
 
          FieldD<T,DESCRIPTOR,descriptors::VELOCITY> u;
          blockLattice.get(iX,iY,iZ).computeU(u.data());
          tPartner->get(iX,iY,iZ).template setField<descriptors::VELOCITY>(u);
        }
      }
    }
 
    for (int iX=newX0; iX<=newX1; ++iX) {
      for (int iY=newY0; iY<=newY1; ++iY) {
        for (int iZ=newZ0; iZ<=newZ1; ++iZ) {
 
          auto u_ppp = tPartner->get(iX+1, iY+1, iZ+1).template getField<descriptors::VELOCITY>();
          auto u_0pp = tPartner->get(iX, iY+1, iZ+1).template getField<descriptors::VELOCITY>();
          auto u_npp = tPartner->get(iX-1, iY+1, iZ+1).template getField<descriptors::VELOCITY>();
          auto u_p0p = tPartner->get(iX+1, iY, iZ+1).template getField<descriptors::VELOCITY>();
          auto u_00p = tPartner->get(iX, iY, iZ+1).template getField<descriptors::VELOCITY>();
          auto u_n0p = tPartner->get(iX-1, iY, iZ+1).template getField<descriptors::VELOCITY>();
          auto u_pnp = tPartner->get(iX+1, iY-1, iZ+1).template getField<descriptors::VELOCITY>();
          auto u_0np = tPartner->get(iX,   iY-1, iZ+1).template getField<descriptors::VELOCITY>();
          auto u_nnp = tPartner->get(iX-1, iY-1, iZ+1).template getField<descriptors::VELOCITY>();
 
          auto u_pp0 = tPartner->get(iX+1, iY+1, iZ  ).template getField<descriptors::VELOCITY>();
          auto u_0p0 = tPartner->get(iX, iY+1, iZ  ).template getField<descriptors::VELOCITY>();
          auto u_np0 = tPartner->get(iX-1, iY+1, iZ  ).template getField<descriptors::VELOCITY>();
          auto u_p00 = tPartner->get(iX+1, iY, iZ  ).template getField<descriptors::VELOCITY>();
          auto u_000 = tPartner->get(iX, iY, iZ  ).template getField<descriptors::VELOCITY>();
          auto u_n00 = tPartner->get(iX-1, iY, iZ  ).template getField<descriptors::VELOCITY>();
          auto u_pn0 = tPartner->get(iX+1, iY-1, iZ  ).template getField<descriptors::VELOCITY>();
          auto u_0n0 = tPartner->get(iX,   iY-1, iZ  ).template getField<descriptors::VELOCITY>();
          auto u_nn0 = tPartner->get(iX-1, iY-1, iZ  ).template getField<descriptors::VELOCITY>();
 
          auto u_ppn = tPartner->get(iX+1, iY+1, iZ-1).template getField<descriptors::VELOCITY>();
          auto u_0pn = tPartner->get(iX, iY+1, iZ-1).template getField<descriptors::VELOCITY>();
          auto u_npn = tPartner->get(iX-1, iY+1, iZ-1).template getField<descriptors::VELOCITY>();
          auto u_p0n = tPartner->get(iX+1, iY, iZ-1).template getField<descriptors::VELOCITY>();
          auto u_00n = tPartner->get(iX, iY, iZ-1).template getField<descriptors::VELOCITY>();
          auto u_n0n = tPartner->get(iX-1, iY, iZ-1).template getField<descriptors::VELOCITY>();
          auto u_pnn = tPartner->get(iX+1, iY-1, iZ-1).template getField<descriptors::VELOCITY>();
          auto u_0nn = tPartner->get(iX,   iY-1, iZ-1).template getField<descriptors::VELOCITY>();
          auto u_nnn = tPartner->get(iX-1, iY-1, iZ-1).template getField<descriptors::VELOCITY>();
 
          const T *h_ppp = & heatFluxCache[blockLattice.getCellId(iX+1, iY+1, iZ+1)];
          const T *h_0pp = & heatFluxCache[blockLattice.getCellId(iX, iY+1, iZ+1)];
          const T *h_npp = & heatFluxCache[blockLattice.getCellId(iX-1, iY+1, iZ+1)];
          const T *h_p0p = & heatFluxCache[blockLattice.getCellId(iX+1, iY, iZ+1)];
          const T *h_pnp = & heatFluxCache[blockLattice.getCellId(iX+1, iY-1, iZ+1)];
          const T *h_00p = & heatFluxCache[blockLattice.getCellId(iX, iY, iZ+1)];
          const T *h_0np = & heatFluxCache[blockLattice.getCellId(iX, iY-1, iZ+1)];
          const T *h_n0p = & heatFluxCache[blockLattice.getCellId(iX-1, iY, iZ+1)];
          const T *h_nnp = & heatFluxCache[blockLattice.getCellId(iX-1, iY-1, iZ+1)];
 
          const T *h_pp0 = & heatFluxCache[blockLattice.getCellId(iX+1, iY+1, iZ  )];
          const T *h_0p0 = & heatFluxCache[blockLattice.getCellId(iX, iY+1, iZ  )];
          const T *h_np0 = & heatFluxCache[blockLattice.getCellId(iX-1, iY+1, iZ  )];
          const T *h_p00 = & heatFluxCache[blockLattice.getCellId(iX+1, iY, iZ  )];
          const T *h_pn0 = & heatFluxCache[blockLattice.getCellId(iX+1, iY-1, iZ  )];
          const T *h_000 = & heatFluxCache[blockLattice.getCellId(iX, iY, iZ  )];
          const T *h_0n0 = & heatFluxCache[blockLattice.getCellId(iX, iY-1, iZ  )];
          const T *h_n00 = & heatFluxCache[blockLattice.getCellId(iX-1, iY, iZ  )];
          const T *h_nn0 = & heatFluxCache[blockLattice.getCellId(iX-1, iY-1, iZ  )];
 
          const T *h_ppn = & heatFluxCache[blockLattice.getCellId(iX+1, iY+1, iZ-1)];
          const T *h_0pn = & heatFluxCache[blockLattice.getCellId(iX, iY+1, iZ-1)];
          const T *h_npn = & heatFluxCache[blockLattice.getCellId(iX-1, iY+1, iZ-1)];
          const T *h_p0n = & heatFluxCache[blockLattice.getCellId(iX+1, iY, iZ-1)];
          const T *h_pnn = & heatFluxCache[blockLattice.getCellId(iX+1, iY-1, iZ-1)];
          const T *h_00n = & heatFluxCache[blockLattice.getCellId(iX, iY, iZ-1)];
          const T *h_0nn = & heatFluxCache[blockLattice.getCellId(iX, iY-1, iZ-1)];
          const T *h_n0n = & heatFluxCache[blockLattice.getCellId(iX-1, iY, iZ-1)];
          const T *h_nnn = & heatFluxCache[blockLattice.getCellId(iX-1, iY-1, iZ-1)];
 
//        cout<<"h_ppn= "<<h_ppp[0]<<endl;
//        cout<<"h_0pn= "<<h_0pp[0]<<endl;
//        cout<<"h_npn= "<<h_npp[0]<<endl;
//        cout<<"h_p0n= "<<h_p0p[0]<<endl;
//        cout<<"h_00n= "<<h_00p[0]<<endl;
//        cout<<"h_n0n= "<<h_n0p[0]<<endl;
//        cout<<"h_pnn= "<<h_pnp[0]<<endl;
//        cout<<"h_0nn= "<<h_0np[0]<<endl;
//        cout<<"h_nnn= "<<h_nnp[0]<<endl;
 
          //Testfilter h 3D
          Vector<T,3> filtered_u;
          T filtered_h = (
                           (
                             ( h_ppp[0] + 2. * h_0pp[0] + h_npp[0]) + 2.*
                             ( h_p0p[0] + 2. * h_00p[0] + h_n0p[0]) +
                             ( h_pnp[0] + 2. * h_0np[0] + h_nnp[0])
                           ) * 0.25 * 0.25 + 2. *
 
                           (
                             ( h_pp0[0] + 2. * h_0p0[0] + h_np0[0]) + 2.*
                             ( h_p00[0] + 2. * h_000[0] + h_n00[0]) +
                             ( h_pn0[0] + 2. * h_0n0[0] + h_nn0[0])
                           ) * 0.25 * 0.25  +
 
                           (
                             ( h_ppn[0] + 2. * h_0pn[0] + h_npn[0]) + 2.*
                             ( h_p0n[0] + 2. * h_00n[0] + h_n0n[0]) +
                             ( h_pnn[0] + 2. * h_0nn[0] + h_nnn[0])
                           ) * 0.25 * 0.25
                         ) * 0.25;
 
          //Testfilter u 3D
          filtered_u = (
                         (
                           ( u_ppp + 2. * u_0pp + u_npp) + 2.*
                           ( u_p0p + 2. * u_00p + u_n0p) +
                           ( u_pnp + 2. * u_0np + u_nnp)
                         ) * 0.25 * 0.25 + 2. *
 
                         (
                           ( u_pp0 + 2. * u_0p0 + u_np0) + 2.*
                           ( u_p00 + 2. * u_000 + u_n00) +
                           ( u_pn0 + 2. * u_0n0 + u_nn0)
                         ) * 0.25 * 0.25 +
 
                         (
                           ( u_ppn + 2. * u_0pn + u_npn) + 2.*
                           ( u_p0n + 2. * u_00n + u_n0n) +
                           ( u_pnn + 2. * u_0nn + u_nnn)
                         ) * 0.25 * 0.25
                       ) * 0.25;
 
          //Ende Testfilter u
          // filtered_u[i] -= u_00[i];
 
//        cout << "u["<<i<<"]= "<< filtered_u[i] << std::endl;
//        cout << "h[0]= "<< filtered_h[0] << std::endl;
          const Vector<T,3> filtered_u_reduced = u_000 - filtered_u;
          const Vector<T,3> filtered_heatFlux_reduced = u_000*h_000[0] - filtered_u*filtered_h;
          T cutoffKinEnergy = filtered_u_reduced[0] * filtered_u_reduced[0]
                              + filtered_u_reduced[1] * filtered_u_reduced[1]
                              + filtered_u_reduced[2] * filtered_u_reduced[2];
          T cutoffHeatFlux = filtered_heatFlux_reduced[0] * filtered_heatFlux_reduced[0]
                             + filtered_heatFlux_reduced[1] * filtered_heatFlux_reduced[1]
                             + filtered_heatFlux_reduced[2] * filtered_heatFlux_reduced[2];
//          cout << "filtered_u_reduced= " << filtered_u_reduced << std::endl;
//          cout << "filtered_heatFlux_reduced= " << filtered_heatFlux_reduced << std::endl;
//          cout << "cutoffKinEnergy= " << cutoffKinEnergy << std::endl;
//          cout << "cutoffHeatFlux= " << cutoffHeatFlux << std::endl;
 
          blockLattice.get(iX,iY,iZ).template setField<descriptors::CUTOFF_KIN_ENERGY>(util::pow(0.5*cutoffKinEnergy, 0.25));
          tPartner->get(iX,iY,iZ).template setField<descriptors::CUTOFF_HEAT_FLUX>(util::pow(0.5*cutoffHeatFlux, 0.25));
 
//        std::cout<<"cutoffKinEnergy_14= "<<cutoffKinEnergy_14[0]<<std::endl;
//          std::cout<<"cutoffHeatFlux_14= "<<cutoffHeatFlux_14[0]<<std::endl;
 
 
          // tau coupling
          auto tauNS = blockLattice.get(iX,iY,iZ).template getField<descriptors::TAU_EFF>();
          auto tauAD = tPartner->get(iX,iY,iZ).template getField<descriptors::TAU_EFF>();
 
 
//        std::cout<<"tau_NS= "<<tauNS[0]<<std::endl;
//        std::cout<<"tau_AD= "<<tauAD[0]<<std::endl;
 
          T tau_mol_NS = 1. / blockLattice.get(iX, iY, iZ).getDynamics()->getParameters(blockLattice).template getOrFallback<descriptors::OMEGA>(0);
          T tau_mol_AD = 1. / tPartner->get(iX, iY, iZ).getDynamics()->getParameters(*tPartner).template getOrFallback<descriptors::OMEGA>(0);
 
          const T temperature = h_000[0];
 
          // computation of the bousinessq force
          auto force = blockLattice.get(iX,iY,iZ).template getField<descriptors::FORCE>();
          // T temperatureDifference = temperature - T0;
          // for (unsigned iD = 0; iD < L::d; ++iD) {
          //   force[iD] = forcePrefactor[iD] * temperatureDifference;
          // }
 
          auto u = tPartner->get(iX,iY,iZ).template getFieldPointer<descriptors::VELOCITY>();
          T rho, pi[util::TensorVal<DESCRIPTOR>::n], j[DESCRIPTOR::d];
          // blockLattice.get(iX,iY).computeAllMomenta(rho, u, pi);
          rho = blockLattice.get(iX,iY,iZ).computeRho();
          blockLattice.get(iX,iY,iZ).computeStress(pi);
 
          int iPi = 0;
          const int dim = DESCRIPTOR::d;
          for (int Alpha=0; Alpha<dim; ++Alpha) {
            for (int Beta=Alpha; Beta<dim; ++Beta) {
              pi[iPi] += rho/2.*(force[Alpha]*u[Beta] + u[Alpha]*force[Beta]);
              ++iPi;
            }
          }
          const T piSqr[6] = {pi[0]*pi[0], pi[1]*pi[1], pi[2]*pi[2], pi[3]*pi[3], pi[4]*pi[4], pi[5]*pi[5]};
          const T PiNeqNormSqr = piSqr[0] + 2.0 * piSqr[1] + 2.0 * piSqr[2] + piSqr[3] + 2.0 * piSqr[4] + piSqr[5];
          const T PiNeqNorm    = util::sqrt(PiNeqNormSqr);
 
          tPartner->get(iX, iY, iZ).computeJ(j);
 
          const Vector<T,3> jNeq = {(j[0] - temperature * u[0]), (j[1] - temperature * u[1]), (j[2] - temperature * u[2])};
 
          const Vector<T,6> t = { 2.0 * jNeq[0] * pi[0], (jNeq[0] + jNeq[1]) * pi[1], (jNeq[0] + jNeq[2]) * pi[2], 2.0 * jNeq[1] * pi[3], (jNeq[1] + jNeq[2]) * pi[4], 2.0 * jNeq[2] * pi[5]};
          const Vector<T,6> tSqr = {t[0]*t[0], t[1]*t[1], t[2]*t[2], t[3]*t[3], t[4]*t[4], t[5]*t[5]};
          const T TnormSqr = tSqr[0] + 2.0 * tSqr[1] + 2.0 * tSqr[2] + tSqr[3] + 2.0 * tSqr[4] + tSqr[5];
          const T Tnorm    = util::sqrt(2.0 * 0.25 * TnormSqr);
 
          auto cutoffKinEnergy_14 = blockLattice.get(iX,iY,iZ).template getField<descriptors::CUTOFF_KIN_ENERGY>();
          auto cutoffHeatFlux_14 = tPartner->get(iX,iY,iZ).template getField<descriptors::CUTOFF_HEAT_FLUX>();
 
          const T tau_turb_NS = C_nu * util::sqrt(util::sqrt(2.)/2.) * invCs2 * invCs2 * util::sqrt(PiNeqNorm / rho / tauNS) * cutoffKinEnergy_14;
          const T tau_turb_AD = C_alpha * invCs2 / rho * util::sqrt(Tnorm * invCs2_g / tauNS / tauAD) * cutoffHeatFlux_14;
 
          blockLattice.get(iX,iY,iZ).template setField<descriptors::TAU_EFF>(tau_mol_NS + tau_turb_NS);
          tPartner->get(iX,iY,iZ).template setField<descriptors::TAU_EFF>(tau_mol_AD + tau_turb_AD);
 
        }
      }
    }
  }
}

References olb::Vector< T, D >::data(), olb::BlockLattice< T, DESCRIPTOR >::get(), olb::BlockStructureD< D >::getCellId(), olb::util::intersect(), olb::util::pow(), and olb::util::sqrt().

Here is the call graph for this function:

---

The documentation for this class was generated from the following files:

• src/dynamics/navierStokesAdvectionDiffusionCouplingPostProcessor3D.h
• src/dynamics/navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh