Claculate Force Term in OpenLB (rayleigh benard problem)

[navidkyo]

Hello All,

I am not an advanced user. I faced problems while trying to calculate the body force term value in each step of the simulation using SuperLatticeGuoZhaoPhysBodyForce2D.

to calculate the body force term in each step: I have implemented these two lines into my code:

>>SuperLatticeGuoZhaoPhysBodyForce2D<T, NSDESCRIPTOR> bodyForce(NSlattice, converter);
>>vtkWriter.addFunctor( bodyForce );

but I get the following error while trying to execute the “make” command:

——————————————————————————————————————————————————————-
/home/user/OlbRoot/olb-1.4r0//src/functors/lattice/latticeGuoZhaoPhysBodyForce2D.hh:70:46: error: no match for ‘operator*’ (operand type is ‘olb::FieldPtr<double, olb::descriptors::D2Q9<olb::descriptors::FORCE>, olb::descriptors::BODY_FORCE>’ {aka ‘olb::VectorPtr<double, 2>’})
70 | output[0] = this->_converter.getPhysForce( *this->_blockLattice.get( input[0], input[1] )
| ^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
71 | .template getFieldPointer<descriptors::BODY_FORCE>() );
| ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
In file included from /home/user/OlbRoot/olb-1.4r0//src/utilities/vectorHelpers.h:35,
——————————————————————————————————————————————————————-

I would appreciate it if you could help me with this problem, or you can suggest to me a more suitable way to calculate the body force term.

Bests,
Navid

[jjessberger]

Hello Navid,

there seems to be a bug in the implementation of the method BlockLatticeGuoZhaoPhysBodyForce2D<T,DESCRIPTOR>::operator().

If you change the field access in the BlockLatticeGuoZhaoPhysBodyForce2D<T,DESCRIPTOR>::operator() method (src/functors/lattice/latticeGuoZhaoPhysBodyForce2D.hh, ll. 67-75) to

template <typename T, typename DESCRIPTOR>
bool BlockLatticeGuoZhaoPhysBodyForce2D<T,DESCRIPTOR>::operator() (T output[], const int input[])
{
output[0] = this->_converter.getPhysForce( this->_blockLattice.get( input[0], input[1] )
.template getFieldPointer<descriptors::BODY_FORCE>()[0] );
output[1] = this->_converter.getPhysForce( this->_blockLattice.get( input[0], input[1] )
.template getFieldPointer<descriptors::BODY_FORCE>()[1] );
return true;
}

, does it work for you?

Yours,
Julius

[navidkyo]

Dear Julius,

Thank you for your time and response. The problem with compiling went away but the simulation still crashes. Here is the error:

