# Isolated Rotor Blade Shape Sensitivity for Aeroacoustic Optimization Using a Discrete Adjoint Framework

Yacine Mohammedi, Majd Daroukh, Martin Buszyk, Antoine Hajczak, Itham

Salah El Din, Marc Bonnet

► To cite this version:

Yacine Mohammedi, Majd Daroukh, Martin Buszyk, Antoine Hajczak, Itham Salah El Din, et al.. Isolated Rotor Blade Shape Sensitivity for Aeroacoustic Optimization Using a Discrete Adjoint Framework. AIAA AVIATION FORUM AND ASCEND 2025, Jul 2025, Las Vegas, United States. 10.2514/6.2025-3367. hal-05190189

HAL Id: hal-05190189

https://hal.science/hal-05190189v1

Submitted on 28 Jul 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Isolated rotor blade shape sensitivity for aeroacoustic optimization using a discrete adjoint framework

Yacine Mohammedi*,<fnref n="1" />,<fnref n="4" />, Majd Daroukh†,<fnref n="1" />, Martin Buszyk‡,<fnref n="1" />, Antoine Hajczak§,<fnref n="2" />, Itham Salah El-Din¶,<fnref n="3" />, and Marc Bonnet||,<fnref n="4" />

A discrete adjoint framework is developed to optimize rotor self-noise from steady fluid simulations in the rotating frame. To this end, a simplified expression of the off-body frequency-domain Ffowcs-Williams and Hawkings (FW-H) equation is derived for far-field observers, following the model of Hanson and Parzych (1993) originally written for on-body surfaces. The latter is implemented and compared against the results given by an established time-domain FW-H solver. Far-field acoustic pressure sensitivities are derived analytically and validated by comparison with second-order accurate finite differences. The sensitivities of any objective function expressed in terms of the acoustic pressure can therefore be reconstructed. Then the discrete adjoint of a Reynolds-averaged Navier-Stokes solver provides the objective function gradients with respect to the blade shape parameters. The complete workflow is validated against finite difference evaluations on an isolated open rotor in cruise conditions.

## Nomenclature

### Roman letters

$\mathcal{A}_\alpha$ = set of admissible designs

$B$ = number of blades

$E$ = fluid total energy

$F$ = FW-H dipole term

$\mathcal{F}$ = objective function

$G_c$ = free space convected Green's function

$H$ = Heaviside function

$i$ = imaginary unit

$\mathcal{J}$ = Bessel function

$\mathcal{L}$ = Lagrangian

$M$ = flow Mach number

$M_{\text{flight}}$ = flight Mach number

$m$ = harmonic index

$\mathbf{n}$ = outer normal

$n_\alpha$ = number of design parameters

$p$ = mean flow static pressure

$\mathcal{Q}$ = FW-H monopole term

$\tilde{\mathcal{R}}$ = discrete residual of Navier-Stokes equations

$r_0$ = rotation radius

$T$ = FW-H quadrupole term

$t$ = observer reception time

$\mathbf{u}$ = mean flow velocity in Cartesian coordinates

$\mathbf{W}$ = vector of flow variables

$\mathbf{X}$ = volumic CFD mesh coordinates

$\mathbf{x}$ = observer position

$\mathbf{y}$ = source position

### Greek letters

$\alpha$ = vector of shape parameters

$\beta$ = Prandtl-Glauert factor

$\gamma$ = air adiabatic index

$\delta$ = Dirac distribution

$\delta_{ij}$ = Kronecker symbol

$\kappa_\infty$ = acoustic wavenumber

$\lambda$ = Lagrange multiplier or adjoint vector

$\mu$ = dynamic viscosity

$\rho$ = mean flow density

$\Sigma$ = FW-H integration surface

$\sigma$ = stress tensor

$\tau_{ij}$ = viscous stress tensor

$\tau$ = source emission time

$\phi_0$ = azimuthal angle of the source

$\psi$ = polar angle

$\Omega$ = angular velocity of the rotor

### Superscripts

$\square'$ = fluctuation

$\square^*$ = optimal

$\square$ = polar projection

### Subscripts

$\square_\infty$ = far-field state

### Operators

$\square^T$ = transposition

$\square$ = complex conjugate

$\hat{\square}$ = Fourier Transform

<footnote n="1">¹DAAA, ONERA, Institut Polytechnique de Paris, 92320 Châtillon, France</footnote>

<footnote n="2">²Sorbonne Université, CNRS, Institut Jean Le Rond d'Alembert, F-75005 Paris, France</footnote>

<footnote n="3">³DAAA, ONERA, Institut Polytechnique de Paris, 92190 Meudon, France</footnote>

<footnote n="4">⁴POEMS (CNRS, INRIA, ENSTA), ENSTA, 91120 Palaiseau, France</footnote>

<footnote>*Ph.D. candidate, yacine.mohammedi@onera.fr</footnote>

<footnote>†Research engineer</footnote>

<footnote>‡Research engineer</footnote>

<footnote>§Associate professor</footnote>

<footnote>¶Research engineer</footnote>

<footnote>||CNRS senior researcher</footnote>
# I. Introduction

With the goal of improving aircraft aerodynamic performance and energy efficiency, the aeronautic industry aims to advance propulsion architectures, notably integrating new systems like open fan engines. With their theoretical infinite bypass ratio, significant fuel savings are expected with these architectures. However, particular attention should be paid to the noise they produce because the absence of nacelle prevents from using acoustic liners or benefiting from duct mode cutting off.

Recent studies [1–3] have focused on open fan tonal noise by combining direct near-field Computational Fluid Dynamics (CFD) and the Ffowcs-Williams and Hawkings (FW-H) [4] integral method. Among them, Lewis et al. [3] and Sugiyama et al. [2] demonstrated that rotor shock noise was dominant at cruise speed. It radiates directly towards the fuselage due to the absence of a nacelle and becomes a source of cabin noise affecting, in particular, passengers comfort.

Minimizing the shock noise from an open rotor can be achieved by optimizing the rotor blade shape. Traditional parametric studies become computationally expensive with many shape parameters, especially when CFD is involved, making efficient shape optimization methods essential. With regard to the optimization process, deterministic gradient-based methods [5, 6] are well suited to the local optimization of industrial configurations. Finding an estimate of the gradient of the functions of interest with respect to shape parameters using the finite differences method quickly becomes prohibitive since the cost of calculation is proportional to the number of design parameters. However, the adjoint method has a calculation cost independent of the number of design variables [7, 8], making it ideal for high-dimensional optimization. In this work, the discrete adjoint approach, where the gradient is obtained from the exact differentiation of the approximate governing equations, is preferred for its flexibility, its easier implementation and its efficiency for large-scale problems. Widely used in aerodynamics [9, 10], aero-structure [11], and aeroelasticity [12, 13], the discrete adjoint method remains still recent when the acoustic contribution is taken into account in the objective function. Algorithmic differentiation is widely adopted in discrete adjoint frameworks to compute sound pressure sensitivities [14–17]. While analytical differentiation is less common, it remains valuable for physical insight and for evaluating the impact of individual terms [18–20].

A preliminary study to this work was carried out by Foglia et al. [21] who developed a steady adjoint framework by computing acoustic pressure sensitivities analytically via Hanson and Parzych (HP) formulation [22]. The latter is a simplified expression of the on-body frequency-domain FW-H integral formulation under far-field assumption and written in the rotating domain. However, this formulation is not suitable for transonic flows, as it only captures far-field noise from surfaces coincident with the blades, missing shock noise near the blades. This paper extends that framework to shock noise by employing a variant of Hanson and Parzych’s formulation adapted to off-body control surfaces encompassing both volume and non-linear sources.

The methodology for optimizing rotor self-noise using the discrete adjoint method is presented in Sec. II. Then, the entire methodology is validated on an isolated open rotor in cruise conditions in Sec. III. Finally, conclusions are provided in Sec. IV.

# II. Discrete adjoint framework for estimating the acoustic metric gradient

## A. Overview

Consider the following acoustic optimization problem under a fluid mechanics constraint:

$$ \alpha^* = \arg \min_{\alpha \in \mathcal{A}_{\alpha}} \mathcal{F}(\mathbf{W}(\alpha), \mathbf{X}(\alpha)) \quad \text{subject to} \quad \mathcal{R}(\mathbf{W}(\alpha), \mathbf{X}(\alpha)) = 0, \qquad (1) $$

where $\alpha$ is the vector of design parameter, $\alpha^*$ the optimal vector satisfying the problem, $\mathcal{F}$ the objective function characterizing a particular acoustic metric which will be detailed in Sec. III.A, $\mathcal{A}_{\alpha}$ the set of admissible designs, $\mathbf{X}$ the fluid mesh coordinates, $\mathcal{R}$ the residual of the discretized Reynolds Average Navier-Stokes (RANS) equations and $\mathbf{W}$ the vector of flow variables. The gradient $d\mathcal{F}/d\alpha$ of the objective function with respect to the design variables $\alpha$ is to be evaluated as in Eq. (2) by satisfying the (RANS) constraint $\mathcal{R} = 0$ for any admissible design, i.e. by enforcing $d\mathcal{R}/d\alpha = 0$:

