forward model
for k in range(nt):
u[:,k+1] = RK4(rhs,u[:,k],dt,*args)
if (km<nt_m) and (k+1==ind_m[km]):
tmp = w[:,km] - ObsOp(u[:,k+1])
tmp = tmp.reshape(-1,1)
floss = floss + np.linalg.multi_dot(( tmp.T, Ri , tmp ))
km = km + 1
floss = floss[0,0]/2
return floss
```
#### *Example: Lorenz 63 System*

Similar to the 3DVAR demonstration, we apply the described 4DVAR using the first-order adjoint method on the Lorenz 63 system. We also begin with the same erroneous initial condition of **u**(0)=[2, 3, 4] *<sup>T</sup>* and observations are collected each 0.2 time units, contaminated with a Gaussian noise with diagonal covariance matrix defined as **R** = *σ*<sup>2</sup> *<sup>m</sup>***I**3, where *σ<sup>m</sup>* = 0.15 is the standard deviation for measurement noise. Moreover, we simply define a linear observation operator defined as *h*(**u**) = **u**, with a Jacobian of identity matrix. We utilize the simple gradient descent for minimizing the cost functional, equipped by a Golden search method for learning rate optimization. A maximum number of iterations is set to 1000, but we highlight that this is highly dependent on the adopted minimization algorithm as well as the line-search technique. In practice, the evaluation of each iteration might be too computationally expensive, so the number of iterations need to be as low as possible. We define two criteria for convergence, and iterations stop whenever any one of them is achieved. The first one is based on the change in the value of the cost or loss functional and the second one is based on the magnitude of its gradient. Extra criteria might be supplied as well.

✝ ✆

Results for running Listing 9 is shown in Figure 2, where we can notice the significant improvement of predictions, compared to the background trajectories. Moreover, we highlight the correction to the initial conditions in Figure 2 which resulted in the analysis trajectory. This is opposed to the 3DVAR implementation, where correction is applied locally at measurements instants only as seen in Figure 1.

```
Listing 9. Implementation of 4DVAR using the first-order adjoint method for the Lorenz 63 system. ✞ ☎
```

```
import numpy as np
import matplotlib.pyplot as plt
#%% Application: Lorenz 63