# Application of Adjoint Data Assimilation Method to Atmospheric Aerosol Transport Problems

Authors: Minjie Xu; Kai Fu; Xianqing Lv
Journal: Advances in Mathematical Physics
 (2017)

---

## Abstract

We propose combining the adjoint assimilation method with characteristic finite difference scheme (CFD) to solve the aerosol transport problems, which can predict the distribution of atmospheric aerosols efficiently by using large time steps. Firstly, the characteristic finite difference scheme (CFD) is tested to compute the Gaussian hump using large time step sizes and is compared with the first-order upwind scheme (US1) using small time steps; the US1 method getsE2 error of 0.2887 using Δt=1/450, while CFD method gets a much smaller E2 of 0.2280 using a much larger time step Δt=1/45. Then, the initial distribution of PM2.5 concentration is inverted by the adjoint assimilation method with CFD and US1. The adjoint assimilation method with CFD gets better accuracy than adjoint assimilation method with US1 while adjoint assimilation method with CFD costs much less computational time. Further, a real case of PM2.5 concentration distribution in China during the APEC 2014 is simulated by using adjoint assimilation method with CFD. The simulation results are in good agreement with the observed values. The adjoint assimilation method with CFD can solve large scale aerosol transport problem efficiently.


---

## Body

## 1. Introduction

Atmospheric aerosols are solid and liquid particles suspended in the air with the aerodynamic diameters between 0.001μm and 100 μm. Aerosol particles that enter the atmosphere directly in the form of particles are called primary particles (primary aerosols), and those generated by gas in the atmosphere are called secondary particles (secondary aerosols). Atmospheric aerosol particles include smoke, haze, dust, pollen, and suspended microorganisms, which not only have an impact on environment, but also are a great concern for human health, as small particles can be inhaled into human body and cause disorders of the respiratory system, endocrine system, immune system, and so on [1–3]. Therefore, it is very important to get the prediction of the temporal and spatial distribution of atmospheric aerosols by numerical simulation.Recently, there were many works on the developments of numerical models of atmospheric aerosols. Grell et al. [4] built the WRF/Chem model to simulate the distribution of pollutants in the northeastern United States; Tie et al. [5] used WRF/Chem model in the characterizations of chemical oxidants in Mexico City; Freitas et al. [6] developed the CATT-BRAMS model and studied the amount of PM2.5 in the surface of the southwest Amazon Basin; Han et al. [7] used RAMS-CMAQ to predict the temporal and spatial distribution of nitrate wet deposition in East Asian region; Hudischewskyj and Seigneur [8] used Lagrangian plume models to study aerosols problem; Fu and Liang [9] developed the conservative characteristic finite difference method to predict the distribution of atmospheric aerosols.The spatial and temporal variations of atmospheric aerosols are affected by several physical and chemical processes, like convection, diffusion, deposition, and so on, which make the prediction of the spatial and temporal distribution of atmospheric aerosols difficult. Meanwhile, parameters in model such as the initial condition play an important role in getting good predictions, but they are hard to be decided. Many studies had demonstrated that the adjoint assimilation method is a good way to optimize the uncertain parameters using the observation data [10, 11].The adjoint assimilation method is the combination of optimal control theory and variation principle. In this method, the atmospheric mathematical model is used as constraint condition, and the parameters of model are adjusted by assimilating observational data which can make the simulations of the model be in good agreement with the observations. Panofsky [12] proposed the data assimilation method in 1949 and obtained the objective analysis with two-dimensional global polynomial interpolation method; Constantinescu et al. [13] used chemical data assimilation techniques for improving chemical initial and boundary conditions to quantify the uncertainties of the air quality forecasts; Yumimoto and Uno [14] used data assimilation to inverse modeling of CO emissions; Koo et al. [15] used inverse modeling to predict PM10 in East Asia; Elbern et al. [16] applied 4D-variational data assimilation with an adjoint air quality model to chemical transport model; Wang et al. [17] applied adjoint assimilation method in a Marine Ecosystem Dynamical Model.Due to the stability constraint of tradition methods such as upwind schemes [18, 19], numerical simulations of atmospheric aerosol transport problem need to be carried out by using very small time step sizes, which causes very high computational cost. Douglas Jr. and Russell [20] proposed characteristic method to solve convection-diffusion equations, which has high order accuracy and enables using large time steps; Pironneau and Tabata [21] studied the stability and convergence of the Galerkin-characteristic schemes; Rui and Tabata [22] studied the characteristic finite element scheme for convection-diffusion problems. In this paper, we propose the adjoint assimilation method with characteristic finite difference scheme (CFD) to predict the distribution of atmospheric aerosols. Numerical tests of the moving of a Gaussian hump show that CFD method can get more accurate solutions than the US1 method. Then numerical experiments of the adjoint method with ideal initial condition are carried out, which show the initial distribution of PM2.5 concentration computed by the adjoint assimilation method with CFD is more accurate than that of adjoint assimilation method with US1, and adjoint assimilation method with CFD uses much less time. Finally, we simulate a realistic PM2.5 aerosol transport problem during APEC 2014. The result shows that the adjoint assimilation method with CFD can obtain very good results even using very large time steps.The paper is organized as follows. Section2 presents the adjoint assimilation method. In Section 3, numerical experiments are taken and the results are analyzed. Finally, conclusions are given in Section 4.


## 2. Model and Method

Taking the studied atmospheric aerosol model as a constrain condition, the adjoint assimilation method can optimize the parameters of model by assimilating observation data and make the simulated results close to the observations. The basic idea of the adjoint assimilation method [23] is as follows. Firstly, we define a model by the governing equations and its parameters such as initial conditions. Then by optimizing the parameters, we minimize the cost function which measures the data misfit between the numerical solution of model and the observed data. In this paper, the adjoint assimilation method includes the atmospheric aerosol transport model, the adjoint model, and assimilation processes. The atmospheric aerosol transport model is used to predict the distribution of atmospheric aerosols and the adjoint model is used to get the gradient of the cost function on the initial condition. In assimilation processes, we can optimize the initial condition of model. Then we can use the atmospheric aerosol transport model with the initial condition to predict the distribution of atmospheric aerosols.


### 2.1. The Atmospheric Aerosol Transport Model

