Scalable Heat Transfer Characterization on Film Cooled Geometries Based on Discrete Green’s Functions

Once the wall temperature is characterized, the heat flux can be computed from the temporal evolution of the temperature following the impulse response method as defined by Oldfield [8] or solving the internal heat conduction as presented by Saavedra et al. [9] or Solano and Paniagua [10] for 1D or 2D geometries, respectively. Finally, the adiabatic wall temperature is computed based on heat transfer measurements at various wall temperatures. For a constant external flow reference temperature, the heat flux decreases as the wall temperature rises. Extrapolating the heat flux decay for various wall thermal conditions, the adiabatic wall temperature is defined as the intersect with the abscissa on the adiabatic condition,$q˙=0$⁠, as represented in Fig. 1. Different approaches have been used to extrapolate the trend toward the adiabatic condition, considering linear evolutions, higher order polynomial fits or exponential fitting [11–13].

Yu et al. [14] relied on a transient liquid crystal technique to identify the distribution of film effectiveness and heat transfer coefficient. In this case, the different heat flux levels for various wall thermal conditions are recorded during the transient experiment. Instead of setting the test article at various temperatures prior to the experiment, they take advantage of the solid heat up during the test. As the hot fluid flows over the test article, it gets warmer and the heat and wall temperature differences are directly monitored during the transient.

Bons et al. [15] studied the adiabatic wall effectiveness on a single row of film-cooling holes mounted on a zero pressure gradient flat plate by means of local wall-mounted thermocouples. Instead of performing multiple experiments at different wall temperatures, to obtain the adiabatic wall temperature they run the facility without coolant for a long period of time and registered the wall temperature once the system reached thermal equilibrium. Though, once active, the coolant modifies the near-wall flow region and the actual adiabatic wall temperature distribution differs from the one measured at thermal equilibrium. Friedrichs et al. [16] exploited a different approach to estimate the film-cooling effectiveness distribution. In their case, they used the ammonia and diazo surface coating technique. Unfortunately, the quantitative description of film coolant effectiveness with ammonia and diazo surface coating requires an in situ calibration, which requires the use of one of the previous methods.

The experimental characterization of the adiabatic effectiveness is commonly performed at reduced temperature levels while maintaining similar Mach and Reynolds numbers, temperature/density ratios, blowing ratios (BRs), capacity ratios, etc. [17]. However, when scaled to engine-like conditions, the adiabatic effectiveness magnitude may deviate from the actual behavior, partially due to viscosity effects. Moffat [18] highlighted the influence of the driving reference temperature to define the convective heat transfer coefficient or the adiabatic effectiveness in this case.

The matrix of coefficients g_ij is defined as the DGF sensitivity matrix. This procedure has already been applied to fundamental studies [20], convective heat transfer in microelectronics [18,21], turbine passages [22], and blade tip heat transfer [23]. Andreoli et al. [23] applied the DGF approach on a turbine blade tip and proved its utility for the adiabatic convective heat transfer coefficient on conditions where the derivation of the adiabatic wall temperature or the driving temperature may become ill posed. The DGF approach enables the accurate extrapolation of convective heat transfer measurements performed at lab conditions to engine settings suffering less error than traditional methods, while enabling the extrapolation in configurations or conditions where retrieving the adiabatic wall temperature or recovery temperature may be unsuitable.

Batchelder and Eaton [24] used the DGF to predict the heat flux for arbitrary thermal boundary conditions using liquid crystal imaging and a heated strip. With the DGF measurements, they were able to measure the heat transfer enhancement due to the presence of strong freestream turbulence. The experimental results were compared against established correlations and numerical simulations. Using a similar approach, Booten and Eaton [25,26] presented the application and DGF on internal flows and cooling passages. In this case, Booten and Eaton used locally distributed thermocouples to characterize the wall temperature distribution, demonstrating its advantages over traditional methods on nonuniform thermal boundary conditions. Additionally, Ortega et al. [27] explored the application of the DGF characterizing the convective heat transfer using thermochromic liquid crystals. They heated a small region in the wall to mimic a source heat and derive from the sensitivity matrix.

The present article aims to define a practical way of experimentally retrieving the DGF coefficients with minimal intrusion and establishing the applicability to film-cooled geometries. Initially, 2D conjugate simulations over a flat plate and a slot-cooled flat plate are used to perform parametric analysis on the DGF coefficient retrieval. This parametric analysis is followed by 3D conjugate simulations over a flat plate and a film-cooled flat plate to quantify the accuracy of the method and its extrapolation to different thermal boundary conditions. Finally, an experimental campaign on a flat plate is carried out to assess the applicability of the developed method and to quantify its accuracy.

The individual geometry elements can be excited with a heat pulse, and then the wall temperature increments are measured along the surface. In a sense, the local heat addition introduces a tracer in the domain. Due to the local heat addition, the flow particles that pass over the excited element increase their temperature, and their path can be characterized looking at the increased wall temperature footprints as they follow their streamline. A geometry consisting of 100 elements will require 100 experiments with local heat addition on each one of the partitions of the domain. Once the DGF coefficients are characterized, the heat flux can be retrieved with Eq. (4) for any thermal boundary condition as long as Reynolds number, Mach number, and temperature ratios are maintained.

The local heat addition that individually excites each element may be achieved by surface heaters, heat gun application, or active laser heating. In this article, a burst-mode laser was used as heat source because of its minimal intrusion on the test article with the unique requirement of optical access. The local heat addition is introduced by firing the burst-mode laser focused at a given spatial location. The target region is easily adjusted by tilting the mirrors guiding the laser toward the test article, while the heat pulse magnitude can also be modulated by adjusting the power of the laser source. The surface temperature can be monitored by infrared thermography [28], temperature-sensitive liquid crystals [5], thermographic phosphorescence [29], thin films, or local wall-mounted thermocouples.

