Lattice Boltzmann Method Simulation Research Articles

Turbulent flow boiling simulation based on pseudopotential lattice boltzmann method with developed wall boundary treatments and unit conversion

The lattice Boltzmann method (LBM) is widely studied for complex flow simulations, including liquid–vapor phase change phenomena. Although previous studies extended the LBM simulation capability from single-phase to two-phase boiling simulations, most of them are restricted to no-flow (pool boiling) or low-flow-rate (laminar flow) boiling conditions owing to the stability problem in collision at low viscosity, representing high Reynolds numbers. However, a high Reynolds number flow is essential for investigating the effects of turbulence on bubble dynamics. Recently, the central moment-based collision method is used to enable high Reynolds numbers for bulk flow; however, extending the method for multiphase-flow simulation with high Reynolds numbers is limited owing to the boundary treatment near the wall. In this study, we augment the previous boundary treatment in LBM simulation with a forcing term, which increases the stability of the simulation even under two-phase conditions. With the improved boundary treatment, the flow boiling phenomena at high Reynolds numbers (up to Re = 107,000) are successfully reproduced. Moreover, in this study, unit conversion is reconsidered with the relationship between conversion factors because the poor outcomes of bubble dynamics are attributed to misinterpreted unit conversion. The new unit conversion method satisfies the necessary conditions for mapping the physical bubble behavior to lattice frame, reproducing bubble dynamics effectively following the Reynolds number, with reasonable spatial and temporal resolution. Consequently, the flow boiling simulation can accurately predict the bubble departure diameter from the real experimental results and the heat transfer trend as a function of the Reynolds number.

An inverse analysis of fluid flow through granular media using differentiable lattice Boltzmann method

This study presents a method for the inverse analysis of fluid flow problems. The focus is put on accurately determining boundary conditions and characterizing the physical properties of granular media, such as permeability, and fluid components, like viscosity. The primary aim is to deduce either constant pressure head or pressure profiles, given the known velocity field at a steady-state flow through a conduit containing obstacles, including walls, spheres, and grains. The lattice Boltzmann method (LBM) combined with automatic differentiation (AD) (AD-LBM) is employed, with the help of the GPU-capable Taichi programming language. A lightweight tape is used to generate gradients for the entire LBM simulation, enabling end-to-end backpropagation. Our AD-LBM approach accurately estimates the boundary conditions for complex flow paths in porous media, leading to observed steady-state velocity fields and deriving macro-scale permeability and fluid viscosity. The method demonstrates significant advantages in terms of prediction accuracy and computational efficiency, making it a powerful tool for solving inverse fluid flow problems in various applications.

LBM Simulation of Free Convection Heat Transfer of Cu/Water Nanofluid in Inclined Cavity with Non-Uniform Heating Temperature Distribution

In the present study, natural convection of Cu/water nanofluid in an inclined square cavity has been investigated numerically using lattice Boltzmann method (LBM). The left wall of the cavity is maintained at a constant hot temperature, while the right wall is subjected to non-uniform temperature distribution. The upper and lower walls are insulated. The fluid flow and heat transfer characteristics are investigated over a wide range of parameters, including Rayleigh numbers (Ra = 103, 104 and 105), solid volume fractions (Φ = 0%, 2%, 4% and 6%), amplitude ratio (A = 0, 0.5 and 1) and phase deviations of sinusoidal temperature distribution (Φ = 0, π/4, π/2, 3π/4 and π). The results are presented graphically in the form of streamlines, isotherms and Nusselt numbers for the different combinations of the considered parameters. The heat transfer rate inside the nanofluid increases as the Rayleigh number and the volume fraction of nanoparticles increase. Further, an increase in the amplitude ratio results in a decrease in the heat transfer rate, with a reduction of up to 8.67% at A = 1. On the other hand, changes in the phase deviation of the sinusoidal temperature distribution lead to an increase in the heat transfer rate, with a surge of up to 32.04% observed at Φ = π.

Effect of the tip speed ratio on the wake characteristics of wind turbines using LBM‐LES