$$ \frac{d\mathcal{F}}{d\alpha} = \frac{d\mathcal{F}}{d\alpha} + \lambda^T \frac{d\mathcal{R}}{d\alpha} = \frac{d\mathcal{L}}{d\alpha}, \qquad (2) $$

where $\mathcal{L}$ is the Lagrangian associated to the optimization problem of Eq. (1), given by:

$$ \mathcal{L}(\alpha, \lambda) = \mathcal{F}(\alpha) + \lambda^T \mathcal{R}(\mathbf{W}(\alpha), \mathbf{X}(\alpha)), \qquad (3) $$
and $\lambda \in \mathbb{R}^{6 \times n_{\text{cells}}}$ (5 variables related to the flow equations and 1 extra variable for the one-equation turbulence model), with $n_{\text{cells}}$ the number of cells, is hence the Lagrange multiplier for the constraint $\mathcal{R}(\mathbf{W}(\alpha), \mathbf{X}(\alpha)) = 0$. On differentiating $\mathcal{F}$ and $\mathcal{R}$ with respect to $\alpha$ by means of the chain rule and rearranging the resulting expression, Eq. (2) provides:

$$ \frac{d\mathcal{F}}{d\alpha} = \left( \frac{\partial \mathcal{F}}{\partial \mathbf{W}} + \lambda^T \frac{\partial \mathcal{R}}{\partial \mathbf{W}} \right) \frac{d\mathbf{W}}{d\alpha} + \left( \frac{\partial \mathcal{F}}{\partial \mathbf{X}} + \lambda^T \frac{\partial \mathcal{R}}{\partial \mathbf{X}} \right) \frac{d\mathbf{X}}{d\alpha}. \quad (4) $$

The derivative formula of Eq. (4) a priori involves the design derivatives $d\mathbf{W}/d\alpha$, whose actual computation is potentially prohibitively costly as it requires performing a number of CFD simulations proportional to the number of design variables denoted by $n_{\alpha}$. However, taking advantage of the fact that Eq. (4) holds for any $\lambda \in \mathbb{R}^{6 \times n_{\text{cells}}}$, this problem can be avoided by selecting $\lambda$ so as to cancel the cofactor of $d\mathbf{W}/d\alpha$ in Eq. (4). We thus fix $\lambda$ by requiring that it satisfies the adjoint equation:

$$ \left( \frac{\partial \mathcal{R}}{\partial \mathbf{W}} \right)^T \lambda = - \left( \frac{\partial \mathcal{F}}{\partial \mathbf{W}} \right)^T. \quad (5) $$

The fact is that this equation does not depend on the vector of shape parameters $\alpha$ and only needs to be solved once, hence the interest of the method. Using the (discrete) adjoint state $\lambda$ defined by Eq. (5), the design derivatives of the objective functional are given by:

$$ \frac{d\mathcal{F}}{d\alpha} = \left( \frac{\partial \mathcal{F}}{\partial \mathbf{X}} - \frac{\partial \mathcal{F}}{\partial \mathbf{W}} \left( \frac{\partial \mathcal{R}}{\partial \mathbf{W}} \right)^{-1} \frac{\partial \mathcal{R}}{\partial \mathbf{X}} \right) \frac{d\mathbf{X}}{d\alpha} \quad (6) $$

Now, the gradient of the functional with respect to the shape parameters can be assembled from its sensitivities, those of the discrete residual and the gradient of the volume mesh. The calculation of each of its terms constitutes a step in the process of estimating the gradient of the function of interest using the discrete adjoint method, as illustrated in Fig. 1.

Based on the design parameter values, a mesh deformation tool is used to update the mesh, yielding the sensitivity of mesh points to a change in each shape parameter. A RANS simulation provides the aerodynamic field over the entire fluid domain. This field is extracted on the blade surface, or on a fictitious surface encompassing the blade as well as possible non-negligible volume sources. Based on the knowledge of this field on the control surface, an acoustic analogy derived from the Ffowcs-Williams and Hawkings (FW-H) equation is applied to radiate acoustics up to the far-field. Analogously, the differentiated version of the acoustic analogy is applied to obtain in turn the sensitivities of acoustic pressure with respect to flow variables and with respect to the fluid mesh. By chain rule differentiation, the sensitivities of any acoustic objective function as a function of acoustic pressure can be calculated. These sensitivities serve as inputs to a discrete adjoint RANS solver, which in turn provides the gradients of the functional with respect to the shape parameters. In the next sections, the various steps involved in this workflow are given in more details.
<figure><img src="image_2.png" /><figcaption>Fig. 1 Gradient estimation proposed workflow with the discrete adjoint method</figcaption></figure>

## B. Prediction of rotor self-noise

This section presents the forward problem of predicting tonal rotor noise using a simplified expression of the far-field off-body frequency domain FW-H formulation in the rotating frame.

### 1. Forward CFD simulation

The first step in predicting rotor self noise is to perform a steady CFD calculation in the rotating frame of reference. To do so, the steady Navier-Stokes equations, written here as:

$$ \mathcal{R}_{\text{exp}}(\mathbf{W}(\alpha), \mathbf{X}(\alpha)) = 0 \qquad (7) $$

are solved numerically on the fluid mesh. The flow solution is the vector of conservative flow variables $\mathbf{W} = (\rho, \rho u_1, \rho u_2, \rho u_3, \rho E, \rho v_t)$ with $\rho, \mathbf{u} = (u_1, u_2, u_3)$ and $E$ the flow density, velocity vector and total energy respectively. $\mathcal{R}_{\text{exp}}$ is called the explicit residual of the Navier-Stokes equations, and is defined as:

$$ \mathcal{R}_{\text{exp}} = \left( \frac{\partial \rho u_j}{\partial x_j} \quad \frac{\partial \rho u_i u_j}{\partial x_j} - \frac{\partial \sigma_{ij}}{\partial x_j} \quad \frac{\partial \rho E u_j}{\partial x_j} + \frac{\partial p}{\partial x_j} \delta_{ij} + \frac{\partial q_j}{\partial x_j} - \frac{\partial u_i}{\partial x_j} \tau_{ij} \right) \qquad (8) $$

where $u_j$ stands for the $j^{th}$ velocity component, $v$ is the kinematic viscosity and $q_j$ corresponds to the $j^{th}$ heat flux vector component determined by the Fourier law. Moreover, the fluid stress tensor $\sigma_{ij}$ has the usual decomposition $\sigma_{ij} = -p\delta_{ij} + \tau_{ij}$, where $\delta_{ij}$ is the Kronecker symbol and the pressure $p$ follows the ideal gas law and the viscous stress tensor $\tau_{ij}$ follows the linear Newtonian constitutive model:

$$ p = \frac{(\gamma - 1)}{\rho} \left( \rho E - \frac{1}{2} \rho u_k \rho u_k \right) \qquad (9a) $$

with $\mu$ the dynamic viscosity. Regarding the Reynolds Averaged Navier-Stokes (RANS) equations, the conservative flow variables are decomposed into a time-averaged part and a fluctuating part. This results in an modified explicit residual whose expression is closed to the explicit residual $\mathcal{R}_{\text{exp}}$ of Navier-Stokes equations, but with several key modifications: the inclusion of the Reynolds stress tensor $\tau_{t,ij}$ in addition to the viscous stress tensor $\tau_{ij}$, the addition of
the turbulent heat flux vector $q_t$ to the heat flux $q$, and the introduction of the operator $SA(v_t)$ to ensure closure of the system [23]. This operator corresponds to the equation of the Spalart-Allmaras turbulence model for the turbulent eddy viscosity $v_t$ described in [24].

Employing discretization schemes, which will be detailed in Sec. III.B, for each component of the adapted expression of $R_{\text{exp}}$ to RANS equations leads to the following equation:

$$ \mathcal{R}(\mathbf{W}(\alpha), \mathbf{X}(\alpha)) = 0 \qquad (10) $$

where $\mathcal{R}$ is the residual of the discretized RANS equations.

## 2. Variant of Hanson and Parzych integral formulation

Typically, the acoustic functional $\mathcal{F}$ to be optimized depends on the acoustic field, which is characterized by the pressure fluctuation $p'$. The latter is determined using an integral formulation suitable for a rotating source, where the required source terms are obtained from conservative flow variables extracted on a control surface denoted by $\Sigma$ following the CFD simulation.