2D conjugate numerical simulations were used to perform parametric analysis on the DGF matrix calculation and qualify its accuracy on assessing the heat transfer distribution over a flat plate. 3D conjugate simulations were used to assess the experimental applicability of the method and its accuracy on film-cooled geometries.

The conjugate unsteady Reynolds-averaged Navier–Stokes (URANS) simulations were carried out with ansys fluent. Both 2D and 3D numerical simulations used the kω-SST (shear stress transport) turbulence closure [30,31]. Figure 2 depicts the different numerical domains used for the parametric analysis and practical assessment. Figures 2(a) and 2(b) illustrate the flat plate numerical domain with and without cooling. The inlet boundary condition is modeled as a pressure inlet boundary where total pressure and temperature are prescribed. The top boundary of the numerical domain is treated as an adiabatic slip line. The outlet is modeled as a pressure outlet, where the static pressure level is set, which reproduces the discharge to a constant pressure reservoir. The border between fluid and solid numerical domain is modeled as a full conjugate viscous wall, with the heat transfer continuity enforced at each time-step. Front and back solid boundaries are treated as adiabatic surfaces, while the bottom surface is modeled as a natural convection boundary with a heat transfer coefficient of 20 W/(m² K) and an external driving temperature of 300 K. The numerical domain has a length of 0.5 m, while the fluid and solid sections are 0.15 m and 0.02 m high, respectively.

For the 2D cooled flat plate geometry, a slot is introduced 0.15 m downstream of the leading edge. The slot has an axial gap of 3.175 mm and a length over slot width ratio of 6. The coolant is introduced with a 30-degree injection angle with respect to the plate. The blowing and intensity ratios are controlled by the total temperature and pressure imposed at the bottom boundary of the slot. The solid boundaries between the slot and the plate slab are also treated as viscous full conjugate boundaries.

Figures 2(c) and 2(d) depict the 3D numerical domain for the flat plate and the film-cooled flat plate respectively. The length and height of the numerical domain are identical to the 2D numerical setup. In this line, the boundary conditions for the 3D numerical domains are maintained, and the inlet and outlet are used to prescribe the main flow conditions. The upper boundary is modeled as an adiabatic slip line, and any boundary between fluid and solid is treated as viscous full conjugate boundary. The front, back, and bottom surfaces of the solid boundary condition definition are also preserved with adiabatic and natural convection conditions. However, the lateral boundaries for both fluid and solid are treated as periodic boundary conditions. The width of the numerical domain is limited to 0.05 m. The film-cooled flat plate numerical domain introduces a 30 degree cylindrical orifice 0.15 m downstream of the inlet of the domain. The film-cooling orifice has a diameter of 3.175 mm and a length-to-diameter ratio (L/Φ) of 6. The bottom surface of the orifice was modeled as an inlet pressure boundary condition, where the total pressure and temperature of the coolant were prescribed. Injection blowing and momentum flux ratios are controlled through the total pressure and temperature at the inlet of the orifice.

The working fluid was air, treated as an ideal gas, while the solid domain was modeled as aluminum. The Sutherland law was used to model the temperature effect on the molecular viscosity. Similarly, temperature dependent conductivity models were introduced for solid and fluid domains. Second-order upwind formulations were adopted for flow and turbulent properties. The time-step and number of inner iterations of the transient simulations, 5 × 10⁻⁶ s and 28, respectively, were selected based on benchmark analysis. The transient progression of the simulation was controlled by second-order bounded implicit formulations.

The numerical domains were discretized following a blocking strategy with ansys icem. To guarantee the results independence from the spatial discretization, the approach presented by Celik et al. [32] was followed. The y+ was maintained below 0.1 for all the wetted surfaces. In this line, Figs. 3(a) and 3(b) represent the heat flux distribution along the centerline for both uncooled and cooled 3D numerical domains, respectively. Four different grid levels are compared for each case. There are only minor differences among all the grids around the boundary layer transition region and by the injection orifice for the cooled geometry. Following Celik’s guidelines, the Grid Convergence Index of the fine grids for both cases is below 1e⁻³, guaranteeing proper mesh convergence. A similar procedure was followed for the 2D numerical domains, and the grids with Grid Convergence Index below 5e⁻³ were selected.

A verification case was also evaluated to characterize the performance of the numerical setup on resolving the heat flux downstream of a film-cooled flat plate. The adiabatic wall effectiveness measurements of film-cooled holes by Gritsch et al. [33] were used as a validation case. The numerical setup replicated Gritsch et al. experimental conditions: external flow Mach number 0.6, 2% turbulence level, L/Φ ∼6, 30 deg cylindrical orifice, temperature ratio of 0.54, BR of 0.6, and coolant to main flow momentum flux ratio (I) of 0.2. Additionally, a laminar boundary layer profile of 1 mm height was prescribed at the inlet of the numerical domain, which led to a δ/Φ of 0.5 upstream of the orifice, equivalent to Gritsch’s experimental conditions. Figures 4(a) and 4(b) depict the streamwise velocity and temperature contours along the mid-span of the numerical domain. Similarly, Fig. 4(c) represents the adiabatic wall temperature distribution along the plate. The coolant injection through the orifice and its trace can be clearly distinguished on the contour, identifying high cooling efficiency just downstream of the injection and how it is diluted as the coolant is mixed with the main stream.

