Darcy–Benard–Oldroyd convection in anisotropic porous layer subject to internal heat generation

An anisotropic horizontal porous layer saturated with viscoelastic liquids of the Oldroyd-B type is explored to determine how the internal heat source affects thermal convection. As a momentum equation, a modified Darcy–Oldroyd model is used that takes into account the anisotropy of the porous layer. The energy equation is formulated in such a way that the influence of internal heat sources and anisotropy in thermal diffusivity on the stability criterion may be easily identified. The effects of anisotropy, viscoelasticity, and internal heat generation on the onset of thermal convection are investigated using linear stability analysis. It is understood that convection begins via an oscillatory mode instead of a stationary mode because viscous relaxation, thermal diffusions, and internal heat generation mechanisms compete with one another. Both steady and unsteady finite-amplitude convections are studied using nonlinear stability analysis with the truncated Fourier series method. The effect of different governing parameters on the system’s stability and on convective heat transfer is studied. The present investigation has been significantly validated by the recovery of several prior results as special situations. The findings presented in this work are anticipated to have significant implications for a number of real-world applications, including modeling of oil reservoirs, crude oil extraction, crystal growth, the pharmaceutical and medical industries, and the use of geothermal energy, among others.

research is the study of the viscoelastic properties of the asthenospheric and mantle components of the Earth (see Lowrie, 2020). Snow systems and the rheology of food transport involve viscoelastic liquids saturated in porous media. Less research has been conducted on Rayleigh-Benard convection in liquids with viscoelastic properties than on thermal convection in a Newtonian fluid (refer to Li and Khayat, 2005).
Due to their coagulated viscosity, polymeric liquids are unaffected by flow and turbulence, whereas Newtonian fluids are affected. Moreover, oscillatory convection is believed to exist as viscoelastic fluids are characterized by their elasticity. According to Kim et al. (2003), the system is unstable due to its elasticity but stable due to its porosity. The oscillatory mode is appropriate for studying convection because convection reduced to supercritical and stable bifurcation forms does not vary with elasticity. However, few authors have investigated oscillatory convection in viscoelastic liquids, and many researchers have not addressed a comparable problem in porous media. Thermal convection in viscoelastic fluids was studied by O' Connell and Budiansky (1977), Griffiths (1987), Rudraiah et al. (1989), Yoon et al. (2004), Malashetty et al. (2006a), and Swamy et al. (2012).
When a porous medium exhibits anisotropy in its mechanical and thermal properties, the flow behavior is significantly altered. The assumption underlying Darcy's model is that the fluid flow is sufficiently slow for inertial effects to be negligible, and that the fluid complies with the laws of continuum mechanics. If the porous medium is anisotropic, then the medium's permeability will depend on the flow direction. In other words, the permeability of the porous medium will vary in different orientations. Darcy's law must be modified to account for the anisotropy of the medium in this instance. This can be accomplished by introducing a permeability tensor that describes the medium's permeability in various orientations. Numerous applications in the actual world involve anisotropic porous media. Anisotropy results from the incorrect orientation of a solid matrix or the asymmetries of the natural porous medium; it is also a characteristic of porous synthetic materials, such as the fibrous material utilized for insulation and pelleting, both of which are beneficial to chemical engineering processes. McKibbin et al. (1985) and Storesletten (1998) have exhaustively documented the work pertaining to convection in anisotropic porous media. Govender (2006) described the effects of anisotropy on thermal convection in a porous layer. Malashetty et al. (2006b)presented a study on the effect of a time-periodic modulated temperature field on the stability of a viscoelastic fluid saturated within an anisotropic porous layer. Saravanan and Purusothaman (2009) studied non-Darcian effects in anisotropic porous media. Many pertinent studies on anisotropy have been conducted (Malashetty and Swamy, 2007a;Malashetty and Swamy, 2007b;Malashetty et al., 2009;Sivakumar and Saravanan, 2009;Malashetty and Swamy, 2010;Agarwal et al., 2011;Malashetty et al., 2011;Srivastava et al., 2011;Chandra and Satyamurty, 2012;Swamy et al., 2013;Swamy et al., 2014;Swamy, 2017;Swamy et al., 2019;Capone et al., 2020).
Understanding thermal convection in a porous stratum with internal heat generation is crucial in a variety of natural and artificial systems. The interior of the Earth provides an example of internal heat production. Due to the disintegration of radioactive isotopes and residual heat from the planet's formation, the Earth's core is thought to be the source of significant internal heating. This heat production is responsible for the elevated temperature of the Earth's interior and is one of the driving forces behind plate tectonics and volcanism. Electronic devices, such as computer processors, are another example of internal heat generation. A computer processor's electronic components generate heat due to the passage of electric current through them. If the generated heat is not effectively dissipated, the processor may malfunction or even fail. To prevent this, electronic devices are frequently equipped with cooling systems, such as heat sinks and fans, to dissipate heat and maintain secure component temperatures. Chemical reactions can also generate heat internally. Due to the exothermic nature of the combustion reaction, for example, heat is produced when fuel is consumed. This internal heat production can be utilized to generate electricity in power facilities or to heat industrial processes.
Internal heat generation results in a nonlinear temperature gradient. A temperature gradient is the temperature change over a specified distance. A nonlinear temperature gradient indicates that the change in temperature over a given distance is not constant. In other terms, the change in temperature is not linear. This refers to the fact that the system or material in question generates its own heat rather than receiving it from an external source. Therefore, thermal convection occurs even though the temperature difference between the lower and upper surfaces is insufficient for convection to begin. Consequently, the production of internal heat is an additional essential mechanism for regulating the onset of convection in the porous layer. Numerous researchers have conducted extensive research on thermal convection in absorbent layers with internal heat generation [refer to (Thirlby, 1970;Mahabaleshwar et al., 2017;Ahmed and Rashed, 2019;Yadav et al., 2021;Enagi et al., 2022;Raju et al., 2022;Upadhya et al., 2022)]. Now, let us see the physical mechanisms involved in the present problem. As a result of the internal heating, a temperature gradient is generated within the porous layer, with the temperature being higher near the heat source and decreasing toward the top surface of the layer. This temperature gradient creates density differences within the viscoelastic liquid, causing it to flow from the heat source toward the top surface. The fluid flow of the viscoelastic liquid is influenced by the anisotropy of the porous layer. If the porous layer is more permeable in one direction than in other directions, the fluid flow will be preferentially directed along the more permeable direction, resulting in anisotropic convection. The properties of the porous medium, such as its porosity, permeability, and tortuosity, also affect the rate of fluid flow and the way in which heat is transferred through the medium. For example, a higher porosity of the porous layer can increase the flow rate of the viscoelastic liquid, while a higher tortuosity can hinder the flow, leading to a slower rate of fluid motion. The viscoelasticity of the liquid can also influence the way in which fluid flow and heat transfer occur within the porous layer. A viscoelastic liquid has both viscous and elastic properties, which can lead to nonlinear behavior in fluid flow and heat transfer. For instance, the elasticity of the liquid can cause it to exhibit oscillatory behavior in response to the temperature gradient, resulting in oscillatory fluid flow and heat transfer. The heating of the bottom surface by an external means can also affect the flow of the viscoelastic liquid. For instance, uniformly heating the bottom surface can generate a temperature gradient that interacts with the temperature gradient generated by the internal heat source, resulting in a more complex flow pattern.
A stability analysis of an Oldroyd-B fluid in a porous medium with the combined influence of anisotropy and the internal heat source is rare. Such investigations are still greatly desired. The Frontiers in Materials frontiersin.org 02 primary objective of this study is to determine how the viscoelastic parameters, internal heat generation coefficient, and anisotropy parameters impact the onset criterion of thermal convection and the heat transfer across the layer.

