© B. Amiot et al., Published by EDP Sciences, 2026

1 Introduction

Numerical modelling underpins the three fundamental stages of photovoltaic (PV) projects: resource assessment, power performance assessment, and production forecasting. Therefore, various energy yield models, or irradiance-to-power chains, have been designed to handle the conversion from meteorological inputs to electrical production and they are now implemented in PV-oriented software or tools like PVSyst [1], PVGIS [2] or pvlib [3]. However, the selection of sub-models and appropriate hypotheses within the chain remain open questions and critically influence the accuracy of the power production outputs [4]. This impact is even more significant for the solar forecasting stage [5], especially for large‑scale, grid‑tied PV installations where both the electrical and atmospheric scientific communities converge on this problem with distinct backgrounds and various numerical tools [6]. The irradiance- part of the chain has probably received the most attention, first due to foundational researches on atmospheric composition and radiative transfers [7,8], and then in recent years, when it has become an ideal ground for data-driven and statistical approaches [9]. Recent advances considering the full chain have also leveraged probabilistic and statistical methods using physically-based approaches [10]; highlighting the dominant influence of sub‑models that handle irradiance separation and transposition [11]. Interestingly, solar forecasting remains generally a one way prediction (i.e., predicting how the environment acts on the PV modules), whereas the prediction of environmental conditions arising from PV‑atmosphere interactions remains less addressed [12,13]. Although this lack of interest has not been a major bottleneck for the development of conventional land-based PV sites, it becomes consequential for eco-photovoltaic applications where panels impact both the biotic and abiotic ecosystem components like in agrivoltaism [14,15] and in floating photovoltaics (FPV) [16,17].

In eco‑photovoltaic systems, PV panels generally lower the temperature of the bottom system (soil [18], lake [19]) by attenuating short‑wavelength radiation. Because panel motion is an intrinsic characteristic of many eco‑photovoltaic installations, quantifying the consequences of such movement constitutes a first methodological challenge. In agrivoltaic systems, dynamic tracking is dictated by crop conditions and is controlled by the PV operator [20,21]; consequently, a variety of tracking and anti‑tracking strategies can be implemented [22], each exerting a pronounced influence on agricultural yields [23] and on energy yields. Floating PV technologies experience a comparatively modest environmental impact from panel motion because they are generally positioned closer to the bottom/water surface. However, the resulting non-uniform irradiance field from the PV movement substantially affects the PV panels themselves as it intensifies the electrical mismatch, an effect that depends on both the configuration of the system and the geographic location [24]. Therefore, various system-level control strategies can be deployed to optimise energy output [25]; reflecting these controls within the irradiance-to-power chain being eventually a next challenge. Beyond shading, PV panels also modify the surrounding environment aerodynamically and the thermal radiation forcings [16,17]. These effects are summarised in Figure 1, in which ${\Phi }_{\mathrm{\Sigma }\leftrightarrow m}^{sw}$ (resp. ${\Phi }_{\mathrm{\Sigma }\leftrightarrow m}^{lw}$) is the net radiation from short wavelength (resp. long wavelength) regions of the spectrum that are decomposed in the contribution from the atmosphere and the bottom surface to the PV module and vice versa, and where the aerodynamic resistances R_m, R_a and R_w are the inverse of the transfer rates (heat or mass convection, depending on the nature of the scalar field Y). The energy balance of the bottom system is therefore affected by the characteristics of PV panels (e.g., height, tilt, technology) and bottom system factors (e.g., soil moisture, waterbody depth), with significant variations over seasonal and diurnal periods. Although bottom system temperature is typically lower during the summer season, changes in long-wavelength forcings may reverse this trend during winter periods [26]. Therefore, a numerical analysis of these interactions requires tools capable of resolving atmospheric and photovoltaic processes at a detailed, granular scale. Computational fluid dynamics (CFD) solvers are particularly suited for this purpose, having been employed in studies of PV/environment interactions [27–29], water evaporation processes in FPV contexts [30] and soil temperature/moisture dynamics in agrivoltaics [31]. In addition, the recent availability of high-quality datasets [32] has opened up the potential to study these interactions for realistic large-scale solar farm arrangements with various geometry and array layout [33,34]. However, the computational cost of CFD can be prohibitive for forecasting tasks (but not impossible to deal with, see [35]), and the treatments of the boundary conditions are not straightforward, since the responses of PV panels, water bodies, crops, and soils to atmospheric forcings are highly interdependent.

Insight into boundary conditions can be drawn from the numerical representations of land-surface exchanges adopted in global circulation models (GCM) and their extension to numerical weather prediction (NWP) tools, like the Weather Research and Forecast − WRF [36] or the Applications of Research to Operations at Mesoscale − AROME [37]. In both of these tools, land surface models digitalise surface exchanges [38,39] and therefore integrate a bi-directionnal communication between the atmospheric system and sub-systems represented by various sets of equations. For instance, if the sub-system is a photovoltaic element, then it is classical to use the Evans approximation [40] as in [41]:

$$\begin{array}{ll}{P}_{out}=\hfill & 𝒜{\phi }_{\mathrm{\Sigma }\to m}^{sw}{\eta }_{STC}\left(1-\beta [T_m-25]\right),\hfill \end{array}$$(1)

where 𝒜 is the PV module surface, ${\phi }_{\mathrm{\Sigma }\to m}^{sw}$is the solar irradiation in the plane-of-array (Wm⁻²), η_STC is the photoconversion efficiency obtained under standard test conditions, β is the fractional decrease of cell efficiency per unit temperature increase and T_m is the module temperature (°C). Note that photoconverted energy is seen as a heat sink from a surface exchange point of view, and therefore it enables studying the large-scale effects of PV panels [41–44] or eco-photovoltaic systems [45,46]. However, granular effects are necessarily parameterized due to the atmospheric system size in front of the PVs sizes, implying empirical methods [47] or downscalling techniques using CFD [45] to identify momentum absorption from PVs, as well as convective and radiative heat transfers. Moreover, two other computational arguments limit the general adoption of NWP tool for solar forecasting: the computational cost to solve the atmospheric equations is prohibitive for operational usage and there is generally no option to differentiate numerical time-stepping for atmospheric solver and photovoltaic solver. When some efforts must be mentioned on the topic of reducing computational cost through reducing domain size and switching from CPU to GPU architecture [48], NWP for photovoltaic forecasting applications also benefits from dedicated specialisation to deal with different time steps, such as WRF-Solar specialisation [49]. Furthermore, by implementing adapted aerosol effects, aerosol cloud coupling, and cloud radiation feedback mechanisms in the atmosphere, this NWP framework provides superior performance relative to conventional NWP [50], and ongoing improvements to aerosol parameterisation schemes can still increase forecast accuracy [51].

