Next Article in Journal
Coherent Streamwise Vortex Structure of a Three-Dimensional Wall Jet
Next Article in Special Issue
Assessing the Applicability of the Structure-Based Turbulence Resolution Approach to Nuclear Safety-Related Issues
Previous Article in Journal / Special Issue
First-Order Comprehensive Adjoint Sensitivity Analysis Methodology for Critical Points in Coupled Nonlinear Systems. I: Mathematical Framework
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

First-Order Comprehensive Adjoint Sensitivity Analysis Methodology for Critical Points in Coupled Nonlinear Systems. II: Application to a Nuclear Reactor Thermal-Hydraulics Safety Benchmark

by
Dan Gabriel Cacuci
Department of Mechanical Engineering, University of South Carolina, Columbia, SC 29208, USA
Fluids 2021, 6(1), 34; https://doi.org/10.3390/fluids6010034
Submission received: 9 December 2020 / Revised: 4 January 2021 / Accepted: 8 January 2021 / Published: 10 January 2021

Abstract

:
Responses defined at critical points are particularly important for reactor safety analyses and licensing (e.g., the maximum fuel and/or clad temperature). The novel mathematical framework of the first-order comprehensive adjoint sensitivity analysis methodology for critical points (1st-CASAM-CP) is applied in this work to develop a reactor safety thermal-hydraulics benchmark model which admits exact closed-form expressions for the adjoint functions and for the first-order sensitivities of responses defined at critical points (maxima, minima, saddle points) in physical systems characterized by imprecisely known parameters, external and internal boundaries. This benchmark model is designed for verifying the capabilities and accuracies of computational tools for modeling numerically thermal-hydraulics systems. The unique and extensive capabilities of the 1st-CASAM-CP methodology are demonstrated in this work by considering two responses of paramount importance in reactor safety, namely, (i) the maximum rod surface temperature, which occurs at the imprecisely known interface between the subsystem that models the heat conduction inside the heated rod and the subsystem modeling the heat convection process surrounding the rod; and (ii) the maximum temperature inside the heated rod, which has a critical point with two components, one located at a precisely known boundary of the subsystem that models the heat conduction inside the heated rod, while the other component depends on an imprecisely known boundary (i.e., the rod length). The exact analytical expressions developed in this work for the sensitivities of the maximum internal rod temperature and maximum rod surface temperature, as well as for the sensitivities of the locations where these respective maxima occur, provide exact benchmarks for verifying the accuracy of thermal-hydraulics computational tools. The sensitivities of such responses and of their critical points with respect to model parameters enable the quantification of uncertainties induced by uncertainties stemming from the system’s parameters and boundaries in the respective responses and their underlying critical points.

1. Introduction

The transfer of the heat generated in reactor fuel rods to the reactor’s coolant and the subsequent transport of this heat by the coolant to the primary heat exchanger are processes generally modeled using three-dimensional thermal-hydraulics computational models. Quantities of particular interest in nuclear reactor safety are the peak (maximum) temperature within the “hottest” fuel rod and the peak cladding temperature (i.e., the maximum temperature on the rod’s surface) of the “hottest” fuel rod within the reactor’s “hottest channel”. The maximum admissible temperature within the rod must remain below a regulatory set limit, in order to avoid structural damage (e.g., melting). The maximum admissible temperature at the surface of the rod must also remain below a regulatory set limit, in order to avoid the onset of chemical and structural interactions between the rod and the surrounding coolant. In general, meeting the safety-imposed criterion for the maximum rod surface temperature is the more important of these two safety criteria since the maximum temperature within the rod usually remains below the limiting safety criterion when the maximum rod surface temperature remains below its safety-imposed limiting temperature.
The thermal-hydraulics computational models comprise many imperfectly known parameters; the uncertainties associated with such parameters induce uncertainties in the computed results. To verify, within known uncertainty bands, the results produced by computational models, it is important to develop accurate benchmarks that admit exact solutions. The “solution verification” process should include the verification of the accuracy of the sensitivities of the results computed by such codes to the uncertain parameters underlying the respective codes. Previous work [1,2,3] has presented a heat transport benchmark model that simulates the steady-state radial conduction in a fuel rod coupled to axial heat convection in a coolant surrounding the rod and flowing along it, as it occurs within a channel in a nuclear reactor and/or within a heated test section of an experimental facility. This benchmark model admits exact analytical solutions for the spatially dependent temperature distributions within the rod and the surrounding coolant, as well as for the adjoint functions needed to obtain exactly sensitivities of the temperature distribution in the coupled rod/coolant system to this system’s uncertain parameters, interfaces, and external boundaries. The development of this benchmark model was motivated by the need to verify the numerical results produced by the commercially developed “computational fluid dynamics” (CFD) software “FLUENT Adjoint Solver” [4], which was used for computing thermal-hydraulics processes within the G4M Reactor [5], an innovative small modular fast reactor cooled by lead-bismuth eutectic.
In previous works [1,2,3], it has been shown that a direct “solution verification” of the “FLUENT Adjoint Solver” is currently not possible because the current “FLUENT Adjoint Solver” does not provide user-access to the adjoint functions it computes. Therefore, the results produced by the “FLUENT Adjoint Solver” could be verified only indirectly, by comparing the sensitivities (which are 1st-order only) computed with the “FLUENT Adjoint Solver” with the exact results obtained from the analytical expression of the corresponding benchmark sensitivities. Furthermore, the current “FLUENT Adjoint Solver” cannot compute sensitivities of the temperature distribution within the solid rod, and most of the important sensitivities of the coolant temperature (including sensitivities to the boundary heat transfer coefficient and sensitivities to material properties such as thermal conductivity, specific heat) were not obtainable from the current post processing output provided by the “FLUENT Adjoint Solver”.
The developments to be presented in this work are motivated, on the one hand, by need to illustrate the application of the “first-order comprehensive sensitivity analysis of critical points” (1st-CASAM-CP) presented in [6] by using a model that admits closed-form exact expressions for the adjoint functions needed for computing efficiently and exactly first-order sensitivities of responses defined at critical points (maxima, minima, saddle points) to uncertain parameters, interfaces, and boundaries. On the other hand, the analytical results presented in this work also provide unique benchmark solutions for the verification of the accuracy of adjoint state functions and sensitivities that would be produced by future developments of the “FLUENT Adjoint Solver” and similar future software that would certainly be used in nuclear reactor safety work. By applying the general methodology developed in [6], the present work provides the added capability of performing sensitivity analysis of the locations (in the phase-space of independent variables) of critical points (maxima, minima, etc.) of such responses.
This work is structured as follows: Section 2 presents the mathematical model of a heated rod surrounded by coolant, which simulates flow in a reactor channel or in an experimental thermal-hydraulics (TH) experimental facility. Also presented in Section 2 are the mathematical definitions of various responses of fundamental importance for reactor design and safety (e.g., maximum temperature inside the heated rod, maximum rod-surface temperature) located at critical points in the phase-space of the independent variables underlying two coupled generic nonlinear physical systems comprising imprecisely known parameters, interfaces, and boundaries. Section 3 illustrates the application of the 1st-CASAM-CP [6] mathematical framework to the obtain the expressions of the first-order sensitivities of the maximum rod surface temperature and of its critical point with respect to the physical systems’ imprecisely known parameters, interfaces, and boundaries. The critical point of the maximum rod surface temperature is located at the imprecisely known interface between the subsystem that models the heat conduction inside the heated rod and the subsystem modeling the heat convection process surrounding the rod. Hence, both components of this critical point are subject to uncertainties. Section 4 presents the application of the 1st-CASAM-CP [6] to obtain the expressions of the first-order sensitivities of the maximum temperature rod surface and of its critical point with respect to the physical systems’ imprecisely known parameters, interfaces, and boundaries on an interface for computing exactly and efficiently the magnitude of the response and of the phase-space location of its critical point. The critical point of the maximum rod surface temperature is located at a precisely known boundary of the subsystem that models the heat conduction inside the heated rod. Therefore, only one of the components of this critical point is subject to uncertainties. The discussion in Section 5 highlights the significance and possible future applications of the exact mathematical expressions derived for the TH-benchmark presented in this work.

2. Mathematical Model of a Heated Rod Surrounded by Coolant

