OpenLB – Open Source Lattice Boltzmann Code Forums on OpenLB General Topics WARNING: no discreteNormal is found

Viewing 4 posts - 1 through 4 (of 4 total)
• Author
Posts
• #6241
Anand
Participant

Hello,
I was preparing the Pipe Geometry and Lattice.
If I choose the pipe diameter 0.04 m and length according to the entrance length for fully developed turbulent flow, I get a “WARNING: no discreteNormal is found” error.

If I slightly increase the diameter (0.0405 m), it works.
I do not understand why. Can you please explain it?

Thank you

#6242
Keymaster

This very likely is a geometry setup artifact (some boundary setter is not able to reconstruct a discrete normal for one or more cell locations). Does the warning also disappear if you keep the diameter at 0.04 and instead tweak the discretization factor (depending on how this is set up in your case, e.g. N+1 instead of N cells per unit length)?

I hope this helps, hard to tell more without further details of your setup.

#6243
Anand
Participant

Here I have copy past the code, if you wish, please have alook….

#include “olb3D.h”
#ifndef OLB_PRECOMPILED // Unless precompiled version is used,
#include “olb3D.hh” // include full template code
#endif
#include <vector>
#include <cmath>
#include <iostream>
#include <iomanip>
#include <fstream>

using namespace olb;
using namespace olb::descriptors;
using namespace olb::graphics;
using namespace olb::util;
using namespace std;
typedef double T;
const T maxPhysT = 10.; // max. simulation time in s, SI unit
int N = 60; // resolution of the model
const T density = 1168; // Liquid Density (kg/m3)
const T viscosity = 0.0044/density; // Liquid kinematic viscosity (m2/s)
const T inletvelocity = 0.774; // Inlet Velocity (m/s)
const T pipediameter = 0.0405; // Pipe Diameter (m)
const T Re = (inletvelocity*pipediameter)/(viscosity); // defined as 1/kinematic viscosity
const T pipelength = 1.5*(1.359*pipediameter*std::pow(Re, 0.25)); // Pipe Length (m) calculated using turbulence entrace length formula plus 1.5 time higher

// Stores geometry information in form of material numbers
void prepareGeometry( UnitConverter<T,DESCRIPTOR> const& converter,
SuperGeometry3D<T>& superGeometry )
{

OstreamManager clout(std::cout, “prepareGeometry”);

clout << “Prepare Geometry …” << std::endl;

superGeometry.rename(0, 2);

superGeometry.rename(2, 1, pipe);

Vector<T, 3> extend = origin;

// Set material number for inflow
origin = -converter.getConversionFactorLength() * 2;
extend = converter.getConversionFactorLength() * 2;
superGeometry.rename(2, 3, 1, inflow);

// Set material number for outflow
origin = pipelength – 2 * converter.getConversionFactorLength();
extend = pipelength + 2 * converter.getConversionFactorLength();
superGeometry.rename(2, 4, 1, outflow);

// Removes all not needed boundary voxels outside the surface
superGeometry.clean();
// Removes all not needed boundary voxels inside the surface
superGeometry.innerClean();
superGeometry.checkForErrors();

superGeometry.print();

clout << “Prepare Geometry … OK” << std::endl;
}

// Set up the geometry of the simulation
void prepareLattice(SuperLattice3D<T, DESCRIPTOR>& sLattice,
UnitConverter<T, DESCRIPTOR>const& converter,
Dynamics<T, DESCRIPTOR>& bulkDynamics,
SuperGeometry3D<T>& superGeometry)
{

OstreamManager clout( std::cout,”prepareLattice” );
clout << “Prepare Lattice …” << std::endl;

const T omega = converter.getLatticeRelaxationFrequency();

// Material=0 –>do nothing
sLattice.defineDynamics( superGeometry, 0, &instances::getNoDynamics<T, DESCRIPTOR>() );

// Material=1 –>bulk dynamics
sLattice.defineDynamics( superGeometry, 1, &bulkDynamics );

// material=2 –> depends on boundary type

sLattice.defineDynamics( superGeometry, 2, &instances::getBounceBack<T, DESCRIPTOR>() );

// Material=3 –>bulk dynamics
sLattice.defineDynamics( superGeometry, 3, &bulkDynamics );
setInterpolatedVelocityBoundary<T,DESCRIPTOR>(sLattice, omega, superGeometry, 3);

// Material=4 –>bulk dynamics
sLattice.defineDynamics( superGeometry, 4, &bulkDynamics );
setInterpolatedPressureBoundary<T,DESCRIPTOR>(sLattice, omega, superGeometry, 4);
// Initial conditions
AnalyticalConst3D<T,T> ux( 0. );
AnalyticalConst3D<T,T> uy( 0. );
AnalyticalConst3D<T,T> uz( 0. );
AnalyticalConst3D<T,T> rho( 1. );
AnalyticalComposed3D<T,T> u( ux,uy,uz );

//Initialize all values of distribution functions to their local equilibrium
sLattice.defineRhoU( superGeometry.getMaterialIndicator({1, 3, 4}), rho, u );
sLattice.iniEquilibrium( superGeometry.getMaterialIndicator({1, 3, 4}), rho, u );

// Make the lattice ready for simulation
sLattice.initialize();

clout << “Prepare Lattice … OK” << std::endl;

}

int main( int argc, char* argv[] )
{
/// === 1st Step: Initialization ===
olbInit(&argc, &argv);
singleton::directories().setOutputDir(“./tmp/”);
OstreamManager clout( std::cout,”main” );
clout.setMultiOutput(true);

UnitConverterFromResolutionAndRelaxationTime<T,DESCRIPTOR> converter(
int {N}, // resolution: number of voxels per charPhysL
(T) 0.5018, // latticeRelaxationTime: relaxation time, have to be greater than 0.5!
(T) pipediameter, // charPhysLength: reference length of simulation geometry
(T) inletvelocity, // charPhysVelocity: maximal/highest expected velocity during simulation in __m / s__
(T) viscosity, // physViscosity: physical kinematic viscosity in __m^2 / s__
(T) density // physDensity: physical density in __kg / m^3__
);

// Prints the converter log as console output
converter.print();
converter.write(“Pipe”);

// === 2nd Step: Prepare Geometry ===

IndicatorLayer3D<T> extendedDomain(pipe, converter.getConversionFactorLength());

// Instantiation of a cuboidGeometry with weights
#ifdef PARALLEL_MODE_MPI
const int noOfCuboids = std::min( 16*N,2*singleton::mpi().getSize() );
#else
const int noOfCuboids = 10;
#endif // ifdef PARALLEL_MODE_MPI
CuboidGeometry3D<T> cuboidGeometry(extendedDomain, converter.getConversionFactorLength(), noOfCuboids);

// Instantiation of a superGeometry

prepareGeometry(converter, superGeometry);

// === 3rd Step: Prepare Lattice ===
SuperLattice3D<T, DESCRIPTOR> sLattice( superGeometry );
const T omega = converter.getLatticeRelaxationFrequency();

SmagorinskyBGKdynamics<T, DESCRIPTOR>bulkDynamics(omega, instances::getBulkMomenta<T, DESCRIPTOR>(),0.1);
}

#6247
Anand
Participant

I managed to resolve this problem by changing the resolutions.

Viewing 4 posts - 1 through 4 (of 4 total)
• You must be logged in to reply to this topic.