AbstractIn this paper, the wake characteristics of Zell 2000 wind turbine under different tip velocity ratios are studied by using the lattice Boltzmann method and large eddy simulation. The adaptive mesh refinement method is performed to capture the fine flow structure and wake characteristics development. In this paper, we mainly focus on the effect of the tip speed ratio on the flow structure and unsteady characteristics of wind turbine wake. The three‐dimensional flow vorticity structures, the section vorticity diagram, the pressure fluctuation of wake and the lift coefficient of wind turbine wake are utilized to explore the effect of the tip speed ratio on the unsteady physics mechanism of wind turbine wake. With the increase of the tip speed ratio, the distance between two adjacent vortex rings along the axial direction gradually decreases, as the position of the broken vortex circles gradually approaches the center of the blade, separated vortexes are rapidly generated, and the coherent structure appears closer to the wind turbine. A relationship is established between the tip speed ratios and the positions of the broken vortex circles. It is further found that the dominant frequency amplitude gradually increases with the increase of tip speed ratio and the pressure amplitude spectra of vortex increases with the decrease of the distance between the wake and the center of the blade axis. The above series of studies can provide significant physical insight into deep understanding the influence of the tip speed ratio on the wake characteristics of wind turbines.

Anisotropic temperatures in multi-layered 2D materials

For multi-layered 2D materials, although its c-axis has a much lower thermal conductivity than the a-axis, its phonon mean free path has been confirmed to be very long, e.g., in the order of 100s nm at room temperature for multi-layered graphene. An anisotropic specific heat concept has been proposed in the past to explain this very long mean free path. This work carries out detailed atomistic modeling to quantify the anisotropic specific heat concept and reports the discovery of anisotropic temperatures in multi-layered 2D materials under ultrafast surface heating. Extremely fast c-phonon energy transport is discovered, and the non-Fourier effect is observed for both a-phonons and c-phonons. The energy coupling factor between these two modes of phonons is determined to be in the order of 1016 W K−1 m−3, with the specific number depending on the structure location. The anisotropic temperature concept is also quantitatively confirmed based on the lattice Boltzmann method simulation. The anisotropic temperature concept does not violate the physics that temperature is a scalar; rather, it is developed to distinguish the temperatures of phonons that travel in different directions. This concept is universally applicable to other 2D materials to describe the heat conduction in the in-plane and out-of-plane directions that feature different interatomic bonds.

Boosting Bulk Oxygen Transport with Accessible Electrode Nanostructure in Low Pt Loading PEMFCs

AbstractDespite the high potential for reducing carbon emissions and contributing to the future of energy utilization, polymer electrolyte membrane fuel cells (PEMFCs) face challenges such as high costs and sluggish oxygen transport in cathode catalyst layers (CCLs). In this study, the impact of pore size distribution on bulk oxygen transport behavior is explored by introducing nano calcium carbonate of varying particle sizes for pore‐forming. Physicochemical characterizations for are employed to examine the electrode structure, while in situ electrochemical measurements are used to scrutinize bulk oxygen transport resistance, effective oxygen diffusivity () and fuel cell performance. Additionally, the CCLs are constructed with aid of Lattice Boltzmann method (LBM) simulations and for CCLs with different pore size distribution are calculated. The findings reveal that initially increases and then decreases as the most probable pore size increases. A “sphere‐pipe” model is proposed to describe practical bulk oxygen transport in CCLs, highlighting the significant role of not only the pore size of secondary pores but also the number of primary pores in bulk oxygen transport.

Lattice Boltzmann Method Simulation of Nanofluid Natural Convection Heat Transfer in a Square Cavity with Constant Heat Flux at Walls

AbstractIn this study, for the first time, a nanofluid's natural convection heat transfer in a two-dimensional square cavity has been numerically investigated by use of the lattice Boltzmann method with the constant heat flux boundary condition. The horizontal walls of the cavity are insulated, and the vertical walls are kept at a constant heat flux. The diameters of the nanoparticles inside the cavity are the same and have a homogeneous distribution, and there is no chemical reaction between the particles. The flow is also assumed to be the steady state and two-dimensional. Constant temperature, streamlines, velocity, and average Nusselt have been investigated for different nanoparticle volume fractions and Rayleigh numbers. The results showed that the lattice Boltzmann method efficiently analyzes the natural heat transfer of nanofluids; moreover, by use of nanofluid in the cavity increases the heat transfer rate. With the increase in the nanoparticle volume fraction, the average Nusselt number on the right wall of the cavity increased. For a volume fraction of 20% with Grashof number 105, the average Nusselt number increased by almost 50% compared to the base fluid at the same Grashof number. It has been observed that as the volume fraction of nanoparticles in the fluid increases, the fluid’s viscosity also increases; consequently, the velocity of the fluid is found to decrease.