The benchmark model considered in this work simulates the steady-state heat transfer processes in an idealized reactor channel or in the test section of thermal-hydraulics experimental facilities. The heat transfer process in this benchmark model comprises two coupled “subsystems” which are defined as follows:
  • “Subsystem I” models the steady-state heat conduction in a cylindrical rod of radius a and length (height) , with a , so that the heat conduction in the axial direction can be neglected by comparison to the heat conduction in the radial direction. The rod is heated by an internal volumetric source of the form q cos ( π z / ) , which simulates the axial power distribution in a nuclear reactor; q [ W m 3 ] denotes a constant volumetric source, while z denotes the coordinate along the rod’s axial (vertical) direction. The rod’s conductivity, k     [ W m 1 K 1 ] , is considered to be a temperature-independent constant. Thus, temperature distribution within the rod, T ( r , z ) , is governed by the following heat conduction equation:
    k r r [ r T ( r , z ) r ] = q cos π z   , 0 r < a ,       2 z 2 ,  
    The rod’s surface is cooled by forced convection to a surrounding liquid flowing along the rod’s length, from the rod’s lower end, taken to be located at z = / 2 , towards the rod’s upper end, located at z = / 2 .
  • “Subsystem II” models the distribution of the temperature, denoted as T f l ( z ) , in the coolant, using the following energy conservation equation:
    d T f l ( z ) d z = π a 2 q W c p cos π z   ,   2 z 2 ,
  • The interface (coupling) relation between the temperature distribution in the rod and the temperature distribution in the coolant is provided by the relation
    k T ( r , z ) r = h [ T ( r , z ) T f l ( z ) ] ,             at       r = a ,
    where the heat transfer coefficient, h     [ W m 2 K 1 ] , from the rod’s surface to the coolant is considered to be an imprecisely known constant.
  • The boundary conditions for T ( r , z ) and T f l ( z ) are as follows:
    T ( r , z ) r = 0 ,         at       r = 0 ,
    T f l ( z ) =   T i n l e t , at         z = / 2
    where T i n l e t [ K ] denotes the inlet temperature.
The system of Equation (1) through Equation (5) can be solved exactly to obtain the following closed form expressions for the temperature distributions within the rod and coolant, respectively:
T ( r , z ) = q ( a 2 r 2 4 k + a 2 h ) cos π z   + T f l ( z ) , 0 r < a ,       2 z 2 ,  
T f l ( z ) = a 2 q W c p ( sin π z + 1 )   + T i n l e t , 2 z 2 .
The imprecisely known parameters underlying the paradigm heat transfer benchmark modeled by Equation (1) through Equation (5) are as follows:
(a)
the model parameters q ,   k ,   h ,   W ,   c p ,   T i n l e t ;
(b)
the interface location   a between the solid rod and the fluid coolant;
(c)
the external boundaries defined by the imprecisely known length, , of the heated rod.
The nominal values of these parameters are considered to be known and will be denoted by using the superscript “zero,” i.e., q 0 ,     k 0 ,     h 0 ,     W 0 ,     c p 0 ,     T i n l e t 0 ,       a 0 ,     0 . Variations in the model parameters, interface, and boundaries around the respective nominal values will be denoted as follows: δ q ,     δ k ,     δ h ,       δ W ,     δ c p ,     δ T i n l e t ,       δ a ,       δ . Such variations will induce variations in the rod and coolant temperatures.
The expression of the rod surface temperature is obtained by setting r = a in Equation (6) to obtain
T ( r , z ) = a q 2 h cos π z   + T f l ( z ) ,       2 z 2
Evidently, the rod surface temperature is defined on the imprecisely known interface r = a between the subsystem that models the heat conduction in the rod and the subsystem that models the heat convection in the surrounding coolant. Therefore, the maximum temperature on the rod’s surface, which will be denoted as T s max , will also occur on the interface r = a between the aforementioned subsystems, at a location denoted as ( r s , z s ) , which is defined by the following conditions:
r s = a   ,       { T ( r , z ) z } ( r s , z s ) = 0
As Equation (9) indicates, both components of the critical point ( r s , z s ) are affected by uncertainties in the imprecisely known model parameters, interface, and boundaries. Applying the conditions provided in Equation (9) to the expression provided in Equation (6) yields the following expression for the location z s :
z s = π arctan ( 2 a h W c p ) .
The closed-form expression of the maximum rod surface temperature will be denoted as T s max ; it is obtained by setting r = a and z = z s in Equation (6) to obtain
T s max = q a 2 h cos π z s   + T f l ( z s )
The maximum temperature within the rod will be denoted as T max ; it occurs at the point ( r max , z max ) where the partial derivatives of T ( r , z ) with respect to the independent variables vanish, i.e.,
T ( r , z ) r = 0   ,       T ( r , z ) z = 0 ,       a t       r = r max ,     z = z max
Applying the conditions provided in Equation (12) to the expression of T ( r , z ) yields the following expressions for the coordinates ( r max , z max ) of the maximum rod temperature:
r max = 0 ,         z max = π arctan 4 a k h W c p ( a h + 2 k )
As indicated by Equation (13), only one component, namely z max , of the critical point ( r max , z max ) is affected by uncertainties in the imprecisely known model parameters, interface, and boundaries. The other component of the critical point ( r max , z max ) , namely r max = 0 , is not affected by any uncertainties since it is located at the exactly known boundary point r = 0 .
Using the results obtained in Equation (13) in Equation (6) yields the following closed-form expression for T max :
T max = q ( a 2 4 k + a 2 h ) cos π z max   + T f l ( z max )
It is evident from the closed-form expressions obtained in Equation (6) through Equation (11) that both T s max and T max , as well as the components of the respective critical points, depend on the uncertain parameters. The closed-form expressions obtained in Equation (6) through Equation (11) will be used in the remainder of this work for verifying and validating the expressions that will be obtained in Section 3 and Section 4, below, for the sensitivities of T max , T s max and of their respective critical points (maxima) with respect to the model’s uncertain parameters, interface, and boundaries.

3. Maximum Rod Surface Temperature: Critical Point Located on Interface

Section 3.1 presents the derivation of the exact, closed-form expressions of the first-order sensitivities of the maximum rod surface temperature with respect to the uncertain model parameters, interface, and boundaries. Section 3.2 presents the derivation of the exact, closed form expressions of the first-order sensitivities of the location (in the space of independent variables) of the critical point of the maximum rod surface temperature with respect to the uncertain model parameters, interface, and boundaries.

3.1. First-Order Sensitivities of the Maximum Rod Surface Temperature

