Yohan Jordan

Dirk Von Wissel

Adrian Dolha

Jakob Mauss

Full Virtualization of Renault's Engine Management Software and Application to System Development

Keywords: Critical Transportation Systems Topic: Processes, methods and tools, in particular: virtual engineering and simulation powertrain control, software development, frontloading, validation and test, virtual ECU

Virtualization allows the simulation of automotive ECUs on a Windows PC executing in closed-loop with a vehicle simulation model. This approach enables to move certain development tasks from road or test rigs and HiL (Hardware in the loop) to PCs, where they can often be performed faster and cheaper. Renault has recently established such a virtualization process for powertrain control software based on Simulink models. If the number of runnables exceeds a threshold (about 1500) the execution of the virtual ECU is no longer straight forward and specific techniques are required. This paper describes the motivation behind a Simulink model based process, the virtualization process and applications of the resulting virtual ECUs.

Motivation

Since 2010, Renault has established a framework to develop engine control software for Diesel and Gasoline engines [START_REF] Dressler | A Walk through EMS 2010 Modular Software Development, 4 th European Congress ERTS[END_REF]. The framework is heavily based on MATLAB/Simulink and the idea of modelbased development, which facilitates the carry-over and carry-across of application software between software projects. In the Renault EMS architecture software is composed in to about 20 functions, such as Air System, Combustion etc. A function consists of modules. A module is the smallest testable software unit and contains runnables to be scheduled and executed by the Operating System (Os) of the ECU. The Renault EMS development process includes basically the following steps [START_REF] Dirk Von Wissel | Renault Model-Based Design -Powertrain control development process[END_REF].

1. Specification of about 200 generic configurable modules per ECU using MATLAB/Simulink. 2. Generation of C code (EMS application software) from all module specifications using MATLAB/Simulink Embedded Coder. 3. MiL (Model in the Loop) test and validation of the resulting executable specifications at module level in a simulated system environment, considering only essential interactions with other modules and system environment. This is essentially a back-to-back test to make sure that the Simulink model of a module and the corresponding production C code show equivalent and intended behaviour. To insure software quality, this step is repeatedly performed with steps 1 and 2, based on the simulation capabilities of MATLAB/Simulink. 4. Configuration of modules to fit to the specific needs of a software project, such as absence or presence of certain components. 5. Integration of generated configured C code and hand-coded platform software (basic software) on supplied target hardware, a real ECU that communicates with other controllers via CAN and other busses. 6. Validation and test of all modules on system level using the real ECU. In contrast to step 3, the interactions of all modules and interactions with the system environment are visible then and subject to testing. For example, the Os runs all scheduled runnables then, not just those of the modules considered to be 'essential' for a module under test.

Critical assessment of the above process shows that there is a considerable delay between delivery of a set of specifications to the software project team (at the end of step 3) and system-level tests based on an ECU that runs entire software (step 6). Typical delays are weeks or months.

Consequently, module developers get feedback about the behaviour of their modules late, i.e., weeks or months after releasing a specification (cf. Figure 1).

Fig. 1: Virtual ECU used to frontload system-level validation and test.

The cycle can be considerable shortened by frontloading system-level validation and test. Ideally, system-level test and validation should be performed interleaved with steps 1, 2 and 3 replacing the actual MiL, SiL validation process by a full ECU validation which would provide the software project team with the possibility to validate integrated software within minutes. Technically, such a frontloading of test activities can be achieved using a virtual ECU, shown as "New" in Figure 1: specifications are used to build a virtual ECU for a model in the loop validation, before production code is generated and integrated with real ECU hardware. This way, the software integration test cycle can then be shortened from weeks and months to minutes.

Virtualization of the Renault EMS

Motivated by the above analysis, Renault decided in 2014 to complement the established process for EMS development by the ability to generate virtual ECUs early, as shown in Fig. 2. This new feature enables the validation and test of all embedded functions before actual generation of industrial C code. Code generation takes long because it includes several time consuming validations. Right from the start, the intention was to build a virtual ECU based on executable specifications (Simulink models) of a software project to eliminate the delay introduced by the code generation step. To implement this idea, Renault invited leading providers of tools for the virtualization of ECUs to participate in a benchmark. Participants were asked to use their tool to setup an incremental process for building a virtual ECU from given specifications. Main selection criteria were: time needed to configure an initial setup, time needed for incremental rebuild of a virtual ECU, execution speed of the resulting virtual ECU. In this process, the virtual ECU tool Silver (by QTronic, see [START_REF] Junghanns | Faster Development of AUTOSAR compliant ECUs through simulation[END_REF][START_REF] Mauss | Chip simulation used to run automotive software on PC[END_REF][START_REF] Linssen | Simulation of Networked ECUs for Drivability Calibration[END_REF]) was selected to extend the Renault's existing code developing processes by the ability to work with virtual ECUs."

Incremental build process for a virtual ECU

The Renault EMS is decomposed into about 200 re-usable configurable modules. A module contains so-called runnables. A runnable is specified by a Simulink block that is executed by the operating system (Os) at certain events: periodically, at given crankshaft positions or in response to driver actions such as ignition switched on. The Os specification maps each event to a list of runnables to be executed in response to the event. This coding follows patterns presented in [START_REF] Burke | Templates for Modeling Scheduling Behavior with Simulink and Stateflow[END_REF]. To handle variants, the configuration of a module for a software project also allows to eliminate some runnables in the C code during the linking phase. The virtualization process is structured as shown in Fig. 2. Starting point is a consistent set of modules, i.e. executable specifications given as Simulink models and MATLAB m scripts, a dictionary that defines data types for all variables exchanged between modules, and an Os specification. Based on this, an incremental, fully automated build process for a virtual ECU is implemented as follows: first, time-stamp information available from the file system is used to identify all modules whose specification has been updated since the previous build. For each such module, a wrapper model is generated that implements the interface of the module to other modules and the Os. The wrapper is a Simulink model and contains all runnables of the corresponding module plus some wiring to enable the Os to call the runnables and to establish signal flow with other modules. In a next step the generated wrapper and the wrapped module is turned into C code using Simulink coder. For this step, the General Realtime Target (grt.tlc) shipped with Simulink Coder is used without modification. The C code generated for the module is then post-processed to support on-line calibration. Finally, the C code of all modules of the EMS is compiled together with an Os generated from the given Os specification. The resulting binary runs on Windows PC and is available either as MATLAB/Simulink Sfunction (.mexw32 or .mexw64), as FMU for Co-Simulation [START_REF]Functional Mock-up Interface (FMI)[END_REF] or as a binary (Windows DLL) for execution in Silver.

The above build process is controlled by three configuration files to be provided to the build script  Module interface data types: A dictionary that defines a data type such as float, int16, uint8 etc. for each input and output of the wrapped module. To generate C code, Simulink actually determines data types of signals dynamically via propagation through the block diagram. This can lead to conflicting data types when generating C code for each module separately, as we do here. To prevent this, we need to define signal data types statically. One task of the wrapper is to fix the data types of inputs and outputs to the ones given in the dictionary. A user does not need to list all inputs and outputs of a module. This would be cumbersome, because interfaces are frequently changing. Instead, the wrapper generator derives the list of inputs and output automatically by inspecting the wrapped module.  vECU config options: This selects configuration options to adapt the vECU to specific use cases. A user can specify which additional inputs and outputs the vECU should expose at toplevel. For example, he can mark a module as 'by-passable', to support co-simulation of the vECU with Simulink, with the bypassed module running in Simulink as described in section 4.2. He can also decide to turn all tuneable parameters of a certain module into inputs, to support online-tuning of these parameters, see section 4.3.  Os specification: this is a table that maps each Os event to a list of runnables. The table is used to generate the Os of the virtual ECU. When the named event occurs, the Os runs the corresponding runnables in the given order. We distinguish periodic events such as "10 ms clock tick" from aperiodic events, such as "ignition on". Periodic events are generated by the Os based on simulation time, while aperiodic events are inputs of the Os, to be generated by the engine-or driver model.

Code generation takes about 40 seconds per module, while compilation of the entire C code for a virtual ECU takes about 90 seconds. This means: when only one module is updated, and incremental rebuild of the entire vECU takes less than 3 minutes. This is fast enough to test and validate an updated module in system context immediately after editing the specification of the module. At runtime, a vECU reads a configuration file that contains calibration data. This means, calibration data can be varied and validated without recompiling the virtual ECU. This way, calibration data becomes part of the validation and test loop as well.

The incremental strategy used to build a vECU allows us to generate C code for a single module independent from all other modules. We also tried the opposite strategy: load all modules of the ECU and a generated model of the Os into Simulink and generate C code in one pass. This proved to be infeasible because initialization of the resulting model takes very long (> 10 hours) and the code generator tends to run out of memory during subsequent code generation, even with 64-bit Simulink.

Differences between real and virtual ECU

A virtual ECU is a model of the real ECU. Consequently, not all properties of the real ECU are visible in the virtual ECU. This means that not all tests that are relevant can also be moved to the PC. For some tests, access to real ECU hardware is needed.