Interactional aerodynamics and acoustics of a rotor with an airframe in hover

The demand for the development of urban air mobility (UAM) powered by electric systems has been steadily rising across private and public sectors. Most UAM flights incorporate a distributed electric propulsion system to enhance aircraft safety and reliability, which entails an increase in the number of rotors or propellers. Consequently, aerodynamics and aeroacoustics are significantly influenced by strong interactions between the rotor and the airframe. In this study, we conducted a computational investigation to examine the impact of rotor–airframe interaction on aerodynamic and aeroacoustic characteristics. This examination considered variations in airframe shape and the distance between the rotor and airframe. The aerodynamic analysis was executed using the lattice-Boltzmann method simulation, in which acoustic predictions were made using the Ffowcs Williams–Hawkings(FW–H) acoustic analogy with a permeable surface. The airframe consists of two geometries: a cylinder and a cone. Tip vortex breakdown and the transition into the turbulent wake state were captured in both airframes, and a fountain flow was affected by the downwash circulation generated under certain proximity of airframe cases. The acoustic prediction results showed that high-intensity noise radiated over the broad surface of the airframe in the conical airframe case. Significant thrust force fluctuations and an increase in noise level were observed at the smallest rotor tip clearance, S/R=−0.1, compared to the isolated rotor. Furthermore, the noise contribution of the rotor and airframe was compared, revealing that the airframe noise level was even higher than the rotor noise at S/R=−0.1.

Rebound dynamics of two droplets asymmetrically impacting single-ridge superhydrophobic surfaces

On a single-ridge superhydrophobic surface, the rebound dynamics of double droplet impacts is investigated via lattice Boltzmann method (LBM) simulations, where one droplet is on the left of macro-ridges, whereas the other droplet being on the top of macro-ridges, leading to the asymmetrical impact of double droplets. The simulation results show that under different impact conditions of various Weber numbers, that is, at low, moderate and high Weber numbers, the rebound regime of double droplet impacts always consists of the complete-coalescence-rebound (CCR), partial-coalescence-rebound (PCR) and coalescence-split-rebound (CSR) regimes. Under these three rebound regimes, the contact time of double droplet impacts is always reduced by the increased We, as the center-to-center distances in the horizontal and vertical directions, Lx and Lz, are constant. At constant Lx and We, the contact time of double droplet impacts on single-ridge superhydrophobic surfaces is increased by the increased Lz because of the earlier formation of liquid bridges at larger Lz. But under different impact conditions of various center-to-center distances in the horizontal direction, that is, at the increased Lx, the rebound regime of double droplet impacts is not always composed of the CCR, PCR and CSR regimes at low, moderate and high Weber numbers. At constant We and Lz, both coalescence strength and rebound regime are significantly affected by the increased Lx for double droplet impacts. So that the contact time of double droplet impacts on singe-ridge superhydrophobic surfaces depends on the center-to-center distance in the horizontal direction, that is, in the x- direction.

Viscous coupling effect on hydraulic conductance in a square capillary tube

Viscous coupling effect on pore-scale flow capacity and its impact on macroscopic flow properties in porous media remain subjects of ongoing debate. This study employs the Stokes equation in conjunction with the zero-velocity interface and continuity interface to investigate viscous coupling effect during two-phase flow through a square capillary tube. Our findings substantiate the significance of the viscous coupling effect in angular pore multiphase flows. Enhanced fluid conductance is linked to fluid viscosity ratios, contact angles, and wetting-film lengths at the pore scale. Proposing a novel upscaling approach, we formulate hydraulic conductance in a square capillary tube for multiphase flows by incorporating a viscous coupling term, and its validation is accomplished through comparison with results from lattice Boltzmann method simulations. Our scaling model predicts hydraulic conductance with a mean relative error of less than 1%, outperforming prior viscous coupling models with errors reaching up to 19%. The derived scaling model, incorporating the viscous coupling effect, holds potential for integration within pore-network models, offering an efficient and precise simulation method for characterizing two-phase flow through porous media at a representative scale.

Hybrid outflow boundary condition for the pseudopotential LBM simulation of flow boiling

