Numerical study of hydrothermal and mass aspects in MHD driven Sisko-nanofluid flow including optimization analysis using response surface method

A steady, incompressible, two-dimensional Sisko-nanofluid flow towards the horizontal direction with no movement in the vertical direction is considered on a stretching/shrinking surface. The power law component (Sisko model) is incorporated under the regime of the porous medium. A magnetic impact is included coming from the MHD in the surface normal direction. In addition, thermal radiation, Brownian diffusion, and thermophoresis are involved in the governing system of equations obtained from the Navier–Stokes model in two-dimensional flow systems. The PDEs are converted into the one-dimensional system using suitable transformations and solved by Galerkin weighted residual method validated with the spectral collocation method. The optimization analysis is performed on heat transfer and skin-friction factors using response surface methodology. The impact of the parameters involved in the model has been testified and is provided in graphical forms. The outcomes indicate that for the values of the porosity factor fluctuating between [0, 2.5], the velocity profile and corresponding boundary layer thickness are lesser towards the maximum value of the parameter, and the results are opposite as the parameter approaches zero. The optimization and sensitivity analysis shows that the transport of heat sensitivity towards thermal radiation, Brownian diffusion, and thermophoresis declined whenever the Nt and Nb increased from low to high and at the medium level of thermal radiation. An increment in the Forchheimer parameter increases the sensitivity of the rate of friction factor, whereas increasing the Sisk-fluid parameter has the reverse effect. Elongation processes like those of pseudopods and bubbles make use of such models. The idea is also widely used in other sectors, such as the textile industry, glass fiber production, cooling baths, paper manufacture, and many more.

www.nature.com/scientificreports/ Governing equations and physical model. The present study incorporates Buongiorno's model utilizing the Sisko-fluid properties adopting a laminar, steady, two-dimensional flow in the positive plane of xy-coordinate system subject to zero movement in vertical direction. The only fluid motion is taken in the x-direction for both shrinking and stretching cases of the surface depending on the stretching/shrinking rate of the surface taken as u = u w . The power law component is represented by j, j = 1 . For j = 1 the model reverts to the simple Newtonian models. A magnetic impact is included in terms of B 0 coming from the MHD in the surface normal direction, as shown in the geometry. The temperature and concentration are taken as constant values (T w , C w ) , respectively, at the solid surface, whereas (T_∞, C_∞) represent their ambient states. The system represents proper boundary layer formulations for both ways (stretching or shrinking). The geometry of the stretching surface has been described in Fig. 1. Mentioning all the boundary conditions. The detailed study of available published literature regarding nanofluid flow shows that porous structures are highly applicable in many industrial applications (see, for example 23,27,29 ). Hence, in the present model, we have emphasized the development of fluid flow using a porous medium. Referring to Khan et al. 40 and Prasannakumara et al. 41 , we can write the complete governing model as follows: The relevant governing boundary conditions (see for example, Mahmood et al. 42 ), are listed below: The following transformations, highly suitable to power law models, are introduced in the governing PDEs to resolve them into ODEs for subsequent solutions via numerical scheme (see, for example 23,39 ): u = εu w , v = 0, T, C = T w , C w aty = 0 , www.nature.com/scientificreports/ Accordingly, the following ODEs are obtained together with the converted boundary conditions given in Eqs. (8)- (12), respectively. A significant number of parameters are involved in the systems (1)- (12); therefore, a complete nomenclature of the parameters is necessarily provided below in the Table 1.
Significantly, the dimensionless expressions of the total viscous frictional factor C fr (i.e., the skin friction coefficient) and the wall thermal transfer rate Nu x (i.e., Nusselt number) are given locally by: Moreover, the involved shear stress components τ wr , τ wϕ , along with the heat flux, q T have the following expressions: The finaly simplified version of system (15)-(16) is given in (17)-(18), as follows: www.nature.com/scientificreports/ Solution methodology and comparative analysis. GWRM, also known as Galerkin weighted residual method, is used to solve the final ODEs. The method can be called a modification case of the Least Squares method. It prioritizes the approximating function derivative over the residual function derivative for the relevant unknown a i . The weight function can be reproduced as follows: GWRM is an effective technique for finding BVP solutions for the following reasons: 1. In the governing differential equations, the initial guesses are assumed to be combinations of trial functions with unknown coefficients. 2. These solutions are incorporated in the equations which contain residuals. 3. The errors are restricted to a low ratio. In addition, a. The easiness of handling BVPs, b. High precision, fast convergence, and efficiency, c. The minimized range within 0 to ∞ , are excellent feature of this method.
The approximate solutions are sought out as follows: where χ(x) is an unknown dependent variable corresponding to the function f (x) within the domain D 0 with the operator L. The approximation of the solution comes from the sum of the initial value with the summation of approximated results for any index 1 ≤ k ≤ n . In general, one can write, Table 2 presents data that explains the arguments and corresponding coefficient values at the respective argument. At approximately 30th argument, the coefficient becomes negligibly small.
A comparative analysis was performed for the spectral collocation method (SCM) and Galerkin weighted residual method (GWRM). Figure 2a,b represent the comparative analysis of velocity profiles for shrinking and stretching surfaces, respectively. An elegant agreement has been found in both cases. Figure 2c,d represent the comparative analysis of SCM and GWRM for temperature and concentration profiles, respectively.

