A 4 . 5 km resolution Arctic Ocean simulation with the global multi-resolution model FESOM 1 . 4

In the framework of developing a global modeling system which can facilitate modeling studies on Arctic Ocean and highto midlatitude linkage, we evaluate the Arctic Ocean simulated by the multi-resolution Finite Element Sea ice-Ocean Model (FESOM). To explore the value of using high horizontal resolution for Arctic Ocean modeling, we use two global meshes differing in the horizontal resolution only in the Arctic Ocean (24 km vs. 4.5 km). The high resolution significantly improves the model’s representation of the Arctic Ocean. The most pronounced improvement is in the Arctic intermediate layer, in terms of both Atlantic Water (AW) mean state and variability. The deepening and thickening bias of the AW layer, a common issue found in coarse-resolution simulations, is significantly alleviated by using higher resolution. The topographic steering of the AW is stronger and the seasonal and interannual temperature variability along the ocean bottom topography is enhanced in the high-resolution simulation. The high resolution also improves the ocean surface circulation, mainly through a better representation of the narrow straits in the Canadian Arctic Archipelago (CAA). The representation of CAA throughflow not only influences the release of water masses through the other gateways but also the circulation pathways inside the Arctic Ocean. However, the mean state and variability of Arctic freshwater content and the variability of freshwater transport through the Arctic gateways appear not to be very sensitive to the increase in resolution employed here. By highlighting the issues that are independent of model resolution, we address that other efforts including the improvement of parameterizations are still required.