[MpiManager] Sucessfully initialized, numThreads=1
[ThermalUnitConv] —————– UnitConverter information —————–
[ThermalUnitConv] — Parameters:
[ThermalUnitConv] Resolution: N= 10
[ThermalUnitConv] Lattice velocity: latticeU= 0.1
[ThermalUnitConv] Lattice relaxation frequency: omega= 1.90375
[ThermalUnitConv] Lattice relaxation time: tau= 0.525278
[ThermalUnitConv] Thermal Lattice relaxation frequency: omega_AD= 1.86705
[ThermalUnitConv] Thermal Lattice relaxation time: tau_AD= 0.535603
[ThermalUnitConv] Characteristical length(m): charL= 0.1
[ThermalUnitConv] Characteristical speed(m/s): charU= 0.0118678
[ThermalUnitConv] Phys. kinematic viscosity(m^2/s): charNu= 1e-05
[ThermalUnitConv] Phys. density(kg/m^d): charRho= 1
[ThermalUnitConv] Characteristical pressure(N/m^2): charPressure= 0
[ThermalUnitConv] Reynolds number: reynoldsNumber= 118.678
[ThermalUnitConv] ————————————————————-
[ThermalUnitConv] —————– ThermalUnitConverter information —————–
[ThermalUnitConv] — Parameters:
[ThermalUnitConv] Phys. Delta X(m): physDeltaX= 0.01
[ThermalUnitConv] Phys. Delta T(s): physDeltaT= 0.0842615
[ThermalUnitConv] Characteristical pressure(N/m^2): charPressure= 0
[ThermalUnitConv] Phys. Thermal Conductivity(W/m/K): physThermalCondcticity= 0.03
[ThermalUnitConv] Phys. specific Heat Capacity(J/kg/K): physSpecificHeatCapacity= 2130
[ThermalUnitConv] Phys. Thermal Expasion Coefficent(K^-1): physThermalExpansionCoefficent= 0.000143573
[ThermalUnitConv] Characteristical Phys. low Temperature(K): charPhysLowTemperature= 273.15
[ThermalUnitConv] Characteristical Phys. high Temperature(K): charPhysHighTemperature= 274.15
[ThermalUnitConv] Prandtl number: prandtlNumber= 0.71
[ThermalUnitConv] Rayleigh number: rayleighNumber= 10000
[ThermalUnitConv] ————————————————————-
[ThermalUnitConv] —————– Conversion factors:—————–
[ThermalUnitConv] Voxel length(m): physDeltaX= 0.01
[ThermalUnitConv] Time step(s): physDeltaT= 0.0842615
[ThermalUnitConv] Velocity factor(m/s): physVelocity= 0.118678
[ThermalUnitConv] Density factor(kg/m^3): physDensity= 1
[ThermalUnitConv] Mass factor(kg): physMass= 1e-06
[ThermalUnitConv] Viscosity factor(m^2/s): physViscosity= 0.00118678
[ThermalUnitConv] Force factor(N): physForce= 1.40845e-06
[ThermalUnitConv] Pressure factor(N/m^2): physPressure= 0.0140845
[ThermalUnitConv] ————————————————————-
[ThermalUnitConv] —————– ThermalConversion factors:—————–
[ThermalUnitConv] Temperature(K): temperature= 1
[ThermalUnitConv] Thermal Duffusity(m^2/s): physThermalDiffusity= 0.00118678
[ThermalUnitConv] specific Heat Capacity(J/kg): physSpecificHeatCapacity= 0.0140845
[ThermalUnitConv] Thermal Coductivity(W/m/K): physThermalCondcticity= 1.67152e-05
[ThermalUnitConv] HeatFlux(W): physHeatFlux= 0.00167152
[ThermalUnitConv] ————————————————————-
[prepareGeometry] Prepare Geometry …
[SuperGeometry2D] cleaned 0 outer boundary voxel(s)
[SuperGeometry2D] cleaned 1 inner boundary voxel(s)
[SuperGeometry2D] the model is correct!
[CuboidGeometry2D] —Cuboid Stucture Statistics—
[CuboidGeometry2D] Number of Cuboids: 1
[CuboidGeometry2D] Delta (min): 0.01
[CuboidGeometry2D] (max): 0.01
[CuboidGeometry2D] Ratio (min): 1
[CuboidGeometry2D] (max): 1.9901
[CuboidGeometry2D] Nodes (min): 20301
[CuboidGeometry2D] (max): 20301
[CuboidGeometry2D] ——————————–
[SuperGeometryStatistics2D] materialNumber=1; count=19898; minPhysR=(0,0.01); maxPhysR=(2,0.99)
[SuperGeometryStatistics2D] materialNumber=2; count=201; minPhysR=(0,0); maxPhysR=(2,0)
[SuperGeometryStatistics2D] materialNumber=3; count=201; minPhysR=(0,1); maxPhysR=(2,1)
[SuperGeometryStatistics2D] materialNumber=4; count=1; minPhysR=(1,0.01); maxPhysR=(1,0.01)
[prepareGeometry] Prepare Geometry … OK
[prepareLattice] defining dynamics
[prepareLattice] Prepare Lattice … OK
[Timer] step=0; percent=0; passedTime=0.049; remTime=581.483; MLUPs=0
[LatticeStatistics] step=0; t=0; uMax=1.49167e-154; avEnergy=0; avRho=1
terminate called after throwing an instance of ‘std::invalid_argument’
what(): This MultiFieldArrayD does not provide FIELD.
[fedora:17569] *** Process received signal ***
[fedora:17569] Signal: Aborted (6)
[fedora:17569] Signal code: (-6)
[fedora:17569] [ 0] /lib64/libpthread.so.0(+0x13a20)[0x7fd0dea0da20]
[fedora:17569] [ 1] /lib64/libc.so.6(gsignal+0x142)[0x7fd0de8662a2]
[fedora:17569] [ 2] /lib64/libc.so.6(abort+0x116)[0x7fd0de84f8a4]
[fedora:17569] [ 3] /lib64/libstdc++.so.6(+0xa1a26)[0x7fd0dec1ba26]
[fedora:17569] [ 4] /lib64/libstdc++.so.6(+0xad24c)[0x7fd0dec2724c]
[fedora:17569] [ 5] /lib64/libstdc++.so.6(+0xad2b7)[0x7fd0dec272b7]
[fedora:17569] [ 6] /lib64/libstdc++.so.6(+0xad519)[0x7fd0dec27519]
[fedora:17569] [ 7] ./rayleighBenard2d[0x4880b1]
[fedora:17569] [ 8] ./rayleighBenard2d[0x485654]
[fedora:17569] [ 9] ./rayleighBenard2d[0x481c9e]
[fedora:17569] [10] ./rayleighBenard2d[0x47e656]
[fedora:17569] [11] ./rayleighBenard2d[0x474f16]
[fedora:17569] [12] ./rayleighBenard2d[0x47a5f7]
[fedora:17569] [13] ./rayleighBenard2d[0x421741]
[fedora:17569] [14] ./rayleighBenard2d[0x41741c]
[fedora:17569] [15] ./rayleighBenard2d[0x40ae14]
[fedora:17569] [16] ./rayleighBenard2d[0x40ba22]
[fedora:17569] [17] /lib64/libc.so.6(__libc_start_main+0xd5)[0x7fd0de850b75]
[fedora:17569] [18] ./rayleighBenard2d[0x4092ae]
[fedora:17569] *** End of error message ***
Aborted (core dumped)

Bests,
Navid

[mathias]

Material 4, there is just one voxel. Numerically, that is not enough to have an outflow or inflow with only one voxel. Best Mathias

[navidkyo]

Hello Mathias,

This is the RayleighBernard2D example of OpenLB V1.4. Material 4 is defined to implement a perturbation; it is not inflow or outflow.

I will redefine material 4 and let you know if the error is resolved.

[Author]

Posts