Large Time Steps Research Articles

An implicit unified gas-kinetic particle method with large time steps for gray radiation transport

For a long time, efficient algorithms for high-dimensional equations, represented by photon radiation transport, have been one important topic in the development of computational methods for particle transport processes. In this paper, we present an implicit unified gas-kinetic particle (IUGKP) method for multiscale gray radiative transfer. Based on the integral solution of the radiative transfer equation, the photon transport processes are categorized into non-equilibrium transport processes with a large photon free path and equilibrium transport processes with a small photon free path. The long-path processes are solved by an implicit Monte Carlo (IMC) method, and the short-path processes are solved by an implicit diffusion system. The closure formulation of photon distribution is derived from the local integral solution of the radiative transfer equation to couple the IMC and diffusion system. The improvement of the proposed IUGKP method over UGKP method is that particles can be tracked continuously instead of just until the first collision, making simulation with large time steps possible. The IUGKP method has the properties of asymptotic-preserving (AP) and regime-adaptive (RA). The AP property states that the IUGKP method converges to the consistent numerical methods for the asymptotic limiting equations of RTE in the limiting regimes. The RA property states that the computational accuracy of the IUGKP method adapts to the regimes. In this paper, the mathematical proof of the AP and RA properties is presented, and the multiscale numerical tests are performed to demonstrate the accuracy and efficiency of the IUGKP method.

Dynamically regularized Lagrange multiplier schemes with energy dissipation for the incompressible Navier-Stokes equations

In this paper, we present efficient numerical schemes based on the Lagrange multiplier approach for the Navier-Stokes equations. By introducing a dynamic equation (involving the kinetic energy, the Lagrange multiplier, and a regularization parameter), we form a new system which incorporates the energy evolution process but is still equivalent to the original equations. Such nonlinear system is then discretized in time based on the backward differentiation formulas, resulting in a dynamically regularized Lagrange multiplier (DRLM) method. First- and second-order DRLM schemes are derived and shown to be unconditionally energy stable with respect to the original variables. The proposed schemes require only the solutions of two linear Stokes systems and a scalar quadratic equation at each time step. Moreover, with the introduction of the regularization parameter, the Lagrange multiplier can be uniquely determined from the quadratic equation, even with large time step sizes, without affecting accuracy and stability of the numerical solutions. Fully discrete energy stability is also proved with the Marker-and-Cell (MAC) discretization in space. Various numerical experiments in two and three dimensions verify the convergence and energy dissipation as well as demonstrate the accuracy and robustness of the proposed DRLM schemes.

High-order accurate implicit-explicit time-stepping schemes for wave equations on overset grids

New implicit and implicit-explicit time-stepping methods for the wave equation in second-order form are described with application to two and three-dimensional problems discretized on overset grids. The implicit schemes are single step, three levels in time, and based on the modified equation approach. Second and fourth-order accurate schemes are developed and they incorporate upwind dissipation for stability on overset grids. The fully implicit schemes are useful for certain applications such as the WaveHoltz algorithm for solving Helmholtz problems where very large time-steps are desired. Some wave propagation problems are geometrically stiff due to localized regions of small grid cells, such as grids needed to resolve fine geometric features, and for these situations the implicit time-stepping scheme is combined with an explicit scheme: the implicit scheme is used for component grids containing small cells while the explicit scheme is used on the other grids such as background Cartesian grids. The resulting partitioned implicit-explicit scheme can be many times faster than using an explicit scheme everywhere. The accuracy and stability of the schemes are studied through analysis and numerical computations.

Phase Field Coupled Finite Deformation Plasticity Formulation of Ductile Fracture With Nonlinear Kinematic Hardening and Modified Energy Release Function

ABSTRACTA ductile damage theory is presented by coupling the covariant formulation of finite deformation plasticity with the phase field modeling of fracture, including kinematic hardening for the ductile response of the materials. A phase field coupled nonlinear kinematic hardening equation is proposed in the reference configuration having equivalent representation through the Lie derivative of the kinematic hardening tensor by push‐forward operation in the spatial configuration, thus ensuring the satisfaction of frame invariance. To capture the correct physical response of the material by the phase field evolution equation, the fracture driving function, that is, the difference between the sum of elastic and plastic energies and threshold energy, after the damage initiation is modified by an energy release controlling function, which is an empirical relation of equivalent plastic strain. In defining the energy release controlling function, well‐defined points with physical significance in the experimental load versus displacement curve are used. To simulate the response of the material in relatively large time steps, a modified staggered scheme is presented, evaluating the fracture driving and energy release controlling functions from the previous converged step and using the updated phase field variable in the weak form of the momentum balance equation. To quantify different material parameters from available experimental results in the literature, the developed phase field coupled elasto‐plastic model uses a neural network optimization procedure consisting of neural network training together with optimization in MATLAB and finite element model evaluation in Abaqus user element subroutine UEL. Model capabilities are demonstrated by simulating the crack propagation in complex 3D geometries such as the second and third Sandia Fracture Challenges.