The most crucial differences between real and virtual ECU as presented here are  Zero-time execution: the virtual ECU behaves like a device with unlimited computing speed. As a consequence, the Os of a virtual ECU cannot interrupt a running task (no pre-emptive multi-tasking, since a task takes no time to execute). All tasks run exactly as scheduled (for example, cyclically with 1, 2 or 10 ms period), no matter how long its execution may take on the real target hardware.  Missing basic software: Hand-coded basic software of the ECU platform (e.g. drivers for Can, Lin, FlexRay) is not part of the virtual ECU and can therefore not be validated. This is not a fundamental restriction of virtual ECUs, but a consequence of the EMS development workflow shown in Fig. 1.  Different production code: C code representation of calibration data and signals as generated with Simulink Coder is close to but not fully equivalent to the corresponding representations of the production code. For example, on both sides, the same scalar data types are used: float, uint8, int16 etc. However, C code generated by Simulink Coder for grt represents all signals and tuneable parameters as members of a C struct, while in production C code, these are represented as global variables. Consequently, certain defects cannot be found using the virtual ECU. Many of these differences can be avoided by switching to the SiL ([1], [START_REF] Linssen | Simulation of Networked ECUs for Drivability Calibration[END_REF]) or vPiL (virtual processor, [START_REF] Mauss | Chip simulation used to run automotive software on PC[END_REF], [START_REF] Linssen | Simulation of Networked ECUs for Drivability Calibration[END_REF]) type of virtual ECU. However, as explained above, this was not an option right now, due to the time constant for code generation and workflow issues, but might become feasible in the future.

Closed-loop simulation with engine model

In order to study the behaviour of different feedback loops of the EMS, the vECU created by the above process are typically executed in closed-loop with an engine model. Currently, real-time capable engine models developed by Renault using the multi-physics modelling tool LMS Imagine.Lab Amesim (LMS Amesim) for HiL applications are re-used for this purpose. However, an appealing property of virtual ECUs is its complete independence from real time. A virtual ECU can easily be coupled with plant models that run slower or faster than real time. This enables the use of computational expensive plant models, for example in the domain of combustion simulation or exhaust after treatment.

Such models cannot be used in a real-time environment, but without problems in conjunction with virtual ECUs.

Runtime performance of the virtual ECU

A virtual ECU for a typical Renault EMS loads and initializes in less than 5 seconds on a Windows PC. In our experiments, ECU behaviour was sampled with 1 ms steps. Execution speed of the vECU on PC depends on the number of variables (outputs of the vECU) to be recorded during simulation.

 When measuring 170 variables per ms, execution of the vECU is 4 times faster than real time.

 When measuring 20.000 variables per ms, execution speed is 3 times slower than real time. The engine model used has been developed for HiL applications and is therefore very fast. In a closed-loop simulation with a virtual ECU on PC, the engine model consumes less than 10% of the computing time.

Virtual ECUs: Related Work

In this section, we review existing work on the virtualization of ECUs and motivate the approach chosen here. As shown in Fig. 4, we distinguish three options to create a virtual ECU for execution on a PC. which instructs the MCU to take the 32-bit integer currently stored in data register d1, decrement the value by one and store the result in data register d5. A chip simulator reads instructions like this from simulated memory, decodes these and performs the corresponding computation, which may update register and memory content. Chip simulators come in different flavours: A chip simulator is instruction accurate, if it correctly simulates the entire instruction set of the target MCU. A chip simulator is cycle accurate, if the simulator also correctly predicts how many clock ticks it takes the MCU to run an instruction. Instruction accurate simulators at typically an order of magnitude faster than cycle accurate simulators because a cycle accurate simulator requires a far more detailed model of the underlying MCU, to take into account effects of multi-core processing, instruction pipelines, caching strategy, and speed of memory access on processor speed. Besides, chip simulators differ in the coverage of their chip model: An instruction set simulator just simulates the instruction set of the MCU, while a platform simulator covers also on-chip peripherals integrated on the chip, and maybe even ASICs and other MCU external chips found in the ECU. A pure instruction set simulator can be evolved into a platform simulator by integrating additional models, typically using some kind of plugin API provided by the simulator.

Cycle accurate platform simulators can be implemented using hardware description languages such as SystemC or VHDL, often based on design data (IP) owned by the chip maker, while instruction set simulators are typically implemented based on a description of the corresponding instruction set.

In the context of engine control (and beyond), there are two major use cases for chip simulation  Suppliers of ECU hardware and MCUs use chip simulators to validate their ECU and chip designs, for example w.r.t computing performance. Such analysis typically requires a cycle accurate platform simulator as e.g. provided by Synopsys.  OEMs that need to integrate and calibrate supplied EMS software use chip simulators to run the EMS software on PC in cases where they have access to the hex file, but no access to the C code of the EMS. Fast instruction set simulators are often used for this type of application. For example, [START_REF] Mauss | Chip simulation used to run automotive software on PC[END_REF][START_REF] Linssen | Simulation of Networked ECUs for Drivability Calibration[END_REF] reports about simulation of an engine controller using Silver's chip simulator. On a typical PC, a modern ECU runs then approximately in real-time.

SiL: Compiling the C Code for PC

In many application scenarios, the C code of the ECU to virtualize is available. A virtual ECU can then be created by compiling the C code for PC. Compared to chip simulation, the resulting virtual ECU runs typically faster and provides better support for C level debugging, if compiled with debug option. This is important, because testing and debugging of EMS code on PC is a mayor use case for virtual ECUs. Typically, only the modules of the application software and optionally selected modules of the basic software are compiled for PC, but not the operating system (Os) and hardware drivers, such as the AUTOSAR MCAL layer. Instead, the services provided by Os and drivers are re-implemented for the PC environment. In the virtual ECU, these PC-specific implementations replace their corresponding MCU-specific counterparts.

Some tools used to create virtual ECUs based on given C code are EVE (ETAS), Silver (QTronic, [START_REF] Junghanns | Faster Development of AUTOSAR compliant ECUs through simulation[END_REF][START_REF] Mauss | Chip simulation used to run automotive software on PC[END_REF][START_REF] Linssen | Simulation of Networked ECUs for Drivability Calibration[END_REF]) and VEOS (dSPACE). Available tools differ w.r.t. the build process for vECUs, their handling of AUTOSAR and non-AUTOSAR projects, their ability to process standardized ECU and bus descriptions (such as ASAP2 .a2l, AUTOSAR .arxml, CAN .dbc, LIN .ldf, and FlexRay/Fibex .xml), treatment of Os and drivers and runtime support for debugging, scripting, measurement and calibration. Good support for the above descriptions enables to automate the creation of virtual ECUs.

In principle, a ECU can also be virtualized using a hypervisor as discussed in [START_REF] Heiser | The Role of Virtualization in Embedded Systems[END_REF]. In our case, this would require to run the Os of the EMS ('guest Os') on Windows ('host Os').

MiL: Run the models used to generate C Code

Today, model-based development (MBD) on PC with tools like MATLAB/Simulink (MathWorks), ASCET (ETAS), SCADE (Ansys) or Rational Rhapsody (IBM) is the established method of choice for the development of control software. For automotive applications, Simulink seems to dominate the market.

Surprisingly, the dominance of MBD does not mean that developers are typically able to simulate their ECU on PC, even if they have full access to all models of the ECU. For example, loading and initializing all 200 modules of a typical Renault EMS into 64-bit Simulink takes about 10 hours. As far as we can see, the delay is mostly caused by the run-time propagation of data types in the model initialization phase. Interactive simulation of a module in full ECU context is clearly out of reach then. The only way to achieve reasonable execution times is to apply compilation techniques, either within the Simulink environment, or (as we did here), in a Simulink external integration environment, based on exported C code. As the case study presented here shows, such a compiled environment is not straight forward to set up. In our case, it took many months. We could also not find reports about similar attempts in the literature. It seems that developers have widely accepted the fact, that they cannot run the ECU model inside their model-based IDE. Instead, they have to wait until the C code generated by the models has been integrated with the target ECU hardware, and the ECU software can be executed on a HiL system. Of course, this practice fundamentally contradicts the idea of MBD which is about executable models whose behaviour can be explored and assessed during design.

Motivation for choosing the MiL type of virtual ECU

The objective of the project reported here has been to improve the quality of executable module specifications by providing a virtual ECU that runs all ECU modules simultaneously, only minutes after a module has been edited. It was therefore a natural choice to build the vECU from the module source (Simulink .mdl files), which lead us to a MiL type of vECU. The production C code generator is a shared and limited resource here, not easily available in the working environment of module developers. This currently rules out the SiL type of vECU. Hex files for the target MCU become only available weeks after a module edit (see Fig. 1), which rules out the vPiL type of vECU.

Applications of virtual ECUs at Renault engine development

In 2016, Renault created the first fully functional virtual EMS using the process described above. The process has been repeated for about 6 releases (updates and different platforms) of the EMS software since then. Besides, we focus work on setting up the tool chain for the most promising applications of virtual ECUs. These applications are briefly described here, ordered by their location on the time-line of the development process. For each such application, we also sketch the current state of the implementation, if applicable.

Definition of system requirements

In the requirement definition phase of a new function, a virtual ECU can be used to complete an environmental model when interactions between systems through control actions need to be taken into account to reflect the real system behaviour. This way, simulation can be used to derive (even quantify) requirements. We have not yet implemented this idea based on virtual ECUs as described above.

Module development in system context

As described above, a main use case for virtual ECUs is validation and test of module specifications in system context before the specifications are forwarded to the software project team. Such tests are best performed from within the development environment (IDE) for modules, which is Simulink here.

To support this, virtual ECUs can be configured at compile time for co-simulation with Simulink. In the vECU configuration dialog, a module developer can select his or her module for execution in Simulink, bypassing (replacing) the implementation of the module provided by the virtual ECU. The module developer can then edit his module in Simulink, for example modify the Simulink block that implements a certain runnable. After such an edit, the developer can immediately run the virtual ECU, without rebuilding it. This runs all runnables of all modules of the vECU, except those of the bypassed module, which run in Simulink then. This way, the effect of an edit on behaviour can immediately be seen by the developer in the context of all other modules of the ECU. In such a co-simulation scenario, rebuild of the vECU is only required when the Os specification or the signal interface for the edited module changes, or to get more recent implementations of other modules into the validation loop. Thanks to incremental build, such a rebuild takes only minutes.