The Arctic air temperature increased more strongly than the global mean temperature in the recent decades (the "Arctic Amplification", e.g., Serreze and Barry, 2011). In the meantime the Arctic Ocean is undergoing unprecedented changes, with a 10 freshening trend in the surface layer McPhee et al., 2009;Rabe et al., 2011;Polyakov et al., 2013b;Haine et al., 2015), warming events (Dmitrenko et al., 2008;Polyakov et al., 2012Polyakov et al., , 2013b, and significant sea ice decline Comiso, 2012;Cavalieri and Parkinson, 2012;Stroeve et al., 2012;Laxon et al., 2013). These changes are accompanied not only by a shift in ocean circulation regimes and physical conditions, but also by substantial changes in biogeochemical processes (Arrigo and van Dijken, 2015;Tremblay et al., 2015). The on-going and future Arctic changes could 15 have a large influence on lower latitude ocean and climate with potential societal impact, although it is a subject that remains under debate (e.g., Vihma, 2014;Wallace et al., 2014). Despite a lot of success in modeling studies of the Arctic Ocean, the state-of-the-art ocean general circulation models still show non-negligible model biases, as illustrated by different model intercomparison studies (Jahn et al., 2012;Johnson et al., 2012;Aksenov et al., 2016;Wang et al., 2016a, b;Ilicak et al., 2016). For example, in the earlier Arctic Ocean Intercomparison 20 Project (AOMIP, Proshutinsky and Kowalik, 2007;Proshutinsky et al., 2011), it was identified that too thick Atlantic Water layers in the Arctic Ocean were simulated in the models, very possibly due to spurious numerical mixing (Holloway et al., 2007;Karcher et al., 2007). In the recent Coordinated Ocean-ice Reference Experiments, phase II (CORE-II, Griffies et al., 2012) project, it was found that this issue still remains one decade later (Ilicak et al., 2016). Large model biases in the upper Arctic Ocean are another common issue in many ocean general circulation models as shown by Wang et al. (2016a, b). They 25 found that the mean state of the liquid freshwater and sea ice simulated in the CORE-II models, including their storage and the Arctic-Subarctic fluxes, shows a very pronounced spread among models, although the temporal variability is more consistently represented. They also showed that all the CORE-II models experienced a dramatic increase in their Arctic liquid freshwater content during the first few decades of model simulations, and the spread of simulated liquid freshwater transport through Fram and Davis Straits amounts to as much as 50% of the mean transport values. The large model uncertainty identified in previous 30 studies calls for further model development efforts in the community.
There are many narrow straits in the world's oceans, which play important roles in connecting different ocean basins but are difficult to resolve with resolution typically used in large scale ocean models. The Arctic Ocean is enclosed by continents and connected to lower latitude oceans via narrow straits. Especially, the three main Canadian Arctic Archipelago (CAA) channels have widths of about 10 km, 30 km and 50 km at their narrowest locations (Melling, 2000). It was shown that the ocean fluxes through these narrow channels can be reasonably resolved when using very high horizontal resolution in model simulations (about 4 km, Wekerle et al., 2013). High resolution is also required to resolve small scale dynamics which could have impact on larger scale circulation and water mass properties. As the first baroclinic Rossby radius is small (about 1 to 15 km) in the Arctic Ocean (Nurser and Bacon, 2014), mesoscale-eddy resolving is difficult to achieve for long simulations even in regional Arctic Ocean models. For process studies, however, simulations with 1-2 km resolution have been used to resolve 5 mesoscale dynamics and ocean circulation in Fram Strait (Kawasaki and Hasumi, 2015;Hattermann et al., 2016;Wekerle et al., 2017b). As computational resources grow with time, the modeling community tends to use higher and higher model resolution. Certainly there is a need in the modeling community to evaluate high resolution models with respect to the common model issues identified in previous model studies as mentioned above.
In the framework of our own model development, we aim to develop a coupled model system that can facilitate to carry 10 out climate research with focus on the Arctic Ocean and Arctic lower-latitude linkage. We use the global Finite Element Sea ice Ocean Model (FESOM, Wang et al., 2014) as its ocean/sea ice component. This model employs unstructured meshes and allows for variable resolution without traditional nesting. With it we are able to allocate finer resolution in the northern high latitudes than in many other parts of the global ocean. In practice, however, optimal ocean resolution in the Arctic Ocean needs to be decided. The finally chosen resolution should help to adequately simulate the key ocean dynamics with confined model 15 biases. At the same time the model system should not be too costly as it will be used in long climate simulations. As one of the first steps towards designing such a system, in this paper we evaluate the simulated Arctic Ocean by FESOM on a global mesh with 4.5 km resolution inside the Arctic Ocean (see Fig. 2). This resolution is higher than typical resolutions used in current climate models (one fourth to one degree) while we still obtain a reasonably high model throughput of about 8 model years per day. 20 For our purpose we carried out ocean simulations driven by prescribed CORE-II atmospheric forcing. A coarse resolution setup of FESOM (with 24 km in the Arctic Ocean) has been used in previous CORE-II studies. We will compare the 4.5 km model results with those from this coarse setup to understand the impact of model resolution. The most climate-relevant metrics of the Arctic Ocean, that is, the Atlantic Water property, the Arctic freshwater budget and sea ice state were used to evaluate the state-of-the-art ocean climate models in the CORE-II Arctic studies (Ilicak et al., 2016;Wang et al., 2016a, b). These studies 25 provided background knowledge on the Arctic Ocean representation in those models and identified their common issues. In this paper we will mainly focus on the key diagnostics used in these studies for evaluating our simulations.
The model setups will be described in Section 2. The results about Atlantic Water and freshwater budget of the Arctic Ocean will be presented in Sections 3 and 4, respectively, followed by discussions (Section 5) and summary (Section 6).
2 Model setup 30 The latest version of FESOM  is used in this study. The ocean dynamical core stems from the early study of Danilov et al. (2004) and Wang et al. (2008). It works with unstructured triangular meshes, so variable grid resolution can be conveniently applied without the necessity of using traditional nesting. It is coupled to a dynamic-thermodynamic sea ice model 4 (a) Mesh resolution (b) Ratio (Rossby radius/resolution) Figure 2. (a) The horizontal grid size of a mesh with 4.5 km in the Arctic Ocean (referred to as mesh HIGH in this paper). (b) Ratio between the first baroclinic Rossby radius and grid size shown with the log2 scale for mesh HIGH. With resolution finer than two grid cells per Rossby radius models may start to resolve mesoscale eddies depending on numerical mixing in the model. The Rossby radius is calculated for each season and the local minimum is used for (b). (Timmermann et al., 2009;Danilov et al., 2015), which is based on the Parkinson and Washington (1979) thermodynamics and uses an updated version of the elastic-viscous-plastic (EVP, Hunke and Dukowicz, 1997) rheology. The sea ice model is discretized on the same surface mesh as the ocean model by using an unstructured-mesh method too.
A blend of two bottom topography data sets is used. North of 69 o N the 2 km resolution version of the International Bathymetric Chart of the Arctic Oceans (IBCAO, Jakobsson et al., 2008) is used, while south of 64 o N the 1 min resolution version 5 of the General Bathymetric Chart of the Oceans (GEBCO) is used. Between 64 o N and 69 o N the topography is taken as a linear combination of the two data sets. An explicit second-order flux-corrected-transport (FCT) scheme (Löhner et al., 1987) is employed in the tracer equations. It helps to preserve monotonicity and eliminate overshoots.
The diapycnal mixing is parameterized with the k-profile scheme proposed by Large et al. (1994). In case of static instability the vertical mixing coefficients are increased as a parametrization for unresolved vertical overturning processes. We apply 10 biharmonic friction with a Smagorinsky (1963) viscosity, which is flow-dependent. The Redi (1982) isoneutral diffusion with small slope approximation and the Gent and McWilliams (1990, GM) parameterization in a skew diffusion form (Griffies, 1998) are used. A reference value is determined for neutral diffusivity and GM thickness diffusivity at each surface grid location by considering the local horizontal resolution. It is then scaled by the squared buoyancy frequency to obtain 3D diffusivity . Two global meshes are compared in this study. The first one (LOW) has 1 degree nominal horizontal resolution in most part of the world's ocean. In the equatorial band the resolution is tripled, and north of 45 o N the resolution is set to about 24 km. On the second mesh (HIGH) the horizontal resolution is further increased to 4.5 km inside the Arctic Ocean (defined by the Arctic gateways of Bering Strait,CAA,Fram Strait and Barents Sea Opening,Fig. 2a). Mesh LOW has been used in the CORE-II model intercomparison studies, and mesh HIGH has been used in a recent study on Arctic sea ice leads (Wang et al., 5 2016c). Judged by comparing the Rossby radius and grid size (Fig. 2b) and inspecting the simulation result, mesh HIGH is not well eddy resolving in the Arctic Ocean, while it permits eddies in the Eurasian and Canadian Basins. To explain the impact of model resolution on the spatial distribution of FW in the Arctic basins, we carried out one sensitivity experiment on an additional mesh (HIGH-CAA). It has 24 km resolution in the Arctic Ocean, except inside the CAA straits where the resolution is increased to 4.5 km. This mesh has been used in the CAA throughflow study by Wekerle et al. (2013). 10 The reference value of the neutral and GM thickness diffusivity is 50 m 2 /s at 24 km resolution, and 4 m 2 /s at 4.5 km. In the vertical 47 z-levels are used with resolution of 10 m in the top 100 m and gradually decreased downwards. The background vertical diffusivity is set to 10 −5 m 2 /s. Three passive tracers are used to illustrate the pathways of different inflow water masses from Fram Strait, Barents Sea Opening (BSO) and Bering Strait. Initially the concentration of these tracers is set to zero. During the simulations they are restored to one in the whole ocean column within these straits. Adding these passive tracers reduces the model throughput by about 20%.
As discussed by Griffies et al. (2009), ocean climate models without a coupled active atmospheric model lack many of the feedbacks present in a fully coupled system, which necessitates restoring of model sea surface salinity (SSS) to observed climatological SSS in global ocean models. In addition, SSS restoring helps to avoid unbounded local salinity trends that can 5 occur in response to inaccuracies in, for example, precipitation. The strength of SSS restoring (defined by a piston velocity) in our simulations is 50 m over 300 days, a value used in many CORE-II models (Danabasoglu et al., 2014). The impact of SSS restoring will be discussed in the discussion section.
The model is forced by the CORE-II interannual atmospheric data set (Large and Yeager, 2009) from 1950 to 2009. The ocean is initialized with temperature and salinity fields from the Polar Science Center Hydrographic Climatology v.3 (PHC3, 10 Steele et al., 2001) and starts from a steady state, and sea ice is initialized with climatological fields obtained from a previous simulation. Interannual monthly mean river runoff is taken from the data provided by Dai et al. (2009), and in the model the river water is spread over a range of 300 km near river mouths to count for unresolved processes .