Figure 4(d) compares the film-cooling effectiveness along the mid-span against the experimental results of Gritsch et al. [33]. The experimental results for equivalent flowing and momentum flux ratio of Schmidt et al. [4] and Schmidt and Bogard [34] and Sinha et al. [35] are also introduced in the graph. The external flow Mach number of Schmidt’s and Sinha’s experimental setup was below 0.1. The influence of the external Mach number is attributed to the difference on the effectiveness close to the injection. However, the dissimilarities between the various Mach number cases are reduced further downstream. Overall, the numerical prediction follows the trend and magnitude characterized by Gritsch et al., which verifies the applicability of the current numerical setup to flat plate film-cooled analysis.

An experimental campaign over an uncooled flat plate was carried out at the Purdue Experimental Turbine Aerothermal Laboratory (PETAL) linear wind tunnel [36] to verify the applicability and accuracy of the developed DGF approach. The Linear Experimental Aerothermal Facility (LEAF) layout is presented in Fig. 5(a). The facility is operated in a blowdown mode, with pressurized dry air stored in reservoirs at 150 bars. The air is withdrawn from the tanks following two different lines. A cooled stream is directly delivered to the test cell, and a second line is diverted through a heat exchanger to provide nonvitiated hot air. The streams are mixed upstream of a critical flow Venturi and delivered to the settling chamber upstream of the test section. While setting up the flow conditions, the air is vented outside through an auxiliary purge line. When the mass flow and temperature are stable, the air is re-directed through the settling chamber and test section. Figure 5(b) illustrates the LEAF test section for the uncooled flat plate experiments. The aluminum flat plate was coated with a layer of black ceramic paint, Cerakote^®, to enhance the temperature rise following the laser heat addition and to allow accurate infrared thermography. The left side of the test section was covered with an aluminum blank while the upper surface accommodated three different inserts. Upstream and downstream inserts were used to introduce total flow measurement probes, Kiel pressure probes, and exposed K-type thermocouples. A fused quartz window was used in the central insert to deliver the laser beam. The right side of the test section was sealed with a large quartz window to allow optical access onto the test article. Downstream of the test section, a sonic valve allows the independent control of Mach and Reynolds numbers when choked. Finally, the air is discharged to a large reservoir that can either be in vacuum or atmospheric conditions.

Figure 5(c) illustrates the test section instrumented with wall-mounted thermocouples to monitor the plate temperature along the experiment. Figure 5(d) depicts the location of the flush mounted thermocouples as well as the streamwise and spanwise location of the total pressure and temperature sensors. Thirty-five exposed welded K-type thermocouples were used to monitor the wall temperature evolution along the plate. The thermocouples were attached to the surface with a thermally conductive paste. The three thermocouples close to the leading edge of the plate were used to measure the wall temperature evolution on the laminar boundary layer portion on the plate, while 22 additional thermocouples were used close to the left side to track the temperature evolution along the transitional and turbulent region. On the opposite side, closer to the right side of the plate, ten additional thermocouples were installed to discretize the central part of the plate for the DGF sensitivity matrix characterization. The center of the plate is divided in ten different regions, as shown in Fig. 5(d). Each one of those regions is excited with a laser heat pulse, and the temperature rise is monitored through the thermocouple array. Additionally, 15 more thermocouples are installed on the opposed side of the plate at the bottom surface of the test article.

The mean flow conditions of the working section were monitored through the mass flow measurements at the critical flow Venturi and by the local Kiel probes and thermocouples measuring the total flow conditions in the test section. The wall temperature evolution is measured with the surface mounted thermocouples at a rate of a 1000 samples per second. Additionally, the heat flux is derived from the temporal signature of the surface thermocouples solving the solid heat conduction problem as presented by Saavedra et al. [9].

Figure 6(a) illustrates the numerical domain used to retrieve the heat transfer based on surface temperature measurements. Both solids are modeled with a perfect junction without loses at the interface. The upper boundary condition on the heat flux calculation is retrieved from a local thermocouple signal as depicted in Fig. 6(b). On the other side, the opposite boundary condition is imposed based on the interpolation of the thermocouple readings distributed along the bottom surface of the plate. A procedure similar to the one adopted by Cuadrado et al. [37] was used to obtain the optimum location of the 15 thermocouples that will allow an accurate reconstruction of the bottom surface temperature. Figure 6(c) displays the temperature contour of the bottom surface at a given time-step taking advantage of Kriging interpolation among the distributed thermocouples.

A frequency-doubled burst-mode Nd-YAG laser was used to impart the impulse heat flux on the flat plate needed to identify the G-matrix,. The burst-mode laser (Spectral Energies, Quasimodo) provides ∼100 pulses with a 100 ns pulse width at 532 nm with a laser beam diameter (FWHM) of ∼8 mm. The output beam from the laser diverges slowly, until reaching the surface of the test article and expanding to a final diameter of ∼13 mm (FWHM). The burst duration is set to 10 ms with a separation duration of 100 µs between each pulse, corresponding to a 10 kHz repetition rate. The measured total energy available at the surface was ∼6 J over the burst duration of 10 ms, which corresponds to a fluence of 3.7 J/cm² at the surface of the test article. Multiple 532 nm coated mirrors with high reflectance (over 99%) were used to guide the laser beam to the test article. The last mirror (M2) was placed above the flat plate, at a distance of 1000 mm, and was used to position the laser beam between each of the ten surface mounted thermocouples. Figure 7, shows the schematic representation of the laser delivery. Focusing the laser beam on the plate reduces the fluid energy absorption, minimizing the influence of the beam on the air molecules. Due to local wall heating, the fluid properties above the active element are altered which affects the behavior of those particles. However, this phenomenon is localized and considering the reduced exposure time of the fluid particles over the heated element its effects are minor.