Pre calibration of the EMS

At least 50% of the development time for engine controllers is spent for engine calibration, i.e. tuning of software parameters of the EMS in order to simultaneously meet customer demands (for example, 'fun to drive'), as well as legal requirements concerning emissions and fuel consumption. With virtual ECUs, Renault started to frontload calibration related activities as well. This is called pre-calibration. The objective is to develop a better starting point for calibration early, in order to gain more time for fine tuning in later phases of development. As pointed out before, a Silver virtual ECU can be configured to load calibration data at runtime from a human readable text file. This means that calibration data can be varied without rebuilding the vECU, which speeds up the calibration process.

To support interactive on-line tuning of parameters, a virtual ECU can also be configured to expose all tuneable parameters of selected modules as permanent (dynamic) inputs of the vECU. This way, a calibration engineer can vary parameters even during simulation, typically using a slider provided by Silver's user interface.

Virtual integration of modules before production C code is generated

Virtual ECUs enable us to virtually integrate all modules of an EMS and to test the resulting virtual ECU in closed loop with an engine simulation before turning the modules into production C code. This helps to shorten development cycles. Integration problems can now be detected more early, before integrating production C code on real ECU hardware. To exploit this idea, we are currently comparing simulation results recorded on a HiL test system with simulation results computed on PC using a virtual ECU for sufficiently similar versions of EMS software and calibration data. Unexpected differences quickly point us to potential problems.

Move selected tests from HiL to MiL

Many of the engine tests that run today on limited resources like HiL test systems, engine test rigs and real vehicles could actually be moved to highly available PCs, because they do not depend on properties of the underlying hardware platform. In such a process, only certain releases of the EMS software (in particular the final release) are tested using real target hardware, while all other intermediate releases are tested on virtual platforms. This can for example be used to increase test coverage. Virtual ECUs are a key tool to implement such a mixed strategy. This idea has not yet been explored using the virtual ECUs described here.

Support joint development within the Renault-Nissan partnership

Renault and Nissan jointly develop certain EMS modules. The development of a shared set of reusable EMS modules requires rework of interfaces and methodologies on both sides. This is a complex process. In this context, Silver virtual ECUs are currently used as a platform for virtual integration of modules. This helps to shorten development cycles and to discover and eliminate problems more early.

Conclusion

Renault started to use virtual ECUs to frontload test and calibration related activities during the development of engine management software. First results of these activities are reported in this paper and have been encouraging so far. In the long run, the existing MiL-based virtual ECUs will be complemented by SiL-based virtual ECUs. The latter are based on production C code and will contain basic software as well, which further minimizes the behavioural gap between real and virtual ECUs. This will enable us to move even more development steps to the virtual platform.

Fig. 2 :

 2 Fig. 2: Process to build a virtual ECU from given specifications

Fig. 3 :

 3 Fig. 3: Virtual ECU of Renault EMS connected to engine model

Fig. 4 :

 4 Fig. 4: Virtual ECU for MiL, SiL and virtual PiLDepending on the use case to be supported, a virtual ECU is created from either  virtual PiL: the hex file resulting from compiling C code for the target processor  SiL: C code compiled for PC  MiL: models that generate the C code 3.1 Virtual PiL: Chip simulation to run a hex file on PC At the heart of an engine control unit (ECU) is a microcontroller unit (MCU), i.e. a digital processor with additional peripherals such as CAN controllers, analog digital converters (ADCs) etc. integrated on a single chip. Wide spread MCU families for engine control are Tricore (Infineon), PowerPC (Freescale/NXP, STM), and v850 (Renesas). MCUs for engine control often have two or more cores for fast concurrent computing. A virtual ECU can be created by simulating the corresponding MCU. The hex file resulting from compiling the EMS software for the target MCU contains binary program code and calibration data. For example, for a Tricore MCU, the four bytes 0x50, 0x1F, 0xF1, 0x8B (given in hexadecimal notation) encode the assembler instruction add %d5, %d1 -1

Fig. 5 :

 5 Fig. 5: D-EIPF running a full vehicle simulation 4.7 Improve existing vehicle-level simulation Renault and Nissan use the Digital Electronic Integration PlatForm (D-EIPF, see [8]) to validate networked ECUs in vehicle context. With D-EIPF, the communication behaviour of networked ECUs of an entire vehicle can be simulated without the need to access real ECU or other vehicle hardware. Since 2010, D-EIPF has been developed in-house on top of the RT-LAB Orchestra API (OPAL RT). Network simulation is based on CANoe (Vector). Silver has been integrated into the D-EIPF environment. This way much more realistic engine models (i.e. a virtual ECU running in closed loop with an engine plant model) can now be used to validate overall vehicle behaviour, which increases the scope of tests that can be executed on D-EIPF.


Christophe Marchetto
email: christophe.marchetto@usherbrooke.ca

Laurent Maxit

Olivier Robin

Alain Berry

Vibroacoustic response of panels under diffuse acoustic

Keywords: numbers: : 43, 40, At, 43, 40, Dx

INTRODUCTION

The experimental vibroacoustic characterization of panels under a diffuse acoustic field (DAF) excitation is of great interest for the industry. This excitation is commonly used to determine the sound reduction index of panels as described in several standards using coupled reverberant-reverberant room [START_REF]Acoustics -Laboratory measurement of sound insulation of building elements -Part 2: Measurement of airborne sound insulation[END_REF] [START_REF]ASTM E90-09 Standard Test Method for Laboratory Measurement of Airborne Sound Transmission Loss of Building Partitions and Elements[END_REF] or reverberant-anechoic room [START_REF]Acoustics -Measurement of sound insulation in buildings and of building elements using sound intensity -Part 1: Laboratory measurements[END_REF] [4] laboratory facilities. Theoretically, a DAF is defined as an infinite set of uncorrelated plane waves with uniformly distributed incidence angles. In standard laboratory measurements, this excitation is reproduced using a reverberant room and only partially corresponds to its theoretical definition, especially below the Schroeder frequency of the room where the sound field is dominated by well-defined acoustic cavity modes. Even above the Schroeder frequency, the pressure field is not perfectly homogeneous and the lack of grazing incidence plane waves has been pointed out in the literature. Inter-laboratory variations of vibroacoustic measurements in reverberant rooms can be attributed to these phenomena, but other parameters are involved such as room dimensions, niche effects, panel mounting conditions, aperture size and measuring protocols [START_REF] Bravo | Variability of low frequency sound transmission measurements[END_REF] [6] [START_REF] Garg | Measurement uncertainty in airborne sound insulation and single-number quantities using sound pressure and sound intensity approaches[END_REF].

In this context, the aim of this study is to investigate an alternative and robust approach to experimentally characterize a panel's response to a DAF excitation by using only the theoretical model of this excitation to overcome the limitations of a reverberant room measurement discussed above. Indeed, the mathematical formulation of a panel's vibro-acoustic response when submitted to random excitations in the wavenumber domain allows estimating the system's response, at any point on the structure or in the acoustic medium, from wall-pressure cross spectral density (CSD) functions (characterizing the excitation) and from so-called 'sensitivity functions', which were introduced in [START_REF] Bonness | Low-wavenumber turbulent boundary layer wall-pressure measurements from vibration data on a cylinder in pipe flow[END_REF] [START_REF] Hambric | Vibrations of plates with clamped and free edges excited by low-speed turbulent boundary layer flow[END_REF] for the analogous problem of panels excited by a turbulent flow. The latter are defined as the panel's response to wallpressure plane waves and characterize the panel's vibroacoustic behavior. The estimation of the panel's response submitted to a DAF excitation therefore only requires the experimental measurement of sensitivity functions in the acoustic wavenumber domain.

A method is proposed for estimating the sensitivity functions experimentally. Whereas the direct interpretation of the sensitivity functions would require exciting the panel by sets of wall plane waves, which is not easy from an experimental point of view, an alternative method based on a reciprocity principle is proposed. The reciprocity principle states that the sensitivity functions at any point on the structure or in the acoustic medium are equivalent to the panel's velocity response expressed in the wavenumber domain when the system is excited by a vibration or acoustic source at the same point. From an experimental point of view, it is then only necessary to excite the system with a vibration or acoustic source and to apply a wavenumber transform to the measured transfer function between the panel velocity and the source magnitude to obtain the sensitivity functions for a wide range of wavenumbers. Globally, the proposed experimental process consists in exciting the panel with the source of given magnitude at the point of interest. The spatial vibratory response of the panel is then measured with a scanning laser vibrometer. In a subsequent post-processing phase, a discrete 2-D wavenumber transform of the measured vibratory field is performed to deduce the sensitivity functions. Finally, using the wall-pressure model of a DAF and the previously estimated sensitivity functions, the response when the panel is excited by a DAF can be deduced at any point on the structure or in the acoustic medium. To evaluate the sound transmission loss, this process is repeated for a series of points belonging to a surface surrounding the panel to estimate the acoustic intensity at these points when the panel is excited by a DAF.

The remainder of the paper is organized as follows: the considered vibro-acoustic problem and the quantities characterizing a panel under a DAF are described in Sec. II. Then, the mathematical formulation of the vibro-acoustic problem is presented in Sec. III where the sensitivity functions appearing in the formulation are defined. An alternative interpretation of these functions based on the reciprocity principle is proposed in Sec. IV. This interpretation suggests a simple implementation for measuring the sensitivity functions. The proposed methodology for characterizing the panel response under a diffuse acoustic field is summarized in Sec. V. Numerical and experimental validations are provided in Sec. VI.

Finally, a comparison with measurements performed following standards [START_REF]Standard Test Method for Laboratory Measurement of Airborne Sound Transmission Loss of Building Partitions and Elements Using Sound Intensity[END_REF] is proposed in Sec. VII.

