Current data usage and short-range forecast performance Incremental 4D-Var, description

1. Developments of ECMWF’s Data Assimilation SystemWith Respect to Higher-density Observations and Higher ResolutionErik Andersson + coworkers Current data usage and short-range forecast performance Incremental 4D-Var, description Accuracy and efficiency The upgrade of 14 Jan 2003. Current ongoing developments Conclusions

2. Satellite sensors available for NWP(from R.Saunders, Met Office) In the incremental formulation the cost function J is expressed in terms of increments xwith respect to the background state, x=x-xb, at initial time. Hi and Miare the TL of Hi and Mi , linearized around x(ti)=Mixb(t0).

3. Number of data used per 12-hour cycle,in the ECMWF operational system (*106)

4. SYNOP Surf.Press, Wind-10m, RH-2m AIREP Wind, Temperature SATOB AMVs Meteosat, GMS, GOES, MODIS DRIBU Surf.Press, Wind-10m TEMP Wind, Temp, Humidity profiles DROPSONDE Wind and Temp profiles PILOT, Am+Eu Profilers Wind profiles PAOB Surface pressure proxy NOAA-15/16/17 HIRS and AMSU-A radiances DMSP-13/14/15 SSMI radiances Meteosat-7, GOES-8/10 Water Vapour radiances QuikScat Ambiguouswinds SBUV Layer ozone GOME Total ozone In preparation: AIRS, AMSU-B, MSG, SSMI/S, Cloud and precipitation data… Used Data, Since Jan 2003 Conventional Satellite

5. Time series of forecast scores, N Hem.

6. Time series of forecast scores, S Hem.

7. Time series of +48h fc scores,Northern Hemisphere, 500 hPa Z WASHN BRAKL ECMWF ECMWF Upgrade

8. Time series of +48h fc scores,Southern Hemisphere, 500 hPa Z BRAKL WASHN ECMWF ECMWF Upgrade

9. Observing SystemExperiments(G. Kelly et al.) 500Z, N.Hem, 89 cases NoSAT= no satellite radiances or winds Control= like operations NoUpper=no radiosondes, no pilot winds, no wind profilers 500Z, S.Hem, 89 cases

10. New BG-errorsfrom an Ensemble of 4D-Var assimilations (M. Fisher) Temperature (K) U-component (m/s) Geopotential (m) (Andersson et al. 2000: Diagnosing Bg-errors for observable quantities, QJRMS.)

11. New BG-error correlations (temperature)from an Ensemble of 4D-Var assimilations

12. Forecast model at T511 (40 km) resolution Observation minus background departures are computed using the full model at full resolution at the observation time. Analysis increments are computed at coarser T159 resolution (125 km), using a tangent linear forecast model and its adjoint. All observations are analysed simultaneously. 12 hours worth of global obser-vations are used in one go. Around 1 500 000 data are used, in total, per 12-hour cycle. Satellite radiances are the most numerous data source The current operational 4D-Var system

13. A few 4D-Var Characteristics • 4D-Var finds the 12-hour forecast evolution that best fits the available observations • It does so by adjusting 1) surface pressure, and the upper-air fields of 2) temperature, 3) wind, 4) specifichumidity and 5) ozone • The control vector has 7,900,000 elements (T159). All observations within a 12-hour period (~1,500,000) are used simultaneously, in one global (iterative) estimation problem

14. 4D-Var incremental formulationCourtier et al. 1994 In the incremental formulation the cost function J is expressed in terms of increments xwith respect to the background state, x=x-xb, at initial time. and are the TL of and , linearized around . The i-summation is over N=25 ½h-longsub-divisions (or time slots) of the 12-hour assimilation period. The innovations di are calculated using the non-linear operators, Hi and Mi . This ensures the highest possible accuracy for the calculation of the innovations di, which are the primary input to the assimilation!

15. The inner iterations 1) The tangent-linear approximations: and 2) Approximations to reduce the cost: this involves degrading the tangent-linear (and its adjoint) with respect to the full model. • Lower resolution (T159 instead of T511), • Simplified physics (some processes ignored), • Simpler dynamics (e.g. spectral instead of grid-point humidity). This results in cheaper TL and AD model integrations during the minimisation - i.e. the inner iterations.

16. The outer iterations After each minimisation at inner level; • is updated: , • and are re-linearized around . • Innovations are re-calculated using the full non-linear model : • Superscript represents the outer iterations. • The full model remains at T511 throughout.

17. Test of incremental approximation(Y. Trémolet 2003, TM 399) • In 4D-VAR the perturbation is not any vector, it is an analysis increment. It is not random and it is the result of a algorithm which involves the linear model. • The linear and non-linear models are used at different resolutions (T511/T159) • The non-linear model uses more physics. • Humidity is represented in spectral space in the linear model, in grid point space in the non-linear model. • Relative error: vs.

18. Evolution of TL model error Operational configuration (T511/T159) The error is large. It grows very rapidly in the first hours. This is not the case in the adiabatic test.

19. TL model resolution T511 outer loop, 12hour. Varying inner loop resolution. The resolution of the inner loop may have reached a limit, at T159?

20. TL model resolution. The adiabatic case. Adiabatic test. Better TL physics is needed at high resolution This is expensive both in development work and CPU.

21. Spatial scales • TL at T255 • 12h integrations • Initial-time increment truncated at varying resolutions • The relative error grows fastest at the smallest scales

22. Hessian eigenvector preconditioningM. Fisher, (Fisher and Andersson, TM 347) The optimal pre-conditioner for the 4D-Var minimisation problem is the Hessian of the cost function, . The full 4D-Var Hessian is not known. So farhas been used as an approximate preconditioner, neglecting the observation term. The consequence is that patches of very dense or particularly accurate observations may deteriorate the conditioning and slow down the rate of convergence.