In this section, the computational domains will be used to analyze the influence of the control parameters involved on the sensitivity matrix calculation. 2D and 3D conjugate URANS simulations are used to develop and establish the procedures to guarantee a successful application of the DGF approach on an experimental campaign. Equation (6) indicates the formulation to retrieve the DGF matrix computation. The dependent parameters are: the spatial discretization, heat pulse magnitude, and duration of the heat pulse. The 2D uncooled numerical domain was used to evaluate the impact of each one of those parameters on the accuracy of the heat flux calculation with the DGF matrix.

The baseline case has a mean flow Mach number of 0.27 imposed through a total pressure of 105 kPa at the inlet and 100 kPa at the outlet. The inlet total flow temperature was set 500 K while the initial solid temperature was set 300 K. To obtain the baseline condition, the fluid domain was initiated at stagnant conditions with an air temperature identical to the solid temperature of 300 K. The flow was then started and maintained for 5 s. The result after 5 s of flow time simulation was then stored and used as the baseline, as depicted in the temperature contour of Fig. 8(a).

Each one of the actively heated cases used to characterize the components of the DGF matrix was started with the baseline data. A user-defined function was introduced to prescribe additional heat flux on each domain partition, mimicking active heating with a laser beam. The domain was initially discretized in 25 different elements, equally spaced Δs ∼ 0.02 m, as depicted in Fig. 8(a). Figures 8(b) and 8(c) represent the surface heat flux and wall temperature distribution for two independent cases. The blue curve represents the results for the case where element 5 is actively heated. The heat flux on the 5th partition was increased by 10 kW/m², and its impact on the downstream and upstream elements heat flux and temperature was registered. Similarly, the black curve represents the heat flux and wall temperature while actively heating the 20th element. The temperature difference and heat flux for each one of the individual cases was then used to retrieve the components of the DGF matrix taking advantage of Eq. (6).

Figures 9(a) and 9(b) depict the heat flux distribution reconstruction and average relative error of the DGF method, respectively, for various heat flux magnitudes. The domain was discretized in 50 different elements equally spaced for eight different heat flux magnitudes. The 10 s duration of the heat flux pulse was also identical, and the data was processed on the last time-step while the active heating was still active. In the heat flux distribution along the plate in Fig. 9(a), there are various DGF predictions compared with the baseline result. The reconstructed heat fluxes are obtained through Eq. (4) and using the baseline temperature difference. The sensitivity matrix obtained with the heat flux magnitude of 12.5 kW introduces a wavy pattern in the heat flux prediction for x ∼ 0.275 m and x ∼ 0.425 m. In fact, all the cases represented in Fig. 9(a) exhibit oscillations around the baseline in those locations. Figure 9 (b) summarizes the average relative error with respect to the baseline case for all the magnitude levels analyzed.

Figures 9(c) and 9(d) represent the heat flux prediction based on the DGF matrix for various heat pulse durations and the relative error, respectively, when compared with the baseline case. The heat flux magnitude was maintained at 10 kW/m² for four different durations such that the energy delivered ranged from 250 to 480 J per unit width. The prediction for the heat pulse durations above 6.25 s displayed again oscillatory or wavy patterns, as can also be identified on the area-average relative error representation. The results for both heat flux magnitude and pulse duration indicate that there is a range of energy required for the DGF method to perform adequately. A low energy input may not be sufficient to excite a given element and track its influence on the downstream or adjacent components. However, an excessive amount of energy may lead to erroneous predictions driven by excessive conduction, which affects the accuracy of the DGF approach. Once an element is actively heated, its surface temperature increases, which then heats the air flow, allowing the tracking of the fluid around the adjacent elements.

Figures 9(e) and 9(f) display the heat flux and average error, respectively, for various spatial discretization with uniform spacing. The pulse heat flux magnitude and duration were constant for the different spatial resolutions at 10 kW/m² and 10 s, respectively. The different distances between partitions for similar heat flux pulses affect the energy delivered to each element. The results indicate that for the given heat flux source terms, a spatial resolution beyond 23 elements exceeds the energy threshold on the convective to conductive heat transfer. On the other hand, the slight error increase from 23 to 18 elements indicates that the energy delivered on the 18 elements case may not be sufficient to fully characterize the convective process.

Based on the lessons learned on the parametric analysis, a cooled 2D flat case was evaluated to assess the applicability and accuracy of the DGF method on cooled configurations. Similar external flow conditions were set: 0.27 Mach number, inlet total temperature of 500 K, and inlet total pressure of 105 kPa. The injection is introduced with a total temperature ratio of 0.7 and a velocity magnitude of 25 m/s. The domain was discretized in 29 different partitions. Fourteen elements were equally spaced upstream of the slot and 15 additional sections downstream of the slot. The first element downstream of the injection had similar spacing to the upstream partitions, and the downstream divisions were set with the spacing increasing at a growing ratio of 1.1. Figure 10(a) depicts the temperature contour for the baseline case and spatial distribution of the different elements for the DGF method.

