From 7a609b2c764f51548ce88dedcdd4c2575a6ae200 Mon Sep 17 00:00:00 2001 From: Frank14f <1515444314@qq.com> Date: Sat, 20 Jun 2026 18:48:23 +0800 Subject: [PATCH] chore: remove legacy, sync generated config headers, fix inlet kernel - Remove legacy/ directory (superseded by current architecture). - Sync auto-generated config headers (config_grid.h, config_objects.h, config_method.h, config_physics.h) for current LBMConfig defaults. - Sync zou_he_local.cuh inlet kernel changes. Co-authored-by: Cursor --- legacy/README.md | 19 - legacy/common_utils.py | 364 -------------- legacy/cuda_compiler_v1.py | 132 ------ legacy/lbm_configs/config_cuda.json | 9 - legacy/lbm_configs/config_flowfield.json | 13 - legacy/lbm_driver.py | 445 ------------------ legacy/macros.h | 108 ----- .../kernels/boundary/inlet/zou_he_local.cuh | 5 +- .../lbm/kernels/config/config_grid.h | 4 +- .../lbm/kernels/config/config_method.h | 8 +- .../lbm/kernels/config/config_objects.h | 2 +- .../lbm/kernels/config/config_physics.h | 2 +- 12 files changed, 10 insertions(+), 1101 deletions(-) delete mode 100644 legacy/README.md delete mode 100644 legacy/common_utils.py delete mode 100644 legacy/cuda_compiler_v1.py delete mode 100644 legacy/lbm_configs/config_cuda.json delete mode 100644 legacy/lbm_configs/config_flowfield.json delete mode 100644 legacy/lbm_driver.py delete mode 100644 legacy/macros.h diff --git a/legacy/README.md b/legacy/README.md deleted file mode 100644 index 6fbca5d..0000000 --- a/legacy/README.md +++ /dev/null @@ -1,19 +0,0 @@ -# Legacy Code Archive - -This directory contains code that has been superseded by the current architecture but is kept for reference. - -## Contents - -| File / Dir | Replaced By | Reason | -|---|---|---| -| `lbm_driver.py` | `src/CelerisLab/simulation.py` + `lbm/field.py` + `lbm/stepper.py` | Monolithic FlowField class. New Simulation API separates concerns: CudaContext / LBMField / LBMStepper / ObjectManager. | -| `cuda_compiler_v1.py` | `src/CelerisLab/cuda/compiler_v2.py` | macros.h-based build system. New compiler writes typed config/*.h headers per architectural layer. | -| `macros.h` | `src/CelerisLab/lbm/kernels/config/*.h` | Single flat macro file. Now split into config_grid.h / config_physics.h / config_method.h / config_objects.h matching the Global/Method/Case/Debug parameter hierarchy. | -| `common_utils.py` | `src/CelerisLab/config.py` + `src/CelerisLab/cuda/context.py` | FlowFieldConfig / CudaConfig NamedTuples and their JSON loaders. Replaced by LBMConfig / BodyConfig dataclasses (config.py) and CudaContext (cuda/context.py). | -| `lbm_configs/` | `src/CelerisLab/configs/` | Old JSON config format used by FlowField / compiler_v1. | - -## Notes - -- None of these files is imported by any active module. -- `lbm_driver.py` (FlowField) depended on `cuda_compiler_v1.py` and `common_utils.py`; all three were removed from src together. -- `macros.h` was the old single-file configuration for `kernel_v2.cu`; kernel_v2.cu now includes `config.h` which aggregates `config/*.h`. diff --git a/legacy/common_utils.py b/legacy/common_utils.py deleted file mode 100644 index a087da7..0000000 --- a/legacy/common_utils.py +++ /dev/null @@ -1,364 +0,0 @@ -# CelerisLab/utils.py - -import pycuda.driver as cuda -import subprocess -import json -import os - -from typing import NamedTuple, Optional, List, Tuple, Union - - -class CudaDeviceInfo(NamedTuple): - name: str - compute_capability: str - multiprocessors: int - total_global_memory: int - max_shared_memory_per_block: int - max_threads_per_block: int - max_blocks_per_multiprocessor: int - device_interconnect: Optional[str] = None - - -class FlowFieldConfig(NamedTuple): - data_type: str - dimensionality: int - lattice: int - field_dim_in_U: Tuple[int, int, int] - viscosity: float - velocity: float - boundary_conditions: Tuple[str, str, str, str, str, str] - - -class CudaConfig(NamedTuple): - multi_gpu: bool - gpu_connection: str - required_cuda_capability: str - threads_per_block: int - unit_dimensions: Tuple[int, int, int] - - -def check_cuda_device_availability(device_id=0): - if cuda.Device.count() == 0: - raise RuntimeError("No CUDA device is available.") - - if device_id < 0 or device_id >= cuda.Device.count(): - raise ValueError( - f"Invalid device_id {device_id}. Must be between 0 and {cuda.Device.count() - 1}." - ) - - try: - subprocess.check_output(["nvidia-smi", "--version"]) - except subprocess.CalledProcessError: - raise RuntimeError("nvidia-smi is not available or not installed correctly.") - - -def query_cuda_device_info(device_id=0) -> CudaDeviceInfo: - check_cuda_device_availability(device_id) - - try: - output = subprocess.check_output( - ["nvidia-smi", "-q", "-d", "TOPOLOGY", "-i", str(device_id)], text=True - ) - if "NVLink" in output: - device_interconnect = "NVLink" - elif "PCIe" in output: - device_interconnect = "PCIe" - else: - device_interconnect = "Unknown" - except Exception as e: - device_interconnect = None - - device = cuda.Device(device_id) - - return CudaDeviceInfo( - name=device.name(), - compute_capability=f"{device.compute_capability()[0]}.{device.compute_capability()[1]}", - multiprocessors=device.get_attribute( - cuda.device_attribute.MULTIPROCESSOR_COUNT - ), - total_global_memory=device.total_memory(), - max_shared_memory_per_block=device.get_attribute( - cuda.device_attribute.MAX_SHARED_MEMORY_PER_BLOCK - ), - max_threads_per_block=device.get_attribute( - cuda.device_attribute.MAX_THREADS_PER_BLOCK - ), - max_blocks_per_multiprocessor=device.get_attribute( - cuda.device_attribute.MAX_BLOCKS_PER_MULTIPROCESSOR - ), - device_interconnect=device_interconnect, - ) - - -def find_config_file(config_filename: str, config_path: Optional[str] = None) -> str: - """ - Find configuration file by searching in multiple locations. - - Search priority: - 1. Provided config_path (if given) - 2. Environment variable CELERISLAB_CONFIG_DIR - 3. Current working directory ./configs/ - 4. Package installation location (relative to this utils.py file) - - Args: - config_filename: Name of the config file (e.g., 'config_cuda.json') - config_path: Optional explicit path to config file - - Returns: - Absolute path to the config file - - Raises: - FileNotFoundError: If config file cannot be found in any location - """ - search_paths = [] - - # Priority 1: Explicit path provided - if config_path: - search_paths.append(config_path) - - # Priority 2: Environment variable - env_config_dir = os.environ.get('CELERISLAB_CONFIG_DIR') - if env_config_dir: - search_paths.append(os.path.join(env_config_dir, config_filename)) - - # Priority 3: Current working directory - search_paths.append(os.path.join(os.getcwd(), 'configs', config_filename)) - - # Priority 4: Package installation location (relative to this utils.py) - # configs are in lbm/configs/ - package_root = os.path.dirname(os.path.dirname(os.path.abspath(__file__))) - search_paths.append(os.path.join(package_root, 'lbm', 'configs', config_filename)) - - # Search for the file - for path in search_paths: - if os.path.isfile(path): - return os.path.abspath(path) - - # File not found, provide helpful error message - error_msg = f"Configuration file '{config_filename}' not found. Searched in:\n" - for path in search_paths: - error_msg += f" - {path}\n" - error_msg += "\nTo fix this, you can:\n" - error_msg += " 1. Set CELERISLAB_CONFIG_DIR environment variable\n" - error_msg += " 2. Place config files in ./configs/ directory\n" - error_msg += " 3. Provide explicit config_path parameter" - raise FileNotFoundError(error_msg) - - -def load_flow_field_config(config_path: Optional[str] = None) -> FlowFieldConfig: - """ - Load flow field configuration from JSON file. - - Args: - config_path: Optional path to config file. If None, searches in standard locations. - Can be relative path like 'configs/config_flowfield.json' or just filename. - - Returns: - FlowFieldConfig object - """ - # Determine config filename and full path - if config_path: - # Check if it's just a filename or a path - if os.path.basename(config_path) == config_path: - # Just a filename, search for it - config_file = find_config_file(config_path, None) - else: - # It's a path, use it if exists, otherwise try to find the basename - if os.path.isfile(config_path): - config_file = config_path - else: - config_file = find_config_file(os.path.basename(config_path), None) - else: - # No path provided, search for default filename - config_file = find_config_file('config_flowfield.json', None) - - try: - with open(config_file, "r") as file: - config = json.load(file) - - required_keys = [ - "data_type", - "dimensionality", - "lattice", - "field_dim_in_U", - "viscosity", - "boundary_conditions", - ] - if not all(key in config for key in required_keys): - raise ValueError("Missing required configuration items.") - - if config["data_type"] not in ["FP32", "FP64"]: - raise ValueError("Data type must be either FP32 or FP64.") - - if config["dimensionality"] not in [2, 3]: - raise ValueError("Dimensionality must be either 2 or 3.") - - if config["dimensionality"] == 2 and config["field_dim_in_U"][2] != 1: - raise ValueError( - "Field dimensions must be 1 in the third dimension for 2D simulations." - ) - - if config["lattice"] not in [9]: - raise ValueError("Lattice must be either 9 or 19.") - - boundary_conditions = tuple( - condition - for key in ["x", "y", "z"] - for condition in config["boundary_conditions"].get(key, []) - ) - if len(boundary_conditions) != 6: - raise ValueError("Boundary conditions must contain exactly six elements.") - - return FlowFieldConfig( - data_type=config["data_type"], - dimensionality=config["dimensionality"], - lattice=config["lattice"], - field_dim_in_U=tuple(config["field_dim_in_U"]), - viscosity=config["viscosity"], - velocity=config["velocity"], - boundary_conditions=boundary_conditions, - ) - except Exception as e: - raise RuntimeError(f"Failed to load or parse the flow field configuration: {e}") - - -def load_cuda_config(config_path: Optional[str] = None) -> CudaConfig: - """ - Load CUDA configuration from JSON file. - - Args: - config_path: Optional path to config file. If None, searches in standard locations. - Can be relative path like 'configs/config_cuda.json' or just filename. - - Returns: - CudaConfig object - """ - # Determine config filename and full path - if config_path: - # Check if it's just a filename or a path - if os.path.basename(config_path) == config_path: - # Just a filename, search for it - config_file = find_config_file(config_path, None) - else: - # It's a path, use it if exists, otherwise try to find the basename - if os.path.isfile(config_path): - config_file = config_path - else: - config_file = find_config_file(os.path.basename(config_path), None) - else: - # No path provided, search for default filename - config_file = find_config_file('config_cuda.json', None) - - try: - with open(config_file, "r") as file: - config = json.load(file) - - required_keys = [ - "multi_gpu", - "gpu_connection", - "required_cuda_capability", - "threads_per_block", - "X_1U", - "Y_1U", - "Z_1U", - ] - - if not all(key in config for key in required_keys): - raise ValueError("Missing required configuration items.") - - return CudaConfig( - multi_gpu=config["multi_gpu"], - gpu_connection=config["gpu_connection"], - required_cuda_capability=config["required_cuda_capability"], - threads_per_block=config["threads_per_block"], - unit_dimensions=(config["X_1U"], config["Y_1U"], config["Z_1U"]), - ) - except Exception as e: - raise RuntimeError(f"Failed to load or parse the CUDA configuration: {e}") - - -def check_cuda_capability( - field_config: FlowFieldConfig, - cuda_config: CudaConfig, - device_id: Union[int, List[int]] = None, -): - SAFE_FACTOR = 0.8 - - if cuda_config.multi_gpu: - if device_id is None or isinstance(device_id, int): - raise ValueError("Multi-GPU support requires a list of device IDs.") - raise NotImplementedError("Multi-GPU support is not implemented yet.") - else: - if isinstance(device_id, list): - if len(device_id) > 1: - raise ValueError( - "Single-GPU mode does not support multiple device IDs." - ) - device_id = device_id[0] - elif device_id is None: - device_id = 0 - device_info = query_cuda_device_info(device_id) - - if device_info.compute_capability != cuda_config.required_cuda_capability: - raise ValueError( - f"Device {device_info.name} has compute capability {device_info.compute_capability}, but {cuda_config.required_cuda_capability} is required." - ) - - field_size = sum( - size * unit - for size, unit in zip( - field_config.field_dim_in_U, cuda_config.unit_dimensions - ) - ) - if ( - device_info.total_global_memory * SAFE_FACTOR - < calc_field_memory_consumption( - field_size, - field_config.dimensionality, - field_config.lattice, - field_config.data_type, - ) - ): - raise ValueError( - f"Device {device_info.name} does not have enough memory to store the flow field." - ) - - if ( - device_info.max_threads_per_block * SAFE_FACTOR - < cuda_config.threads_per_block - ): - raise ValueError( - f"Device {device_info.name} does not have enough threads per block to run the simulation." - ) - - block_size = cuda_config.threads_per_block - if ( - device_info.max_shared_memory_per_block * SAFE_FACTOR - < 2 - * calc_field_memory_consumption( - block_size, - field_config.dimensionality, - field_config.lattice, - field_config.data_type, - ) - ): - raise ValueError( - f"Device {device_info.name} does not have enough shared memory per block to run the simulation." - ) - - -def calc_field_memory_consumption( - field_size: int, dimensionality: int, directions: int, data_type: str -) -> int: - if data_type == "FP32": - data_size = 4 - elif data_type == "FP64": - data_size = 8 - else: - raise ValueError(f"Unsupported data type {data_type}.") - - return ( - field_size * directions * data_size * 2 - + field_size * dimensionality * data_size - + field_size - ) diff --git a/legacy/cuda_compiler_v1.py b/legacy/cuda_compiler_v1.py deleted file mode 100644 index 18eb528..0000000 --- a/legacy/cuda_compiler_v1.py +++ /dev/null @@ -1,132 +0,0 @@ -# CelerisLab/cuda/compiler.py - -import subprocess -import re -import os - -from ..common.utils import FlowFieldConfig, CudaConfig - - -def kernel_path(file_name: str) -> str: - # kernels are in lbm/kernels/ - current_dir = os.path.dirname(os.path.dirname(os.path.abspath(__file__))) - return os.path.join(current_dir, "lbm", "kernels", file_name) - - -def read_lines(file_path): - with open(file_path, "r") as file: - lines = file.readlines() - return lines - - -def write_lines(file_path, lines): - with open(file_path, "w") as file: - file.writelines(lines) - - -def modify_macro(lines, macro_name, new_value): - pattern = re.compile(rf"(#define\s+{macro_name}\s+)(\S+)") - for i, line in enumerate(lines): - if pattern.match(line): - lines[i] = pattern.sub(rf"\g<1>{new_value}", line) - break - return lines - -def modify_const(lines, const_name, new_type, new_value): - pattern = re.compile(rf"(__constant__\s+)(\S+\s+{const_name}\s*=\s*)([^;]+)(;)") - for i, line in enumerate(lines): - if pattern.match(line): - lines[i] = pattern.sub(rf"\g<1>{new_type} {const_name} = {new_value}\4", line) - break - return lines - -def compile_kernel(): - subprocess.run( - [ - "nvcc", - "-ptx", - kernel_path("kernel.cu"), - "-o", - kernel_path("kernel.ptx"), - ] - ) - - -def compile_kernel_v2(): - """Compile the new modular kernel (kernel_v2.cu → kernel_v2.ptx).""" - subprocess.run( - [ - "nvcc", - "-ptx", - kernel_path("kernel_v2.cu"), - "-o", - kernel_path("kernel_v2.ptx"), - ] - ) - -def config_kernal(config_cuda: CudaConfig, config_field: FlowFieldConfig): - lines = read_lines(kernel_path("macros.h")) - lines = modify_macro(lines, "MULT_GPU", config_cuda.multi_gpu) - lines = modify_macro(lines, "NT", config_cuda.threads_per_block) - lines = modify_macro(lines, "X_1U", config_cuda.unit_dimensions[0]) - lines = modify_macro(lines, "Y_1U", config_cuda.unit_dimensions[1]) - lines = modify_macro(lines, "Z_1U", config_cuda.unit_dimensions[2]) - - if config_field.data_type == "FP32": - lb_type = "float" - else: - raise ValueError(f"Unsupported data type {config_field.data_type}.") - lines = modify_macro(lines, "LBtype", lb_type) - lines = modify_macro(lines, "UX", config_field.field_dim_in_U[0]) - lines = modify_macro(lines, "UY", config_field.field_dim_in_U[1]) - lines = modify_macro(lines, "UZ", config_field.field_dim_in_U[2]) - lines = modify_macro(lines, "NX", config_field.field_dim_in_U[0] * config_cuda.unit_dimensions[0]) - lines = modify_macro(lines, "NY", config_field.field_dim_in_U[1] * config_cuda.unit_dimensions[1]) - lines = modify_macro(lines, "NZ", config_field.field_dim_in_U[2] * config_cuda.unit_dimensions[2]) - lines = modify_macro(lines, "DIM", config_field.dimensionality) - lines = modify_macro(lines, "NQ", config_field.lattice) - lines = modify_macro(lines, "VIS", config_field.viscosity) - lines = modify_macro(lines, "U0", config_field.velocity) - - write_lines(kernel_path("macros.h"), lines) - - -def config_kernal_v2(config_cuda: CudaConfig, config_field: FlowFieldConfig, - collision_model: int = 2, - streaming_model: int = 0, - store_precision: int = 0, - use_ddf_shifting: int = 0, - use_les: int = 0, - les_cs: float = 0.16): - """Configure macros.h for the new modular kernel architecture. - - Args: - collision_model: 0=SRT, 1=TRT, 2=MRT (default) - streaming_model: 0=double-buffer (default), 1=Esoteric-Pull - store_precision: 0=FP32 (default), 1=FP16S, 2=FP16C - use_ddf_shifting: 0=off (default), 1=on - use_les: 0=off (default), 1=Smagorinsky LES - les_cs: Smagorinsky constant C_s - """ - # First apply legacy config - config_kernal(config_cuda, config_field) - - # Then apply new architecture macros - lines = read_lines(kernel_path("macros.h")) - lines = modify_macro(lines, "COLLISION_MODEL", collision_model) - lines = modify_macro(lines, "STREAMING_MODEL", streaming_model) - lines = modify_macro(lines, "STORE_PRECISION", store_precision) - lines = modify_macro(lines, "USE_DDF_SHIFTING", use_ddf_shifting) - lines = modify_macro(lines, "USE_LES", use_les) - lines = modify_macro(lines, "LES_CS", f"{les_cs:.6f}f") - write_lines(kernel_path("macros.h"), lines) - -def config_object(n_obj: int): - lines = read_lines(kernel_path("macros.h")) - lines = modify_macro(lines, "N_OBJS", n_obj) - write_lines(kernel_path("macros.h"), lines) - -def config_sensor(n_sen: int): - lines = read_lines(kernel_path("macros.h")) - lines = modify_macro(lines, "N_SENS", n_sen) - write_lines(kernel_path("macros.h"), lines) \ No newline at end of file diff --git a/legacy/lbm_configs/config_cuda.json b/legacy/lbm_configs/config_cuda.json deleted file mode 100644 index 3d4c10b..0000000 --- a/legacy/lbm_configs/config_cuda.json +++ /dev/null @@ -1,9 +0,0 @@ -{ - "multi_gpu": false, - "gpu_connection": "NVLink", - "required_cuda_capability": "7.0", - "threads_per_block": 128, - "X_1U": 128, - "Y_1U": 32, - "Z_1U": 1 -} \ No newline at end of file diff --git a/legacy/lbm_configs/config_flowfield.json b/legacy/lbm_configs/config_flowfield.json deleted file mode 100644 index 3f2f42d..0000000 --- a/legacy/lbm_configs/config_flowfield.json +++ /dev/null @@ -1,13 +0,0 @@ -{ - "data_type": "FP32", - "dimensionality": 2, - "lattice": 9, - "field_dim_in_U": [10, 16, 1], - "viscosity": 0.002, - "velocity": 0.01, - "boundary_conditions": { - "x": ["parabolic", "outflow"], - "y": ["noslip", "noslip"], - "z": ["none", "none"] - } -} \ No newline at end of file diff --git a/legacy/lbm_driver.py b/legacy/lbm_driver.py deleted file mode 100644 index 8f844b7..0000000 --- a/legacy/lbm_driver.py +++ /dev/null @@ -1,445 +0,0 @@ -# CelerisLab/lbm/driver.py - -import pycuda.driver as cuda -import numpy as np -import struct -from scipy.special import jv, expi -from typing import List, Tuple, Union, Optional - -from ..common import utils -from ..common import preprocess as preproc -from ..cuda import compiler - -FLUID = 0b00000001 -SOLID = 0b00000010 -GAS = 0b00000100 -INTERFACE = 0b00001000 -SENSOR = 0b00010000 -V_TAYLOR = np.int32(1) - -class FlowField: - def __init__( - self, - field_config: utils.FlowFieldConfig, - cuda_config: utils.CudaConfig, - device_id: Union[int, List[int]] = None, - use_kernel_v2: bool = True, - collision_model: int = 0, - streaming_model: int = 0, - store_precision: int = 0, - use_ddf_shifting: int = 0, - use_les: int = 0, - les_cs: float = 0.16, - ): - self.field_config = field_config - self.cuda_config = cuda_config - cuda.init() - - # Sanity checks - if cuda_config.multi_gpu: - if device_id is None or isinstance(device_id, int): - raise ValueError("Multi-GPU support requires a list of device IDs.") - # self.devices = [cuda.Device(id) for id in device_id] - raise NotImplementedError("Multi-GPU support is not implemented yet.") - else: - if isinstance(device_id, list): - if len(device_id) > 1: - raise ValueError( - "Single-GPU mode does not support multiple device IDs." - ) - device_id = device_id[0] - elif device_id is None: - device_id = 0 - utils.check_cuda_device_availability(device_id) - self.device = cuda.Device(device_id) - self.context = self.device.make_context() - - utils.check_cuda_capability(field_config, cuda_config, device_id) - - self.use_kernel_v2 = bool(use_kernel_v2) - self.collision_model = int(collision_model) - self.streaming_model = int(streaming_model) - self.store_precision = int(store_precision) - self.use_ddf_shifting = int(use_ddf_shifting) - self.use_les = int(use_les) - self.les_cs = float(les_cs) - - if self.collision_model not in (0, 1, 2): - raise ValueError("collision_model must be 0(SRT), 1(TRT), or 2(MRT).") - if self.streaming_model not in (0, 1): - raise ValueError("streaming_model must be 0(double-buffer) or 1(esopull).") - if self.store_precision not in (0, 1, 2): - raise ValueError("store_precision must be 0(FP32), 1(FP16S), or 2(FP16C).") - if self.use_ddf_shifting not in (0, 1): - raise ValueError("use_ddf_shifting must be 0 or 1.") - if self.use_les not in (0, 1): - raise ValueError("use_les must be 0 or 1.") - if not (0.0 < self.les_cs < 1.0): - raise ValueError("les_cs must be in (0, 1).") - - # Set constants - if field_config.data_type == "FP32": - self.DATA_TYPE = np.float32 - else: - raise ValueError(f"Unsupported data type {field_config.data_type}.") - - self.FIELD_SHAPE = tuple( - size * unit - for size, unit in zip( - field_config.field_dim_in_U, cuda_config.unit_dimensions - ) - ) - self.FIELD_SIZE = np.prod(self.FIELD_SHAPE) - self.LATTICE = field_config.lattice - self.DIM = field_config.dimensionality - if field_config.lattice == 9 and field_config.dimensionality == 2: - self.E = np.array( - [0, 0, 1, 0, 0, 1, -1, 0, 0, -1, 1, 1, -1, 1, -1, -1, 1, -1], - dtype=np.int32, - ).reshape(9, 2) - self.OPP = np.array([0, 3, 4, 1, 2, 7, 8, 5, 6], dtype=np.int32) - self.WW = np.array( - [4 / 9, 1 / 9, 1 / 9, 1 / 9, 1 / 9, 1 / 36, 1 / 36, 1 / 36, 1 / 36], - dtype=self.DATA_TYPE, - ) - else: - raise NotImplementedError( - f"Unsupported lattice type {field_config.lattice} in {field_config.dimensionality} dimensions." - ) - - self.objects = {} - - # Compile and load kernel - self._rebuild_kernel() - - # Initialize memory - self.ddf = np.zeros(self.FIELD_SIZE * self.LATTICE, dtype=self.DATA_TYPE) - self.ddf_save = np.zeros(self.FIELD_SIZE * self.LATTICE, dtype=self.DATA_TYPE) - self.flag = np.ones(self.FIELD_SIZE, dtype=np.uint8) - self.indx = np.zeros(self.FIELD_SIZE, dtype=np.int32) - self.delta_curve = np.zeros(0, dtype=self.DATA_TYPE) - self.vortex_config = np.zeros(7, dtype=float) - - self.ddf_gpu = cuda.mem_alloc(self.ddf.nbytes) - self.temp_gpu = cuda.mem_alloc(self.ddf.nbytes) - self.flag_gpu = cuda.mem_alloc(self.flag.nbytes) - self.indx_gpu = cuda.mem_alloc(self.indx.nbytes) - self.delta_gpu = cuda.mem_alloc(1) - self.vortex_gpu = cuda.mem_alloc(self.vortex_config.nbytes) - - self.action = np.zeros(0, dtype=self.DATA_TYPE) - self.obs = np.zeros(0, dtype=self.DATA_TYPE) - - self.initflow( - self.flag_gpu, - self.ddf_gpu, - block=(self.cuda_config.threads_per_block, 1, 1), - grid=( - int(self.FIELD_SHAPE[0] / self.cuda_config.threads_per_block), - int(self.FIELD_SHAPE[1]), - int(self.FIELD_SHAPE[2]), - ), - ) - cuda.memcpy_dtoh(self.flag, self.flag_gpu) - cuda.memcpy_dtoh(self.ddf, self.ddf_gpu) - - def _configure_kernel(self): - if self.use_kernel_v2: - compiler.config_kernal_v2( - self.cuda_config, - self.field_config, - collision_model=self.collision_model, - streaming_model=self.streaming_model, - store_precision=self.store_precision, - use_ddf_shifting=self.use_ddf_shifting, - use_les=self.use_les, - les_cs=self.les_cs, - ) - else: - compiler.config_kernal(self.cuda_config, self.field_config) - - def _compile_and_load_kernel(self): - if self.use_kernel_v2: - compiler.compile_kernel_v2() - self.ptx = cuda.module_from_file(compiler.kernel_path("kernel_v2.ptx")) - self.step = self.ptx.get_function("OneStep") - self.initflow = self.ptx.get_function("InitTubeFlow_v2") - else: - compiler.compile_kernel() - self.ptx = cuda.module_from_file(compiler.kernel_path("kernel.ptx")) - self.step = self.ptx.get_function("OneStep") - self.initflow = self.ptx.get_function("InitTubeFlow") - - def _rebuild_kernel(self): - self._configure_kernel() - compiler.config_object(len(self.objects)) - self._compile_and_load_kernel() - - def add_cylinder(self, center: Tuple[float, float, float], radius: float, id_obj: Optional[int] = None): - x_c, y_c, z_c = center - - if ( - x_c - radius <= 0 - or x_c + radius >= self.FIELD_SHAPE[0] - 1 - or y_c - radius <= 0 - or y_c + radius >= self.FIELD_SHAPE[1] - 1 - ): - raise ValueError("Cylinder is out of bounds.") - - index = self.delta_curve.size if self.delta_curve.size > 0 else 0 - - if self.DATA_TYPE == np.float32: - id_object = np.int32(len(self.objects)) - # max_id = max(self.objects.keys()) - else: - raise ValueError(f"Unsupported data type {self.DATA_TYPE}.") - - # Ensure host-side DDF mirrors current device state before local edits. - cuda.memcpy_dtoh(self.ddf, self.ddf_gpu) - - for x in range(int(x_c - radius) - 1, int(x_c + radius) + 1): - for y in range(int(y_c - radius) - 1, int(y_c + radius) + 1): - if (x - x_c) ** 2 + (y - y_c) ** 2 < radius**2: - k = x + y * self.FIELD_SHAPE[0] - self.flag[k] = SOLID - for i in range(self.LATTICE): - self.ddf[k + i * self.FIELD_SIZE] = self.WW[i] - delta_temp = np.zeros(11, dtype=self.DATA_TYPE) - delta_temp[0] = id_object.view(self.DATA_TYPE) - for i in range(self.LATTICE): - x_neb = x + self.E[i][0] - y_neb = y + self.E[i][1] - if (x_neb - x_c) ** 2 + (y_neb - y_c) ** 2 >= radius**2: - self.flag[k] |= INTERFACE - x_i, y_i = preproc.find_circle_intersection( - x, y, x_neb, y_neb, x_c, y_c, radius - ) - d_neb = np.sqrt((x_i - x_neb) ** 2 + (y_i - y_neb) ** 2) - delta_temp[i] = d_neb / np.sqrt( - self.E[i][0] ** 2 + self.E[i][1] ** 2 - ) - if self.flag[k] & INTERFACE: - delta_temp[9] = (y_c - y) / radius - delta_temp[10] = (x - x_c) / radius - self.delta_curve = np.concatenate( - (self.delta_curve, delta_temp) - ) - self.indx[k] = index - index += delta_temp.size - - self.objects[id_object] = { - "type": "cylinder", - "center": center, - "radius": radius, - } - - if hasattr(self, "delta_gpu"): - self.delta_gpu.free() - self.delta_gpu = cuda.mem_alloc(self.delta_curve.nbytes) - - self.action = np.zeros(len(self.objects), dtype=self.DATA_TYPE) - if hasattr(self, "action_gpu"): - self.action_gpu.free() - self.action_gpu = cuda.mem_alloc(self.action.nbytes) - - self.obs = np.zeros(len(self.objects) * self.DIM, dtype=self.DATA_TYPE) - if hasattr(self, "obs_gpu"): - self.obs_gpu.free() - self.obs_gpu = cuda.mem_alloc(self.obs.nbytes) - - cuda.memcpy_htod(self.delta_gpu, self.delta_curve) - cuda.memcpy_htod(self.flag_gpu, self.flag) - cuda.memcpy_htod(self.indx_gpu, self.indx) - cuda.memcpy_htod(self.ddf_gpu, self.ddf) - cuda.memcpy_htod(self.temp_gpu, self.ddf) - - self._rebuild_kernel() - - def add_sensor(self, center: Tuple[float, float, float], radius: float): - x_c, y_c, z_c = center - - if ( - x_c - radius <= 0 - or x_c + radius >= self.FIELD_SHAPE[0] - 1 - or y_c - radius <= 0 - or y_c + radius >= self.FIELD_SHAPE[1] - 1 - ): - raise ValueError("Sensor is out of bounds.") - - id_object = len(self.objects) - for x in range(int(x_c - radius) - 1, int(x_c + radius) + 1): - for y in range(int(y_c - radius) - 1, int(y_c + radius) + 1): - if (x - x_c) ** 2 + (y - y_c) ** 2 < radius**2: - k = x + y * self.FIELD_SHAPE[0] - self.flag[k] |= SENSOR - self.indx[k] = id_object - - self.objects[id_object] = { - "type": "sensor", - "center": center, - } - - self.action = np.zeros(len(self.objects), dtype=self.DATA_TYPE) - if hasattr(self, "action_gpu"): - self.action_gpu.free() - self.action_gpu = cuda.mem_alloc(self.action.nbytes) - - self.obs = np.zeros(len(self.objects) * self.DIM, dtype=self.DATA_TYPE) - if hasattr(self, "force_gpu"): - self.obs_gpu.free() - self.obs_gpu = cuda.mem_alloc(self.obs.nbytes) - - cuda.memcpy_htod(self.flag_gpu, self.flag) - cuda.memcpy_htod(self.indx_gpu, self.indx) - - self._rebuild_kernel() - - def add_vortex(self, center: Tuple[float, float, float], radius: float, strength: float, direction: float, type: str): - x_c, y_c, z_c = center - - if ( - x_c - radius <= 0 - or x_c + radius >= self.FIELD_SHAPE[0] - 1 - or y_c - radius <= 0 - or y_c + radius >= self.FIELD_SHAPE[1] - 1 - ): - raise ValueError("Vortex is out of bounds.") - - if type not in ["lamb", "oseen", "taylor"]: - raise ValueError("Vortex type" + type + " is not supported.") - - x = np.linspace(-x_c, self.FIELD_SHAPE[0] - 1 - x_c, self.FIELD_SHAPE[0]) - y = np.linspace(-y_c, self.FIELD_SHAPE[1] - 1 - y_c, self.FIELD_SHAPE[1]) - X, Y = np.meshgrid(x, y) - r = np.sqrt(X**2 + Y**2) - nu = self.field_config.viscosity - theta = np.arctan2(Y, X) - psi = np.zeros_like(r) - - if type == "lamb": - b = 3.831705970207512 - n = b / radius - u0 = strength - inside = r <= radius - outside = r > radius - - psi[inside] = (2 * u0 / n / jv(0, b) * jv(1, n * r[inside]) - u0 * r[inside]) * np.sin(theta[inside]) - psi[outside] = -u0 * radius**2 / r[outside] * np.sin(theta[outside]) - - u_vor = np.gradient(psi, axis=0) - v_vor = -np.gradient(psi, axis=1) - p_vor = -2 * (np.gradient(v_vor, axis=1) - np.gradient(u_vor, axis=0)) * psi - (u_vor**2 + v_vor**2) / 2 - elif type == "oseen": - # 4 nu t = radius^2 / 4 - kappa = 2 * np.pi * radius **2 * strength - u_vor = - kappa / (2 * np.pi * r) * (1 - np.exp(-4 * r**2 / radius**2)) * np.sin(theta) - v_vor = kappa / (2 * np.pi * r) * (1 - np.exp(-4 * r**2 / radius**2)) * np.cos(theta) - zeta = 4 * r**2 / radius**2 - p_vor = -kappa**2 / 8 / np.pi**2 / r**2 * (-2 * zeta * (expi(-zeta) - expi(-2 * zeta)) + (1 - np.exp(-zeta))**2) - elif type == "taylor": - # 4 nu t = radius^2 - M = strength * np.pi * radius**4 / 8 / nu - u_vor = - M * r * 4 * nu / radius**4 * np.exp(-r**2 / radius**2) * np.sin(theta) - v_vor = M * r * 4 * nu / radius**4 * np.exp(-r**2 / radius**2) * np.cos(theta) - p_vor = -4 * M**2 * nu**2 * np.exp(-2 * r**2 / radius**2) / np.pi**2 / radius**6 - - cuda.memcpy_dtoh(self.ddf, self.ddf_gpu) - ddf_temp = self.ddf.copy().reshape((self.LATTICE, self.FIELD_SHAPE[1], self.FIELD_SHAPE[0])).transpose(2, 1, 0) - u_ddf = ddf_temp[:, :, 1] + ddf_temp[:, :, 5] + ddf_temp[:, :, 8] - ddf_temp[:, :, 3] - ddf_temp[:, :, 6] - ddf_temp[:, :, 7] - v_ddf = ddf_temp[:, :, 2] + ddf_temp[:, :, 5] + ddf_temp[:, :, 6] - ddf_temp[:, :, 4] - ddf_temp[:, :, 7] - ddf_temp[:, :, 8] - p_ddf = np.sum(ddf_temp, axis=2) / 3 - - for i in range(self.FIELD_SHAPE[0]): - for j in range(self.FIELD_SHAPE[1]): - k = i + j * self.FIELD_SHAPE[0] - if (j == 0 or j == self.FIELD_SHAPE[1] - 1) or (i == 0 or i == self.FIELD_SHAPE[0] - 1): - continue - else: - for e in range(self.LATTICE): - u = u_ddf[i, j] + u_vor[j, i] - v = v_ddf[i, j] + v_vor[j, i] - p = p_ddf[i, j] + p_vor[j, i] - eu = self.E[e][0] * u + self.E[e][1] * v - u2 = u ** 2 + v ** 2 - self.ddf[k + e * self.FIELD_SIZE] = self.WW[e] * (3 * p + 3 * eu + 4.5 * eu ** 2 - 1.5 * u2) - - cuda.memcpy_htod(self.ddf_gpu, self.ddf) - - # def add_vortex_gpu(self, center: Tuple[float, float, float], radius: float, strength: float, direction: float, type: str): - # x_c, y_c, z_c = center - - # if ( - # x_c - radius <= 0 - # or x_c + radius >= self.FIELD_SHAPE[0] - 1 - # or y_c - radius <= 0 - # or y_c + radius >= self.FIELD_SHAPE[1] - 1 - # ): - # raise ValueError("Vortex is out of bounds.") - - # if type not in ["lamb", "oseen", "taylor"]: - # raise ValueError("Vortex type" + type + " is not supported.") - - # add_vortex = self.ptx.get_function("AddVortex") - - # self.vortex_config[0:3] = np.array(center, dtype=float) - # self.vortex_config[3] = radius - # self.vortex_config[4] = strength - # self.vortex_config[5] = direction - # if type == "taylor": - # self.vortex_config[6] = - - def run(self, num_steps: int, action_target: np.ndarray): - if ( - action_target.size != len(self.objects) - or action_target.dtype != self.DATA_TYPE - ): - raise ValueError("action data type or size does not match the objects.") - elif len(self.objects) == 0: - raise ValueError("No objects have been added to the flow field.") - - weight = 0.1 - stream = cuda.Stream() - action_pinned = cuda.pagelocked_empty_like(self.action) - action_pinned[:] = self.action - obs_pinned = cuda.pagelocked_empty_like(self.obs) - self.obs[:] = 0 - for i in range(num_steps): - action_pinned = (1 - weight) * action_pinned + weight * action_target - cuda.memcpy_htod_async(self.action_gpu, action_pinned, stream) - self.step( - self.flag_gpu, - self.ddf_gpu, - self.temp_gpu, - self.indx_gpu, - self.delta_gpu, - self.action_gpu, - self.obs_gpu, - block=(self.cuda_config.threads_per_block, 1, 1), - grid=( - int(self.FIELD_SHAPE[0] / self.cuda_config.threads_per_block), - int(self.FIELD_SHAPE[1]), - int(self.FIELD_SHAPE[2]), - ), - stream=stream, - ) - self.ddf_gpu, self.temp_gpu = self.temp_gpu, self.ddf_gpu - cuda.memcpy_dtoh_async(obs_pinned, self.obs_gpu, stream) - cuda.memset_d32_async(self.obs_gpu, 0, self.obs.size, stream) - self.obs += obs_pinned - stream.synchronize() - self.obs = (self.obs / num_steps).astype(self.DATA_TYPE) - - def apply_ddf(self): - cuda.memcpy_htod(self.ddf_gpu, self.ddf) - - def get_ddf(self): - cuda.memcpy_dtoh(self.ddf, self.ddf_gpu) - - def save_ddf(self): - self.ddf_save = self.ddf.copy() - - def restore_ddf(self): - self.ddf = self.ddf_save.copy() - - def __del__(self): - self.context.pop() \ No newline at end of file diff --git a/legacy/macros.h b/legacy/macros.h deleted file mode 100644 index 7962b26..0000000 --- a/legacy/macros.h +++ /dev/null @@ -1,108 +0,0 @@ -// CelerisLab/kernels/macros.h - -// cuda parameters -#define MULT_GPU False -#define NT 128 -#define X_1U 384 -#define Y_1U 192 -#define Z_1U 1 - -// flow parameters -#define LBtype float -#define UX 1 -#define UY 1 -#define UZ 1 -#define NX 384 -#define NY 192 -#define NZ 1 -#define DIM 2 -#define NQ 9 -#define VIS 0.0144000000 -#define RHO 1.0 -#define U0 0.04 - -// constants -#define PI 3.141592653589793238 -#define FLUID 0b00000001 -#define SOLID 0b00000010 -#define GAS 0b00000100 -#define INTERFACE 0b00001000 -#define SENSOR 0b00010000 - -// vortex type -#define V_TAYLOR 0b00000001 - -// variables -#define N_OBJS 0 -// #define N_SENS 2 - -// ============================================================================ -// New architecture configuration (Stage 1) -// These defaults are safe for backward compatibility. -// compiler.py can override any of them via modify_macro(). -// ============================================================================ - -// Collision model: 0=SRT, 1=TRT, 2=MRT -#ifndef COLLISION_MODEL -#define COLLISION_MODEL 0 -#endif - -// Streaming model: 0=double-buffer, 1=esoteric-pull -#ifndef STREAMING_MODEL -#define STREAMING_MODEL 0 -#endif - -// Storage precision: 0=FP32, 1=FP16S, 2=FP16C -#ifndef STORE_PRECISION -#define STORE_PRECISION 0 -#endif - -// DDF-shifting: 0=off, 1=on -#ifndef USE_DDF_SHIFTING -#define USE_DDF_SHIFTING 0 -#endif - -// LES model: 0=off, 1=Smagorinsky -#ifndef USE_LES -#define USE_LES 0 -#endif - -// Smagorinsky constant C_s -#ifndef LES_CS -#define LES_CS 0.160000f -#endif - -// Inlet profile: 1=parabolic (channel), 0=uniform (external flow) -#ifndef INLET_PROFILE -#define INLET_PROFILE 1 -#endif - -// Outlet mode: 0=non-equilibrium extrapolation, 1=zero-gradient copy (more dissipative) -#ifndef OUTLET_MODE -#define OUTLET_MODE 0 -#endif - -// Outlet blend factor for damped outlet mode (OUTLET_MODE=2): -// f_out = a*(non-eq extrapolation) + (1-a)*(zero-gradient copy) -#ifndef OUTLET_BLEND_ALPHA -#define OUTLET_BLEND_ALPHA 0.700f -#endif - -// Outlet backflow clamp: 0=off, 1=force non-negative streamwise velocity at outlet target -#ifndef OUTLET_BACKFLOW_CLAMP -#define OUTLET_BACKFLOW_CLAMP 1 -#endif - -// Global collision omega guardrails -#ifndef OMEGA_COLLISION_MIN -#define OMEGA_COLLISION_MIN 0.01f -#endif - -#ifndef OMEGA_COLLISION_MAX -#define OMEGA_COLLISION_MAX 1.999f -#endif - -// TRT magic parameter Lambda used to map omega+ -> omega- -#ifndef TRT_MAGIC_PARAM -#define TRT_MAGIC_PARAM 0.187500f -#endif \ No newline at end of file diff --git a/src/CelerisLab/lbm/kernels/boundary/inlet/zou_he_local.cuh b/src/CelerisLab/lbm/kernels/boundary/inlet/zou_he_local.cuh index 5c3d837..ec33219 100644 --- a/src/CelerisLab/lbm/kernels/boundary/inlet/zou_he_local.cuh +++ b/src/CelerisLab/lbm/kernels/boundary/inlet/zou_he_local.cuh @@ -17,11 +17,10 @@ // Free-slip y-walls: at inlet rows y=1 and y=NY-2, pull can source wall nodes for // some known directions. Copy those from stored DDF at (x=1, same y) only. // -// NOTE: This helper is NOT Zou-He-specific. All west inlet schemes that use +// NOTE: This helper is not Zou-He-specific. All west inlet schemes that use // west_velocity_rho_closure_d2q9() need clean known-direction values. The // free-slip wall interferes with these at the top/bottom inlet corners. -// Renamed from repair_zou_he_west_knowns_d2q9 for clarity. The old name is -// kept for backward compatibility during the transition. +// The legacy name is kept for now (not renamed yet). __device__ inline void repair_zou_he_west_knowns_d2q9( float* __restrict__ f, const fpxx* __restrict__ fi_in, diff --git a/src/CelerisLab/lbm/kernels/config/config_grid.h b/src/CelerisLab/lbm/kernels/config/config_grid.h index d0937eb..77da63b 100644 --- a/src/CelerisLab/lbm/kernels/config/config_grid.h +++ b/src/CelerisLab/lbm/kernels/config/config_grid.h @@ -6,8 +6,8 @@ #define NT 256 #define MULT_GPU 0 -#define NX 361 -#define NY 161 +#define NX 512 +#define NY 256 #define NZ 1 // ---- Lattice model (single source of truth) ---- diff --git a/src/CelerisLab/lbm/kernels/config/config_method.h b/src/CelerisLab/lbm/kernels/config/config_method.h index 8f5b2b9..3b82402 100644 --- a/src/CelerisLab/lbm/kernels/config/config_method.h +++ b/src/CelerisLab/lbm/kernels/config/config_method.h @@ -6,18 +6,18 @@ #define COLLISION_MODEL 0 #define STREAMING_MODEL 0 #define STORE_PRECISION 0 -#define USE_DDF_SHIFTING 1 +#define USE_DDF_SHIFTING 0 #define USE_LES 0 #define LES_CS 0.160000f #define LES_CLOSED_FORM 1 -#define INLET_PROFILE 0 -#define INLET_SCHEME 3 +#define INLET_PROFILE 1 +#define INLET_SCHEME 0 #define OUTLET_MODE 0 #define OUTLET_BLEND_ALPHA 0.700f #define OUTLET_BACKFLOW_CLAMP 1 -#define Y_WALL_BC 1 +#define Y_WALL_BC 0 #define OMEGA_COLLISION_MIN 0.01f #define OMEGA_COLLISION_MAX 1.990f diff --git a/src/CelerisLab/lbm/kernels/config/config_objects.h b/src/CelerisLab/lbm/kernels/config/config_objects.h index 8357ada..e7cad65 100644 --- a/src/CelerisLab/lbm/kernels/config/config_objects.h +++ b/src/CelerisLab/lbm/kernels/config/config_objects.h @@ -3,6 +3,6 @@ #ifndef CELERIS_CONFIG_OBJECTS_H #define CELERIS_CONFIG_OBJECTS_H -#define N_OBJS 1 +#define N_OBJS 0 #endif diff --git a/src/CelerisLab/lbm/kernels/config/config_physics.h b/src/CelerisLab/lbm/kernels/config/config_physics.h index a08cc9b..3d6d171 100644 --- a/src/CelerisLab/lbm/kernels/config/config_physics.h +++ b/src/CelerisLab/lbm/kernels/config/config_physics.h @@ -4,7 +4,7 @@ #define CELERIS_CONFIG_PHYSICS_H #define LBtype float -#define VIS 0.0090000000 +#define VIS 0.0035000000 #define RHO 1.0 #define U0 0.03