A comprehensive theoretical formulation has been developed to examine global and regional stability characteristics of a PHWR. The dynamic stability of a large Pressurized Heavy Water Reactor (PHWR) has been investigated through an integrated modelling framework based on the multipoint reactor kinetics formulation along with appropriate reactivity feedbacks. The model incorporates fuel and coolant temperature coefficients, xenon dynamics, and control system representation to capture the dominant neutronic–thermal interactions governing reactor core behavior. The resulting closed-loop representation of the PHWR explicitly accounts for the interaction between neutronics and the associated thermal hydraulics through reactivity feedbacks. An eigenvalue-based stability analysis has been performed to quantify the system response under varying reactivity feedback conditions. The resulting stability maps distinctly delineate stable and unstable operating domains in terms of key reactivity coefficients, identifying regimes associated with in-phase and out-of-phase oscillatory modes. Furthermore, a methodology has been proposed to characterize the onset of regional instabilities across different spatial modes. Analytical results demonstrate the ordered emergence of front-to-back, side-to-side, and top-to-bottom power tilts as the fuel temperature coefficient of reactivityαf tends to cross critical thresholds. The Multipoint reactor kinetics framework provides insight into the existence of diverse regional instability modes. Since, space time kinetics is not only computationally expensive but also lacks a systematic approach for identifying the ordered emergence of these instabilities, the proposed methodology provides a computationally tractable and systematic basis for determining stability margins and facilitating -oriented operational planning in large nuclear reactors.
This paper outlines an approach for carrying out a linear stability analysis of spatial xenon-induced power oscillations in a large pressurized heavy water reactor (PHWR) using the modal synthesis technique. The results obtained from the stability analysis are compared with those obtained by a solution using the improved quasi-static method. The present approach offers a computationally efficient linear stability assessment rather than obtaining the same using the improved quasi-static approach. In this approach, the space-time-dependent neutron flux is expressed as a product of the time-dependent amplitude function and the space-dependent eigenfunctions known as modes. This technique is used to solve time-dependent neutron diffusion equations, incorporating delayed neutrons, iodine, and xenon feedback effects. An eigenvalue-based linear stability analysis is performed for the first three harmonic modes, considered one at a time. Stability boundaries are generated in the parameter space defined by the operating power fraction and the power coefficient of reactivity, identifying both the stable and unstable regimes. Numerical simulations validated the stability boundary, showing growing, sinusoidal, and damped oscillations in the unstable boundary and stable regions, respectively. The model is verified against the space-time neutronics code IQS-3D, which uses the improved quasi-static method for the solution to the top-to-bottom and side-to-side modes for different operating power levels. The results indicated a stability threshold of similar to 46% full power (FP) for nominal reactivity device configurations, with some configurations stabilizing the xenon oscillations even at 100%FP. These findings highlight the key role of xenon dynamics in the determination of reactor stability and present a computationally efficient tool for mapping the operational stability boundaries in PHWRs. This model can be further used to perform nonlinear stability analyses.
Pulsating Heat Pipes (PHPs) exhibit chaotic thermo-hydrodynamic behavior arising from the coupling of phase change, pressure oscillations, and frictional effects. This study explores the route to chaos in a simplified threedimensional autonomous model of a PHP, governed by phase change and friction coefficients. A period-doubling cascade leading to chaos is identified through bifurcation analysis and validated by Fast Fourier Transform (FFT) and Lyapunov spectrum. The transition is traced through successive oscillatory regimes-first doubling from two to four, then to eight, and eventually becoming chaotic-capturing the classical period-doubling route. At each bifurcation, the second Lyapunov exponent touches zero, revealing a one-to-one correspondence between bifurcation points and dynamical instability. The dissipative behavior of the system is characterized through its divergence, analytically proportional to a friction parameter representing viscous-like damping. This damping term, embedded in the momentum equation, acts as a continuous energy sink that extracts kinetic energy from oscillatory motion and dissipates it as heat. The sum of Lyapunov exponents, calculated using the Wolf algorithm, matches the divergence of the vector field, highlighting the intrinsic link between energy dissipation and trajectory divergence in phase space. The universality of the transition is supported by the calculated Feigenbaum constant, which converges to approximate to 4.6692. These findings offer a deeper understanding of nonlinear dynamics and dissipation in PHPs, providing a theoretical foundation for future control and optimization of chaotic oscillations in thermal systems.
This study performs a weakly nonlinear analysis of a light-induced bioconvection model for a phototactic algal suspension confined between two horizontal plates. The lower boundary is rigid, while the upper is considered under both rigid and stress-free conditions, illuminated by collimated irradiation from above. The primary objective is to investigate the weakly nonlinear effects through bifurcation analysis and convective rolls, which significantly influence the observable structure and dynamics of the bioconvective system. The weakly nonlinear analysis of phototactic bioconvection is essential due to its critical role in optimizing biomass concentration and pattern formation in photobioreactors for advanced biofuel production. The analysis begins with linear stability theory to establish the critical conditions for the onset of convection from an infinitesimal perturbation. Progressing to a finite-amplitude perturbation, a weakly nonlinear analysis is employed, leading to the derivation of a Stuart-Landau equation that governs the amplitude of the convective mode. The sign of the Landau constant reveals whether the bifurcation is supercritical or subcritical. Our findings indicate that the bifurcation type, whether pitchfork or Hopf, as well as its supercritical or subcritical nature, are parametrically dependent on factors including critical light intensity, cell swimming speed, optical depth, and boundary conditions. Although supercritical bifurcations are common, a transition to subcritical bifurcations is observed at higher cell swimming speeds and optical depths. The dynamics are further elucidated using phase-space trajectories, time-evolution plots, and flow pattern visualizations.
This study presents a comprehensive numerical investigation of a four-dimensional non-autonomous memristive FitzHugh–Nagumo neuron model incorporating a generalized charge-dependent memristor (GCDM) and a magnetic flux-controlled memristor (MFCM) under periodic external forcing. Following a physically motivated circuit formulation and nondimensionalisation, the system is transformed into an autonomous representation through auxiliary phase variables, enabling systematic dynamical analysis. High-accuracy numerical integration combined with QR-based tangent-space propagation is employed to compute the complete Lyapunov spectrum and construct bifurcation diagrams through systematic variation of forcing amplitude and frequency. The dynamical analysis is further supported by vector-field divergence calculations and phase-space volume evolution obtained from ensemble propagation and four-dimensional convex-hull reconstruction. Numerical comparisons between phase-volume dynamics, divergence, and Lyapunov spectra demonstrate consistency with established results from dynamical-systems theory. In particular, the temporal mean divergence is found to closely match the sum of the Lyapunov exponents across the investigated parameter ranges, confirming the dissipative nature of the flow. The combined bifurcation–Lyapunov analysis reveals Feigenbaum-type period-doubling cascades, intermittent chaotic windows, reverse period-doubling transitions, and entrainment under coherent forcing. Computation of the complete Lyapunov spectrum enables identification of expanding, contracting, and neutral phase-space directions, providing a complete characterization of stability transitions than analyses based solely on the largest Lyapunov exponent. The results highlight the rich nonlinear dynamics arising from the interplay between neuronal excitability, memristive memory, and periodic stimulation.
Pulsating heat pipes have demonstrated great potential in advanced thermal management owing to their ability to passively transfer heat through self-sustained oscillations of the working fluid. These oscillations, primarily driven by phase-change phenomena at the evaporator, are inherently nonlinear and depend on the dynamic coupling between heat transfer, wall temperature distribution, and fluid motion. Unlike earlier simplified approaches that modeled phase change phenomena using an overall coefficient, the present study develops a nonlinear mathematical model that explicitly accounts for three dominant factors governing phase change, namely, phase-change limit, wall temperature gradient, and evaporation rate. By analyzing the steady-state behavior and temporal evolution of the meniscus position, the system is shown to undergo Hopf bifurcations for each contributing factor, indicating the onset of oscillatory motion. Through numerical continuation and time-series analysis, the combined variation of two factors revealed well-defined stability boundaries separating stable, unstable, and constant-amplitude oscillatory regimes. The results demonstrate that any combination involving the wall temperature gradient and phase-change limit promotes growing, self-sustained oscillations, whereas the evaporation rate alone fails to sustain periodic motion near the instability boundary. This systematic analysis provides a new theoretical framework to identify the operational regions required for continuous and thermally stable performance of PHPs. The study highlights that proper control of wall thickness and heat input is crucial for maintaining the phase-change limit and wall temperature gradient within the oscillatory domain, thereby ensuring safe, efficient, and self-sustained operation under transient conditions.
Recent experiments in MagnetoHydroDynamics (MHD) flow in heated duct reveal the presence of low frequency high amplitude temperature fluctuation, termed as magneto-convective fluctuations (MCF). To better understand the behaviour of MCF, numerical simulations have been carried out with various direction of applied magnetic field as well as temperature difference applied at different walls. The numerical methodology is validated against the experimental data available in the literature. The validated approach is then used for detailed numerical simulation. The different flow conditions lead to different type of convection rolls leading to two distinctive flow features. First, when the magnetic field is transverse to gravity, the convection rolls get aligned along the applied magnetic field direction, which would have aligned along the axis of the duct without the magnetic field. This results in observation of magneto-convective fluctuations. Second, when the magnetic field is parallel to gravity, the magnetic field slices the convection roll into multiple convection rolls aligned along the duct axis. This results in a wavy temperature pattern in the duct vertical cross section. Such phenomenon has not been reported in literature. The probable reason for such transformation of convection rolls have been discussed in detail.
The performance of Low-Temperature Proton Exchange Membrane Fuel Cells (LT-PEMFCs) is not constant, particularly in applications like vehicles that experience frequent load changes. Transient models explicitly capture these dynamics, unlike steady-state models, allowing for better design and control for optimal cell performance. This study introduces a 3D, multiphase, non-isothermal physics-based model of the LT-PEMFCs to comprehensively explore its transient behavior under various current profiles, temperatures, pressures, and humidity levels. Unlike the earlier studies, this model provides spatio-temporal insights into cell performance by incorporating liquid water dynamics for various flow configurations. Here, the primary objective is to evaluate the design and operating conditions that lead to the cumulative accumulation of liquid water inside the cell. Therefore, the performance characteristics of the multi-pass serpentine and pin flow configurations are compared under real-world driving scenarios such as the World harmonized Light vehicle Test Procedure. A comparison of experimental and simulated results at 75 degrees C and 80 % RH yields a root-mean-square error of 22 mV, which demonstrates the model's accuracy. At low RH and temperature (60 %, 65 degrees C), the serpentine configuration has delivered nearly a 10 % higher cell potential than the pin flow configuration. However, its average liquid water saturation has increased by 30 % due to sluggish liquid water dynamics, leading to higher possibilities of cumulative water accumulation inside the cell. This work offers valuable insights for identifying and mitigating potential issues such as water flooding and cell starvation.
Bottom-heated viscoelastic fluids in a cavity transit from conduction to convection through periodic oscillations or steady-states of flow patterns, depending on the Rayleigh number and other fluid parameters. It is in contrast to the Newtonian fluids, where the transition is always to steady-state convection. The trapezoidal cavities filled with Oldroyd-B fluid have been explored with side walls inclined from to using OPENFOAM-based RheoTool simulations. The obtuse angle trapezoidal cavities have five different types of solutions. Four of these solutions consist of one-roll and two-roll solutions (TRS) with and without oscillations. The fifth solution is the conduction-dominated solution with low flow and heat transfer. However, since it lacks two-roll solutions, there are only three kinds of solutions for Rayleigh-B & eacute;nard convection (RBC) in an acute angle trapezoidal cavity (AATZC). The bifurcation maps and heat transfer characteristics are presented for various types of rolls. One-roll periodic solution exists beyond the viscosity ratio of 0.5 in AATZCs, while it is up to in square and obtuse angle trapezoidal cavity (TZCs). Moreover, two-roll periodic solutions are observed beyond the viscosity ratio of in obtuse angle TZCs. Isotherms and flow patterns are presented to illustrate the dynamics of each flow regime. The effect of sidewall angles on the stabilization of the flow is also investigated. The periodic oscillations are observed up to higher Rayleigh numbers for trapezoidal cavities with smaller cavity angles compared to those with higher cavity angles.
Natural circulation is commonly accounted for and designed into thermal fluid systems such as nuclear reactors and concentrated solar thermal plants as part of safety systems and normal operation. This work presents the development of a Fourier-based model of single-phase natural circulation loops with application-relevant boundary conditions to study their performance, stability, and response to geometry-induced turbulence fluctuations. This results in a reduced order model consisting of eight ordinary differential equations that reproduces the phenomena observed in experiments and numerical simulations. This model enables a robust study of the stability and dynamics of the system. The model results compare well against past experiments for steady-state conditions and instability predictions. Supercritical and subcritical Hopf bifurcations are obtained, and a chaotic attractor is observed in line with previous predictions and studies. The effect of transient fluctuations arising from complex geometry in the system is modeled using a data-driven statistical emulator with the help of a modified minor loss coefficient or form friction parameter. The effect of the stochastic friction parameter or stochastic forcing on the model predictions was analyzed. Under certain conditions, the stochastic forcing considerably shrinks the region of stability as the system is perturbed from a metastable state. The stochastic forcing is also found to accelerate the transition to chaos. For cases that do not escape to the chaotic attractor, the reaction of the system to the stochastic parameter depends on the response time of the attractor and its limiting behavior.
Nodal integral methods (NIMs) have been proven effective in solving a wide range of scientific and engineering problems by providing accurate solutions with coarser grids. Despite notable advantages, these methods have encountered limited acceptance within the fluid flow community, primarily due to the lack of robust and efficient nonlinear solvers for the algebraic equations arising from discretization using NIM. A preconditioned Jacobian-free Newton-Krylov approach has been recently developed to solve Navier-Stokes equations to overcome this limitation. The developed approach has extended the acceptability of NIM and demonstrated considerable gains in computational time. However, a challenge persists in the efficiency of the proposed approach, particularly in solving the pressure Poisson equation. Addressing this, we offer novel strategies and algorithms to solve the pressure Poisson equation. These strategies aim to improve the computational efficiency of NIMs, making them more effective in solving complex problems in scientific and engineering applications.
An efficient coarse-mesh nodal integral method (NIM), based on cell-centered variables and termed the cell-centered NIM (CCNIM), is developed and applied to solve multi-dimensional, time-dependent, nonlinear Burgers equations, extending the applicability of CCNIM to nonlinear problems. To overcome the existing limitation of CCNIM to linear problems, the convective velocity in the nonlinear convection term is approximated using two different approaches, both demonstrating accuracy comparable to or better than traditional NIM for nonlinear Burgers problems. Unlike traditional NIM, which utilizes surface-averaged variables as discrete unknowns, this innovative approach formulates the final expression of the numerical scheme using discrete unknowns represented by cell-centered (or node-averaged) variables. Using these cell centroids, the proposed CCNIM approach presents several advantages compared to traditional NIM. These include a simplified implementation process in terms of local coordinate systems, enhanced flexibility regarding the higher order of accuracy in time, straightforward formulation for higher-degree temporal derivatives, and offering a viable option for coupling with other physics. The multi-dimensional time-dependent Burgers problems (propagating shock, propagation, and diffusion of an initial sinusoidal wave, shock-like formation) with known analytical solutions are solved in order to validate the developed scheme. Furthermore, a detailed comparison between the proposed CCNIM approach and other traditional NIM schemes is conducted to demonstrate its effectiveness. The proposed approach has shown quadratic convergence in both space and time, i.e., O[(Δx)^2, (Δt)^2], for the considered test problems. The simplicity and robustness of the approach provide a strong foundation for its seamless extension to more complex fluid flow problems.
Nodal integral methods (NIMs) are a class of highly efficient coarse-mesh techniques for the numerical solution of partial differential equations (PDEs). Among these, the cell-centered NIM (CCNIM) stands out as an enhanced version, proving its effectiveness in addressing fluid flow problems. However, it faces certain limitations, such as its unsuitability for one-dimensional problems and the necessity of using a complex system of differential–algebraic equations (DAEs) to represent discrete unknowns per node in transient problems. To address these challenges, this study introduces a modified version of CCNIM for the solution of the transient diffusion equation. We discretize the spatial and temporal coordinates following the principles of the nodal method, achieving second-order accuracy in both spatial and temporal variables. Unlike the previous version of CCNIM, our proposed scheme uses algebraic equations to represent discrete variables at each node, eliminating the need for complex DAE systems. To validate the proposed method, we solved several transient diffusion problems in one dimension with known analytical solutions. The simplification of the proposed scheme offers a strong basis for its smooth expansion to more complicated fluid flow problems.
We report experiments on the dynamics and characteristics of surface waves generated by the low Weber number impact of immiscible silicone oil droplets on a water pool. The study employed high-speed videography and background-oriented schlieren (BOS), in tandem, to capture the droplet impact and spatiotemporally resolved interfacial topography and surface waves, respectively. In addition to the applications of fundamental research interests, the importance of this study lies in the generation of whole-field experimental data, which can be used to develop and validate models for surface waves. The investigation explored effects of droplet viscosity, Weber number, and pool height on the characteristics of generated waves. The impact parameters led to various modes of droplet-pool interactions, including bouncing, coalescence, secondary jet formation, droplet pinning, and lens formation. The study revealed that the wave amplitude increases with the Weber number, while the wave velocity is independent of the wave amplitude and droplet Weber number. The wave velocity was found to be lower for thin films and higher for shallow and deep pools, and it also increases with the droplet viscosity for shallow and deep pools. The study identified three distinct types of waves: primary, leading, and trailing waves, each with different velocities and amplitudes. In comparing the experimental wave velocities with the theoretical model, it was observed that the experimentally measured wave velocity corresponded to the different critical velocities observed at corresponding pool heights. The velocity of the primary wave aligns with the mean phase velocity, while the velocities of the leading and trailing waves correspond to the group velocity predicted by the wave dispersion equations.
Rayleigh-Benard convection in square closed cavities filled with Oldroyd-B fluid was studied using OpenFOAMbased RheoTool. For the RBC in Newtonian fluids, the transition always occurs from conduction to steady state convection with increasing Rayleigh number (Ra). On the other hand, the viscoelastic fluids may also show the transition from conduction to oscillatory convection. Further increase in Ra may result in a steady state convective solutions. It is further noted that the behavior is similar to Newtonian fluids for larger values of viscosity ratio (B). Considering the abovementioned different flow behavior at different values of the parameters, it is noted that there are five different types of solutions possible for the viscoelastic fluids viz. pure conduction (PC), one roll periodic oscillations (ORPO), one roll steady state (ORSS) convection, two roll periodic oscillations (TRPO), simultaneous one and two roll steady state convection. Therefore, a bifurcation diagram in the parametric space of Ra and B is presented, depicting these five regions corresponding to each type of solution. The boundaries of these regions have been identified by numerical simulation. Note that all these regions exist in the laminar flow regime, and the transition to turbulence is not considered here. Interestingly, at low values of B, as one increases Ra, it is seen that the ORSS region is sandwiched between ORPO and TRPO. The likely reason for this interesting behavior is explained. Moreover, representative solutions in each region in terms of isotherms, streamlines, and vector plots have been included to demonstrate the dynamics of each delineated region.
The application of Nodal Integral Method (NIM) is mostly limited to the low Reynolds number problems; due to lack of suitable nonlinear solvers. This necessitates the development of advanced nonlinear solvers for higher nonlinearity. Newton’s method is a well-established solver for such equations. However, the rate of convergence in Newton methods is dependent on the initial guess. Additionally, the condition number of Jacobian matrix dictates the time required for matrix inversion. In this work, for the first time a preconditioner is developed for fluid flow problem using Modified-NIM (MNIM). This preconditioned Jacobian free Newton Krylov algorithm uses the linearized Modified-MNIM (M2NIM) as the preconditioner, resulting in better eigenvalue clustering, reducing Krylov iterations. The effectiveness of proposed method is demonstrated using 1-D Burgers' equation for larger time steps and higher Reynolds number up to 2500. The proposed preconditioner drastically reduces the spectral radius, CPU run-time and Krylov iterations.