Influence of Temperature Gradient on the Stability of a Viscoelastic Film with Gradually Fading Memory over a Topography with Ridges and Furrows

A long-wave evolution equation is derived using an asymptotic analysis, and the linear stability of a viscoelastic film flowing along the direction of parallel grooves over a uniformly heated topography is studied. A numerical approach adopting spectral collocation technique is used to demarcate the stable and unstable flow regimes. The combined influence of thermocapillarity and viscoelasticity on the films stability is analyzed. By accounting the bottom topography comprising longitudinal gutters, the possibilities of regulating the film dynamics under iso- and nonisothermal conditions and/or optimizing design structure of an apparatus for desirable flow outcomes have been focussed.


Introduction
Falling films are relevant to a broad class of interfacial instability problems over a wide range of length and time scales in various technological setups since they offer small thermal resistance and large contact area at small specific flow rates. This serves as an advantage in many processes involving cooling, condensation, absorption, and evaporation, thereby regulating the transport of mass, momentum, and heat across the liquid-gas and liquid-solid interfaces. For a detailed review of falling film problems, refer to Nepomnyashchy et al. [1], Kalliadasis et al. [2], and the references therein.
Although many studies available in literature have focused attention on Newtonian fluids, many fluids used in industrial applications are rheologically complex and non-Newtonian in nature. This has led to the generalization of Navier-Stokes model to satisfy highly nonlinear constitutive laws to arrive at complex partial differential equations, which are one order higher than the Navier-Stokes equations [3][4][5][6][7][8].
Unlike Newtonian fluids, which respond virtually instantaneously to an imposed deformation rate, viscoelastic fluids respond on a macroscopically large time scale, known as the relaxation time. The viscoelastic fluids, a subclass of microstructure flows, display both elastic (for deformation rates larger than the inverse relaxation time) and viscous (for deformation rates smaller than the magnitude of the inverse relaxation time) characteristics. The stress in this liquid is neither directly proportional to the strain nor to the rate of strain but displays a complex relationship [9]. A lot of work on the flow and heat transfer characteristics of non-Newtonian fluids has also been done in order to control the quality of the end product in many manufacturing and processing industries ( [10] and the references therein). This area of fluid dynamics can simulate accurately many complex polymeric, biotechnological, and tribological fluids and has attracted special attention with the advent and growth of synthetic material industries. It finds applications in clothing, food, and pharmaceutical industry, and in analyzing the lubrication behavior of bearings, gears, and cams.
There are various constitutive models which address the elasticoviscous aspect of a viscoelastic fluid. The Walters-B liquid model [11] can simulate accurately many complex polymeric, biotechnological, and tribological fluids and forms an important basis for the manufacture of plastic engineering equipments, contact lenses, and so forth. The mixture of polymethyl methacrylate and pyridine at 25 ∘ C containing 2 Chinese Journal of Engineering 30.5 g of polymer per liter behaves very close to this model. The physical properties of such a fluid constitute its density = 0.98 × 10 3 Kg/m 3 , limiting viscosity = 0.79 Ns/m 2 , kinematic viscosity ] = 806 × 10 −6 Nsm/Kg, surface tension = 40 × 10 −3 N/m, and viscoelastic coefficient Γ 0 = 0.04 Ns 2 /m 2 . Hydrodynamic stability studies for the mixture of above polymer have been reported by Andersson and Dahl [12], Sadiq and Usha [13], Lin and Chen [14], and Sirwah and Zakaria [15] for thin films.
Chang [16] in his report categorized different types of waves while observing a falling film. He identified interesting regimes where the wave behavior changes depending upon the magnitude of Reynolds number, Re. In the low-moderate Reynolds number regime, long gravity-capillary waves dominate. This region is interesting both from a practical point of view (since many industrial equipments operate under this condition) as well as from a theoretical point of consideration (due to regular structures, it is possible to describe the aspects of the flow by using relatively simple models).
Tackling the complete set of mass, momentum, and energy equations poses a formidable challenge. Therefore, with lubrication approximation at the forefront, dimension reduction techniques that take advantage of thin film geometry have been successfully developed in the majority of theoretical work that has appeared to date. This assumption neglects relatively unimportant terms while scaling complex equations by considering the depth-length ratio to be small and is an essential analytical tool in modeling thin film low-Reynolds number flows (a subregime of 1 < Re < 300). Under this assumption, the physical variables are expanded in terms of the small aspect ratio to yield a more tractable yet highly nonlinear degenerate fourth-order PDE for the local layer, which is still capable of retaining the dominant physics [17]. For nonisothermal problems, some notable pieces of research include Joo et al. [18], Miladinova et al. [19], Scheid [20], and Sadiq et al. [21].
Fluid actuation and transport are core functionalities and challenges in many fluidic systems. The shape and property of the substrate (its curvature and/or flexibility, patterns and/or obstacles) play an important role in the development of flow instabilities and affect the film dynamics and heat transfer. In particular, the interplay between the physical surface and the local fluid flow affects the entire system, where the viscosity, surface tension, thermocapillarity, and the local interface, for instance decides the flow mechanism. A better understanding of these properties and how they can be harnessed and controlled is extremely important. Their impact in relation to optimizing resources and manufacturing cost in the industrial sector have yet to be fully realized, with their use in relation to forming novel functional coatings for applications in many engineering disciplines still at an early stage.
Literature on non-Newtonian fluids, their stability, and dynamics on substrates with topography is not abundant. For a power-law based non-Newtonian model, the effect of substrate topography, inertia, and non-Newtonian rheology was examined for a two-dimensional thin Carreau-Bird model based liquid which has direct relevance in problems connected to the early stages of coating process [40]. The effect of local microtopography with steps, hills, and periodicity on the steady flow was considered by Mogilevskii and Shkadov [41]. The effect of inertia and surface tension on the flow of a power-law film over an undulated surface has been investigated by Heining and Aksel [30]. They identified that the shear-thinning rheology of the film leads to possible resonance. Usha and Uma [42] extended the study carried by Tougou [43] to viscoelastic Walters-B liquid on a wavy incline and found the impact of the viscoelastic parameter on the films stability and long waves generated. Using a convected Maxwell model, Saprykin et al. [44] investigated the influence of inertia and viscoelasticity on a film flowing on a stepdown topography. The effect of viscoelastic properties under creeping flow conditions was investigated by Pavlidis et al. [45] considering a Phan-Thien and Tanner viscoelastic model for fluid flowing on a horizontal trench. By incorporating an asymptotic expansion of the upper-convected Maxwell model for elastic fluids, recently, Bouchut and Boyaval [46] proposed a new model for gravity-driven free-surface flows of shallow viscoelastic fluids across a topography. An important feature of their model is the convexity with respect to the physically relevant pseudoconservative variables.
Shaqfeh et al. [47] assessed the stability of gravity driven viscoelastic flow of an Oldroyd-B model in a broad spectrum of dimensionless groups categorizing the flow into thick dilute, thin dilute and thick concentrated or melt film regimes. In the present problem, the influence of nonisothermal effects and bottom topography is considered on a viscoelastic fluid's stability. For this purpose, a Walters-B liquid is considered to model viscoelastic effects. Though this investigation does not cover the entire class of viscoelastic flows (thin to thick films), yet it gives a qualitative understanding and a first insight into the mechanism occurring in thin film regime, where the surface-tension forces become very important. Also, it forms as the first study in modeling such fluids on a corrugated geometry comprising longitudinal grooves. To cater non-isothermal effects, a uniform temperature (larger than the constant ambient air temperature) is prescribed on the bottom of the topography, and a long-wave equation is derived. The stability of the primary flow to infinitesimal disturbances is discussed under non-isothermal effects.  and the geometry extends infinitely along the -anddirections. The bottom topography, ( , ) = cos(2 / ) undulates weakly along the -axis such that ( , , ) = ( , ) + ℎ( , , ) measures the total thickness. The parameters and represent the amplitude (the peak deviation of from its center position) and the wavelength of a typical periodic groove. The air above the liquid is passive such that its temperature, ∞ , and pressure, ∞ , are constants sufficiently far from the interface. A uniform temperature distribution is prescribed at the bottom of the topography assuming the substrate to be a perfect conductor of heat of infinite heat capacity. Boussinesq approximation is used to describe the bulk movement of the liquid. In the framework of this approximation, density ( ), kinematic viscosity (]), thermal conductivity ( ) and thermal diffusivity ( ) are assumed constants [1,20,48]. The surface tension decreases linearly as a function of temperature and is given by = ∞ − ( − ∞ ), where ∞ is the mean surface tension at gas temperature and = − / is a positive constant as for most common liquids.