Outcomes and discussions
The impact of various parameters named as magnetic field, thermal radiation, Forchheimer number, porosity number, Sisko-fluid parameter, and power law index represented here with notation j, Brownian diffusion, thermophoresis & Prandtl number on all the three main profiles of fluid flow are analyzed by finding numerical data which is subsequently plotted in graphical forms. A dual analysis is performed within the same frame of reference for each parameter for shrinking and stretching surfaces. Solid lines in the graph represent shrinking surfaces, www.nature.com/scientificreports/ while dashed lines represent the stretching surface. Specifically, the diagram in Fig. 3 depicts the importance of magnetic field parameters and the consequences of incremental trends in its values on the flow profile. In both cases, the profile undergoes the same declining trend for shrinking and stretching. The velocity profile and the corresponding boundary layer thickness are much smaller for larger values of the parameters. In contrast, when the values of the parameters are closer to zero or approaching zero, the profile shows more thickness and a higher temperature profile. The rise of a Lorentz force due to magnetic parameters is solely responsible for this case. The surface normal direction of the field lines directly interrupts the fluid movement, and therefore, the profile shows a decline when the magnetic field impact is enhanced. The influence of the involvement of porous medium in fluid flow analysis generally appears in terms of the Forchheimer and porosity numbers. Figure 4 reflects the traits of the velocity profile governed by the elevated values of the Forchheimer number, which consequently results in an inevitable decline in fluid flow. Physically, the inertial influence occurring due to the Fr results in the deterioration of the fluid velocity. Similar trends appear in the case of the porosity number given  www.nature.com/scientificreports/ in Fig. 5. For the values fluctuating between 0 ≤ ≤ 2.5 , the given profile and corresponding boundary layer thickness are lesser towards = 2.5, and the results are opposite as = 0 is approached. This trend appears due to the enhancement of frictional force by adding the porosity factor. At = 0 , the fluid flow medium returns to a simple medium having no porosity and no additional frictional force. The impact of the material parameter, also known as the sisko-fluid parameter, on the velocity distribution is given in Fig. 6, which reflects the rise in velocity profile for larger parameter values. Looking back, the constituency term of the material parameter is shown in Table 1, the material parameter (Sisko-fluid parameter) comprises the consistency index and the shear rate viscosity. The inverse relations depict that a rise in material parameter diminishes the viscosity, and therefore, the fluid movement catches a jump.   www.nature.com/scientificreports/ Figure 7 reflects the trends in temperature distribution uprising due to the thermal radiation factor involved in the model subject to Rossland's radiative process. For both cases, stretching & shrinking, a significant rise is noticed in thermal distribution because thermal radiation predominates significantly over the fluid's thermal conduction property, resulting in increased temperature distribution and corresponding thickness of the boundary layer. Figure 8 exhibits the development in the temperature distribution due to elevated values of the generalized Prandtl number. Since a higher Prandtl number appears to keep the lower thermal conductivity of the material under consideration, known from the constituency term of the Prandtl number given in Table 1, therefore, the conduction declines. Ultimately the situation gives rise to the boundary layer thickness, and a declining trend is noticed in the temperature distribution profile. Figure 9 reflects the temperature distribution profile for rising values of the Brownian diffusion parameter. The activation of Brownian diffusion creates a hassle for the nanoparticles to find their colder region, thus colliding more rapidly and resulting in the rise of the thermal state. The boundary layer thickness diminishes accordingly. Figure 10 reflects the trends of temperature distribution for elevated values of the thermophoresis parameter. The temperature distribution results in an inevitable rise due to the incident in predictive force activated within the fluid, resulting in rapid collisions and a heat transfer process. For both the cases, shrinking and stretching, we have noticed that the concentration of the particles shows particular rising behavior for the range 0 ≤ N Th ≤ 1, but the impact of thermophoresis becomes negligible close to 0, as portrayed in Fig. 11; however, the results for both cases are opposite for the Brownian diffusion parameter for the range of parametric values for the concentration profile. For both shrinking & stretching, we notice the decline in concentration profile is prominent, even close to the lowest possible value for the Brownian diffusion parameter given in Fig. 12. The impact of Lewis number on the concentration profile is portrayed in Fig. 13 reflects a certain decline in the concentration of the nanoparticles. Since the constituency of Lewis number is given in Table 1. Shows a ratio of thermal diffusion & species diffusion rate within the given geometry. Thus, incremental values of the Lewis number reduce the thickness of the thermal boundary layer, resulting in a cumulative trend in the concentration profile. Numerical findings of skin friction and Nusselt number are present in Table 3. The corresponding value for SF and Nu is calculated for different parameters for each altered value. Nusselt shows an augment for higher values of radiation factor and Brownian diffusion. Skin friction rises for porosity factors.  www.nature.com/scientificreports/ Optimization analysis. To determine the optimal response of heat transfer rate and the skin friction coefficients for the efficiency of specific components, a recently developed statistical method known as Response Surface Methodology (RSM) is used. The current simulation is dependent on thermal radiation, Brownian motion, and thermophoresis variables for the local Nusselt number, and the simulation is dependent on the magnetic field parameter, Sisk-fluid parameter, and Forchheimer parameter for the skin friction coefficient, respectively. The simulation process is carried out within the allowed range of these parameters. These factors are considered when determining the heat transfer rate as 0.
To enforce the suggested second-order model,     www.nature.com/scientificreports/ we use the CCD's in-built linear, quadratic, and interaction terms with low, medium, and high values, as the face-centered design method recommended. These ranges with their corresponding levels for the components  Nusselt number and the skin friction coefficient are arrayed through Tables 4 and 5, respectively. The distribution of the corresponding answers obtained from 20 separate runs employing the various parameters is shown in Tables 6 and 7, which are given one after the other. In addition, the multivariate models for the response function in terms of the relevant elements are provided by the entire quadratic polynomial as We then reach the following response functions by removing the insignificant components from the response function and calculating the regression coefficients.  Table 5. Range of independent factors and their levels to assess Skin friction.