An Orbital Basis Set for Double Photoionization of Atoms and Molecules.

The ab initio theoretical treatment of one-photon double photoionization processes has been limited to atoms and diatomic molecules by the challenges posed by large grid-based representations of the double ionized continuum wave function. To provide a path for extensions to polyatomics, an energy-adapted orbital basis approach is demonstrated that reduces the dimensions of such representations and simultaneously allows larger time steps in time-dependent computational descriptions of double ionization. Additionally, an algorithm that exploits the diagonal nature of the two-electron integrals in the grid basis and dramatically accelerates the transformation between grid and orbital representations is presented. Excellent agreement between the present results and benchmark theoretical calculations is found for H- and Be atoms, as well as the hydrogen molecule, including for the triply differential cross sections that relate the angular distribution and energy sharing of all of the particles in the molecular frame.

A multi-scale IMEX second order Runge-Kutta method for 3D hydrodynamic ocean models

Understanding complex physical phenomena often involves dealing with partial differential equations (PDEs) where different phenomena exhibit distinct timescales. Fast terms, associated with short characteristic times, coexist with slower ones requiring relatively longer time steps for resolution. The challenge becomes more manageable when, despite the varying characteristic times of fast and slow terms, the computational cost associated with faster terms is significantly lower than that of slower terms. Additionally, slower terms can also exhibit two distinct longer characteristic times, adding complexity to the system and resulting in a total of three characteristic timescales. In this paper, an innovative split second-order IMEX (IMplicit-EXplicit) temporal scheme is introduced to address this temporal complexity. It is used to solve the primitive equation ocean model. Extremely short times are handled explicitly with small time steps, while longer timescales are managed explicitly and semi-implicitly using larger time steps. The decision to solve a portion of the slower terms semi-implicitly is due to the fact that it does not significantly increase the total computational cost, allowing for greater flexibility in the time step without imposing a substantial burden on the overall computational efficiency. This strategy enables efficient management of the various temporal scales present in the equations, thereby optimizing computational resources. The proposed scheme is applied to solve 3D hydrodynamics equations encompassing three time scale: fast terms representing wave phenomena, slow terms describing horizontal aspects and stiff terms for vertical ones. Furthermore, the scheme is designed to respect crucial physical properties, namely global and local conservation. The obtained results on different test cases demonstrate the robustness and efficiency of the IMEX approach in simulating these complex systems.

A Novel Approach Integrating the Method of Characteristics with Large Time-Step Scheme (MOC-LTS) for Efficient Transient Flow Simulation in Liquid Pipelines

Water hammer incidents pose a significant risk to the safety and stability of pipeline operations. Therefore, rapid and precise transient flow simulation is essential for efficiently developing scientific water hammer control strategies. Nevertheless, the prevailing transient flow simulation methods for liquid pipelines predominantly employ explicit schemes to solve transient flow control equations, necessitating adherence to the stability criterion of Courant-Friedrichs-Lewy (CFL) ≤ 1. This results in limited time step sizes, which in turn constrains computational efficiency, particularly when simulating hydraulic behavior in large-scale pipeline networks. In this paper, a novel approach integrating the method of characteristics (MOC) with a large time-step scheme (LTS) is proposed to enable rapid and accurate simulation of transient flow in liquid pipelines. The proposed approach discretizes the computational domain into contiguous control volumes, ategorizing them as either boundary or internal control volumes depending on whether they are affected by boundary conditions within a large time step. For internal control volumes, an LTS scheme based on the first-order Godunov format is employed to improve computational efficiency. For boundary control volumes, MOC is applied iteratively to accurately capture boundary dynamics until the simulated time matches that of the internal control volumes. Quantitative analysis is conducted to verify the performance of the proposed approach. The results confirm that the proposed approach outperforms the MOC in computational efficiency while maintaining minimal accuracy loss.

