summaryrefslogtreecommitdiff
path: root/examples/cavity_2d
diff options
context:
space:
mode:
Diffstat (limited to 'examples/cavity_2d')
-rw-r--r--examples/cavity_2d/.nix/derivation.nix36
-rw-r--r--examples/cavity_2d/SConscript41
-rw-r--r--examples/cavity_2d/SConstruct80
-rw-r--r--examples/cavity_2d/cavity_2d.cpp436
4 files changed, 0 insertions, 593 deletions
diff --git a/examples/cavity_2d/.nix/derivation.nix b/examples/cavity_2d/.nix/derivation.nix
deleted file mode 100644
index a98795b..0000000
--- a/examples/cavity_2d/.nix/derivation.nix
+++ /dev/null
@@ -1,36 +0,0 @@
-{ lib
-, stdenv
-, scons
-, clang-tools
-, forstio
-, pname
-, version
-}:
-
-stdenv.mkDerivation {
- inherit pname version;
- src = ./..;
-
- nativeBuildInputs = [
- scons
- clang-tools
- ];
-
- buildInputs = [
- forstio.core
- forstio.async
- forstio.codec
- forstio.codec-unit
- forstio.codec-json
- ];
-
- doCheck = true;
- checkPhase = ''
- scons test
- ./bin/tests
- '';
-
- preferLocalBuild = true;
-
- outputs = [ "out" "dev" ];
-}
diff --git a/examples/cavity_2d/SConscript b/examples/cavity_2d/SConscript
deleted file mode 100644
index 077fb99..0000000
--- a/examples/cavity_2d/SConscript
+++ /dev/null
@@ -1,41 +0,0 @@
-#!/bin/false
-
-import os
-import os.path
-import glob
-
-
-Import('env')
-
-dir_path = Dir('.').abspath
-
-# Environment for base library
-cavity2d_env = examples_env.Clone();
-
-examples_env.sources = sorted(glob.glob(dir_path + "/*.cpp"))
-examples_env.headers = sorted(glob.glob(dir_path + "/*.hpp"))
-
-env.sources += examples_env.sources;
-env.headers += examples_env.headers;
-
-# Cavity2D
-examples_objects = [];
-examples_env.add_source_files(examples_objects, ['cavity_2d.cpp'], shared=False);
-examples_env.cavity_2d = examples_env.Program('#bin/cavity_2d', [env.library_static, examples_objects]);
-
-# Set Alias
-env.examples = [
- examples_env.meta_2d,
-# examples_env.cavity_2d,
-# examples_env.cavity_3d,
-# examples_env.particle_ibm_2d
- examples_env.poiseulle_2d,
- examples_env.poiseulle_channel_2d,
- examples_env.poiseulle_particles_channel_2d
-];
-env.Alias('examples', env.examples);
-
-if env["build_examples"]:
- env.targets += ['examples'];
- env.Install('$prefix/bin/', env.examples);
-#endif
diff --git a/examples/cavity_2d/SConstruct b/examples/cavity_2d/SConstruct
deleted file mode 100644
index fc60882..0000000
--- a/examples/cavity_2d/SConstruct
+++ /dev/null
@@ -1,80 +0,0 @@
-#!/usr/bin/env python3
-
-import sys
-import os
-import os.path
-import glob
-import re
-
-
-if sys.version_info < (3,):
- def isbasestring(s):
- return isinstance(s,basestring)
-else:
- def isbasestring(s):
- return isinstance(s, (str,bytes))
-
-def add_kel_source_files(self, sources, filetype, lib_env=None, shared=False, target_post=""):
-
- if isbasestring(filetype):
- dir_path = self.Dir('.').abspath
- filetype = sorted(glob.glob(dir_path+"/"+filetype))
-
- for path in filetype:
- target_name = re.sub( r'(.*?)(\.cpp|\.c\+\+)', r'\1' + target_post, path )
- if shared:
- target_name+='.os'
- sources.append( self.SharedObject( target=target_name, source=path ) )
- else:
- target_name+='.o'
- sources.append( self.StaticObject( target=target_name, source=path ) )
- pass
-
-def isAbsolutePath(key, dirname, env):
- assert os.path.isabs(dirname), "%r must have absolute path syntax" % (key,)
-
-env_vars = Variables(
- args=ARGUMENTS
-)
-
-env_vars.Add('prefix',
- help='Installation target location of build results and headers',
- default='/usr/local/',
- validator=isAbsolutePath
-)
-
-env_vars.Add('build_examples',
- help='If examples should be built',
- default="true"
-)
-
-env=Environment(ENV=os.environ, variables=env_vars, CPPPATH=[],
- CPPDEFINES=['SAW_UNIX'],
- CXXFLAGS=[
- '-std=c++20',
- '-g',
- '-Wall',
- '-Wextra'
- ],
- LIBS=[
- 'forstio-core',
- 'forstio-codec'
- ]
-);
-env.__class__.add_source_files = add_kel_source_files
-env.Tool('compilation_db');
-env.cdb = env.CompilationDatabase('compile_commands.json');
-
-env.objects = [];
-env.sources = [];
-env.headers = [];
-env.targets = [];
-
-Export('env')
-SConscript('SConscript')
-
-env.Alias('cdb', env.cdb);
-env.Alias('all', [env.targets]);
-env.Default('all');
-
-env.Alias('install', '$prefix')
diff --git a/examples/cavity_2d/cavity_2d.cpp b/examples/cavity_2d/cavity_2d.cpp
deleted file mode 100644
index e228a06..0000000
--- a/examples/cavity_2d/cavity_2d.cpp
+++ /dev/null
@@ -1,436 +0,0 @@
-#include "../c++/descriptor.hpp"
-#include "../c++/macroscopic.hpp"
-#include "../c++/lbm.hpp"
-#include "../c++/component.hpp"
-#include "../c++/collision.hpp"
-#include "../c++/boundary.hpp"
-
-#include <forstio/codec/data.hpp>
-
-#include <iostream>
-#include <fstream>
-#include <cmath>
-
-namespace kel {
-namespace lbm {
-namespace sch {
-using namespace saw::schema;
-
-/**
- * Basic distribution function
- * Base type
- * D
- * Q
- * Scalar factor
- * D factor
- * Q factor
- */
-using T = Float32;
-using D2Q5 = Descriptor<2u,5u>;
-using D2Q9 = Descriptor<2u,9u>;
-
-template<typename Desc>
-using DfCell = Cell<T, Desc, 0u, 0u, 1u>;
-
-template<typename Desc>
-using CellInfo = Cell<UInt8, D2Q9, 1u, 0u, 0u>;
-
-/**
- * Basic type for simulation
- */
-template<typename Desc>
-using CellStruct = Struct<
- Member<DfCell<Desc>, "dfs">,
- Member<DfCell<Desc>, "dfs_old">,
- Member<CellInfo<Desc>, "info">
->;
-
-template<typename T, uint64_t D>
-using MacroStruct = Struct<
- Member<FixedArray<T,D>, "velocity">,
- Member<T, "pressure">
->;
-
-using CavityFieldD2Q9 = CellField<D2Q9, CellStruct<D2Q9>>;
-}
-
-
-/*
-template<typename T, typename Encode>
-class df_cell_view;
-*/
-/**
- * Minor helper for the AA-Pull Pattern, so I can use only one lattice
- *
- * Am I sure I want to use AA this way?
- * Esoteric Twist technically reduces the needed memory access footprint
- */
-/*
-template<typename Desc, size_t SN, size_t DN, size_t QN, typename Encode>
-class df_cell_view<sch::Cell<sch::T, Desc, SN, DN, QN>, Encode> {
-public:
- using Schema = sch::Cell<sch::T,Desc,SN,DN,QN>;
-private:
- std::array<std::decay_t<typename saw::native_data_type<sch::T>::type>*, QN> view_;
-public:
- df_cell_view(const std::array<std::decay_t<typename saw::native_data_type<sch::T>::type>*, QN>& view):
- view_{view}
- {}
-};
-*/
-namespace cmpt {
-struct MovingWall {};
-}
-
-/**
- * Full-Way moving wall Bounce back, something is not right here.
- * Technically it should reflect properly.
- */
-template<typename Desc>
-class component<sch::T, Desc, cmpt::MovingWall> {
-public:
- std::array<typename saw::native_data_type<sch::T>::type, Desc::D> lid_vel;
-
-public:
- void apply(
- saw::data<sch::DfCell<Desc>>& dfs
- ){
- using dfi = df_info<sch::T,Desc>;
-
- // Technically use .copy()
- /*
- auto dfs_cpy = dfs;
-
- for(uint64_t i = 0u; i < Desc::Q; ++i){
- dfs({dfi::opposite_index.at(i)}) = dfs_cpy({i}) - 2.0 * dfi::weights[i] * 1.0 * ( lid_vel[0] * dfi::directions[i][0] + lid_vel[1] * dfi::directions[i][1]) * dfi::inv_cs2;
- }
- */
- }
-};
-}
-}
-
-constexpr size_t dim_size = 2;
-constexpr size_t dim_x = 128;
-constexpr size_t dim_y = 128;
-
-
-template<typename Func>
-void apply_for_cells(Func&& func, saw::data<kel::lbm::sch::CavityFieldD2Q9>& dat){
- for(std::size_t i = 0; i < dat.meta().at({1u}).get(); ++i){
- for(std::size_t j = 0; j < dat.meta().at({0u}).get(); ++j){
- saw::data<saw::schema::UInt64> di{i};
- saw::data<saw::schema::UInt64> dj{j};
- auto& cell_v = dat({{dj,di}});
- func(cell_v, j, i);
- }
- }
-}
-
-void set_geometry(saw::data<kel::lbm::sch::CavityFieldD2Q9>& latt){
- using namespace kel::lbm;
- apply_for_cells([](auto& cell, std::size_t i, std::size_t j){
- uint8_t val = 0;
- if(j == 1){
- val = 2u;
- }
- if(i == 1 || (i+2) == dim_x || (j+2) == dim_y){
- val = 3u;
- }
- if(i > 1 && (i+2) < dim_x && j > 1 && (j+2) < dim_y){
- val = 1u;
- }
- if(i == 0 || j == 0 || (i+1) == dim_x || (j+1) == dim_y){
- val = 0u;
- }
- cell.template get<"info">()(0u).set(val);
- }, latt);
-}
-
-void set_initial_conditions(saw::data<kel::lbm::sch::CavityFieldD2Q9>& latt){
- using namespace kel::lbm;
-
- saw::data<sch::T> rho{1.0};
- {
- saw::data<sch::FixedArray<sch::T,sch::D2Q9::D>> vel{{0.0,0.0}};
- auto eq = equilibrium<sch::T,sch::D2Q9>(rho, vel);
-
- apply_for_cells([&eq](auto& cell, std::size_t i, std::size_t j){
- (void) i;
- (void) j;
- auto& dfs = cell.template get<"dfs">();
- auto& dfs_old = cell.template get<"dfs_old">();
- auto info = cell.template get<"info">()(0u).get();
- for(uint64_t k = 0; k < sch::D2Q9::Q; ++k){
- dfs(k) = eq.at({k});
- dfs_old(k) = eq.at({k});
- }
- }, latt);
- }
-
- {
- saw::data<sch::FixedArray<sch::T,sch::D2Q9::D>> vel{{0.1,0.0}};
- auto eq = equilibrium<sch::T,sch::D2Q9>(rho, vel);
-
- apply_for_cells([&eq](auto& cell, std::size_t i, std::size_t j){
- (void) i;
- (void) j;
- auto& dfs = cell.template get<"dfs">();
- auto& dfs_old = cell.template get<"dfs_old">();
- auto info = cell.template get<"info">()(0u).get();
- if(info == 2u){
- for(uint64_t k = 0; k < sch::D2Q9::Q; ++k){
- dfs(k) = eq.at({k});
- dfs_old(k) = eq.at({k});
- }
- }
- }, latt);
- }
-}
-
-void lbm_step(
- saw::data<kel::lbm::sch::CavityFieldD2Q9>& latt,
- bool even_step,
- uint64_t time_step
-){
- using namespace kel::lbm;
- using dfi = df_info<sch::T,sch::D2Q9>;
-
- /**
- * 1. Relaxation parameter \tau
- */
- component<sch::T, sch::D2Q9, cmpt::BGK> coll{0.5384};
- component<sch::T, sch::D2Q9, cmpt::BounceBack> bb;
-
- component<sch::T, sch::D2Q9, cmpt::MovingWall> bb_lid;
- bb_lid.lid_vel = {0.1,0.0};
-
- auto meta = latt.meta();
-
- // Collide
- for(saw::data<sch::UInt64> i{0u}; i < meta.at(0u); ++i){
- for(saw::data<sch::UInt64> j{0u}; j < meta.at(1u); ++j ){
- auto& cell = latt({{i,j}});
- auto& info = cell.template get<"info">();
-
- switch(info({0u}).get()){
- case 1u: {
- coll.apply(latt, {{i,j}}, time_step);
- break;
- }
- case 2u: {
- auto& df = even_step ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
- bb_lid.apply(df);
- break;
- }
- case 3u: {
- bb.apply(latt, {{i,j}}, time_step);
- break;
- }
- default:
- break;
- }
- }
- }
-
- // Stream
- for(uint64_t i = 1u; (i+1u) < latt.template get_dim_size<0>().get(); ++i){
- for(uint64_t j = 1u; (j+1u) < latt.template get_dim_size<1>().get(); ++j){
- auto& cell = latt({{i,j}});
- auto& df_new = even_step ? cell.template get<"dfs">() : cell.template get<"dfs_old">();
- auto& info_new = cell.template get<"info">();
-
- if(info_new({0u}).get() > 0u && info_new({0u}).get() != 2u){
- for(uint64_t k = 0u; k < sch::D2Q9::Q; ++k){
- auto dir = dfi::directions[dfi::opposite_index[k]];
- auto& cell_dir_old = latt({{i+dir[0],j+dir[1]}});
-
- auto& df_old = even_step ? cell_dir_old.template get<"dfs_old">() : cell_dir_old.template get<"dfs">();
- auto& info_old = cell_dir_old.template get<"info">();
-
- if( info_old({0}).get() == 2u ){
- auto& df_old_loc = even_step ? latt({{i,j}}).template get<"dfs_old">() : latt({{i,j}}).template get<"dfs">();
- df_new({k}) = df_old_loc({dfi::opposite_index.at(k)}) - 2.0 * dfi::inv_cs2 * dfi::weights.at(k) * 1.0 * ( bb_lid.lid_vel[0] * dir[0] + bb_lid.lid_vel[1] * dir[1]);
- // dfs({dfi::opposite_index.at(i)}) = dfs_cpy({i}) - 2.0 * dfi::weights[i] * 1.0 * ( lid_vel[0] * dfi::directions[i][0] + lid_vel[1] * dfi::directions[i][1]) * dfi::inv_cs2;
- } else {
- df_new({k}) = df_old({k});
- }
- }
- }
- }
- }
-}
-
-int main(){
- using namespace kel::lbm;
-
- saw::data<sch::FixedArray<sch::UInt64,sch::D2Q9::D>> dim{{dim_x, dim_y}};
-
- saw::data<sch::CavityFieldD2Q9, saw::encode::Native> lattice{dim};
-
- converter<sch::T> conv{
- {0.1},
- {0.1}
- };
-
- print_lbm_meta<sch::T, sch::D2Q9>(conv, {1e-3});
-
- //auto& df_field = lattices.at(0).template get<"dfs">();
- //for(uint64_t i = 0; i < df_field.get_dim_size<0u>(); ++i){
- // lattices.at(i) = {dim_x, dim_y};
- //}
-
- /**
- * Set meta information describing what this cell is
- */
- set_geometry(lattice);
-
- /**
- *
- */
- set_initial_conditions(lattice);
-
- /**
- * Timeloop
- */
-
- /**
- * Print basic setup info
- */
- apply_for_cells([](auto& cell, std::size_t i, std::size_t j){
- // Not needed
- (void) j;
- std::cout<<(static_cast<uint32_t>(cell.template get<"info">()({0}).get()));
- if( (i+1) < dim_x){
- std::cout<<" ";
- }else{
- std::cout<<"\n";
- }
- }, lattice);
-
- uint64_t lattice_steps = 512000u;
- bool even_step = true;
-
- uint64_t print_every = 256u;
- uint64_t file_no = 0u;
-
- saw::data<sch::Array<sch::MacroStruct<sch::T,sch::D2Q9::D>,sch::D2Q9::D>> macros{dim};
-
- for(uint64_t step = 0; step < lattice_steps; ++step){
-
- if(step % print_every == 0u){
- std::cout<<"\n";
- typename saw::native_data_type<sch::T>::type sum = 0.0;
-
-
- apply_for_cells([&](auto& cell, std::size_t i, std::size_t j){
- auto& dfs = cell.template get<"dfs">();
- typename saw::native_data_type<sch::T>::type rho;
- std::array<typename saw::native_data_type<sch::T>::type, sch::D2Q9::D> vel;
- compute_rho_u<sch::T,sch::D2Q9>(dfs,rho,vel);
-
- if(i > 1 && (i+2) < dim_x && j > 1 && (j+2) < dim_y){
- sum += rho;
- }
- auto angle = std::atan2(vel[1],vel[0]);
- double vel_mag = vel[1] * vel[1] + vel[0] * vel[0];
-
- auto dir_char = [&]() -> std::string_view {
- /**
- * Flipped y due to print order
- */
- constexpr auto pi = M_PI;
- return "■";
- if(vel_mag < 1e-4){
- return "•";
- }
- if(angle > 7.0 / 8.0 * pi){
- return "←";
- }
- if(angle > 5.0 / 8.0 * pi){
- return "↙";
- }
- if(angle > 3.0 / 8.0 * pi){
- return "↓";
- }
- if(angle > 1.0 / 8.0 * pi){
- return "↘";
- }
- if(angle > -1.0 / 8.0 * pi){
- return "→";
- }
- if(angle > -3.0 / 8.0 * pi){
- return "↗";
- }
- if(angle > -5.0 / 8.0 * pi){
- return "↑";
- }
- if(angle > -7.0 / 8.0 * pi){
- return "↖";
- }
- return "←";
- }();
-
- std::array<uint32_t,3u> rgb_start{64,64,255};
- std::array<uint32_t,3u> rgb_stop{255,64,64};
- std::array<uint32_t,3u> rgb_middle{255,255,255};
-
-
- auto col_interpol = [&](){
- std::array<uint32_t, 3u> rgb_interpol = rgb_start;
- double vel_mag_cut = std::min(1.0,std::max(vel_mag/(0.07*0.07),0.0));
-
- if(vel_mag_cut < 0.5){
- uint32_t vel_mag_i = static_cast<uint32_t>(2.0 * vel_mag_cut * 256);
- for(uint8_t i = 0u; i < 3u; ++i){
- rgb_interpol[i] = (rgb_middle[i] * vel_mag_i + (256-vel_mag_i) * rgb_start[i]) / 256;
- }
- }else{
- uint32_t vel_mag_i = static_cast<uint32_t>((2.0*(vel_mag_cut-0.5)) * 256);
- for(uint8_t i = 0u; i < 3u; ++i){
- rgb_interpol[i] = (rgb_stop[i] * vel_mag_i + (256-vel_mag_i) * rgb_middle[i]) / 256;
- }
- }
- return rgb_interpol;
- };
- auto rgb_interpol = col_interpol();
-
- std::cout<<"\x1b[38;2;"<<rgb_interpol[0]<<";"<<rgb_interpol[1]<<";"<<rgb_interpol[2]<<"m"<<dir_char;
- if( (i+1) < dim_x){
- std::cout<<" ";
- }else{
- std::cout<<"\n";
- }
- }, lattice);
- std::cout<<"\x1b[38;2;255;255;255m";
-
- /// std::cout<<"Average rho: "<<(sum / ((dim_x-4)*(dim_y-4)))<<std::endl;
- std::cout.flush();
-
- apply_for_cells([&](auto& cell, std::size_t i, std::size_t j){
- auto dfs = even_step ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
-
- auto& rho = macros.at({{i,j}}).template get<"pressure">();
- auto& vel = macros.at({{i,j}}).template get<"velocity">();
- compute_rho_u<sch::T,sch::D2Q9>(dfs,rho,vel);
- }, lattice);
-
- std::string vtk_f_name{"tmp/cavity_2d_"};
- vtk_f_name += std::to_string(file_no) + ".vtk";
- write_vtk_file(vtk_f_name, macros);
-
- ++file_no;
- }
-
- lbm_step(lattice, even_step, step);
-
- even_step = not even_step;
- }
-
- /**
- * Flush cout
- */
- std::cout<<"\n\n";
- std::cout.flush();
- return 0;
-}
generated by cgit v1.2.3 (git 2.25.1) at 2025-09-17 17:40:57 +0000