The maximum temperature of the rod’s surface, T s max , can be represented in the following form:
T s max T ( a , z s ) = 0 a r d r / 2 / 2 d z   T ( r , z ) δ ( r a ) r δ ( z z s )
where z s is implicitly defined by the relation in Equation (9). The first-order total differential, δ T s max , of T s max depends, in principle, on all of the first-order sensitivities of T s max with respect to the imprecisely known model and boundary parameters through the following relation:
δ T s max = T s max q   δ q + T s max h δ k + T s max h   δ h + T s max W   δ W           + T s max c p     δ c p + T s max   T i n l e t δ T i n l e t + T s max a δ a + T s max   δ .
The total differential δ T s max is obtained by applying the definition of the G-differential to Equation (15), which yields:
δ T s max =   { d d ε 0 a 0 + ε δ a r d r ( 0 + ε δ ) / 2 ( 0 + ε δ ) / 2   [ T ( r , z ) + ε δ T ( r , z ) ] δ ( r a 0 ε δ a ) r   × δ ( z z s 0 ε δ z s )   d z     } ε = 0 = { δ T s max } d i r + { δ T s max } i n d .
The indirect-effect term, { δ T s max } i n d , in Equation (17) depends only on the variation in the respective state function, namely δ T ( r , z ) , and is defined as follows:
{ δ T s max } i n d 0 a 0 r d r 0 / 2 0 / 2   δ T ( r , z ) δ ( r a 0 ) r δ ( z z s 0 )   d z  
Thus, the indirect-effect term depends on the parameter variations indirectly, through the variation δ T ( r , z ) in the rod temperature. In contradistinction, the direct-effect term { δ T s max } d i r in Equation (17) depends directly on parameter variations and is defined as follows:
{ δ T s max } d i r = ( δ ) { T ( r , z ) } ( r = a 0 , z = z s 0 ) + ( δ z s )   { T ( r , z ) z s } ( r = a 0 , z = z s 0 )                                                       + ( δ a ) { T ( r , z ) r } ( r = a 0 , z = z s 0 ) .
Thus, the direct-effect term { δ T s max } d i r e c t can already be computed at this stage by using Equation (6) to obtain:
{ δ T s max } d i r e c t = ( δ a ) ( a 0 q 0 2 k 0 cos π z s 0 0 )
The variation δ T ( r , z ) is the solution of the “first-level forward sensitivity system” (1st-LFSS) which is obtained by G-differentiating the original system defined by Equation (1) through Equation (5). Applying the definition of the G-differential to Equation (1) through Equation (5) yields the following relations:
{ d d ε k 0 + ε δ k r r [ r ( T 0 + ε δ T ) r ] } ε = 0 = { d d ε ( q 0 + ε δ q ) cos π z 0 + ε δ } ε = 0   ,
{ d d ε ( T 0 + ε δ T ) r } ε = 0 = 0 ,         at       r = 0 ,
{ d d ε ( k 0 + ε δ k ) [ T 0 ( a 0 + ε δ a , z ) + ε δ T ( a 0 + ε δ a , z ) ] ( a 0 + ε δ a ) } ε = 0 = { d d ε ( h 0 + ε δ h ) [ T 0 ( a 0 + ε δ a , z ) + ε δ T ( a 0 + ε δ a ) T f l 0 ( z ) ε δ T f l ( z ) ] } ε = 0 ,       at       r = a 0 ,
{ d d ε d [ T f l 0 ( z ) + ε δ T f l ( z ) ] d z } ε = 0 = π { d d ε ( a 0 + ε δ a ) 2 ( q 0 + ε δ q ) ( W 0 + ε δ W ) ( c p 0 + ε δ c p ) cos π z 0 + ε δ } ε = 0 ,
{ d d ε [ T f l 0 ( 0 + ε δ 2 ) + ε δ T f l ( 0 + ε δ 2 ) ] } ε = 0 =   { d d ε ( T i m l e t 0 + ε δ T i m l e t ) } ε = 0 .
Carrying out in Equation (21) through Equation (25) the differentiations with respect to ε and setting ε = 0 in the resulting expressions yields the following set of equations, which constitute the 1st-LFSS:
k 0 r r { r r [ δ T ( r , z ) ] } = ( δ k ) r r [ r T 0 ( r , z ) r ] ( δ q ) cos π z 0 + ( δ ) q 0 π z ( 0 ) 2 sin π z 0 ,  
r [ δ T ( r , z ) ] = 0 ,             at       r = 0 ,
( δ k ) { T 0 ( r , z ) r } r = a 0 k 0 ( δ a ) { 2 T 0 ( r , z ) r 2 } r = a 0 k 0 { r [ δ T ( r , z ) ] } r = a 0 = ( δ h ) { [ T 0 ( r , z ) T f l 0 ( z ) ] } r = a 0 + h 0 [ δ T ( r ) δ T f l ( z ) ] + ( δ a ) h 0 { T 0 ( r , z ) r } r = a 0   ,
d d z [ δ T f l ( z ) ] = π [ 2 a 0 q 0 ( δ a ) W 0 c p 0 + ( a 0 ) 2 ( δ q ) W 0 c p 0 ( a 0 ) 2 q 0 ( δ W ) ( W 0 ) 2 c p 0 ( a 0 ) 2 q 0 ( δ c p ) W 0 ( c p 0 ) 2 ] cos π z 0                                                                     + ( δ ) π 2 z ( a 0 ) 2 q 0 W 0 c p 0 ( 0 ) 2 sin π z 0 Q f l ( z )   ,         0 2 z 0 2 ,  
δ T f l ( z ) =   d T f l 0 ( z ) d z δ 2 + δ T i n l e t = δ T i n l e t , at         z = 0 / 2 .  
The first term on the right-side of Equation (26) can be simplified by using Equation (1) to obtain the following equation:
k 0 r r { r r [ δ T ( r , z ) ] } = [ ( δ k ) q 0 k 0 ( δ q ) ] cos π z 0 + ( δ ) q 0 π z ( 0 ) 2 sin π z 0 Q ( z ) .
The terms containing derivatives of T ( r , z ) in Equation (28) can also be simplified using Equations (1) and (3) to obtain the following equation:
{ k 0 r [ δ T ( r , z ) ] h 0 [ δ T ( r , z ) δ T f l ( z ) ] } r = a 0 =   ( δ h ) a 0 q 0 2 h 0 cos π z 0             ( δ k ) a 0 q 0 2 k 0 cos π z 0   ( δ a ) q 0 2 ( h 0 a 0 k 0 + 1 )   cos π z 0   .
Since the equations underlying the 1st-LFSS, cf. Equations (27) and (29) through Equation (32), depend on the parameter variations, it is computationally expensive to repeatedly solve the 1st-LFSS for all possible parameter variations. The need for repeatedly solving the 1st-LFSS can be circumvented by expressing the indirect-effect term defined in Equation (18) in terms of the solution of a “first-level adjoint sensitivity system” (1st-LASS), which will be constructed next by applying the general principles of the 1st-CASAM-CP presented in [6].
The Hilbert space appropriate for the heat transport benchmark under consideration comprises the space of all square-integrable two-component vector functions of the form u ( x )   [ u 1 ( r , z ) , u 2 ( z ) ] , endowed with an inner product u ( x ) , ψ ( x ) of the form
u ( x ) , ψ ( x ) 0 a 0 r d r 0 / 2 0 / 2 d z   [ u 1 ( r , z ) ψ 1 ( r , z ) + u 2 ( z ) ψ 2 ( z ) ]   .  
Using the definition provided in Equation (33), construct the inner product of a square integrable vector function ψ ( x ) = [ ψ ( r , z ) ,   ψ f l ( z ) ] , where ψ ( r , z ) and ψ f l ( z ) denote the adjoint sensitivity functions that correspond to the forward functions δ T ( r , z ) and δ T f l ( z ) , with Equations (29) and (31), respectively, to obtain the following relation:
0 a 0 r d r 0 / 2 0 / 2 d z { ψ ( r , z ) k 0 r r [ r r δ T ( r , z ) ] + ψ f l ( z ) d [ δ T f l ( z ) ] d z }         = 0 a 0 r d r 0 / 2 0 / 2 d z [ ψ ( r , z ) Q ( z ) + ψ f l ( z ) Q f l ( z ) ]     .
The left-side of Equation (34) is now integrated by parts (twice over the variable r and once over the variable z ) to obtain
0 a 0 r d r 0 / 2 0 / 2 d z { ψ ( r , z ) k 0 r r [ r r δ T ( r , z ) ] + ψ f l ( z ) d [ δ T f l ( z ) ] d z }   = 0 a 0 r d r 0 / 2 0 / 2 d z [ δ T f l ( z ) ] [ d ψ f l ( z ) d z ] + 0 a 0 r d r 0 / 2 0 / 2 d z { [ δ T ( r , z ) ] k 0 r r [ r ψ ( r , z ) r ] } + 0 a 0 r d r [ ψ f l ( 0 2 ) δ T f l ( 0 2 ) ψ f l ( 0 2 ) δ T f l ( 0 2 ) ] + 0 / 2 0 / 2 d z [ ψ ( r , z ) r k 0 r δ T ( r , z ) δ T ( r , z ) r k 0 ψ ( r , z ) r ] r = a 0 0 / 2 0 / 2 d z [ ψ ( r , z ) r k 0 r δ T ( r , z ) δ T ( r , z ) r k 0 ψ ( r , z ) r ] r = 0 .
Using the boundary condition given in Equation (27) and imposing the boundary condition
r ψ ( r , z ) r = 0 ,       at       r = 0 ,
eliminates the last term on the right-side of Equation (35). Imposing the boundary condition
ψ f l ( z ) =   0 , at         z = 0 / 2
eliminates the unknown function δ T f l ( z = 0 / 2 ) on the right-side of Equation (35). Using the boundary condition given in Equation (30) to replace the term δ T f l ( z = / 2 ) on the right side of Equation (35) and replacing the left-side of Equation (35) by the right-side of Equation (34) yields the following expression equivalent to Equation (35):
0 a 0 r d r 0 / 2 0 / 2 d z [ ψ ( r , z ) Q ( z ) + ψ f l ( z ) Q f l ( z ) ]   = 0 a 0 r d r 0 / 2 0 / 2 d z { [ δ T ( r , z ) ] k 0 r r [ r ψ ( r , z ) r ] + [ δ T f l ( z ) ] [ d ψ f l ( z ) d z ] } 0 a 0 r d r     ψ f l ( 2 ) δ T i n l e t + 0 / 2 0 / 2 d z [ ψ ( r , z ) r k 0 r δ T ( r , z ) δ T ( r , z ) r k 0 ψ ( r , z ) r ] r = a 0 .
The unknown quantity { [ δ T ( r , z ) ] / r } r = a 0 , which appears in the last term on the right-side of Equation (38) is eliminated by using the boundary condition given in Equation (32); this operation transforms Equation (38) into the following form:
0 a 0 r d r 0 / 2 0 / 2 d z [ ψ ( r , z ) Q ( z ) + ψ f l ( z ) Q f l ( z ) ]   = 0 a 0 r d r 0 / 2 0 / 2 d z { [ δ T ( r , z ) ] k 0 r r [ r ψ   ( r , z ) r ] + [ δ T f l ( z ) ] [ d ψ f l ( z ) d z ] } 0 a 0 r d r   ψ f l ( 2 ) ( δ T i n l e t ) 0 / 2 0 / 2 d z [ δ T ( r , z ) r k 0 ψ ( r , z ) r ] r = a 0 0 / 2 0 / 2 d z { ψ ( a 0 , z ) h 0 a 0 [ δ T ( r , z ) δ T f l ( z ) ] } r = a 0 0 / 2 0 / 2 d z     ψ ( a 0 , z ) a 0 { ( δ h ) a 0 q 0 2 h 0 cos π z 0   ( δ k ) a 0 q 0 2 k 0 cos π z 0   ( δ a ) q 0 2 ( h 0 a 0 k 0 + 1 )   cos π z 0 } .
The unknown quantity δ T ( a 0 , z ) , which appears in third and fourth terms on the right-side of Equation (39), is eliminated by imposing the following interface condition on the (adjoint) function Ψ ( r , z ) :
k 0 ψ ( r , z ) r = h 0   ψ ( r , z ) ,             at       r = a 0
Inserting Equation (40) into the right-side of Equation (39) reduces it to the following form:
0 a 0 r d r 0 / 2 0 / 2 d z [ ψ ( r , z ) Q ( z ) + ψ f l ( z ) Q f l ( z ) ]   = 0 a 0 r d r 0 / 2 0 / 2 d z { [ δ T ( r , z ) ] k 0 r r [ r ψ ( r , z ) r ] + [ δ T f l ( z ) ] [ d ψ f l ( z ) d z ] } ( a 0 ) 2 2 ψ f l ( 2 ) ( δ T i n l e t ) + h 0 a 0 0 / 2 0 / 2 d z   ψ ( a 0 , z ) [ δ T f l ( z ) ] a 0 [ ( δ h ) a 0 q 0 2 h 0 ( δ k ) a 0 q 0 2 k 0 ( δ a ) q 0 2 ( h 0 a 0 k 0 + 1 ) ] 0 / 2 0 / 2 d z   ψ ( a 0 , z ) cos π z 0 .
The two terms that contain the unknown function δ T f l ( z ) in Equation (41) are grouped together, transforming Equation (41) into the following form
0 a 0 r d r 0 / 2 0 / 2 d z [ ψ ( r , z ) Q ( z ) + ψ f l ( z ) Q f l ( z ) ]   = ( a 0 ) 2 2 ψ f l ( 2 ) ( δ T i n l e t ) + 0 a 0 r d r 0 / 2 0 / 2 d z { [ δ T ( r , z ) ] k 0 r r [ r ψ ( r , z ) r ] + [ δ T f l ( z ) ] [ d ψ f l ( z ) d z + 2 h 0 a 0 ψ ( a 0 , z ) ] } a 0 [ ( δ h ) a 0 q 0 2 h 0 ( δ k ) a 0 q 0 2 k 0 ( δ a ) q 0 2 ( h 0 a 0 k 0 + 1 ) ] 0 / 2 0 / 2 d z   ψ ( a 0 , z ) cos π z 0 .
The second term on the right-side of Equation (42) will represent the indirect-effect term { δ T s max } i n d defined in Equation (18) by requiring that the following equations be satisfied:
k 0 r [ r ψ ( r , z ) r ] = δ ( r a 0 ) δ ( z z s 0 )   , 0 r < a 0 ,       0 2 z 0 2 ,    
ψ f l ( z ) z + 2 h 0 a 0 ψ ( a 0 , z )   = 0 ,   0 2 z 0 2 ,  
Inserting the relations provided in Equations (43) and (44) into Equation (18) and re-arranging the resulting equation yields the following expression for the indirect-effect term { δ T s max } i n d :
{ δ T s max } i n d = 0 a 0 r d r 0 / 2 0 / 2 d z [ ψ ( r , z ) Q ( z ) + ψ f l ( z ) Q f l ( z ) ]   + a 0 [ ( δ h ) a 0 q 0 2 h 0 ( δ k ) a 0 q 0 2 k 0 ( δ a ) q 0 2 ( h 0 a 0 k 0 + 1 ) ] 0 / 2 0 / 2 d z   ψ ( a 0 , z ) cos π z 0         + ( a 0 ) 2 2 ψ f l ( 2 ) ( δ T i n l e t ) .
Inserting the definitions provided for Q ( z ) and Q f l ( z ) in Equations (29) and (31), respectively, into Equation (45) yields the following expression:
{ δ T s max } i n d = [ ( δ k ) q 0 k 0 ( δ q ) ] 0 a 0 r d r   0 / 2 0 / 2 ψ ( r , z ) cos π z 0 d z + + ( δ ) q 0 π ( 0 ) 2 0 a 0 r d r   0 / 2 0 / 2 ψ ( r , z ) z sin π z 0 d z + ( δ ) π 2 ( a 0 ) 2 q 0 W 0 c p 0 ( 0 ) 2 ( a 0 ) 2 2 0 / 2 0 / 2 ψ f l ( z )     z sin π z 0 d z + ( a 0 ) 2 2 π [ 2 a 0 q 0 ( δ a ) W 0 c p 0 + ( a 0 ) 2 ( δ q ) W 0 c p 0 ( a 0 ) 2 q 0 ( δ W ) ( W 0 ) 2 c p 0   ( a 0 ) 2 q 0 ( δ c p ) W 0 ( c p 0 ) 2 ] 0 / 2 0 / 2 ψ f l ( z ) cos π z 0 d z   + ( a 0 ) 2 2 ψ f l ( 2 ) ( δ T i n l e t ) + a 0 [ ( δ h ) a 0 q 0 2 h 0 ( δ k ) a 0 q 0 2 k 0 ( δ a ) q 0 2 ( h 0 a 0 k 0 + 1 ) ] 0 / 2 0 / 2   ψ ( a 0 , z ) cos π z 0 d z   .
Inserting the results obtained in Equations (20) and (46) into Equation (17) and identifying the expressions that multiply the various parameter variations by comparing Equation (46) to Equation (16) yields the following results for the partial sensitivities of T s max with respect to the uncertain parameters, interface, and boundaries:
T s max q = 0 a 0 r d r   0 / 2 0 / 2 ψ ( r , z ) cos π z 0 d z + ( a 0 ) 4 W 0 c p 0 π 2 0 / 2 0 / 2 ψ f l ( z ) cos π z 0 d z
T s max k = q 0 k 0 0 a 0 r d r   0 / 2 0 / 2 ψ ( r , z ) cos π z 0 d z ( a 0 ) 2 q 0 2 k 0 0 / 2 0 / 2   ψ ( a 0 , z ) cos π z 0 d z
T s max h = ( a 0 ) 2 q 0 2 h 0 0 / 2 0 / 2   ψ ( a 0 , z ) cos π z 0 d z
T s max W = π q 0 ( a 0 ) 4 2 ( W 0 ) 2 c p 0 0 / 2 0 / 2 ψ f l ( z ) cos π z 0 d z  
T s max c p = π ( a 0 ) 4 q 0 2 W 0 ( c p 0 ) 2 ] 0 / 2 0 / 2 ψ f l ( z ) cos π z 0 d z
T s max T i n l e t = ( a 0 ) 2 2 ψ f l ( 2 )
T s max = q 0 π ( 0 ) 2 0 a 0 r d r   0 / 2 0 / 2 ψ ( r , z ) z sin π z 0 d z   + π 2 ( a 0 ) 4 q 0 2 W 0 c p 0 ( 0 ) 2 0 / 2 0 / 2 ψ f l ( z ) z     sin π z 0 d z ,
T s max a = a 0 q 0 2 ( h 0 a 0 k 0 + 1 ) 0 / 2 0 / 2   ψ ( a 0 , z ) cos π z 0 d z                                   + π ( a 0 ) 3 q 0 W 0 c p 0 0 / 2 0 / 2 ψ f l ( z ) cos π z 0 d z a 0 q 0 2 k 0 cos π z s 0 0 .
The sensitivities obtained in Equation (47) through Equation (54) can be computed efficiently, using just quadrature formulas, after the adjoint sensitivity functions ψ ( r , z ) and ψ f l ( z ) are obtained by solving the “first-level adjoint sensitivity system” (1st-LASS), which comprises Equations (43) and (44) together with the interface condition provided in Equation (40) and the boundary conditions given in Equations (36) and (37). The 1st-LASS is solved in a manner that is “reverse/backwards” by comparison to the way in which the solution proceeds for solving the 1st-LFSS and/or for the original heat transport model. Thus, while the 1st-LFSS and the original heat transport model are solved by starting with the fluid flow equation (which is solved from the inlet to the outlet of the fluid flow) and subsequently solving the heat conduction equation in the rod, the solution of the 1st-LASS proceeds in the reverse manner, by first solving the heat conduction in the rod, followed by solving the fluid flow equation from the outlet to the inlet. Notably, the 1st-LASS, is independent of parameter variations and needs to be solved just once to obtain the adjoint functions ψ ( r , z ) and ψ f l ( z ) .
Solving the 1st-LASS yields the following expressions for the adjoint functions ψ ( r , z ) and ψ f l ( z ) :
ψ ( r , z ) = δ ( z z s 0 ) [ 1 a 0 h 0 + H ( r a 0 ) 1 k 0 ln r a 0 ] ,     0 < r < a 0 ,       0 2 z , z s 0 0 2 ,
ψ f l ( z ) =   2 ( a 0 ) 2 H ( z s 0 z ) ,         0 2 z , z s 0 0 2 .  
Using Equations (55) and (56) in Equation (47) through Equation (54) and carrying out the respective integrations yields the following expressions for the 1st-order sensitivities of the maximum rod temperature, T s max :
T s max q = a 0 2 h 0 cos π z s 0 0 + ( a 0 ) 2 0 W 0 c p 0 ( sin π z s 0 0 + 1 ) ,
T s max k = 0
T s max h = a 0 q 0 2 ( h 0 ) 2 cos π z s 0 0
T s max W = q 0 0 c p 0 ( a 0 W 0 ) 2 ( sin π z s 0 0 + 1 )
T s max c p = q 0 0 W 0 ( a 0 c p 0 ) 2 ( sin π z s 0 0 + 1 )
T s max T i n l e t = 1
T s max = q 0 a 0 2 h 0 π z ( 0 ) 2 sin π z p 0 + ( a 0 ) 2 q 0 W 0 c p 0 ( π z p 0 cos π z p 0 + sin π z p 0 + 1 )
T s max a = q 0 2 h 0 cos π z s 0 + 2 a 0 q 0 0 W 0 c p 0 ( sin π z p 0 + 1 )

