The formulation in musoc
The restoring algorithm follows
D.E. Dietrich et al (2004), where the fluxes are decomposed into two terms:
instantaneous nudging terms and monthly ensembled non-restoring terms. However, for fresh water flux, we use
the same approach with heat flux. Such approach means we add 'salt', not the 'fresh water' into the model,
which is somehow not realistic but rather effective ways to keep model on track with climatology.
For simplicity, we will just show temperature alogorithm. In a given model time step from `t` to `t+dt`,
ommiting advection and diffusion terms, the surface temperature `T` is updated as:
`T_{t+dt}=T_{t}+Q_{m+1/2,n}+T_{t+dt/2}^{n\udg\e}`
where, capital `Q` represents the monthly ensembled flux; `T^{n\udg\e}` is the instantaneous nudging toward
surface climatology; n and m are the current modeling year and month. `Q_{m+1/2,n}` reflects the temperature change
during one time-step due to non-restoring, model-determined ensembled flux (Dierich, 2004); The `T_{t+dt/2}^{n\udg\e}`
term stands for instantaneous physical damping of SST anomalies to the atmosphere (Haney, 1971). For `T^{n\udg\e}`
term, the value is computed at each time-step (and thus we call 'Instantaneous nudging terms'); while the `Q`
term is computed once each month, and thus it is a monthly-fixed value at each time-step.
Instantaneous nudging terms `T^{n\udg\e}`
As previous stated, `T^{n\udg\e}` represents the instantaneous nudging toward surface climatology at each model
time-step. The formulation follows:
`T_{t+dt/2}^{n\udg\e}=(dt/tau_T)*(T^{cli}-T_{t+dt/2})`
Here, `dt` is the model time increment; `tau_T` is the time scale for nudging toward surface climatology, and
`T^{cli}` means the surface climatology interpolated to the present time; `T_{t+dt/2}` is the present temperature
field. Notes that `tau_T` should be chosen to emulate a physically reasonable scales for coupling model simulation.
(ie: In Dietrich 2004, they use `tau_T=30days`.) For testing here, we use a rather long time scale, `tau_T=365days`
to examine the model responses (see figure below).
Monthly ensembled non-restoring terms `Q`
In a similar vein, monthly ensembled terms `Q` is computed from the deviation of current state and surface
climatology. However, it is more complex and involved in the ensmble average to help model retain its variabilities
while tracking with climatology. Conceptually, the `Q` term is computed once each month, and the value `Q` is
derived from averaging previous years' month value(so we call Monthly ensembled flux). The derivation of `Q`
follows two steps. For step one, we need to find the intermediate variable `q` as follows:
`q_{m+1/2,n+1}=Q_{m+1/2,n}+a_{m+1/2,n}`
whereas, `a_{m+1/2,n}` is derived by:
`a_{m+1/2,n}=[T_{m+1}^{cli}-T_{t+dt/2}+sumT^{n\udg\e}]/N`
Here, `N` is the number of time-steps per month; `sumT^{n\udg\e}` is the sum of the instantaneous nudings within
the month. To explain further, the intermediate variable `q` is computed from last years' month value `Q`,
the deviation between current state and climatology and summation of instantaneous nudgings. This is the key
feature of the restoring algorithm, and we will shed more light on it later.
In last step, `Q` is computed from ensemble average `q`:
`Q_{m+1/2,n+1}=(1/{n+1})sum_{j=1}^{n+1}q_{m+1/2,j}`
Now, we have the updated `Q` for current month and complete the entire restoring algorithm for surface fluxes.
Initialization
During the interation of `q`, we need first year's value to initialize the algorithm. Namely, when `n=0`:
`q_{m+1/2,0}=Q_{m+1/2,0}+a_{m+1/2,0}`
For first year, we simply assume `a_{m+1/2,0}=0` in the absence of a reasonable estimation. As for the `Q_{m+1/2,0}`
, a reasonable guess could be:
`Q_{m+1/2,0}=(T_{m+1}^{cli}-T_{m}^{cli})/N`