Atlantic Water in the Arctic Ocean
At the beginning of this and the next sections we will briefly introduce the present understanding of the Arctic Ocean dynamics 15 and changes, and the major issues to be discussed. Then the model results will be presented. The AW below the halocline circulates mainly cyclonically along the continental slope and mid-ocean ridges as topographically steered boundary currents (Rudels et al., 1994;Karcher et al., 2003). The warmer Fram Strait branch and colder BSO 30 branch converge north of the Kara Sea (Schauer et al., 2002;Karcher and Oberhuber, 2002;Maslowski et al., 2004) and continue eastward along the Eurasian slope. After passing the Laptev Sea slope, the boundary current bifurcates into one branch 7 following the Lomonosov Ridge and another following the continental slope (Woodgate et al., 2001). The former brings the AW toward the Fram Strait, while the latter continues into the Canadian Basin. Interannual changes in AW temperature can propagate into the Arctic Ocean via the Fram Strait inflow, leading to temperature variability along the AW boundary current . Pronounced warming events in the Arctic AW layer have been observed in recent decades (Polyakov et al., 2012(Polyakov et al., , 2013b. These recent unprecedented warming implies that the Arctic deep basins are undergoing significant changes. 5 In previous model intercomparison studies with focus on Arctic AW (Holloway et al., 2007;Karcher et al., 2007), it was found that one outstanding issue in most ocean models is the unrealistic deepening and thickening of the AW layer. Numerical mixing associated with the advection operator was suggested to be the major cause (Holloway et al., 2007). The recent CORE-II study indicates that the state-of-the-art ocean general circulation models which are currently used in climate studies still suffer from the deepening of the AW layer (Ilicak et al., 2016). In the following we will explore whether and to what extent 10 this problem can be alleviated by increasing horizontal resolution. Besides the mean state of the AW, we will investigate the model representation of decadal warming and variability on seasonal and interannual time scales.

Spin-up of the AW in Arctic basins
The annual mean temperature horizontally averaged in Eurasian and Canadian basins is plotted as a function of time and depth in Fig. 3. The two basins are defined as the Arctic region where the ocean bottom is deeper than 500 m, and separated by afterwards. The longer simulation presented in the CORE-II model intercomparison work indicates that the depth of the AW layer temperature maxima in the Canadian Basin continues deepening and stays at around 600 m depth after 300 model years (Ilicak et al., 2016). In HIGH no thickening and deepening trend is found in the Canadian Basin. In both simulations the 25 Eurasian Basin is featured with decadal warming events, and the Canadian Basin shows more pronounced warming in recent years. Besides the mean state, the two simulations are also different in their representation of variability and decadal changes, which will be assessed below.

Mean state of AW
To assess the spatial distribution of the AW in the Arctic Ocean, we show the Atlantic Water core temperature (AWCT) derived 30 from 30 years mean model results in Fig. 4a. The AWCT is defined as the maximum temperature over the depth at each location.
The typical spatial pattern of AW is shown by the climatology. The WSC brings warm AW into the Fram Strait, with a fraction recirculating southwards and the remaining part entering the Arctic Ocean. The latter passes the northern slope of Svalbard  The same as ( while it is horizontally more spread in LOW. The observed AWCT is located above 300 m depth in most part of the Eurasian Basin, and deepens towards the Beaufort Sea (Fig. 4b). Simulation HIGH largely reproduces the spatial change of the AWCT depth, and the depth in both the Eurasian and Canadian Basins is well represented. In simulation LOW the AWCT is deeper than the observation in most of the Arctic regions, and the contrast between the Eurasian and Canadian Basins is not as obvious as in the observation and simulation HIGH.   The AW circulation pattern was examined by comparing the topostrophy (Holloway et al., 2007) in the two simulations.
Cyclonic circulation dominates the AW layer boundary currents in both ocean basins similarly in the two simulations (not shown). The fact that simulation LOW also has the correct circulation direction in the Canadian Basin is very possibly just because its resolution is already fine compared to those with problems reported in previous studies (e.g., Holloway et al., 2007).
Indeed, the Arctic Ocean hydrography obtained on mesh LOW is well simulated when compared to the suite of state-of-the-art 5 ocean climate models analyzed in Ilicak et al. (2016).

Variability of AW
Although the Arctic Ocean is at the far end of the North Atlantic Current northern limb, strong warming events have been observed in the Arctic AW layer at the end of the 20th century and beginning of the 21st century Karcher et al., 2003;Polyakov et al., 2012). Despite very limited observations in the remote Arctic deep basins, averaged over decadal 10 time scales the warming events are outstanding and the compiled datasets can be used to assess the model representation of the AW warming (Polyakov et al., 2012). In Fig. 6a the warming in the 1990s relative to the 1970s in the two model simulations and observation is shown. The observation indicates a basinwise warming by about 1 o C in the Eurasian Basin, which propagates into the Canadian Basin crossing the Lomonosov Ridge along the continental slope. Simulation HIGH similarly shows a basinwise warming in the Eurasian Basin, and slightly weaker penetration of the warming signal into the Canadian Basin 5 compared to the observation. In this simulation the boundary current along the continental slope and Lomonosov Ridge shows stronger warming than the basin interior, which is not seen in the observation. This could be partly due to the sparseness of hydrography observations. Simulation LOW also obtains a warming signal in the 1990s, but mainly in the eastern Eurasian Basin and over a large part of the Canadian Basin. As shown by the time series of the basin mean temperature in Fig

15
In the following we will focus on the resolution dependency of the interannual variability of AWCT in the two simulations.
We use the standard deviation (std) of annual mean AWCT as an indicator of the interannual variability. As shown in Fig. 7, the interannual variability is stronger in the Eurasian Basin and weakens along the AW advection pathway in both simulations.
In simulation HIGH, the std is more than 0. simulations. However, the spatial pattern of the std is very different from any of our simulations. In their simulation the highest interannual variability is found starting from the Laptev Sea coast toward the Lomonosov Ridge directly crossing the Eurasian Basin (figure 13 of Lique and Steele, 2012). In this respect the difference in the AW interannual variability induced by different resolutions, although significant, is less pronounced than the difference between two different models.
The magnitude of mean seasonal cycle of the AWCT is also compared in Fig. 7. Both simulations show that the seasonal another penetrating into the Canadian Basin, the seasonal variability also propagates further along these branches. However, the magnitude becomes smaller with distance, which is less than 0.05 o C in HIGH and even lower in LOW. When compared to the strength of interannual variability, in both simulations the seasonal variability is negligible except in the region north of Svalbard and within the narrow boundary current. The most significant seasonal variability is found within the narrow boundary current between the St. Anna Trough and the Laptev Sea continental slope in simulation HIGH, where the ratio between the 5 magnitude of the AWCT seasonal cycle and the std of the annual mean AWCT is about 0.8. Therefore, the BSO branch supplies a large part of the seasonal variability shown in this slope region. The spatial pattern of the AWCT seasonal variability in the Arctic Ocean in our simulations is similar to that derived from a different model at 10 km resolution shown by Lique and Steele (2012). However, the strength of the AWCT seasonal variability in their model simulation is similar to that obtained in our simulation LOW and lower than in HIGH. 10 13 4 Salinity and freshwater budget

Background
A schematic of freshwater (FW) circulation in the pan-Arctic Ocean is shown in Fig. 1. The Arctic Ocean receives a large amount of FW from river runoff, net precipitation, and Pacific Water through Bering Strait (Serreze et al., 2006;Dickson et al., 2007;Haine et al., 2015;Carmack et al., 2016). Liquid FW is stored in the Arctic Ocean with a very non-uniform spatial 5 distribution. We will discuss the FW content in this section, which is defined as the amount of pure FW that could be taken out of the upper ocean so that the ocean salinity is changed to 34.8, the Arctic Ocean reference salinity (see Aagaard and Carmack, 1989). In the calculation of the FW content in our model simulations, the integration is taken from ocean surface to the depth where salinity is equal to the reference salinity. The Canadian Basin is characterized by the largest FW content. Especially in the Beaufort Gyre the FW amounts to about 20 m, whereas it is about 5-10 m in the Eurasian Basin (e.g., Rabe et al., 2011). 10 The anticyclonic Beaufort Gyre is driven by the Beaufort Sea High in atmospheric pressure, which changes the FW content in Beaufort Gyre and the FW distribution between the ocean basins by modulating convergence/divergence of Ekman transport (e.g., Proshutinsky et al., 2002Proshutinsky et al., , 2009Giles et al., 2012). Wind variability over continental shelves can locally induce more significant changes in FW content than the variability from river fluxes, and the variation in large scale atmospheric circulation (Arctic Oscillation) can modify the pathway of river runoff, thus changing the FW distribution between the Arctic basins 15 (Dmitrenko et al., 2008;Morison et al., 2012).
Both liquid FW and sea ice are drained by the Transpolar Drift and released through the Fram Strait. The liquid FW exported through the Fram Strait is slightly larger than sea ice export (Serreze et al., 2006), but the difference has increased during the last decade (Haine et al., 2015). Arctic liquid FW is also released to lower latitudes through the CAA and then the Davis Strait, with an amount similar to that released through the Fram Strait (Serreze et al., 2006;Curry et al., 2014). Sea ice export through 20 the Davis Strait is much less than that from the Fram Strait. The possible climate relevance of the FW cycle in the Arctic Ocean and FW release to the North Atlantic is one of the main reasons for continued research on the Arctic FW budget in both the observation and modeling communities (Carmack et al., 2016).
The liquid FW stored in the Arctic Ocean has been increasing starting from the mid-1990s as shown by observations McPhee et al., 2009;Giles et al., 2012;Polyakov et al., 2013b;Rabe et al., 2014), while sea ice 25 has a persistent declining trend in thickness and volume Laxon et al., 2013). The liquid FW export through the Davis Strait was observed to be lower in the 2000s than in the 1990s (Curry et al., 2014), while the Fram Strait liquid FW export has slightly increased in the 2000s compared to the climatological value (Haine et al., 2015). In recent CORE-II model studies using a suite of global ocean-sea ice models, it was shown that the recent increase in Arctic liquid FW content is caused by both sea ice melting and reduction of total liquid FW export, with the former being more significant in most of the models 30 (Wang et al., 2016b). However, current observations, especially those of liquid FW budget, are still too sparse for the purpose of quantitative verification of the finding based on models.
In the CORE-II model intercomparison project, it was found that the simulated mean state of Arctic FW (FW content and its spatial distribution, and FW transport through Arctic gateways) has significantly large model spreads, even though the same atmospheric forcing was used (Wang et al., 2016a, b). Interannual variability of FW export via the Fram Strait is the least consistently simulated among different Arctic gateways, in both AOMIP and CORE-II models (Jahn et al., 2012;Wang et al., 2016b). Besides, all the CORE-II models show a dramatic increase in the simulated Arctic liquid FW content during the model spin-up phase; afterwards, the FW content stays at the overestimated level (unless overestimated AW salt inflow causes it to drop in one particular model, Wang et al., 2016b). There are indications in some studies that higher model resolution might 5 improve the pathway and spatial distribution of liquid FW (Koldunov et al., 2014;Aksenov et al., 2016). In the following we will compare the Arctic FW budget simulated with FESOM using two different horizontal resolutions. The focus will be on the impact on model spin-up, mean state, and interannual to decadal variability of the FW budget.