3.2. First-Order Sensitivities of the Critical Point (Maximum) of the Rod Surface Temperature

As has been shown in Equation (9), the maximum value, T s max , of the rod surface temperature T ( a , z ) is attained at the critical point ( r s , z s ) . As Equation (10) indicates, the components of the critical point ( r s , z s ) are subject to uncertainties since they depend on imprecisely known parameters. The sensitivities of the radial component r s = a are evident: the only non-zero sensitivity is r s / a = 1 .
On the other hand, the component z s of the critical point ( r s , z s ) of T s max depends, in principle, on all of the uncertain parameters, which means that the total differential δ z s will have the following expression:
δ z s = z s q   δ q + z s h δ k + z s h   δ h + z s W   δ W   + z s c p     δ c p + z s   T i n l e t δ T i n l e t + z s a δ a + z s   δ .
Since the closed-form expression provided in Equation (10) would not be available in general, the expression of δ z s will be obtained by applying the general methodology presented in [6]. This application commences by writing the relations provided in Equation (9), which implicitly define the location, z s , in the following form:
0 a r d r / 2 / 2 T ( r , z ) z δ ( r a ) r δ ( z z s ) d z   = 0
Taking the G-differential of Equation (66) yields the following relation:
0 =   { d d ε 0 a 0 + ε δ a d r ( 0 + ε δ ) / 2 ( 0 + ε δ ) / 2   [ T ( r , z ) + ε δ T ( r , z ) ] z δ ( r a 0 ε δ a ) δ ( z z s 0 ε δ z s ) d z } ε = 0 = ( δ a ) 0 a 0 d r 0 / 2 0 / 2   [ T ( r , z ) ] z δ ( r a 0 ) δ ( z z s 0 ) d z         ( δ z s ) 0 a 0 d r 0 / 2 0 / 2   [ T ( r , z ) ] z δ ( r a 0 ) δ ( z z s 0 ) d z     + 0 a 0 d r 0 / 2 0 / 2   [ δ T ( r , z ) ] z δ ( r a 0 ) δ ( z z s 0 ) d z   .
The relation in Equation (67) can be re-written in the following form:
( δ z s ) {   2 [ T ( r , z ) ] z 2 } ( r = a 0 , z = z s 0 )                             = 0 a 0 d r 0 / 2 0 / 2   δ T ( r , z ) δ ( r a 0 ) δ ( z z s 0 ) d z   ( δ a ) { 2 [ T ( r , z ) ] z r } ( r = a 0 , z = z s 0 )   ,
where the quantities {   2 [ T ( r , z ) ] / z 2 } ( r = a 0 , z = z s 0 ) and { 2 [ T ( r , z ) ] / z r } ( r = a 0 , z = z s 0 ) are computed using Equation (6) for this benchmark model (or are computed numerically when the closed form of T ( r , z ) is unavailable) to obtain
{   2 [ T ( r , z ) ] z 2 } ( r = a 0 , z = z s 0 ) = { ( π ) 2 q a [ 1 2 h cos π z   + a W c p sin π z ] } ( r = a 0 , z = z s 0 )
{ 2 [ T ( r , z ) ] z r } ( r = a 0 , z = z s 0 )   = a 0 q 0 2 k 0 π 0 sin π z s 0 0
The term 0 a 0 d r 0 / 2 0 / 2   δ T ( r , z ) δ ( r a 0 ) δ ( z z s 0 ) d z in Equation (68) can be expressed in terms of adjoint functions by applying the procedure outlined in [6] to obtain the following result:
0 a 0 d r 0 / 2 0 / 2   δ T ( r , z ) δ ( r a 0 ) δ ( z z s 0 ) d z   = [ ( δ k ) q 0 k 0 ( δ q ) ] 0 a 0 r d r   0 / 2 0 / 2 Φ ( r , z ) cos π z 0 d z + + ( δ ) q 0 π ( 0 ) 2 0 a 0 r d r   0 / 2 0 / 2 Φ ( r , z ) z sin π z 0 d z + ( δ ) π 2 ( a 0 ) 2 q 0 W 0 c p 0 ( 0 ) 2 ( a 0 ) 2 2 0 / 2 0 / 2 Φ f l ( z )     z sin π z 0 d z + ( a 0 ) 2 2 π [ 2 a 0 q 0 ( δ a ) W 0 c p 0 + ( a 0 ) 2 ( δ q ) W 0 c p 0 ( a 0 ) 2 q 0 ( δ W ) ( W 0 ) 2 c p 0   ( a 0 ) 2 q 0 ( δ c p ) W 0 ( c p 0 ) 2 ] 0 / 2 0 / 2 Φ f l ( z ) cos π z 0 d z   + ( a 0 ) 2 2 Φ f l ( 2 ) ( δ T i n l e t ) + a 0 [ ( δ h ) a 0 q 0 2 h 0 ( δ k ) a 0 q 0 2 k 0 ( δ a ) q 0 2 ( h 0 a 0 k 0 + 1 ) ] 0 / 2 0 / 2   Φ ( a 0 , z ) cos π z 0 d z .
where the adjoint functions are the solutions of the following 1st-LASS:
k 0 r [ r Φ ( r , z ) r ] = δ ( r a 0 ) δ ( z z s 0 )   ,
r Φ ( r , z ) r = 0 ,       at       r = 0 ,
k 0 Φ ( r , z ) r = h 0 Φ ( r , z ) ,             at       r = a 0
Φ f l ( z ) z + 2 h 0 a 0 Φ ( a , z )   = 0 ,   0 2 z 0 2 ,  
Φ f l ( z ) =   0 , at         z = 0 / 2
Introducing the results obtained in Equation (69) through Equation (76) into Equation (68) and equating the coefficients corresponding to the same parameter variations yields the following expressions for the sensitivities of the location z s to the imprecisely known model, interface, and boundary parameters:
z s q = { 0 a 0 r d r   0 / 2 0 / 2 Φ ( r , z ) cos π z 0 d z ( a 0 ) 2 2 π ( a 0 ) 2 W 0 c p 0 0 / 2 0 / 2 Φ f l ( z ) cos π z 0 d z ( π ) 2 q a [ 1 2 h cos π z   + a W c p sin π z ] } z = z s
z s k = { q 0 k 0 0 a 0 r d r   0 / 2 0 / 2 Φ ( r , z ) cos π z 0 d z + ( a 0 ) 2 q 0 2 k 0 0 / 2 0 / 2   Φ ( a 0 , z ) cos π z 0 d z ( π ) 2 q a [ 1 2 h cos π z   + a W c p sin π z ] } z = z s
z s h = {   ( a 0 ) 2 q 0 2 h 0 0 / 2 0 / 2   Φ ( a 0 , z ) cos π z 0 d z . ( π ) 2 q a [ 1 2 h cos π z   + a W c p sin π z ] } z = z s
z s W   = { ( a 0 ) 2 2 π ( a 0 ) 2 q 0 ( W 0 ) 2 c p 0   0 / 2 0 / 2 Φ f l ( z ) cos π z 0 d z ( π ) 2 q a [ 1 2 h cos π z   + a W c p sin π z ] } z = z s
z s c p   = { ( a 0 ) 2 2 π ( a 0 ) 2 q 0 W 0 ( c p 0 ) 2   0 / 2 0 / 2 Φ f l ( z ) cos π z 0 d z ( π ) 2 q a [ 1 2 h cos π z   + a W c p sin π z ] } z = z s
z s   = { q 0 π ( 0 ) 2 0 a 0 r d r   0 / 2 0 / 2 Φ ( r , z ) z sin π z 0 d z + π 2 ( a 0 ) 4 q 0 2 W 0 c p 0 ( 0 ) 2 0 / 2 0 / 2 Φ f l ( z )     z sin π z 0 d z ( π 0 ) 2 q 0 a 0 [ 1 2 h 0 cos π z 0   + a 0 0 W 0 c p 0 sin π z 0 ] } z = z s
z s a = { ( π ) 2 q a [ 1 2 h cos π z s 0   + a W c p sin π z s 0 ] } 1 { a 0 q 0 2 k 0 π 0 sin π z s 0 0 + ( a 0 ) 2 2 π 2 a 0 q 0 W 0 c p 0 0 / 2 0 / 2 Φ f l ( z ) cos π z 0 d z a 0 q 0 2 ( h 0 a 0 k 0 + 1 ) 0 / 2 0 / 2   Φ ( a 0 , z ) cos π z 0 d z } z = z s ,
z s T i n l e t = { ( a 0 ) 2 2 Φ f l ( 2 ) ( π ) 2 q a [ 1 2 h cos π z   + a W c p sin π z ] } z = z s
The expressions obtained in Equation (77) through Equation (84) can be evaluated after solving (numerically or analytically) the 1st-LASS represented by Equation (72) through Equation (76) to determine the adjoint functions Φ ( r , z ) and Φ f l ( z ) . Since Φ ( r , z ) and Φ f l ( z ) do not depend on any parameter variations, the 1st-LASS needs to be solved only once to obtain the respective adjoint functions, which are subsequently used in Equation (77) through Equation (84) to compute, using simple quadrature formulas, all of the sensitivities of z s to the imprecisely known model, interface, and boundary parameters. Solving the 1st-LASS analytically yields the following expressions for the adjoint functions Φ ( r , z ) and Φ f l ( z ) :
Φ ( r , z ) = δ ( z z s 0 ) [ 1 a 0 h 0 + 1 k 0 H ( r a 0 ) ln r a 0 ] ,       0 2 z , z s 0 0 2 ,
Φ f l ( z ) = 2 δ ( z z s 0 ) ( a 0 ) 2 ,         0 2 z , z s 0 0 2
Using the expressions for Φ ( r , z ) and Φ f l ( z ) provided in Equations (85) and (86) in Equation (77) through Equation (84) and performing the respective integrations yields the following closed-form expressions for the sensitivities of z s to the imprecisely known model, interface, and boundary parameters:
z s q = 0
z s k = 0
z s h = 2 a 0 ( 0 ) 2 π W 0 c p 0   [ 1   + ( 2 h 0 a 0 0 W 0 c p 0 ) 2 ] 1
z s W = 2 h 0 a 0 ( 0 ) 2 π ( W 0 ) 2 c p 0   [ 1   + ( 2 h 0 a 0 0 W 0 c p 0 ) 2 ] 1
z s c p   = 2 ( 0 ) 2 h 0 a 0 π W 0 ( c p 0 ) 2   [ 1   + ( 2 h 0 a 0 0 W 0 c p 0 ) 2 ] 1
z s   = z s 0 0 + 2 a 0 0 h 0 π W 0 c p 0   [ 1   + ( 2 h 0 a 0 0 W 0 c p 0 ) 2 ] 1 ,
z s a = 2 ( 0 ) 2 h 0 π W 0 c p 0   [ 1   + ( 2 h 0 a 0 0 W 0 c p 0 ) 2 ] 1
z s T i n l e t = 0
For this benchmark model, the closed-form expression for z s is available (which would not be the case for large-scale models) in Equation (10). Consequently, Equation (10) can be used to verify that the expressions obtained using the general 1st-CASAM-CP (which has specifically yielded the 1st-LASS comprising Equation (72) through Equation (76)) has produced exact expressions, in Equation (87) through Equation (94), for the sensitivities of the location z s of the maximum rod surface temperature to all of the system’s imprecisely known parameters.