The baseline condition was retrieved after 5 s of simulation starting at stagnant conditions with the solid and fluid temperatures at 300 K. In this case, local active heating was used to mimic the firing of the burst-mode laser. The burst-mode laser can fire with a frequency rate of 10 kHz during each 10 ms burst. To model the laser discharge, the heat addition process was limited to 10 ms and the power was increased to 125 kW/m². Assuming continuous firing during 10 ms, the total energy delivered to each element is approximately 50 J/m. To process the surface heat flux and temperature toward the DGF matrix identification, the surface was monitored at a frequency of 1 kHz. The datasets from each one of those time instances can then be processed to retrieve the g_ij matrix components.

As identified on the parametric analysis, depending on the amount of energy deposited on each element, the accuracy on estimating convective heat flux may vary. In this sense, Fig. 10(b) depicts the heat flux prediction with the datasets of various time-steps compared with the baseline result. Excluding the results based on the 12 ms surface measurements, the predictions accurately follow the baseline results. Figure 10(c) represents the relative error of the heat flux prediction. The error of the DGF method obtained while the laser is being fired is below 6%, identifying the applicability and accuracy of the method. The data set coinciding with the end of the active laser heating actually deviates from the baseline beyond 27%. Afterward, during the next 4 ms, the energy delivered to the local elements is still high enough to derive the DGF matrix and the error is gradually reduced to 6%. However, as the heat flux impulse effects diffuse, the error slightly increases for processing times greater than 20 ms.

The surface roughness plays a relevant role on the convective heat transfer augmentation. The current computational fluid dynamics (CFD) simulations were performed assuming smooth wall conditions. However, in real applications the surface finish may not grant a smooth wall. To enable proper convective heat transfer characterization, the experimental or modeled geometry must replicate the surface roughness of the actual target configuration. Under identical or scaled roughness conditions the DGF approach will be able properly scale the convective heat flux. However, if the actual geometry roughness is neglected, an error associated with the roughness heat transfer augmentation will be present on the DGF extrapolation.

After assessing the method accuracy and applicability on 2D geometries, the DGF approach was extended to 3D configurations. Initially, the DGF characterization was applied to a 3D uncooled flat plate. The main flow conditions are identical to the 2D uncooled flat plate: Mach number of 0.27, T₀ of 500 K, and static pressure of 100kPa. Figure 11(a) displays the temperature contour at mid-span plane, outlet, conjugate plate, and back wall of the solid. The baseline condition was obtained after 5 s of transient simulation, initiated at stagnant conditions and with a fluid and solid temperature of 300 K.

Figure 11(b) displays the temperature contour at the conjugate plate for the baseline condition. The spatial discretization used for the DGF method is displayed on top of the contour, with five divisions on the spanwise direction and 22 in the streamwise direction. A total of 110 individual simulations are required to find the components of the DGF matrix in this case with a size of 110 × 110. Each one of the separate heat pulse cases is initiated based on the baseline condition. The heat flux pulse is uniformly imposed on each element taking advantage of a user-defined function. The power of the source term is equivalent to 125 kW/m² and maintained during 10 ms, delivering a total energy of 20 J on the largest elements.

Figure 12(a) represents the temperature contour along the top surface of the plate while actively heating element 7 at time t = 5 ms. The temperature in the element 7 is greatly increased as a consequence of the heat pulse. The air flowing above the element 7 is heated and leaves a trace of increased wall temperature as it moves downstream. To process the DGF method, the conjugate wall data are reduced to the 22 × 5 spatial discretization following an area-average process. In this line, Fig. 12(b) represents the surface temperature reduction at t = 5 ms while actively heating element 7. The wall-to-inlet total flow temperature difference and the heat flux are introduced in Eq. (6) to retrieve the g_ij components. Figure 12(c) depicts the relative error of the heat flux prediction for various processing times against the baseline result. The data set after 1 ms of heat pulse has a 10% error driven by the lack of convective energy transport to adjacent elements. While for the data sets processed between 2 and 5 ms the relative error is below 3%. For later processing times the conduction starts to affect the g_ij matrix characterization and the error slightly increases to ∼4.5%. The choice of processing time, after the heat input, is driven by an optimization of the time at which there is influence of maximum convective effect and minimum conductive effect.

Figure 13(a) depicts the DGF matrix magnitude contour for the dataset based on the 3 ms processing time. The largest absolute magnitude of the matrix is present in the diagonal, displaying the self-influence of each element, which implies that each domain partition is most sensitive to its own local temperature gradient. Figure 13(b) represents the heat flux distribution along the plate for the baseline case. As the boundary layer develops, the heat flux level is reduced in the first portion of the plate at x < 0.05. The boundary layer transition process starts at x ∼ 0.05 and develops until x ∼ 0.11. The heat flux magnitude raises along the transitional region 0.05 < x < 0.11 and then decays as the turbulent boundary layer grows, x > 0.12. Figure 13(c) depicts the heat flux prediction obtained with the baseline wall-to-inlet total temperature difference and the g_ij matrix displayed in Fig. 13(a). In this line, Fig. 13(d) depicts the relative error of the prediction when compared with the baseline case. The area-average relative error is ∼1.8%, highlighting the accuracy of the DGF method on predicting the heat flux level. The largest sources of error are present around the transitional region driven by the lack of spatial resolution. Additionally, the error on the first streamwise components is also larger due to the lack of any upstream element to retrieve their sensitivity. Using the DGF method, the influence of heat addition on adjacent elements is used to derive the sensitivity and hence predict the heat transfer. Since the first elements of the domain lack of any upstream flow sensitivity, the heat flux prediction on the first axial elements is constrained.