Mathematical formulation
A shallow horizontal anisotropic porous layer saturated with Oldroyd-B liquid is considered. The surfaces held at z 0 and z h are regarded as being stress-free and isothermal. The gravitational force g ≡ (0, 0, −g) is acting downward in the direction of the z-axis. The adverse temperature gradient ΔT between the two surfaces is maintained by heating the lower surface uniformly. Internal heat generation is considered as an additional source. The temperatures of the solid and liquid phases are assumed to have reached equilibrium. The conservation law of linear momentum is represented by a modified Darcy-Oldroyd model incorporating local time derivatives, Boussinesq approximation, and anisotropy. The convection velocities are expected to be negligible. Thus, the effects of Forchheimer inertia and advection are disregarded. Consequently, the pertinent mathematical model is where denotes the velocity vector, v D εv is the seepage velocity, Λ 1 is the stress-relaxation time, Λ 2 is the strainretardation time, ε is the porosity, μ is the viscosity, ρ is the density, α is
Frontiers in Materials frontiersin.org the thermal expansion coefficient, γ (ε(ρc) f + (1 − ε)(ρc) s )/(ρc) f is the ratio of specific heat capacities, K (e 1 e 1 + e 2 e 2 )K −1 1 + (e 3 e 3 )K −1 3 and κ (e 1 e 1 + e 2 e 2 )κ 1 + (e 3 e 3 )κ 3 denote the anisotropy-induced permeability and thermal diffusivity tensor, respectively, with e 1 , e 2 , e 3 being the unit vectors along the x, y, and z-axes, respectively. The last term in Eq. 3 is due to the presence of internal heat generation, where Q I quantifies the amount of heat generation within the bulk of the porous layer Because the fluid is at rest in the basic state, we can determine its mass, pressure, and temperature by  Frontiers in Materials frontiersin.org 05 Swamy et al. 10.3389/fmats.2023.1158644 The stability of the system is studied by imposing the perturbations on the basic state.
where prime represents the quantity in the perturbed state. We use Eqs. 5-8 in the governing Eqs. 1-4 and eliminate p′. Furthermore, assume the flow to be two-dimensional and thus incorporate the stream function such that (u′, w′) (z/z z , −z/zx)Ψ′. Then, use (x, z) h(x*, z*), t (h 2 /k 3 )t*, Ψ k 3 Ψ* and T′ (ΔT)T* to express the equations in the nondimensional form as follows: where f(z) Ra I √ cos( Ra I √ (1 − z))/sin( Ra I √ ), ∇ 2 1 z 2 /zx 2 + (1/ξ)(z 2 /zz 2 ), and ∇ 2 2 η(z 2 /zx 2 ) + z 2 /zz 2 . Ra ρ 0 αgΔTh 3 /μk 3 denotes the thermal Rayleigh number, Ra I Qh 2 /k 3 is the internal Rayleigh number, Pr μh 2 /ερ 0 k 3 K 3 is the Darcy-Prandtl number, λ 1 Λ 1 k 3 /h 2 is the Deborah number or stress-relaxation parameter, λ 2 Λ 2 k 3 /h 2 is the strain-retardation parameter, and ξ K 1 /K 3 and η k 1 /k 3 denote mechanical and thermal anisotropy parameters, respectively. Because the boundaries are assumed to be stress-free and isothermal, the relevant boundary conditions are 3 Galerkin weighted-residual technique In the linear stability theory, the imposed perturbations are anticipated to be infinitesimal, and hence, the nonlinear terms in Eqs. 9-11 can be ignored. The normal mode analysis is used to solve the resulting eigenvalue problem. Thus, it is supposed that the infinitesimal perturbations are the periodic waves of form where the real variable k is the wavenumber, and the complex variable ω is the temporal growth rate. It decides whether these tiny-recurrent-perturbations either enlarge or degenerate. In thermal convection, wavenumber is a term used to describe the spatial variation of temperature and fluid flow patterns that arise due to temperature differences within the fluid. Wavenumber is a measure of the number of waves that occur in a given distance. The concept of wavenumber is closely related to the concept of wavelength, which is the distance between successive peaks or troughs of a wave. The wavelength and wavenumber are related by the formula k = 2π/λ, where λ is the wavelength of the wave. This shows that the wavenumber is inversely proportional to the wavelength; that is, waves with shorter wavelengths have larger wavenumbers, and waves with longer wavelengths have smaller wavenumbers.
In natural convection, the spatial variation of temperature and fluid flow patterns can take on different wavelengths, depending on factors such as the geometry of the system, the temperature differences within the fluid, and the properties of the fluid itself. The wavelength of these patterns can be quantified in terms of the wavenumber, and the behavior of natural convection can be analyzed in terms of the relationship between the wavenumber and other parameters such as the Rayleigh number.
According to the Galerkin method, we choose where A 1 and B 1 are constants, and Ψ 1 (z), Θ 1 (z) are so designated that they satisfy the boundary conditions (15). On multiplying Eqs. 13, 14, respectively, by Ψ 1 (z), Θ 1 (z) and integrating the resultant equations w. r. t. z between the limits 0 and 1, we acquire where X 4 〈f(z)Θ 1 Ψ 1 〉, X 5 〈Θ 2 1 〉, and X 6 〈Θ 1 (D 2 − ηk 2 )Θ 1 〉. The angular brackets denote the integral w. r. t. z between the limits 0 and 1. The requirement for the existence of a non-trivial solution of Eqs. 17 and 18 yield the expression for Rayleigh number: Because ω ω r + iω i , the case ω r < 0 specifies a stable state, while ω r > 0 corresponds to the unstable mode. A neutral state is attained for ω r 0. The steady-marginal stability can be observed for ω 0 (i.e., ω r ω i 0), with which Eq. 19 abridges to the expression of stationary Rayleigh number: The trial functions that satisfy the boundary conditions (15) are obviously On using these, one can estimate X 1 , X 2 , X 4 , X 5 and X 6 value and then substituting into Eq. 20 to get where δ 2 1 k 2 + π 2 /ξ and δ 2 2 ηk 2 + π 2 . Eq. 21 is free from viscoelasticity, so it resembles the equation obtained for a viscous Newtonian fluid. The validity of our work can be ascertained through Table 1, wherein we recovered the previous classical results as a special case of Eq. 21 Now, we discuss the behavior of Eq. 19 with the nonzero growth rate, that is, ω ≠ 0. As Ra portrays a physical phenomenon, it should not be imaginary and hence, Eq. 19 admits N i 0 as ω ≠ 0. This affords the expression forω 2 : The real part of Eq. 19 then symbolizes the expression for the oscillatory Rayleigh number: To estimate Ra Osc c , we minimize (23) w. r. t. k, after substituting ω 2 ( > 0) from Eq. 22.