Generally speaking, the spatial and temporal variations of atmospheric aerosols are affected by several physical and chemical processes [24, 25], such as convection, diffusion, deposition, chemical reaction, and emissions, which make the prediction of the spatial and temporal distribution of atmospheric aerosols very difficult. In this paper, taking the physical and chemical processes as source term without considering the specific details, we consider the following atmospheric aerosol transport model:(1)∂C∂t+u∂C∂x+v∂C∂y-∂∂xAH∂C∂x-∂∂yAH∂C∂y=S,(2)Ct0,x,y=C0x,y,where t is time and x,y are components of the Cartesian coordinate system; C is the mass concentration of atmospheric aerosol; u,v is the horizontal wind velocity; AH is the horizontal diffusivity coefficient; S is the source term; C0x,y is the given initial condition. Constant boundary conditions are used at the inflow boundary ΓIN, and nongradient boundary conditions are used at the outflow boundary ΓOUT:(3)∂C∂t=0,on⁡ΓIN,∂C∂n=0,on⁡ΓOUT.LetΔx and Δy be the space step of x- and y-directions and Δt the time step. Let Ci,jn denote the atmospheric aerosol concentration at iΔx,jΔy at tn=nΔt.Governing equation (1) can be solved by several numerical schemes, such as the first-order upwind difference scheme, which is usually adopted in the adjoint assimilation method for its simplicity [26–28] and is given as follows:(4)Ci,jn+1-Ci,jnΔt+ui,jnCuupn-CudnnΔxj+vi,jnCvupn-CvdnnΔy=1ΔxjAHCi+1,jn-Ci,jnΔxj-AHCi,jn-Ci-1,jnΔxj+1ΔyAHCi,j+1n-Ci,jnΔy-AHCi,jn-Ci,j-1nΔy+Si,jn,where the upwind scheme is used in the advection term:(5)Cuupn=Ci,jn,Cudnn=Ci-1,jnifui,jn>0,Cuupn=Ci+1,jn,Cudnn=Ci,jnifui,jn<0,Cvupn=Ci,jn,Cvdnn=Ci,j-1nifvi,jn>0,Cvupn=Ci,j+1n,Cvdnn=Ci,jnifvi,jn<0.The stabilization condition of scheme (4) is(6)max⁡AH+ui,jnΔxj2ΔtΔxj2,AH+vi,jnΔy2ΔtΔy2≤12.Due to stabilization condition (6), very small time step has to be used in the computation, which causes long simulation time. Therefore, we propose using the adjoint assimilation method with characteristic finite difference scheme. As the variation of atmospheric aerosol mass concentration is small along the characteristic curve, then, by computing (1) along the characteristic direction, this method can get more accurate solutions even using large time step size.As shown in Figure1, let Ci,jn be the aerosol concentration at x→=iΔx,jΔy at tn=nΔt. We assume that the atmospheric aerosol mass concentration at each grid point at t=tn is known, and we want to know the atmospheric aerosol mass concentration at t=tn+1. Let the characteristic direction be denoted by τ, and let Xτ;x→,tn+1 be the characteristic curve [29]:(7)dXτ;x→,tn+1dτ=UXτ;x→,tn+1,τ,(8)Xtn+1;x→,tn+1=x→.Denote the intersection point of Xτ;x→,tn+1 with the time level tn by x→∗ (point D in Figure 1). We solve x→∗ from (7)-(8) by (9)x→∗=Xtn;x→,tn+1=x→-Ux→,tn+1Δt.Figure 1
The process of constructing characteristic finite difference schemes.The atmospheric aerosol concentrationC¯i,jn at x→∗ can be determined by the interpolation of the concentrations at the points Ck,ln,Ck+1,ln,Ck,l+1n,Ck+1,l+1n surrounding x→∗. Then the characteristic finite difference scheme is given as(10)Ci,jn+1-C¯i,jnΔt=1ΔxjAHCi+1,jn+1-Ci,jn+1Δxj-AHCi,jn+1-Ci-1,jn+1Δxj+1ΔyAHCi,j+1n+1-Ci,jn+1Δy-AHCi,jn+1-Ci,j-1n+1Δy+Si,jn,where(11)C¯i,jn=a00Ck,ln+a01Ck,l+1n+a10Ck+1,ln+a11Ck+1,l+1n,xk<x¯i<xk+1,yl<y¯j<yl+1,a00=xk+1-x¯iΔxyl+1-y¯jΔy,a01=xk+1-x¯iΔxy¯j-ylΔy,a10=x¯i-xkΔxyl+1-y¯jΔy,a11=x¯i-xkΔxy¯j-ylΔy.Then the atmospheric aerosol transport problem can be solved by scheme (10) in N time steps.


### 2.2. The Adjoint Model

The adjoint method, which is derived by using the Lagrange multiplier method and the adjoint operator in functional analysis, plays an important role in the estimation of model parameters, such as initial condition. Nguyen et al. [30] used the adjoint method to estimate the state and parameter in 1D hyperbolic PDEs; Zhao and Lu [31] estimated the parameter in ecosystem model.We define the cost function as(12)JC=12∫ΩC-C′TKCC-C′dΩ,where C is the solution of atmospheric aerosol transport problem; C′ is the observed data of atmospheric aerosol; T denotes matrix transposition; Ω denotes model domain; KC is the weighting matrix of C′ and is defined as(13)KC=1,if the observations are available,0,otherwise.The problem is then transformed into an unconstrained minimization problem. Define the Lagrange function as(14)Lλ,C,C′=JC+λ,FC,S=JC+∫Ωλ∂C∂t+u∂C∂x+v∂C∂y-∂∂xAH∂C∂x-∂∂yAH∂C∂y-SdΩ,where λ is the Lagrangian multiplier, which is a function of t,C,S. To get the minimum of the cost function, based on Lagrange multiplier theory, we make the first-order derivatives of Lagrange function be zero, as follows:(15)∂L∂λ=0,(16)∂L∂C=0,(17)∂L∂S=0.In fact, (15) is (1), and the adjoint equation can be derived from (16).Firstly, we simplify∂L/∂C by subsection integration(18)∂L∂C=-∫Ω∂λ∂t+∂uλ∂x+∂vλ∂y+∂∂xAH∂λ∂x+∂∂yAH∂λ∂ydΩ+∫ΩKCC-C′dΩ.Then the adjoint equation becomes(19)∂λ∂t+∂uλ∂x+∂vλ∂y+∂∂xAH∂λ∂x+∂∂yAH∂λ∂y=KCC-C′.For the solution of adjoint equation (19), the first-order upwind scheme is given similar to (4) as(20)λi,jn-λi,jn+1Δt=ui,jn+1λuupn+1-λudnn+1Δxj+vi,jn+1λvupn+1-λvdnn+1Δy+1ΔxjAHλi+1,jn+1-λi,jn+1Δxj-AHλi,jn+1-λi-1,jn+1Δxj+1ΔyAHλi,j+1n+1-λi,jn+1Δy-AHλi,jn+1-λi,j-1n+1Δy-KCCi,jn+1-C′i,jn+1,where(21)λuupn+1=λi,jn+1,λudnn+1=λi-1,jn+1,ifui,jn+1<0,λuupn+1=λi+1,jn+1,λudnn+1=λi,jn+1,ifui,jn+1>0,λvupn+1=λi,jn+1,λvdnn+1=λi,j-1n+1,ifvi,jn+1<0,λvupn+1=λi,j+1n+1,λvdnn+1=λi,jn+1,ifvi,jn+1>0.We propose the characteristic finite difference schemes for (19)(22)λi,jn-λ¯i,jn+1Δt=1ΔxjAHλi+1,jn-λi,jnΔxj-AHλi,jn-λi-1,jnΔxj+1ΔyAHλi,j+1n-λi,jnΔy-AHλi,jn-λi,j-1nΔy-KCCi,jn+1-C′i,jn+1,where(23)λ¯i,jn+1=a00λk,ln+1+a01λk,l+1n+1+a10λk+1,ln+1+a11λk+1,l+1n+1,xk<x¯i<xk+1,yl<y¯j<yl+1,a00=xk+1-x¯iΔxyl+1-y¯jΔy,a01=xk+1-x¯iΔxy¯j-ylΔy,a10=x¯i-xkΔxyl+1-y¯jΔy,a11=x¯i-xkΔxy¯j-ylΔy.Based on (14), we get the gradient of the cost function on the initial conditions of the aerosol mass concentration Ci,j0:(24)∂J∂C0=∂λ∂t+∂uλ∂x+∂vλ∂y+∂∂xAH∂λ∂x+∂∂yAH∂λ∂y0.Then the optimization of the initial condition can be obtained by the steepest descent method. The relationship betweenC0 and gradient is as follows:(25)C0new=C0old-α∂J∂C0,where α is the step size of steepest descent method.With the initial condition obtained in the adjoint method, we can get more accurate simulation results by the aerosol transport model.


