From 527e8aa9c33d7c82b548b2280c979f653a7131ae Mon Sep 17 00:00:00 2001 From: Adrian Kummerlaender Date: Tue, 5 Nov 2019 20:14:29 +0100 Subject: Add basic AA pattern support for C++ example --- lid_driven_cavity/cpp/CMakeLists.txt | 5 +- lid_driven_cavity/cpp/config.py | 4 +- lid_driven_cavity/cpp/generate.py | 3 +- lid_driven_cavity/cpp/ldc.cpp.AA.mako | 202 ++++++++++++++++++++++++++++++++++ lid_driven_cavity/cpp/ldc.cpp.AB.mako | 178 ++++++++++++++++++++++++++++++ lid_driven_cavity/cpp/ldc.cpp.mako | 178 ------------------------------ 6 files changed, 389 insertions(+), 181 deletions(-) create mode 100644 lid_driven_cavity/cpp/ldc.cpp.AA.mako create mode 100644 lid_driven_cavity/cpp/ldc.cpp.AB.mako delete mode 100644 lid_driven_cavity/cpp/ldc.cpp.mako (limited to 'lid_driven_cavity/cpp') diff --git a/lid_driven_cavity/cpp/CMakeLists.txt b/lid_driven_cavity/cpp/CMakeLists.txt index 5b5fb90..cc526e3 100644 --- a/lid_driven_cavity/cpp/CMakeLists.txt +++ b/lid_driven_cavity/cpp/CMakeLists.txt @@ -16,7 +16,10 @@ add_custom_command( WORKING_DIRECTORY ${CMAKE_CURRENT_SOURCE_DIR} DEPENDS - generate.py config.py ldc.cpp.mako + generate.py + config.py + ldc.cpp.AA.mako + ldc.cpp.AB.mako ) include_directories( diff --git a/lid_driven_cavity/cpp/config.py b/lid_driven_cavity/cpp/config.py index 2759ebf..4ea30b4 100644 --- a/lid_driven_cavity/cpp/config.py +++ b/lid_driven_cavity/cpp/config.py @@ -4,10 +4,12 @@ from boltzgen import Geometry descriptor = D2Q9 geometry = Geometry(256, 256) tau = 0.52 -precision = 'single' +precision = 'double' +streaming = 'AA' ## 3D LDC #descriptor = D3Q19 #geometry = Geometry(64, 64, 64) #tau = 0.52 #precision = 'single' +#streaming = 'AA' diff --git a/lid_driven_cavity/cpp/generate.py b/lid_driven_cavity/cpp/generate.py index 27c7997..dae2514 100755 --- a/lid_driven_cavity/cpp/generate.py +++ b/lid_driven_cavity/cpp/generate.py @@ -18,6 +18,7 @@ generator = Generator( model = BGK(config.descriptor, tau = config.tau), target = 'cpp', precision = config.precision, + streaming = config.streaming, index = 'XYZ', layout = 'AOS') @@ -30,7 +31,7 @@ with open('%s/kernel.h' % args.output, 'w') as kernel: kernel.write(generator.kernel(config.geometry, functions)) ldc_src = '' -with open('ldc.cpp.mako', 'r') as template: +with open('ldc.cpp.%s.mako' % config.streaming, 'r') as template: ldc_src = template.read() with open('%s/ldc.cpp' % args.output, 'w') as app: diff --git a/lid_driven_cavity/cpp/ldc.cpp.AA.mako b/lid_driven_cavity/cpp/ldc.cpp.AA.mako new file mode 100644 index 0000000..a35e8f0 --- /dev/null +++ b/lid_driven_cavity/cpp/ldc.cpp.AA.mako @@ -0,0 +1,202 @@ +#include +#include +#include +#include +#include +#include + +#include "kernel.h" + +void collect_moments_to_vtk(const std::string& path, ${float_type}* f) { + std::ofstream fout; + fout.open(path.c_str()); + + fout << "# vtk DataFile Version 3.0\n"; + fout << "lbm_output\n"; + fout << "ASCII\n"; + fout << "DATASET RECTILINEAR_GRID\n"; +% if descriptor.d == 2: + fout << "DIMENSIONS " << ${geometry.size_x-2} << " " << ${geometry.size_y-2} << " 1" << "\n"; +% else: + fout << "DIMENSIONS " << ${geometry.size_x-2} << " " << ${geometry.size_y-2} << " " << ${geometry.size_z-2} << "\n"; +% endif + + fout << "X_COORDINATES " << ${geometry.size_x-2} << " float\n"; + for( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + fout << x << " "; + } + + fout << "\nY_COORDINATES " << ${geometry.size_y-2} << " float\n"; + for( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + fout << y << " "; + } + +% if descriptor.d == 2: + fout << "\nZ_COORDINATES " << 1 << " float\n"; + fout << 0 << "\n"; + fout << "POINT_DATA " << ${(geometry.size_x-2) * (geometry.size_y-2)} << "\n"; +% else: + fout << "\nZ_COORDINATES " << ${geometry.size_z-2} << " float\n"; + for( std::size_t z = 1; z < ${geometry.size_z-1}; ++z ) { + fout << z << " "; + } + fout << "\nPOINT_DATA " << ${(geometry.size_x-2) * (geometry.size_y-2) * (geometry.size_z-2)} << "\n"; +% endif + + ${float_type} rho; + ${float_type} u[${descriptor.d}]; + + fout << "VECTORS velocity float\n"; +% if descriptor.d == 2: + for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + collect_moments_tock(f, x*${geometry.size_y}+y, rho, u); + fout << u[0] << " " << u[1] << " 0\n"; + } + } +% else: + for ( std::size_t z = 1; z < ${geometry.size_z-1}; ++z ) { + for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + collect_moments_tock(f, x*${geometry.size_y*geometry.size_z}+y*${geometry.size_z}+z, rho, u); + fout << u[0] << " " << u[1] << " " << u[2] << "\n"; + } + } + } +% endif + + fout << "SCALARS density float 1\n"; + fout << "LOOKUP_TABLE default\n"; +% if descriptor.d == 2: + for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + collect_moments_tock(f, x*${geometry.size_y}+y, rho, u); + fout << rho << "\n"; + } + } +% else: + for ( std::size_t z = 1; z < ${geometry.size_z-1}; ++z ) { + for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + collect_moments_tock(f, x*${geometry.size_y*geometry.size_z}+y*${geometry.size_z}+z, rho, u); + fout << rho << "\n"; + } + } + } +% endif + + fout.close(); +} + +void simulate(std::size_t nStep) +{ + auto f_aa = std::make_unique<${float_type}[]>(${geometry.volume*descriptor.q}); + + ${float_type}* f = f_aa.get(); + + std::vector ghost; + std::vector bulk; + std::vector lid_bc; + std::vector box_bc; + + for (int iX = 0; iX < ${geometry.size_x}; ++iX) { + for (int iY = 0; iY < ${geometry.size_y}; ++iY) { +% if descriptor.d == 2: + const std::size_t iCell = iX*${geometry.size_y} + iY; + if (iX == 0 || iY == 0 || iX == ${geometry.size_x-1} || iY == ${geometry.size_y-1}) { + ghost.emplace_back(iCell); + } else if (iY == ${geometry.size_y-2}) { + lid_bc.emplace_back(iCell); + } else if (iX == 1 || iX == ${geometry.size_x-2} || iY == 1) { + box_bc.emplace_back(iCell); + } else { + bulk.emplace_back(iCell); + } +% elif descriptor.d == 3: + for (int iZ = 0; iZ < ${geometry.size_z}; ++iZ) { + const std::size_t iCell = iX*${geometry.size_y*geometry.size_z} + iY*${geometry.size_z} + iZ; + if ( iX == 0 || iY == 0 || iZ == 0 + || iX == ${geometry.size_x-1} + || iY == ${geometry.size_y-1} + || iZ == ${geometry.size_z-1}) { + ghost.emplace_back(iCell); + } else if (iZ == ${geometry.size_z-2}) { + lid_bc.emplace_back(iCell); + } else if ( iX == 1 || iX == ${geometry.size_x-2} + || iY == 1 || iY == ${geometry.size_y-2} + || iZ == 1) { + box_bc.emplace_back(iCell); + } else { + bulk.emplace_back(iCell); + } + } +% endif + } + } + + std::cout << "#ghost : " << ghost.size() << std::endl; + std::cout << "#bulk : " << bulk.size() << std::endl; + std::cout << "#lid : " << lid_bc.size() << std::endl; + std::cout << "#wall : " << box_bc.size() << std::endl; + std::cout << std::endl; + + #pragma omp parallel for + for (std::size_t iCell = 0; iCell < ${geometry.volume}; ++iCell) { + equilibrilize_tick(f, iCell); + } + + auto start = std::chrono::high_resolution_clock::now(); + + for (std::size_t iStep = 1; iStep <= nStep; ++iStep) { + if (iStep % 2 == 0) { + #pragma omp parallel for + for (std::size_t i = 0; i < bulk.size(); ++i) { + collide_and_stream_tick(f, bulk[i]); + } + ${float_type} u[${descriptor.d}] { 0. }; + #pragma omp parallel for + for (std::size_t i = 0; i < box_bc.size(); ++i) { + velocity_momenta_boundary_tick(f, box_bc[i], u); + } + u[0] = 0.05; + #pragma omp parallel for + for (std::size_t i = 0; i < lid_bc.size(); ++i) { + velocity_momenta_boundary_tick(f, lid_bc[i], u); + } + + } else { + #pragma omp parallel for + for (std::size_t i = 0; i < ghost.size(); ++i) { + equilibrilize_tick(f, ghost[i]); + } + #pragma omp parallel for + for (std::size_t i = 0; i < bulk.size(); ++i) { + collide_and_stream_tock(f, bulk[i]); + } + ${float_type} u[${descriptor.d}] { 0. }; + #pragma omp parallel for + for (std::size_t i = 0; i < box_bc.size(); ++i) { + velocity_momenta_boundary_tock(f, box_bc[i], u); + } + u[0] = 0.05; + #pragma omp parallel for + for (std::size_t i = 0; i < lid_bc.size(); ++i) { + velocity_momenta_boundary_tock(f, lid_bc[i], u); + } + } + + if (iStep % 1000 == 0) { + auto duration = std::chrono::duration_cast>( + std::chrono::high_resolution_clock::now() - start); + std::cout << "iStep = " << iStep << "; ~" << 1000*${geometry.volume}/(1e6*duration.count()) << " MLUPS" << std::endl; + + collect_moments_to_vtk("result/ldc_" + std::to_string(iStep) + ".vtk", f); + + start = std::chrono::high_resolution_clock::now(); + } + } +} + +int main() { + simulate(20000); +} diff --git a/lid_driven_cavity/cpp/ldc.cpp.AB.mako b/lid_driven_cavity/cpp/ldc.cpp.AB.mako new file mode 100644 index 0000000..c66679f --- /dev/null +++ b/lid_driven_cavity/cpp/ldc.cpp.AB.mako @@ -0,0 +1,178 @@ +#include +#include +#include +#include +#include +#include + +#include "kernel.h" + +void collect_moments_to_vtk(const std::string& path, ${float_type}* f) { + std::ofstream fout; + fout.open(path.c_str()); + + fout << "# vtk DataFile Version 3.0\n"; + fout << "lbm_output\n"; + fout << "ASCII\n"; + fout << "DATASET RECTILINEAR_GRID\n"; +% if descriptor.d == 2: + fout << "DIMENSIONS " << ${geometry.size_x-2} << " " << ${geometry.size_y-2} << " 1" << "\n"; +% else: + fout << "DIMENSIONS " << ${geometry.size_x-2} << " " << ${geometry.size_y-2} << " " << ${geometry.size_z-2} << "\n"; +% endif + + fout << "X_COORDINATES " << ${geometry.size_x-2} << " float\n"; + for( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + fout << x << " "; + } + + fout << "\nY_COORDINATES " << ${geometry.size_y-2} << " float\n"; + for( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + fout << y << " "; + } + +% if descriptor.d == 2: + fout << "\nZ_COORDINATES " << 1 << " float\n"; + fout << 0 << "\n"; + fout << "POINT_DATA " << ${(geometry.size_x-2) * (geometry.size_y-2)} << "\n"; +% else: + fout << "\nZ_COORDINATES " << ${geometry.size_z-2} << " float\n"; + for( std::size_t z = 1; z < ${geometry.size_z-1}; ++z ) { + fout << z << " "; + } + fout << "\nPOINT_DATA " << ${(geometry.size_x-2) * (geometry.size_y-2) * (geometry.size_z-2)} << "\n"; +% endif + + ${float_type} rho; + ${float_type} u[${descriptor.d}]; + + fout << "VECTORS velocity float\n"; +% if descriptor.d == 2: + for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + collect_moments(f, x*${geometry.size_y}+y, rho, u); + fout << u[0] << " " << u[1] << " 0\n"; + } + } +% else: + for ( std::size_t z = 1; z < ${geometry.size_z-1}; ++z ) { + for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + collect_moments(f, x*${geometry.size_y*geometry.size_z}+y*${geometry.size_z}+z, rho, u); + fout << u[0] << " " << u[1] << " " << u[2] << "\n"; + } + } + } +% endif + + fout << "SCALARS density float 1\n"; + fout << "LOOKUP_TABLE default\n"; +% if descriptor.d == 2: + for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + collect_moments(f, x*${geometry.size_y}+y, rho, u); + fout << rho << "\n"; + } + } +% else: + for ( std::size_t z = 1; z < ${geometry.size_z-1}; ++z ) { + for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { + for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { + collect_moments(f, x*${geometry.size_y*geometry.size_z}+y*${geometry.size_z}+z, rho, u); + fout << rho << "\n"; + } + } + } +% endif + + fout.close(); +} + +void simulate(std::size_t nStep) +{ + auto f_a = std::make_unique<${float_type}[]>(${geometry.volume*descriptor.q}); + auto f_b = std::make_unique<${float_type}[]>(${geometry.volume*descriptor.q}); + + ${float_type}* f_prev = f_a.get(); + ${float_type}* f_next = f_b.get(); + + std::vector bulk; + std::vector lid_bc; + std::vector box_bc; + + for (int iX = 1; iX < ${geometry.size_x-1}; ++iX) { + for (int iY = 1; iY < ${geometry.size_y-1}; ++iY) { +% if descriptor.d == 2: + const std::size_t iCell = iX*${geometry.size_y} + iY; + if (iY == ${geometry.size_y-2}) { + lid_bc.emplace_back(iCell); + } else if (iX == 1 || iX == ${geometry.size_x-2} || iY == 1) { + box_bc.emplace_back(iCell); + } else { + bulk.emplace_back(iCell); + } +% elif descriptor.d == 3: + for (int iZ = 0; iZ < ${geometry.size_z}; ++iZ) { + const std::size_t iCell = iX*${geometry.size_y*geometry.size_z} + iY*${geometry.size_z} + iZ; + if (iZ == ${geometry.size_z-2}) { + lid_bc.emplace_back(iCell); + } else if (iX == 1 || iX == ${geometry.size_x-2} || iY == 1 || iY == ${geometry.size_y-2} || iZ == 1) { + box_bc.emplace_back(iCell); + } else { + bulk.emplace_back(iCell); + } + } +% endif + } + } + + std::cout << "#bulk : " << bulk.size() << std::endl; + std::cout << "#lid : " << lid_bc.size() << std::endl; + std::cout << "#wall : " << box_bc.size() << std::endl; + std::cout << "#steps : " << nStep << std::endl; + std::cout << std::endl; + +#pragma omp parallel for + for (std::size_t iCell = 0; iCell < ${geometry.volume}; ++iCell) { + equilibrilize(f_prev, f_next, iCell); + equilibrilize(f_next, f_prev, iCell); + } + + const auto start = std::chrono::high_resolution_clock::now(); + + for (std::size_t iStep = 0; iStep < nStep; ++iStep) { + if (iStep % 2 == 0) { + f_next = f_a.get(); + f_prev = f_b.get(); + } else { + f_next = f_b.get(); + f_prev = f_a.get(); + } + +#pragma omp parallel for + for (std::size_t i = 0; i < bulk.size(); ++i) { + collide_and_stream(f_next, f_prev, bulk[i]); + } + ${float_type} u[${descriptor.d}] { 0. }; +#pragma omp parallel for + for (std::size_t i = 0; i < box_bc.size(); ++i) { + velocity_momenta_boundary(f_next, f_prev, box_bc[i], u); + } + u[0] = 0.05; +#pragma omp parallel for + for (std::size_t i = 0; i < lid_bc.size(); ++i) { + velocity_momenta_boundary(f_next, f_prev, lid_bc[i], u); + } + } + + auto duration = std::chrono::duration_cast>( + std::chrono::high_resolution_clock::now() - start); + + std::cout << "MLUPS : " << nStep*${geometry.volume}/(1e6*duration.count()) << std::endl; + + collect_moments_to_vtk("test.vtk", f_next); +} + +int main() { + simulate(10000); +} diff --git a/lid_driven_cavity/cpp/ldc.cpp.mako b/lid_driven_cavity/cpp/ldc.cpp.mako deleted file mode 100644 index c66679f..0000000 --- a/lid_driven_cavity/cpp/ldc.cpp.mako +++ /dev/null @@ -1,178 +0,0 @@ -#include -#include -#include -#include -#include -#include - -#include "kernel.h" - -void collect_moments_to_vtk(const std::string& path, ${float_type}* f) { - std::ofstream fout; - fout.open(path.c_str()); - - fout << "# vtk DataFile Version 3.0\n"; - fout << "lbm_output\n"; - fout << "ASCII\n"; - fout << "DATASET RECTILINEAR_GRID\n"; -% if descriptor.d == 2: - fout << "DIMENSIONS " << ${geometry.size_x-2} << " " << ${geometry.size_y-2} << " 1" << "\n"; -% else: - fout << "DIMENSIONS " << ${geometry.size_x-2} << " " << ${geometry.size_y-2} << " " << ${geometry.size_z-2} << "\n"; -% endif - - fout << "X_COORDINATES " << ${geometry.size_x-2} << " float\n"; - for( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { - fout << x << " "; - } - - fout << "\nY_COORDINATES " << ${geometry.size_y-2} << " float\n"; - for( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { - fout << y << " "; - } - -% if descriptor.d == 2: - fout << "\nZ_COORDINATES " << 1 << " float\n"; - fout << 0 << "\n"; - fout << "POINT_DATA " << ${(geometry.size_x-2) * (geometry.size_y-2)} << "\n"; -% else: - fout << "\nZ_COORDINATES " << ${geometry.size_z-2} << " float\n"; - for( std::size_t z = 1; z < ${geometry.size_z-1}; ++z ) { - fout << z << " "; - } - fout << "\nPOINT_DATA " << ${(geometry.size_x-2) * (geometry.size_y-2) * (geometry.size_z-2)} << "\n"; -% endif - - ${float_type} rho; - ${float_type} u[${descriptor.d}]; - - fout << "VECTORS velocity float\n"; -% if descriptor.d == 2: - for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { - for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { - collect_moments(f, x*${geometry.size_y}+y, rho, u); - fout << u[0] << " " << u[1] << " 0\n"; - } - } -% else: - for ( std::size_t z = 1; z < ${geometry.size_z-1}; ++z ) { - for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { - for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { - collect_moments(f, x*${geometry.size_y*geometry.size_z}+y*${geometry.size_z}+z, rho, u); - fout << u[0] << " " << u[1] << " " << u[2] << "\n"; - } - } - } -% endif - - fout << "SCALARS density float 1\n"; - fout << "LOOKUP_TABLE default\n"; -% if descriptor.d == 2: - for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { - for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { - collect_moments(f, x*${geometry.size_y}+y, rho, u); - fout << rho << "\n"; - } - } -% else: - for ( std::size_t z = 1; z < ${geometry.size_z-1}; ++z ) { - for ( std::size_t y = 1; y < ${geometry.size_y-1}; ++y ) { - for ( std::size_t x = 1; x < ${geometry.size_x-1}; ++x ) { - collect_moments(f, x*${geometry.size_y*geometry.size_z}+y*${geometry.size_z}+z, rho, u); - fout << rho << "\n"; - } - } - } -% endif - - fout.close(); -} - -void simulate(std::size_t nStep) -{ - auto f_a = std::make_unique<${float_type}[]>(${geometry.volume*descriptor.q}); - auto f_b = std::make_unique<${float_type}[]>(${geometry.volume*descriptor.q}); - - ${float_type}* f_prev = f_a.get(); - ${float_type}* f_next = f_b.get(); - - std::vector bulk; - std::vector lid_bc; - std::vector box_bc; - - for (int iX = 1; iX < ${geometry.size_x-1}; ++iX) { - for (int iY = 1; iY < ${geometry.size_y-1}; ++iY) { -% if descriptor.d == 2: - const std::size_t iCell = iX*${geometry.size_y} + iY; - if (iY == ${geometry.size_y-2}) { - lid_bc.emplace_back(iCell); - } else if (iX == 1 || iX == ${geometry.size_x-2} || iY == 1) { - box_bc.emplace_back(iCell); - } else { - bulk.emplace_back(iCell); - } -% elif descriptor.d == 3: - for (int iZ = 0; iZ < ${geometry.size_z}; ++iZ) { - const std::size_t iCell = iX*${geometry.size_y*geometry.size_z} + iY*${geometry.size_z} + iZ; - if (iZ == ${geometry.size_z-2}) { - lid_bc.emplace_back(iCell); - } else if (iX == 1 || iX == ${geometry.size_x-2} || iY == 1 || iY == ${geometry.size_y-2} || iZ == 1) { - box_bc.emplace_back(iCell); - } else { - bulk.emplace_back(iCell); - } - } -% endif - } - } - - std::cout << "#bulk : " << bulk.size() << std::endl; - std::cout << "#lid : " << lid_bc.size() << std::endl; - std::cout << "#wall : " << box_bc.size() << std::endl; - std::cout << "#steps : " << nStep << std::endl; - std::cout << std::endl; - -#pragma omp parallel for - for (std::size_t iCell = 0; iCell < ${geometry.volume}; ++iCell) { - equilibrilize(f_prev, f_next, iCell); - equilibrilize(f_next, f_prev, iCell); - } - - const auto start = std::chrono::high_resolution_clock::now(); - - for (std::size_t iStep = 0; iStep < nStep; ++iStep) { - if (iStep % 2 == 0) { - f_next = f_a.get(); - f_prev = f_b.get(); - } else { - f_next = f_b.get(); - f_prev = f_a.get(); - } - -#pragma omp parallel for - for (std::size_t i = 0; i < bulk.size(); ++i) { - collide_and_stream(f_next, f_prev, bulk[i]); - } - ${float_type} u[${descriptor.d}] { 0. }; -#pragma omp parallel for - for (std::size_t i = 0; i < box_bc.size(); ++i) { - velocity_momenta_boundary(f_next, f_prev, box_bc[i], u); - } - u[0] = 0.05; -#pragma omp parallel for - for (std::size_t i = 0; i < lid_bc.size(); ++i) { - velocity_momenta_boundary(f_next, f_prev, lid_bc[i], u); - } - } - - auto duration = std::chrono::duration_cast>( - std::chrono::high_resolution_clock::now() - start); - - std::cout << "MLUPS : " << nStep*${geometry.volume}/(1e6*duration.count()) << std::endl; - - collect_moments_to_vtk("test.vtk", f_next); -} - -int main() { - simulate(10000); -} -- cgit v1.2.3