II. VIBROACOUSTIC CHARACTERIZATION OF PANELS UNDER DIFFUSE ACOUSTIC FIELD

Let us consider a baffled panel of surface Σ p with arbitrary boundary conditions separating two semi-infinite acoustic domains. Each of these acoustic domains is characterized by a To experimentally characterize the vibroacoustic behavior of this panel, two quantities are considered:

-the one-sided normal velocity frequency response v at a given point on the panel, which for random excitations is given by the auto spectral density (ASD) function of the normal velocity G vv (x, f ), -the transmission loss (TL) defined by

T L (f ) = 10 log 10 Π inc (f ) Π rad (f ) , ( 1
)
where f is the frequency and is considered positive.

The incident acoustic power and the radiated acoustic power are denoted Π inc and Π rad , respectively. For a DAF exciting a panel of area Σ p , the incident acoustic power is theoretically

given by [START_REF] Robin | Experimental vibroacoustic testing of plane panels using synthesized random pressure fields[END_REF] Π

inc (f ) = G p b p b (f ) Σ p 8 ρ 0 c 0 , (2)
where G p b p b (f ) is the one-sided frequency wall-pressure ASD function and the factor 8 accounts for pressure doubling at the panel surface [START_REF] Chazot | Diffuse Acoustic Field Produced in Reverberant Rooms: A Boundary Diffuse Field Index[END_REF]. The radiated power of the panel into the receiving medium is obtained by integrating the normal active sound intensity flow passing through a virtual surface Σ v surrounding the panel

Π rad (f ) = Σv I act (x, f ) n dx, (3)
n being the unit exterior vector normal to Σ v , dx the surface element and I act (x, f ) the active sound intensity at point x. The latter is directly related to the one-sided frequency CSD function G pv 0 (x, f ) between the sound pressure p and the particle velocity v 0 at point

x with [START_REF] Fahy | Sound Intensity[END_REF]

I act (x, f ) = Re{G pv 0 (x, f )}. ( 4
)
Theoretically, the radiated power is obtained by solving the formal integral in Eq. (3).

For numerical and experimental applications, the integral of Eq. ( 3) may be approximated by discretizing the surface Σ v and using the rectangular integration rule:

Π rad (f ) ≈ x∈σv I act (x, f ) n δx, (5)
where σ v represents the set of points defined on Σ v and δx is the elementary point area.

To characterize the vibration response of a panel to a DAF, it is then necessary to evaluate G vv at the considered point on the panel while the evaluation of G pv 0 for the set of points σ v is required to estimate the TL. An approach for evaluating these quantities based on deterministic transfer functions and using a reciprocity principle is presented in the following sections.

III. MATHEMATICAL FORMULATION OF THE VIBROACOUSTIC RESPONSE OF PANELS UNDER DIFFUSE ACOUSTIC FIELD

Let us consider the blocked wall-pressure field p b (x, t) exerted on the panel by a DAF excitation at point x as a function of time. The response of the panel at point x when the panel is excited by p b (x, t) is denoted α(x, t). If x is on the panel, α stands for v whereas it stands for p or v 0 if x is in the acoustic domain. This response can be expressed by the convolution product [START_REF] Maury | A wavenumber approach to modelling the response of a randomly excited panel, part 1: general theory[END_REF] α

(x, t) = Σp ∞ -∞ h α/Fn (x, x, t -τ ) p b (x, τ ) dτ dx, (6)
where h α/Fn (x, x, t) is the impulse response (structural velocity, acoustic pressure or particle velocity-wise) at point x for a normal unit force applied at point x. Assuming that the random process is ergodic, the cross-correlation function R αα (x, t) is defined by

R αα (x, t) = ∞ -∞ α (x, τ ) α (x, t + τ ) dτ. ( 7
)
where α also designates v, p or v 0 . Introducing Eq. ( 6) in Eq. ( 7) and performing a time Fourier transform of the resulting expression of the cross-correlation function gives the spacefrequency spectrum S αα (x, ω), which after some manipulations (see [START_REF] Maury | A wavenumber approach to modelling the response of a randomly excited panel, part 1: general theory[END_REF] for details) can be written as

S αα (x, ω) = Σp Σp H α/Fn (x, x, ω) H * α /Fn x, x, ω S p b p b x, x, ω dxd x, (8)
where H α/Fn (x, x, ω) is the time Fourier transform of h α/Fn (x, x, t) and corresponds to the panel frequency response function (velocity, pressure or particle velocity-wise) at point x when it is excited by a normal force F n applied at point x; S p b p b x, x, ω is the time Fourier transform of the cross-correlation function of the blocked wall-pressure; finally * is the complex conjugate. Defining the wavenumber-frequency spectrum of the wall-pressure S p b p b (k, ω) as the wavenumber transform of the space-frequency spectrum S p b p b x, x, ω , one has

S p b p b x, x, ω = 1 4π 2 ∞ -∞ S p b p b (k, ω) e jk( x-x) dk. ( 9
)
where k = (k x , k y ) is the wavevector defined in the plane (x, y) and dk is the two-dimensional wavenumber element. By introducing Eq. ( 9) in Eq. ( 8) and rearranging the terms, one obtains

S αα (x, ω) = 1 4π 2 ∞ -∞ H α (x, k, ω) H * α (x, k, ω) S p b p b (k, ω) dk, (10)
where

H α (x, k, ω) = Σp H α/Fn (x, x, ω) e -jkx dx. (11)
The H α (x, k, ω) functions are called the sensitivity functions [START_REF] Lin | Probabilistic theory of structural dynamics[END_REF] and characterize the vibroacoustic behavior of the panel. The function H α/Fn (x, x, ω) is the time Fourier transform h α/Fn (x, x, t) and, therefore, corresponds to the transfer function between the panel velocity frequency response and the frequency spectrum of the applied effort.

The wall-pressure CSD function in the space-frequency domain of a DAF can be expressed by [START_REF] Cook | Measurement of correlation coefficients in reverberant sound fields[END_REF]

S p b p b (r, ω) = S p b p b (ω) sin (k 0 r) k 0 r , (12)
where r = |x -x|, k 0 = ω/c 0 is the acoustic wavenumber and S p b p b (ω) is the wall-pressure ASD function. The space-wavenumber transform of Eq. ( 12) gives the wall-pressure CSD function in the wavenumber-frequency space

S p b p b (k, ω) = S p b p b (ω) Φ p b p b (k, ω) , (13)
where

Φ p b p b (k, ω) =        2π k 0 1 √ k 2 0 -|k| 2 if |k| < k 0 0 if |k| ≥ k 0 . ( 14
)
As the wall-pressure CSD function of a DAF is null for wavenumbers larger than the acoustic wavenumber, the integration domain involved in Eq. ( 10) can be restricted to the wavenumbers contained in the acoustic domain (i.e., |k| < k 0 ). Moreover, in practice, this integral is approximated considering a set of wavevectors in the acoustic domain Ω k and using the rectangular integration rule. It should be stressed here that S αα (x, ω) is a two-sided spectrum as a function of the angular frequency. It can be related to the one-sided spectrum as a function of the frequency G αα (x, f ) by

G αα (x, f ) = 4π S αα (x, ω) . (15)
For the sake of coherence with experiments, one-sided frequency spectra will be considered in the remainder of the article.

According to Eq. ( 10) and ( 15), the one-sided frequency ASD function of the velocity v at point x (z = 0) of a panel excited by a DAF can be estimated with

G vv (x, f ) ≈ 1 4π 2 k∈Ω k |H v (x, k, ω)| 2 G p b p b (f ) Φ p b p b (k, ω) δ 2 k, (16)
whereas the one-sided frequency CSD function between the pressure p and the particle velocity v 0 at a given point x into the acoustic domain can be estimated with

G pv 0 (x, f ) ≈ 1 4π 2 k∈Ω k H p (x, k, ω) H * v 0 (x, k, ω) G p b p b (f ) Φ p b p b (k, ω) δ 2 k, (17)
where δk represents the wavenumber resolution and

G p b p b (f ) = 4π S p b p b (ω) is the one-sided
frequency ASD function of the blocked wall-pressure.

To evaluate these two quantities, the sensitivity functions H v , H p and H v 0 for wavenumbers belonging to Ω k are thus to be determined. A direct interpretation of these sensitivity functions can be deduced from Eq. [START_REF] Chazot | Diffuse Acoustic Field Produced in Reverberant Rooms: A Boundary Diffuse Field Index[END_REF]. Since H α/Fn (x, x, ω) is the response α at point x for a unit normal force at point x, H α represents the frequency response α at point x due to a wall-pressure plane wave of wavevector -k (i.e., due to the pressure field e -jkx ). This direct interpretation is depicted in Figs. 2(a), 2(c) and 2(e) for H v , H p and H v 0 , respectively.

The sensitivity functions must, therefore, be estimated only at the point of interest x and for the set of wavevectors Ω k . A large number of waves should be considered to entirely cover the acoustic wavenumber domain. Moreover, from an experimental point of view, wall-pressure plane waves cannot be easily reproduced. To circumvent these issues, another interpretation of these sensitivity functions based on the reciprocity principle is given in the next section.

IV. ALTERNATIVE INTERPRETATION OF THE SENSITIVITY FUNCTIONS

In order to propose another interpretation of the sensitivity functions, let us consider the standard reciprocity principle which states that the response of a system is invariant with respect to the exchange of points of excitation and observed response [START_REF] Fahy | Some Applications of the Reciprocity Principle in Experimental Vibroacoustics[END_REF]. Following the previous notation, it can be translated into