Consider an aircraft with uniform translational motion at speed $c_{\infty}M_{\text{flight}}$ along the $\mathbf{e}_1$ axis. $M_{\text{flight}}$ is the flight Mach number of opposite sign to the fluid flow Mach number $M > 0$, as illustrated in Fig. 2 and $c_{\infty}$ is the far-field sound velocity. The fan is animated not only by this uniform translational motion, but also by a rotational motion with angular velocity $\Omega$. As illustrated in Fig. 2, $(\mathbf{e}_1, \mathbf{e}_2, \mathbf{e}_3)$ is the uniform translation frame of reference with speed $c_{\infty}M_{\text{flight}}$ while $(\tilde{\mathbf{e}}_1, \tilde{\mathbf{e}}_2, \tilde{\mathbf{e}}_3)$ is the frame of reference associated with the rotor blade. $\mathbf{x} = (x_1, x_2, x_3)$ denotes the vector of coordinates of an observer placed in the far-field and $\mathbf{y} = (y_1, y_2, y_3)$ represents the vector of coordinates of a source point.

<figure><img src="image_3.png" /><figcaption>Fig. 2 Definition of the translational reference frame and the reference frame attached to the blade</figcaption></figure>

In order to evaluate the fluctuating pressure $p'$ generated aerodynamically by a flow bounded by surfaces in arbitrary motion, Ffowcs Williams and Hawkings [4] established a generalization of Lighthill's acoustic analogy, which was later adapted by Goldstein [25] to account for mean flow effects, as expressed in equation (11):

$$ \frac{D^2}{Dt^2} \left[ p' H(f) \right] - c_{\infty}^2 \frac{\partial^2}{\partial x_i^2} \left[ p' H(f) \right] = \frac{\partial^2}{\partial x_i \partial x_j} \left[ T_{ij} H(f) \right] - \frac{\partial}{\partial x_i} \left[ F_{ij} n_j \delta(f) \right] + \frac{D}{Dt} \left[ Q_j n_j \delta(f) \right] \quad (11) $$

where $n_j$ is the $j^{th}$ component of the outer normal to the surface $\Sigma$, $\delta$ the Dirac distribution and $H$ the Heaviside function. $f$ is defined as a function of distance from the interface. This function is respectively negative, positive or zero if the observer is respectively inside, outside or on the surface. $D/Dt$ is the convected derivative operator defined as:

$$ \frac{D\square}{Dt} = \frac{\partial\square}{\partial t} - c_{\infty} M_{\text{flight}} \frac{\partial\square}{\partial y_1} \qquad (12) $$
The quadrupolar term $T_{ij}$ corresponds to the shear noise generated by the turbulence within the fluid, the dipolar term $F_{ij}$ to the loading noise associated with the pressure forces applied to the blade and the monopolar term $Q_j$ to the thickness noise due to the volume flow of air expelled by the rotating blade:

$$T_{ij} = \rho u_i u_j - \tau_{ij} + (p' - c_{\infty}^2 \rho') \delta_{ij} \quad (13a)$$

$$F_{ij} = p' \delta_{ij} - \tau_{ij} + \rho u_i (u_j - u_j^\Sigma) \quad (13b)$$

$$Q_j = \rho_{\infty} u_j^{\Sigma} + \rho(u_j - u_j^{\Sigma}) \quad (13c)$$

with $u_1 = u_1 + c_{\infty}M_{\text{flight}}$ since CFD simulation is carried out in the translational relative frame of reference and $u^\Sigma$ represents the control surface velocity. The solution of Eq. (11) is obtained by convolving each of the source terms with the 3D free space convected Green function $G_c$ expressed as:

$$G_c(\mathbf{x}, t|\mathbf{y}, \tau) = \frac{\delta(t - \tau - \frac{r_{\beta} - M_{\text{flight}}(x_1 - y_1)}{c_{\infty}\beta^2})}{4\pi r_{\beta}} \quad (14)$$

with $r_{\beta} = \sqrt{(x_1 - y_1)^2 + \beta^2 ((x_2 - y_2)^2 + (x_3 - y_3)^2)}$ the observer-source distance with convection correction where $\beta = \sqrt{1 - M_{\text{flight}}^2}$ is the Prandtl-Glauert factor.

In practice, the quadrupole contribution can be neglected when $MS_t < 1$, where $S_t$ is the Strouhal number [26]. Doing so, the formulation for calculating far-field acoustic pressure is reduced to two surface integrals, the first, $p'_F$, refers to far-field loading noise and the second, $p'_Q$, refers to far-field thickness noise:

$$H(f) p'(\mathbf{x}, t) \approx - \underbrace{\int_{-\infty}^{t} \int_{\Sigma(\tau)} \frac{\partial}{\partial x_i} \left[ F_{ij}(\mathbf{y}, \tau) n_j \right] G_c(\mathbf{x}, t|\mathbf{y}, \tau) d\Sigma(\mathbf{y}) d\tau +}_{p'_F} \underbrace{\int_{-\infty}^{t} \int_{\Sigma(\tau)} \frac{D}{Dt} \left[ Q_j(\mathbf{y}, \tau) n_j \right] G_c(\mathbf{x}, t|\mathbf{y}, \tau) d\Sigma(\mathbf{y}) d\tau}_{p'_Q} \quad (15)$$

The variable $\mathbf{x}$ refers to the observer position, $\mathbf{y}$ to the source position, $t$ to the reception time and $\tau$ to the source emission time. $\Sigma$ is the control surface which depends on the time of emission $\tau$ since it is subjected to any motion, which in the case of the rotor is a rotation at angular speed $\Omega$. Using convolving product properties, the spatial and convected derivative operator are applied to $G_c$ leading to:

$$H(f) p'(\mathbf{x}, t)|_{f>0} \approx - \int_{-\infty}^{t} \int_{\Sigma(\tau)} \left[ F_{ij}(\mathbf{y}, \tau) n_j \right] \frac{\partial G_c(\mathbf{x}, t|\mathbf{y}, \tau)}{\partial x_i} d\Sigma(\mathbf{y}) d\tau + \int_{-\infty}^{t} \int_{\Sigma(\tau)} \left[ Q_j(\mathbf{y}, \tau) n_j \right] \frac{D G_c(\mathbf{x}, t|\mathbf{y}, \tau)}{Dt} d\Sigma(\mathbf{y}) d\tau \quad (16)$$

In open rotor applications, the flow typically has a high Reynolds number, so the viscous term present in $F_{ij}$ can be neglected [27], which yields:

$$F_{ij} \approx p' \delta_{ij} + \rho u_i (u_j - u_j^\Sigma) \quad (17)$$

where $p' = p - p_s^\infty$ with $p$ calculated in Eq. (9a) and $p_s^\infty$ the mean flow static pressure. The control surface $\Sigma$ is the surface where the fields of the conservative variables are extracted after performing a CFD simulation so that an acoustic analogy can be applied. $\Sigma$ can be of two kinds: either it coincides with the object walls and is called an on-body surface, or it does not and is referred to as an off-body surface.

In the case of an on-body surface $\Sigma_{on}$ that coincides with the blade skin, the monopole and dipole source terms are reduced to the volume flow expelled by the blade and the pressure fluctuations that apply to it:

$$F_{ij}^{\text{on}} = p' \delta_{ij} \quad (18)$$
$$Q_j^{\text{on}} = \rho_{\infty} u_j^{\Sigma} \qquad (19)$$

From the formulation of Eq. (16) fed with the source terms of Eqs. (18) and (19), Hanson and Parzych (HP) established a simplified frequency formulation valid in the far-field [22]. The transition to the frequency domain is initially made by applying a Fourier transform to the variable $t$. Since $\Sigma$ depends on $\tau$, this complicates the integration. For this reason, a transition is made to the rotating frame of reference, where the surface is steady and the source terms are assumed to be known. The frequency-convected free space Green's function is replaced by its far-field approximation. Finally, HP developments lead to the far-field expressions of rotor loading and thickness noise at the blade passing frequencies (BPFs) given by equations (56) and (48) of reference [22]. With the current paper's notations, it writes, for the far-field loading noise at $m^{\text{th}}$ BPF:

$$\begin{aligned} \widehat{p}'_F(\mathbf{x}, mB\Omega) = & -\frac{iB\kappa_{\infty}e^{i\kappa_{\infty}\bar{r}_{\beta}}}{4\pi\bar{r}_{\beta}}e^{imB(\arctan(\frac{x_3}{x_2}) - \frac{\pi}{2})} \left[ \begin{aligned} & S_c \int_{\Sigma_{\text{on}}} J_{mB}(\kappa_{\infty}r_0S_s)e^{-imB\phi_0}e^{-i\kappa_{\infty}S_c y_1}F_1^{\text{on}}(\tilde{\mathbf{y}})d\Sigma_{\text{on}}(\tilde{\mathbf{y}}) + \\ & iS_s \int_{\Sigma_{\text{on}}} J'_{mB}(\kappa_{\infty}r_0S_s)e^{-imB\phi_0}e^{-i\kappa_{\infty}S_c y_1}F_2^{\text{on}}(\tilde{\mathbf{y}})d\Sigma_{\text{on}}(\tilde{\mathbf{y}}) + \\ & \int_{\Sigma_{\text{on}}} J_{mB}(\kappa_{\infty}r_0S_s)\frac{e^{-imB\phi_0}e^{-i\kappa_{\infty}S_c y_1}}{M_r}F_3^{\text{on}}(\tilde{\mathbf{y}})d\Sigma_{\text{on}}(\tilde{\mathbf{y}}) \end{aligned} \right] \end{aligned} \qquad (20)$$