Validity of the model. A statistical instrument that is reliable. Analysis of Variance (ANOVA) is used
to assess the regression models and the different statistical tests with F-values, p-values, lack of fitting, error, and total error for the simulated outcomes of Nusselt number and the skin friction coefficient is deployed via Tables 8 and 9. The results of the tests on the suggested data set are produced for either a significance level of 5% or a degree of confidence of 95%. The previous research demonstrates that to get a superior accuracy model, it is always preferable to choose the F-values that are more than one and the p-values that are lower than 0.05. In addition, based on the computed result that is displayed in Table 8 (which was obtained from the ANOVA), and depending on the range of the factors with respective p-values, the proposed model of the Nusselt number is eliminated along with the interaction terms, and the expression for this is shown in the body of the text. Similarly, the model of skin friction is provided by excluding the interaction variables in the manner directed in Table 9 (obtained from the ANOVA). As a result, a best-fitting model for the responses has been projected. This confidence in the model's ability to optimize and forecast has been established by the computed and adjusted values of 100%.
The normal plot for the distributions of the data with the fitting values is shown in Figs. 14 and 15, as well as the histogram for the same distribution and the order of the data distributions for the responses of Nusselt number and the skin friction coefficients, respectively. Due to the fact that the residuals are located close to the straight line, it is evident from these numbers that the model is appropriate.    Fig. 17a,b that the significant enhancement in the skin friction coefficient is marked with higher levels of, lower levels of, and medium levels of. Figure 17c,d exhibits the impact of M 1 and Fr affects the skin friction coefficient rate. It is witnessed that the enhancement in the response function is observed at higher levels ofand lower level of magnetic field parameters. Figure 17e,f portrays the role of and Fr on the skin friction coefficient enhancement while M 1 maintaining the middle level. The contours and surface depict that the skin friction coefficient is enhanced significantly for the increasing increase.  www.nature.com/scientificreports/  The present study investigates the sensitivity of the heat transfer rate and the shear force rate to modest changes in parameters under controlled settings. It is defined as the partial derivative of the responses, which are the local Nusselt number and the shear rate, concerning the characterizing parameters (Rd, Nb, Nt) , which are referred to as coded variables. Calculating the partial derivative concerning the effective parameters looks like this: www.nature.com/scientificreports/ . Furthermore, at the medium level of thermal radiation, the transport of heat sensitivity (Rd, Nb, Nt) declined whenever the Nb and Nt increased from low-level to high-level values. It is also observed that the magnetic field parameter and Forchheimer parameter controls the rise in the sensitivity of the rate of friction factor M 1 and Fr . In contrast, an opposite trend is celebrated with the increase in the Sisk-fluid parameter.