4. Maximum Rod Temperature: Critical Point Located on Boundary

Section 4.1 presents the derivation of the exact, closed-form expressions of the first-order sensitivities of the maximum temperature within the rod with respect to the uncertain model parameters, interface, and boundaries. Section 4.2 presents the derivation of the exact, closed-form expressions of the first-order sensitivities of the location (in the space of independent variables) of the critical point of the maximum rod temperature with respect to the uncertain model parameters, interface, and boundaries.

4.1. First-Order Sensitivities of the Maximum Temperature Inside the Rod

The maximum rod temperature, T max , can be represented in the following form:
T max = 0 a r d r / 2 / 2 d z   T ( r , z ) δ ( r ) r δ ( z z max )
where the z max is implicitly defined by Equation (12). The first-order total differential, δ T max , of T max depends, in principle, on all of the first-order sensitivities of T s max with respect to the imprecisely known model and boundary parameters through the following relation:
δ T max = T max q   δ q + T max h δ k + T max h   δ h + T max W   δ W           + T max c p     δ c p + T max   T i n l e t δ T i n l e t + T max a δ a + T max   δ .
The total differential of T max , which will be denoted as δ T max , is obtained by taking the first-order Gateaux- (G-)-differential of Equation (95), which is computed by definition as follows:
δ T max { d d ε 0 a 0 + ε δ a r d r ( 0 + ε δ ) / 2 ( 0 + ε δ ) / 2   [ T 0 ( r , z ) + ε δ T ( r , z ) ] δ ( r ) r δ ( z z max 0 ε δ z max )   d z } ε = 0                         = 0 a 0 r d r 0 / 2 0 / 2   δ T ( r , z ) δ ( r ) r δ ( z z max 0 )   d z .
It is noteworthy that the entire contribution to δ T max stems just from the indirect-effect term; the direct-effect term vanishes for this particular response because of the delta-functionals in its definition. Thus, δ T max can be expressed alternatively in terms of adjoint sensitivity functions which are determined by following the general methodology presented in [1], as has already been illustrated in Section 3. Comparing the expression in Equation (97) to the expression in Equation (46) and following the same procedure as in Section 3 leads to the following expression for δ T max :
δ T max = [ ( δ k ) q 0 k 0 ( δ q ) ] 0 a 0 r d r   0 / 2 0 / 2 Ψ max ( r , z ) cos π z 0 d z + + ( δ ) q 0 π ( 0 ) 2 0 a 0 r d r   0 / 2 0 / 2 Ψ max ( r , z ) z sin π z 0 d z + ( δ ) π 2 ( a 0 ) 2 q 0 W 0 c p 0 ( 0 ) 2 ( a 0 ) 2 2 0 / 2 0 / 2 Ψ f l max ( z )     z sin π z 0 d z + ( a 0 ) 2 2 π [ 2 a 0 q 0 ( δ a ) W 0 c p 0 + ( a 0 ) 2 ( δ q ) W 0 c p 0 ( a 0 ) 2 q 0 ( δ W ) ( W 0 ) 2 c p 0   ( a 0 ) 2 q 0 ( δ c p ) W 0 ( c p 0 ) 2 ] 0 / 2 0 / 2 Ψ f l max ( z ) cos π z 0 d z   + ( a 0 ) 2 2 Ψ f l ( 2 ) ( δ T i n l e t ) + a 0 [ ( δ h ) a 0 q 0 2 h 0 ( δ k ) a 0 q 0 2 k 0 ( δ a ) q 0 2 ( h 0 a 0 k 0 + 1 ) ] 0 / 2 0 / 2   Ψ max ( a 0 , z ) cos π z 0 d z   ,
where the adjoint sensitivity functions Ψ max ( r , z ) and Ψ f l max ( z ) are the solutions of the following 1st-LASS:
k 0 r [ r Ψ max ( r , z ) r ] = δ ( r ) δ ( z z max )   , 0 r < a 0 ,       0 2 z 0 2 ,  
r Ψ max ( r , z ) r = 0 ,       at       r = 0 ,
k 0 Ψ max ( r , z ) r = h 0 Ψ max ( r , z ) ,             at       r = a 0
Ψ f l max ( z ) z + 2 h 0 a 0 Ψ max ( a 0 , z )   = 0 ,   0 2 z 0 2 ,  
Ψ f l max ( z ) =   0 , a t         z = 0 / 2
The expressions of the sensitivities of δ T max in terms of the adjoint sensitivity functions Ψ max ( r , z ) and Ψ f l max ( z ) are formally similar to those in Equation (47) through Equation (54), namely:
T max q = 0 a 0 r d r   0 / 2 0 / 2 Ψ max ( r , z ) cos π z 0 d z + ( a 0 ) 4 W 0 c p 0 π 2 0 / 2 0 / 2 Ψ f l max ( z ) cos π z 0 d z
T max k = q 0 k 0 0 a 0 r d r   0 / 2 0 / 2 Ψ max ( r , z ) cos π z 0 d z ( a 0 ) 2 q 0 2 k 0 0 / 2 0 / 2   Ψ max ( a 0 , z ) cos π z 0 d z
T max h = ( a 0 ) 2 q 0 2 h 0 0 / 2 0 / 2   Ψ max ( a 0 , z ) cos π z 0 d z
T max W = π q 0 ( a 0 ) 4 2 ( W 0 ) 2 c p 0 0 / 2 0 / 2 Ψ f l max ( z ) cos π z 0 d z  
T max c p = π ( a 0 ) 4 q 0 2 W 0 ( c p 0 ) 2 0 / 2 0 / 2 Ψ f l max ( z ) cos π z 0 d z
T max T i n l e t = ( a 0 ) 2 2 Ψ f l max ( 2 )
T max = q 0 π ( 0 ) 2 0 a 0 r d r   0 / 2 0 / 2 Ψ max ( r , z ) z sin π z 0 d z   + π 2 ( a 0 ) 4 q 0 2 W 0 c p 0 ( 0 ) 2 0 / 2 0 / 2 Ψ f l max ( z ) z     sin π z 0 d z ,
T max a = a 0 q 0 2 ( h 0 a 0 k 0 + 1 ) 0 / 2 0 / 2 Ψ max ( a 0 , z ) cos π z 0 d z + π ( a 0 ) 3 q 0 W 0 c p 0 0 / 2 0 / 2 Ψ f l max ( z ) cos π z 0 d z   .    
Solving the 1st-LASS represented by Equations (99)–(103) yields the following expressions for the sensitivity functions Ψ max ( r , z ) and Ψ f l max ( z ) :
Ψ max ( r , z ) = δ ( z z max ) [ 1 a 0 h 0 + H ( r ) ln r k 0 ln a 0 k 0 ]   ,       0 r a 0 ,       0 2 z , z max 0 2
Ψ f l max ( z ) =   2 ( a 0 ) 2 H ( z max z ) ,         2 z , z max 2  
Inserting the expressions obtained in Equations (112) and (113) into Equation (104) through Equation (111) and carrying out the respective integrations yields the following expressions for the first-order sensitivities of the maximum rod temperature, T max :
T max q = [ ( a 0 ) 2 4 k 0 + a 0 2 h 0 ] cos π z max 0 + ( a 0 ) 2 0 W 0 c p 0 ( sin   π z max 0 + 1 )
T max k = q 0 4 ( a 0 k 0 ) 2 cos π z max 0
T max h = q 0 a 0 2 ( h 0 ) 2 cos π z max 0
T max W = q 0 0 c p 0 ( a 0 W 0 ) 2 ( 1 + sin   π z max 0 )
T max c p = q 0 0 W 0 ( a 0 c p 0 ) 2 ( 1 + sin   π z max 0 )
T max T i n l e t = 1
T max = q 0 π 2 ( a 0 0 ) 2 ( 1 a 0 h 0 + 1 2 k 0 ) z max sin π z max 0         + q 0 ( a 0 ) 2 W 0 c p 0 ( π 0 z max cos π z max 0 + sin π z max 0 + 1 ) ,
T max a = q 0 2 ( a 0 k 0 + 1 h 0 ) cos π z max 0 + 2 a 0 q 0 0 W 0 c p 0 ( 1 + sin   π z max 0 )