and for the far-field thickness noise at $m^{\text{th}}$ BPF:

$$\begin{aligned} \widehat{p}'_Q(\mathbf{x}, mB\Omega) = & -\frac{iB\kappa_{\infty}c_{\infty}e^{i\kappa_{\infty}\bar{r}_{\beta}^{+}}}{4\pi\beta^2\bar{r}_{\beta}} \left(1 + M_{\text{flight}}\frac{x_1}{\bar{r}_{\beta}}\right) e^{imB\left(\arctan\left(\frac{x_3}{x_2}\right) - \frac{\pi}{2}\right)} \\ & \times \int_{\Sigma_{\text{on}}} J_{mB}(\kappa_{\infty}r_0S_s y_1)e^{-imB\phi_0}e^{-i\kappa_{\infty}S_c y_1}Q^{\text{on}}(\tilde{\mathbf{y}})d\Sigma_{\text{on}}(\tilde{\mathbf{y}}) \end{aligned} \qquad (21)$$

where $B$ is the number of rotor blades and $\kappa_{\infty} = mB\Omega/c_{\infty}$ denotes the acoustic wavenumber. The observer's position is specified by Cartesian coordinates $x_i$, while the sources are described in cylindrical coordinates $\tilde{y}_i$, with $i \in \{1, 2, 3\}$ and $i = 1$ corresponding to the rotation axis. The source location in cylindrical coordinates is given by $\tilde{\mathbf{y}} = (\tilde{y}_1, r_0, \phi_0)$, where $\tilde{y}_1$ is the axial coordinate, $r_0$ is the rotation radius, and $\phi_0$ is the azimuthal angle of the source. The axial $F_1^{\text{on}}(\tilde{\mathbf{y}})$, radial $F_2^{\text{on}}(\tilde{\mathbf{y}})$, and tangential $F_3^{\text{on}}(\tilde{\mathbf{y}})$ components of the loading source term as well as the thickness source term $Q^{\text{on}}(\tilde{\mathbf{y}})$ for the use of an on-body surface, are all represented in a normalized polar basis, reflecting their dependence on the cylindrical coordinates $\tilde{\mathbf{y}}$ of a source point $y$ on the surface:

$$\left( \begin{matrix} F_1(\tilde{\mathbf{y}}) \\ F_2(\tilde{\mathbf{y}}) \\ F_3(\tilde{\mathbf{y}}) \end{matrix} \right) = \underbrace{\left( \begin{matrix} 1 & 0 & 0 \\ 0 & \cos \phi_0 & \sin \phi_0 \\ 0 & -\sin \phi_0 & \cos \phi_0 \end{matrix} \right)}_{R(\phi_0)} \left( \begin{matrix} F_1(\mathbf{y}) \\ F_2(\mathbf{y}) \\ F_3(\mathbf{y}) \end{matrix} \right) \qquad (22)$$

where $R(\phi_0)$ is the rotation matrix that transforms coordinates from cylindrical to Cartesian. The functions $J_{mB}$ are Bessel functions of order $mB$ and $J'_{mB}$ denotes their derivatives. The rotational Mach number of the source is defined as $M_r = \kappa_{\infty}r_0$. The terms $\bar{r}_{\beta}, \bar{r}_{\beta}^{\dagger}, S_c$ and $S_s$ depend on the observer's coordinates $(x_1, x_2, x_3)$ as defined in [21]. The total far-field acoustic pressure in the frequency domain is nothing more than the sum of Eqs. (20) and (21). Throughout the paper, the latter will be referred to as the HP formulation.

In the case of an off-body surface $\Sigma_{\text{off}}$, a similar procedure can be applied with the differences hidden in the source terms $F_{i,j}$ and $Q_j$ given by Eqs. (17) and (13c) respectively. This leads to similar expressions for the far-field loading noise in the frequency domain:

$$\begin{aligned} \widehat{p}'_F(\mathbf{x}, mB\Omega) = & -\frac{iB\kappa_{\infty}e^{i\kappa_{\infty}\bar{r}_{\beta}}}{4\pi\bar{r}_{\beta}}e^{imB(\arctan(\frac{x_3}{x_2}) - \frac{\pi}{2})} \left[ \begin{aligned} & S_c \int_{\Sigma_{\text{off}}} J_{mB}(\kappa_{\infty}r_0S_s)e^{-imB\phi_0}e^{-i\kappa_{\infty}S_c y_1}F_1(\tilde{\mathbf{y}})d\Sigma_{\text{off}}(\tilde{\mathbf{y}}) + \\ & iS_s \int_{\Sigma_{\text{off}}} J'_{mB}(\kappa_{\infty}r_0S_s)e^{-imB\phi_0}e^{-i\kappa_{\infty}S_c y_1}F_2(\tilde{\mathbf{y}})d\Sigma_{\text{off}}(\tilde{\mathbf{y}}) + \\ & \int_{\Sigma_{\text{off}}} J_{mB}(\kappa_{\infty}r_0S_s)\frac{e^{-imB\phi_0}e^{-i\kappa_{\infty}S_c y_1}}{M_r}F_3(\tilde{\mathbf{y}})d\Sigma_{\text{off}}(\tilde{\mathbf{y}}) \end{aligned} \right] \end{aligned} \qquad (23)$$
and for the far-field thickness noise in the frequency domain:

$$ \hat{p}'_Q(\mathbf{x}, mB\Omega) = -\frac{iB\kappa_\infty c_\infty e^{i\kappa_\infty \bar{r}_\beta}}{4\pi\beta^2 \bar{r}_\beta} \left(1 + M_{\text{flight}} \frac{x_1}{\bar{r}_\beta}\right) e^{imB\left(\arctan\left(\frac{x_3}{x_2}\right) - \frac{\pi}{2}\right)} \\ \times \int_{\Sigma_{\text{off}}} J_{mB}(\kappa_\infty r_0 S_s y_1) e^{-imB\phi_0} e^{-i\kappa_\infty S_c y_1} Q(\tilde{\mathbf{y}}) d\Sigma_{\text{off}}(\tilde{\mathbf{y}}) \quad (24) $$

The sum of Eqs. (23) and (24) will be referred to as the HP formulation variant used in this paper to capture shock noise.

### C. Evaluation of rotor self-noise gradient

In order to assemble the gradient of the objective function with respect to the shape parameters given by Eq. (6), it is essential to know the sensitivities of the functional with respect to flow state variables $\partial \mathcal{F}/\partial \mathbf{W}$ and with respect to coordinates $\partial \mathcal{F}/\partial \mathbf{X}$.

#### 1. Sensitivities with respect to flow state variables

The sensitivities with respect to flow state variables are obtained here by analytically differentiating the variant of the Hanson and Parzych formulation in Eqs. (23) and (24), so that any functional expressed in terms of $\hat{p}'$ can be obtained by applying the chain rule. The sum of Eqs. (23) and (24) can be re-written as:

$$ \hat{p}'(\mathbf{x}, mB\Omega) = \xi_F \int_{\Sigma_{\text{off}}} I_i F_i(\tilde{\mathbf{y}}) \, d\Sigma_{\text{off}}(\tilde{\mathbf{y}}) + \int_{\Sigma_{\text{off}}} I Q(\tilde{\mathbf{y}}) \, d\Sigma_{\text{off}}(\tilde{\mathbf{y}}) \quad (25) $$

where the constant $\xi_F$ given by:

$$ \xi_F = -\frac{iB\kappa_\infty e^{i\kappa_\infty \bar{r}_\beta}}{4\pi\bar{r}_\beta} e^{imB\left(\arctan\left(\frac{x_3}{x_2}\right) - \frac{\pi}{2}\right)} \quad (26) $$

only depends on the observer location, and the terms $I$ and $I_i$ write:

$$ I_1 = S_c J_{mB}(\kappa_\infty r_0 S_s) e^{-imB\phi_0} e^{-i\kappa_\infty S_c y_1} \quad (27) $$