Conclusions
The current study provides numerical approximations for a steady, incompressible, two-dimensional Siskonanofluid flow on a stretching/shrinking surface equipped with power-law component (Sisko model), magnetic impact, thermal radiation, Brownian diffusion, and thermophoresis. In addition to boundary layers formulations, optimization analysis of the heat and mass transfer attributes is performed using Response Surface Methodology. The followings are the key outcomes of the study: 1. The influence of the involvement of porous medium in fluid flow analysis is prominent. Forchheimer number results in an inevitable decline in fluid flow. And so is the case with values of porosity parameter fluctuating between 0 ≤ ≤ 2.5 . At = 0 , the fluid flow medium returns to a simple medium having no porosity and no additional frictional force. 2. Temperature distribution receives enhancement for larger thermal radiation, Brownian diffusion, and thermophoresis, which is quite logical from their constituency; however, a more significant Prandtl number results in a low thermal profile. 3. For both the cases, shrinking & stretching, we have noticed that the concentration of the particles shows specific rising behavior for the range 0 ≤ N Th ≤ 1, but the impact of thermophoresis becomes negligible close to 0. However, the results for both cases are the opposite of Brownian diffusion. 4. The involvement of Darcy-medium is quite helpful in controlling the fluid movement and balancing the fluid's thermal state for various manufacturing procedures while maintaining the nanofluid's low viscosity and higher thermal conductivity. 5. The coefficient of heat transfer is significantly achieved at higher levels of thermophoresis and lower levels of thermal radiation and Brownian motion factors. The surface friction rate at the plate boundary is dominant at higher values of Sisk-fluid parameter and Forchheimer parameter factors and lower levels of magnetic field parameter. 6. The transport of heat has the highest sensitivity value at the levels of (A = 0, B = −1, C = −1) Rd = 0.2, Nb = 0.4 and Nt = 0.1 while the highest sensitivity value of the rate of friction factor is observed at the levels of (A = 0, B = 1, C = −1)M 1 = 0.2, = 0.5 and Fr = 0.2. 7. The transport of heat sensitivity (Rd, Nb, Nt) declined whenever the Nb and Nt increased from low to high and at the medium level of thermal radiation. Increasing the magnetic field parameter and the Forchheimer parameter increases the sensitivity of the rate of friction factor, whereas increasing the Sisk-fluid parameter has the reverse effect. www.nature.com/scientificreports/

Data availability
All data generated or analyzed during this study are included in this published article.