Weak nonlinear theory
Nonlinear stability analysis is preferred to measure convection amplitudes and heat transfer. This facilitates comprehension of the physical mechanism with a little mathematical labor. This is a basic step toward grasping the full nonlinearity of the problem. Because the perturbations are assumed to be of finite amplitude, it is reasonable to represent them in the form of a limited Fourier series, as follows: The finite amplitudes of A and B subscripts for unsteady nonlinear convection are to be assessed by the dynamics of the system. Using Eqs. 24 and 25 in Eqs. 9 and 10 and comparing the coefficients of like terms, the following fourth-order Lorenz system of autonomous nonlinear differential equations is obtained: Frontiers in Materials frontiersin.org where D 1 −Pr δ −2 λ −1 1 δ 2 1 A 11 + λ 2 G 1 ( )+ δ 2 Pr −1 G 1 + RakB 11 + RakG 2 , G 2 −γ −1 kf z ( )A 11 + δ 2 2 B 11 − Ra I B 11 + 2πkA 11 B 02 , D 2 −γ −1 4π 2 − Ra I B 02 − πkA 11 B 11 2 .
There is no suitable analytical method to obtain a closed-form solution of the aforementioned system. Thus, a competent numerical technique is recommended. Although it may not be possible to make precise quantitative predictions, there are several qualitative insights that can be gleaned from the available data or theoretical models. As system of Eq. 26 is homogeneously bounded with time, it retains numerous features of the full problem. For Ra I < (0.5)δ 2 2 + 2π 2 , the velocity field possesses a constant negative divergence; that is, This implies that the system is constrained. In dynamical systems theory, a system is said to be constrained if it is subject to certain limitations or conditions. These constraints can arise due to physical, mathematical, or other reasons. For example, a mechanical system may be constrained by rigid walls or other physical barriers that restrict the motion of its components. When a system is constrained, the possible states that it can occupy are limited to a subset of its phase space. The phase space is the space of all possible states of the system, and it is often represented as a high-dimensional space in which each dimension corresponds to a particular variable or parameter of the system. Because the system is constrained, the phase space paths are drawn toward a set of measure zero or a fixed point. A set of measure zero is a subset of the phase space that has zero volume or probability. In other words, it is a set of states that is extremely unlikely to be reached by the system. A fixed point, on the other hand, is a state of the system that does not change over time. When the phase space paths are drawn toward a set of measure zero or a fixed point, this can result in volume shrinkage in the phase space. This means that the volume of the phase space that is accessible to the system becomes smaller over time as the system is forced to occupy a smaller subset of the phase space due to the constraints. This is revealed by Eq. 28 through the fact that if a set of preliminary points in space fills a volume V(0) at t = 0, then after time t, the endpoints of the corresponding paths will occupy a region The larger values of the Darcy-Prandtl number and strainretardation number and smaller values of the Deborah number and heat-capacities ratio are used to emphasize the exponential deterioration of volume with time.
Upon switching from qualitative exploration, we now look into the existence of an analytical solution. As the finite amplitude instability can be well explored analytically using the truncated model shown in Eqs. 24 and 25, a closed-form solution of Eq. 26 is used for the steady case. The following expression is obtained by Eq. 26 after removing all terms except A 11 : δ 2 1 − Rak 2 f z ( ) δ 2 2 − Ra I + 2π 2 k 2 π 2 − Ra I 4 A 2 11 8 −1 A 11 0.
The solution A 11 0 symbolizes the pure conduction state. Thus, the second option guarantees the likelihood of finite amplitude steady convection by offering the value of the finiteamplitude A 2 11 /8 of convective motions in the form A 2 11 8 π 2 − Ra I 4 2π 2 k 2 Rak 2 f z ( )δ −2 1 − δ 2 2 − Ra I . (29) Rather than merely defining the onset criterion, the impact of the Rayleigh number can be swiftly documented by analyzing its effect on heat transport. In the study of convection, determining the quantity of heat transported past the layer is of utmost importance. Because the basic state is immobile, heat transfer in this state is limited to conduction. However, as the Rayleigh number exceeds the threshold, convection develops. The Nusselt number is used to describe the convection-heat transport throughout the stratum.
This analysis is valid for Rayleigh numbers close to their threshold value. By including more terms in the Fourier expansion, one can anticipate better outcomes. The Runge-Kutta-Gill method is used to clarify the unsteady Eq. 26. The calculated amplitude values for various time intervals are then used to estimate Nu as a function of t.

