[kolotinsky]

What does “collides without updating u” mean? Does it mean that equilibrium distribution function in collision BGK operator which is f0 = f0(u) does not change during a simulation. And for example I can set u to zero?

[kolotinsky]

I have found interesting method named “staticCollide” in doxygen. I suppose it can help me for my purpose. But I do not really understand what this method do. Could you explain it for me?

[kolotinsky]

I understand this. I mean, could you advice me which file should I rewrite to change equilibrium distribution function for MRT or BGK forced dynamic?

[kolotinsky]

My aim is to simulate the following equation
df/dt + v*df/dr + (F/m)*df/dv = w*(f-f0)
As I understand in OpenLb algorithm for this case f0 is shifted maxwell distribution.
I suppose to simulate collisionless ion gas in neutral gas it is enough to replace shifted maxwell distribution function of ions with maxwell distribution function of neutral gas. How can I apply my idea using OpenLB? Is it possible to change some source file to implement my idea?

[kolotinsky]

So, if I understand correctly, then using OpenLB it is not possible to simulate a collisionless gas, isn’t it?

[kolotinsky]

If I add viscosity to external force term to describe collisions with neutral gas can I set zero relaxation frequency in MRT operator or simply use BGK operator with zero relaxation frequency?

[kolotinsky]

Oh, thank you very much, such a stupid mistake!

[kolotinsky]

Dear Mathias,
Of course, velocity must go to infinity, I missed up it yesterday. In my code I used forcedMRT dynamics.
Here is my code

//////////////////////////////////////////////////////////////////////////
#include “olb3D.h”
#include “olb3D.hh” // use only generic version!
#include <cstdlib>
#include <iostream>

using namespace olb;
using namespace olb::descriptors;
using namespace olb::graphics;
using namespace std;

typedef double T;
#define DESCRIPTOR ForcedMRTD3Q19Descriptor

// Parameters for the simulation setup
const int maxIter = 2000;
const int nx = 30;
const int ny = 30;
const int nz = 30;

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

OstreamManager clout( std::cout,”prepareGeometry” );
clout << “Prepare Geometry …” << std::endl;

// Sets material number for fluid
superGeometry.rename( 0,1 );

// 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,
Dynamics<T, DESCRIPTOR>& bulkDynamics1,
SuperGeometry3D<T>& superGeometry ) {

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

// Initial conditions

std::vector<T> v( 3,T() );
AnalyticalConst3D<T,T> zeroVelocity( v );
AnalyticalConst3D<T,T> Rho( 1. );

// Initialize force

std::vector<T> force(3, T());
force[0] = 0;
force[1] = 0;
force[2] = 0.00075;

//olb::BlockVTIreader3D< T,BaseType > datareader(“forcefield.vti”,”ForceField”);
//BlockDataF3D<T, BaseType> blockfield(datareader.getBlockData());
//AnalyticalF3D<T,T>* field;
//field = new SpecialAnalyticalFfromBlockF3D<T,T>(blockfield, datareader.getCuboid(),spacing);
AnalyticalConst3D<T,T> constant(force);
sLattice.defineExternalField( superGeometry, 1,
DESCRIPTOR<T>::ExternalField::forceBeginsAt,
DESCRIPTOR<T>::ExternalField::sizeOfForce, constant );

// Initialize all values of distribution functions to their local equilibrium
sLattice.defineRhoU( superGeometry, 1, Rho, zeroVelocity );
sLattice.iniEquilibrium( superGeometry, 1, Rho, zeroVelocity );

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

//delete field;
//field = 0;
}

// Output to console and files
void getResults( SuperLattice3D<T, DESCRIPTOR>& sLattice, int iT,
SuperGeometry3D<T>& superGeometry, Timer<T>& timer ) {

OstreamManager clout( std::cout,”getResults” );

SuperVTMwriter3D<T> vtmWriter( “ionwake3d” );
SuperLatticeDensity3D<T, DESCRIPTOR> density( sLattice );
vtmWriter.addFunctor( density );

const int statIter = 1;
vtmWriter.createMasterFile();
vtmWriter.write();

// Writes output on the console
if ( iT%statIter==0 ) {
// Timer console output
timer.update( iT );
timer.printStep();

// Lattice statistics console output
sLattice.getStatistics().print( iT,iT );
}
}