## 2.1. The Atmospheric Aerosol Transport Model

Generally speaking, the spatial and temporal variations of atmospheric aerosols are affected by several physical and chemical processes [24, 25], such as convection, diffusion, deposition, chemical reaction, and emissions, which make the prediction of the spatial and temporal distribution of atmospheric aerosols very difficult. In this paper, taking the physical and chemical processes as source term without considering the specific details, we consider the following atmospheric aerosol transport model:(1)∂C∂t+u∂C∂x+v∂C∂y-∂∂xAH∂C∂x-∂∂yAH∂C∂y=S,(2)Ct0,x,y=C0x,y,where t is time and x,y are components of the Cartesian coordinate system; C is the mass concentration of atmospheric aerosol; u,v is the horizontal wind velocity; AH is the horizontal diffusivity coefficient; S is the source term; C0x,y is the given initial condition. Constant boundary conditions are used at the inflow boundary ΓIN, and nongradient boundary conditions are used at the outflow boundary ΓOUT:(3)∂C∂t=0,on⁡ΓIN,∂C∂n=0,on⁡ΓOUT.LetΔx and Δy be the space step of x- and y-directions and Δt the time step. Let Ci,jn denote the atmospheric aerosol concentration at iΔx,jΔy at tn=nΔt.Governing equation (1) can be solved by several numerical schemes, such as the first-order upwind difference scheme, which is usually adopted in the adjoint assimilation method for its simplicity [26–28] and is given as follows:(4)Ci,jn+1-Ci,jnΔt+ui,jnCuupn-CudnnΔxj+vi,jnCvupn-CvdnnΔy=1ΔxjAHCi+1,jn-Ci,jnΔxj-AHCi,jn-Ci-1,jnΔxj+1ΔyAHCi,j+1n-Ci,jnΔy-AHCi,jn-Ci,j-1nΔy+Si,jn,where the upwind scheme is used in the advection term:(5)Cuupn=Ci,jn,Cudnn=Ci-1,jnifui,jn>0,Cuupn=Ci+1,jn,Cudnn=Ci,jnifui,jn<0,Cvupn=Ci,jn,Cvdnn=Ci,j-1nifvi,jn>0,Cvupn=Ci,j+1n,Cvdnn=Ci,jnifvi,jn<0.The stabilization condition of scheme (4) is(6)max⁡AH+ui,jnΔxj2ΔtΔxj2,AH+vi,jnΔy2ΔtΔy2≤12.Due to stabilization condition (6), very small time step has to be used in the computation, which causes long simulation time. Therefore, we propose using the adjoint assimilation method with characteristic finite difference scheme. As the variation of atmospheric aerosol mass concentration is small along the characteristic curve, then, by computing (1) along the characteristic direction, this method can get more accurate solutions even using large time step size.As shown in Figure1, let Ci,jn be the aerosol concentration at x→=iΔx,jΔy at tn=nΔt. We assume that the atmospheric aerosol mass concentration at each grid point at t=tn is known, and we want to know the atmospheric aerosol mass concentration at t=tn+1. Let the characteristic direction be denoted by τ, and let Xτ;x→,tn+1 be the characteristic curve [29]:(7)dXτ;x→,tn+1dτ=UXτ;x→,tn+1,τ,(8)Xtn+1;x→,tn+1=x→.Denote the intersection point of Xτ;x→,tn+1 with the time level tn by x→∗ (point D in Figure 1). We solve x→∗ from (7)-(8) by (9)x→∗=Xtn;x→,tn+1=x→-Ux→,tn+1Δt.Figure 1
The process of constructing characteristic finite difference schemes.The atmospheric aerosol concentrationC¯i,jn at x→∗ can be determined by the interpolation of the concentrations at the points Ck,ln,Ck+1,ln,Ck,l+1n,Ck+1,l+1n surrounding x→∗. Then the characteristic finite difference scheme is given as(10)Ci,jn+1-C¯i,jnΔt=1ΔxjAHCi+1,jn+1-Ci,jn+1Δxj-AHCi,jn+1-Ci-1,jn+1Δxj+1ΔyAHCi,j+1n+1-Ci,jn+1Δy-AHCi,jn+1-Ci,j-1n+1Δy+Si,jn,where(11)C¯i,jn=a00Ck,ln+a01Ck,l+1n+a10Ck+1,ln+a11Ck+1,l+1n,xk<x¯i<xk+1,yl<y¯j<yl+1,a00=xk+1-x¯iΔxyl+1-y¯jΔy,a01=xk+1-x¯iΔxy¯j-ylΔy,a10=x¯i-xkΔxyl+1-y¯jΔy,a11=x¯i-xkΔxy¯j-ylΔy.Then the atmospheric aerosol transport problem can be solved by scheme (10) in N time steps.


## 2.2. The Adjoint Model

