Introduction
The Arctic Ocean is a hotspot in the changing Earth System and is experiencing rapid warming of surface air temperature \cite{amap2019}, changes in the freshwater content of the Canadian Basin \cite{Proshutinsky_2019}, in the mass exchange with the Atlantic Ocean \cite{Dmitrenko_2008} and the Pacific Ocean \cite{Woodgate_2012} and a dramatic decline in sea ice cover \cite{Kwok_2018}. The sea ice coverage at the ocean-atmosphere interface modulates the heat, freshwater, and momentum fluxes and is potentially changing the heat and freshwater budgets of the Arctic system. Improving our understanding of the Arctic Ocean circulation and its changes, its interactions with sea ice and the atmosphere, as well as its exchanges with the Pacific and Atlantic Oceans, is crucial for making predictions and projections of Arctic changes.
Due to the harsh environmental conditions and diplomatic constraints, the Arctic Ocean remains one of the most under-sampled regions of the global Oceans. Although sea ice concentration (SIC) is frequently observed by satellites, the sea ice cover limits hydrographic observations and degrades the accuracy of altimetry-based sea level observations \cite{Armitage_2016,Rose_2019}. The sparseness of observations limits their interpretation in terms of mechanisms and feedbacks in the Arctic Ocean. Numerical model simulations, which provide spatiotemporal-varying ocean states, are therefore used to complement in-situ and satellite observations, and to understand the Arctic Ocean variability and its causes. However, model simulations also suffer from model deficiencies, such as biases in transports \cite{Aksenov_2016,Fieg_2010} and misrepresentations of the ice edge.
To further improve state-of-the-art numerical simulations of the Arctic Ocean circulation and, at the same time, to improve our interpretation of the sparse observations, models are being constrained by ocean observations through data assimilation techniques. Resulting ocean-sea ice reanalyses (e.g., \citealt{Fenty_2013,Fenty_2015,Koldunov_2017}) provide an invaluable basis for assessing variability and trends of the Arctic sea ice cover \cite{Chevallier_2016}, or oceanic transports and their variability \cite{Uotila_2018}. However, analyzing the dynamics of the sea ice changes and the full freshwater budget of the Arctic remains difficult from reanalyses.
A crucial factor in this context is that most of the coupled ocean-sea ice reanalyses are produced with sequential methods (so-called filters), which introduce unphysical mass and energy increments during each analysis step (see \citealt{Stammer_2016} for a review on methods). For instance, the TOPAZ4 reanalysis \cite{Xie_2017} uses an ensemble Kalman Filter to assimilate near-real-time observations, resulting in potential discontinuities in the reanalysis during each analysis step. In contrast, the PIOMAS reanalysis \cite{Lindsay_2006} uses a nudging scheme to assimilate SIC observations and sea surface temperature (SST), thereby continuously adding source terms to the governing equations. Discontinuities introduced by the data assimilation are thereby alleviated. From both reanalyses, we cannot expect to obtain detailed dynamical insights into mechanism leading to observed ocean-sea ice changes.
In contrast, the adjoint assimilation technique preserves model dynamics by bringing the model simulation close to observations through changing control parameters, and the resulting reanalysis is dynamically consistent over a long assimilation window (years to decades) \cite{Stammer_2016}. Applications of the adjoint method excluding the assimilation of sea ice data have been a mature field \cite{Forget_2015,K_hl_2014} but remain challenging when sea ice data is to be incorporated.
Great efforts have been made to implement the adjoint method to coupled ocean-sea ice model (e.g., \citealt{Fenty_2013,Heimbach_2010}) with encouraging results. However, the adjoint model may still suffer from strong unphysical sensitivities, which degrades the usefulness of the adjoint sensitivity, limits the number of iterations performed, and stalls the optimization process. For instance, \citet{Koldunov_2017} reported that only five iterations could be performed in their Arctic ocean-sea ice assimilation system within a 1-year assimilation window because locally large unphysical sensitivities arose, rendering the adjoint sensitivity useless for the optimization algorithms. Causes for the problem still have to be understood in detail, but may likely arise from linearizing nonlinear processes in the sea ice model.
Here, we use the adjoint method to constrain a coupled ocean-sea ice model by ocean observations. Building on the work of \cite{Koldunov_2017} , we aim to improve the previous setup of the Arctic reanalysis by modifying the adjoint model to enable more iterations over longer assimilation windows covering the period 2007-2016. In this work, we will examine the changes imposed in the model through the assimilation procedure and compare the resulting reanalysis against the TOPAZ4 \cite{Xie_2017} and PIOMAS \cite{Lindsay_2006} reanalyses.
The structure of the remaining paper is as follows. The data assimilation system, observations, and the TOPAZ4 and PIOMAS reanalyses are described in Section 2. In Section 3, we assess the details of the model-data misfit reduction. Improvements in SIC and sea ice thickness (SIT) are evaluated and compared against TOPAZ4 and PIOMAS reanalyses in Section 4. Section 5 focuses on ocean state changes concerning SST, freshwater content (FWC), and transports through Fram Strait, Davis Strait, and the Barents Sea Opening. A comparison with the TOPAZ4 reanalysis is also shown in that section. Adjustments of the control variables are examined in Section 6. Section 7 provides concluding remarks.
The data assimilation system and observations
This study builds on the coupled ocean-sea ice data assimilation system presented in the study of
\citet{Koldunov_2017} . The system uses a pan-Arctic configuration of the MITgcm
\cite{Marshall_1997} coupled with a dynamic-thermodynamic sea ice model
\citep[see][]{Losch_2010}. The adjoint model is generated by the Transformation of Algorithms in Fortran
\cite{Giering_1998} (
TAF, Giering & Kaminski, 1998) and is then modified to stabilize the adjoint when calculated over a 4-year assimilation window following the approach in
\citet{Stammer_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 \cite{SERRA_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/NCAR-RA1 reanalysis \cite{dennis1996}, 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 \citet{Large_1994} is used to parameterize unresolved vertical mixing effects.