Ultimately, the DGF method is evaluated based on film-cooled flat plate simulations. The external flow was imposed with a total pressure of 105 kPa and total temperature of 500 K, while prescribing a static pressure of 100 kPa at the outlet. The density ratio and blowing ratio of the injection were kept at 1.4 and 0.3, respectively, for a coolant to external flow total temperature ratio of 0.7. The outer total flow to initial wall temperature ratio was set to 0.6.

The baseline simulation was obtained after 5 s of flow time simulation starting from stagnant conditions at a fluid and wall temperature of 300 K. Figure 15(a) represents the temperature contour at mid-span, conjugate plate, solid bottom wall, and outlet for baseline conditions. Figure 15(b) represents the heat flux contour along the wetted surface of the plate. The spatial discretization used for the DGF approach is superimposed on the heat flux contour, with refined partitions around the injection orifice. Eighteen different divisions are considered along the axial direction and five across the spanwise direction for a total of 90 different elements. The heat flux pulse is introduced by a user-defined function. In this case, instead of assuming uniform heating on each element, a 3 mm radius Gaussian distribution is introduced centered on each partition. The Gaussian distribution with a 3 mm radius simulates the diffused laser beam impact on the surface. The heat addition is maintained during a 10 ms pulse, and the magnitude of the heat pulse decays from a maximum of 125 kW/m² at the center according to a Gaussian profile. The total energy delivered per element is 8 J. The heat pulse magnitude and the spatial distribution of the heat source combined with the pulse duration define the total energy delivered per element. There needs to be a compromise between the energy delivered per element and the size of the partitions. In this sense, heat pulse magnitude, duration, spatial distribution of the heat pulse, and spatial discretization need to be defined accordingly. If an excessive amount of energy is delivered to a refined discretization, the conduction to the adjacent elements in the domain may become dominant, preventing the correct characterization of the sensitivity matrix. On the other hand, if the total amount of energy delivered results minimal, the temperature rise of the local element will be insufficient to set the trace to characterize the sensitivity. In terms of heat pulse spatial resolution, it must fit within the element size and its aimed area must be set in agreement with the maximum heat pulse magnitude to prevent excessive energy delivery.

Figures 15(c) and 15(d) depict the baseline and predicted heat flux contour on the DGF domain. The prediction is calculated based on the g_ij matrix obtained 4 ms after initiating the heat pulse and with the baseline total to wall temperature difference. Figure 15(e) summarizes the relative error of the prediction along the discretized region when compared with the baseline results. The area-average relative error is below 1%. The strongest sources of error are clustered around the injection orifice, driven by the lack of sufficient spatial resolution. The satisfactory results in this test case prove the applicability and accuracy of the DGF method on characterizing the convective heat flux on film-cooled geometries.

Two additional baseline simulations are performed to test the scalability of the DGF method under different thermal boundary conditions.

Table 1 summarizes the main flow parameters of the three different cases evaluated. Reynolds and Mach numbers were maintained at 1.54 × 10⁶ and 0.266, respectively. The same geometry was studied for the three different cases and the Reynolds scaling was achieved by modifying the static pressure level. Similarly, the coolant density ratio and blowing ratio were kept constant at 1.4, and 0.3, respectively. The external flow total temperature to coolant total temperature was also preserved at 0.7. Similarly, the external total flow temperature to initial wall temperature ratio was matched to 0.6. Both fluid and solid were also identical for the different thermal boundary conditions explored, using air and aluminum, respectively. Figures 16(a) and 16(b) depict the temperature contour at mid-span and the heat flux contour over the wet surface of the plate for the two additional cases evaluated: inlet T₀ = 700 K, and inlet T₀ = 900 K, respectively.

Figures 17(a) and 17(d) display the baseline heat flux contours for both extended thermal boundary conditions on the discretized DGF domain. While Figs. 17(b) and 17(e) represent the heat flux prediction for both thermal boundary conditions using the g_ij matrix obtained with the dataset at 4 ms on the T₀ = 500 K thermal boundary conditions.

Both predictions seem to deviate from the baseline results, and Figs. 17(c) and 17(f) illustrate the relative error distribution along the discretized domain. Figure 17(g) summarizes the area-average relative error for the three thermal boundary conditions using the available processing times. The relative error of the DGF method when applied on the thermal boundary conditions that it was obtained remains below 4% for most of the processing times. However, when used at different thermal boundary conditions the prediction fails on describing the actual heat flux magnitude. The spatial distribution seems accurate, given that the main flow parameters describing the external flow and the coolant are identical. However, the area-average error rises to 18% for the inlet total temperature of 700 K, and it is worsened for the 900 K case with an average error around 28%.

For comparison purposes, Fig. 18(a) displays the Nusselt number obtained from the adiabatic wall heat transfer coefficient at baseline conditions, T₀ = 500 K. The approach outlined in Fig. 18(b) was used to retrieve the adiabatic wall temperature and the adiabatic wall heat transfer. Based on two different time-steps, once the simulation was converged the heat flux decay for larger wall temperatures was extrapolated to the adiabatic condition. Figures 18(c) and 18(d) represent the heat flux prediction for the 700 and 900 K boundary conditions using the Nusselt number displayed in Fig. 18(a). The heat flux prediction seems to be inaccurate on both cases. In this sense, Figs. 18(e) and 18(f) display the relative error contour along the discretized region. Besides the overall magnitude deviation, the heat flux prediction based on the heat transfer coefficient seems to miss also some of the physical patterns around the injection orifice. The area-average relative error for the 700 K case extrapolation is 9.2%, while it rises to 19.1% at 900 K. These results also highlight the shortcomings of using the heat transfer coefficient on scaled thermal boundary conditions.