4.2. First-Order Sensitivities of the Critical Point (Maximum) of the Rod Temperature

Following the principles presented in the accompanying work [6], the relations provided in Equation (12), which implicitly define the location z max , are written in the following form:
0 a r d r / 2 / 2 T ( r , z ) z δ ( r ) r δ ( z z max ) d z   = 0
The first-order total differential, δ z max , of z max depends, in principle, on all of the first-order sensitivities of z max with respect to the imprecisely known model and boundary parameters through the following relation:
δ z max = z max q   δ q + z max h δ k + z max h   δ h + z max W   δ W           + z max c p     δ c p + z max   T i n l e t δ T i n l e t + z max a δ a + z max   δ .
Taking the G-differential of Equation (122) yields the following relation:
0 =   { d d ε 0 a 0 + ε δ a d r ( 0 + ε δ ) / 2 ( 0 + ε δ ) / 2   [ T ( r , z ) + ε δ T ( r , z ) ] z δ ( r ) δ ( z z max 0 ε δ z max ) d z } ε = 0 =   ( δ z max ) 0 a 0 d r 0 / 2 0 / 2   [ T ( r , z ) ] z δ ( r ) δ ( z z max 0 ) d z     + 0 a 0 d r 0 / 2 0 / 2   [ δ T ( r , z ) ] z δ ( r ) δ ( z z max 0 ) d z   .
The relation in Equation (124) can be re-written in the following form:
( δ z max ) {   2 [ T ( r , z ) ] z 2 } ( r = 0 , z = z max 0 ) = 0 a 0 d r 0 / 2 0 / 2   δ T ( r , z ) δ ( r ) δ ( z z max 0 ) d z
where the quantity {   2 [ T ( r , z ) ] / z 2 } ( r = 0 , z = z max 0 ) is computed using Equation (6) for this benchmark model (or is computed numerically when the closed form of T ( r , z ) is unavailable) to obtain
{   2 [ T ( r , z ) ] z 2 } ( r = 0 , z = z max 0 ) = { q a ( π ) 2 [ a h + 2 k 4 k h cos π z + a W c p sin π z ] } ( z = z max 0 )
The term 0 a 0 d r 0 / 2 0 / 2   δ T ( r , z ) δ ( r ) δ ( z z max 0 ) d z in Equation (125) can be expressed in terms of adjoint functions by applying the procedure outlined in [6], to obtain the following result:
0 a 0 d r 0 / 2 0 / 2   δ T ( r , z ) δ ( r ) δ ( z z max 0 ) d z   = [ ( δ k ) q 0 k 0 ( δ q ) ] 0 a 0 r d r   0 / 2 0 / 2 u ( r , z ) cos π z 0 d z + + ( δ ) q 0 π ( 0 ) 2 0 a 0 r d r   0 / 2 0 / 2 u ( r , z ) z sin π z 0 d z + ( δ ) π 2 ( a 0 ) 2 q 0 W 0 c p 0 ( 0 ) 2 ( a 0 ) 2 2 0 / 2 0 / 2 u f l ( z )     z sin π z 0 d z + ( a 0 ) 2 2 π [ 2 a 0 q 0 ( δ a ) W 0 c p 0 + ( a 0 ) 2 ( δ q ) W 0 c p 0 ( a 0 ) 2 q 0 ( δ W ) ( W 0 ) 2 c p 0   ( a 0 ) 2 q 0 ( δ c p ) W 0 ( c p 0 ) 2 ] 0 / 2 0 / 2 u f l ( z ) cos π z 0 d z   + ( a 0 ) 2 2 u f l ( 2 ) ( δ T i n l e t ) + a 0 [ ( δ h ) a 0 q 0 2 h 0 ( δ k ) a 0 q 0 2 k 0 ( δ a ) q 0 2 ( h 0 a 0 k 0 + 1 ) ] 0 / 2 0 / 2   u ( a 0 , z ) cos π z 0 d z .
where the adjoint functions u ( r , z ) and u f l ( z ) are the solutions of the following 1st-LASS:
k 0 r [ r u ( r , z ) r ] = δ ( r ) δ ( z z max 0 )   ,
r u ( r , z ) r = 0 ,       at       r = 0 ,
k 0 u ( r , z ) r = h 0 u ( r , z ) ,             at       r = a 0
u f l ( z ) z + 2 h 0 a 0 u ( a , z )   = 0 ,   0 2 z 0 2 ,  
u f l ( z ) =   0 , at         z = 0 / 2
Introducing the results obtained in Equation (128) through Equation (132) into Equation (127) and equating the coefficients corresponding to the same parameter variations yields the following expressions for the sensitivities of the location z max to the imprecisely known model, interface, and boundary parameters:
z max q = { 0 a 0 r d r   0 / 2 0 / 2 u ( r , z ) cos π z 0 d z ( a 0 ) 2 2 π ( a 0 ) 2 W 0 c p 0 0 / 2 0 / 2 u f l ( z ) cos π z 0 d z ( π ) 2 q a ( a h + 2 k 4 k h cos π z + a W c p sin π z ) } z = z max
z max k = { q 0 k 0 0 a 0 r d r   0 / 2 0 / 2 u ( r , z ) cos π z 0 d z + ( a 0 ) 2 q 0 2 k 0 0 / 2 0 / 2   u ( a 0 , z ) cos π z 0 d z ( π ) 2 q a ( a h + 2 k 4 k h cos π z   + a W c p sin π z ) } z = z max
z max h = {   ( a 0 ) 2 q 0 2 h 0 0 / 2 0 / 2   u ( a 0 , z ) cos π z 0 d z . ( π ) 2 q a ( a h + 2 k 4 k h cos π z   + a W c p sin π z ) } z = z max
z max W   = { ( a 0 ) 2 2 π ( a 0 ) 2 q 0 ( W 0 ) 2 c p 0   0 / 2 0 / 2 u f l ( z ) cos π z 0 d z ( π ) 2 q a ( a h + 2 k 4 k h cos π z   + a W c p sin π z ) } z = z max
z max c p   = { ( a 0 ) 2 2 π ( a 0 ) 2 q 0 W 0 ( c p 0 ) 2   0 / 2 0 / 2 u f l ( z ) cos π z 0 d z ( π ) 2 q a ( a h + 2 k 4 k h cos π z   + a W c p sin π z ) } z = z max
z max   = { q 0 π ( 0 ) 2 0 a 0 r d r   0 / 2 0 / 2 u ( r , z ) z sin π z 0 d z + π 2 ( a 0 ) 4 q 0 2 W 0 c p 0 ( 0 ) 2 0 / 2 0 / 2 u f l ( z )     z sin π z 0 d z ( π 0 ) 2 q 0 a 0 ( a h + 2 k 4 k h cos π z 0   + a 0 0 W 0 c p 0 sin π z 0 ) } z = z max
z max a = { ( π ) 2 q a ( a h + 2 k 4 k h cos π z max 0   + a W c p sin π z max 0 ) } 1 { a 0 q 0 2 k 0 π 0 sin π z s 0 0 + ( a 0 ) 2 2 π 2 a 0 q 0 W 0 c p 0 0 / 2 0 / 2 u f l ( z ) cos π z 0 d z a 0 q 0 2 ( h 0 a 0 k 0 + 1 ) 0 / 2 0 / 2   u ( a 0 , z ) cos π z 0 d z } z = z max
z max T i n l e t = { ( a 0 ) 2 2 u f l ( 2 ) ( π ) 2 q a ( a h + 2 k 4 k h cos π z   + a W c p sin π z ) } z = z max
The expressions obtained in Equation (133) through Equation (140) can be evaluated after solving (numerically or analytically) the 1st-LASS represented by Equation (128) through Equation (132) to determine the adjoint functions u ( r , z ) and u f l ( z ) . Since u ( r , z ) and u f l ( z ) do not depend on any parameter variations, this 1st-LASS needs to be solved only once to obtain the respective adjoint functions. Solving the 1st-LASS analytically yields the following expressions for the adjoint functions u ( r , z ) and u f l ( z ) :
u ( r , z ) = δ ( z z max 0 ) [ 1 a 0 h 0 + H ( r ) ln r k 0 ln a 0 k 0 ] ,       0 2 z , z max 0 0 2 ,
u f l ( z ) = 2 δ ( z z max 0 ) ( a 0 ) 2 ,         0 2 z , z max 0 0 2
Inserting the expressions for u ( r , z ) and u f l ( z ) provided in Equations (141) and (142) into Equation (133) through Equation (140) and performing the respective integrations yields the following closed-form expressions for the sensitivities of the location z max to the imprecisely known model, interface, and boundary parameters:
z max q = 0
z max k = 4 a 2 2 h 2 π W c p ( a h + 2 k ) 2 [ 1   + x 2 ] 1
z max h = 8 a 2 k 2 π W c p ( a h + 2 k ) 2 [ 1   + x 2 ] 1
z max W   = 4 a 2 k h π W 2 c p ( a h + 2 k ) [ 1   + x 2 ] 1
z max c p   = 4 a 2 k h π W c p 2 ( a h + 2 k )   [ 1   + x 2 ] 1
z max   = z max + 4 a k h W c p ( a h + 2 k )   [ 1   + x 2 ] 1 ,
z max a = 8 2 k 2 a h π W c p ( a h + 2 k ) 2   [ 1   + x 2 ] 1
z max T i n l e t = 0
where
x 4 a 2 h π W c p ( a h + 2 k )
The subscript “zero” (which indicates “nominal values”) has been omitted, for simplicity, in writing the expressions given in Equation (143) through Equation (151), but all of the quantities which appear in the respective expressions are to be evaluated/computed at the nominal values of the respective quantities.