Designing a Fully‐Tunable and Versatile TKE‐l Turbulence Parameterization for the Simulation of Stable Boundary Layers

AbstractThis study presents the development of a so‐called Turbulent Kinetic Energy (TKE)‐l, or TKE‐l, parameterization of the diffusion coefficients for the representation of turbulent diffusion in neutral and stable conditions in large‐scale atmospheric models. The parameterization has been carefully designed to be completely tunable in the sense that all adjustable parameters have been clearly identified and the number of parameters has been minimized as much as possible to help the calibration and to thoroughly assess the parametric sensitivity. We choose a mixing length formulation that depends on both static stability and wind shear to cover the different regimes of stable boundary layers. We follow a heuristic approach for expressing the stability functions and turbulent Prandlt number in order to guarantee the versatility of the scheme and its applicability for planetary atmospheres composed of an ideal and perfect gas such as that of Earth and Mars. Particular attention has been paid to the numerical stability and convergence of the TKE equation at large time steps, an essential prerequisite for capturing stable boundary layers in General Circulation Models (GCMs). Tests, parametric sensitivity assessments and preliminary tuning are performed on single‐column idealized simulations of the weakly stable boundary layer. The robustness and versatility of the scheme are assessed through its implementation in the Laboratoire de Météorologie Dynamique Zoom GCM and the Mars Planetary Climate Model and by running simulations of the Antarctic and Martian nocturnal boundary layers.

An novel finite difference dispersion error elimination mechanism in the Lax–Wendroff high‐order time discretization

AbstractTime domain finite difference methods have been widely used for wave‐equation modelling in exploration geophysics over many decades. When using time domain finite difference methods, it is desirable to use a larger time step so as to save numerical simulation time. The Lax–Wendroff method is one of the well‐known methods to allow larger time step without increasing the time grid dispersion. However, the Lax–Wendroff method suffers from more time consumption because there are more spatial derivatives required to be approximated by the finite difference operators. We propose a new finite difference scheme for the Lax–Wendroff method so as to reduce the numerical simulation time. Then we determine the finite difference operator coefficients and analyse the dispersion error of the proposed finite difference scheme for the Lax–Wendroff method. At last, we apply the proposed finite difference scheme for the Lax–Wendroff method to different velocity models. The numerical simulation results indicate that the proposed finite difference scheme for the Lax–Wendroff method can effectively suppress time grid dispersion and is more efficient compared to the traditional finite difference scheme for the Lax–Wendroff method.

Spatially and temporally high-order dynamic nonlinear variational multiscale methods for generalized Newtonian flows

Non-Newtonian fluids are of interest in industrial sectors, biological problems and other natural phenomena. This work proposes rheologically-dependent, spatially and temporally high-order non-residual stabilized finite element formulations. The accuracy of the methods is assessed by tackling highly-convective time-dependent power-law flows. The spatial approximation uses Lagrangian finite elements up to fourth order. The temporal integration is done via backward differentiation formulas of order one, two and three. A key aspect of our work is using a non-residual orthogonal variational multiscale (VMS) formulation to stabilize dominant convection and to allow equal-order interpolation of velocity and pressure. Our VMS method uses dynamic nonlinear subscales, which have not been tested so far for generalized Newtonian fluids. In this work, the use of high-order temporal discretizations for the subscale components is systematically evaluated. Numerical experiments consider the flow over a confined cylinder for Reynolds numbers between 40 and 400 and power-law indices between 0.4 and 1.8. Numerical testing demonstrates the method to be stable in all combinations of spatial and temporal orders. Our results show that using high-order spatial discretizations more accurately approximates boundary layers and viscosity fields. Moreover, higher temporal orders allow using larger time steps while still capturing highly dynamic behaviors with better resolution in frequency spectra.

Second-Order Modified Nonstandard Explicit Euler and Explicit Runge–Kutta Methods for n-Dimensional Autonomous Differential Equations

Nonstandard finite-difference (NSFD) methods, pioneered by R. E. Mickens, offer accurate and efficient solutions to various differential equation models in science and engineering. NSFD methods avoid numerical instabilities for large time steps, while numerically preserving important properties of exact solutions. However, most NSFD methods are only first-order accurate. This paper introduces two new classes of explicit second-order modified NSFD methods for solving n-dimensional autonomous dynamical systems. These explicit methods extend previous work by incorporating novel denominator functions to ensure both elementary stability and second-order accuracy. This paper also provides a detailed mathematical analysis and validates the methods through numerical simulations on various biological systems.