Problem Description
Constitutive equations which result from a retardedmotion expansion represent the correct constitutive equations for flows in which the rate-of-strain tensor and its time derivatives are small [49]. The second order fluid is derived through a retarded motion expansion and shows two different normal stress functions [49] in simple shear flows, being valid for slow and weak flows [50]. Though the shear thinning seen in experiments on polymeric fluids is not described by this model, yet it finds usage in many applications, and it also predicts most non-Newtonian effects qualitatively, if not always quantitatively [51].
Dandapat and Samanta [52] assessed the constitutive equations of the second order incompressible fluid [49] framed on the postulate of gradually fading memory [53]. Such fluids show normal stress differences in the shear flow (which is a characteristic property of a viscoelastic fluid), and the constitutive equation is applicable to the flow of dilute polymer solutions (e.g., polyethylene oxide in aqueous solution (POLYOX), polyisobutylene in cetane), which are only slightly viscoelastic and may have a bearing on the study of fluids with drag-reducing properties (Toms phenomenon).
Regarding the applicability of their model, if the disturbance time scale is large compared with the characteristic time scale (relaxation time) of the fluid, then the second order fluid model is internally consistent with the stress-relaxing fluid due to Oldroyd [3]. As pointed out by Porteous and Denn [54], this would happen if the viscoelastic parameter is small.
The components of the symmetric Cauchy stress tensor, T, in a Walters-B fluid are given by [55][56][57] where is the stress tensor, represents the isotropic pressure, ( ) is the total viscosity of the Maxwell elements whose relaxation time lies between and , ( ) represents the coordinates of a fluid element at time , ( ) designates the coordinates of the same fluid element at (> ), is the metric tensor, and (1) indicates the rate of strain tensor. Neglecting ∫ ∞ 0 ( ) , ≥ 2 for liquids with short relaxation times, (1) reduces to a model involving only one non-Newtonian parameter [12,57]: where Γ 0 = ∫ ∞ 0 ( ) is the viscoelastic parameter exhibiting low relaxation times, = ∫ ∞ 0 ( ) is the limiting dynamic viscosity at small rates of shear, / evaluates the upper convected time derivative or Oldroyd derivative of a tensor [58], and are the components of the velocity u = ( , V, ) in the direction ( , , ) of a fixed coordinate system. Refer to Appendix A for the complete expressions of the components of the stress tensor.
Since the model proposed by Dandapat and Samanta [52] is for weakly viscoelastic fluids, it yields the same evolution equation for the free surface which results by adopting (2a)-(2b) [59]. This is true since 0 [52], the coefficient of crossviscosity [60] defined by Dandapat and Samanta [52], is ignored by assuming the fluid to be weakly elastic, which occurs exactly by neglecting ∫ ∞ 0 ( ) , ≥ 2. Such an idealized model represents an approximation to first order in elasticity [12] and is a valid approximation of Walters-B liquid taking very short memory [61]. Barnes et al. [62] have introduced a general canonical form of constitutive equations for viscoelastic fluids. The rheological model given by (2a) and (2b) has been recovered as a special case of this general form and is classified as a second-order model according to their nomenclature. Also, the justifications proposed by Dandapat and Samanta [52], Gupta [63] and Dandapat and Gupta [64] qualify and validate (2a)-(2b).
The equations governing the flow are those for an incompressible fluid such that the compressibility effects and viscous heat generation are neglected in the heat equation: where F = ( sin , 0, − cos ) is the body force, is the standard gravity constant, and is the temperature. To complete the mathematical description of = ℎ + u ⋅ ∇ on = ( , , ) , (kinematic boundary condition at the free surface) .
In (5)-(7) above, n = ∇( − )/|∇( − )| is the unit outward normal at any point on the free surface, and t and t are the unit tangential vectors such that n⋅t = 0, = , . The heat exchange coefficient between the liquid and the passive air is described by .
In order to nondimensionalize systems (3)-(8), the problem is related to the flow of a falling film on an infinitely extending rigid plane [28,38,39]. For a stationary flow of uniform thickness ℎ 0 , the solution of the Navier-Stokes and Fourier (expressing energy balance) equations is a parallel flow dependent only upon the cross-stream coordinate such that viscosity balances gravity and the heat is propagated by pure conduction. The primary flow variables [65] are (neglecting the capillary boundary layers at the side walls The depth averaged velocity of the Nusselt flow is defined

and it serves
to scale the streamwise component of the velocity, for which the constant flux is given by 0 = 0 ℎ 0 . The Nusselt thickness is considered as the characteristic length scale in thedirection. The characteristic wavelength on the free surface is considered to be of the same order as the wavelength of the periodic undulation and serves as a scale in the -and -directions. Mass balance condition serves to choose the scales for the velocity components along the spanwise and transverse directions. The pressure is scaled using 2 0 and the viscoelastic parameter [13] by ℎ 2 0 . Lastly, the temperature scale is chosen to allow the largest possible temperature difference [21,33]. In short, the scales can be summarized as After nondimensionalizing, the asterisks are dropped from the scaled variables (cf. Appendix B for the complete set of nondimensional equations). The set of dimensionless parameters arising during the nondimensionalization procedure is Re = 0 ℎ 0 / (Reynolds number), Γ = Γ 0 / ℎ 2 0 (viscoelastic parameter), We = ∞ /( 2 0 ℎ 0 ) (Weber number), Ca = ( − ∞ )/ ∞ (measures the ratio of temperature dependent part of the surface tension over its basic value; note that it is defined sometimes as the Capillary number [18,19]), Pr = ]/ (Prandtl number), Pe = Re Pr (Péclet number), = ( − ∞ )/ 0 (Marangoni number), and Bi = ℎ 0 / (Biot number). Referring to Joo et al. [18] and the references therein, the orders of magnitude assigned to the dimensionless parameters are Re ≃ O(1), 2 where is the small aspect ratio. Note that the assumption ≃ O( ) allows to construct a series solution analytically as mentioned later [26,30]. However, when ≃ O(1), the base solution is found through a numerical procedure [26,30].
The physical variables , V, , , and are asymptotically expanded in powers of the small parameter and solved sequentially at each order [21,67] using the symbolic computation toolbox in MATLAB to obtain the long-wave model (cf. Appendix B for the equations and solutions at O(1) and O( )): It is to be noted that the above long-wave model is obtained without estimating the time derivatives from the leading order solutions as traditionally done [35,36]. For small Reynolds numbers, the long-wave model described by (11) captures the dynamics efficiently since the temperature and velocity fields are slaved to the kinematics of the free surface. The second term in (11) describes wave propagation and steepening effect. The remainder of the terms in (11) within the curly braces decides the growth or decay of waves of which the first two terms are due to the inertial contribution and correspond to the destabilizing effects of the mean shear flow, the third and the fourth terms measure jointly the viscoelastic effects, and the fifth term is responsible for the effects of mean surface-tension. While the sixth term accounts for the stabilizing effects of the hydrostatic pressure, the seventh term is due to thermocapillarity and represents the perturbation of the interface temperature induced by variations of ℎ, when there is heat transfer to the gas phase.