$$ I_2 = i S_s J'_{mB}(\kappa_\infty r_0 S_s) e^{-imB\phi_0} e^{-i\kappa_\infty S_c y_1} \quad (28) $$

$$ I_3 = J_{mB}(\kappa_\infty r_0 S_s) \frac{e^{-imB\phi_0} e^{-i\kappa_\infty S_c y_1}}{M_r} \quad (29) $$

$$ I = -\frac{iB\kappa_\infty c_\infty e^{i\kappa_\infty \bar{r}_\beta}}{4\pi\beta^2 \bar{r}_\beta} \left(1 + M_{\text{flight}} \frac{x_1}{\bar{r}_\beta}\right) e^{imB\left(\arctan\left(\frac{x_3}{x_2}\right) - \frac{\pi}{2}\right)} e^{-imB\phi_0} e^{-i\kappa_\infty S_c y_1} \quad (30) $$

Equation (25) is expressed in the rotating frame of reference with $F_i$ and $Q$ expressed in polar coordinates as they depend on $\tilde{\mathbf{y}}$. As the differentiation of the source terms of Eq (25) is more suitable when expressed in the Cartesian reference frame, the transition matrix from cylindrical to Cartesian coordinates $R(\phi_0)$ is introduced so that Eq. (25) is re-written in Cartesian coordinates ($Q$ is a scalar):

$$ \hat{p}'(\mathbf{x}, mB\Omega) = \xi_F \int_{\Sigma_{\text{off}}} I_i R_{ij} F_j(\mathbf{y}) \, d\Sigma_{\text{off}}(\mathbf{y}) + \int_{\Sigma_{\text{off}}} I Q(\mathbf{y}) \, d\Sigma_{\text{off}}(\mathbf{y}) \quad (31) $$

where only source terms $F_j$ and $Q$ depend on flow variables. By noting $(\lambda, \mu)$ the structure indices of surface cells, integration can be numerically approximated following the midpoint rule and allowing to rewrite the expression with sums:

$$ \hat{p}'(\mathbf{x}, mB\Omega) \approx \sum_{\lambda} \sum_{\mu} \left[ \xi_F I_i R_{ij} F_{jl}(\mathbf{y}) + I Q_l(\mathbf{y}) \right]_{\lambda\mu} (\Delta\Sigma_{\text{off},l})_{\lambda\mu} \quad (32) $$

The general expression for the derivative of sound pressure in the cell of indices $(\lambda, \mu)$ with respect to a flow variable denoted $W_{k,\lambda\mu}$ is written as:

$$ \frac{\partial \hat{p}'(\mathbf{x}, mB\Omega)}{\partial W_{k,\lambda\mu}} = \left[ \xi_F I_i R_{ij} \frac{\partial F_{jl}}{\partial W_k} + I \frac{\partial Q_l}{\partial W_k} \right]_{\lambda\mu} (\Delta\Sigma_{\text{off},l})_{\lambda\mu} \quad (33) $$
where, in the case of an off-body surface as acoustic source, the source terms are expressed as in Eqs. 17 and 13c. After some calculations, one can obtain the final sensitivities of monopolar and dipolar terms with respect to flow variables as:

$$ \frac{\partial Q_l}{\partial \rho} = \frac{\rho u_l}{\rho} - u_l^{\Sigma} \quad (34) $$
$$ \frac{\partial F_{jl}}{\partial \rho} = \frac{1}{\rho^2} \left[ \frac{(\gamma - 1)}{2} \rho u_k \rho u_k \delta_{jl} - \rho u_j \rho u_l \right] \quad (37) $$

$$ \frac{\partial Q_l}{\partial \rho u_k} = \delta_{lk} \quad (35) $$
$$ \frac{\partial F_{jl}}{\partial \rho u_k} = \frac{1}{\rho} \left[ -(\gamma - 1) \rho u_k \delta_{jl} + (\rho u_l - \rho u_l^{\Sigma}) \delta_{jk} + \rho u_j \delta_{lk} \right] \quad (38) $$

$$ \frac{\partial Q_l}{\partial \rho E} = 0 \quad (36) $$
$$ \frac{\partial F_{jl}}{\partial \rho E} = (\gamma - 1) \delta_{jl} \quad (39) $$

where as before $\rho c_{\infty} M_{\text{flight}}$ is added to $\rho u$ since the CFD simulation is performed in the translational relative frame of reference.

## 2. Sensitivities with respect to mesh coordinates

A similar procedure can be applied to evaluate the sensitivities with respect to coordinates, this time taking care to derive all terms since they all depend on the source coordinates y and paying attention to node/cell data conversion. The calculation steps and results are not presented here to keep conciseness.

# III. Preliminary assessment on a realistic open rotor blade

In this part, the entire discrete adjoint framework described in the previous section is applied to the calculation of the sound power level gradient of a high-speed open rotor, where the design variable is a Gaussian deformation parameter.

## A. Description of the configuration

The chosen configuration is the realistic rotor of an open fan configuration, referred to as Rotor ONERA Supercritical Experiment (ROSE), which is detailed by Dumont in [28] and illustrated in Fig. 3. This thirteen-bladed open rotor is considered in an isolated configuration and operating at cruise conditions. The geometrical characteristics and physical parameters of operation are described in Tab. 1.

<figure><img src="image_4.png" /><figcaption>Fig. 3 ROSE geometry</figcaption></figure>

<table><caption>Table 1 Physical parameters of the considered open rotor: number of blades, blade radius, flight Mach number, angular velocity, static temperature, static pressure and viscosity ratio</caption><thead><tr><th>B</th><th>R<sub>blade</sub></th><th>M</th><th>Ω</th><th>T<sub>s</sub><sup>∞</sup></th><th>p<sub>s</sub><sup>∞</sup></th><th>μ<sub>t</sub>/μ</th></tr></thead><tbody><tr><td>13</td><td>0.8 m</td><td>0.75</td><td>850 RPM</td><td>218.808 K</td><td>23842.3 Pa</td><td>0.1</td></tr></tbody></table>

The considered acoustic objective function for the assessment of the framework is the acoustic power level (PWL) obtained by integrating the sound intensity over the surface of a sphere $S$, which can be written in spherical coordinates as:

$$ \mathcal{F} = \text{PWL} = 10 \log_{10} \left( \int_S \frac{\hat{p}' \hat{p}'}{2\rho_{\infty} c_{\infty}} F_M(\psi, M_{\text{flight}}) dS \right) - 10 \log_{10} (1.10^{-12}) \quad (40) $$

with $\rho_{\infty}$ the mean flow density, $\hat{p}'$ the frequency-domain acoustic pressure and $\overline{\hat{p}'}$ its complex conjugate. $F_M$ is a modulation function that depends on the polar angle $\psi$ and the flight Mach number. It can be obtained by simplifying the formula for sound intensity given by Cantrell and Hart [29] (as done by Blandeau [30] for example):

$$ F_M(\psi, M) = \sqrt{1 - M^2 \sin^2 \psi} \left( \sqrt{1 - M^2 \sin^2 \psi} + M \cos \psi \right)^2 \quad (41) $$
## B. Numerical setup for forward evaluation

The forward simulation is carried out using a structured multi-block Navier-Stokes solver available in the simulation tool elsA [31] where governing equations are discretized in space using a cell-centered finite-volume approach. The steady-state RANS equations are discretized and solved so that their residual noted $R$ (cf. Eq. (10)) converges to zero. The equations are closed using the Spalart-Allmaras one-equation turbulence model. Convective terms are discretized using a Roe upwind scheme with Van Albada limiting function and Harten entropy correction. Viscous terms are discretized using a centered scheme. Time integration is performed using a backward Euler scheme. The Courant Friedrichs Lewy number (CFL) increases linearly from 0.5 to 10 in the first 1000 iterations.

The simulation is carried out considering an angular periodicity of $2\pi/B$. Material surfaces as the hub, the blade and the spinner are treated as no-slip surfaces, while a wall-slip condition is applied upstream of the spinner that is extended with a small radius in order to ease the meshing process. The azimuthal boundaries are treated as periodic while the rest of external boundaries are defined as far-field conditions. The boundary conditions of the problem are summarized in Figs. 4a and 4b.

<figure><img src="image_5.png" /><figcaption>Fig. 4 Boundary conditions</figcaption></figure>

The CFD mesh is a structured grid of nearly 25 millions cells, including 71000 cells over the blade surface. At the walls, the first layer has a height such that $y^+$ is unity. In the remaining part of the domain, the grid is stretched with an expansion ratio strictly less than 1.2. Figs. 5a and 5b show the volume grid of the fluid domain as well as the blade surface mesh.