The pseudopotential lattice Boltzmann method (LBM) has been successful in simulating the pool boiling with a density ratio up to 1000. However, the flow boiling simulation of using the same LBM is rare and the density ratio is limited to about 17. The current study investigated the challenge of flow boiling by first simulating the isothermal bubble outflow, where large pressure deviation from the initial saturation pressure is observed due to perturbation of the bubble outflow. Pressure drifts continuously as there is no constraint at both inlet and outlet when using the usual convective outflow boundary condition (OBC). The subsequent flow boiling simulation shows that boiling will be suppressed by the drifting/increasing pressure. A hybrid OBC is proposed to stabilize the pressure by imposing a fixed pressure either on the whole outlet when the bubbles are still far away from the outlet, or on a part of the outlet when the bubbles are crossing the other part of the outlet where the convective OBC is still used. With the improvement by using the hybrid OBC, the heat transfer coefficient of saturated nucleate boiling in a channel flow is studied and the results are in good agreement with the experiments. More realistic problems with structured wall surface, various surface properties and density ratios are simulated to show the algorithm stability and robustness. The current study has extended the density ratio to about 125, which is much larger than the previous record of about 17.

Lattice Boltzmann study of the effect of catalyst layer structure on oxygen reduction reaction within a PEMFC

The catalyst layer (CL) structure has a great impact on proton exchange membrane fuel cell (PEMFC) performance. In this study, the microporous structure of a cathode CL is reconstructed via a stochastic algorithm. Subsequently, the lattice Boltzmann method (LBM) is used to investigate the influence of the CL structure regarding carbon carrier and platinum (Pt) particle sizes and ionomer thickness on oxygen reduction reaction (ORR) processes within the cathode CL. The LBM simulation results indicate that an increased carbon carrier size promotes oxygen transport through the CL pores owing to the decreased tortuosity of the porous CL but reduces the overall ORR rate. A larger Pt particle size leads to a reduction in the number of reaction sites but shortens the ionomer thickness for oxygen diffusion to the reaction sites. The oxygen transport resistance through the ionomer increases with increasing ionomer thickness and thus degrades the overall ORR rate. However, an excessively small ionomer thickness leads to the exposure and shedding of Pt, thereby reducing the number of reaction sites. The findings here demonstrate that proper design of the carbon carrier and Pt particle sizes and ionomer thickness are essential for achieving a CL with high oxygen transport performance and overall ORR rate.

LBM simulation for combined thermal radiation and natural convection in 2D enclosure with multiple solid blocks

A numerical investigation of conjugate natural convection in a 2D square enclosure with multiple solid blocks, considering the impact of volumetric thermal radiation, was conducted using the lattice Boltzmann method (LBM). The effect of Planck number wall surface emissivity extinction coefficient solid-to-fluid volume fraction and the number of solid blocks on the fluid-structure and heat transfer within the enclosure were analyzed. The numerical model was validated using benchmark problems from the literature, showing good agreement. The results revealed that at and increasing wall surface emissivity significantly increased the Nusselt number, ranging from approximately 51-116% for different wall emissivity values, as compared to scenarios without surface radiation (). On the other hand, reducing χ, Pl, and β improved the heat transfer rate. The results also demonstrated that subdividing the initial block reduced the heat transfer rate depending on the number of blocks and the value of the Planck number. For instance, when 64 solid blocks were present in the cavity and the Nusselt number decreased by 59.36% compared to the case with a lower solid-to-fluid volume fraction (χ = 15%). The parametric analysis provided useful insights for engineering design in the modern construction industry, with possible applications in building block design and energy efficiency enhancements.

Unrestricted component count in multiphase lattice Boltzmann: A fugacity-based approach.

Studies of multiphase fluids utilizing the lattice Boltzmann method (LBM) are typically severely restricted by the number of components or chemical species being modeled. This restriction is particularly pronounced for multiphase systems exhibiting partial miscibility and significant interfacial mass exchange, which is a common occurrence in realistic multiphase systems. Modeling such systems becomes increasingly complex as the number of chemical species increases due to the increased role of molecular interactions and the types of thermodynamic behavior that become possible. The recently introduced fugacity-based LBM [Soomro et al., Phys. Rev. E 107, 015304 (2023)2470-004510.1103/PhysRevE.107.015304] has provided a thermodynamically consistent modeling platform for multicomponent, partially miscible LBM simulations. However, until now, this fugacity-based LB model had lacked a comprehensive demonstration of its ability to accurately reproduce thermodynamic behavior beyond binary mixtures and to remove any restrictions in a number of components for multiphase LBM. In this paper we closely explore these fugacity-based LBM capabilities by showcasing comprehensive, thermodynamically consistent simulations of multiphase mixtures of up to ten chemical components. The paper begins by validating the model against the Young-Laplace equationfor a droplet composed of three components. The model is then applied to study mixtures with a range of component numbers from one to six, showing agreement with rigorous thermodynamic predictions and demonstrating linear scaling of computational time with the number of components. We further investigate ternary systems in detail by exploring a wide range of temperature, pressure, and overall composition conditions to produce various characteristic ternary diagrams. In addition, the model is shown to be unrestricted in the number of phases as demonstrated through simulations of a three-component three-phase equilibrium case. The paper concludes by demonstrating simulations of a ten-component, realistic hydrocarbon mixture, achieving excellent agreement with thermodynamics for both flat interface vapor-liquid equilibrium and curved interface spinodal decomposition cases. This study represents a significant expansion of the scope and capabilities of multiphase LBM simulations that encompass multiphase systems of keen interest in engineering.