5. Discussion

This work has applied the 1st-CASAM-CP methodology pioneered in [6] to develop a thermal-hydraulics benchmark model for the verification of the adjoint sensitivity functions and sensitivities of response defined at critical points for thermal-hydraulics computational tools such as the “Adjoint FLUENT Solver” [4]. It has been shown that the computation of the first-order response sensitivities of the magnitude of a response defined at a critical point of a function in the phase-space of the systems’ independent variables requires a single large-scale computation of the 1st-LASS which corresponds to the magnitude of the critical point. Solving the 1st-LASS represents a “large-scale” computation, which is not more extensive, however, than solving the original coupled systems. This is because the 1st-LASS is linear in the dependent variables, whereas the original coupled systems are nonlinear in the dependent variables. Furthermore, the computation of the sensitivities of the location in phase-space of each critical point requires solving one 1st-LASS for each of the components in the phase-space of independent variables of the respective critical point. The same operators appear on the right-sides of the 1st-level adjoint sensitivity systems needed for computing the adjoint functions corresponding to the magnitude and the location (in phase-space) of a response defined at a critical point. Only the sources on the left-sides of the respective 1st-LASS differ from each other. Therefore, the same solver can be used (albeit with different sources) to compute the adjoint functions. These unique characteristics of the general 1st-CASAM-CP methodology [6] have been demonstrated in this work by considering two responses of paramount importance in reactor safety, namely, (i) the maximum rod surface temperature, which occurs at the imprecisely known interface between the subsystem that models the heat conduction inside the heated rod and the subsystem modeling the heat convection process surrounding the rod; and (ii) the maximum temperature inside the heated rod, which has two components, one located at a precisely known boundary of the subsystem that models the heat conduction inside the heated rod, while the other component depends on an imprecisely known boundary (i.e., the rod length). The exact analytical expressions developed in this work for the sensitivities of the maximum internal rod temperature and maximum rod surface temperature, as well as the sensitivities of the respective locations in phase-space where these respective maxima occur, provide accurate benchmarks for verifying the accuracy of computational tools for modeling thermal-hydraulics systems. This “solution verification” process should include the verification of the accuracy of the sensitivities of the results computed by such codes to the uncertain parameters underlying the respective codes.
Even the most advanced current thermal-hydraulics and/or CFD-computational tools, such as the Adjoint-FLUENT solver [4], are insufficiently developed and currently lack capabilities for computing sensitivities of the locations of critical points (maxima, minima, saddle points) with respect to model parameters, interfaces, and boundaries. The thermal-hydraulics benchmark presented in this work serves for the future verification/validation of the future developments/extensions of the “adjoint sensitivity analysis” capabilities of future thermal-hydraulics and/or CFD computational tools, particularly for the verification of the adjoint functions needed to compute the crucially-important (for reactor design and licensing) “maximum fuel temperature” and “maximum clad temperature” in the reactor’s hot channels. Applications of the 1st-CASAM-PC [6] and use of the benchmark presented in this work for adjoint sensitivity analysis of “responses at critical points” must eventually come in the future, since such responses are crucially important for the design, optimization, and licensing of engineering systems, see e.g., Refs. [7,8]