H α/Fn (x, x, ω) = H v/ ᾱ(x, x, ω), (18)
where H v/ ᾱ(x, x, ω) is the frequency response function between the panel velocity and a source ᾱ, dual of α. As the normal force F n is applied on the panel, the exchanged observation point is also on the panel which explains why the right-hand side of Eq. ( 18) remains the velocity response of the panel, regardless of α. However the type of excitation source ᾱ depends on α and it will be detailed below for each quantity considered for α.

Sticking to the general case, by introducing Eq. ( 18) in Eq. ( 11) one obtains

H α (x, k, ω) = Σp H v/ ᾱ(x, x, ω)e -jkx dx. (19)
The right hand side of Eq. ( 19) can be interpreted as the space-wavenumber transform of

H v/ ᾱ(x,
x, ω) with respect to the spatial variable x. The points x become observation points on the panel surface Σ p , which means that the space-wavenumber transform is performed over the vibration velocity field of the panel. To sum up, the sensitivity function

H α (x, k, ω)
may be obtained by exciting the panel with a source ᾱ at point x and by calculating the space-wavenumber transform of the panel velocity frequency response normalized by the source frequency spectrum. This second interpretation of the sensitivity functions is now detailed for the three cases involved in the evaluation of the panel response excited by a DAF.

Case of plate velocity (α = v):

The reciprocity principle states [START_REF] Maxit | Prediction of flow induced sound and vibration of periodically stiffened plates[END_REF] that the ratio of the normal velocity of the panel at point x over the applied normal force at point x is equal to the ratio of the normal velocity of the panel at point x over the normal force applied at point x. Eq. ( 18) becomes

H v/Fn (x, x, ω) = H v/Fn (x, x, ω). (20)
In this case, ᾱ is a normal force and thus, the sensitivity function

H v (x, k, ω) is obtained
by exciting the panel with a normal force at point x and by performing a space-wavenumber transform of the transfer function between the panel vibration velocity response and the force frequency spectrum, as illustrated in Fig. 2(b).

Case of radiated pressure (α = p): Lyamshev reciprocity relations for elastic structures excited by point forces [START_REF] Fahy | Some Applications of the Reciprocity Principle in Experimental Vibroacoustics[END_REF] indicate that the ratio of the pressure at point x over the applied normal force at point x is equal to the ratio of the normal velocity of the panel at point x over the volume velocity Q v of a monopole source placed at point x. The demonstration of this particular reciprocity relation according to Lyamshev [START_REF] Lyamshev | A method for solving the proble; of sound radiation by thin elastic shells and plates[END_REF] is provided in Appendix A 1 (this classical demonstration is useful to introduce reciprocity in terms of particle velocity and dipole source strength, which is demonstrated in appendix A 2). In this case ᾱ is a monopole source of volume velocity Q v and Eq. ( 18) becomes

H p/Fn (x, x, ω) = H v/Qv (x, x, ω). (21)
The sensitivity function H p (x, k, ω) is therefore obtained by exciting the panel with a monopole source at point x and by performing a space-wavenumber transform of the panel vibration velocity response normalized by the frequency spectrum of the monopole's volume velocity, as illustrated in Fig. 2(d).

Case of particle velocity (α = v 0 ): Fahy stated that the reciprocity relationship could be extended to acoustic dipoles and particle velocities [START_REF] Fahy | Some Applications of the Reciprocity Principle in Experimental Vibroacoustics[END_REF], however no demonstration could be found in the literature. As a point dipole can be represented by a point force injected in the fluid, the ratio of the particle velocity at point x over the applied normal force at point

x is equal to the ratio of the normal velocity of the panel at point x over the force injected in the fluid F 0 at point x. This latter reciprocity relation is demonstrated in Appendix A 2.

The dual source ᾱ being a dipole source of force F 0 injected in the fluid, Eq. ( 18) now becomes

H v 0 /Fn (x, x, ω) = H v/F 0 (x, x, ω) (22)
As demonstrated in Appendix A 2, the force F 0 is injected in the same direction n as the desired direction of the particle velocity v 0 . One can obtain the sensitivity function To sum up, the sensitivity functions can be obtained by exciting the system at the point of interest x and by performing a space-wavenumber transform of the panel velocity field.