int main( int argc, char *argv[] ) {

// === 1st Step: Initialization ===
olbInit( &argc, &argv );
singleton::directories().setOutputDir( “./tmp/” );
OstreamManager clout( std::cout,”main” );

const T omega1 = 0.0041;

// === 2rd Step: Prepare Geometry ===

const int noOfCuboids = 1;

CuboidGeometry3D<T> cuboidGeometry( 0, 0, 0, 1, nx, ny, nz, noOfCuboids );

// Periodic boundaries in x- and y- and z-direction
cuboidGeometry.setPeriodicity( true, true, true );

// Instantiation of a loadBalancer
HeuristicLoadBalancer<T> loadBalancer( cuboidGeometry );

// Instantiation of a superGeometry
SuperGeometry3D<T> superGeometry( cuboidGeometry,loadBalancer,2 );

prepareGeometry( superGeometry );

// === 3rd Step: Prepare Lattice ===
SuperLattice3D<T, DESCRIPTOR> sLattice( superGeometry );

ForcedMRTdynamics<T, DESCRIPTOR> bulkDynamics1 ( omega1, instances::getBulkMomenta<T, DESCRIPTOR>() );

prepareLattice( sLattice, bulkDynamics1, superGeometry );

// === 4th Step: Main Loop ===
int iT = 0;
clout << “starting simulation…” << endl;
Timer<T> timer( maxIter, superGeometry.getStatistics().getNvoxel() );
timer.start();

//std::vector<T> spacing(3, T());
//spacing[0] = 1;
//spacing[1] = 1;
//spacing[2] = 1;
for ( iT = 0; iT < maxIter; ++iT ) {

// === 5th Step: Collide and Stream Execution ===
//olb::BlockVTIreader3D< T,BaseType > datareader(“forcefield.vti”,”ForceField”);
//BlockDataF3D<T, BaseType> blockfield(datareader.getBlockData());
//AnalyticalF3D<T,T>* field;
//field = new SpecialAnalyticalFfromBlockF3D<T,T>(blockfield, datareader.getCuboid(),spacing);
//sLattice.defineExternalField( superGeometry, 1,
// DESCRIPTOR<T>::ExternalField::forceBeginsAt,
// DESCRIPTOR<T>::ExternalField::sizeOfForce, *field );
//sLattice.collideAndStream();
//delete field;
//field = 0;

// === 7th Step: Computation and Output of the Results ===
getResults( sLattice, iT, superGeometry, timer );
}

timer.stop();
timer.printSummary();
}
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

Best regards,
Kolotinskiy

[kolotinsky]

Dear Fabian,
I have changed my vti data from Appended to Ascii:

/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
<VTKFile type=”ImageData” version=”1.0″ byte_order=”LittleEndian” header_type=”UInt64″>
<ImageData WholeExtent=”-1 50 -1 50 -1 50″ Origin=”0 0 0″ Spacing=”1 1 1″>
<Piece Extent=”-1 50 -1 50 -1 50″>
<PointData>
<DataArray type=”Float64″ Name=”ForceField” NumberOfComponents=”3″ format=”ascii” RangeMin=”0.00014142135266500606″ RangeMax=”0.00014142135266500606″>
0.000099999997474 0.000099999997474 0 0.000099999997474 0.000099999997474 0
0.000099999997474 0.000099999997474 0 0.000099999997474 0.000099999997474 0
……………………………………………………………………………………………………………………………………………
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

Now vtiReader reads my data correctly, but on the step “define external field” I have segmentation fault in case of such
definition:

//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
std::vector<T> spacing(3, T());
spacing[0] = 1;
spacing[1] = 1;
spacing[2] = 1;
olb::BlockVTIreader3D< T, BaseType > datareader(“forcefield.vti”,”ForceField”);
BlockDataF3D<T, BaseType> blockfield(datareader.getBlockData());
AnalyticalF3D<T,T>* field;
field = new SpecialAnalyticalFfromBlockF3D<T,BaseType>(blockfield, datareader.getCuboid(), spacing);
sLattice.defineExternalField( superGeometry, 1,
DESCRIPTOR<T>::ExternalField::forceBeginsAt,
DESCRIPTOR<T>::ExternalField::sizeOfForce, *field );
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

I also tried to define force field such way:

//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
olb::BlockVTIreader3D< T, BaseType > datareader(“forcefield.vti”,”ForceField”);
BlockDataF3D<T, BaseType> blockfield(datareader.getBlockData());
AnalyticalF3D<T,T>* field;
field = new AnalyticalFfromBlockF3D<T,BaseType>(blockfield, datareader.getCuboid(), 0);
sLattice.defineExternalField( superGeometry, 1,
DESCRIPTOR<T>::ExternalField::forceBeginsAt,
DESCRIPTOR<T>::ExternalField::sizeOfForce, *field );
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