23. Diagnosing slow rate of convergenceAdding dense ATOVS data (Andersson et al. 2000, QJ) Four possible solutions: 1) Reducing b in the stratosphere 2) Thinning the ATOVS data 3) Modifying R to account for horizontal correlation of observation error 4) Pre-conditioning Temperature. Leading Hessian EV with dense ATOVS. Xsection (top). Lev10 (below) Lev60, reduced ATOVS density

24. Computing leading Hessian eigenvectors Conjugate-gradient/Lanczos method – M.Fisher (pers com.) The close connection between Conjugate Gradients and the Lanczos algorithm allows us to simultaneously: Minimize the costfunction Calculate the K leading eigenvectors kand eigenvalues kof the Hessian. Spectrum of k in ECMWFs 4D-Var, in a ‘multi-resolution’ incremental test-suite - kand k are calculated at low resolution, and used to pre-condition at higher resolutions Note C<400 !!

25. Preconditioning Convergence is roughly twice as fast with Hessian preconditioning, The Conjugate-Gradient minimisation algorithm is used.

26. The revised 4D-Var algorithm: Specification • Quadratic inner iterations.Variational quality control and SCAT ambiguity removal moved to outer-loop level. • Conjugate Gradient minimisation.With objective stopping-criterion based on the gradient-norm reduction. • Hessian eigenvector pre-conditioning.Updated after each inner minimisation. • Multi-resolution incremental, T95/T159.With some tests at T255. • Interpolation of the trajectory.From T511 to T95/T159. Inspired by discussions with J. Nocedal.

27. Conjugate-gradient: Reduction of Norm of gradient 0.05 T95 T159 T42 With C.G. minimisation the gradient norm reduces nearly monotonically with iteration. It is therefore possible to introduce an objective stopping-criterion based on its ratio. We have chosen a value =0.05.

28. Multi-incremental: RMS of T analysis increments Most of the total An-increment is formed at T42. There is a clear scale-separation between successive minimisation. The rapid decrease beyond ~T100 is due to the filtering properties of Jb, and the lack of observational information on smallest scales.

29. Testing at T255. Surface pressure increments. T159 T255 Scores over a 34-day period were very slightly positive, nut not enough to warrant the extra expense.

30. Implementation on the 14 Jan. 2003: Conj. Gradient minimisation Hessian pre-conditioning Inner/outer iteration algorithm Improved TL approximations Multi-incremental T95/T159 These developments will help facilitate: Use of higher density data Higher resolution (T255) Enhanced use of (relatively costly) TL physics Cloud and rain assimilation The revised 4D-Var solution algorithm More work is needed to improve the representation of the smallest scales in the inner loop.

31. Data Assimilation: Use of the Omega-equation and Non-linear balance in Jb Jb statistics based on 4D-Var ensemble More selective Jc (DFI for Div only) Direct assimilation of SSMI radiances Model: Improved cloud-scheme numerics Revised cloud physics Revised convection scheme – solved the “North America Problem” Other ingredients in the 14-Jan-03 update On the 4 March 2003, the operational system was successfully migrated from the Fujitsu to the IBM clusters.

32. The ‘North America problem’ solved This energetic behaviour caused serious problem in 4D-Var, and it affected FC performance in spring and summer

33. Results from parallel testingof the 14-January upgrade, 255 cases, 500 hPa Z Europe N.Hem N.Atl. S.Hem

34. Results from parallel testingof the 14-January update, Tropics, 255 cases 200 hPa Temp 200 hPa wind 850 hPa Temp 850 hPa wind

35. Ongoing developments in data assimilation • Assist in implementation of data from new satellite systems (e.g. AIRS, MSG, AMSU-B, SSMI/S …) • New humidity analysis formulation (Holm et al. TM383) • A more accurate description of humidity Bg errors • It removes our problematic spin-down in tropical rainfall • OSEs with all main types of satellite humidity data • Model Error (weak constraint 4D-Var) research (Y. Tremolet) • O3 and CO2 assimilation • OSEs (Cardinali et al. TM371) • Preparation for cloud and rain assimilation (TM 383) • Implement more accurate TL physics • Shorter assimilation window? • RRKF, further flow dependence in Jb (Fisher and Andersson TM 347) • Diagnostics • Sensitivity w.r.t. observations, • Realism of SVs and Forecast sensitivity patterns • Information content studies (Fisher 2003, TM 397) • Improved surface analysis (snow, emissivity)

36. The new humidity analysis (E. Holm)Implied BG-error for Relative Humidity 540 hPa The new humidity variable is a ‘normalized relative humidity’, with asymmetric p.d.f. at zero and saturation, see TM383 Holm et al. 2003.

37. The New Humidity AnalysisRemoves spin-down in tropical precipitation Hydrological budget, Sea, New Hydrological budget, Sea, Control

38. Difference in RMS of FC-error, Z 500 hPa31 days of August 2002 NoHIRS NoSSMI NoGEOS AddAMSUB

39. The revised 4D-Var solution algorithm has been implemented, for higher accuracy and efficiency: • Higher density data (preparing for AIRS) • Higher resolution increments (towards T255) Conclusions Our current emphasis is • 1) The humidity analysis: • A new formulation has been developed. • Spin down in tropical precipitation has been cured • A large effort on cloud and rain assimilation is underway – TL moist physics + ‘RTrain’ • 2) Jb developments • Flow dependence • Regional variation (wavelet) • Ensemble techniques • 3) Model Error and Biases