Dynamics of an acoustically levitated fluid droplet captured by a low-order immersed boundary method

We present a variant of the immersed boundary (IB) method that implements acoustic perturbation theory to model acoustically levitated fluid droplets. Instead of resolving sound waves numerically, our hybrid method solves acoustic scattering semi-analytically to obtain the corresponding time-averaged acoustic forces on the droplet. This framework allows the droplet to be simulated on inertial timescales of interest, and therefore works with much larger time steps than traditional compressible flow solvers. To benchmark this technique and demonstrate its utility, we implement the hybrid IB method for a single droplet in a standing wave. Simulated droplet shape deformations and streaming profiles agree with available theoretical predictions. Our simulations also yield insights into the streaming profiles for elliptical droplets, for which a comprehensive analytic solution does not yet exist.

Exponential time propagators for elastodynamics

We propose a computationally efficient and systematically convergent approach for elastodynamics simulations. We recast the second-order dynamical equation of elastodynamics into an equivalent first-order system of coupled equations, so as to express the solution in the form of a Magnus expansion. With any spatial discretization, it entails computing the exponential of a matrix acting upon a vector. We employ an adaptive Krylov subspace approach to inexpensively and accurately evaluate the action of the exponential matrix on a vector. In particular, we use an apriori error estimate to predict the optimal Krylov subspace size required for each time-step size. We show that the Magnus expansion truncated after its first term provides quadratic and superquadratic convergence in the time-step for nonlinear and linear elastodynamics, respectively. We demonstrate the accuracy and efficiency of the proposed method for one linear (linear cantilever beam) and three nonlinear (nonlinear cantilever beam, soft tissue elastomer, and hyperelastic rubber) benchmark systems. For a desired accuracy in energy, displacement, and velocity, our method allows for 10−100× larger time-steps than conventional time-marching schemes such as Newmark-β method. Computationally, it translates to a ∼1000× and ∼10−100× speed-up over conventional time-marching schemes for linear and nonlinear elastodynamics, respectively.

Exponential time difference methods with spatial exponential approximations for solving boundary layer problems

AbstractIn this work, an efficient and stable exponential time difference method is presented for solving boundary layer problems. By combining exponential time difference schemes with spatial direct discontinuous Galerkin discretization based on exponential boundary layer approximations, the proposed algorithm not only may admit large time step sizes but also could provide good spatial approximations even if on rather coarse spatial grids in the boundary layer. Some energy stabilities of the numerical scheme are rigorously derived. Numerical examples illustrate the accuracy, stability and efficiency of the algorithm.

An efficient numerical solver for highly compliant slender structures in waves: Application to marine vegetation

This paper presents a fully explicit coupled wave–vegetation interaction model capable of efficiently solving the coupled wave dynamics and flexible vegetation motion with large deflections. The flow model is formulated using the continuity equation and linearized momentum equations of an incompressible fluid, with additional terms within the canopy region accounting for the presence of vegetation. This linearized flow solver is unconditionally stable and second-order accurate. The flow model is validated and verified against experimental measurements and analytical solutions for waves over a rigid canopy, demonstrating its capability to accurately capture the wave dissipation and flow velocity profiles, even with a relatively coarse grid. A truss-spring model is proposed to capture vegetation motion with substantial deflections, and is proven to be mathematically consistent with the governing equation for the flexible vegetation motion. It allows for explicit time integration with large time steps when dealing with highly compliant vegetation. The truss-spring model is validated and verified by experimental and numerical results for large-amplitude motions of a single elastic blade subjected to waves and sinusoidal oscillatory flows. The coupled model, combining the linearized flow solver and the truss-spring model, is applied to investigate wave propagating over a heterogeneous, suspended, and flexible canopy, showing high efficiency and good agreement with the experiments concerning wave attenuation and the hydrodynamic loads on the vegetation.

(Invited) On the Definition of Velocity in Statistically Accurate, Discrete-Time Simulations of Thermodynamics