Spin-up of salinity and freshwater content
The annual mean salinity horizontally averaged over the Eurasian and Canadian Basins is plotted as a function of time and  (Fig. 10d). In the simulation LOW, more FW takes the release route through the Fram Strait, because the CAA straits are poorly resolved with the coarse resolution and the CAA outflow is restricted (see more details in the section of mean state). This increases the FW content in the western Eurasian Basin. Therefore, it is mainly resolving the narrow straits   in simulation HIGH that leads to the difference of FW content spatial distribution from simulation LOW, rather than the high resolution inside the Arctic Ocean.

Mean state of liquid freshwater
As a consequence of salinity drift during the model spin-up, the basin mean salinity shows biases in the halocline in both Arctic basins (Fig. 5b). The biases are largest at the mid-depth of the halocline, as the salinity is restored to the climatology at the ocean surface, and below the halocline the salinity is determined by that of the AW. As the Eurasian Basin bias in simulation LOW is larger than in simulation HIGH, the overestimation of Arctic FW content is more significant in LOW ( 26% compared 5 to 18%, Table 1). As mentioned above (Section 4.2), the spatial distribution of liquid FW content is better reproduced in HIGH than in LOW (Fig. 10a-d), albeit with overestimation in both simulations, because the high resolution more faithfully represents the narrow channels in CAA. The variety of FW content distributions simulated in different ocean models shown by Wang et al.
(2016b) presumably can be partly attributed to different model representations of the CAA region.
The spatial pattern of FW content is manifested in the simulated sea surface height (SSH, see Fig. 11), since the steric height is dominated by the halosteric component in the Arctic Ocean (e.g., Griffies et al., 2014). In simulation HIGH the SSH field shows a better represented Beaufort Gyre. The CAA resolution not only impacts the FW content pattern, but also the circulation and export pathways of water masses. For example, as illustrated by passive tracers (Fig. 11), in simulation LOW the Pacific Water penetrates more into the Canadian Basin, and has a higher concentration at the Fram Strait than in HIGH. And the better 5 resolved CAA channels in HIGH allow more Atlantic Water from BSO to be released through the CAA.
To access the simulated mean state of FW transport through main Arctic gateways, we  Table 1). Observations suggested that more FW is released through the CAA than through the Fram Strait. This is reproduced in simulation HIGH, while the FW transports through the two export gateways are nearly the same in LOW. Although the simulated CAA FW export in both simulations 10 is lower than the synthesized value, the CAA FW export in HIGH is significantly higher than in LOW, and still within the observational uncertainty range. At the Fram Strait both the ocean volume and FW transports in LOW are higher than in HIGH, as expected from the impact of resolution in the CAA discussed above. Although using higher resolution reduces the Fram Strait FW export, the mean value is still close to the lower bound of the observational range. At the Bering Strait, the FW import is underestimated in the two simulations, with simulation HIGH obtaining a slightly higher value, very close to 15 the lower bound of the observational range. As the Bering Strait ocean volume transports in the two simulations are within the range suggested by observations, the underestimation of FW transports is due to biases in the Pacific Water salinity, which could be still in a phase of large scale spin-up within the model integration period.