<figure><img src="image_6.png" /><figcaption>Fig. 5 3D view of volume mesh and blade surface mesh</figcaption></figure>
Since the accuracy of the adjoint calculation is very sensitive to the convergence of the direct calculation [7], its convergence is checked in Figs. 6a and 6b where the evolution of the normalized residuals of flow and turbulent variables and that of normalized aerodynamic quantities are given. Fig. 6a demonstrates the convergence of the direct CFD calculation. Indeed, for each conservative variable, the associated normalized residual converges to a low value, losing around 3 orders of magnitude between the 1st and the 100000th iteration. Conversely, the normalized residual associated with the Spalart-Allmaras equation increases until it reaches a convergence plateau at 10². In fact, the initial guess for the turbulent quantity $v_t$ is far from the converged solution, and as the mean flow field adjusts, the turbulence model equation can experience rapid changes, leading to a temporary growth in the associated residual.

As shown in Fig. 6b, monitoring of integrated aerodynamic quantities reveals convergence of the normalized drag, lift and torque forces associated with blade rotation, which all reach a plateau after 5000 iterations.

<figure><img src="image_7.png" /><figcaption>Fig. 6 Validation of forward CFD simulation convergence</figcaption></figure>

### C. Validation of the HP formulation variant

At the end of the CFD calculation, the aerodynamic field is extracted on an off-body surface and feeds the monopolar and dipolar source terms of the integral formulation, enabling acoustic radiation to be performed. As explained in Sec. III.B, the solver adopts a cell-centered approach for which the CFD data is stored at the centre of the cells. To improve the accuracy of sensitivity calculations, the off-body surface was defined by extracting a single layer of cells from the fluid mesh and performing interpolation on nodes, rather than extracting two cell layers and then directly interpolating cell-center data which could compromise the accuracy of the off-body surface data.

Therefore, in order to ensure better accuracy of the sensitivity calculation, the off-body surface has been defined by extracting a layer of cells from the fluid mesh and performing interpolation on the nodes avoiding extraction of two layers of cells and direct interpolation on data located at cell-centres which could affect the accuracy of the resulting data on the off-body surface.

Since the aim is to capture the shock noise correctly, the ideal off-body surface must be chosen so that it lies in the linear region where the acoustic power is conserved. As the upstream and downstream axial positions of the surface have no influence on the radiated noise, these have been set so that the surface is positioned 0.5 m upstream and 0.5 m downstream of the rotor. In addition, the surface is only partially closed downstream of the rotor in order to avoid to be crossed by the blade tip vortices that can create spurious noise due to source term truncation errors [32]. To define this ideal surface, the acoustic power is computed for different off-body surfaces characterized by the radial distance between their top part and the blade tip as illustrated in Fig. 7a. To this end, one hundred observers are uniformly positioned in the far-field on a semi-circle with centre (0, 0, 0) and radius $r_{obs}$ such as $\kappa_{\infty}r_{obs} \approx 30$, as illustrated in Fig. 7b.
<figure><img src="image_8.png" /><figcaption>Fig. 7 Acoustic post-processing for power level calculation</figcaption></figure>

Acoustic pressure is obtained on the arc from Eqs. 23 and 24 at the blade passing frequency (1st BPF) and its first (2nd BPF) and second (3rd BPF) harmonics. Using Eq. 40, acoustic power levels are obtained for these frequencies by integrating the sound intensity over a sphere obtained by revolution of the semicircle for each off-body surface, leading to the results of Fig. 8a. It can be noticed that it is important not to place the surface too close from the rotor (radial distance to the blade tip less than 0.5 m) in the non-linear region, or too far away where the mesh becomes coarser (radial distance to blade tip greater than 1.0 m). The results show that the ideal off-body surface, for which linear region is reached, is placed at a radial distance of 1.0 m to the blade tip. The pressure field over this surface is shown in Fig. 8b for illustration.

<figure><img src="image_9.png" /><figcaption>Fig. 8 Ideal surface study</figcaption></figure>

In order to validate the proposed variant HP formulation, the sound pressure levels due to the loading and to the thickness noise sources are evaluated separately. As a first step, the blade surface is chosen as the integration surface in Eq. (16) and the directivity diagram of loading and thickness noise at 1st BPF obtained with the classic Hanson and Parzych formulation (Eqs. (20) and (21)) is compared in Fig. 9a with that given by a well-validated time-domain FW-H solver available at ONERA (Eq. 2 with $u_n = v_n$ from [33]), denoted by FW-H in the figure. The results obtained with the variant of HP formulation at 1st BPF for off-body surfaces (Eqs. (23) and (24)) are also compared in Fig. 9b with that provided by a similar calculation with Eq. 3 from [34]. Figs. 9a and 9b show an excellent agreement between the
far-field sound pressure level predicted by both formulations and with both types of control surface. The thickness noise generated by the displacement of air due to the blade's movement is dominant in the upstream direction, while the loading noise, generated by pressure fluctuations over the blades, is more pronounced in the downstream direction. Fig. 9b clearly shows the contribution of shock noise (difference between on-body and off-body surfaces) that radiates mainly upstream of the rotor. Some discrepancies that can be attributed to numerical errors in the integration techniques of the time-domain FW-H formulation can be noticed close to the axis of rotation (0 – 30° and 150 – 180°), where the absolute levels are expected to be very low.

<figure><img src="image_10.png" /><figcaption>Fig. 9 Validation of the variant of HP formulation at 1st BPF</figcaption></figure>

In order to illustrate the sound radiation pattern, the fields of the sound pressure real part as well as the sound pressure level obtained via the variant of the HP formulation using the ideal off-body surface of Fig. 8b are plotted in Figs. 10a and 10b for 1st BPF.

In Fig. 10a, the sound pressure fluctuation real part for 1st BPF is visualized on an axial plane and an horizontal plane. In the axial plane, 13 lobes, characteristic of the rotor-locked mode ($B = 13$), can be observed. In the plane of rotation of the rotor, the spiral shape of the modes reveals the rotating nature of the acoustic field in the direction of rotation of the rotor blade. This spiral shape of the modes can be explained by the fact that the analytical acoustic field in the variant of the HP formulation in Eqs. (23) and (24) is written as a linear combination of Bessel functions of the first kind, which are representative of a rotating source. Moreover, the amplitude of the real part of the acoustic pressure fluctuations decreases in $1/r_β$.
<figure><img src="image_11.png" /><figcaption>Fig. 10 Pressure field in case of off-body surface at 1st BPF</figcaption></figure>

Fig. 10b shows that most of the noise is radiated in the plane of rotation of the rotor and at 90 ± 30° in the horizontal plane where the sound pressure level reaches nearly 170 dB in the direct vicinity of the rotor. These preferred directions for noise propagation provide important information for optimizing the shape of the rotor in order to target these directions in priority.

## D. Validation of the acoustic metric sensitivities

Sensitivities of the sound power level (Eq. (40)) with respect to the flow variables W and the mesh X are obtained analytically using the analytic sensitivities of acoustic pressure developed in Sec. II.C, and by 2nd order finite differences (FD) for validation purposes. On-blade sensitivities of the sound power level with respect to density using both methods are shown in Fig. 11.

<figure><img src="image_12.png" /><figcaption>Fig. 11 Sensitivity of acoustic power level with respect to density on pressure side (PS) 11a, 11b and on suction side (SS) 11c, 11d</figcaption></figure>

Fig. 11 illustrates the importance of visualizing sensitivities to determine areas of importance for blade parametrization prior to optimization. In the interest of noise reduction and prior to calculation, it is expected to unload the blade on the pressure side and load it on the suction side. Fig. 11 shows that the absolute levels of sensitivities are higher close to the blade tip and near leading and trailing edges, which is expected since the velocity is highest there. Also, the development of the shock on the blade explains the discontinuity observed in the sensitivities since between the leading edge and the trailing edge, the density is increased, making the blade more loaded. Unsurprisingly, on Figs. 11c and 11a, the sensitivity values on the pressure side are positive close to the blade tip near leading edge and negative close to the blade tip near trailing edge. This suggests that an increase in density near trailing edge and a decrease of it near leading
edge on the pressure side would reduce the far-field acoustic power level and hence the radiated noise. The opposite observation can be made for the suction side. Figs. 11d and 11b only confirm and validate the analytical derivation method of the HP formulation variant to obtain the sound pressure sensitivities. Indeed, on both figures, the absolute deviation between the analytical derivation and the finite-difference approach averages $10^{-9}$ over the entire blade, with a few localized areas where the deviation is of the order of 15. The limits of computer rounding have thus been reached, and the numerical error predominates, hence the noisy effect.

On-blade sensitivities of the sound power level with respect to the axial coordinate using both methods and their deviation are now shown in Figs. 12c, 12a, 12d and 12b:

<figure><img src="image_13.png" /><figcaption>Fig. 12 Sensitivity of acoustic power level with respect to axial coordinate on pressure side (PS) 12a, 12b and suction side (SS) 12c, 12d</figcaption></figure>