One of the most basic struggles in simulations of statistical and dynamical systems is how to appropriately balance the desire for simulation accuracy, obtained in the small time step limit, with simulation efficiency, obtained for large time steps. Obviously, we would like to have the benefits of both limits of accuracy and efficiency. Thus, understanding the influence of discrete time on the behavior of equations of motion is crucial for the understanding and optimization of numerical simulations in physical science and engineering. There are two different types of crucial discretization errors. One is that the simulated, discrete-time trajectory increasingly deviates from the true, continuous-time trajectory as the time step is increased. Another is that the simulated discrete-time velocity increasingly deviates from the expected velocity of the corresponding simulated trajectory as the time step is increased. While the former error is a signature that the simulated trajectory is not exactly the intended one, the simulation may still be representative of the true system. In contrast, the latter error is a signature of a fundamental internal inconsistency in the simulation, which may greatly complicate physical interpretation of the simulation results, especially when both configurational and kinetic measures are required. It is therefore crucially important to understand the relationship between a discrete-time coordinate and its corresponding velocity.It has previously been shown that, unlike in continuous-time, discrete-time configurational dynamics does not need a discrete-time velocity variable. This realization allowed for the derivation of the complete configurational GJ set of stochastic Verlet-type algorithms for simulating the Langevin equation. This set has the key feature of simultaneously giving exact, time-step independent statistical sampling of damped, noisy harmonic oscillators, and giving exact spatial transport values, drift and Einstein diffusion, on flat surfaces. The apparent fact that a discrete-time velocity is not required to simulate a coordinate trajectory implies that discrete-time velocity is an ambiguous concept, which, in turn, implies that multiple velocity definitions can be entertained in order to avoid internal simulation inconsistencies between configurational and kinetic measures. From the derivation of the stochastic GJ methods it was shown that a statistically correct velocity obtained at the same time as the coordinate (an on-site velocity) cannot be obtained. In contrast, a velocity at the half-step between consecutive configurational time steps (a half-step velocity) was identified to produce exact statistical kinetic response for harmonic oscillators regardless of the applied time step. However, only for a single GJ method does this velocity also give the correct drift value.We here derive new, practical velocity definitions that are optimized for statistical kinetic measures, such that, e.g., measures of kinetic energy and correlations with the configurational coordinate are exact for the noisy harmonic oscillator, while also producing the correct kinetic measure of transport on flat surfaces. We identify and analyze three half-step velocity definitions that fulfill these criteria, and we outline simple and practical Verlet-type algorithms that incorporate these velocities within the statistically desirable GJ format that provides the optimized configurational statistics. The application of the methods are demonstrated through characteristic Molecular Dynamics simulations of both simple and complex systems, illustrating how the resulting methods are able to efficiently and simultaneously simulate configurational and kinetic statistical ensembles in different situations using very large time steps.

Adaptive conservative time integration for unsteady compressible flow

Unsteady flow computations typically rely on time integration schemes which employ the same step size everywhere. However, in cases with strong variations of wave speed and/or spatial resolution, local stability criteria and truncation errors would allow for much larger time steps in parts of the domain. To exploit this potential of improving the computational efficiency, adaptive time-stepping methods have been developed. Recently, an adaptive conservative time integration (ACTI) scheme for explicit finite volume methods was devised. Opposed to previous methods, ACTI guarantees exact conservation and periodic synchronization. Both are achieved by splitting a global time step into 2L (L≥0 is a cell dependent level) sub-steps, whereas the order of cell updates is critical. It has been demonstrated for flows in fractured porous media and for the Euler equations that ACTI can reduce the computational cost dramatically while preserving second order accuracy in space and time. In this paper, the ACTI scheme is extended to the compressible Navier-Stokes-Fourier system, and special attention is required for the spatio-temporal discretization of the convective and viscous fluxes. Numerical studies of unsteady flows involving shock boundary layer interaction demonstrate that the same accuracy is achieved with the new compressible ACTI flow solver as with classical time integration, but at much lower cost. Moreover, since the method is explicit, it has the potential for efficient computations on multiple CPUs and GPUs.

Improved discrete unified gas-kinetic scheme for interface capturing.