Figure 19(e) summarizes the area-average relative error for the g_ij matrix for various processing times on the three different thermal conditions. After applying the near-wall thermal conduction correction, the DGF prediction can be accurately scaled to different thermal boundary conditions. At higher temperatures, the radiative heat flux magnitude will be noticeable. Under those circumstances, the DGF approach will fail to predict the total heat flux magnitude, missing the radiative term. However, the convective heat flux prediction will still be accurate. At temperatures where radiation plays a relevant role, the heat flux prediction should combine the use of DGF to predict the convective heat flux and an additional model to predict the radiative term. The DGF approach is limited to purely linear convective heat transfer processes. Any nonlinear convective process will not be properly extrapolated by the DGF approach.

An additional set of simulations was carried out on a steel substrate to check the DGF accuracy when a different solid is used. The selected steel has a conductivity of 16.7 (W/(mK)), which is 1/10th of the original aluminum conductivity. The convective process on both geometries is identical, but the different material conductivities may affect to the DGF discretization.

Figure 20 depicts the baseline and DGF heat flux predictions for both solids. Figures 20(d) and 20(e) illustrate the prediction error for aluminum and steel, respectively. The overall area weighted error for the aluminum substrate is 2.5% of the baseline heat flux, while the error on the steel slab is reduced to 1.28%. The lower conductivity on the steel slows down the heat conduction and allows a more accurate characterization of the sensitivity matrix. Based on the conductive time scales estimated by Eq. (8), the lower the heat conduction of the solid, the higher the chances to obtain a proper representation of the convection process with the sensitivity-based method.

An experimental campaign on an uncooled flat plate geometry was performed to study the applicability and accuracy of the DGF method. Table 2 summarizes the main flow parameters for the test conditions that were run. Both tests started with the plate at room temperature, approximately 295 K. The massflow was measured based on the pressure and temperature upstream of the critical Venturi, while the dynamic pressure was measured as the total-to-static pressure difference at the test section core.

Figure 21(a) depicts the massflow stability during the experiments. Throughout the test duration, the pressure upstream of the Venturi was accurately controlled and the massflow was held constant. Similarly, because of the proper performance of the flow conditioning system, the Mach number at mid-span and mid-height of the test section was constant throughout the test, as displayed in Fig. 21(b). Driven by the conduction loses along the upstream piping and the settling chamber, the total flow temperature throughout the test suffers minor variations.

Figure 21(c) represents the total flow temperature at the core of the test section for both experiments. The heater was switched on and set at 500 K for the experiment FP 2.1. As the piping and settling chamber warmed up, the conduction losses were reduced, resulting in a gradual total temperature increase throughout the experiment. When the heater was switched off for the experiment 2.2 to run a lower inlet total temperature, this gradual temperature rise is absent. Experiment 2.2 was run shortly after 2.1, so the thermal inertia on the piping and hardware was enough to increase the flow temperature to approximately 330 K. However, in this case as the hardware is cooled with the airstream, the test section total temperature is gradually reduced.

The agreement between CFD, experimental results, and Sieder-Tate’s correlation is excellent in the turbulent region of the plate, x > 0.22. Both the CFD and experimental results deviate from the correlation upstream on the plate in the laminar and transitional region, where the CFD seems to predict a later transition to turbulent boundary layer. The numerical domain used for these cases is the uncooled flat plate also used in previous sections of the article. Without modeling the contraction area or any of the components of the flow conditioning system, those simulations assume that the boundary layer starts at the leading edge of the plate. While in the test section, the plate ingests the boundary layer that has grown along the contraction area. The disagreement on the laminar and transitional region of the plate at x < 0.22 m is attributed to the differences in the boundary layer at the leading edge of the plate.

The laser heating pulse was used to individually excite each one of the partitions on the DGF domain presented in Fig. 5(d). A line close to the plate mid-span was discretized in ten different elements labeled 1 to 10. Thermocouples were wall mounted on the center of each one of those elements to measure the transient wall temperature evolution and derive the heat flux. The characterization of the DGF matrix for this case required ten different experiments to calculate the 10 × 10 components of the sensitivity matrix. The laser could be fired every 10 s and delivered 6 J on each pulse burst. The laser beam was centered upstream of the thermocouple at each one of the partitions, providing the individual heat flux pulse required to derive the g_ij matrix coefficients.

Figure 23(a) represents the temporal evolution of the temperature at the thermocouples being used to characterize the DGF matrix. There is a clear temperature peak on the fourth thermocouple when the laser is fired in the 4th element. Similarly, some of the downstream thermocouples seem to depict also a temperature increase. Figure 23(b) presents a zoom in during the 4th firing of the laser over the fourth element of the domain. The 4th thermocouple shows a sudden temperature increase, followed by a slightly delayed temperature increase on the fifth thermocouple of similar magnitude. Further downstream, the sixth thermocouple shows a temperature increase of half magnitude shortly after.

The temporal signal of each one of the thermocouples is processed to retrieve the heat flux. Ultimately, the temperature and heat flux values at each location for the ten different experiments are processed toward obtaining the g_ij coefficients. Figure 24(a) represents the magnitude of the sensitivity matrix. Figures 24(b) and 24(c) depict the heat flux predictions compared against the baseline direct measurements and the CFD predictions obtained with similar boundary conditions.

The direct measurement heat flux is acquired with the wall temperature evolution on the thermocouples prior to laser heating on the first element of the domain. The DGF prediction is obtained with the temperature difference between inlet total flow and wall measurements and applying the g_ij matrix represented in Fig. 24(a). The agreement between the DGF prediction and the direct measurement falls within the uncertainty of the measurements and is in close agreement with the CFD calculation. These results verify the experimental applicability of the described DGF method and its accuracy on representing the forced convection heat flux over a flat plate.