Funding

This research received no external funding.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Not applicable.

Conflicts of Interest

The author declares no conflict of interest.

References

  1. Cacuci, D.G.; Fang, R.; Ilic, M.; Badea, M.C. A Heat Conduction and Convection Analytical Benchmark for Adjoint Solution Verification of CFD Codes Used in Reactor Design. Nucl. Sci. Eng. 2015, 182, 452–480. [Google Scholar] [CrossRef]
  2. Cacuci, D.G. Second-Order Adjoint Sensitivity and Uncertainty Analysis of a Heat Transport Benchmark Model—I: Analytical Results. Nucl. Sci. Eng. 2016, 183, 1–21. [Google Scholar] [CrossRef]
  3. Cacuci, D.G.; Ilic, M.; Badea, M.C.; Fang, R. Second-Order Adjoint Sensitivity and Uncertainty Analysis of a Heat Transport Benchmark Model—II: Computational Results Using G4M Reactor Thermal-Hydraulic Parameters. Nucl. Sci. Eng. 2016, 183, 22–38. [Google Scholar] [CrossRef]
  4. ANSYS® Academic Research. Release 16.0, FLUENT Adjoint Solver; ANSYS, Inc.: Canonsburg, PA, USA, 2015. [Google Scholar]
  5. Gen4 Energy, Inc. G4M Reactor Core Design; Gen4 Energy, Inc.: Denver, Colorado, USA, 2012. [Google Scholar]
  6. Cacuci, D.G. First-Order Comprehensive Adjoint Sensitivity Analysis Methodology for Critical Points in Coupled Nonlinear Systems. I: Mathematical Framework. Fluids 2020. submitted. [Google Scholar]
  7. Cacuci, D.G. Sensitivity Theory for Nonlinear Systems: II. Extensions to Additional Classes of Responses. J. Math. Phys. 1981, 22, 2803–2812. [Google Scholar] [CrossRef]
  8. Cacuci, D.G.; Maudlin, P.J.; Parks, C.V. Adjoint Sensitivity Analysis of Extremum-Type Responses in Reactor Safety. Nucl. Sci. Eng. 1983, 83, 112–135. [Google Scholar] [CrossRef]
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Cacuci, D.G. First-Order Comprehensive Adjoint Sensitivity Analysis Methodology for Critical Points in Coupled Nonlinear Systems. II: Application to a Nuclear Reactor Thermal-Hydraulics Safety Benchmark. Fluids 2021, 6, 34. https://doi.org/10.3390/fluids6010034

AMA Style

Cacuci DG. First-Order Comprehensive Adjoint Sensitivity Analysis Methodology for Critical Points in Coupled Nonlinear Systems. II: Application to a Nuclear Reactor Thermal-Hydraulics Safety Benchmark. Fluids. 2021; 6(1):34. https://doi.org/10.3390/fluids6010034

Chicago/Turabian Style

Cacuci, Dan Gabriel. 2021. "First-Order Comprehensive Adjoint Sensitivity Analysis Methodology for Critical Points in Coupled Nonlinear Systems. II: Application to a Nuclear Reactor Thermal-Hydraulics Safety Benchmark" Fluids 6, no. 1: 34. https://doi.org/10.3390/fluids6010034

APA Style

Cacuci, D. G. (2021). First-Order Comprehensive Adjoint Sensitivity Analysis Methodology for Critical Points in Coupled Nonlinear Systems. II: Application to a Nuclear Reactor Thermal-Hydraulics Safety Benchmark. Fluids, 6(1), 34. https://doi.org/10.3390/fluids6010034

Article Metrics

Back to TopTop