Variability of liquid freshwater
The simulated liquid FW contents do not show significant interannual variability, but rather large decadal changes (Fig. 12a). 20 In both simulations the FW content decreases from the beginning of 1980s until the mid-1990s, and then increases afterwards.
The descending trend of observed FW content (Polyakov et al., 2013b) before the mid-1990s is much lower (Fig. 12a). Most of the models used in the CORE-II model intercomparison obtained a significant descending trend before the mid-1990s (  FW distribution in the two Arctic basins described by Morison et al. (2012). They explained that the changes were due to a cyclonic shift in the ocean pathway of Eurasian runoff associated with an increased Arctic Oscillation index.
We are also interested in the model representation of temporal variation of FW content spatial distribution. In Fig. 10e The spatial pattern of positive changes in FW content in HIGH is very similar to that obtained in a model with about 12 10 km resolution in the Arctic Ocean shown by Wang et al. (2016b). In their study most other models show a quite different pattern because of too coarse model resolution used (about 1 degree resolution). Besides the difference in FW content between the two periods, we also calculated the linear trend of vertically integrated FW content from 1996 to 2009 (Fig. 10i-k). The two methods of diagnosing the temporal variation of the FW content provide similar conclusions on the impact of model resolution (compare Fig. 10e-g with Fig. 10i-k). As mentioned above, the resolution inside the CAA plays an important role in 15 representing the mean state of the Arctic Ocean FW content. Here the additional sensitivity experiment, where high resolution is only applied in the CAA channels, helps to illustrate that the high resolution inside the Arctic Ocean does have some impact on the representation of FW content spatial variation, for example, in the Beaufort Gyre and the central Eurasian Basin. The interannual variability of FW transport through the Arctic gateways shows large similarity between the two simulations after the spin-up phase (Fig. 13a-c). The correlation coefficients between the FW transports from the two simulations are similar at the Davis and Fram Straits (0.75 and 0.78, respectively for the period of 1980-2009, Table 1). The correlation is lower than the correlation for ocean volume transports, indicating that the simulated interannual variability of salt transport changes between the two simulations and leads to reduced inter-simulation correlation for FW transports. The current model 5 results are largely similar to the multi-model mean result analyzed by Wang et al. (2016b, also plotted in Fig. 13). The most significant difference is in the Fram Strait FW transport. For example, the changes of FW transport from the mid-1990s to the beginning of 2000s is more pronounced in our two simulations (Fig. 13a). The variability of FW transport at the Fram Strait has been found to be the least consistently simulated among both AOMIP and CORE-II models (Jahn et al., 2012;Wang et al., 2016b). At the Bering Strait the variability is nearly not distinguishable between the two simulations and the multi-model mean obtained in the past model study (Fig. 13c). 5 On decadal time scales, the observed FW export through the Davis Strait in the 2000s is about 10% lower than the climatology of 1980-2000 (Haine et al., 2015). Both simulations reproduce the reduction in the Davis Strait FW export, but the magnitude of reduction is less significant than the observed (Table 1). In simulation HIGH the reduction (about 5%) is larger than in LOW. At the Fram Strait the FW export is suggested to be slightly higher in the 2000s than in the period of 1980-2000 (Haine et al., 2015), while the two simulations similarly show an opposite result, obtaining a reduction of ∼ 300 km 3 /yr in the 10 2000s. The Bering Strait FW transport remains nearly at the same level after the 2000s, which is reproduced by the simulations.
Note that the uncertainty in observations is large due to the sparseness of measurements, and both the observed and simulated changes in FW transports through the Arctic gateways between the two periods are smaller than the magnitude of respective observational uncertainty.

15
The sea ice volume (and corresponding solid FW content) in the two simulations is nearly the same (Table 1), because both the sea ice thickness and concentration are not significantly influenced by the model resolution (Fig. 14). At 4.5 km the sea ice model starts to capture some small scale features (sea ice leads) with reasonable spatial and temporal variability (Wang et al., 2016c). However, the mean sea ice thickness and concentration is not impacted by whether those small scale features are represented or not in the model. Note that much higher model resolution is required in order to simulate sea ice leads with 20 realistic width, because they are typically narrower than 1 km in reality (Tschudi et al., 2002).
The summer sea ice area along the sea ice edge on the Eurasian side is slightly overestimated in both simulations, and the simulated sea ice thickness is about half a meter thicker than the satellite observation in the last few model years (Fig. 14).
Because of lacking sufficient long-term sea ice thickness observations, we compare our simulated solid FW content with the estimate from the PIOMAS Arctic sea ice volume reanalysis (Schweiger et al., 2011). The simulated mean solid FW content 25 in the period of 1980-2000 is about 20% higher than the PIOMAS estimate (Table 1).
The time series of annual mean solid FW content show that the two simulations obtain a descending trend very similar to that from the PIOMAS estimate (Fig. 12b). Compared to the mean value before 2000, the solid FW content decreases by about 4000 km 3 averaged over the 2000s in the simulations, similar to the PIOMAS result (Table 1). Because our simulated sea ice thickness is underestimated before 2000 compared to the submarine observations (as shown in figure 13 of Wang et al., 2016a) 30 and overestimated in later years compared to satellite observations (Fig. 14a-c), the descending trend of solid FW content over the last 3 decades in our simulation and in the PIOMAS estimate as well might be lower than reality.
Arctic sea ice is mainly exported through the Fram Strait. The two simulations produced very similar solid FW transports through the Fram Strait, well representing the observed value (Table 1). Although the sea ice area export through the Fram Strait has been increasing in recent decades due to increasing sea ice drift (e.g., Smedsrud et al., 2017), sea ice volume and thus solid FW export has been decreasing due to the thinning of Arctic sea ice. Compared to the estimate of the 1980-2000 period, the solid FW export flux decreased by 400 km 3 /yr in the 2000s (Haine et al., 2015). The two simulations similarly produce a decrease in the Fram Strait solid FW export of about 300 km 3 /yr between the two periods (Table 1). On interannual time scales the two simulated solid FW transports are well correlated (Table 1 and Fig. 13d). As shown in Wang et al. (2016a, b), ocean 5 climate models can more consistently simulate the interannual variability of solid FW transports through Arctic gateways than the liquid FW transports.

Atlantic Water (a) Heat content and water mass sources
We found that at the end of the simulations the Arctic heat content is higher than the climatology in both simulations, but it is about 4 × 10 21 J higher in simulation LOW than in HIGH. This difference in heat content requires an additional heat flux 5 of 2 TW over 60 years. Due to inaccuracy in diagnosing heat budget terms (e.g., caused by interpolation) and missing heat diffusion terms in our model output, the mismatch between the ocean heat content changing rate and Arctic net heat flux can have the same order of magnitude as this value. Therefore it is hard to carry out analysis of closed heat budget in this and previous modeling studies (e.g., Lique and Steele, 2013). In the following we try to better understand the difference of ocean heat content between the two simulations by analyzing AW passive tracers. 10 The temperature and heat content in the AW layer is influenced by both the warm Fram Strait and the cold BSO AW branches, the latter of which joins the former mainly through the St. Anna Trough (Schauer et al., 2002). Using passive tracers we can obtain the spatial distribution of the two water sources (Fig. 15). The locations of the maxima of the Fram Strait passive tracer coincide with the maxima of temperature in both basins (cf. Fig. 15a and Fig. 5a). The maxima of the Fram Strait passive tracer are located deeper in simulation LOW than in HIGH, consistent to the deepening of the AW layer shown by its temperature 15 maxima. Below about 350 m depth the concentration of the Fram Strait passive tracer in LOW is higher in both Arctic basins than in HIGH (Fig. 15a,c). The BSO passive tracer indicates that cold AW (lower than 0 o C) from the BSO has a lower concentration in the Eurasian Basin in simulation LOW than in HIGH, and the situation is opposite in the Canadian Basin (Fig. 15b,d). AW from both branches have replenished the Canadian Basin more intensively in simulation LOW. BSO AW has the effect to reduce the temperature of the AW layer, so the higher temperature and heat content in the Canadian Basin in LOW should be attributed to 25 the larger amount of warm Fram Strait AW. Because the temperature of the BSO branch is similar after the atmospheric cooling over the continental shelves, the slightly lower volume transport through BSO in LOW (Table 1) has a positive contribution to the overall AW layer heat content. southwards in different paths (Quadfasel et al., 1987;Gascard et al., 1988;Saloranta and Haugan, 2001;Marnela et al., 2013;de Steur et al., 2014). Strong variability associated with mesoscale eddies was observed in the Fram Strait (von Appen et al., 2016), which may play an important role in setting the AW recirculation . The first baroclinic Rossby radius in the Fram Strait is very small (about 2 km in winter), thus our high resolution (4.5 km grid size) simulation cannot resolve mesoscale eddies. At this resolution the simulated warm AW is confined to the strong boundary current and does not 5 reach the central Fram Strait, presenting a cold bias in the center of the strait (Wekerle et al., 2017). As in other high resolution, but not eddy-resolving models (e.g., Fieg et al., 2010), our simulated AW temperature in the boundary current is too high at Fram Strait and north of Svalbard (Fig. 4a). The deficiency indicates a clear requirement for eddy resolving resolution in the Fram Strait region in order to faithfully simulate the amount and property of AW that enters the Arctic basins through the Fram Strait (as shown by Wekerle et al. (2017b)). In long climate simulations, however, it is hardly possible to afford 1 km model 10 resolution in the near future. Accordingly, further effort on parameterizing mesoscale eddy effects is required to represent AW circulation in the Fram Strait.
The AW is located at intermediate depths in the Arctic Ocean and is separated from surface water and sea ice by a strong halocline. However, recent pan-Arctic microstructure measurements of turbulent kinetic energy dissipation reveal that tides can significantly enhance vertical mixing and bring up substantial heat in some areas (Rippeth et al., 2015), implying an impact 5 of AW heat on Arctic sea ice. It was shown that tides can explain a non-negligible part of the sea ice volume reduction in numerical simulations (Luneva et al., 2015). Tides are not simulated in our model, so their potential impact on sea ice and AW characteristics is not explicitly considered. If tides were present in the simulations, and indeed have significant impact on heat uptaken, the influence of AW on sea ice would be different in the two simulations, because the temperature and depth of the AW layer are different between them. Dedicated studies are required to investigate such effect. 10 After the AW warming in the Arctic basins in the 1990s, unprecedented warming has been observed in the 2000s (Polyakov et al., 2013b). However, no warming as strong as observed was obtained in the latter period in the two model simulations (Fig.   3). The AW transport calculated in the northern Fram Strait was found to decrease in the 2000s in the simulations. As the warming in 1990s is reasonably represented in the model, the discrepancy between the observed and simulated temperature variation in recent years could be attributed to model deficiency in representing ocean processes under the condition of sea ice 15 decline, or to the quality of the atmospheric forcing data used. Furthermore, the AW layer temperature in the Arctic interior is not only determined by the amount of warm AW through the Fram Strait and cold AW from the Barents Sea, but also the circulation details of the two branches inside the basins. Research on these subjects is required in future work.

Freshwater
(a) Freshwater content drift and sea surface salinity restoring 20 In both simulations the Arctic liquid FW content increases rapidly during the first 20-30 years, the same as in other ocean climate models participating in the CORE-II intercomparison project analyzed by Wang et al. (2016b). They showed that the source of excessive FW is sea surface salinity (SSS) restoring. We repeated the low resolution simulation with SSS restoring switched off. In this simulation, the salinity in the Canadian Basin has a positive bias instead of a negative one, most pronounced at the surface (Fig. 5c). In the Eurasian Basin the salinity bias is still negative, but becomes smaller. The spin-up in this 25 sensitivity run also takes about 20 to 30 years (Fig. 9). The FW content in Canadian Basin decreases in the spin-up phase, with a magnitude similar to that of FW content increase in the two simulations with SSS restoring (Fig. 9c). In the Eurasian Basin the FW content remains lower than in LOW by nearly a constant offset after 30 model years (Fig. 9b). The total Arctic FW content does not have a significant model drift (Fig. 9a), because the opposite drifts in the two basins largely cancel each other. In the last 30 model years, the variability of FW content in both basins in the sensitivity simulation is similar to that in 30 simulation LOW and HIGH.
In the sensitivity simulation without SSS restoring, the salinity has a positive bias at the surface and negative bias in the lower halocline in the Canadian Basin (Fig. 5c). This implies that too much vertical mixing has taken place, which could be linked to the fact that brine rejection induced convection on very small spatial scales is neither resolved nor properly parameterized in the model: If salt rejected during ice formation is added to the ocean surface, the static instability on the model grid may initialize strong vertical mixing and weaken the vertical salinity gradient, resulting in negative salinity anomaly in the halocline and positive salinity anomaly near the ocean surface. The ocean temperature profile in this depth range is also smoothed out. This issue was discussed by, for example, Duffy et al. (1999) and Nguyen et al. (2009), who proposed to distribute rejected salt in 5 the ocean column with some vertical distribution function, thus preventing static instability. By doing so they got significantly improved salinity profiles. We have implemented this parameterization for brine rejection in the model and are able to achieve improvement on the salinity representation in the Canadian Basin. However, it is not easy to define one particular salt vertical distribution function that can satisfy different Arctic basins and the Southern Ocean at the same time. Some research is required before we can suggest a default scheme for brine rejection in our global model simulations. The background vertical diffusivity 10 was suggested to be one of the key parameters controlling the simulated Arctic Ocean hydrography and circulation, especially in the Canadian Basin Nguyen et al., 2009). In our next model tuning phase, FESOM sensitivity to such model parameters should be carefully examined.
The salinity bias and overestimated FW content in the Eurasian Basin is very possibly caused by inaccurate representation of the pathways of upper ocean circulation (Fig. 8e,f). The Transpolar Drift carrying fresh Pacific Water and river water is located 15 too much to the Eurasian side of the Lomonosov Ridge, and the anticyclonic surface circulation in the Canadian Basin occupies a too large spatial range compared to the observation (Fig. 10a-d). The low resolution inside CAA in simulation LOW causes more FW to release through the Fram Strait, which further increases the FW content in the Eurasian Basin. In the sensitivity simulation without SSS restoring, the SSS is still well represented in the Eurasian Basin. The Eurasian Basin salinity bias in the halocline becomes smaller in this sensitivity simulation, because the FW content is lower and the anticyclonic circulation 20 shrinks in the Canadian Basin, with less FW penetrating into the Eurasian Basin. The upper ocean circulations are mainly driven by surface wind stress, so it is required to investigate the wind forcing fields and the impact of sea ice on the ocean surface stress in order to better understand the Eurasian Basin salinity drift.