The adjoint method, which is derived by using the Lagrange multiplier method and the adjoint operator in functional analysis, plays an important role in the estimation of model parameters, such as initial condition. Nguyen et al. [30] used the adjoint method to estimate the state and parameter in 1D hyperbolic PDEs; Zhao and Lu [31] estimated the parameter in ecosystem model.We define the cost function as(12)JC=12∫ΩC-C′TKCC-C′dΩ,where C is the solution of atmospheric aerosol transport problem; C′ is the observed data of atmospheric aerosol; T denotes matrix transposition; Ω denotes model domain; KC is the weighting matrix of C′ and is defined as(13)KC=1,if the observations are available,0,otherwise.The problem is then transformed into an unconstrained minimization problem. Define the Lagrange function as(14)Lλ,C,C′=JC+λ,FC,S=JC+∫Ωλ∂C∂t+u∂C∂x+v∂C∂y-∂∂xAH∂C∂x-∂∂yAH∂C∂y-SdΩ,where λ is the Lagrangian multiplier, which is a function of t,C,S. To get the minimum of the cost function, based on Lagrange multiplier theory, we make the first-order derivatives of Lagrange function be zero, as follows:(15)∂L∂λ=0,(16)∂L∂C=0,(17)∂L∂S=0.In fact, (15) is (1), and the adjoint equation can be derived from (16).Firstly, we simplify∂L/∂C by subsection integration(18)∂L∂C=-∫Ω∂λ∂t+∂uλ∂x+∂vλ∂y+∂∂xAH∂λ∂x+∂∂yAH∂λ∂ydΩ+∫ΩKCC-C′dΩ.Then the adjoint equation becomes(19)∂λ∂t+∂uλ∂x+∂vλ∂y+∂∂xAH∂λ∂x+∂∂yAH∂λ∂y=KCC-C′.For the solution of adjoint equation (19), the first-order upwind scheme is given similar to (4) as(20)λi,jn-λi,jn+1Δt=ui,jn+1λuupn+1-λudnn+1Δxj+vi,jn+1λvupn+1-λvdnn+1Δy+1ΔxjAHλi+1,jn+1-λi,jn+1Δxj-AHλi,jn+1-λi-1,jn+1Δxj+1ΔyAHλi,j+1n+1-λi,jn+1Δy-AHλi,jn+1-λi,j-1n+1Δy-KCCi,jn+1-C′i,jn+1,where(21)λuupn+1=λi,jn+1,λudnn+1=λi-1,jn+1,ifui,jn+1<0,λuupn+1=λi+1,jn+1,λudnn+1=λi,jn+1,ifui,jn+1>0,λvupn+1=λi,jn+1,λvdnn+1=λi,j-1n+1,ifvi,jn+1<0,λvupn+1=λi,j+1n+1,λvdnn+1=λi,jn+1,ifvi,jn+1>0.We propose the characteristic finite difference schemes for (19)(22)λi,jn-λ¯i,jn+1Δt=1ΔxjAHλi+1,jn-λi,jnΔxj-AHλi,jn-λi-1,jnΔxj+1ΔyAHλi,j+1n-λi,jnΔy-AHλi,jn-λi,j-1nΔy-KCCi,jn+1-C′i,jn+1,where(23)λ¯i,jn+1=a00λk,ln+1+a01λk,l+1n+1+a10λk+1,ln+1+a11λk+1,l+1n+1,xk<x¯i<xk+1,yl<y¯j<yl+1,a00=xk+1-x¯iΔxyl+1-y¯jΔy,a01=xk+1-x¯iΔxy¯j-ylΔy,a10=x¯i-xkΔxyl+1-y¯jΔy,a11=x¯i-xkΔxy¯j-ylΔy.Based on (14), we get the gradient of the cost function on the initial conditions of the aerosol mass concentration Ci,j0:(24)∂J∂C0=∂λ∂t+∂uλ∂x+∂vλ∂y+∂∂xAH∂λ∂x+∂∂yAH∂λ∂y0.Then the optimization of the initial condition can be obtained by the steepest descent method. The relationship betweenC0 and gradient is as follows:(25)C0new=C0old-α∂J∂C0,where α is the step size of steepest descent method.With the initial condition obtained in the adjoint method, we can get more accurate simulation results by the aerosol transport model.


## 3. Numerical Simulation

### 3.1. Numerical Experiment of Aerosol Transport Model

In this section, we first consider a transport of a Gaussian hump. The characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1) are used to solve atmospheric aerosols model (1), and the results of two schemes are compared.The initial condition is given as(26)cx,y,0=exp⁡-x-x02+y-y022σ02,x,y∈Ω,where the initial center is x0,y0=-0.4,0 and 2σ02=0.02. The horizontal diffusivity coefficient is AH=0.001, the spatial domain is Ω=-1,1×-1,1, and t∈0,T=0,π/4. The velocity is u,v=-4y,4x.The exact solution of the problem with the given initial condition is(27)cx,y,t=exp⁡-x~-x02+y~-y022σ02+4AH,where(28)x~=xcos⁡4t+ysin⁡4t,y~=-xsin⁡4t+ycos⁡4t.LetCn(x,y) denote the approximate solution. The errors L∞ and L2 are defined as follows:(29)E∞n=maxi,j⁡cxi,yj,tn-Cnxi,yj,E2n=∑i,jΔxΔycxi,yj,tn-Cnxi,yj2.We choose different time gridsNt= 45, 50, 55, 60, 75 and small space step of Δx=Δy=h=1/200 to compute the errors and ratios in time of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1). The results are shown in Table 1. We can find the ratio in time of CFD method is first order, but US1 method can not get stable results using these large time step sizes. Therefore, due to the limit of stability (6), we choose different time grids Nt= 450, 500, 550, 600, 750 and the same space steps of h=1/90,1/100,1/110,1/120,1/150 which are proportional to time steps. The experimental results are shown in Figure 2. Comparing to the maximum value 0.8649 of the exact solution, the US1 gets only 0.6680 when Δt=T/750, while CFD method gets a better result of 0.8134 using a much larger Δt of T/75. The errors and ratios are shown in Table 2. We can find the ratios of US1 method and CFD method are both first order, while CFD method converges faster than US1 method and gets better results. Even the CFD method uses large time steps; the errors with different Δt are much smaller than those of the first-order upwind scheme (US1). For example, when Nt=450, E∞ and E2 of US1 are 2.8873×10-1 and 4.3042×10-2, respectively; when Nt=45, E∞ and E2 of CFD are 2.2796×10-1 and 4.1029×10-2, respectively.Table 1
Errors and ratios in time of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1) using small space step ofh=1/200.
N
t
45
50
55
60
75
US1
E
∞
1.5720
e
+
51
3.2108
e
+
56
4.0470
e
+
61
3.5640
e
+
66
2.2081
e
+
80
Ratio
—
-
1.1605
e
+
2
-
1.2322
e
+
2
-
1.3085
e
+
2
-
1.4232
e
+
2
E
2
3.4389
e
+
49
5.9727
e
+
54
8.1163
e
+
59
7.2351
e
+
64
5.2907
e
+
78
Ratio
—
-
1.1451
e
+
2
-
1.2401
e
+
2
-
1.3099
e
+
2
-
1.4306
e
+
2
CFD
E
∞
2.1983
e
-
1
1.9882
e
-
1
1.8218
e
-
1
1.6750
e
-
1
1.3652
e
-
1
Ratio
—
0.9535
0.9169
0.9660
0.9163
E
2
4.0790
e
-
2
3.7005
e
-
2
3.3867
e
-
2
3.1225
e
-
2
2.5355
e
-
2
Ratio
—
0.9243
0.9297
0.9336
0.9331Table 2
Errors and ratios in time of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1) using different time steps and space step ofh=1/90, 1/100, 1/110, 1/120, 1/150.
(a)
N
t
450
500
550
600
750
US1
E
∞
2.8873
e
-
1
2.6926
e
-
1
2.5213
e
-
1
2.3696
e
-
1
2.0070
e
-
1
Ratio
—
0.6624
0.6899
0.7131
0.7443
E
2
4.3042
e
-
2
3.9833
e
-
2
3.7073
e
-
2
3.4672
e
-
2
2.9038
e
-
2
Ratio
—
0.7352
0.7535
0.7694
0.7947
(b)
N
t
45
50
55
60
75
CFD
E
∞
2.2796
e
-
1
2.0689
e
-
1
1.8933
e
-
1
1.7449
e
-
1
1.4119
e
-
1
Ratio
—
0.9204
0.9306
0.9386
0.9488
E
2
4.1029
e
-
2
3.7005
e
-
2
3.3867
e
-
2
3.1225
e
-
2
2.5355
e
-
2
Ratio
—
0.9019
0.9125
0.9214
0.9333Figure 2
The results of (b) the US1 method and (c) the CFD method atT=π/4 using different time steps and proportional space steps comparing with (a) the exact solution.
(a)
The exact solution
(b)
US1 method
(c)
CFD methodThen we choose different space gridsΔx=Δy=h=1/40,1/50,1/60,1/70,1/80 and small time step of Nt=400 to compute the errors and ratios in space of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1). As exhibited in Table 3, both US1 method and CFD method get first-order accuracy in space, while the convergence rate of US1 method is less than CFD method. Besides, the errors with different h are smaller than those of the first-order upwind scheme (US1). For example, when h=1/60, E∞ and E2 of US1 are 3.7864×10-1 and 5.8315×10-2, respectively. While E∞ and E2 of CFD are 3.7034×10-1 and 5.7361×10-2, respectively.Table 3
Errors and ratios in space of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1) using small time step ofNt=400.
h
1
/
40
1
/
50
1
/
60
1
/
70
1
/
80
US1
E
∞
4.7095
e
-
1
4.2104
e
-
1
3.7863
e
-
1
3.4243
e
-
1
3.1126
e
-
1
Ratio
—
0.5020
0.5823
0.6520
0.7147
E
2
7.5471
e
-
2
6.5945
e
-
2
5.8315
e
-
2
5.2054
e
-
2
4.6817
e
-
2
Ratio
—
0.6047
0.6744
0.7368
0.7940
CFD
E
∞
4.6345
e
-
1
4.1301
e
-
1
3.7034
e
-
1
3.3406
e
-
1
3.0294
e
-
1
Ratio
—
0.5353
0.5907
0.6303
0.6546
E
2
7.4605
e
-
2
6.5022
e
-
2
5.7361
e
-
2
5.1086
e
-
2
4.5847
e
-
2
Ratio
—
0.6161
0.6876
0.7516
0.8102The results clearly show that the characteristic finite difference scheme (CFD) can get better solutions of the atmospheric model (1) than the first-order upwind scheme while greatly saving the computational time by using large time step size.