Sensitivities of sound power levels to axial coordinates behave quite differently from sensitivities to flow variables. From Figs 12c and 12a, it can be noticed that the error in mesh sensitivities is larger at the edges of the zones which is an artifact that is counterbalanced by an opposite sign contribution in the neighboring zones. In both case, the logarithmic errors are very small, and sometimes even converge towards zero machine.

Similar analysis can be done for off-blade sensitivities of the sound power level. Only the results for the sensitivities with respect to the density and to the axial coordinate using both analytical derivation and finite differences along with their deviation are shown in Figs. 13a, 13b, 13c and 13d. Again, the results validate the analytical derivation of the sensitivities.

<figure><img src="image_14.png" /><figcaption>Fig. 13 Sensitivity of acoustic power level with respect to density 13a, 13b and axial coordinate at the upstream of the off-body surface 13c, 13d</figcaption></figure>

## E. Noise discrete gradient validation in case of an analytical mesh deformation

Now that the sound power sensitivities have been validated, the whole workflow is evaluated here, this time by comparing the final gradient $dF/d\alpha$ obtained after adjoint calculation with that obtained by finite differences. For the
purpose of validation, an analytical mesh deformation is performed with the following Gaussian form:

$$X_2(\alpha) = X_{2,0} + \alpha \exp\left(-\frac{(\mathbf{X} - \mathbf{X}_0)^2}{2\sigma^2}\right) \qquad (42)$$

where $X_2$ is the second coordinate of the surface mesh in the direction transverse to the blade surface, $X_{2,0}$ its initial value, $\sigma = 0.01$ is the standard deviation monitoring the bump width, $\alpha$ is the design parameter controlling the Gaussian height and $\mathbf{X}_0 = (-1.25, 0, 1.5)$ is the point of application of the Gaussian deformation. In this case, the analytical derivative of the mesh with respect to the Gaussian-shaped parameter $\alpha$ is easily established in Eq. (43):

$$\frac{dX_2(\alpha)}{d\alpha} = \exp\left(-\frac{(\mathbf{X} - \mathbf{X}_0)^2}{2\sigma^2}\right) \qquad (43)$$

On the other hand, acoustic power sensitivities with respect to flow variables and mesh coordinates $\partial\mathcal{F}/\partial\mathbf{W}$ and $\partial\mathcal{F}/\partial\mathbf{X}$ obtained in the previous section serve as inputs to a discrete adjoint RANS solver which provides the gradient of the functional with respect to shape parameters. The system of adjoint equations is solved using a Flexible Generalized Minimal Residual (FGMRES) where GMRES itself is used as preconditioner. The main numerical parameters of this adjoint simulation are summarized in Tab. 2.

<table><caption>**Table 2** GMRES method parameters for adjoint computation</caption><tr><td>Maximum dimension of the outer Krylov subspace</td><td>150</td></tr><tr><td>preconditioner</td><td>GMRES</td></tr><tr><td>Maximum dimension of the inner Krylov subspace</td><td>15</td></tr><tr><td>Inner residual decrease ratio</td><td>0.6</td></tr><tr><td>GMRES tolerance</td><td>10<sup>-3</sup></td></tr><tr><td>Tolerance for RHS deflation in block GMRES</td><td>10<sup>-12</sup></td></tr><tr><td>Number of deflated vectors</td><td>60</td></tr></table>

Due to a lack of convergence of the right-hand member residual for the adjoint calculation, a coarser mesh, adapted to the 1st BPF studied here, was used for gradient validation. Fig. 14 plots the convergence history of GMRES for the adjoint equation, showing that it takes almost 350 iterations to reach the GMRES tolerance set to $10^{-3}$. The value of the gradient obtained using an on-body surface is about $-1.54 \times 10^{-2}\text{dB.m}^{-1}$, a higher value than when using an off-body surface where the gradient is $-6.81 \times 10^{-3}\text{dB.m}^{-1}$.
<figure><img src="image_15.png" /><figcaption>Fig. 14 Normalized residual of the GMRES Right Hand Side for the simulation using off-body surface</figcaption></figure>

For validation purposes, the power level gradient is now computed using finite differences. To this end, the objective function of Eq. (40), that depends implicitly on the shape parameters via the flow variable and the mesh, and can be written as:

$$ \mathcal{F}(\alpha) = F(\mathbf{W}(\alpha), \mathbf{X}(\alpha)) \qquad (44) $$

where F is a function symbolizing the dependency of F on W and X. The order-2 central finite difference approximate evaluation of its gradient is defined by:

$$ \frac{d\mathcal{F}}{d\alpha} \approx \frac{1}{2\Delta\alpha} \left[ F(\mathbf{W}(\alpha + \Delta\alpha), \mathbf{X}(\alpha + \Delta\alpha)) - F(\mathbf{W}(\alpha - \Delta\alpha), \mathbf{X}(\alpha - \Delta\alpha)) \right] \qquad (45) $$

A small ±Δα disturbance is therefore applied locally to the blade as shown in Figs. 15a and 15b. Forward simulations are then performed for the positive and negative deformations in order to compute the resulting power levels, and finally its gradient.

<figure><img src="image_16.png" /><figcaption>Fig. 15 Local Gaussian blade deformation for Δα = ±10⁻²</figcaption></figure>
The procedure is applied for several $\Delta\alpha$ values and the results are compared against the adjoint evaluation of the gradient in Fig. 16.

<figure><img src="image_17.png" /><figcaption>Fig. 16 Relative deviation (finite differences - adjoint) on the acoustic power level gradient for local Gaussian deformation</figcaption></figure>