Nevertheless, altering the governing equations and parameterization of PV systems in GCM/NWP models remains a complex task that calls for deep knowledge with these meteorological tools. Furthermore, the electrical and thermal behavior of grid‑tied PV installations depends on the electrical load, thereby adding difficulty to the generation of reliable energy forecasts [5] as it would require a numerical support for grid systems, exacerbating the computational bottleneck. One solution to deal with this degree of numerical adaptability and multiphysics integration is to use a co‑simulation architecture, where specialized models are connected and run simultaneously. When co-simulations found their place among system-based models (i.e., lumped models implementing ODEs or DAEs), they now echo in the atmospheric community as urban meteorology reunites both building systems and microclimate models, requiring to adapt co-simulations to distributed parameter systems (i.e., PDEs).' Where and how communications between models must take place?' is in fact a legitimate double-sided question; and both custom [52,53] and standard frameworks [54] have been proposed as an answer. The main feature of standardized model exchange interfaces is the frame they provide, enabling modellers to focus on high-level codes rather than low-level ones. The Functional Mock-Up Interface (FMI, [55]) is the most common standard which has been applied to various multiphysics applications (e.g., automotive [56], hydraulic [57], indoor ventilation [58]) and seems promising in the context of energy systems [59,60]. Noting the recent effort on the numerical tools supporting distributed parameter systems to incorporate this standard using boundary conditions or source terms in the governing equations, some examples of coupled systems between a lumped model and a fluid model have been achieved while maintaining the original high-fidelity fluid flow computation [61,62]. Even though these examples demonstrate that the frame is compatible with atmospheric-like models, where variables should be exchanged remains application-dependant.

The objective of the work is to propose a co-simulation method based on the FMI standard in order to evaluate the impact of photovoltaics on the environment, along with energy yield and PV module temperatures, in a deterministic and granular fashion. The ambition is to address eco-photovoltaics applications, such as agrivoltaism and floating PV, where factors beyond maximizing energy yield influence PV system control. Based on an irradiance-to-power chain and a 3-D atmospheric CFD model, the exchange of data between these models becomes generalizable owing to the regular geometry of the PV systems and therefore allows one to rely on a standardized protocol instead of in-house data-exchange schemes. The numerical architecture is tested in the context of performance assessment and energy forecasting, where it is of principal interest because of the small size of the CFD model, and specifically for a floating photovoltaic installation. Here, the meteorological forcings are obtained from large-scale weather forecasts provided by third-party entities, which are then downscalled using the atmospheric model based on the CFD implementation. Finally, the particular implementation of the irradiance-to-power model with a simple battery system allows to evaluate the influence of the basic curtailment scenario over microclimate forcings. The paper is organized as follows: Section 2 introduces the FMI and the code architecture that includes the communication structure between the irradiance-to-power chain and the CFD model; Section 3 evaluates the applicability of the FMI strategy in the given forecast context along with a simple performance assessment scenario, while Section 4 explores the possibilities offered by co-simulating PV energy and the environment, especially considering the change in ambient temperature and long-wavelength irradiation when PV energy is curtailed. Future usage and potential limitations of the proposed FMI strategy for solar energy meteorology are also discussed.

Fig. 1 Illustration of the physical components involved in land surface model for floating photovoltaics: radiative components (Φsw, Φlw), aerodynamic resistances R and the electrical production Pout. Impacts from floats are neglected here.

2 Method and material

In this section, a generalizable co-simulation structure for PV/atmosphere interaction is proposed as illustrated in Figure 2a. It includes an orchestration system (Sect. 2.1) to coordinate data exchange between irradiance-to-power (Sect. 2.2) and environmental models (Sect. 2.3). (Sect. 2.4) and (Sect. 2.5) describe a simplified version of the co-simulation to PV forecasting at a specific floating photovoltaic location. The framework is illustrated in Figure 2b and implies large-scale forecasts from a third-party provider, with the co-simulation kernel enhancing local meteorological-PV interaction predictions.

Fig. 2 Illustration of the co-simulation structure for forecast purposes.

2.1 Orchestrating data exchanges using a co-simulation environment and the functional mock-up interface

The Functional Mock‑up Interface is a tool‑independent standard that defines how dynamic system models are packaged and exchanged for co-simulation purposes. The functional mock-up unit − FMU, is the minimal element of the co-simulation and it bundles the compiled binary for a given numerical model, together with an XML file that describes its variables, interfaces, and metadata. Exporting a numerical model as an FMU can be performed following the standard guidelines [55], and many simulation software implement this operation (e.g. Simulink^® [63], Openmodelica [64]). After exporting one or more numerical models as FMU(s), a general-purpose simulation environment can be parameterized to orchestrate the FMU(s) instantiation, the data exchange, and the coupling scheme of the co-simulation (e.g., Python using PyFMI [65], Julia using FMI.jl [66]). Two principal coupling options are available: the FMU for Co-Simulation − CS-FMU for which each FMU integrates its own numerical solver and the FMU for Model Exchange − ME-FMU which implies a tight coupling of numerical solvers through the simulation environment. Note that the coupling scheme must be defined both at the simulation environment level and during the model export as FMU stage. For this proof-of-concept, only two CS-FMUs are considered and form the core of PV/environment co-simulation: a power-to-irradiance model and an atmospheric model. To this end, and following the schematic diagram of heat transfers illustrated in Figure 1, a minimal set of exchanges can be numerically reproduced. Note that thermal energy conservation is assumed to be solved in both models, and the choices made about this matter are shown in Section 2.2 and Section 2.3.

The heat transfer between systems is managed by the simulation environment; this is described by the structure of the orchestrator presented in Figure 2a. In order to obtain a generalizable framework for this orchestrator, a first requirement is to identify the FMUs boundaries over which data exchange will take place. These boundaries can be assumed from the physical and numerical frontier of the given numerical models. Considering two models for microclimate and PV systems with distinct computational domains, it is logical to treat the PV surfaces as the intersection of these domains. Given that the tightness of the PV module measure approximately 5 × 3 m (excluding module frames), while the characteristic lengths for the front and rear surfaces are greater than 1, it is reasonable to consider the intersection of computational domains only on the front and rear surfaces. As a result, a generalisable framework for co-simulation involves establishing fourteen endpoints matching the terms indicated in Figure 1, and relative to the PV surfaces. Figure 3 illustrates this fundamental model-pairing stage, where fr and re denote, respectively, the front and rear module surfaces, a is the atmosphere, w is the bottom surface and CHTC is the convective heat transfer coefficient. Note that even if the connection between these endpoints is managed by the orchestrator using in-memory storage, each FMU must declare the memory location for its own endpoints. Another important element for generalising the method holds in the orchestrator potential to interact with the transmitted data, which allows one to make the bridge among models that do not expose all variables. For example, it may be more affordable from the microclimate case to output a spatially averaged surface temperature instead of a radiative heat flux, which is then processed as a radiative heat flux using a radiosity method during data transmission.

A second requirement for the orchestrator is to arrange the temporality of data exchanges. Given the two computational domains where specific constitutive equations are solved, Figure 4 shows a simple coupling scheme that allows one to compute quantities of interest for both microclimate and PV models. After loading the external large-scale meteorological data in both models (i.e., initial conditions), the serial explicit bidirectional scheme consists in:

• Running the PV model for a single time iteration.

• Transferring the computed data to the microscale model, with potential data adaptation at the simulation environment level.

• Running the microclimate model for a single time iteration.

• Transferring the computed data to the PV model, with potential data adaptation at the simulation environment level.

These four sub-iteration stages are applied at each global time step t_n of the co-simulation where n stands for the global iterator.

When these two requirements are met, it is still possible to simplify the data-mapping and coupling schemes based on:

• The physics solved by the microclimate solver (e.g., convection and radiation, convection only).

• The type of explicit (i.e., geometrically resolved) or implicit approach to represent the PV panels into the microclimate domain.

• The structure of the PV model (i.e., lumped capacitance, 1-D thermal model, single-diode model) and the associated hypotheses like empirical coefficients.

• The temporal averaged or instantaneous nature of the co-simulation under investigation.