### 3.2. Numerical Experiment of the Adjoint Model

In this subsection, we consider getting initial field ofPM2.5 aerosol mass concentration of atmospheric transport model (1) by using the adjoint method. Initial conditions have important effects on the simulation results of the aerosol transport model. However, in most cases, instead of getting initial distribution of all the simulated mesh grids, we can only get observations for a few locations, which leads to big error of the final results. Therefore, we use the adjoint method to obtain reasonable initial fields that can get good simulation results.In this experiment, we first give an ideal initial distribution ofPM2.5 aerosol mass concentration and solve aerosol transport model (1). Taking the solution as the observation, the experiments are carried out in the following steps.Step 1.
Give a priori distribution of initial condition, solve the aerosol transport model with this initial condition, and get simulated results.Step 2.
Solve cost function (12) using the simulated results and the observation. If the solution decreases to a very small value or the number of iterations exceeds the given iterative steps, then stop; otherwise, continue to run Step 3.Step 3.
Calculate the gradient of the cost function by the adjoint model and adjust initial field with the gradient. Then, new initial value is obtained and run Step1.The simulation domain is 70°E~140°E, 20°N~55°N with0.5∘×0.5∘ spatial resolution that covers China mainland, and the simulation time is 7 days. The simulation time steps of the US1 method and the CFD method are 600 s and 7200 s, respectively. The horizontal diffusion coefficient is AH=100 m2/s.Since thePM2.5 pollution is more serious in the north of China than other areas [32, 33], in EX 1 the mass concentration of PM2.5 is given as follows:(30)Ci,j=-loni-113.62+latj-34.72+130,if-loni-113.62+latj-34.72+130>35,35,if-loni-113.62+latj-34.72+130≤35.First-order upwind scheme (20) and CFD scheme (22) are used to solve adjoint model (19). And we can compute the initial distribution of PM2.5 by the adjoint assimilation method. The results of experiment are shown in Figures 3 and 4. Comparing to the results of adjoint assimilation method with US1, the adjoint assimilation method with CFD gets better agreements with the ideal initial distribution even using large time steps.Figure 3
The 3D images of the results inverted by the adjoint assimilation method with different schemes.
(a)
The ideal initial distribution
(b)
The result of the adjoint assimilation method with CFD
(c)
The result of the adjoint assimilation method with US1Figure 4
The 2D images of the results inverted by the adjoint assimilation method with different schemes.
(a)
The ideal initial distribution
(b)
The result of the adjoint assimilation method with CFD
(c)
The result of the adjoint assimilation method with US1The errors and running time are shown in Table4. J50 is the final value of cost function and J1 is the initial value of cost function. We find the adjoint assimilation method with US1 gets J50/J1 of 4.1174×10-3, mean absolute error (MAE) reducing to 3.6432 μg/m3, E∞ of 2.9095 μg/m3, and E2 of 0.5355 μg/m3 using computational time of 257.11 s, while the adjoint assimilation method with CFD gets better results with J50/J1 of 5.8143×10-3, mean absolute error (MAE) reducing to 5.3705 μg/m3, E∞ of 11.1768 μg/m3, and E2 of 0.6618 μg/m3 using a much less computational time of 64.86 s.Table 4
J
50
/
J
1, errors, and run time inverted by the adjoint assimilation method with different schemes.
Method
J
50
J
1
MAE (μg/m3)
E
∞ (μg/m3)
E
2 (μg/m3)
Time (s)
Before assimilation
After assimilation
EX1
CFD
4.1174
e
-
3
33.8274
3.6432
2.9095
0.5355
64.86
US1
5.8143
e
-
3
33.8274
5.3705
11.1768
0.6618
257.11
EX2
CFD
2.8311
e
-
3
57.7805
3.6802
8.8527
0.4689
53.63
US1
2.9662
e
-
3
57.7805
6.5276
10.8295
0.7204
253.06We then carry out an experiment with another ideal mass concentration ofPM2.5 of folding line distribution:(31)Ci,j=loni-70.0×3-25,90<loni<117,114.5-loni-116.5×5,117≤loni≤132,35.0,others.The results of experiment are shown in Figure 5. Similarly, we can see that even using large time steps, the inverted initial distribution of the adjoint assimilation method with CFD is more accurate than the adjoint assimilation method with US1 comparing to the ideal initial distribution.Figure 5
The 2D images of the results inverted by the adjoint assimilation method with different schemes.
(a)
The ideal initial distribution
(b)
The result of the adjoint assimilation method with CFD
(c)
The result of the adjoint assimilation method with US1Table4 shows the errors and running time. We find the adjoint assimilation method with US1 gets J50/J1 of 2.8311×10-3, mean absolute error (MAE) reducing to 3.6802 μg/m3, E∞ of 8.8527 μg/m3, and E2 of 0.4689 μg/m3 using 600 s time step, while the adjoint assimilation method with CFD gets smaller J50/J1 of 2.9662×10-3, mean absolute error (MAE) reducing to 6.5276 μg/m3, E∞ of 10.8295 μg/m3, and E2 of 0.7204 μg/m3 using a much larger time step of 7200 s.Running the adjoint assimilation method with US1 costs the computational time of 253.06 s, while it costs only 53.63 s using CFD.The experiments show that the adjoint assimilation method with CFD can invert ideal initial distribution of aerosol concentration very well using large time steps.