(b) Basinwise and Beaufort Gyre freshwater content variability
In this work we have assessed the total Arctic FW content and its distribution between the Eurasian and Canadian Basins. It was 25 found that the increase of FW storage in the Canadian Basin in recent years behaves nearly identically in different simulations ( Fig. 9c). On the contrary, the trend of FW content in the Beaufort Gyre region indicates difference among the simulations ( Fig.   10i-k). Recent research indicates that mesoscale eddy fluxes counteract Ekman pumping, thus playing a crucial role in Beaufort Gyre FW content variability (e.g., Manucharyan et al., 2016;Yang et al., 2016). In model simulations eddy parameterization (applied on coarse meshes) and the effect of implicit numerical mixing will certainly influence the dynamical balance and the 30 Beaufort Gyre FW content. Further effort is required to investigate the model representation of Beaufort Gyre FW content and more importantly its relationship to Arctic FW release to the North Atlantic.

Unstructured-mesh modeling
The variable-resolution functionality provided by unstructured-mesh models offers new possibility in ocean modeling. One can increase model resolution locally where research interest is located, without the necessity of using traditional nesting. On the mesh the resolution can vary in space conveniently according to given functions chosen for particular applications. Many ocean process studies have been carried out taking use of FESOM in global and regional simulations, for example, with focus 5 on overflows , ice shelf cavities (Timmermann et al., 2012), deep water formation (Scholz et al., 2013), polynyas (Haid and Timmermann, 2013), and Arctic sea ice and ocean dynamics (Wekerle et al., 2013;Wang et al., 2016c;Wekerle et al., 2017b). In global ocean climate simulations, the value of unstructured meshes can be more outstanding. One can design meshes with resolution varying continuously in space according to the strength of ocean variability, for example, by considering observed sea surface height variability (Sein et al., 2016) and/or Rossby radius (Sein et al., 2017), to permit or 10 resolve mesoscale eddies in mid to low latitudes. It would be interesting to use this kind of global meshes together with specific mesh refinement in the Arctic Ocean for the purpose of Arctic Ocean studies, as the lower latitude ocean will be better resolved with acceptable increase of computational cost and provide more faithful oceanic linkage with the Arctic Ocean. Developing such a model configuration is aligned with our strategic plan for Arctic Ocean modeling using FESOM and the coupled climate model. It will facilitate us to study and predict not only Arctic changes, but also large scale linkage between high and lower 15 latitudes. Towards this goal, we need to understand, for example, the impact of regional resolution in the Arctic region, using economy configurations as reported in this paper.
With an unstructured-mesh model like FESOM, one can locally increase model resolution to accurately resolve the narrow channels in the CAA and faithfully simulate the FW export (Wekerle et al., 2013). However, if the finest grid size is just used in narrow straits, the model time step and the overall model throughput can be constrained by this grid size (the 20 Courant-Friedrichs-Lewy (CFL) constraint). In ocean climate simulations, therefore, it is not preferable to design resolution in narrow straits to be much higher than the highest resolution used in large ocean basins. Table 2 shows the computational performance of the three simulations studied in this work. The number of grid points in HIGH-CAA is similar to that in LOW, but its time step is one third of LOW. The consequence is that the throughput in HIGH-CAA is three times lower than LOW and the CPU cost is three times higher. Therefore, meshes like HIGH-CAA are mainly used in process studies (one example 25 is this work where we use it to isolate and understand the role of better simulating the CAA throughflow) 1 . Because of good scalability of FESOM (Biastoch et al., 2018), simulation HIGH has a throughput similar to that of HIGH-CAA (Table 2). We usually try to use as many CPUs as possible until the computational performance (in terms of simulated years per day, SYPD) does not further increase effectively 2 . In this case, the model throughput is mainly determined by the time step.  , therefore developing such a functionality has not been put to a high priority.
In the structured-mesh model community, global and near-global ocean models with mesoscale-eddy resolving resolutions 5 have been developed in many groups (e.g., Chassignet et al., 2009;Storch et al., 2012;Oke et al., 2013;Metzger et al., 2014;Dupont et al., 2015;Iovino et al., 2016), and coupled climate models with eddy resolving ocean have also been used in practice (e.g., Griffies et al., 2015). Most of the models analyzed in past CORE-II model intercomparison studies have relatively coarse resolution. For developing our unstructured-mesh model system with regional focus, it would be helpful to communicate experience with the large structured-mesh model community in future high resolution climate model intercomparison projects 10 (for example, through the future CMIP projects where increasing model resolution will be pursued (Haarsma et al., 2016)).