In this paper, we extend the improved discrete unified gas-kinetic scheme (DUGKS) from solving the hydrodynamic equationsto addressing the phase field equations, building upon our prior work [Wang et al., Phys. Fluids 35, 017106 (2023)10.1063/5.0128912]. The conservative Allen-Cahn equationand its modified form are presented first, followed by the construction of two improved DUGKS methods for interface capturing, based on the corresponding kinetic equations. The improved DUGKS for interface capturing utilizes the node distribution function instead of the interface center distribution function for evaluating the interface flux. The improved DUGKS enhances the numerical stability of the original DUGKS, and the good stability allows the calculations to be performed using large time steps, reducing the cumulative error from which more accurate predictions can be obtained. To verify the validity of the scheme, a series of numerical experiments were further carried out, including the diagonal translation, Zalesak's disk rotation, reversed single vortex, and deformation field. The comparison with the benchmark data shows that the improved DUGKS can simply and effectively capture the sharp interface and complex deformation interface of the two-phase flow interface.

Do not forget the electrons: Extending moderately-sized nuclear networks for multidimensional hydrodynamic codes

Nuclear networks are widely used coupled with hydrodynamical simulations of explosive scenarios to account for the change of nuclear species and energy generation rate due to nuclear reactions. In this way, there is a feedback mechanism between the hydrodynamical state and the nuclear processes. Unfortunately, the timescale of nuclear reactions is orders of magnitude smaller than the dynamical timescale that drives hydrodynamical simulations. Therefore, these nuclear networks are usually very small, reduced in most cases to a dozen elements, especially when simulations are carried out in more than one dimension. We present here an extended nuclear network, with 90 species, designed for being coupled with hydrodynamic simulations, which includes neutrons, protons, electrons, positrons, and the corresponding neutrino and anti-neutrino emission. This network is also coupled with temperature, making it extremely robust and, together with its size, unique of its kind. The inclusion of electron captures on free protons makes the network very appropriate for multidimensional studies of Type Ia supernova explosions, especially when the exploding object is a massive white dwarf. We perform several tests that are relevant to simulate explosive scenarios, such as Type Ia supernovae and core-collapse supernovae. We compare the results of the 90 nuclei network with a standard alpha -chain network with 14 elements to evaluate the differences in the energy generation rate. We also evaluate the relevance of including the electrons in the network in terms of generated yields and how it affects the pressure of a degenerate fluid such as that of white dwarfs. The results obtained with the 90-nuclei network have been verified with a much larger 2000-nuclei network built from REACLIB (WinNet), in terms of nuclear energy generation rate, pressure, and produced yields. The results obtained with the proposed medium-sized network compare fairly well, to a few percent, with those computed with in scenarios reproducing the gross physical conditions of current Type Ia supernova explosion models. In those cases where the carbon and oxygen fuel ignites at high density, the high-temperature plateau typical of the nuclear statistical equilibrium regime is well defined and stable, allowing large integration time steps. We show that the inclusion of electron captures on free protons substantially improves the estimation of the electron fraction of the mixture. Therefore, the pressure is better determined than in networks where electron captures are excluded, which will ultimately lead to more reliable hydrodynamic models. Explosive combustion of helium at low density, occurring near the surface layer of a white dwarf, is also better described with the proposed network, which gives nuclear energy generation rates much closer to than typical reduced alpha networks. A nuclear network with N=90 species, including electrons, aimed at multidimensional calculations of supernova explosions is described and verified. The proposed network is suitable for the study of Type Ia supernova explosions because it provides better values of pressure and electron abundance than other existing networks with smaller or even a similar size but without including electron capture processes.

Transient Friction Analysis of Pressure Waves Propagating in Power-Law Non-Newtonian Fluids

Modulated pressure waves propagating in the drilling fluids inside the drill string are a reliable real-time communication technology that transmit data from downhole to the surface during oil and gas drilling. In the analysis of pressure waves’ propagation characteristics, the modeling of transient friction in non-Newtonian fluids remains a great challenge. This paper establishes a numerical model for transient pipe flow of power-law non-Newtonian fluids by using the weighted residual collocation method. Then, the Newton–Raphson method is applied to solve the nonlinear equations. The numerical method is validated by using the theoretical solution of Newtonian fluids and is proven to converge reliably with larger time steps. Finally, the influencing factors of the wall shear stress are analyzed using this numerical method. For shear-thinning fluids, the friction loss of periodic flow decreases with the increase in flow rate, which is opposite to the variation law of friction with the flow rate for stable pipe flow. Keeping the amplitude of pressure pulsation unchanged, an increase in frequency leads to a decrease in velocity fluctuations; therefore, the friction loss decreases with the increase in frequency.