In Fig. 16, the relative deviation between 2-order finite difference and discrete adjoint evaluation of the gradient is evaluated with both the on-body surface HP formulation and its extension for an off-body surface. For very large deformations, such as when ($\Delta\alpha = 10^{-1}$ and $\Delta\alpha = 1$, the resulting deviations reached extremely high values (on the order of $10^5$), which were omitted from the graph to maintain clarity. This occurs because the step size selected for the finite-difference method is too large, resulting in the function of interest changing too much over the interval $[-\Delta\alpha, \Delta\alpha]$ for an accurate gradient approximation. Specifically, the finite-difference formula is based on a truncated Taylor series expansion. When the step size is large, the higher-order terms that are omitted become important, causing the approximation to differ from the actual derivative. Similarly, a very small step size (lower than $\Delta\alpha = 10^{-8}$) leads to round-off error (not represented in Fig. 16 because subtracting $F(\alpha - \Delta\alpha)$ and $F(\alpha + \Delta\alpha)$ nearly equal numbers results in loss of significant digits.

For values of $\alpha \in [10^{-8}, 10^{-3}]$, the relative deviation is around 6.2% for the use of an on-body control surface and 6.0% for the use of an off-body surface in the far-field sound power evaluation process. One reason for the smaller deviation obtained with the off-body configuration could be that the control surface does not correspond to the surface where the deformation is applied. Consequently, the effects of the blade deformation would be attenuated when propagating to the fictitious surface in the CFD calculation. The extracted aerodynamic field would be much less perturbed than in the on-body case, where the surface coincides with the blade, which would have an impact on the evaluation of the function of interest and the calculation of its gradient.

## IV. Conclusion

This study has successfully established a discrete adjoint methodology for the evaluation of the rotor tonal noise gradient using a far-field off-body FW-H formulation in the rotating domain. Adapting HP formulation to the use of an off-body surface, where the source terms differ from the use of an on-body surface, enables innovative open-rotor architectures to be addressed with a view to assessing and reducing shock noise, the dominant noise source for these engines at cruising speeds.

The proposed approach accurately predicts far-field noise levels, as validated by close agreement with results from a well-established acoustic solver. Analytical differentiation of the variant of HP integral formulation allowed calculation
of sensitivities for any acoustic objective function depending on far-field acoustic pressure, establishing a robust basis for gradient-based optimization. These analytical results were then thoroughly compared with those obtained via a second-order centered finite-difference method, showing strong agreement and thereby confirming the reliability of the analytical approach. Furthermore, the integrity and robustness of the entire computational chain were rigorously verified by comparing the final gradients from the CFD adjoint calculations with those from the finite-difference method, and a satisfactory agreement is obtained. Overall, the developed methodology not only offers a reliable and efficient tool for rotor noise optimization but also establishes a solid framework for future research on aeroacoustic optimization in complex rotating systems. This work paves the way for further advancements in reducing tonal noise emissions, contributing to quieter and more environmentally friendly rotorcraft designs.

The advances presented in this paper provide a sound basis for extending the use of Hanson and Parzych's formulation to more complex and realistic operational regimes, particularly those involving transonic blade tip motion. By addressing these challenging conditions, the methodology can be effectively leveraged for aeroacoustic shape optimization, targeting the reduction of shock-induced noise, a critical concern in modern turbomachinery and propeller design. A multi-point optimization strategy would be explored, enabling the identification of optimal design compromises that balance the often-conflicting acoustic requirements of both take-off and cruise flight conditions. This approach ensures that noise reduction benefits are realized across a broader operational envelope, ultimately contributing to quieter and more efficient rotorcraft operation. Moreover, with the aim of accounting in the future for the rotor's unsteady loading noise-which is not considered in the present work-the integration of a time spectral method would offer a promising pathway for efficiently handling unsteady, periodic configurations. This approach circumvents the computational expense typically associated with time-domain adjoint methods, thereby enhancing the practicality and scalability of the discrete adjoint framework for real-world engineering applications.

## Acknowledgments

This work is supported by a thesis grant through the project CLEOPATRA (CLEaner future through open fan inTegRAtion) which has received funding from the French DGAC (Direction Générale de l'Aviation Civile).

## References

[1] Daydé-Thomas, V., Polacsek, C., Fauqueux, S., Gloerfelt, X., and Marjono, J., "Tone noise characterization of an open-fan engine using source-mode integral formulations," *30th AIAA/CEAS Aeroacoustics Conference (2024)*, 2024, p. 3152.

[2] Sugiyama, M., Falissard, F., Reboul, G., Jaouani, N., and Gloerfelt, X., "Noise prediction of an open fan using acoustic analogy approaches," *INTER-NOISE and NOISE-CON Congress and Conference Proceedings*, Vol. 268, Institute of Noise Control Engineering, 2023, pp. 2921-2932.

[3] Lewis, D., Perrin, C., and Falissard, F., "Tonal Noise Predictions of a Generic Open-Fan Based on a Hybrid URANS Approach," *30th AIAA/CEAS Aeroacoustics Conference (2024)*, 2024, p. 3264.

[4] Ffowcs Williams, J. E., and Hawkings, D. L., "Sound generation by turbulence and surfaces in arbitrary motion," *Philosophical Transactions of the Royal Society of London. Series A, Mathematical and Physical Sciences*, Vol. 264, No. 1151, 1969, pp. 321-342.

[5] Pacini, B., Yildirim, A., Davoudi, B., Martins, J. R., and Duraisamy, K., "Towards efficient aerodynamic and aeroacoustic optimization for urban air mobility vehicle design," *AIAA AVIATION 2021 FORUM*, 2021, p. 3026.

[6] Groom, M. J., Zhou, B. Y., and Lopes, L. V., "Aeroacoustic Optimization of Installed Propeller Configurations with Scattering and Shielding Effects," *30th AIAA/CEAS Aeroacoustics Conference (2024)*, 2024, p. 3039.

[7] Peter, J., "Contributions to discrete adjoint method in aerodynamics for shape optimization and goal-oriented mesh-adaptation," Ph.D. thesis, UNIVERSITE DE NANTES, 2020.

[8] Jameson, A., "Aerodynamic design via control theory," *Journal of scientific computing*, Vol. 3, 1988, pp. 233-260.

[9] "Calculs de gradients pour l'optimisation des performances aérodynamiques d'un rotor d'hélicoptère en vol stationnaire," Ph.D. thesis, Poitiers, 2010.

[10] Al-Garni, A. Z., Abdallah, A. M., et al., "Aerodynamic-Shape Optimization of Supersonic-Missiles Using Monte-Carlo," *International Review of Aerospace Engineering (I. RE. AS. E)*, 2008.
[11] Ghazlane, I., "Adjoint-based aerostructural sensitivity analysis for wing design," Ph.D. thesis, Université Nice Sophia Antipolis, 2012.

[12] Blondeau, C., and Liauzun, C., "A modular implementation of the time spectral method for aeroelastic analysis and optimization on structured meshes," *IFASD 2019*, 2019.

[13] Liauzun, C., and Blondeau, C., "A modular TSM solver for aeroelastic analysis and optimization," *IFASD 2024*, 2024.

[14] Zhou, B. Y., Albring, T. A., Gauger, N. R., Economon, T. D., Palacios, F., and Alonso, J. J., "A discrete adjoint framework for unsteady aerodynamic and aeroacoustic optimization," *16th AIAA/ISSMO multidisciplinary analysis and optimization conference*, 2015, p. 3355.

[15] Abergo, L., Morelli, M., Pullin, S. F., Zhou, B. Y., and Guardone, A., "Adjoint-Based Aeroacoustic Optimization of Propeller Blades in Rotating Reference Frame," *AIAA Aviation 2023 Forum*, 2023, p. 3836.

[16] Zhou, B. Y., Albring, T. A., Gauger, N. R., Economon, T. D., and Alonso, J. J., "An efficient unsteady aerodynamic and aeroacoustic design framework using discrete adjoint," *17th AIAA/ISSMO multidisciplinary analysis and optimization conference*, 2016, p. 3369.

[17] Zhou, B. Y., Albring, T., Gauger, N. R., da Silva, C. R., Economon, T. D., and Alonso, J. J., "Efficient airframe noise reduction framework via adjoint-based shape optimization," *AIAA Journal*, Vol. 59, No. 2, 2021, pp. 580-595.

[18] Abergo, L., Guardone, A., Pullin, S., and Zhou, B. Y., "Adjoint-Based Tonal and Broadband Aeroacoustic Optimization of Propeller Blades With Amiet Model," *30th AIAA/CEAS Aeroacoustics Conference (2024)*, 2024, p. 3130.

[19] Hajczak, A., and Buszyk, M., "Shape Sensitivities of 2D Airfoils for Broadband Noise Reduction Using the Adjoint Method and Semi-Analytical Models," *30th AIAA/CEAS Aeroacoustics Conference (2024)*, 2024, p. 3182.

[20] Chelius, A., El Din, I. S., and Koch, R., "Toward an Adjoint Based Aeroacoustic Optimisation for Propeller Noise Reduction," *ICSV26*, 2019.

[21] Foglia, E., Daroukh, M., Buszyk, M., Salah El-Din, I., and Moreau, S., "Towards the Optimization of Propeller Tonal Noise in the Frequency Domain Using the Discrete Adjoint Method," *30th AIAA/CEAS Aeroacoustics Conference (2024)*, 2024, p. 3153.

[22] Hanson, D. B., and Parzych, D. J., "Theory for noise of propellers in angular inflow with parametric studies and experimental verification," Tech. rep., 1993.

[23] Parente, E., Tarsia Morisco, C., Alauzet, F., and Spalart, P., "Operating the Spalart-Allmaras turbulence model under high free stream turbulence," *AIAA SCITECH 2024 Forum*, 2024, p. 2740.

[24] Spalart, P., and Allmaras, S., "A one-equation turbulence model for aerodynamic flows," *30th aerospace sciences meeting and exhibit*, 1992, p. 439.

[25] Goldstein, M., "Unified approach to aerodynamic sound generation in the presence of solid boundaries," *Acoustical Society of America*, 1974.

[26] Hajczak, A., Sanders, L., Vuillot, F., and Druault, P., "A comparison between off and on-body control surfaces for the FW-H equation: Application to a non-compact landing gear wheel," *Journal of Sound and Vibration*, Vol. 490, 2021, p. 115730.

[27] Glegg, S., and Devenport, W., *Aeroacoustics of low Mach number flows: fundamentals, analysis, and measurement*, Academic Press, 2017.

[28] Dumont, A., "Design of the open-rotor engine ONERA-SE for research application in aero-acoustic optimization," 2024.

[29] Cantrell, R., and Hart, R., "Interaction between sound and flow in acoustic cavities: Mass, momentum, and energy considerations," *The Journal of the Acoustical Society of America*, Vol. 36, No. 4, 1964, pp. 697-706.

[30] Blandeau, V., "Aerodynamic broadband noise from contra-rotating open rotors," Ph.D. thesis, University of Southampton, 2011.

[31] Cambier, L., Heib, S., and Plot, S., "The Onera elsA CFD software: input from research and feedback from industry," *Mechanics & Industry*, Vol. 14, No. 3, 2013, pp. 159-174.

[32] Spalart, P. R., Belyaev, K. V., Shur, M. L., Kh Strelets, M., and Travin, A. K., "On the differences in noise predictions based on solid and permeable surface Ffowcs Williams-Hawkings integral solutions," *International Journal of Aeroacoustics*, Vol. 18, No. 6-7, 2019, pp. 621-646.
[33] Prieur, J., and Rahier, G., “Aeroacoustic integral methods, formulation and efficient numerical implementation,” *Aerospace Science and Technology*, Vol. 5, No. 7, 2001, pp. 457–468.

[34] Rahier, G., Huet, M., and Prieur, J., “Additional terms for the use of Ffowcs Williams and Hawkings surface integrals in turbulent flows,” *Computers & Fluids*, Vol. 120, 2015, pp. 158–172.