### 3.3. Practical Experiments

In this part, a real case of thePM2.5 concentration during APEC 2014 in China is carried out by our adjoint assimilation method. We compare the results of adjoint assimilation method with CFD using large time steps with the results of US1 method using small time steps. We take the experiment for the period from November 5 to November 11, 2014. The studied area (70°E~140°E, 20°N~55°N) covers China and is divided into 140×70 grid cells with the horizontal resolution of 0.5∘×0.5∘. The temporal resolutions of adjoint assimilation method with US1 and adjoint assimilation method with CFD are 600 s and 7200 s, respectively. We obtain the distribution of PM2.5 aerosol mass concentration in November from of the historical database of the air quality (https://wat.epmap.org/). Besides, we obtain wind data in November from the National Centers for Environmental Prediction (NCEP). The background wind is determined by the interpolation of these wind data. The horizontal diffusion coefficient is taken as AH=100 m2/s. When iteration times are over 300, the computing stops.By the adjoint method, we get the initial distribution ofPM2.5 aerosol mass concentration in Nov. 5. The spatial distribution of inverted PM2.5 concentration in Nov. 5 is shown in Figure 6. Comparing to the observation of PM2.5 in Nov. 5 in Figure 6(a), the initial distribution of PM2.5 concentration inverted by the adjoint assimilation method with the characteristic finite difference scheme (CFD) using large time steps is more accurate than the results of first upwind scheme (US1) using small time steps. J300 is the final value of cost function and J1 is the initial value of cost function. Table 5 shows J300/J1, errors between inverted initial distribution and observations, and run time of the adjoint assimilation method with two difference schemes. When we use adjoint assimilation method with CFD with Δt=7200 s, we find J300/J1 reduces to 9.8045×10-2, mean absolute error (MAE) reduces to 10.8625 μg/m3, and E2 is 2.4111 μg/m3. While adjoint assimilation method with US1 gets J300/J1 of 1.0901×10-1, mean absolute error (MAE) and E2 are 12.3036 μg/m3 and 2.7571 μg/m3 using small time step size of Δt=600 s.Table 5
J
300
/
J
1, errors, and run time inverted by the adjoint assimilation method with different schemes.
Method
J
300
J
1
MAE (μg/m3)
E
∞ (μg/m3)
E
2 (μg/m3)
Time (s)
Before assimilation
After assimilation CFD
CFD
9.8045
e
-
2
37.8672
10.8625
41.0113
2.4111
2950.09
US1
1.0901
e
-
1
37.8672
12.3036
53.7658
2.7571
45877.78Figure 6
The 2D images of the initial distribution inverted by the adjoint assimilation method with different schemes.
(a)
The distribution of PM2.5 on Nov. 5
(b)
The initial distribution of PM2.5 inverted by the adjoint assimilation method with CFD
(c)
The initial distribution of PM2.5 inverted by the adjoint assimilation method with US1Further, taking the inverted initial distribution ofPM2.5 concentration as the initial condition, we simulate the distribution of PM2.5 concentration from Nov. 5 to Nov. 11 by the aerosol transport model. Time series of PM2.5 concentration simulated by adjoint assimilation method with CFD and US1 in Beijing, Harbin, Shenyang, Xian, Yuxi, and Xiamen during APEC 2014 are shown in Figure 7. It is clear that the simulated PM2.5 concentration by our adjoint assimilation method with CFD using large time steps is in good agreement with the observed concentration. Figure 8 gives the comparison of the time-varying of the average PM2.5 concentration in China from Nov. 5 to Nov. 11 simulated by two methods. We can see that the adjoint assimilation method with CFD can simulate aerosol concentration using large time steps accurately.Figure 7
Comparisons of the time series of simulatedPM2.5 concentrations in 6 cities from Nov. 5 to Nov. 11.
(a)
(b)
(c)Figure 8
The averagedPM2.5 concentrations in China from Nov. 5 to Nov. 11.


## 3.1. Numerical Experiment of Aerosol Transport Model

In this section, we first consider a transport of a Gaussian hump. The characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1) are used to solve atmospheric aerosols model (1), and the results of two schemes are compared.The initial condition is given as(26)cx,y,0=exp⁡-x-x02+y-y022σ02,x,y∈Ω,where the initial center is x0,y0=-0.4,0 and 2σ02=0.02. The horizontal diffusivity coefficient is AH=0.001, the spatial domain is Ω=-1,1×-1,1, and t∈0,T=0,π/4. The velocity is u,v=-4y,4x.The exact solution of the problem with the given initial condition is(27)cx,y,t=exp⁡-x~-x02+y~-y022σ02+4AH,where(28)x~=xcos⁡4t+ysin⁡4t,y~=-xsin⁡4t+ycos⁡4t.LetCn(x,y) denote the approximate solution. The errors L∞ and L2 are defined as follows:(29)E∞n=maxi,j⁡cxi,yj,tn-Cnxi,yj,E2n=∑i,jΔxΔycxi,yj,tn-Cnxi,yj2.We choose different time gridsNt= 45, 50, 55, 60, 75 and small space step of Δx=Δy=h=1/200 to compute the errors and ratios in time of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1). The results are shown in Table 1. We can find the ratio in time of CFD method is first order, but US1 method can not get stable results using these large time step sizes. Therefore, due to the limit of stability (6), we choose different time grids Nt= 450, 500, 550, 600, 750 and the same space steps of h=1/90,1/100,1/110,1/120,1/150 which are proportional to time steps. The experimental results are shown in Figure 2. Comparing to the maximum value 0.8649 of the exact solution, the US1 gets only 0.6680 when Δt=T/750, while CFD method gets a better result of 0.8134 using a much larger Δt of T/75. The errors and ratios are shown in Table 2. We can find the ratios of US1 method and CFD method are both first order, while CFD method converges faster than US1 method and gets better results. Even the CFD method uses large time steps; the errors with different Δt are much smaller than those of the first-order upwind scheme (US1). For example, when Nt=450, E∞ and E2 of US1 are 2.8873×10-1 and 4.3042×10-2, respectively; when Nt=45, E∞ and E2 of CFD are 2.2796×10-1 and 4.1029×10-2, respectively.Table 1
Errors and ratios in time of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1) using small space step ofh=1/200.
N
t
45
50
55
60
75
US1
E
∞
1.5720
e
+
51
3.2108
e
+
56
4.0470
e
+
61
3.5640
e
+
66
2.2081
e
+
80
Ratio
—
-
1.1605
e
+
2
-
1.2322
e
+
2
-
1.3085
e
+
2
-
1.4232
e
+
2
E
2
3.4389
e
+
49
5.9727
e
+
54
8.1163
e
+
59
7.2351
e
+
64
5.2907
e
+
78
Ratio
—
-
1.1451
e
+
2
-
1.2401
e
+
2
-
1.3099
e
+
2
-
1.4306
e
+
2
CFD
E
∞
2.1983
e
-
1
1.9882
e
-
1
1.8218
e
-
1
1.6750
e
-
1
1.3652
e
-
1
Ratio
—
0.9535
0.9169
0.9660
0.9163
E
2
4.0790
e
-
2
3.7005
e
-
2
3.3867
e
-
2
3.1225
e
-
2
2.5355
e
-
2
Ratio
—
0.9243
0.9297
0.9336
0.9331Table 2
Errors and ratios in time of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1) using different time steps and space step ofh=1/90, 1/100, 1/110, 1/120, 1/150.
(a)
N
t
450
500
550
600
750
US1
E
∞
2.8873
e
-
1
2.6926
e
-
1
2.5213
e
-
1
2.3696
e
-
1
2.0070
e
-
1
Ratio
—
0.6624
0.6899
0.7131
0.7443
E
2
4.3042
e
-
2
3.9833
e
-
2
3.7073
e
-
2
3.4672
e
-
2
2.9038
e
-
2
Ratio
—
0.7352
0.7535
0.7694
0.7947
(b)
N
t
45
50
55
60
75
CFD
E
∞
2.2796
e
-
1
2.0689
e
-
1
1.8933
e
-
1
1.7449
e
-
1
1.4119
e
-
1
Ratio
—
0.9204
0.9306
0.9386
0.9488
E
2
4.1029
e
-
2
3.7005
e
-
2
3.3867
e
-
2
3.1225
e
-
2
2.5355
e
-
2
Ratio
—
0.9019
0.9125
0.9214
0.9333Figure 2
The results of (b) the US1 method and (c) the CFD method atT=π/4 using different time steps and proportional space steps comparing with (a) the exact solution.
(a)
The exact solution
(b)
US1 method
(c)
CFD methodThen we choose different space gridsΔx=Δy=h=1/40,1/50,1/60,1/70,1/80 and small time step of Nt=400 to compute the errors and ratios in space of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1). As exhibited in Table 3, both US1 method and CFD method get first-order accuracy in space, while the convergence rate of US1 method is less than CFD method. Besides, the errors with different h are smaller than those of the first-order upwind scheme (US1). For example, when h=1/60, E∞ and E2 of US1 are 3.7864×10-1 and 5.8315×10-2, respectively. While E∞ and E2 of CFD are 3.7034×10-1 and 5.7361×10-2, respectively.Table 3
Errors and ratios in space of the characteristic finite difference scheme (CFD) and the first-order upwind scheme (US1) using small time step ofNt=400.
h
1
/
40
1
/
50
1
/
60
1
/
70
1
/
80
US1
E
∞
4.7095
e
-
1
4.2104
e
-
1
3.7863
e
-
1
3.4243
e
-
1
3.1126
e
-
1
Ratio
—
0.5020
0.5823
0.6520
0.7147
E
2
7.5471
e
-
2
6.5945
e
-
2
5.8315
e
-
2
5.2054
e
-
2
4.6817
e
-
2
Ratio
—
0.6047
0.6744
0.7368
0.7940
CFD
E
∞
4.6345
e
-
1
4.1301
e
-
1
3.7034
e
-
1
3.3406
e
-
1
3.0294
e
-
1
Ratio
—
0.5353
0.5907
0.6303
0.6546
E
2
7.4605
e
-
2
6.5022
e
-
2
5.7361
e
-
2
5.1086
e
-
2
4.5847
e
-
2
Ratio
—
0.6161
0.6876
0.7516
0.8102The results clearly show that the characteristic finite difference scheme (CFD) can get better solutions of the atmospheric model (1) than the first-order upwind scheme while greatly saving the computational time by using large time step size.