For forecasting application, a simplified framework for co-simulation is applied (Sect. 2.5) due to the specific irradiance-to-power and microclimate models used in this work.

Fig. 3 Generalisable framework for the co-simulation using fourteen endpoints connecting the microclimate (gray colored) and the PV (orange colored) models.

Fig. 4 Description of the four sub-iteration process for the serial explicit bidirectional scheme.

2.2 Irradiance-to-power chain

The irradiance-to-power chain has two roles in the co-simulation process and the forecasting framework:

• Computing the quantities of interest for the PV module monitoring like the produced electricity and the module temperature.

• Providing entry-points, as heat fluxes, for the microclimate model.

Many alternatives exist to elaborate this chain, providing various amounts of precision; see [11]. As the photoconversion efficiency of PV cells is intimately connected to the cell temperature, it is recommended to construct the chain around a tight coupling between the thermal and electrical models describing the PV systems and the balance of the system. This allows one to close the heat balance at the PV module level, and therefore reduce the uncertainty on the computation of boundary conditions later on. Instead of using two FMUs (one for the thermal part, one for the electrical part), it is suggested to materialize the tight coupling by considering a single FMU that gathers the two distinct physical types, and therefore the resolution is best optimized by the numerical solver. Note that the conversion from global horizontal to plane-of-array irradiation is distinct from the irradiance-to-power chain and is dealt with at the orchestrator level; see Section 2.5.

The irradiance-to-power model is generated on the Openmodelica software in this study because it eases the implementation of multiphysics components and allows to export models to FMU with the CVODE solver. The open-source library Photovoltaics [67] is used to facilitate the implementation of the electrical chain as it provides high-level components (i.e., PV module based on ideal-diode model, a discrete power tracker model, a DC/DC conversion model and a battery storage − BES), and adapted structure to specify the current and voltage output of any PV module. Assuming the PV array as a single big-PV module, the proposed electrical chain is illustrated in Figure 5a which is typical for an off-grid system. In order to mimic a grid-tied system, the battery system can be parameterized as an infinite storage (allowing one to operate at a realistic V_mp), or as a high state of charge level (allowing one to mimic a total-curtailment scenario). It is assumed that the grid voltage and frequency dynamics are not simulated using such electrical approximation, and therefore the co-simulation is not adapted to identify very short phenomena (sub-second scales). However, it should not influence phenomena on larger time scales as well as temporally averaged quantities if a clear sky dominates during the time-averaged window. Note that a post-processing stage can be applied to remove model bias or extrapolate behavior if a state of reference is known. This strategy is used to extrapolate the AC electrical power based on the computed P_out which stands for the DC side only.

Then in order to match the boundary requirement set in Section 2.1, the thermal model is specifically developed to let heat flux terms appear for both the rear and the front sides of the PV module. Here, a photovoltaic system composed of three layers (glass, silicon, and glass) is implemented, taking into account thermal resistance R^cd and thermal inertia C in each layer, as indicated in Figure 5b. The electrical production of the module P_out and the absorbed short wavelength radiation, expressed as ${\mathrm{\Phi }}_{\mathrm{\Sigma }\leftrightarrow m}^{sw}=𝒜\times {\alpha }_{si}{\tau }_{gl}\times {\phi }_{\mathrm{\Sigma }\to m}^{sw}$ where 𝒜 = 2m² and the product of short wavelength absorptivity and glass transmission α_siτ_gl = 0.9, are implemented at the cell layer level. They are prescribed as a thermal heat sink and a thermal heat source, respectively. Then, convection and radiation are prescribed to the front and rear surfaces, they are parameterized based on the inputs from the microclimate model and the orchestrator stages. For the front surfaces, convective heat transfer is expressed as Φ_cv,fr = 𝒜 × CHTC_fr × (T_a − T_fr), while long wavelength radiative energy from the sky is converted to an apparent sky temperature (${\mathrm{\Phi }}_{a\leftrightarrow m}^{lw}=𝒜\times ϵ\sigma \left(T_{sky}^4-T^4\right)$, where ϵ = 0.95 and σ = wm⁻²K⁻²5.67 × 10⁻⁸K⁻⁴). Note that the CHTC and the sky temperature are converted from the microclimate FMU endpoints by the orchestrator; see Section 2.5. For the rear module side, the convective heat transfer coefficient is set constant to 1 Wm⁻² K⁻¹ (see [68]) and the ambient temperature is used to determine the convection (Φ_cv,re = 𝒜 × CHTC_re × (T_a − T_re)) and the radiative heat transfer (${\mathrm{\Phi }}_{w\leftrightarrow m}^{lw}=𝒜\times ϵ\sigma \left(T_a^4-T_{re}^4\right)$).

When the thermal and electrical models are combined into a single FMU, running a simulation with typical meteorological data becomes straightforward in any general‑purpose numerical environment. Listing 1 illustrates how this can be implemented with the FMI.jl library in Julia for a specific weather forcing (e.g., a positive trend in ${\mathrm{\Phi }}_{\mathrm{\Sigma }\leftrightarrow m}^{sw}$ from 0 s and 200 s followed by a constant value until 2000 s). Figure 6 then displays the time evolution of the computed rear‑side module temperature T_m. The calculation is completed in less than 5 × 10⁻¹ s on a single core laptop processor, delivering tight coupling with minimal time-overhead.

using FMI, Plots

