(Marshall et al., 1997) coupled with a dynamic-thermodynamic sea ice model (see Losch et al., 2010). The adjoint model is generated by the Transformation of Algorithms in Fortran (TAF, Giering & Kaminski, 1998) and is then modified to stabilize the adjoint when calculated over a 4-year assimilation window following the approach in Stammer et al. (2018).

The coupled ocean-sea ice model 

The model domain covers the entire Arctic Ocean north of the Bering Strait and the Atlantic Ocean north of 44°N (see Figure 1). The open boundaries are nested laterally into a 16-km Atlantic-Arctic configuration (Serra et al., 2010). The system has 50 vertical z-levels ranging from 10 m at the surface and 456 m in the deep ocean.  In the horizontal, a curvilinear grid with a resolution ~16 km is used.
Atmospheric states from the 6-hourly NCEP-RA1 reanalysis (Kalnay et al., 1996), including 2-m air temperature, 10-m wind vectors, precipitation rate, 2-m specific humidity, downward longwave radiation, and net shortwave radiation, and bulk formulae are used to compute the surface momentum, heat, and freshwater fluxes. The river runoff is applied near the river mouth with seasonal-varying river discharge. A virtual salt flux parameterization is used to simulate the effects of freshwater input on salinity changes. The K-profiles scheme of Large et al. (1994) is used to parameterize unresolved vertical mixing effects.
 The thermodynamic sea ice model is based on the zero-layer formulation of Semtner (1976) and Hibler (1980). And the dynamic sea ice model is based on Hibler (1979) and is implemented following Zhang and Hibler (1997). The thermodynamic-dynamic sea ice model is then modified for the application of adjoint data assimilation (Fenty & Heimbach, 2013; Losch et al., 2010).