## 3.2. Numerical Experiment of the Adjoint Model

In this subsection, we consider getting initial field ofPM2.5 aerosol mass concentration of atmospheric transport model (1) by using the adjoint method. Initial conditions have important effects on the simulation results of the aerosol transport model. However, in most cases, instead of getting initial distribution of all the simulated mesh grids, we can only get observations for a few locations, which leads to big error of the final results. Therefore, we use the adjoint method to obtain reasonable initial fields that can get good simulation results.In this experiment, we first give an ideal initial distribution ofPM2.5 aerosol mass concentration and solve aerosol transport model (1). Taking the solution as the observation, the experiments are carried out in the following steps.Step 1.
Give a priori distribution of initial condition, solve the aerosol transport model with this initial condition, and get simulated results.Step 2.
Solve cost function (12) using the simulated results and the observation. If the solution decreases to a very small value or the number of iterations exceeds the given iterative steps, then stop; otherwise, continue to run Step 3.Step 3.
Calculate the gradient of the cost function by the adjoint model and adjust initial field with the gradient. Then, new initial value is obtained and run Step1.The simulation domain is 70°E~140°E, 20°N~55°N with0.5∘×0.5∘ spatial resolution that covers China mainland, and the simulation time is 7 days. The simulation time steps of the US1 method and the CFD method are 600 s and 7200 s, respectively. The horizontal diffusion coefficient is AH=100 m2/s.Since thePM2.5 pollution is more serious in the north of China than other areas [32, 33], in EX 1 the mass concentration of PM2.5 is given as follows:(30)Ci,j=-loni-113.62+latj-34.72+130,if-loni-113.62+latj-34.72+130>35,35,if-loni-113.62+latj-34.72+130≤35.First-order upwind scheme (20) and CFD scheme (22) are used to solve adjoint model (19). And we can compute the initial distribution of PM2.5 by the adjoint assimilation method. The results of experiment are shown in Figures 3 and 4. Comparing to the results of adjoint assimilation method with US1, the adjoint assimilation method with CFD gets better agreements with the ideal initial distribution even using large time steps.Figure 3
The 3D images of the results inverted by the adjoint assimilation method with different schemes.
(a)
The ideal initial distribution
(b)
The result of the adjoint assimilation method with CFD
(c)
The result of the adjoint assimilation method with US1Figure 4
The 2D images of the results inverted by the adjoint assimilation method with different schemes.
(a)
The ideal initial distribution
(b)
The result of the adjoint assimilation method with CFD
(c)
The result of the adjoint assimilation method with US1The errors and running time are shown in Table4. J50 is the final value of cost function and J1 is the initial value of cost function. We find the adjoint assimilation method with US1 gets J50/J1 of 4.1174×10-3, mean absolute error (MAE) reducing to 3.6432 μg/m3, E∞ of 2.9095 μg/m3, and E2 of 0.5355 μg/m3 using computational time of 257.11 s, while the adjoint assimilation method with CFD gets better results with J50/J1 of 5.8143×10-3, mean absolute error (MAE) reducing to 5.3705 μg/m3, E∞ of 11.1768 μg/m3, and E2 of 0.6618 μg/m3 using a much less computational time of 64.86 s.Table 4
J
50
/
J
1, errors, and run time inverted by the adjoint assimilation method with different schemes.
Method
J
50
J
1
MAE (μg/m3)
E
∞ (μg/m3)
E
2 (μg/m3)
Time (s)
Before assimilation
After assimilation
EX1
CFD
4.1174
e
-
3
33.8274
3.6432
2.9095
0.5355
64.86
US1
5.8143
e
-
3
33.8274
5.3705
11.1768
0.6618
257.11
EX2
CFD
2.8311
e
-
3
57.7805
3.6802
8.8527
0.4689
53.63
US1
2.9662
e
-
3
57.7805
6.5276
10.8295
0.7204
253.06We then carry out an experiment with another ideal mass concentration ofPM2.5 of folding line distribution:(31)Ci,j=loni-70.0×3-25,90<loni<117,114.5-loni-116.5×5,117≤loni≤132,35.0,others.The results of experiment are shown in Figure 5. Similarly, we can see that even using large time steps, the inverted initial distribution of the adjoint assimilation method with CFD is more accurate than the adjoint assimilation method with US1 comparing to the ideal initial distribution.Figure 5
The 2D images of the results inverted by the adjoint assimilation method with different schemes.
(a)
The ideal initial distribution
(b)
The result of the adjoint assimilation method with CFD
(c)
The result of the adjoint assimilation method with US1Table4 shows the errors and running time. We find the adjoint assimilation method with US1 gets J50/J1 of 2.8311×10-3, mean absolute error (MAE) reducing to 3.6802 μg/m3, E∞ of 8.8527 μg/m3, and E2 of 0.4689 μg/m3 using 600 s time step, while the adjoint assimilation method with CFD gets smaller J50/J1 of 2.9662×10-3, mean absolute error (MAE) reducing to 6.5276 μg/m3, E∞ of 10.8295 μg/m3, and E2 of 0.7204 μg/m3 using a much larger time step of 7200 s.Running the adjoint assimilation method with US1 costs the computational time of 253.06 s, while it costs only 53.63 s using CFD.The experiments show that the adjoint assimilation method with CFD can invert ideal initial distribution of aerosol concentration very well using large time steps.