Results and discussion
The onset of convection refers to the point at which fluid motion due to buoyancy forces begins to occur in a fluid that is initially at rest. This occurs when a fluid is subjected to a temperature difference that is large enough to cause density variations within the fluid, which in turn generate buoyancy forces that drive fluid motion. The onset of convection can be characterized by a critical value of a dimensionless parameter called the Rayleigh number. The Rayleigh number represents the ratio of buoyancy forces to viscous forces in the fluid. For a fluid layer that is initially at rest, the onset of convection occurs when the Rayleigh number exceeds a critical value that is specific to the geometry and boundary conditions of the system. Above this critical value, the buoyancy forces overcome the viscous forces and initiate fluid motion. The critical Rayleigh number can be determined through theoretical analysis or experimental observation. The onset of convection is an Frontiers in Materials frontiersin.org 08 important phenomenon in many natural and industrial systems, including geophysical flows, crystal growth, and industrial heat transfer processes. Understanding the onset of convection is important for predicting the behavior of these systems and for designing efficient heat transfer systems.
The primary objective of investigating a convection problem is to determine the smallest possible Rayleigh number that demonstrates convection. Fine-tuning the parameters that define the Rayleigh number is advantageous for deferring or accelerating convective motions. This study examines the interaction between variable permeability, thermal diffusivity, and internal heat generation. The point (k c , Ra c ) at which the marginal stability curve reaches the least signifies the convection threshold. The detailed behavior of this critical value is deliberated as a function of the strain-retardation number. Figures 1A-D display the dependence of critical Ra on the strain-retardation parameter λ 2 . It has already been mentioned that Ra St c is independent of viscoelasticity. This fact is justified by a horizontal dotted line (corresponding to Ra St c ) in Figure 1A, which is invariant with respect to λ 2 , λ 1 and located at the higher level. All the Ra Osc curves lie at the lower level, indicating that the onset of convection is through oscillatory mode. It is observed that Ra Osc c increases with λ 2 . Thus, the strain-retardation parameter makes the system more stable, but this stabilizing effect is retarded by λ 1 because there is a significant decrease in Ra Osc c w.r.t. λ 1 . We note that the influence of λ 2 on Ra Osc c is constrained to a specific range depending on the value of λ 1 . Beyond this range, the oscillatory convection ceases to occur.
In Figure 1B, the effect of the mechanical anisotropy parameter, which signifies heterogeneity in the permeability of the porous layer, is expressed. The values of both Ra St c and Ra Osc c decline with ξ. This indicates that the onset of convection can be advanced by increasing the anisotropy in permeability. Figure 1C exhibits the impact of anisotropy in thermal diffusivity on the threshold of convection. It portrays that by choosing larger values for η, one can enhance the values of Ra St c and Ra Osc c . Hence, stabilization can be achieved through increasing η. Figure 1D depicts the variation of Ra c with the internal Rayleigh number. Convective motions in both steady and oscillatory modes vary considerably with the aid of internal heating. This is revealed through the fact that when Ra I is increased, the Ra c curves are shifted downward in both stationary and oscillatory cases. This confirms that internal heating favors the onset of convection.
One of the main objectives of the present paper is to analyze the significance of controlling the convection by the anisotropic nature of the porous layer. The detailed behavior of Ra Osc c as a function of anisotropy parameters is demonstrated via Figure 2A. The black solid curve with regards to the left-side axis shows that the value of Ra Osc c decreases with ξ. Therefore, the anisotropy in permeability causes the oscillatory convective motions to occur at the earlier stage. However, note that this destabilization is more sensitive for the small and moderate values of ξ. Furthermore, the red curve drawn with reference to the right-side axis makes us aware of the stabilizing role of anisotropic thermal diffusivity. The value of Ra Osc c increases almost linearly with η.
The impact of Ra I , which signifies the strength of the internal heat source, is expressed through the curves in the (Ra c , Ra I ) plane (see Figure 2B). Both Ra St c and Ra Osc c plummet with Ra I . This indicates that onset of convection can be brought forward by increasing the rate of internal heating. The figure also exhibits that these critical curves are shifted upward when we increase the thermal anisotropy parameter. So, the heterogeneity of thermal diffusivity retards the destabilization caused by an internal heat source.
To explore the activity of Ra c with respect to the varying Prandtl number, we plot the critical curves in the (Ra c , Pr) plane in Figure 2C. For smaller values of Pr, there is a swift decrease in the value of Ra Osc c . This trend continues up to some specific value of Pr, beyond which Ra Osc c becomes independent of Pr. Thus, a destabilizing effect is reported for smaller values of Pr. The dotted horizontal line, which corresponds to the stationary case, shows that Ra St c is not influenced by the varying Prandtl number. Furthermore, the value of Ra St c is much larger than Ra Osc c . Another fact that can be noticed through this figure is the enhancement in the values of Ra Osc c with respect to the increasing heat capacities ratio. Therefore, the destabilization caused by increasing Pr can be suppressed by γ. Linear stability analysis, which provided us with a glimpse of the convection threshold, is followed by weak nonlinear stability analysis.
This study helped us to measure the amplitude of convective motions and the amount of heat transfer. The Nusselt number (Nu) that signifies the extent of convective heat transfer is calculated. The control of Rayleigh number over Nu is presented in Figures 3A-C. The value of Nu is found to be 1 at the onset, while as the Ra c increases to about three times the value of the critical Rayleigh number, the Nusselt number also increases. Thereafter, Nu becomes less sensitive to Ra c . Thus, one can conclude that in the vicinity of the onset of convection, the enhancement of heat transfer occurs, and the same magnitude is maintained even after the further increase in Ra c .
From Figures 3A,B, Nu is found to upsurge for the higher Ra I and ξ. Therefore, heat transport is amplified by introducing an internal source within the porous layer and by choosing the anisotropy in permeability. Furthermore, through Figure 3C, we notice that there is a significant reduction in the value of Nu w.r.t. the thermal anisotropy parameter. The heat transfer across the layer decreases considerably with increasing anisotropy in thermal diffusivity.
In unsteady finite-amplitude analysis, the fourth-order Lorenz model has been solved numerically using the Runge-Kutta-Gill method. The amplitudes are obtained as the functions of t and then substituted into the expression of Nu. In general, the Nusselt number is a function of several parameters, including the fluid velocity, temperature, and physical properties of the fluid, and it can vary with time during transient heat transfer processes. The behavior of the heat transfer coefficient with respect to time is displayed through the curves in the (Nu, t) plane, as shown in Figures 4A-D. One can observe from these figures that at the onset of natural convection, heat transfer initially occurs through conduction alone, and Nu has a value of 1, which corresponds to purely conductive heat transfer. As the fluid flow starts to develop, convective heat transfer becomes dominant, and the Nusselt number becomes sensitive to time. The behavior of the Nusselt number during this transient phase is oscillatory, with the value of Nu fluctuating about a mean value. This means that the heat transfer coefficient, which is related to the Nusselt number, varies periodically in time, with some oscillations. However, as time progresses, the oscillatory behavior of the Nusselt Frontiers in Materials frontiersin.org 09 number starts to decay, and after a short period of time, the heat transfer process reaches a steady state. At this point, the Nusselt number reaches a mean value that is similar to the value of Nu obtained in the steady-state, finite-amplitude case. Thus, in general, these figures display that the behavior of the Nusselt number during a transient heat transfer process is initially oscillatory and sensitive to time but eventually reaches a steady-state mean value. This behavior is a result of the interplay between conduction and convection heat transfer mechanisms during the transient phase of the process. Figures 4A, B depict a considerable enhancement in the value of the heat transfer coefficient with the Prandtl number and the heat capacities ratio. It also shows that the sensitivity of Nu with respect to t increases with increasing Pr and γ. Figure 4C shows that the amount of convective heat transfer rises with Deborah's number. However, the strain-retardation parameter shows a decrease in the heat-transfer coefficient (see Figure 4D). This justifies the fact that the influence of retardation time is to inhibit the heat transfer.