fmu_pv = loadFMU("/path/to/
 irradiance_to_power_model.fmu
 ")

solution = simulate(
 fmu_pv,
 (0., 2000.);
 dt = 1.,
 recordValues = ["sw", "T_m"
 ],
 freeInstance=true,
 terminate=true
)

plot(solution)

Fig. 5 Irradiance-to-power chain.

Fig. 6 Result of the plot(solution) command from Listing 1.

2.3 Microclimate model

The microclimate model details how atmospheric flow and scalar quantities interact with the PV system, and should align with the data-exchange process outlined in Section 2.1. Fulfilling the entry points of the irradiance-to-power model also requires a sufficient numerical flexibility from a software perspective as it must support the determination of an arbitrary location for the computation of the CHTCs along with a temperature of reference (T_a), as well as providing a radiative solver for short wavelength and long wavelength irradiation. The open source computational fluid dynamics solver code_saturne [69] is used in this work due to its ability to solve the governing equations for atmospheric flows and for radiative transfers as well as bringing an adapter to export the microclimate model to CS-FMU 2.0. Moreover, the definition of variable causality is simplified as it relies on high-level parameter declaration of the numerical environment.

The code_saturne solver typically serves as a high-fidelity solver for various applications related to atmospheric sciences (e.g., [70–72]) and therefore it is possible to create a microclimate model with a thin spatial discretization allowing one to simulate the flow around PV modules that are explicitly represented in the domain. We call this the geometrically resolved model. Although the generalizable co-simulation would be compatible with this approach as the irradiance-to-power model implements both front- and rear-side boundary conditions, meshing the domain and running the airflow computation will have a significant cost that precludes using the co-simulation in operational context. Therefore, the proposed proof-of-concept implements an implicit representation of PV modules, allowing to simplify the microclimate system while computing the atmospheric forcings for the irradiance-to-power chain. This simplification is materialized by the use of a 1D Cartesian mesh of 10 height that represents the near-surface atmosphere integrating the PV array and the air below and above it. The momentum and energy equations are solved in this domain via the Reynolds-Averaged Navier-Stokes (RANS) method (denoting the averaging operator by (̅.) and the fluctuating part by (¯. ^′)), using a centred numerical scheme for advection and an implicit first-order Euler scheme for time-stepping.

The effect of the PV array on the atmosphere is mainly seen in the energy equation solved in the domain, it reads:

$${\text{c}}_{p,d}\left({\partial }_t\left({\rho }_0{\displaystyle \overline{\theta }}\right)+{\partial }_j({\rho }_0{{\displaystyle \overline{u}}}_j{\displaystyle \overline{\theta }}\right))={\partial }_j\left(\text{k}{\partial }_j{\displaystyle \overline{\theta }}-c_{p,d}{\rho }_0{\displaystyle \overline{u^{\prime }{\hspace{0.17em}}_j\theta \text{'}}}\right)+{\mathrm{\Gamma }}_{\theta }\text{,}$$(2)

where c_p,d = 1006JK⁻¹ kg⁻¹ is the specific heat capacity of air, ∂_t and ∂_j denote the temporal derivative and the spatial derivative in the j direction respectively, ρ₀ is the mass density, u_j is the velocity in the j direction, θ is the potential temperature defined by $\theta =T\times {\left(\frac{p_s}{p}\right)}^{\frac{{\text{R}}_d}{c_{p,d}}}$, where T is the temperature field, p_s = 1 × 10⁻⁵ Pa is the reference pressure, p is the pressure field, R_d = 287.058 JK⁻¹ kg⁻¹ is the gas constant ; k is the thermal conductivity of air, and Γ_θ is a source term. A simplified gradient hypothesis is applied to model the turbulent heat flux ${\displaystyle \overline{{u^{\prime }}_i{\theta }^{\prime }}}=-\frac{{\nu }_t}{P{r}_t}{\partial }_i{\displaystyle \overline{\theta }}$ in which ν_t is the turbulent viscosity and Pr_t is the turbulent Prandtl number set to 0.7. The turbulent viscosity is determined using the k − ϵ LP turbulence model [73]. The coupling of the PV model with the microclimate model relies on the Γ_θ term (Wm⁻³), representing heat lost from the PV module to the air via convection, as illustrated in Figure 7. It is implemented in the numerical cells of the domain where the PV panels would be located if the PV systems were geometrically resolved. Importantly, since convection is averaged per unit of surface, the thermal source term must be adjusted by the global covering ratio (GCR) for the particular PV array under investigation, which is the surface of the PV module divided by the spacing between the modules. The application of cyclic conditions for the lateral surfaces of the domain also requires one to modify the source term applied to each numerical cell in the domain. Denoted ${\mathrm{\Gamma }}_{\theta }^{\text{'}}$, it is calculated to close the energy balance in the numerical domain (i.e., the sum of each ${\mathrm{\Gamma }}_{\theta }^{\text{'}}$ per cell has to be equal to the Γ_θ term). In doing so, the model can reach a permanent state, where the wind velocity and the temperature profiles are constant in time (i.e., infinite PV array situation) as a no-flux condition is imposed at the bottom surface and a Dirichlet condition is prescribed to the top. Finally, the absorption of momentum energy by the PV array is mimicked using a wind direction-dependent roughness length applied to the bottom surface, whereas Dirichlet conditions are prescribed at the top of the fluid-domain using large-scale meteorological dataset.

For the radiative solver, a two-stream approximation based on the Lacis and Hansen parameterization [74] is employed to resolve solar radiation from the top of the atmosphere (z = 11000m) to the ground surface. The approximation consists in separating the ultraviolet-visible band (300 nm to 700 nm) and the solar infrared band (700 nm to 3000 nm), while assuming homogeneous atmosphere layers that take into account various optical properties as described in [75]. The code_saturne allows one to determine the number of radiative layers as well as computing the flow and the scalar quantities which are needed to estimate the absorption and scattering processes from tropospheric gases. When this could also allows to determine the impact from clouds and microphysics phenomena, the proposed microclimate model does not implement this option as computational cost would have been prohibitive for operational purpose. Therefore, atmospheric quantities in the troposphere are assumed from a large-scale meteorological forecast. This explains why the fluid domain only covers a 10 m atmospheric segment near the surface.

Similarly to Listing 1, the launch of the microclimate model as an FMU instance requires as little programming effort as launching the power-to-irradiance model. However, the time-to-solution is larger for the microclimate model and reaches approximately 115 s on a single CPU and for the same 2000 s of simulation. In fact, the code_saturne solver only require around 15 s for computing atmospheric flow and radiative transfers and the latency due to data exchange through the interface becomes a computational bottleneck. When further code optimization could reduce the latency in the communication protocol, some simulation optimization can be explored to hide the latency if a more granular solution of environmental forcing is required (e.g., increasing the numerical complexity/size of the domain of the microclimate model, reducing the number of calls to the microclimate model through the interface).

Fig. 7 Thermal and radiative conditions applied to the fluid domain.

2.4 Integrating weather forecasts into the co‑simulation framework for a specific FPV array

The co-simulation framework is adapted to compute the PV/environment interaction for a single powerplant, but requires large-scale weather data to initialize and force the boundary conditions of each FMU. They can be provided by retrospective weather data analysis, large-scale measurements, or numerical weather prediction models. Here, NWP forecasts issued by the Météo-France meteorological agency for the French area with ≃1.3 km resolution are used (AROME-France [76,77]). In practice, the meteorological fields from 01:00 UTC to 24:00 UTC for day J+1 are collected from the simulation run launched at 09:00 UTC on day J, as shown in Figure 2b. Thus, 24 samples are gathered on a daily basis and before utilizing the co-simulation framework with these extensive forecasts. Retrieving these high‑quality data is however an operational challenge due to the dataset size and network constraints; therefore, a trade-off is made to limit the number of atmospheric layers to be imported. First, the wind velocity components at heights 10 m and 100 m are imported so that the velocity in the fluid domain can be initialized following a neutral atmosphere profile (Monin-Obukhov theory) and the boundary condition is imposed at the top of the fluid domain. Then, the ambient temperature, relative humidity, and pressure profiles must be defined for the fluid and radiative domains. Two or more heights can be retrieved depending on the desired accuracy for the radiative computation. Eight levels are chosen from 10 m to 3000 m to balance optimal accuracy and timely retrieval of large-scale data. Note that more than eight levels would be needed in order to capture clouds, therefore requiring more computational effort; however, our interest here was primarily to test the validity of the model under clear-sky conditions. Finally, a persistence model for aerosol quantities was applied to refine the atmospheric solver. The corrected measurements of optical properties from the open source AERONET database were used considering two wavelengths per band for the spectral integration (namely 440 nm, 870 nm, 1020 nm, and 1640 nm), and selecting the closest data source from the power plant of interest.

A specific FPV site is selected to test the proposed co-simulation approach and the forecasting application. Located in the Hautes-Alpes region in France (44.34 N, 5.86 E and 630 m at altitude), the selected FPV power plant has a nominal capacity of ≃20MWp and a GCR of ≃70%. Figure 8 shows an overview of the power plant that is composed of classical mono-float systems with a narrow space between the water level and the modules (≃30cm). Several quantities of interest were monitored for a period of 2 months from 1st in January 2025 to 28th in February 2025: electrical production at the array scale P_out (industrial operator courtesy, AC-side, 10 min sampling) as well as averaged module temperature (T_m, 1 min sampling) and meteorological conditions (mainly the ambient temperature T_m at 2.5 height, 1 min sampling). Note that module temperatures are spatially retrieved for three positions onto a single PV module. Once curated and resampled over a period of 1 h centered in the middle of hourly intervals, those measurements are used as experimental references to estimate the co-simulation accuracy. Note that the wind speed U_w (2.5 m height), the irradiance in the plane of the array ${\phi }_{\mathrm{\Sigma }\to m}^{sw}$, and the downward thermal irradiance ${\phi }_{\mathrm{\Sigma }\perp w}^{lw}$ were also collected to test other framework scenarios with various FMU arrangements.

Fig. 8 Picture of mono-float system at the FPV site.

2.5 Adopted co-simulation framework for forecasting PV/environment interactions

The general framework presented in Section 2.1 is adapted according to the irradiance-to-power model and the features of the microclimate model. First, the simulation environment is set in the Julia language, and the example provided in Listing 1 is extended to two FMU instances. Then, the fourteen endpoints are reduced to five, which are: T_sky, ${\phi }_{\mathrm{\Sigma }\to m}^{sw}$, CHTC_fr, T_a and Γ_θ. This reduction is applied because of the implicit representation of PV modules in the CFD model and the small tilt angulation of the PV modules coupled with the large GCR of the PV array, hence:

• The radiative fluxes between the PV module and the bottom surface are neglected for the radiant energy conservation in the CFD model (i.e., adiabatic condition and black body assumption at the bottom surface of the CFD domain).

• Reference temperature T_a is shared across both convective heat transfers and rear-side upward thermal radiative flux.

• Rear-side convective heat transfer coefficient has a low magnitude in front of front-side convective heat transfer coefficient.

Figure 9 shows the data communication scheme and the modification brought to the endpoints to maintain the compatibility between the models. With the exception of the ambient temperature, which is directly transmitted from the microclimate model to the irradiance-to-power model, all remaining endpoints are modified in one way or another by the orchestrator (i.e., differences in data exposure). More precisely, the wind velocity U_w and ambient temperature T_a at 2.5 m height are used, as well as the global horizontal irradiance (visible ${\phi }_{\mathrm{\Sigma }\perp w}^{sw}$ and ${\phi }_{\mathrm{\Sigma }\perp w}^{lw}$ infrared) at the surface level. Importantly, decomposition and transposition models are applied to translate the global horizontal irradiance into the plane of array irradiation. The selected models, namely the Erbs decomposition and the isotropic transposition, are implemented by the simulation environment, which is advantageous because it eliminates the need for developing intricate codes in the CFD solver or specific routines in the PV model. Similarly, the orchestrator implements the empirical CHTC relation, selected from previous experiments [45]. Characteristics of the PV module for the irradiance-to-power chain are based on the manufacturer datasheet, see Section 6.1.

Finally, one of the advantages of the co-simulation framework is that it can be applied to simulate the quantities of interest along a given timeline either considering transient simulations or a series of independent steady-state runs. Choosing the appropriate simulation workflow is motivated by the nature of the model encapsulated in FMUs, the large-scale data forcings, and eventual computation acceleration. In the case of day-ahead forecasting where the co-simulation is fed by hourly runs for the J+1 time window, it is more efficient to perform a series of independent steady-state runs as it allows one to parallelize the resolution across several CPUs. Therefore, each simulation is dedicated to a single timestamp along the J+1 timeline, considering a fixed time step of 1 s per iteration for a total time of 2000 s. Note that each run is performed in the middle of the preceding hourly interval, so that the solar irradiation is representative of the averaged solar irradiation for the interval, considering that only clear-sky situations are evaluated in this proof-of-concept.

Fig. 9 Specialised framework for the co-simulation using a restrictive set of endpoints connecting the microclimate (gray colored) and the PV (orange colored) models.

3 Results

This section evaluates the co‑simulation framework's capability to model electrical production, module temperature, and ambient temperature. As illustrated in Figure 10b, three simulation scenarios are considered: the baseline local scenario (BL), adapted for performance assessment which employs the irradiance‑to‑power model driven by on‑site measurements to identify model corrections for subsequent forecast runs; the baseline forecast scenario (BF), wherein the same irradiance‑to‑power model is fed exclusively with forecast large-scale data, omitting any PV‑environment feedback; and the proposed co‑simulation scenario (CS), which couples the irradiance‑to‑power model with a local micro‑climate model, thereby incorporating the PV/atmosphere interaction. Note that in all scenarios, including those in which only the irradiance-to-power chain is launched alone, the FMU standard is used, see Listing 1. For the sake of reducing the uncertainties on the model performance related to cloud passages, hereafter only the clear-sky conditions are considered. Clear-sky conditions are first assessed from the predicted clearness index of the large-scale forecast, and then they are confirmed using the irradiation measurement on the monitored site. Therefore, for the two-month forecast period, 60 hr of the data met these criteria and were retained for analysis.

Fig. 10 Scenarios of simulations and co-simulation, including various external forcings (large scale prediction or local measurements at the industrial site) and FMU instances.

3.1 Validating the FMU of PV system in performance assessment use case

The ability of the FMU irradiance-to-power chain to simulate PV operational metrics is first assessed using optimal input data (Fig. 10b, baseline local scenario (BL)). In doing so, it is possible to identify the shortcomings of the irradiance-to-power chain before introducing uncertainty into the meteorological conditions (NWP forecasts). In addition, and because the irradiance-to-power chain was deliberately made simple, two reference simulations are elaborated using the classical Evans model from (1) assuming $\beta \simeq {\beta }_{V_{OC}}$: a single layer model (SL) and the Ross model (RO). These supplementary models do not involve the FMU standard and are assumed to be representative of the state-of-the-art models for PV surfaces implemented in NWP/GCM. In the SL model, the multilayer thermal model is reduced to only one layer with a homogeneous temperature, whereas the Ross model uses the approximation:

$$T_m=T_a+k\times {\phi }_{\mathrm{\Sigma }\to m}^{sw}\text{,}$$(3)

where the Ross coefficient k = 0.05Km²W⁻¹, as in [41].

First, electrical productions for the baseline local, the single layer and the Ross model are shown in Figure 11a for a clear sky day. Note that all the electrical values are rescaled considering the maximum electrical production measured on-site (AC side), for the current day. First, it is clearly observed that the computed DC output follows the monitored AC output on site, considering a proportional bias. Applying a constant derating ratio of 16% corrects the bias and the corrected/computed AC output overlaps the monitored AC output. As the measured plane-of-array irradiance, which is inputted to the irradiance-to-power chain, takes into account the complete FPV scene (i.e., far-away shading, motions due to the floating nature of the PV system, varying water albedo), the uncertainty in the photoconversion modelling process is reduced significantly for almost all hourly simulations. This correction can be attributed to external parameters that are not simulated within the irradiance-to-power chain FMU: series and parallel resistances, PV module soiling, AC and DC wire losses, inverter conversion rate, ageing; the absentia of series and parallel resistances being hypothesized as the largest provider of bias. Only the first hourly simulation of the day fails in reproducing the AC production level by exhibiting a particular and systematic numerical overestimation that cannot be solely attributed to the list of potential effects included in the constant bias correction. As these discrepancies appear at low irradiance levels, it is hypothesized that the enforced electrical schemes are not adapted to deal with these conditions (i.e., ideal-diode model, ideal battery), and more operational data are needed to identify the root causes. For instance, as monitored inverters require a minimal voltage level to switch on, whereas there are no cut-off limits on the ideal battery system, the implemented conversion stage may consider a higher production level compared to the real scenario. With the exception of these short periods of little level of solar irradiance, the irradiance-to-power FMU is assessed to be a sufficiently robust solution considering the nominal electrical range and slowly evolving irradiation levels over time. The use of the Evans electrical model (SL and RO) also demonstrates that it can be a good alternative compared to the ideal-diode model as the computed DC electrical is close to the measured AC curve with a little bias as well. However, as it computes directly the electrical power instead of the IV curve, it is less adaptable to scenario where the maximum power point is constrained by external parameters (e.g., partial shading, inverter curtailment).

The second operational quantity of interest is the module temperature T_m and the results for the three models are provided in Figure 11b that the module temperature is set to be the rear-side temperature for the FMU-based model so that it can be compared to the measured module temperature. In this context, it is found that the baseline scenario achieves good performance during peak sun hours (T_m = 34.8 °C compared to T_m = 35.9 ° C), considering that PV modules exhibit a spatial uncertainty of around 2.5° C. The largest errors are observed at 08:30 and 15:30 (ΔT_m ≃ 5 ° C), they are attributed to varying wind direction that pollutes the measured module temperature. As the probe is placed around a corner of the floating PV array, the effect of wind is intensified and cannot be reproduced using a big PV module approach (i.e., spatial heterogeneity of temperature). It is also evident that using a thermal model which takes into account the action from convection and radiation in a separate fashion dramatically reduces the prediction bias, as demonstrated by the larger errors produced by the Ross model. In fact, it is simpler to change the Ross coefficient to match the desired temperature output compared to changing the parameters for the multi- a single-layer models; however, in doing so, one assumes to force the magnitude of the photovoltaic/atmosphere interactions whereas the multi- and single- approaches relax these constraints because an energy balance is solved. Interestingly, it is found that the multilayer model offers a slight enhancement compared to the single-layer model (e.g., 0.4 °C at 12:30), therefore, one may consider encapsulating a simple thermal model as an FMU to save computational time.

To emphasize the opportunity offered by the proposed implementation of the FMU, Figure 12 shows an example of the temperature evolution during a maintenance operation at the floating PV site, during which the PV modules were briefly reconnected to the grid. Observed temperatures and spatial variations are denoted Mean and Deviation, respectively. Here, the FMU is launched in a fully transient mode (1 s time-step) between 12:00 and 14:00 (BL), and a unit commitment signal is applied for a period of 35 min from 12:25 to 13:10. This signal is recomputed and interpolated from existing AC power measurements with a frequency sample of 1 meas per 10 min. When the unit is committed, the photoconversion process begins almost instantly, so the module temperature is decreased from 69.5 °C to approximately 64.5 °C. Conversely, when the system is disconnected around 13:10, the photoconversion is interrupted so that heat is rejected at the PV cell level and the temperature increases to 71 °C. The advantage of the FMU architecture is clearly observed here as the unit commitment signal is reproduced using the load and the conversion stage, whereas the SL and the Ross models cannot reproduce this effect with such simplicity. Moreover, the transient thermal model applied in the multi-layer model allows to capture the rapid evolution of the temperature, whereas more simple models would definitely fail in addressing these time-steps.

Fig. 11 Performance assessments for a clear sky day and for three model implementation: the baseline local BL, the single layer SL and the Ross model RO; compared to quantities measured on-site.

Fig. 12 Evolution of the module temperature during a maintenance operation where the PV modules were is briefly reconnected to the grid. Global horizontal irradiation evolving from 800 Wm−2 to 850 Wm−2.

3.2 Recomputing short-wavelength irradiation to enhance forecasting

Although considering an efficient irradiance-to-power chain allows one to optimize the electrical production forecast when the irradiance inputs are precise, solar forecasting remains a key element for the overall forecast performance. This is why the co-simulation implements a radiative scheme within the microclimate model so that it is used to recompute the short-wavelength downward irradiation when launching the co-simulation framework for energy forecasting purposes. The focus here is on clear-sky days and an example of short-wavelength prediction accuracy is proposed in Figure 13. During peak sun hours, the co-simulation framework experiences a reduced level of short-wavelength irradiation in comparison to the AROME forecast, it is closer to the observed short-wavelength irradiation at the industrial site as indicated in Figure 13a. The gain in precision is up to 16 Wm⁻² in this precise case, whereas the largest fluctuations in solar irradiation observed on-site are not captured by both of the co-simulation and AROME forecasts, for instance the differences in forecast and observed irradiation reaches approximately 50 Wm⁻² around 09:30 UTC. In fact, sub-grid effects are not reproduced in the two-stream approximation for the co-simulation framework which limits the gain in performance from a variability perspective (root mean square error metrics for instance); however, the increase in accuracy around sun peak hours, coming from an optimized parameterization for the atmospheric system may participate in reducing the averaged bias. In Figure 13b, the better performance of the code_saturne radiative model over AROME forecasts is evident across the clear-sky dataset, with notable improvements of approximately ≃6 Wm⁻² and ≃15 Wm⁻² in short-wavelength radiation forecasts. This consistent enhancement across various time frames and irradiation conditions suggests that the atmospheric parameterization implemented is the primary factor driving this improvement. In fact, the persistence model applied to parameterize the aerosol optical depth and the eight-level spatial discretization for the thermodynamical properties in the troposphere allows to refine the radiative transfer scheme by better accounting for the ozone absorption and scattering processes.

Knowing that the short-wavelength radiative transfer has a reduced biased in the co-simulation framework, the increase in prediction performances for the electrical production and the module temperature are indicated in Figures 14a and 14b, respectively. Noting that the photoconversion process is primarily controlled by the short-wavelength irradiation, the reduced bias on the radiative field in the co-simulation directly contributes to the little overestimation in the electrical production forecasts CS compared to the baseline forecast BF. This is evidenced by the linear relationships in Figure 14a, indicating that electrical production of CS is generally greater than BF. The variabilities of predictions are approximately the same (R² values), so that it supports the claim that the largest errors are provoked by the subgrid effects, which are neither considered in the co-simulation CS nor in the baseline forecast BF. When it was anticipated that the co-simulation framework would dramatically enhance module temperature prediction, the results presented in Figure 14b demonstrate that the increase in performance is in fact little impacted by operating temperatures. Overall, both co-simulation CS and baseline forecast BF reach the same linear relationships, only the variabilities are affected, and co-simulation shows a better behavior on this point. The environmental conditions are in fact close between the co-simulation and baseline forecast, and they exhibit smooth temporal evolutions: the external heat fluxes are therefore nearly identical and result in module temperatures of similar magnitudes. Therefore, using the co-simulation instead of the baseline forecast option offers an increase in power prediction around 3% which corresponds to around 6Wh and 12Wh per module, mainly because of the gains obtained by the radiative solver implemented in the co-simulation.

Fig. 13 Evaluation of the short-wavelength forecasts for clear sky days.

Fig. 14 Comparison of simulated PV operational metrics with respect to the observed quantities on the floating PV site of reference (clear-sky situations) in the two forecast scenario (co-simulation CS and baseline forecast BF).

3.3 Predicting the air temperature above the PV panels

The heat rejected by the PV panels depends on the characteristics of the PV module and the operating conditions that evolve as a result of large-scale forcing over the PV system. Figure 15a shows a diurnal evolution of the ambient temperature recorded at 2.5 m height, the forecasted ambient temperature at 10 m height by AROME (AROME-BF) and the simulated ambient temperature at 2.5 m height using the co-simulation setting (CS). It is observed that the simulated temperature largely overestimates the observed temperature before 11:30 and after 14:30. This level of discrepancies, up to 5 °C, is correlated with the numerical assumptions made when constructing the microclimate as a 1D atmospheric domain with significant control from boundary conditions. As the top of the domain is prescribed to a Dirichlet using the 10 m estimation from the large-scale AROME temperature, the temperature profile in the case becomes bounded to this assumption. As the AROME forecast is largely overestimating the real temperature during this period, the case amplifies the existing error as the PV modules heat up the atmosphere. Between 12:30 and 14:30, the AROME forecast becomes more relevant, and so does the ambient temperature computed by the co-simulation which is around 1 °C above the temperature at 10 m predicted by AROME. This surplus in temperature is also observed at the industrial site, so it is assumed that the proposed co-simulation method is capable of reproducing the heat effect by using an implicit thermal heat source and a given surface temperature. In contrast, before 12:30 and after 14:30, there are no heat sinks implemented that can stir the temperature down in the microclimate domain assuming that the large-scale ambient temperature can be trusted. Moreover, since the sensible heat transfer from the bottom surface is reduced by the loss in shear stress induced by the PV panels, the temperature of the bottom surface has, in fact, a little influence over the 1-D air column.

Constructing the analysis around the clear-sky periods during the peak sun hours of the days, when the large-scale forecast of the ambient temperature offers a more appropriate set of boundary conditions, Figure 15b shows the capacity of the co-simulation to compute the surplus in air temperature provoked by the PV/atmosphere interaction. From a statistical perspective, the co-simulation reaches a satisfactory linear relationship whereas the AROME prediction underestimates the observed air temperature, demonstrating that the proposed co-simulation approach can descend in scale, closer to the PV system acting as a heat source, from 10 m genuine prediction (AROME), to 2.5 m (CS). However, the relatively low R² value must be signalled, as it appears that the AROME forecast tends to underestimate ambient temperature more when the temperature magnitude is high (around 8 °C in this dataset). This has an impact on the proposed co-simulation framework, limiting the ability of the microclimate model to anticipate the magnitude of the PV heat islands effect, for example. To less rely on the original temperature forecast but keep a 1-D microclimate model, one may consider implementing a more precise land surface model, including radiation balance and potential microphysics phenomena, that are left out in this study, so that the heat sources and sinks may be better reproduced in the microclimate domain.

Fig. 15 Comparison of the simulated ambient temperature with respect to the observation on the floating PV site of reference (clear-sky situations) in the two forecast scenario (co-simulation CS, and baseline forecast AROME-BF).

4 Discussions

4.1 On the effects of energy curtailment on the environment

The co-simulation framework offers a unique opportunity to evaluate the role of the critical operational mode in the environment, for instance, inverter clipping or energy curtailment. Energy curtailment is explored in this section because it is the most impactful situation during which no electrical energy is generated, but heat. To do so, a fictive case where a high state of charge for the BES is imposed and the infinite storage situation is no longer valid (CS-Curtailed case). This forces the PV system to run at V_OC and the day of interest is used to support this numerical exploration. Three quantities of interest are displayed in Figure 16, the module temperature, the air temperature and the long-wavelength irradiation received by the bottom surface (assuming a 100% sky coverage from the PV panels), the co-simulation framework running with infinite storage feature being labelled CS- V_mp. First in Figure 16a, it is found that the temperature surplus is approximately 6 °C around peak sun hours, which echoes with the module temperature variation shown in Figure 12 between 5 °C to 6.5 °C . When the larger global horizontal irradiance during the summer period leads to larger variations of the PV panel temperature, the curtailment effect remains almost as significant during the clear-sky days of the cold season. Note that in these situations of high module temperatures, one may eventually consider to implement the endpoint for the convective heat transfers at the rear side and use a correlation adapted for natural convection. This would lead to reduce the module temperature compared to the proposed fixed CHTC_re = 1Wm⁻²K⁻¹ rate in this study. The effect of energy curtailment on air temperature at height 2.5 m is presented in Figure 16b, along to the temperature level in the nominal case and the genuine forecasted temperature by AROME (at 10 m). Interestingly, the curtailment has little impact on air temperature between 0.2 °C and 0.3 °C around peak sun hours, assuming that the convective heat transfer coefficient is independent of temperature. Therefore, it is likely that the estimated surplus in air temperature due to energy curtailment is a conservative approximation, but would remain in the sub-1 °C range for the nominal magnitude of convection. As a consequence, the heat-up effect on the environment, due to energy curtailment, is relatively weak compared to other contributions such as long-wavelength irradiation, as shown in Figure 16c. In this case, the significant increase in module temperature changes downward thermal irradiation with a variation up to approximately 50 Wm⁻² in this study case. This level may be even greater in cases where the irradiation levels are much higher, and it does not evolve linearly with respect to the temperature variation; however, the rate of the convective heat transfer coefficient would have a converse effect on the thermal irradiation, as the greater magnitude of convection reduces the temperature of the PV module. In addition, the estimation of the thermal radiation effect on the bottom surface assumes a 100% sky-view coverage, which is potentially compatible with floating photovoltaics, but not for agrivoltaics where the covering ratio is generally lower than 40% so the various effects of curtailment would be lower.

It is worth mentioning that the proposed FMU implementation is not restricted to a two-state scenario (i.e., V_mp or V_OC) and one may consider to constrain the MPP tracker instead of the BES state of charge. In doing so, the PV system can be controlled and all situation with various levels of demand-side constraints can be managed. The effects on module and air temperatures, and long-wavelength irradiation would be delimited by the two-state scenario proposed in the study as they are boundary situations.

Fig. 16 Effects from curtailment over three quantities of interest.

4.2 Effect from varying convective heat transfer coefficients

As shown in Figure 11, a limitation of GCM/NWP is that module temperature is generally parameterized (see the RO case) instead of being solved using an energy balance. In contrast, the proposed co-simulation interface offers the opportunity to compute the temperature of PV panels based on distinct parameterization for radiative and convective heat transfers. Hence it is possible to set a convective heat transfer coefficient adapted to the PV geometry, or even determining the appropriate coefficient if one implements a geometrically resolved PV system in the microclimate instance. As the choice was to make the PV module implicit in the computational domain so that it is not possible to use the wall function directly from a PV surface, for example, the only option is to set the CHTC value. Therefore, it becomes a key element in the distribution of heat fluxes in the co-simulation. Initially chosen from past experiments, the CHTC formula cannot be universally applied to all PV installations. Instead, it might be more effective to select or construct a suitable correlation from existing records [78]. In the case of dynamic tracking systems, adjusting the CHTC values based on the installation geometry rather than estimating the module temperature from external relations may allow one to identify the evolution of module temperature and hence electrical output. Therefore, the co-simulation CS has been re-launched assuming two other magnitudes of convective heat transfer coefficients and the results are presented in Figure 17. Note that the case with low coefficient denoted CS-L implements CHTC_fr = 1.5 × U_w + 1.5Wm⁻²K⁻¹ and the high coefficient scenario denoted CS-H implements CHTC_fr = 1.5 × U_w + 7.5Wm⁻²K⁻¹. First, the module temperature is indicated in Figure 17a and the low coefficient simulation provides the highest magnitude of the temperature during the day and the lowest during the night. The maximum difference between the temperature obtained using the lowest CHTC and the highest CHTC reaches 8.2 °C at 11:30, which is of the same order of magnitude as the difference in module temperature between a nominal case and the curtailed case as indicated in Figure 16a. This change in module temperature and the associated magnitude in convective transfer led to an increase in air temperature around the PV installation and Figure 17b shows that the increase in air temperature reaches 0.8 °C when comparing the cases of low and high convection. There is a trade-off between increasing the CHTC, for instance by adapting the geometry of the PV installation, and reducing the impact on the nearby environment. Note that in this case, the CHTC has been de-correlated from the change in momentum absorption that happens when PV geometry is modified (e.g., increase of the tilt angulation). As the PV system protrudes more and more above the ground, increasing the CHTC value, the wind velocity would be reduced due to the larger obstacle sizes. Therefore, the sensible heat rejection may be overestimated in the case of high CHTC, and conversely the sensible heat may be underestimated in the case of low CHTC. Finally, Figure 17c indicates that the short-wavelength irradiation received by the ground surface is principally impacted by the CHTC value, as the difference between the low and the high CHTC values results in 50 Wm⁻². Another time, this change in environmental forcing is approximately equal to the simulated change occurring between curtailed and nominal PV operation. In this two cases, the PV panels increase the long-wavelength irradiation with respect to the sky irradiation so that the ecosystems would still be heated-up through thermal radiative transfers. Although the bottom condition in the proposed CFD model is set constant in temperature due to the floating set up and the hypothesis made, bio-physical model as crop model may be implemented to picture the response from the crops with the modified-by-pv, environmental forcings [31]. In these circumstances, the co-simulation offers an opportunity to anticipate both the energy output due to the tracking configuration, but also this modified-by-pv environmental forcings.

Fig. 17 Effects from curtailment over three quantities of interest.

5 Conclusions

A co‑simulation of photovoltaic arrays integrated into the environment is proposed, allowing us to compute PV operating quantities together with the ambient temperature, in a deterministic and granular way. The co‑simulation follows the Functional Mock‑up Interface standard, and we present a general architecture from data mapping to communication of variables between sub-models, that enables the combination of irradiance-to-power chain and microclimate model, which are containerized as Functional Mock‑up Units to facilitate numerical interoperability. As a proof-of-concept, two basic sub-models are implemented in the co-simulation loop: a simplified 3-D microclimate model (Navier-Stokes equations, radiative transfers and energy conservation) and an irradiance-to-power model represented by a multi-layer thermal scheme combined to an ideal-diode connected to a DC/DC converter with MPP tracking and an ideal battery storage. Appropriate endpoints (i.e., inputs/outputs) for these two sub-models are presented and the management of the co-simulation is held by a simulation environment developed on the Julia language. Following this framework and endpoint scheme should enable adjusting sub-model complexity at the modeller convenience (e.g., increase forecast precision, address various forecast timelines).

Here, a specific implementation for day-ahead forecasting in the context of a floating photovoltaic array is tested. Focusing on clear-sky conditions and using large-scale prevision issued by Météo-France, the microclimate FMU is shown to be appropriate for recomputing the downward short-wavelength irradiation. On average, the solar prediction accuracy was found to be increased by approximately 10 Wm⁻² compared to the baseline scenario. Then, due to photovoltaic / atmosphere feedback, the co-simulation was found to reduce the error on the module temperature compared to state-of-the-art method implemented in comparable meteorological tool, while the prediction of the ambient temperature at 2.5 m above the array was sporadically improved during the peak sun hours. General improvements for these two quantities would require better anticipating the large-scale meteorological fields and implementing a specific land-surface model to the microclimate with respect to the already existing large-scale simulation. In addition, given the small dimensions of the proposed microclimate model, it could ultimately be used to generate probabilistic forecasts with low computational resources.

Using battery storage to mimic a unit commitment profile and a full-curtailment scenario, we identified to what extent the module temperature varies when grid-tied PV systems are disconnected from the grid. Estimated from 5 °C and 7 °C in the proposed PV configuration, this effect can lead to an increase in temperature for systems below the PV panels due to the reduction in net thermal radiative energy, mitigating the effect of shading. However, air temperature is found to be less affected by energy curtailment than by the varying magnitude of the convective heat transfer coefficients that can occur, for instance, when considering dynamic PV tracking. Therefore, using the co-simulation framework with a refined geometry for the PV arrays may offer a good opportunity to explore how the biotic and abiotic components beneath the PV panels are impacted due to the PV operator decisions.

Acknowledgments

Writefull's model (including large-scale language-generation model and language checking) has been used to improve the writing style of this article. B. Amiot reviewed, edited, and revised the generated texts to his own liking and takes the ultimate responsibility for the content of this publication.

We thank G. Bayart and the EDF Power Solutions company for access to the operational data.

We thank Philippe Goloub and his team for the effort in establishing and maintaining the OHP_OBSERVATOIRE site (AERONET Database).

Funding

This research received no external funding.

Conflicts of interest

The authors have nothing to disclose.

Data availability statement

The data that support the findings of this study are available from the corresponding author upon reasonable request.

Author contribution statement

All the authors were involved in the preparation of the manuscript. All the authors have read and approved the final manuscript.

References

All Figures

	Fig. 1 Illustration of the physical components involved in land surface model for floating photovoltaics: radiative components (Φsw, Φlw), aerodynamic resistances R and the electrical production Pout. Impacts from floats are neglected here.
In the text

	Fig. 2 Illustration of the co-simulation structure for forecast purposes.
In the text

	Fig. 3 Generalisable framework for the co-simulation using fourteen endpoints connecting the microclimate (gray colored) and the PV (orange colored) models.
In the text

	Fig. 4 Description of the four sub-iteration process for the serial explicit bidirectional scheme.
In the text

	Fig. 5 Irradiance-to-power chain.
In the text

	Fig. 6 Result of the plot(solution) command from Listing 1.
In the text

	Fig. 7 Thermal and radiative conditions applied to the fluid domain.
In the text

	Fig. 8 Picture of mono-float system at the FPV site.
In the text

	Fig. 9 Specialised framework for the co-simulation using a restrictive set of endpoints connecting the microclimate (gray colored) and the PV (orange colored) models.
In the text

	Fig. 10 Scenarios of simulations and co-simulation, including various external forcings (large scale prediction or local measurements at the industrial site) and FMU instances.
In the text

	Fig. 11 Performance assessments for a clear sky day and for three model implementation: the baseline local BL, the single layer SL and the Ross model RO; compared to quantities measured on-site.
In the text

	Fig. 12 Evolution of the module temperature during a maintenance operation where the PV modules were is briefly reconnected to the grid. Global horizontal irradiation evolving from 800 Wm−2 to 850 Wm−2.
In the text

	Fig. 13 Evaluation of the short-wavelength forecasts for clear sky days.
In the text

	Fig. 14 Comparison of simulated PV operational metrics with respect to the observed quantities on the floating PV site of reference (clear-sky situations) in the two forecast scenario (co-simulation CS and baseline forecast BF).
In the text

	Fig. 15 Comparison of the simulated ambient temperature with respect to the observation on the floating PV site of reference (clear-sky situations) in the two forecast scenario (co-simulation CS, and baseline forecast AROME-BF).
In the text

	Fig. 16 Effects from curtailment over three quantities of interest.
In the text

	Fig. 17 Effects from curtailment over three quantities of interest.
In the text