## 3.3. Practical Experiments

In this part, a real case of thePM2.5 concentration during APEC 2014 in China is carried out by our adjoint assimilation method. We compare the results of adjoint assimilation method with CFD using large time steps with the results of US1 method using small time steps. We take the experiment for the period from November 5 to November 11, 2014. The studied area (70°E~140°E, 20°N~55°N) covers China and is divided into 140×70 grid cells with the horizontal resolution of 0.5∘×0.5∘. The temporal resolutions of adjoint assimilation method with US1 and adjoint assimilation method with CFD are 600 s and 7200 s, respectively. We obtain the distribution of PM2.5 aerosol mass concentration in November from of the historical database of the air quality (https://wat.epmap.org/). Besides, we obtain wind data in November from the National Centers for Environmental Prediction (NCEP). The background wind is determined by the interpolation of these wind data. The horizontal diffusion coefficient is taken as AH=100 m2/s. When iteration times are over 300, the computing stops.By the adjoint method, we get the initial distribution ofPM2.5 aerosol mass concentration in Nov. 5. The spatial distribution of inverted PM2.5 concentration in Nov. 5 is shown in Figure 6. Comparing to the observation of PM2.5 in Nov. 5 in Figure 6(a), the initial distribution of PM2.5 concentration inverted by the adjoint assimilation method with the characteristic finite difference scheme (CFD) using large time steps is more accurate than the results of first upwind scheme (US1) using small time steps. J300 is the final value of cost function and J1 is the initial value of cost function. Table 5 shows J300/J1, errors between inverted initial distribution and observations, and run time of the adjoint assimilation method with two difference schemes. When we use adjoint assimilation method with CFD with Δt=7200 s, we find J300/J1 reduces to 9.8045×10-2, mean absolute error (MAE) reduces to 10.8625 μg/m3, and E2 is 2.4111 μg/m3. While adjoint assimilation method with US1 gets J300/J1 of 1.0901×10-1, mean absolute error (MAE) and E2 are 12.3036 μg/m3 and 2.7571 μg/m3 using small time step size of Δt=600 s.Table 5
J
300
/
J
1, errors, and run time inverted by the adjoint assimilation method with different schemes.
Method
J
300
J
1
MAE (μg/m3)
E
∞ (μg/m3)
E
2 (μg/m3)
Time (s)
Before assimilation
After assimilation CFD
CFD
9.8045
e
-
2
37.8672
10.8625
41.0113
2.4111
2950.09
US1
1.0901
e
-
1
37.8672
12.3036
53.7658
2.7571
45877.78Figure 6
The 2D images of the initial distribution inverted by the adjoint assimilation method with different schemes.
(a)
The distribution of PM2.5 on Nov. 5
(b)
The initial distribution of PM2.5 inverted by the adjoint assimilation method with CFD
(c)
The initial distribution of PM2.5 inverted by the adjoint assimilation method with US1Further, taking the inverted initial distribution ofPM2.5 concentration as the initial condition, we simulate the distribution of PM2.5 concentration from Nov. 5 to Nov. 11 by the aerosol transport model. Time series of PM2.5 concentration simulated by adjoint assimilation method with CFD and US1 in Beijing, Harbin, Shenyang, Xian, Yuxi, and Xiamen during APEC 2014 are shown in Figure 7. It is clear that the simulated PM2.5 concentration by our adjoint assimilation method with CFD using large time steps is in good agreement with the observed concentration. Figure 8 gives the comparison of the time-varying of the average PM2.5 concentration in China from Nov. 5 to Nov. 11 simulated by two methods. We can see that the adjoint assimilation method with CFD can simulate aerosol concentration using large time steps accurately.Figure 7
Comparisons of the time series of simulatedPM2.5 concentrations in 6 cities from Nov. 5 to Nov. 11.
(a)
(b)
(c)Figure 8
The averagedPM2.5 concentrations in China from Nov. 5 to Nov. 11.


## 4. Conclusions

In this paper, we use the adjoint assimilation method with the characteristic finite difference scheme (CFD) to solve the aerosol transport problem. The adjoint assimilation method with CFD can solve the problem effectively by large time steps. The effectiveness of CFD method was first shown in computing a Gaussian hump. Numerical results exhibit that the CFD method can get good solutions using 10 times time step of the US1 method and get better results. Further, the ideal initial distribution was inverted by adjoint assimilation method with CFD and US1. Comparing to the results of adjoint assimilation method with US1, the adjoint assimilation method with CFD gets better agreements with the ideal initial distribution even using large time steps. At last, a real case ofPM2.5 concentration distribution in China during the APEC 2014 was simulated and analyzed by using adjoint assimilation method with CFD. The inverted initial distribution of PM2.5 concentration by adjoint assimilation method with CFD was in good agreement with the observation. Besides, the adjoint assimilation method with CFD with large time steps can obtain vary good simulations. The adjoint assimilation method with characteristic finite difference scheme can solve large scale aerosol transport problem efficiently.



---

*Source: 5865403-2017-03-30.xml*