Lattice Boltzmann modeling of the effective thermal conductivity in plant fiber porous media generated by Quartet Structure Generation Set

This study proposes a novel approach that combines a Quartet Structure Generation Set (QSGS) algorithm with the Lattice Boltzmann Method (LBM) to investigate effective thermal conductivity in plant fiber materials. The QSGS algorithm is utilized to construct a two-dimensional digital structure of the porous media, enabling the determination of effective thermal conductivity through LBM simulations. The simulation results demonstrate good agreement with experimental data, exhibiting a maximum relative error of only 1.7%. The study examines the influence of key factors on effective thermal conductivity, including temperature, porosity, fiber distribution, and fractal dimension. The effective thermal conductivity increased by approximately 6.1% for each 10 °C increase in temperature between 10 °C and 50 °C. Moreover, the effective thermal conductivity showed an average decrease of 16.3% with each 10% increment in porosity. The increase in thermal conductivity becomes more significant with increasing temperature and decreasing porosity. Conversely, an increased fractal dimension correlates with a decline in effective thermal conductivity. The effective thermal conductivity is influenced synergistically by the directional growth probability of the solid phase and the direction of heat flow. The findings offer valuable insights into thermal conductivity behavior and heat transfer mechanisms in plant fiber porous media.

Experimental and numerical study on the two-phase flow inside a cracked gas diffusion layer of PEMFC

The presence of cracks in the gas diffusion layer (GDL) can have a significant impact on the flow of two phases within Proton Exchange Membrane Fuel Cells (PEMFCs). This study focuses on investigating the transport process of two-phases flow in a GDL sample using both experimental and numerical methods. The porous structure of a hydrophobic carbon paper, serving as the macro-porous substrate (MPS) of the GDL, is examined using synchrotron radiation X-ray phase contrast Computed Tomography (CT). A 3-D model with various crack geometries is then created, and a Lattice Boltzmann Method (LBM) simulation is employed to analyze the flow behavior within the GDL model. The findings reveal that the water intrusion mode in the MPS primarily depends on its hydrophobic nature, but the width and depth of the cracks also significantly influence both the rate and path of water intrusion. Increasing the crack width from 7 μm to 28 μm results in a 30% increment in water intrusion rate, while the intrusion path remains unchanged. However, when the crack width reaches 42 μm, the crack itself becomes the primary path for water intrusion, doubling the intrusion rate. Furthermore, an increase in crack depth leads to a higher invasion rate of the liquid phase, consequently increasing the liquid storage volume within the GDL. This research provides a novel approach for analyzing cracked GDLs and demonstrates the potential to optimize the water management process of PEMFCs through GDL microstructure design.

Physical and Numerical Simulation of Tight Gas Flow at the Microscale

The porous media in tight reservoirs are mainly composed of micro- and nanopores, gas seepage through which is complex, making it difficult to study. Physical simulation using micron tubes is an intuitive and effective method to study the seepage mechanism of tight gas. The lattice Boltzmann method (LBM) is the most effective method for the tight gas seepage simulation, and it has been widely used. Microscale gas seepage simulation experiments and LBM simulations of micron tubes with different inner diameters were performed. The results showed that in micron tubes, the gas flow increases nonlinearly with an increasing pressure gradient. Influenced by compression and rarefaction effects, the degree of the nonlinearity of pressure distribution in series micron tubes increases with inlet pressure. The existence of a connecting channel between parallel micron tubes breaks the linear distribution of pressure in the original micron tubes, and the gas forms a raised relative high-pressure area at the connection of the two micron tubes; the wider the channel, the greater the bulge. The average gas flow rate in the whole micron tube increases with the channel width, and the seepage capacity increases instead of decreases. The diameter change of one micron tube has no effect on the gas flow in the other micron tube. Although the two micron tubes are connected, they are still relatively independent individuals. These research results lay a foundation for the correct understanding of the characteristics and laws of tight gas seepage in the pores of reservoirs at the micro- and nanoscales, and they have important theoretical significance for the study of seepage mechanisms in tight gas reservoirs.