Conclusion
The aforementioned results support the subsequent conclusions. Due to a competition between thermal diffusion, viscoelastic relaxation, anisotropy in permeability, thermal diffusivity, and internal heat generation, the conduction state degenerated into convective motions via the oscillatory mode. Viscoelasticity, the heat capacity ratio, and the Prandtl number have no effect on stationary convection. Increasing anisotropy in permeability, the coefficient of internal heat generation, and stress-relaxation parameters are associated with an early onset. It has been discovered that anisotropy in thermal diffusivity and heat capacity ratios delays convection. Strain-retardation time reinforces oscillatory case stability. The range of retardation parameter values within which oscillatory convection occurs is determined by the magnitude of the relaxation time. Beyond this range, oscillatory convection ceases, and stationary mode instability is then established. The Prandtl number indicates the effect of destabilization on oscillatory convection. The coefficient of heat transfer increases with the Rayleigh number, the internal Rayleigh number, the Deborah number, the Prandtl number, anisotropy in permeability, and the heat capacity ratio. Increasing values of the strain-retardation parameter and anisotropy in thermal diffusivity indicate an inverse trend. It is possible to promote or inhibit convection in a given system by adjusting the relevant parameters in accordance with practical application. In other words, the onset and intensity of convection can be manipulated by adjusting the system's controlling parameters, such as temperature gradients, fluid properties, and geometrical factors. Several of the prior results were obtained through the use of special cases. This provides substantial support for the results of the current investigation (Eswaramoorthi et al., 2023).

Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding authors.

Author contributions
MS, BH, AS: conceptualization, methodology, software, formal analysis, and writing-original draft. AH: writing-original draft, data curation, investigation, visualization, and validation. UK: conceptualization, writing-original draft, writing-review and editing, supervision, and resources. RK: validation, investigation, writing-review and editing, and formal analysis. VK: writing-review and editing, data curation, validation, and resources. All authors contributed to the article and approved the submitted version.

Funding
This work was partially funded by the research center of the Future University in Egypt 2023.