290 likes | 481 Views
Numerical Weather Prediction on Linux-clusters – Operational and research aspects Nils Gustafsson SMHI. Weather Prediction as a problem in mechanics and physics (V. Bjerknes,1904).
E N D
Numerical Weather Prediction on Linux-clusters – Operational and research aspects Nils Gustafsson SMHI
Weather Prediction as a problem in mechanics and physics (V. Bjerknes,1904) • 7 differential equations: 3 momentum equations, the thermodynamic equation, the equations for conservation of mass and moisture and the equation of state for gases. • 7 model state variables: 3 velocity components, pressure, temperature, moisture and density • Bjerknes suggested that these differential equations and these variables could be used for forecast the future development of the weather • Bjerknes also stated that the calculation problem was impossible to solve
Historical development of Numerical Weather Prediction (NWP) • 1904: The idea of V. Bjerknes • 1922: L.F. Richardson makes the first trial to integrate a 2D problem by manual calculations • 1950: First NWP calculation with a quasi-geostrophic 2D model on a computer by Charney et al. • 1950’s: First operational NWP on the BESK computer in Sweden • 1960’s: Quasi-geostrophic multi-level models • 1970’s: Primitive equation models (hydrostatic equations), development of physical parameterizations • 1980’s: Global models • 1990 - : Advanced data assimilation techniques, use of satellite data, ensemble prediction systems, non-hydrostatic model (back to Bjerknes non-approximated equations)
Popular numerical techniques for NWP • Finite difference approximations, finite element approximations, spectral transform techniques • Semi-implicit techniques: treat fast processes with implicit time integration • Semi-Lagrangian techniques: Move air packages along backward trajectories that end in gridpoints every time step. Use accurate spatial interpolation at the starting point of the trajectories.
NWP in Europe today • Medium range (2-10 days) forecasting is handled by ECMWF (European Centre for Medium Range Forecasting) • Short range forecasting is handled by national weather services • Collaboration in 4 groups for development of short range forecasting
The international HIRLAM project(HIgh Resolution Limited Area Model) • Started in 1985 • Collaboration between the Nordic countries, Spain, Netherlands, Ireland and France • Development of a complete system for Numerical Weather prediction: Model and Data Assimilation • Present emphasis on data assimilation and on high resolution model (1-3 km)
HIRLAM data assimilation . . . . .
Basic data assimilation problems • Degrees of freedom of the model >> number of observation • Balances between state variables in the atmosphere are important • Forecast errors have preferred spatial scales • Fast growth of small-scale perturbations (non-linear instabilities) • Model have multiple time-scales • Observations are made irregularly in space and time • Many different quantities are observed • Random and systematic observation errors
4D-Var cost function Basic non-incremental formulation: X = Model state to be determined by minimization of J Xb = Background model state (e.g. 6 h forecast) B = Covariance matrix (background state error) y = observations H = observation operator (model state to observed variables) R = Covariance matrix (observation errors) In case observations are distributed over a time interval, H includes forward integration of the atmospheric model and backward integration of an adjoint of the atmospheric model.
Successful operational 4D-Var implementations • ECMWF, global model, 12 h assimilation window, incremental formulation, (3D-Var and) 4D-Var has made it possible to introduce a wide range of satellite observations operationally, ECMWF is at present 12 h better than any other NWP center at all forecast ranges • Japan Meteorological Agency, mesoscale incremental 4D-Var, emphasis on improving short range precipitation forecasting, use of radar and GPS observations
HIRLAM 3D-Var and 4D-Var developments 1995-1997: Tangent linear and adjoint of the Eulerian spectral adiabatic HIRLAM. Sensitivity experiments. ”Poor man´s 4D-Var”. 1997-1998: Tangent linear and adjoints of the full HIRLAM physics. 1996-2000: Development of HIRLAM 3D-Var: Background error constraint and observation processing. 2000: First experiments with ”non-incremental” 4D-Var. 2001-2002: Incremental 4D-Var. Simplified physics packages (Buizza vertical diffusion and Meteo France package). 2002: 4D-Var feasibility study. 2003: Semi-Lagrangian scheme (SETLS), outer loops (spectral or gridpoint HIRLAM) and multi-incremenal minimization.
HIRLAM 3D-Var and 4D-Var area extension • FFT’s are used in the spectral tangent-linear and adjoint models and in the back-ground constraint (utilizing the assumption of homogeneity with respect to horizontal correlations). • Double periodic variations through area extension
HIRLAM 3D-VAR – background error formulation Transform the model state increment vector in such a way that the corresponding transform of a model forecast error state vector could be assumed to have a covariance matrix equal to the identity matrix. The following series of transforms are applied in the reference HIRLAM 3D-VAR: • Normalize with forecast error standard deviations • Horizontal spectral transforms (using an extension zone technique) • Reducing dependencies between the mass field and the wind field increments by subtracting geostrophic wind increments from the full wind increments • Project on eigen-vectors of vertical correlation matrices • Normalize with respect to horizontal spectral densities and vertical eigen-values Structure function are non-separable; different horizontal scales at different levels and different vertical scales for different horizontal scales Forecast error statistics are derived by the NMC-method from differences between 48 h and 24 h forecasts valid at the same time. Forecast error standard deviations are re-scaled.
4D-Var feasibility and impact study at DMI • DMI G area, 202 x 190 x 31 gridpoints, 50 km grid • 4D-Var minimization with 150 km increments • Buizza vertical diffusion • One single outer loop • Eulerian dynamics • 6 h assimilation window, 1 h observation windows • Several test runs: 1-10 Dec 1999, 20-30 Dec 1999, Feb 2002 • Conventional observations + AMSU A for Feb 2002 • On the average neutral impact in comparison with 3D-Var (positive for 06 and 18 UTC, negative for 00 and 12 UTC) DMI operational G area
HIRLAM forecasts for the French Christmas storm, valid 28 Dec 1999 00UTC. 3D-Var (left), 4D-Var (right). c. 24h, d. 12h, e. 6h, f. analysis
Parallelization of HIRLAM 4D-Var • FFTs: (1) For the FFT in x-direction the area is divided according to the y-direction; (2) For FFT in the y-direction the area is divided according to x and z (3) Transpose of all data in between the 1D FFT transforms. • Area decomposition also for vertical transforms of spectral coefficients • Semi-Lagrangian interpolation: Area decomposition is the same as for x-direction FFT; A halo zone includes values communicated from neighbor processors. • Observation processing: Each processor has equally many observation of each type to take care of. Field values are obtained by interpolation in the same area decomposition as for x-direction FFT and communicated. • Use of MPI and OPEN MP as standard (SHMEM as an option).
Operational HIRLAM 4D-Varon a linux-cluster at SMHI? Computer timings, 438 x 310 x 40 grid-points, 22 km horizontal resolution: Computer timing estimates, 306 x 306 x 40 gridpoints, 22 km horizontal resolution, SL scheme, 64 PEs SGI 3000: Non-linear forecast, +48h, 5 min timestep 14 min 4D-Var minimization, 3 x 20 iterations, 40 km increment 36 min
Scaling of the spectral HIRLAM on Monolith Large problem size: 438 x 310 x 40 points, 22 km resolution, Semi-Lagrangian time integration, 7,5 min timestep, calculation time in seconds for 3 h forecast
Scaling of HIRLAM 4D-Var on Monolith Small problem size: 202 x 178 x 31 points, 44 km resolution, 88 km resolution in minimization, 20 iterations over a 6 hour window
HIRLAM work on use of remote sensing datain 3D-Var and 4D-Var ATOVS data: AMSU-A radiances over sea (data thinning and bias-correction). Operational at DMI and DNMI. EUMETSAT re-transmission of ATOVS data from the Atlantic and Arctic areas. Work with AMSU-A over land and ice, HIRS and AMSU-B. Radar radial wind vectors and radar VAD profiles: De-aliasing (wind speed ambiguity). Formation of super-observations. Spatial and temporal filters. Development of observation operators taking radar beam bending and spread into account. Impact studies. Ground-based GPS data: Assimilation of Zenith Total Delay. Bias correction (individual for each station) and possibly a horizontally correlated error. Impact studies. Use of slant delays. Other data: Scatterometer winds, MODIS/MERIS IWV and WV above clouds, wind profilers.
Forecast verification scores from implementation of AMSU-A radiances at DMI. Impact study for January 2003. NOA without AMSU A; WIA with AMSU A.
Ten-day assimilation experiment with radar radial wind data: 1-10 December, 1999 Integration area and radar sites Verification of time-series of +24 h wind forecasts (against observations) Observation fit statistics
GPS Assimilation of ground-based GPS data“One person’s noise is another person’s signal” Advantages: • High resolution • 55 to 60 stations (in Sweden) every 15 minutes • All weather, all the time • Very cheap Disadvantages: • Essentially only Integrated Water Vapour • Possible biases and spatially correlated errors
Challenges • Operational 4D-Var on a computer that a small weather service can afford • Use of mesoscale moisture related observations, for example, radar reflectivity • Improved treatment of non-linear processes (utilize ideas from ensemble Kalman filters) • 4D-Var for nowcasting and very short range forecasting with a NWP model at the km-scale