In this case I don’t get segmentation fault, but my system’s behaviour is wrong. To find out if the external field was written correctly I generated file contained external field using SuperVtmWriter3D:

//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
SuperVTMwriter3D<T> vtmWriter( “ionwake3d” );

SuperLatticeExternalField3D<T, DESCRIPTOR> force( sLattice, DESCRIPTOR<T>::ExternalField::forceBeginsAt,
DESCRIPTOR<T>::ExternalField::sizeOfForce );
vtmWriter.addFunctor( force );
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

This generated file doesn’t math up with uploaded:

/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

<VTKFile type=”ImageData” version=”1.0″ byte_order=”LittleEndian” header_type=”UInt64″>
<ImageData WholeExtent=”-1 50 -1 50 -1 50″ Origin=”0 0 0″ Spacing=”1 1 1″>
<Piece Extent=”-1 50 -1 50 -1 50″>
<PointData>
<DataArray type=”Float32″ Name=”externalField” NumberOfComponents=”3″ format=”ascii” RangeMin=”0″ RangeMax=”18.759683655122647″>
0 0 0 0 0 0
0.26010000705718994 0.26010000705718994 0 0.5202000141143799 0.5202000141143799 0
0.780299961566925 0.780299961566925 0 1.0404000282287598 1.0404000282287598 0
1.3004999160766602 1.3004999160766602 0 1.56059992313385 1.56059992313385 0
……………………………………………………………………………………………………………………………………………………………..

/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

I think it is also important to say a few words about my system. Now I’m using uniform 50*50*50 grid with periodic boundary conditions.

By the way , is it comfortable for you to communicate via email or any social network? It will allow me to describe the problem in more details and attach required files. I suppose it will make our dialog more effective.
If yes this is my mail address: kolotinsky.daniil1998@gmail.com

Sincerely,
Daniil Kolotinsky

[kolotinsky]

Dear Fabian, following your example I have written code for my case:

//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Initialize force

std::vector<T> spacing(3, T());
spacing[0] = 1;
spacing[1] = 1;
spacing[2] = 1;

olb::BlockVTIreader3D< T, BaseType > datareader(“forcefield.vti”,”AppendedData”);
BlockDataF3D<T, BaseType> blockfield(datareader.getBlockData());
AnalyticalF3D<T,T>* field;
field = new SpecialAnalyticalFfromBlockF3D<T,BaseType>(blockfield, datareader.getCuboid(), spacing);
sLattice.defineExternalField( superGeometry, 1,
DESCRIPTOR<T>::ExternalField::forceBeginsAt,
DESCRIPTOR<T>::ExternalField::sizeOfForce, *field );
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

However, during the simulation I am receiving wrong zero or NaN values of Energy and Umax:

////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
[LatticeStatistics] step=190; t=190; uMax=0; avEnergy=nan; avRho=nan
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

I have verified my vti file in Paraview but everything is ok.

I also tried to write the data uploaded from my vti file using vtiWriter provided by OpenLB:

///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
olb::VTIwriter3D< T, BaseType > datawriter;
datawriter.writeData(“field”,”AppendedData”,datareader.getBlockData(),datareader.getCuboid());
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

Generated force.vti file contains only the following data:

/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
<?xml version=”1.0″?>
<VTKFile type=”ImageData” version=”0.1″ byte_order=”LittleEndian”>
<ImageData WholeExtent=” 0 49 0 49 0 49 ” Origin=”0 0 0″ Spacing=”1 1 1″>
<Piece Extent=”0 49 0 49 0 49″>
<PointData>
<DataArray type=”Float32″ Name=”AppendedData” NumberOfComponents=”0″>
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

What is wrong?

[kolotinsky]

As I understand I need to load my data into an AnalyticalF3D object. Initially my force field is located in a vti file. Thus my key question is: how to fill AnalyticalF3D object using the data from the vti file?

[kolotinsky]

Dear mlm,
could you tell me in more detail how can I create my own VTI file and in particular how can I read the information from this file and apply it for my purpose?

[kolotinsky]

Thank you very much, Marc! Your advice has helped!

[kolotinsky]

Quote:

Quote from kolotinsky on September 18, 2018, 11:29
I have recently installed openlb package on mac os. And I have faced the following problem. During the execution of examples I have noticed the warning in console output “”./tmp/imageData/data/heatMapl2physVelocity.p”, line 5: unexpected or unrecognized token”. Then when I have tried to open Jpeg output files system has given the error “could not open the file because it is empty”. Could you give me an advice, what should I do with it?

[Author]

Posts