CovidMIP simulation
The multi-model ensemble simulations from CovidMIP were used to investigate changes in solar PV and wind energy potentials in response to carbon-neutral policies. In CovidMIP, the baseline scenario follows the SSP2-4.5, a medium pathway of future greenhouse gas emissions with a global mean warming of about 2 °C in the 2050s and 2.6 °C in 210055. Accounting for anthropogenic emissions decline during the COVID-19 pandemic and the possible fast renewable transition during the post-pandemic recovery, CovidMIP produced two carbon-neutral pathways from 2020 to 2050: a moderate green recovery (MOD) and a strong green recovery (STR)35. These two carbon-neutral scenarios both include a ‘two-year-blip’ period (2020–2021), which is not the focus of our analysis here, followed by the moderate MOD and faster STR scenarios towards reaching global carbon neutrality by 2060 and 2050, respectively.
Six ESMs participated in CovidMIP, four of which, ACCESS-ESM-5, MIROC-ES2L, MPI-ESM1-2-LR and MRI-ESM2-0 (Supplementary Table 1), are used here, because they made available daily output of surface air temperature (T), surface downwelling shortwave radiation (I) and 10 m wind speed (W, calculated by meridional and zonal winds), which are necessary to quantify solar PV and wind energy, especially their day-to-day variability. All four models provided large ensemble simulations (30-member for ACCESS-ESM-5 and MIROC-ES2L, and 10-member for MPI-ESM1-2-LR and MRI-ESM2-0), which allows a robust quantification of projection uncertainty. In addition, the same number of simulations for each model (as in Supplementary Table 1) are selected from the historical experiment (1995–2014) from CMIP6 to evaluate model skill and correct biases, and to provide a base period to quantify future changes of solar PV and wind energy. In this study, we first calculate renewable energy in individual ensemble members, which is then averaged to each model and the multi-model generation.
Bias correction
To evaluate model performance in simulating solar PV and wind energy during the historical period (1995–2014), we obtain daily (averaged from hourly) T, I and W (calculated by meridional and zonal winds) from the ERA5 reanalysis. ERA5 has been evaluated extensively in earlier studies and is found to be one of the best global reanalysis products56,57, widely used as a benchmark to evaluate and correct model simulations58,59. An early study showed that the ERA5 reanalysis overestimated global I by 4.05 W m−2 (~3%) using site-based solar radiation from the Baseline Surface Radiation Network57. Here, we further evaluate T and W from the ERA5 reanalysis based on global observations at 3,511 weather stations. Compared with observations, the ERA5 reanalysis well captures global T and W with low normalized mean biases of −0.2% and 3.0%, respectively (Supplementary Fig. 9), much smaller than model–ERA5 discrepancy (Extended Data Fig. 6c,d). To maintain the consistency of spatial resolution among the four ESMs, both reanalysis and model outputs are re-gridded to a median resolution of 2 × 2° using the bilinear interpolation method. To improve the robustness of the model projection, we use the multivariate bias correction technique based on the n-dimensional probability density function transform (MBCn) to simultaneously correct daily T, I and W in historical and future model simulations using the ERA5 reanalysis as the benchmark.
MBCn is a multivariate generalization of quantile delta mapping, which considers the dependence among different variables60. In using MBCn, three datasets are included: historical observations (Xobs), historical simulations (Xhist) and projected simulations (Xproj). First, we rotate Xobs, Xhist and Xproj with an N × N uniformly distributed random orthogonal rotation matrix R[j] at the jth iteration:
$$\left\{\begin{array}{c}{\widetilde{X}}_{{{\mathrm{obs}}}}^{[j]}={X}_{{{\mathrm{obs}}}}^{[j]}{R}^{[j]}\,\\ {\widetilde{X}}_{{{\mathrm{hist}}}}^{[j]}={X}_{{{\mathrm{hist}}}}^{[j]}{R}^{[j]}\,\\ {\widetilde{X}}_{{{\mathrm{proj}}}}^{[j]}={X}_{{{\mathrm{proj}}}}^{[j]}{R}^{[j]}\,\end{array}\right.$$
(1)
Second, the quantile delta mapping method uses the same empirical cumulative distribution function (CDF) for simulations (historical and future) and observation, but it preserves the signal of future changes in climate projections61. This method is applied to obtain bias-corrected datasets in historical and projected simulations (\({\hat{X}}_{{{\mathrm{hist}}}}^{[j]}\) and \({\hat{X}}_{{{\mathrm{proj}}}}^{[j]}\)):
$$\left\{\begin{array}{c}{\hat{X}}_{{{\mathrm{hist}}}}^{[j]}={F}_{{{\mathrm{obs}}}}^{-1}\left({F}_{{{\mathrm{hist}}}}\left({\widetilde{X}}_{{{\mathrm{hist}}}}^{[j]}\right)\right)\,\\ {\hat{X}}_{{{\mathrm{proj}}}}^{[j]}={\widetilde{X}}_{{{\mathrm{proj}}}}^{[j]}+{F}_{{{\mathrm{obs}}}}^{-1}\left({F}_{{{\mathrm{proj}}}}\left({\widetilde{X}}_{{{\mathrm{proj}}}}^{[j]}\right)\right)-{F}_{{{\mathrm{hist}}}}^{-1}\left({F}_{{{\mathrm{proj}}}}\left({\widetilde{X}}_{{{\mathrm{proj}}}}^{[j]}\right)\right)\end{array}\right.$$
(2)
where Fhist and Fproj represent the CDFs of \({\widetilde{X}}_{{{\mathrm{hist}}}}^{[j]}\) and \({\widetilde{X}}_{{{\mathrm{proj}}}}^{[j]}\), respectively. \({F}_{{{\mathrm{obs}}}}^{-1}\) and \({F}_{{{\mathrm{hist}}}}^{-1}\) represent the inverse CDFs of \({\widetilde{X}}_{{{\mathrm{obs}}}}^{[j]}\) and \({\widetilde{X}}_{{{\mathrm{hist}}}}^{[j]}\), respectively.
Finally, the bias-corrected datasets are rotated back:
$$\left\{\begin{array}{c}{X}_{{{\mathrm{hist}}}}^{[j+1]}={\hat{X}}_{{{\mathrm{hist}}}}^{[j]}{R}^{{[j]}^{-1}}\,\\ {X}_{{{\mathrm{proj}}}}^{[j+1]}={\hat{X}}_{{{\mathrm{proj}}}}^{[j]}{R}^{{[j]}^{-1}}\end{array}\right.$$
(3)
For historical or projected simulation correction, we repeat the above three steps until the multivariate distribution of \({X}_{{{\mathrm{hist}}}}^{[j+1]}\) or \({X}_{{{\mathrm{proj}}}}^{[j+1]}\) matches that of Xobs.
The MBCn is applied to individual members of each ESM’s simulation, separately.
Calculation of solar PV
Solar PV power yield depends on PV power generation potential (PVPOT) and installed capacity. PVPOT is a dimensionless value, which describes the performance of PV cells relative to the nominal power capacity under actual environmental conditions. Therefore, PVPOT multiplied by the nominal installed watts of PV power capacity is the actual PV power generation. Following previous studies15,43,62, we used daily T, I and W to calculate PVPOT:
$${\mathrm{{{PV}}_{{POT}}}}\,=\,{P}_{{\mathrm{R}}}\frac{I}{{I}_{{{\mathrm{STC}}}}}$$
(4)
where I represents surface downwelling shortwave radiation and ISTC represents shortwave flux on the PV panel under standard test conditions, defined as a constant of 1,000 W m−2. PR is the performance ratio, representing temperature influence on PV efficiency:
$${P}_{{\mathrm{R}}}\,=1+\gamma \left({T}_{{{\mathrm{cell}}}}-{T}_{{{\mathrm{STC}}}}\right)$$
(5)
where γ is defined as −0.005 °C−1 in monocrystalline silicon solar panels, representing the negative impact on conversion efficiency, and TSTC is the cell temperature under standard test conditions (25 °C). Tcell is the actual cell temperature, which is approximated by T, I and W:
$${T}_{{{\mathrm{cell}}}}\,=a1+a2\times T+a3\times I+a4\times W$$
(6)
where a1, a2, a3 and a4 are taken as 4.3 °C, 0.943 (unitless), 0.028 °C (W m−2)−1 and −1.528 °C (m s−1)−1, respectively. These coefficients represent the influence of meteorological conditions on the cell temperature. The ambient T determines the base temperature of the cell, a strong I increases the cell temperature and W decreases cell temperature. These coefficients are found to be fairly independent on site location and cell technology type63, which have been widely used to predict PV cell temperature15,43,64.
Calculation of wind energy
WPD (W m−2) is a typical measure of wind energy potential65, defined as follows:
$${{\mathrm{WPD}}}\,=\,\frac{1}{2}\rho {W}_{{\mathrm{h}}}^{3}$$
(7)
where ρ represents the air density, which is assumed to be a constant value of 1.213 kg m−3 at standard atmospheric conditions, and Wh represents the wind speed at the 100 m hub height.
It is noted that Wh is not available from climate model outputs here. Similar to previous studies22,23,25, Wh is extrapolated from the 10 m wind speed (W) using the wind power law:
$$\frac{{W}_{z}}{{W}_{{z}_{{{\mathrm{ref}}}}}}\,=\,{\left(\frac{z}{{z}_{{{\mathrm{ref}}}}}\right)}^{\alpha }$$
(8)
where Wz represents the wind speed at a height z and \({W}_{{z}_{{{\mathrm{ref}}}}}\) represents the wind speed at a reference height zref. The scaling factor of α, representing how quickly the wind decays towards the ground, is often approximated as a constant of 0.143 over land surface in previous studies22,26,45. As the ERA5 reanalysis provided wind speeds at both 10 m and 100 m, here we estimate α at each location grid (Extended Data Fig. 7) to account for spatial disparity. The higher values of 0.2–0.25 are mainly located in the eastern United States, eastern China, Amazon, India and northern Asia due to large forest coverage, but the lower values of 0.12–0.16 usually occur in flat terrain of desert and steppe. As an improvement to a few previous studies using a constant scaling factor21,22, the wind speed at 100 m here estimated using a spatially variant scaling factor is closer to the benchmark (Supplementary Fig. 5c versus b). The normalized mean bias decreases from −10.1% to −0.4% on global average. In contrast to the large spatial heterogeneity, the scaling factor only shows a small temporal variability (Supplementary Fig. 10 showing seasonal change as an example), resulting in limited benefits on estimation of 100 m wind speed (Supplementary Fig. 5d versus c; −0.4% improved to −0.3%). Therefore, the spatially variant but temporally invariant scaling factor is adopted here to estimate 100 m wind speed from the 10 m wind speed in the model output.
The actual wind power (WP in KW) is sensitive to wind speed and wind turbine. Here, we adopt a typical wind turbine of Sinovel SL3000 ( as an example to describe WP20:
$$\begin{array}{l}{\mathrm{WP}}=\\\left\{\!\begin{array}{lc}0 & {W}_{{\mathrm{h}}}\in [0,\,3]\,({\mathrm{below}}\,{\mathrm{cut}}{\hbox{-}}{\mathrm{in}})\\ & \begin{array}{c}-0.24151\times {{W}_{{\mathrm{h}}}}^{5}+6.9287\times {{W}_{{\mathrm{h}}}}^{4}-74.2354\times {{W}_{{\mathrm{h}}}}^{3}+412.0241\times {{W}_{\mathrm{{h}}}}^{2}\\ -1,049.58726\times {W}_{{\mathrm{h}}}+956.1936\,{W}_{{\mathrm{h}}}\in [3,\,11]\,({\mathrm{ramp}}{\hbox{-}}{\mathrm{up}})\,\end{array}\\ 3,000 & \,{W}_{{\mathrm{h}}}\in [11,\,25]\,({\mathrm{rated}}{\hbox{-}}{\mathrm{power}})\\ 0 & {W}_{{\mathrm{h}}}\in (25,\,+\infty )\,({\mathrm{above}}\,{\mathrm{cut}}{\hbox{-}}{\mathrm{out}})\end{array}\right.\end{array}$$
(9)
where the turbine starts functioning at 3 m s−1, generating the rated-power of 3,000 KW at 11 m s−1 or above, but needs be shut down at 25 m s−1. Therefore, 3 m s−1 and 25 m s−1 are defined as cut-in and cut-out speeds, respectively. The wind speed of 3 to 11 m s−1 is defined as ramp-up condition and the wind speed of 11 to 25 m s−1 is defined as rated-power condition.
Temporal variability of solar PV and wind energy
For both solar and wind energy, we first calculate the daily value using model outputs, which is then averaged to monthly and annual generation over all locations. It is well known that solar PV and wind energy generation are heavily influenced by weather fluctuation, which yields strong variability at various time scales. Understanding the variability of renewable energy is vital for coordinating compensatory energy sources and storage in order to secure a stable energy supply66,67.
Here, we use the metric of NMAD to quantify the day-to-day, month-to-month and year-to-year variability of renewable energy. For any given time series (Ti, i = 1,2,…N), the NMAD (%) is defined as the mean absolute deviation divided by the mean:
$${{\mathrm{NMAD}}}=\,\frac{\mathrm{{{Mean}}}\left[{\rm{|}}{T}_{i}-{{\mathrm{Mean}}}\left(T\right){\rm{|}}\right]}{{{\mathrm{Mean}}}(T)}\times 100 \% .$$
(10)
For day-to-day and month-to-month variability of solar PV and wind energy, we first calculate the corresponding NMAD for each year, and then present the multi-year average of NMAD.
Model evaluation with bias correction
Climate model simulations of solar PV and wind energy resources remain highly uncertain21,36,45. Here, we evaluate the simulated solar PV and wind energy in the historical period (1995–2014) against ERA5 (Extended Data Fig. 6). The observed solar PVPOT shows a smooth global spatial contrast (Extended Data Fig. 6a). In global arid and semi-arid regions, including the western United States, northern Africa, western Asia and Australia, the solar PVPOT shows higher values of >0.24. However, the solar PVPOT is relatively lower in global monsoon regions, including east Asia, south Asia, central Africa, southeastern North America and the Amazon. Such a spatial pattern can be attributed to more clouds in those monsoon regions and possibly denser vegetation cover, leading to less solar radiation available at the ground level68,69,70. Compared with reanalysis, the raw output from model simulation overestimates the solar PVPOT by more than 15% in southeastern Asia, the Amazon and western Africa, but underestimates it by approximately 10% in the western United States, northeastern Asia and western Asia (Extended Data Fig. 6c). Using the MBCn technique to jointly correct T, I and W, the calculated solar PVPOT from the model simulation agrees well with the observation, with a relative bias of less than 1% over all land grids (Extended Data Fig. 6e).
In contrast to the solar PVPOT, the observed annual mean WPD shows very large spatial heterogeneity (Extended Data Fig. 6b). The large values of >160 W m−2 are mainly located in the central United States, Europe, midwestern Australia, northern Africa and northern Asia. The annual mean WPD of 60–100 W m−2 is found in eastern China, India and southern Africa. However, the annual mean WPD is smaller than 20 W m−2 in the tropics, including the Amazon, central Africa and southeastern Asia, where there are light winds due to a small horizontal temperature gradient and large surface friction over vegetated lands. Compared with observation, the raw ensemble mean simulation shows very large biases across the world, with a mean relative bias of 49.4% over global lands (Extended Data Fig. 6d). Regionally, simulated annual mean WPD is overestimated by more than 100% in the western United States, northern Africa, India, Australia and most of Asia, but underestimated by 50–70% in Europe, the Amazon and central Africa. After the MBCn bias correction, the simulated annual mean WPD is close to observation at both global and regional scales (Extended Data Fig. 6f). The global spatial relative bias decreases from 49.4% to −7.5%. In most regions of the world (other than Greenland and southern Australia), the relative bias is smaller than 15%.
The evaluations demonstrate the need of bias correction to improve the fidelity of model simulations. Here, the same bias correction is applied to the ensemble future projection in CovidMIP to investigate the future changes in solar PV and wind energy. This type of multivariant bias correction for climate variables has been done for many climate impact studies71,72,73, but not in previous renewable energy resources analysis. Compared with an early assessment based on traditional bias-corrected climate data24, our study uses the MBCn method, which considers the dependence across climate variables, and represents an advance from the widely used single-variable bias correction.