H v 0 (x,
In practice, the vibratory field has to be measured on a regular grid of points denoted Γ x, using a scanning laser vibrometer, for example. The space-wavenumber transform is therefore approximated by a discrete Fourier transform. In order to avoid aliasing effects, the spatial resolution δx over Γ x should be determined so that the spatial variations of the vibratory field can be correctly represented by the grid of points. For a homogeneous isotropic thin panel, δx should be less than or equal to a quarter of the natural flexural wavelength of the panel λ f at the highest frequency of interest. For a more complex panel, a preliminary study should be carried out to define this parameter (for instance, by using a numerical model of the panel or by using a trial and error procedure).

V. DESCRIPTION OF THE PROPOSED METHODOLOGY

A methodology for experimentally estimating the vibroacoutic response of a panel excited by a DAF is now presented. This methodology is based on Eqs. ( 16) and [START_REF] Maxit | Prediction of flow induced sound and vibration of periodically stiffened plates[END_REF], and the second interpretation of the sensitivity functions, as described in the previous section.

Vibration response of the panel: The methodology for evaluating the velocity ASD function G vv at a given point x of the panel (z = 0) can be summarized as follows:

-Excite the panel with a normal mechanical force at point x (for instance by using a shaker) and measure the normal velocity response of the panel at points x ∈ Γ x to determine H v/Fn (x, x, ω),

-Perform a discrete Fourier transform of the panel velocity response H v/Fn (x, x, ω)

(with respect to x) to obtain the sensitivity functions

H v (x, k, ω) at point x for k ∈ Ω k ,
-Use Eqs. ( 16) and ( 14) to estimate the velocity ASD function G vv at point x under an ideal DAF excitation.

Acoustic response of the panel: The acoustic response of the panel is characterized by the TL as described in Sec. II. It can be obtained by following the next five steps:

-Excite the panel with a monopole source at a given point of interest x and measure the normal velocity response of the panel at points x ∈ Γ x to determine H v/Qv (x, x, ω),

-Excite the panel with a dipole source at a given point of interest x and measure the normal velocity response of the panel at points x ∈ Γ x to determine H v/F 0 (x, x, ω),

-Perform a discrete Fourier transform of the panel velocity responses obtained for both monopole and dipole cases to estimate the sensitivity functions H p (x, k, ω) and

H v 0 (x, k, ω) at point x,
-Calculate the pressure -particle velocity CSD function at point x using Eqs. [START_REF] Maxit | Prediction of flow induced sound and vibration of periodically stiffened plates[END_REF] and ( 14), -Calculate the active sound intensity at point x using Eq. ( 4).

The five previous steps are repeated for points x ∈ σ v (discretizing the whole virtual surface Σ v surrounding the panel) to calculate the radiated power using Eq. ( 5). The TL is finally deduced using Eq. ( 1) while the incident acoustic power is evaluated with Eq. ( 2).

VI. NUMERICAL AND EXPERIMENTAL VALIDATION OF THE PROPOSED APPROACH A. Test case description

For numerical and experimental validation purposes a test case is considered, which consists in a rectangular thin aluminum plate, simply supported on its four edges, baffled, and submitted to a DAF excitation on one side. This baffled plate separates two semi-infinite domains filled with air (ρ 0 = 1.3 kg.m -3 and c 0 = 343 m.s -1 ). The plate's geometrical and mechanical properties are detailed in Table I. The structural loss factor η mn of the (m, n) mode has been experimentally estimated using the -3 dB bandwidth method on the first few resonances of the plate and is taken into account in the numerical simulations. A mean value of η = 0.005 has been measured. Simply-supported boundary conditions have been chosen because they lead to a simple analytical solution of the plate equation of motion.

In addition, the experimental setup proposed by Robin et al. [START_REF] Robin | A plane and thin panel with representative simply supported boundary conditions for laboratory vibroacoustic test[END_REF] for reproducing these boundary conditions has already been validated.

The frequency range of interest is [170, 2000 Hz] with a frequency resolution of 0.625 Hz.

The low frequency limit is set according to the frequency response of the monopole source and the high frequency limit has been chosen arbitrarily. This frequency range is below the critical frequency of the panel, f c , given by

f c = c 2 0 2π ρh D , (23)
where

D = Eh 3 12(1-ν 2 )
is the flexural stiffness. For the considered case, f c = 3867 Hz.

In this section, we will focus on:

-the velocity sensitivity functions H v at point x M of coordinates (x = 0.06 m, y = 0.3 m, z = 0 m) on the plate, -the pressure and particle velocity sensitivity functions H p and H v 0 at point x N of coordinates (x = 0.06 m, y = 0.3 m, z = 0.1 m) into the acoustic medium. The particle velocity v 0 will be determined in direction z,

-the frequency response of the velocity response at point x M and of the active intensity at point x N in the direction z.

To apply the methodology described in Sec. V, the panel velocity field has to be measured or calculated on a grid of points Γ x. In the following, a uniform mesh of 15 × 13 points is considered in directions x and y respectively and a gap of 30 mm along the edges is left for practical reasons. This leads to a spatial resolution of δ x = δ y = 30 mm and ensures at least 4 points per flexural wavelength for all frequencies of interest. The highest wavenumbers k max

x and k max y that can be resolved in directions x and y, respectively, are given by

k max x = k max y = π δ x = π δ y 105 m -1 . ( 24
)
These wavenumbers are well above twice the acoustic wavenumber (related to the Shannon criterion) at the highest frequency of interest (i.e., k 0 = 37 m -1 at 2000 Hz). As a consequence, the considered grid of points provides correct estimation of the sensitivity functions in the acoustic wavenumber domain Ω k . The wavenumber resolutions δk x and δk y in directions x and y respectively, are given by

δk x = 2π L x 13 m -1 ; δk y = 2π L y 15 m -1 . ( 25
)
These wavenumber resolutions are relatively large because of the small dimensions of the panel. In order to improve the wavenumber resolution, zero-padding is used to obtain a wavenumber resolution of 0.5 m -1 along k x and k y .

In order to assess the accuracy of the reciprocity approach for evaluating the panel sensitivity functions, the results obtained with this approach are compared with those obtained by considering the direct interpretation of these functions (as described in Sec. III). This comparison will allow validating the uniform mesh of 15 × 13 points used for the discrete spatial Fourier transform of the panel velocity field. The numerical model used for this study is described in Appendix B.

B. Experimental set-up

Fig. 3 shows the experimental setup used to measure the sensitivity functions. The plate was glued on thin blades and fixed on a frame as described in [START_REF] Robin | A plane and thin panel with representative simply supported boundary conditions for laboratory vibroacoustic test[END_REF] to reproduce simply supported boundary conditions. To determine the velocity sensitivity functions H v , the plate was excited by a normal force at point x M of coordinates (x = 0.06 m, y = 0.3 m, z = 0 m).

This force was applied using a TMS SmartShaker K2007E01 with integrated amplifier, which was fed with a swept sine over the considered frequency range and the force was measured using an impedance head PCB288D01 (as shown in Fig. 3(a)). An adapter was used between the impedance head and the plate reducing the area of mechanical coupling to approximately a 5 mm diameter circle.

For acoustic applications, the plate was baffled in a 1 × 1 m 2 plywood panel of 2 cm thickness. The experiment was performed in a hemi-anechoic room and 10 cm thick sound absorbing foam (Decibel France Polyphone 63 T) was placed on the ground and around the plate (see Figs. given distance r in anechoic conditions for a given input voltage U and using the theoretical model of a monopole in free field for the relation between p and Q v . The effect of the tube on the frequency response has thereby been accounted for.

To estimate the particle velocity sensitivity functions at point x N in the direction z, the response of the plate to a dipole source has been reconstructed by exciting the plate with the same monopole source, still fed with a white signal and moved from the previous position by a distance d = 3 cm in the direction of the injected force (in this case z).

Two monopoles close to each other and out of phase have been thereby reconstructed by subtracting the measured transfer functions H v/Qv . A preliminary experimental study in free field conditions and using the plate was performed to determine an appropriate separation d.

It showed that below a separation of 0.5 cm, the vibration fields induced for both positions of the monopole source were not sufficiently different to be noticeably measured. On the other hand, above a separation of 5 cm, the directivity of the reconstructed dipole did not match that of a theoretical dipole. A value of 3 cm for d appeared to be an optimal value for the present case. It should be noted that the condition k 0 d << 1 does not hold at the highest frequencies. However, the induced vibrations were in accordance with the response of a plate to a theoretical dipole. Since the numerical and experimental sensitivity functions and pressure -particle velocity CSD functions are in good agreement (see next section), these results validate the experimental method for reconstructing a dipole source.

According to the methodology described in Sec. V, the vibratory response of the panel has been measured on the grid of 15 × 13 points with a single point Polytec laser vibrometer for each case of the excitation (i.e., force, monopole, dipole). Also, in each case, the time Fourier transform was directly performed in the post-processing software with ten linear averages.

C. Comparison between numerical and experimental results

Sensitivity functions

Fig. 4 shows the velocity sensitivity functions H v obtained with the direct calculation and the reciprocal approach using numerical and experimental data, respectively. They are provided for three different frequencies, the lowest corresponding to the (2,1) vibration mode frequency (Fig. 4(a)) and the two others being off-resonance cases (Figs. 4(b) and 4(c)).

The product of sensitivity functions H p and H * v 0 , which is involved in the expression of G pv 0 (Eq. ( 17)), is shown in Fig. 5 at the same frequencies. The circles of radius k 0 and k f = 2πfc c 0 k 0 , corresponding to the acoustic and flexural natural wavenumbers respectively, are also indicated in Figs. 4 and5.

In Fig. 4 and5, results obtained by simulating numerically the direct and the reciprocity methods match perfectly. This validates the grid of points considered on the panel and the use of zero-padding to improve the wavenumber resolution without affecting the results.

It also validates the method described above to reconstruct a dipole from two monopole sources.

In Fig. 4 and5, the numerical and experimental results are generally in good agreement. One can observe that the sensitivity functions are slightly overestimated experimentally at the (2,1) vibration mode frequency compared to the numerical results. This can be explained by the fact that the modal damping loss factor has been estimated from the response of the plate to a shaker excitation. The added mass from the shaker possibly had an influence on the evaluation of the damping of the (2,1) mode.

A good agreement is particularly noticed within the acoustic wavenumber circle (delineated by a full line). Again, only values in the acoustic wavenumber domain contribute to the plate's vibroacoutic response to a DAF. However, the sensitivity functions are also correctly estimated experimentally for wavenumbers higher than the acoustic wavenumber.

Plate velocity ASD function

The velocity ASD function G vv (x, f ) at point x M of the panel excited by a DAF with a unit wall-pressure ASD function (G p b p b (f ) = 1 Pa 2 .Hz -1 ) has been estimated using Eq. ( 16) and the three previously described sensitivity functions. Fig. 6(a) compares the results obtained with the direct calculation and the numerical reciprocity approach. The two curves are perfectly superimposed, showing that the sensitivity functions obtained with the two approaches are essentially identical in the acoustic domain and on the whole considered frequency range. The small noticeable peaks between resonance frequencies of the plate on both curves are attributable to the wavenumber sampling in Eq. ( 16) which induces errors especially as |k| approaches k 0 , in which case the DAF wall-pressure CSD goes to infinity (see Eq. ( 14)). reciprocity principle.

The experimental and theoretical pressure -particle velocity CSD functions are compared in Fig. 7(b). Their good agreement shows that the active sound intensity is correctly estimated in this situation using the reciprocity approach. Furthermore, the experimental monopole and dipole excitations correctly reproduce the theoretical conditions. The resonance peaks are better estimated compared to those obtained to estimate the plate's velocity ASD function in Fig. 6(b), particularly at high frequencies. This can be possibly explained by the accuracy in positioning the source which, for the experiment with the shaker was done manually (subject to more errors) whereas for the acoustic applications, it was controlled with a robot allowing a higher accuracy. A more likely explanation is the dynamic influence of the mass added with the shaker, which explains the slight shifts of the resonance peaks in Fig. 6(b).

VII. COMPARISON WITH REVERBERANT ROOM MEASUREMENTS

The proposed approach is finally compared with measurements performed at the University of Sherbrooke transmission loss facility (coupled reverberant-anechoic rooms) using a plate similar to the one used in the previous section (similar dimensions, material and boundary conditions) and following test standard ASTM E2249-02 (2016) [START_REF]Standard Test Method for Laboratory Measurement of Airborne Sound Transmission Loss of Building Partitions and Elements Using Sound Intensity[END_REF]. The reverberant room has a volume of approximately 140 m 3 (7.5 × 6.2 × 3 m 3 ), and the Schroeder frequency of the room, above which the sound field can be considered diffuse, is approximately 410 Hz.

The plate was mounted in an existing niche between the coupled anechoic-reverberant rooms (the panel being flush mounted on the reverberant room side). A double-wall structure with mechanical decoupling was then built around the plate to prevent acoustic leaks and flanking paths, as described in [START_REF] Robin | A plane and thin panel with representative simply supported boundary conditions for laboratory vibroacoustic test[END_REF]. A loudspeaker fed with a white noise signal excited the reverberant chamber. The sensitivity functions of this second plate were estimated experimentally using the reciprocity method. Some differences with the sensitivity functions of the plate considered in Sec. VI (not shown here) indicate that the positioning of the force applied with the shaker is not exactly at the considered point x M . Indeed, the position has a significant influence on the vibratory response, particularly at high frequencies where the mode shapes get more complex. The mounting base of the force sensors also has a finite diameter of approximately 5 mm, which makes the applied force not perfectly punctual. However, the modal frequencies correspond to the theory for both plates. The velocity ASD function was calculated using Eq. ( 16) and ( 14) whilst including the measured wall-pressure ASD function Ḡp b p b (f ).

The plate velocity ASD function measured in the reverberant room at point x M is compared to the result obtained with the proposed method in Fig. 8. The two obtained responses are in good agreement up to 800 Hz. Above this frequency, differences are noticeable and can be explained with:

-the inaccuracy in the positioning of point x M where the velocity ASD function is measured in the reverberant chamber (corresponding to the point force position in the reciprocity method) whose influence increases with the frequency, -the deviation of the pressure field in reverberant room to an ideal DAF. An analysis of the pressure field measured with the microphone array shows spatial variations in contradiction with the assumptions of a perfect DAF. Moreover, it is well known that a reverberant chamber has difficulty creating grazing incidence waves. The absence of grazing incidence waves can hardly be quantified as it varies from one reverberant room to another. Some authors suggest corrections on the theoretical model of the excitation to better represent the actual excitation in a reverberant room [START_REF] Kang | An experimental investigation on the directional distribution of incident energy for the prediction of sound transmission loss[END_REF] [21].

B. Sound transmission loss

A second experiment in the transmission loss facility was conducted for estimating the plate TL. A 1/2 inch Bruel & Kjaer rotating microphone was used to measure the spatiallyaveraged sound pressure level L p in far radiation field. The average sound intensity level L I was measured on the anechoic side by using a Bruel & Kjaer sound intensity probe composed of two 1/2 inch microphones with a 12 mm spacing. The sound intensity probe was manually moved 5 cm away from the plate to scan over a parallel surface identical to the plate area, as described in [START_REF]Standard Test Method for Laboratory Measurement of Airborne Sound Transmission Loss of Building Partitions and Elements Using Sound Intensity[END_REF] for the case of a plate flush mounted on the source side. The transmission loss of the structure is given by T L = L p -L I -6 [START_REF] Guy | Measurement of sound transmission loss by sound intensity[END_REF]. An illustration of this experiment and the considered virtual surface is given in Fig. 9(a).

In addition, the proposed methodology described in Sec. V for estimating the TL was applied. A numerical study on the definition of the virtual surface Σ v (which is not detailed in this paper) showed that considering a virtual surface equal to the plate area and positioned 5 cm away from it would result in a slight overestimation of the TL, particularly at high frequencies. This overestimation is due to sound intensity levels outside the virtual surface Σ v that are thus not taken into account in the calculation of the radiated power. The In both case, the points were positioned at z v /2. Note that in the reciprocal approach, the direction of the active sound intensity is defined by the direction of the force injected by the dipole (see Fig. 2(f)). To determine the active sound intensity at point x using the reciprocity principle, the plate was excited successively by a monopole and dipole source at point x.

To reduce the number of excitation points and the measurement time, the symmetries of the system (with respect to x = L x /2 and y = L y /2) were considered. Only the points belonging to a fourth of Σ v , as illustrated in Fig. 9(c), were considered, leading to a total of 40 positions of excitation as compared to 160 in Fig. 9(b). The experiment was performed in an anechoic room, using two translating robots to automate the process (see Figs. mesh.

As a reference to compare the following numerical and experimental results, a numerical model has been established to accurately evaluate the radiated power from wall pressure and wall velocity values. At the wall, the particle velocity sensitivity function H v 0 normal to the plate is equal to the plate velocity sensitivity function H v . Expressing it in the wavenumber domain allows using the wavenumber relation between the acoustic pressure and the particle velocity [START_REF] Williams | Numerical evaluation of the Rayleigh integral for planar radiators using the FFT[END_REF] and thus determining H p at the wall. Finally, by expressing Eqs. ( 4), ( 5) and [START_REF] Maxit | Prediction of flow induced sound and vibration of periodically stiffened plates[END_REF] in the wavenumber domain and using H p and H v 0 determined at the wall, one obtains the radiated power directly at the wall.

Fig. 10(a) shows experimental results obtained based on the reciprocity method, while considering the virtual surface Σ v described at the beginning of this section and illustrated in Figs. 9(b) and 9(c), versus simulation results for which the radiated power was evaluated at the panel surface. The curves match very well on the whole frequency range. This demonstrates that the proposed reciprocity approach accurately reproduces the theoretical TL for the considered test case. Moreover, it validates the definition of the virtual surface Σ v considered in the experiment.

The experimental results obtained with the reciprocity principle are compared in Fig. 10(b) to the experimental results obtained in the transmission loss facility. The transmission loss facility results are noisier than those derived from the reciprocity principle. A good agreement is however noticed in general. Relatively small differences are noticeable above the acter of the sound field in the reverberant room. On the other hand, below 410 Hz the TL measured in the transmission loss facility is on average 6.5 dB lower than the one obtained with the reciprocal approach. This is explained by the modal behavior of the reverberant room below the Schroeder frequency, which enhances the non diffuse character of the incident sound field.

VIII. CONCLUSION

In this paper, a methodology for characterizing the response of flat panels to a diffuse acoustic field excitation without using a reverberant room was proposed. This approach is based on the mathematical formulation of the random excitation problem in the wavenumber domain. This formulation indicates that the panel's response at point x (on the panel or in the acoustic medium) to a random field depends on two quantities in the wavenumber domain: the wall-pressure cross spectral density function of the excitation and on so-called 'sensitivity functions' at point x which characterize the panel. Using the reciprocity principle, it has been shown that these functions can be determined from the panel velocity field in the wavenumber domain when the system is excited by a source of unit amplitude at the point of interest x. The sensitivity functions can be estimated easily by experiment based on the reciprocal interpretation.

The proposed approach avoids the use of a reverberant room to determine the sound transmission loss factor and vibration response of plane panels under a diffuse field excitation.

As the excitation is represented by an analytical model, this approach can be applied to experimentally characterize the vibroacoustic response of a panel to an ideal diffuse acoustic field. It should however be underlined that the main limitations of the proposed approach rely on the assumptions of the mathematical formulation of the problem: the system should be linear (i.e., elastic material, small deformations) and time invariant, and the condition of a baffled panel in an anechoic environment should be verified (particularly for acoustic applications). It offers however a large field of applications.

From a practical point of view, the vibration response of a panel to a diffuse field excitation can be easily estimated using a mechanical source of effort and a vibration measuring device to determine the vibratory response of the panel (in this study, a shaker and a scanning laser vibrometer were used). In this particular case, the use of a baffle is not crucial, because no acoustic excitation are considered.

For acoustic applications, the condition of baffled panel in an anechoic environment is fundamental to agree with the assumption of blocked pressure. Monopole and dipole sources are required to determine the radiated pressure and the particle velocity, respectively. It can be quite challenging to experimentally reproduce those sources. In this study, a monopolelike source was used and moved from a certain distance to represent two monopoles close to each other and out of phase. A dipole source could thereby be reproduced. For the determination of the sound transmission loss factor, the monopole and dipole excitations have to be applied on several points discretizing a virtual surface, which encloses the whole panel. For each position of the excitation, the vibration response of the entire panel has to be measured. In this study, two translating robots were used to automate the process. As opposed to measurements in a transmission loss facility, this experiment was highly timeconsuming. However, with the recently developed vibration measuring techniques (such as optical measurement), the time of experiment could be largely reduced.

To conclude, the method has been validated numerically and experimentally for the considered test case. Comparisons of numerical and experimental results have shown that the sensitivity functions have been well estimated both inside and outside the acoustic circle in the wavenumber domain. A good agreement between numerical and experimental results has also been obtained whether for the velocity spectrum at a point on the panel or for the sound intensity spectrum at a point in the acoustic domain. An application of the proposed methodology for estimating the sound transmission loss of the plate has been presented and the results have been compared with standard measurements in a coupled room facility. In the near future, the method will be extended to the characterization of panels excited by a turbulent boundary layer.

(1) (x) + k 2 0 p (1) (x) = 0, (A1)
excited by the pressure distribution P (x), x ∈ Σ p is then obtained with:

v (x, ω) = jω m,n F mn φ mn (x) M mn (ω 2 mn -ω 2 + jηωω mn ) , (B4)
where the modal force F mn is defined by

F mn = Σp P (x)φ mn (x) dx. (B5)
1. Calculation of the sensitivity functions with the direct interpretation

The direct interpretation described in Sec. III and in Fig. 2 indicates that the sensitivity functions are equal to the system response at point x when the panel is excited by wall-pressure plane waves of wavenumber -k = (-k x , -k y ). The modal force is therefore calculated by considering the pressure distribution P (x) = e -jkx . In this case, the analytical solution to Eq. (B5) is

F mn = I x m I y n , (B6)
where for ξ = x or ξ = y,

I ξ p =        pπ L ξ (-1) p e -jk ξ L ξ -1 k 2 ξ -pπ L ξ 2 , if |k ξ | = pπ L ξ 1 2 jL ξ , otherwise. (B7)
The sensitivity functions H v for a point x on the panel can be directly estimated using Eqs. (B4) to (B7). Based on the velocity response of the plate to a wall-pressure plane wave, the radiated pressure at a point x in the acoustic domain may be calculated using the Rayleigh integral [START_REF] Williams | Numerical evaluation of the Rayleigh integral for planar radiators using the FFT[END_REF] whereas the particle velocity may be deduced from the Euler equation.

Doing so for a set of wall-pressure plane waves allows calculating the sensitivity functions H p and H v0 at point x in the acoustic domain.

Calculation of the sensitivity functions with the reciprocity interpretation

The calculation of the velocity sensitivity function H v at point x using the reciprocity principle involves exciting the plate with a normal unit force at point x. The modal force is thus simply given by F mn = φ mn (x) . (B8)

According to the process described in Sec. V, the normal velocity of the panel should be calculated with Eq. (B4) for points x ∈ Γ x. A discrete Fourier transform is then applied to deduce the sensitivity functions H v (x, k, ω).

Similarly, the pressure sensitivity function H p at point x in the acoustic domain is obtained by exciting the plate by a monopole at point x. The particle velocity sensitivity function H v 0 at point x is obtained analogously using a dipole excitation in the direction n at point x. The modal force is then obtained by approximating the integral of Eq. (B5)

with the rectangular integration rule for a wall-pressure defined by

P (x) = jωρ 0 Q v e -jk 0 r 2πr (B9)
for the monopole case, with r = |x -x| and Q v = 1 m 3 .s -1 . For the dipole case,

P (x) = F 0 d e -jk 0 r 1 2πr 1 - e -jk 0 r 2 2πr 2 . ( B10
)
where r 1 and r 2 are two positions of monopoles separated by a distance d representing a dipole at a distance r [START_REF] Pierce | Acoustics: An Introduction to Its Physical Principles and Applications[END_REF] and the dipole force will be then considered of unit amplitude (i.e., F 0 = 1 N). (a) direct calculation (bold gray line) vs. numerical reciprocity approach (light black line). (b) numerical reciprocity approach (bold gray line) vs. experimental reciprocity approach (light black line).

FIG. 1 .

 1 FIG. 1. Panel (gray line) and coordinate system. (a) receiving side: semi-infinite domain. (b) source side: theoretical DAF.

FIG. 2 .

 2 FIG. 2. Direct interpretation of the sensitivity functions: (a) H v , (c) H p , (e) H v 0 and corresponding reciprocal interpretation (b) H v , (d) H p , (f) H v 0 . (d) and (f) see appendix for demonstration.

  Fig. 2(f).

FIG. 3 .

 3 FIG. 3. Experimental setup. (a) plate excited by a shaker to determine H v . (b) and (c) baffled plate excited by a monopole source to determine H p and H v 0 . 1 -shaker with impedance head. 2 -plate. 3 -frame. 4 -baffle. 5 -sound absorbing foam. 6 -monopole source mounted on 3-axis robot. 7 -single-point laser vibrometer mounted on 2-axis robot.

  3(b) and 3(c)) to avoid potential reflections and possible influence of the background noise generated by the robot used to estimate the TL (see Figs. 3(c) and section VII B). This allowed approaching fully anechoic conditions and ideal monopole and dipole excitations. To estimate the pressure sensitivity functions at point x N of coordinates (x = 0.06 m, y = 0.3 m, z = 0.1 m), the plate was excited by a Microflown Mid-High frequency monopole-HFM source at point x N fed with a white noise signal on the considered frequency range. The monopole source consists of a high impedance loudspeaker connected to a socket by a tube with an inner diameter of 15 mm. The frequency range over which the source is effective and acts like a monopole is [100, 7000 Hz]. The calibration of the source (volume velocity Q v per unit input voltage U ) was obtained by measuring the radiated sound pressure p at a

FIG. 4 .

 4 FIG. 4. Velocity sensitivity functions at point x M , |H v | 2 (dB, ref. 1 m 2 .s -2 ): direct calculation (left), numerical reciprocity approach (middle), experimental reciprocity approach (right). (a) f = 178 Hz. (b) f = 600 Hz. (c) f = 1710 Hz. ---, circle of radius k f . --, circle of radius k 0 .

FIG. 5 .

 5 FIG. 5. Product of sensitivity functions at point x N , |Re{H p × H * v 0 }| (dB, ref. 1 Pa.m.s -1 ): direct calculation (left), numerical reciprocity approach (middle), experimental reciprocity approach (right). (a) f = 178 Hz. (b) f = 600 Hz. (c) f = 1710 Hz. ---, circle of radius k f . --, circle of radius k 0 .

Fig. 6 ( 3 .FIG. 7 .

 637 Fig.6(b) compares the results obtained numerically and experimentally while considering the reciprocity method to estimate the sensitivity functions. A good agreement is observed between the two results, which experimentally validates the proposed methodology for the considered test case. Slight shifts of the resonance peaks in the high frequency range are noticed. They can be explained by small differences between the experimental and the theoretical boundary conditions of the panel or more likely by the added mass from the shaker.

A.

  Panel velocity responseA first experiment in the transmission loss facility was carried out to evaluate the vibratory response of the panel under a DAF. A Polytec scanning laser vibrometer placed on the anechoic side was used to measure the plate velocity ASD function. A 9 × 9 microphones array (1/4 inch B&K 4957) separated by 10 cm in directions x and y was centered to the plate and used to directly measure wall-pressure fluctuations 1 cm away from the plate. An average sound pressure level over all 81 microphones was then calculated to evaluate a mean wall-pressure ASD function Ḡp b p b (f ) on the reverberant side.

FIG. 8 .

 8 FIG.8. Velocity ASD functions G vv (dB, ref. 1 m 2 .s -2 .Hz -1 ): reverberant room measurements (bold gray line) vs. experimental reciprocity approach (light black line).

FIG. 9 .

 9 FIG. 9. Illustration of the virtual surface Σ v considered for estimating the radiated power. (a) in the transmission loss facility. (b) with the reciprocal approach. (c) with the reciprocal approach considering the symmetry properties of the system.

FIG. 10 .

 10 FIG.10. Transmission loss (dB). (a) numerical approach (bold gray line) vs. experimental reciprocity approach (light black line). (b) transmission loss facility (bold gray line) vs. experimental reciprocity approach (light black line).

FIG. 11 .

 11 FIG. 11. Illustration of the three considered vibroacoustic problems of an elastic structure excited by: (a) a normal force at point xf , (b) a monopole source at point x m (internal or external to the structure), (c) a dipole source at point x d (internal or external to the structure).

Fig. 1 .Fig. 2 .

 12 Fig. 1. Panel (gray line) and coordinate system. (a) receiving side: semi-infinite domain. (b) source side: theoretical DAF.

Fig. 3 .

 3 Fig. 3. Experimental setup. (a) plate excited by a shaker to determine H v . (b) and (c) baffled plate excited by a monopole source to determine H p and H v 0 . 1 -shaker with impedance head. 2 -plate. 3 -frame. 4 -baffle. 5 -sound absorbing foam. 6monopole source mounted on 3-axis robot. 7 -single-point laser vibrometer mounted on 2-axis robot.

Fig. 4 .Fig. 5 .Fig. 6 .

 456 Fig. 4. Velocity sensitivity functions at point x M , |H v | 2 (dB, ref. 1 m 2 .s -2 ): direct calculation (left), numerical reciprocity approach (middle), experimental reciprocity approach (right). (a) f = 178 Hz. (b) f = 600 Hz. (c) f = 1710 Hz. ---, circle of radius k f . --, circle of radius k 0 .

Fig. 7 .

 7 Fig. 7. Pressure -particle velocity CSD functions Re{G pv 0 } (dB, ref. 1 W 2 .m -2 .Hz -1 ).

Fig. 8 .

 8 Fig. 8. Velocity ASD functions G vv (dB, ref. 1 m 2 .s -2 .Hz -1 ): reverberant room measurements (bold gray line) vs. experimental reciprocity approach (light black line).

Fig. 9 .

 9 Fig. 9. Illustration of the virtual surface Σ v considered for estimating the radiated power. (a) in the transmission loss facility. (b) with the reciprocal approach. (c) with the reciprocal approach considering the symmetry properties of the system.

Fig. 10 .

 10 Fig. 10. Transmission loss (dB). (a) numerical approach (bold gray line) vs. experimental reciprocity approach (light black line). (b) transmission loss facility (bold gray line) vs. experimental reciprocity approach (light black line).

Fig. 11 .

 11 Fig. 11. Illustration of the three considered vibroacoustic problems of an elastic structure excited by: (a) a normal force at point xf , (b) a monopole source at point x m (internal or external to the structure), (c) a dipole source at point x d (internal or external to the structure).

  FIG.6. Velocity ASD functions G vv (dB, ref. 1 m 2 .s -2 .Hz -1 ). (a) direct calculation (bold gray line) vs. numerical reciprocity approach (light black line). (b) numerical reciprocity approach (bold gray line) vs. experimental reciprocity approach (light black line).

	10 log 10 (G vv ) [dB]	-100 -80 -60 -40
		200	400	600	800	1000	1200	1400	1600	1800	2000
						Frequency [Hz]
						(a)
	10 log 10 (G vv ) [dB]	-100 -80 -60 -40
		200	400	600	800	1000	1200	1400	1600	1800	2000
						Frequency [Hz]
						(b)

ACKNOWLEDGMENTS

This work was supported by the Labex CeLyA of Université de Lyon, operated by the French National Research Agency (ANR-10-LABX-0060/ANR-11-IDEX-0007).

Special thanks must go to Mr Patrick Blachier for his considerable contribution to all

∆p (2) (x) + k 2 0 p (2) (x) = -jωρ 0 Q (2) v δ (x -x m ) , (A2) ∆p (3) (x) + k 2 0 p (3) (x) = -F (3)

where ω is the angular frequency, Q v the volume velocity flow of the monopole source, F

v the dipole source strength (or dipole force) and k 0 = ω/c 0 the acoustic wavenumber. The equilibrium equation for the structure in each case is

where L is the self-adjoint operator of the structure and v is the structure velocity in direction n normal to the structure. Euler's formula provides a relation between the pressure gradient in the fluid and the velocity along the normal external to the fluid. On the surface of the structure, the normal external to the fluid corresponds to -n. In this particular case, Euler's formula becomes

1. Reciprocity principle for the radiated pressure (monopole source)

Multiplying Eq. (A1) by p (2) (x) and Eq. (A2) by -p (1) (x), adding them and integrating them over the entire acoustic domain, one obtains

The volume integral on the left-hand side of Eq. (A8) can be transformed to a surface integral using Green's theorem. Then, using Euler's formula together with Eqs. (A4) and (A5), one finally obtains [18]

From the property of the Dirac delta function (i.e., V f (x) δ (x -x 0 ) dx = f (x 0 ) for any function f defined on V and any point x 0 ∈ V ) one has:

Using Eq. (A7), Eq. (A10) becomes

Eq. (A11) shows that the pressure radiated by a structure at point x m when it is excited by a normal unit point force at point xf equals the structure normal velocity at point xf when it is excited by a monopole source of unit volume velocity at point x m .

2. Reciprocity principle for the particle velocity (dipole source)

Multiplying Eq. (A1) by p (3) (x) and Eq. (A3) by -p (1) (x), adding them and integrating them over the entire acoustic domain, one obtains

As previously, the volume integral on the left-hand side of Eq. (A12) can be transformed to a surface integral using Green's theorem. Then, using Euler's formula together with Eqs. (A4) and (A6), one finally obtains [18]

Using the previously described property of the Dirac delta function and the property of the distributional derivative of the Dirac delta function (i.e., V f (x

for any function f defined on V and any point

x 0 ∈ V ) Eq. (A13) becomes

Using Euler's formula, the normal velocity of the structure at point xf and the particle velocity v n d in direction n d at point x d are introduced in Eq. (A14). One finally obtains