Pore Structure and Permeability of Tight-Pore Sandstones: Quantitative Test of the Lattice–Boltzmann Method

This paper presents a characterization of the pore structure of tight-pore sandstones of the Achimov suite and examines the application of Lattice–Boltzmann method (LBM) simulations to estimate the permeabilities of rock formations with a single-scale porosity. Porosity is characterized by pore volume distribution, pore throat connectivity, and tortuosity, which are calculated from 3D computer tomography pore network maps. The tight sandstones are poorly permeable, with permeabilities from 0.7 to 13 mD. For comparison, sandstones and carbonates with higher porosity and permeability from the existing database are also considered. For the more permeable reference samples with wider pores (250 µm), LBM simulations show good agreement with the experiments and somewhat outperform the selected state-of-the-art direct simulations from the literature. For samples with the tightest pores and lowest porosity, LBM simulations tend to somewhat overestimate the permeability in comparison with the direct simulation methods, whereas for samples of higher porosity, a slight underestimation is obtained. We explain the inconsistencies by an interplay between the compressibility effects neglected by our LBM simulations in wider pores and the friction at the pore-wall interface, which is underestimated due to the use of the bounce-back conditions. However, the general agreement with experimental and direct simulation methods is very reasonable and suitable for practical use, which means that LBM is fast, highly parallel, and computationally sound even in tight pores.

Isothermal CO2 injection into water-saturated porous media: Lattice-Boltzmann modelling of pulsatile flow with porosity, tortuosity, and optimal frequency characterization

The Lattice Boltzmann Method (LBM) is used to simulate the isothermal injection of CO2 into a water-saturated, homogeneous porous medium. The complex multiphase flow is modeled using the 2D-pseudopotential Shan-Chen multiphase LBM. This work proposes a unique fluid-pulse technology similar to the one applied to enhance oil recovery with CO2 sequestration. This technology is different from classical CO2 sequestration in aquifers and, to our knowledge, has never been used before. The primary purpose is to investigate whether our LBM simulation allows the introduction of fluid-pulse CO2 injection patterns similar to those reported in the literature on hydrocarbon reservoirs and to evaluate the benefit of using such a technique compared to classical uniform injection. First, our simulations consider a transient homogeneous flow with uniform inlet velocity to determine the natural frequency of the instabilities triggered downstream of the interstitial pore array. Then, CO2 is injected at the inlet in a pulsating regime, varying the forcing frequency near and far from the natural frequency to determine its effect on CO2 penetration. Results demonstrate that a pulsatile injection effectively enhances CO2 sequestration in aquifers when forced with a characteristic frequency. However, this frequency resulting after CO2 is injected to displace the water differed from the natural frequency preliminarily obtained from the homogeneous water flow. Injecting at the characteristic frequency has proven an effective technique for boosting CO2 penetration into a water-saturated porous medium, with a significant enhancement estimated between 1.5% and 16% for a range of pulsatile amplitudes between 1% and 10%.

Coupled Lattice Boltzmann and Discrete Element Simulations of Ship-Ice Interactions

Evaluating ice loads acting on ships is essential for the safety of ships navigating in ice-covered seas. In this study, we develop a CFD method to handle ship, ice, and fluid interaction. The lattice Boltzmann method, capable of large-scale calculations, is applied to the simulation of free-surface fluids. The ice motion is computed by solving the equations of motion of a rigid body, and the discrete element method models the ice-ice and ice-ship contact interactions. A momentum exchange scheme couples the lattice Boltzmann method and particle-based rigid body simulation. We introduce tree-based adaptive mesh refinement and multiple GPU computing to improve grid resolution and computational time. The proposed method is applied to model scale simulations of ship navigation in a brash ice channel. Simulations were performed for various conditions with different ice concentrations and ship velocities, and we observed that ice resistance increased with the ice concentration and the ship velocity increased. The ice motions and resistances obtained from our simulations are reasonable compared to model equations of Finnish-Swedish ice class rules (FSICR) and numerical analyses of a previous study.