Temporal Linear Stability and Results
It is assumed that the primary flow, ℎ = ( ), is independent of and and depends only upon the spanwise coordinate . The primary flow satisfies the following ordinary differential equation: For ≪ 1 (a weak undulating topography), a closed form expression for the thickness profile of the primary flow can be computed [26,30,[37][38][39] by expanding ( ) = 1 + 1 ( )+O( 2 ) such that 1 ( ) = 1 sin(2 )+ 2 cos(2 ) is periodic and represents the first order correction to the primary flow. This assumption used in (12) yields 1 = 0 and . Since is small, another alternative to construct a mathematical form for the primary flow would be to set ( ) = 1 + 1 ( ) + O( 2 ). This ansatz, however, imposes a mathematical condition on the Reynolds number, angle of inclination, wall amplitude and the Weber number associated with the flow at the leading order. And at O( ), it involves solving a quartic auxiliary equation associated with the differential equation for determining 1 ( ), which as a result yields complicated complementary functions depending upon the nature of roots. Therefore, this approach is avoided here. Having known the undisturbed thickness of the basic state, the timeindependent properties of the fluid like velocity, temperature, and so forth can be evaluated at the free surface and averaged over a typical groove [33] by referring to the leading and first order solutions in Appendix B. This yields The above equations give a first insight on the flow mechanism associated with the topography comprising parallel grooves. From (13a), it is inferred that the averaged timeindependent stream velocity enhances (which would cause the averaged flow rate,̃, to increase) as the groove amplitude increases, suggesting possible flow destabilization. Also, the temperature at the bottom topography is larger than the averaged temperature at the free surface (cf. (B.2) and (13b), and note that the temperature averaged on a groove at the bottom substrate has the value 1). It is known that when the interfacial flow is driven by local changes in interfacial tension due to heating of the substrate, the modification of the temperature at the free surface due to surface deformation triggers long-wave thermocapillary instability [68,69], sometimes leading to unwetted regions [70]. Therefore, the difference in temperature affects the surface tension and causes the fluid to flow from hotter region to colder regions, thereby increasing the wave amplitude. It should be remarked that contrary to the observations in Gambaryan-Roisman and Stephan [35], Sadiq et al. [37], and Sadiq [38], (13a) proposes possible flow destabilization on a patterned topography as the averaged flow rate over a typical groove increases. However, at this juncture, it is really premature to speculate and conclude on this behavior of the time-independent velocity on a grooved topography, given the fact that (13a) and (13b) are obtained only for small . Adopting a numerical solution of ( ) (as discussed later), a detailed investigation on the equilibria would only reveal the true picture of the basic state getting either stabilized or destabilized. In this regard, a small harmonic disturbance (a function of times a sinusoid, which is a function of and ) is imposed on the primary flow thickness to investigate the temporal stability and to see whether the disturbance dies away in the course of time (stable flow), persists as a disturbance of similar magnitude (neutral stability) or continues to grow in time leading to a different laminar or turbulent flow (instability) in the form, where = − is the temporal pulsation (eigenvalue) such that , , and are real numbers and represent the wavenumber, the linear amplification, and the linear wave velocity, respectively (by fixing to be real, the solution can oscillate but is not allowed to grow/decay in ). The absolute value of the disturbance amplitude function ( ) is small such that max(| ( )|) = O(1) and ≪ 1. Substituting (14) in (11), a linear ordinary differential equation of fourth order in ( ) is found: For the flow on a planar wall, the terms contributing to the real and imaginary parts of the complex eigenvalue are found as ( = 1, = 0, and is a constant) For ≪ 1, the above result reduces to the asymptotic analysis result of the Orr-Sommerfeld system considered in Sadiq et al. [21] for a Newtonian fluid, however on a rigid impermeable plane. The condition for the surface wave instability ( > 0) from (16a) becomes which shows that the critical Reynolds number decreases when the viscoelastic effect increases and also when the thermocapillary effect increases, indicating flow destabilization. However, when increases the critical Reynolds number increases implying stabilization of the flow. Also, note that ignoring surface tension and heating effects at the leading order for a Newtonian film, the well-known hydrodynamic instability condition Re > (5/6)cot is recovered, which shows that a vertically falling film is always unstable.
When is a function of and the flow is on a planar wall, (15) becomes a fourth-order linear differential equation with constant coefficients, which can be straightforwardly solved. On a topographical substrate when Marangoni and surface-tension effects are neglected, (15) reduces to the form ( )( 2 / 2 ) + ( ) = 0. For an eigenvalue , the eigenfunction could be explicitly solved using a series solution technique for the above second-order differential equation at ordinary points ̸ = (1/2 )cos −1 (1/ ). Since the coefficients of (15) are periodic with period 1, according to Floquet theory, can be expressed in the form [71] ( ) = ( ) , where ( ) has period 1. The parameter is taken real by neglecting the exponential growth of disturbances along the spanwise direction. Together with the periodic boundary conditions for ( ), when the topography is uneven, the eigenvalue problem (15) can be solved numerically to determine ( , ) and ( , ). Also, note that (15) remains invariant under the transformation → − , suggesting that (15) is symmetric about the ( , ) plane (remark that (11) is actually symmetric about this plane), and hence, one can concentrate only in the region ∈ [0, 1/2] and demand / = 0, 3 / 3 = 0 at the end points. Another alternative would be the application of symmetry boundary conditions along = 0 and periodic boundary conditions at the end points of a typical groove. However, in the rest of the analysis, only periodic conditions are assumed.
In order to solve the linear stability problem (15) The derivatives / are replaced by a matrix of dimension (0 : , 0 : ), whose elements are where , are distributed according to the cosine law (18). Further, 0 = = 2, = 1, and 1 ≤ ≤ − 1. The vector obtained by multiplying with the matrix is the discrete approximation of / . Higher order derivatives are approximated and obtained by repeated multiplication by the matrix . For the cosine law based on the odd values of and for the theoretical foundation for approximating the derivatives through successive matrix multiplications, the readers are advised to refer to the standard literature on spectral methods [72,73].
Instead of choosing a particular fluid and analyzing its stability (which has only a limited scope), the parameters are changed by varying the components of P = ( , , Re, , Γ, , and Bi) to analyze the equilibria. However, following Joo et al. [18] and Miladinova et al. [19], = 0.2 is fixed. Also, while varying the parameters such as Re and Bi, it is ensured that these parameters do not take unrealistic values. For instance, the long-wave equation on a planar surface displays finite-time blow-up behavior when the Reynolds number exceeds a particular value [20,74], and for liquid films in contact with gases in many realistic experimental situations, Biot number is usually small [20]. Further, as the rheological model used in the study is assumed to be realistic only for weakly elastic fluids, the resulting flow should be considered as a perturbation of a Newtonian viscous flow, and therefore, small values of Γ are chosen. Lastly, appreciable surfacetension effect is included bearing in mind the thinness of the film.
In the numerical calculations performed, the maximal values of correspond to = 0, so that the most unstable modes are those for which ( ) is periodic. Therefore, the dependence of on alone is considered from this point onwards. Note that the numerically generated solution of (12) was considered in the study for plotting all the figures presented here. The reason for doing this is twofold: (i) the accuracy of ( ) theoretically constructed using the imposition ≪ 1 can be checked and (ii) solutions at ∼ O(1) can be computed. This was accounted by integrating (12) numerically by approximating the derivatives using the spectral collocation technique and then solving the resulting nonlinear algebraic system using Newton's method. However, an initial condition is required to initiate integration. For this purpose, the evaluated analytical expression of the primary flow was considered as the initial condition, and periodic conditions were chosen as boundary conditions. Figure 2 displays the free surface obtained through approximated analytical and numerical solutions of ( ). It is observed that the numerical solution of ( ) shows a largest discrepancy of 8% from the analytical solution. This reduces further either when the magnitude of the Weber number increases or when the topography is more inclined, leading to agreement between the numerical and analytical values of ( ), as shown in Figure 2. Interestingly, D' Alessio et al. [33] also reported on the behavior of free surface approaching unity as the magnitude of surface tension increases (cf. Figures 4 and 23 in their research paper). Finally, note that in order to ensure the accuracy of the numerical results yielded by the spectral method, the derivatives were approximated using a central difference scheme of second-order accuracy and solved (only the numerical results obtained through spectral method are however presented in this paper). The numerical results ensuing due to these two different procedures agreed with each other. For this purpose, ≥ 40 was considered in the spectral technique, and ≥ 200 was chosen to implement the finite difference scheme.
Stability characteristics of Newtonian and non-Newtonian isothermal layers are examined first. It is known that the flow of an isothermal liquid film with a smooth free surface down a vertical plane is unstable at all Reynolds numbers [75]. Also, since the amplitude of the long waves on a vertically falling film is larger than that along an inclined substrate case [37,38], in the rest of the investigation, = /2 is set. It is ascertained from the numerical behavior of the linear amplification curves that the inertial force increases the instability threshold for Newtonian films falling along parallel longitudinal grooves ( Figure 3). For weakly undulating topography at small Reynolds numbers, the grooves aligned parallel to the flow direction decrease the growth rate leading to waves with small amplitude, and reduces the range of unstable wavenumbers in relation to the flow on a flat substrate. As the Reynolds number increases, there exists a certain wavenumber below which the growth rate for the flow along grooves is larger than that along a planar substrate, and such a trend changes beyond this wavenumber (Figure 3(b) (i) and Figure 3(c) (i)). Take a note on the maximum amplitude profiles presented in Sadiq et al. [37] and Sadiq [38], where max | >0 > max | =0 for a short period of time before this behavior changes. The probable explanation to this behavior can be gathered from (13a), which suggests increase in average velocity at the free surface on a typical groove. Therefore, when a small disturbance is imposed on the primary flow, the grooves tend to increase the growth rate initially. When the bottom undulations increase, the numerical investigation exhibits interesting behavior on the films stability. There are instances at which the maximal growth rate corresponding to a large Weber number is larger than that related to a small Weber number; however, the range of unstable wavenumbers remains smaller than the case corresponding to small Weber numbers (Figure 3(b) (ii)). For a fixed moderate groove amplitude, when the Weber number increases, the flow destabilizes at small Reynolds numbers (Figure 3(a) (ii) and Figure 3(b) (iii)) and stabilizes when the magnitude of the inertial force increases (Figure 3(c) (iii)). Pointing out on the observations reported through direct [27,28,33] and indirect problems [76] for the flow across wavy undulations on an inclined topography (where the critical Reynolds number decreased on an undulating topography when the effect of surface tension is appreciable), here, at fixed Re and , the variation in We shows either destabilizing or stabilizing tendency on the flow.
When the non-Newtonian effects are included, the infinitesimal disturbance to the primary flow leads the viscoelastic film to destabilize in relation to the Newtonian film. The sketch in Figure 4 displays interesting trends associated with the growth rate curves of isothermal Newtonian and viscoelastic layers on a vertical substrate, when there is a pronounced inertial effect. At Re = 5, the numerical results of the stability analysis predict decrement in the value of growth rate when the surface-tension effect increases (compare the curves in each column in Figure 4). Analogous to the Newtonian flow, the growth rate curves corresponding to the viscoelastic film on a grooved topography show increasing and decreasing behavior before and after a certain wavenumber in connection to the flow on a flat substrate. However, at large Weber numbers, when the viscoelastic parameter increases, the growth rate and the window of unstable wavenumbers for the flow on a topographical substrate are larger than the case when = 0.  the undulating groove becomes larger, the growth rate curves and the range of unstable wavenumbers start decreasing (cf. curve 4 in Figure 4(c) (ii) and Figure 4(c) (iii)). See also Figure 11 presented later displaying the above characteristics on a heated uneven substrate. Figure 5 shows the dispersion curves for linear growth rate of a viscoelastic film at small Reynolds numbers. Just as in the case of a Newtonian film, at fixed moderate groove amplitudes, the topographical substrate destabilizes the small Reynolds number viscoelastic flow as the magnitude of the surface-tension effect increases ( Figure 5 (ii), (iii)). As a note from the numerical study carried, it should be mentioned that when a value as high as Γ = 1 was considered (such a value may be beyond the realistic range of weakly elastic fluids), the longitudinal grooves destabilized the flow completely.
To substantiate the information availed from Figures 3, 4, and 5, neutral stability curves are presented. While the neutral stability curves in Figure 6 clearly display the destabilizing nature of the periodic bottom undulations at small Reynolds and large Weber numbers, Figure 7 distinctly shows the instability at small when the viscoelastic effects increase at Re = 5, as already discussed in relation to Figure 4.
In Figure 8, the phase velocity curves varying with are displayed. The phase speed value is 3 for the flow on a planar substrate initially (this agrees with the theoretical prediction suggested by (16b) at ≪ 1). However, the nonlinear nature of the phase velocity curves can be attributed due to the gross effect of all the dimensionless parameters involved in the problem without truncating or approximating them at any level; for example, see (16b). The linear wave velocity curves first decrease up to a certain wavenumber and then increase beyond it. It is seen that the phase velocity decreases as the viscoelastic parameter increases. This fact can be confirmed by referring to the flow over a flat substrate (cf. (16b)), where it is observed that when > 0, would increase as a function of Γ only when Γ > 5/12. Otherwise, despite > 0, can decrease when Γ increases provided that 0 ≤ Γ ≤ 5/12. Only numerical calculations reveal the influence of an uneven geometry on the wave velocity. Computations show that the phase velocity is larger for the case of a flow on a grooved topography than on a planar substrate. For the case of a topographical substrate subjected to uniform heating, D' Alessio et al. [33] reported on very interesting phenomena associated with Marangoni effect. Their study witnessed that thermocapillarity can either stabilize or destabilize the flow depending upon the strength of the Marangoni number and can also lead to a reversal in flow stability. Figure 9 shows the response of thermocapillary effect on Newtonian and viscoelastic films over a flat and uneven topography. At low Reynolds number, the amplitude of the grooves decides the instability mechanism depending upon the strength of heating. For small groove amplitudes and for a fixed Weber number, increase in Marangoni number leads to film destabilization. Apart from the already observed behavior of flow destabilization at moderate reasonable with increase in Weber number, it is seen that increase in Marangoni effect can promote stability at large Weber numbers (cf. Figure 9(c)). Such an effect is also observed when the Reynolds number increases (figures are not included here in order to avoid repetition). Observe that in the recent research report of D' Alessio et al. [33], they have identified on a semiinclined plane (cf. Figure 8 in their paper) that for a given groove measure with increase in Marangoni number, the critical Reynolds number increases leading to film stabilization, provided the surface-tension effect is large. When the viscoelastic effects are included and when the Reynolds number associated with the flow is very small, the film tends to destabilize in the absence of Marangoni effect. However, at = 0.7 and for pronounced surface-tension effects, increase in heating corresponds to flow stabilization (in Figure 9(d), compare curves 3 with 4 and 5 with 6). The influence of Marangoni stress on the stability threshold at low Reynolds number is presented by means of Figure 10, which shows improved stability at large Weber numbers at certain regions of . However, it is also inferred from the figure that as increases, the trend changes, and once again Marangoni stresses start destabilizing the flow. In conjunction with the results displayed in Figure 4(c) (ii), Figure 11 reveals film destabilization at Re = 5 when the viscoelastic effects increase on a grooved topography (cf. curves for = 0 and 0.4 in Figure 11). The effect of thermocapillarity in this case is to increase the growth rate. However, at moderate = 0.9 increase in heating stabilizes the flow, similar to the case discussed in Figure 9(d). Figure 12 presents linear amplification curves at Re = 5 for different values of Biot number for both Newtonian and non-Newtonian liquids. When the rate of heat transfer from the liquid to the ambient gas phase increases, the flow tends to destabilize (Figures 12(a) and 12(b)). For a non-Newtonian film on a planar substrate with increase in the magnitude of Bi, the growth rate curves dominate the Newtonian case causing more flow destabilization (Figure 12(a)). On a substrate with topography, the groove amplitude plays a crucial role in deciding films stability. Note that increase in the value of viscoelastic parameter for a flow on a topography for which = 0.4 causes flow destabilization compared to the planar case at small values of Bi as inferred from the curves 4 in Figures 12(a) and 12(b) for which the range of unstable wavenumbers, respectively is [0, 0.984] and [0, 0.972], despite the respective maximal growth rates being max = 0.27 and max = 0.303. However, when the heat transfer rate increases, the growth rate curves over a topography with small amplitude increase up to a certain wavenumber in relation to the growth rate curve over a flat substrate, and this trend changes beyond such a wavenumber causing reversal in stability (cf. curves 6(a) and 6(b) in Figure 12(d), where these curves correspond to curve 6 from Figures 12(a) and  12(b), resp.). At moderate values of groove amplitude, it was already observed that the Marangoni effect can stabilize the flow mechanism. This stabilizing mechanism enhances further when the heat transfer rate increases (Figure 12(c)). When the Reynolds number corresponding to the flow is small and when and are large, increase in Biot number stabilizes the flow configuration as seen in Figure 13 (notice that the stabilization process corresponding to increased thermocapillary effect is also seen).
The present investigation, although valid only at the region where the viscous forces dominate inertia due to longwave modeling, has revealed interesting results associated with the linear stability. The results showed the dramatic impact of the bottom topography on the films stability in conjunction with the force of surface-tension and Marangoni stresses. The analytical expansion of the streamwise velocity at the free surface averaged over a typical bottom undulation does reveal only its increasing tendency contributing towards flow destabilization, but the numerical investigation on the films stability exposes the additional effects involved in the problem. A close-up theoretical investigation can be performed by referring to Figure 2 to get an insight into the physical mechanism involved. In this figure, it is observed that the free surface flattens up at = 1 for a large value of and a small value of bottom amplitude. Using this fact and referring to Appendix B, the leading order solutions can lead to An integration of (21) yields the base flow rate as and the averaged flow rate as indicating that for large values of surface tension force the averaged value of the flow rate increases leading to flow destabilization. The undisturbed thickness being a function of the spanwise variable, when perturbed along the -direction, leads to an eigenvalue problem comprising terms arising as a consequence of streamwise perturbation. The net effect of these terms with the aid of numerical approach shows that both flow stabilization and destabilization are possible in various parametric regimes.

Concluding Remarks and Objective
The linear stability of a thin Walters-B fluid due to the linear variation of the surface-tension along a liquid interface induced by temperature gradients along the liquid layer was examined on a topography comprising sinusoidal longitudinal grooves. The nontrivial eigenvalue problem was numerically solved using a Chebyshev spectral collocation technique and the results were compared with those occurring from a central difference scheme of second order accuracy to ensure numerical correctness.
Numerical investigation of the eigenvalue problem corresponding to a falling film on a vertically held substrate revealed interesting results. In principle, when the viscous effects dominate the inertial forces, the flow tends to stabilize the system on a longitudinal groove topography, provided the surface-tension effects are weak. However, when the magnitude of the surface-tension force increases, interesting trends were exhibited by the stability curves. Some of the interesting key observations in the study can be enlisted as follows: (i) increased surface-tension effects contribute towards flow destabilization at low-Reynolds number as the amplitude of the bottom undulations increase (ii) Marangoni stresses can lead the system towards flow stabilization for certain increasing range over , provided the strength of the Weber number is appreciable, (iii) viscoelastic effects destabilize the system compared to a Newtonian flow, but the structured topography tends to stabilize the flow in relation to the flow on a planar substrate for weak viscoelastic effects, (iv) increase in viscoelasticity promotes instability for small groove amplitudes when the Weber and Reynolds number magnitude increases and (v) increase in heat transport rate at large Weber numbers causes flow stabilization at low-Reynolds numbers.
Calculations were performed through numerical approach to analyze the effect of perturbation on the undisturbed state on a longitudinal grooved topography. Theoretically, it was also found that the average flow rate increases when the surface-tension force increases at small groove amplitudes. Although the above facts are findings based on a numerical linear stability study, such a behavior is yet to be confirmed through laboratory investigations. However, the latest research article by D' Alessio et al. [33] reported strikingly on similar characteristics for the flow of a Newtonian film but across the undulations over an inclined topography.
It is not only with a fundamental research point of interest but also due to benefits gained through peer research that the thin film community intensively and actively strives to shed light on understanding intricate and complex mechanisms involved in falling film problems and fill up unknown gaps. It should be remarked that Mazouchi and Homsy [77] compared the solution to two dimensional Stokes flow to their earlier thin-film based solution. This comparison indicated that the thin-film hypothesis remains valid even in the presence of steep topographic variation. In this regard, theoretical understanding and modeling of flows and deformation in Newtonian/non-Newtonian fluids on small scales on a patterned topography assuming the film to be thin would help in the design and fabrication of small scale structures over the topography to regulate heat and mass transfer in soft matter with smart and advanced functionalities to suppress instabilities and in preventing rupture mechanism.
In addition, such studies may help in understanding the unwanted influence of surface irregularities, which result from a particular stage in a manufacturing process thereby affecting the flow mechanism.
Walters-B fluid finds its usage in medicinal diagnosis and clinical applications. For instance, Nadeem et al. [78] have mentioned that viscoelastic fluid can easily flow in intestine because elastic materials strain instantaneously when stretched and just as quickly return to their original state once the stress is removed and that the Walters-B model is quite close to our real life system. The present study may serve as a part of the basis to study the peristaltic flow of Walters-B model in asymmetrical channels or tubes with geometries inside. Also, Walters-B type finds relevance and importance [60] in geophysical fluid dynamics, chemical technology, and petroleum industries. An implication of the study could be possibly in the product quality and patterning of the substrates, which require multilayer coating where one of the purposes of the bottom layer is to planarize the topography, that is, to smooth out the uneven topography. And a nonuniform coating leads to low quality products or to manufacturing failures. Pavlidis et al. [45] reported on the steady viscoelastic film flow over 2D topography with the prime objective to analyze the effect of the viscoelastic properties of polymeric solutions in spreading over and planarizing various topographical features, especially during spin-coating. Such studies are interesting and could be extended for topographies with longitudinal geometries for second-order fluids with or without invoking lubrication approximation.
Mathematical models obtained by combining a systematic gradient expansion with weighted residual techniques using polynomials as test functions may possibly shed more light on the stability thresholds beyond the low-Reynolds number regime [20,79]. Such a study has been accomplished by Amatousse et al. [80] to study traveling waves on falling Walters-B viscoelastic fluid. As a future objective, a Galerkin strategy will be employed to obtain a system of equations coupling the film thickness and flow rates at various orders to investigate the stability [20,38,79,80]. Also, it would be interesting to investigate the combined effect of nonuniform heating and the substrate topography on the equilibria of a film, for the case when the flow rate is high and when the topography is inundated. Of particular interest would also be to analyze the phenomena in case of a partially filled topography occurring due to low flow rate, where the peaks of the topography are exposed and remain dry. Another direction of study could be to use an Oldroyd-B model and investigate the stability if the nonlinearities inherent in more realistic constitutive equations become important [47]. Such studies will be addressed in future investigations.

A. Stress Tensor Components
The expressions for the components of the stress tensor T in the explicit form are Re