Summary
A faithful model representation of the ocean circulation, water mass property and sea ice state in the Arctic Ocean is still challenging, not only for its mean state, but also for the variability of some of the key diagnostics, in state-of-the-art ocean-sea ice models (e.g., Jahn et al., 2012;Wang et al., 2016a, b;Ilicak et al., 2016). With the development of computing resources and  The simulations of the unstructured-mesh ocean-sea ice model FESOM  with two global meshes differing in resolution in the Arctic Ocean are evaluated. The coarse resolution mesh has been used in previous CORE-II model intercomparison studies (e.g., Griffies et al., 2014;Danabasoglu et al., 2014). Its resolution in the Arctic Ocean is 24 km. On the high resolution mesh the Arctic resolution is increased to 4.5 km. As our intention is to provide information for developing model configurations that can be used for ocean climate simulations, a reasonably high model throughput is a prerequisite. With 4.5 km resolution in the Arctic Ocean we can run FESOM for about 8 model years per day. Using further higher resolution, though preferable for the Arctic region due to very small Rossby radius, would prevent us and groups working with other ocean climate models from carrying out long simulations at the current stage. For ocean process studies, we certainly can use the variable resolution functionality of FESOM to even better resolve local dynamics (for example, using 1 km horizontal grid size 5 locally to resolve mesoscale eddies in the Fram Strait, Wekerle et al. (2017b)). This aspect of Arctic Ocean modeling is beyond the scope of this paper. As we kept the same model resolution outside the Arctic region, we are able to attribute the difference in the two simulations to the Arctic Ocean resolution. Note that we did not try to tune the two model setups separately in this paper. We used a model configuration (schemes and parameters for ocean and sea ice) similar to what has been used in the CORE-II studies (e.g., Danabasoglu et al., 2014), except that eddy diffusivity is scaled by the resolution. 10 At 24 km resolution the simulated AW layer is unrealistically deep and thick, which currently is a common issue in coarse resolution models (Ilicak et al., 2016). Such a model bias was found to be caused by numerical mixing in past AOMIP studies (Holloway et al., 2007). When using 4.5 km resolution, the AW in both Arctic basins is located at the observed depth with a very reasonable thickness (Fig. 4). Note that the tracer advection scheme (a second order FCT scheme) used in our simulations is the one suggested for large scale applications in FESOM, because it enforces monotonicity and has decent computational 15 cost. As we kept the vertical resolution the same in the two simulations, which needs separate investigation, the reduction in numerical mixing is only due to the change in horizontal resolution.
With higher resolution the cyclonic AW boundary current becomes narrower and more energetic. Moreover, the topographic steering on the current is stronger, causing more AW to recirculate along the Lomonosov Ridge in the Eurasian Basin. The resulting constrained penetration of AW into the Canadian Basin in the high resolution simulation helps to eliminate the intensive 20 warming and deepening trend of the AW in the Canadian Basin present in the coarse resolution setup. More AW recirculates in the Eurasian Basin and leaves the Arctic Ocean, which can partly explain that the increase in Arctic heat content is much lower in the high resolution simulation. The strength of interannual and seasonal variability of AW temperature, especially in the boundary current along the continental slope and Lomonosov Ridge, becomes significantly higher with increasing resolution.
The impact of horizontal resolution on ocean surface circulation and FW cycle is limited to the spatial pattern of liquid FW 25 content and pathways of different water masses. It mainly stems from the difference in the representation of the Canadian Arctic Archipelago (CAA) channels, not the resolution in the Arctic basins. The CAA channels are often treated very differently in different ocean climate models, for example, for the number of CAA channels and number of active grid points across the channels, as shown in the model intercomparison study by Wang et al. (2016b). They found that the spread in simulated CAA and Fram Strait FW transports is considerably large. Therefore, inspecting and tuning CAA representation is one of the 30 important tasks in future development of ocean climate models.
The mean state and variability of total and basinwise liquid FW content, the variability of liquid FW transports through Arctic gateways and the characteristics of Arctic sea ice volume and export do not change significantly with increasing resolution.
The recent upward trend of FW content in the Beaufort Gyre shows some sensitivity to the resolution inside the Arctic basin.
How well mesoscale eddies are resolved in the Canadian Basin in the high resolution simulation, and how realistic the effect of 35 eddies is parameterized in the low resolution simulation, both need to be assessed in the context of the interplay with Ekman pumping in future studies. Here it is important to note that the variability of both the Arctic FW storage and release to North Atlantic is insensitive to the model resolution applied in our simulations. We also found that only better resolving the CAA channels (in the simulation where only the CAA is resolved with 4.5 km) did not significantly impact the representation of the AW layer. 5 Besides identifying the impact of horizontal resolution on the Arctic Ocean circulation, we also discussed scientific questions and model issues that need to be explored in future work, and some of the illustrated model issues are common in many other ocean-sea ice models. Overall, increasing model resolution does considerably improve the performance of the Arctic Ocean simulation, while further efforts are necessary to solve remaining issues that are not linked to applied model resolution, and to develop/improve parameterizations that are still required even with best resolution affordable now. Competing interests. The authors declare that they have no conflict of interest.
Acknowledgements. The public availability of different observational data sets and reanalysis data used in this work is a great help for model development, so the efforts of respective working groups are appreciated. We would like to thank I. Polyakov and B. Rabe for providing 20 us their data compiled from large data sets. Q. Wang is funded by the Helmholtz Climate Initiative REKLIM (Regional Climate Change) project. C. Wekerle is funded by the FRontiers in Arctic marine Monitoring program (FRAM). The model simulations were performed at the North-German Supercomputing Alliance (HLRN). We thank the two reviewers and the editor for their helpful comments.