A novel engine-scalable method to characterize the convective heat flux on film-cooled geometries based on DGF was investigated. The method is based on the superposition principle and the determination of a sensitivity matrix that relates the heat flux to the temperature gradient between the surface and the external flow total temperature. The sensitivity matrix is identified based on individual experiments where each one of the geometry partitions is excited by a heat pulse. To correctly characterize the convective heat transfer process, the response of the domain to convective heat transfer must be much shorter than the conductive heat transfer. The DGF matrix is obtained at a time interval where convective effects are maximized while still maintaining minimal influence of conduction and radiation heat transfer modes. Thus, the influence of the heat pulse on the adjacent elements defines the components of the matrix and effectively describes the forced convection process. A parametric analysis on simplified 2D geometries identified the range of energy that is required to derive the components of the DGF sensitivity matrix.

The applicability of the method using active laser heating was assessed with 3D conjugate simulations. Similarly, the scalability of the DGF method using a near-wall thermal conductivity correction was also proved by comparing the heat flux prediction based on a matrix obtained at lower inlet total temperatures to warmer thermal boundary conditions. A correction on the near-wall region thermal conductivity was required to account for the air property changes as the thermal conditions are altered. The comparison of the DGF results against the traditional application of the Nusselt number conservation shows higher accuracy on convective heat transfer extrapolation to different thermal boundary conditions. Additionally, the use of DGF allows the description of the convective heat transfer in processes or configurations where the characterization of the adiabatic wall temperature may be unavailable due to geometrical or performance constraints.

Finally, experiments on an uncooled flat plate were completed to characterize the DGF matrix. A burst-mode laser was used to individually excite each of the elements, while the temperature was monitored with wall-mounted thermocouples. The sensitivity matrix was characterized thanks to the laser heat addition and the wall temperature measurements. The DGF heat flux prediction was compared with the measured heat flux, displaying excellent agreement and proving the applicability and accuracy of the DGF method on forced convection experimental analysis.

This work was supported in part by an appointment to the U.S. Department of Energy (DOE) Postgraduate Research Program at the National Energy Technology Laboratory (NETL). This appointment was administered by the Oak Ridge Institute for Science and Education (ORISE).We would also like to thank Paul Clark, Lakshya Bhatnagar, and Francisco Lozano for their help setting up the instrumentation and the test article for the experimental campaign.

There are no conflicts of interest.

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request. The authors attest that all data for this study are included in the paper.

Neither the United States Government nor any agency thereof, nor any of their employees, nor any of their employees, makes any warranty, expressed or implied, or assumes any legal liability or responsibility for the accuracy, completeness, or usefulness of any information, apparatus, product, or process disclosed, or represents that its use would not infringe privately owned rights. Reference herein to any specific commercial product, process, or service by trade name, trademark, manufacturer, or otherwise, does not necessarily constitute or imply its endorsement, recommendation, or favoring by the United States Government or any agency thereof. The views and opinions of authors expressed herein do not necessarily state or reflect those of the United States Government or any agency thereof.

d =

thickness (m)

t =

time (s)

u =

axial flow velocity (m/s)

x =

axial distance (m)

y =

normal distance (m)

z =

spanwise distance (m)

A =

area (m²)

C =

heat capacity (J/K)

H =

heat transfer coefficient (W/(m² K))

I =

coolant-to-flow momentum flux ratio

L =

material width (m)

P =

pressure (Pa)

T =

temperature (K)

$q˙$ =

heat flux (W/m²)

c_p =

specific heat (J/(kg K))

c_k =

Cerakote

$gij−1$ =

inverse discrete Green’s function

FWHM =

full width at half maximum (mm)

Δ =

gradient

Δs =

spacing between elements (m)

η =

film-cooling effectiveness ((T_aw − T_ref)/(T_c − T_ref))

κ =

thermal conductivity (W/(mK))

ρ =

density (kg/m³)

Φ =

orifice diameter (m)

0 =

stagnant

∞ =

freestream

aw =

adiabatic wall

c =

coolant

ref =

reference

rec =

recovery

wall =

wall conditions

In order to quantify the uncertainty of the measurements, all primary variables used in their derivation factors are taken into account. The main sources of uncertainty in the methodology are summarized in Table 3. All error evaluations are given at 95% confidence level. The absolute uncertainty of each derived quantity is estimated based on the variable mean value and impact of the uncertainty of each one of the prime factors used to derive that quantity. For example, Table 3 represents the Venturi massflow uncertainty the massflow through the critical Venturi depends only on three prime variables: total pressure, total temperature, and discharge coefficient (obtained from Venturi calibration). Based on the mean value for each one of these quantities, the mean massflow is obtained. Then, additional evaluations of the massflow computation are repeated including the uncertainty of each one of the prime variables. Considering the mean value of the total pressure plus its uncertainty the massflow is reevaluated keeping the other quantities at its mean value. The outcome of that calculation reflects the massflow uncertainty associated to the total pressure uncertainty reading. The global uncertainty is then computed as the square root of the sum of squares of the individual uncertainties. Table 4 summarizes the uncertainty estimation on the massflow measurements. In this line, Tables 5 and 6 define the uncertainty on the heat flux measurements and the heat transfer coefficient calculation, respectively. Finally, Table 7 identifies the uncertainty associated with the discrete Green’s function heat transfer characterization.