Thank you very much for the comments. We noticed that you left the comments very early. We are very grateful for this because it gives us more time to consider the improvement of the simulation and manuscript. You mentioned that the focus of our manuscript is "the model configuration and the concept of ISOM as a tool, rather than original scientific results". Yes, this is our main objective, which is why we selected the Development and Technical Paper sector in GMD. However, we will incorporate more original scientific results in the new manuscript as you suggested.

In addition, you mentioned that "the configuration (of ISOM), while differing in certain respects, is not so different to that of Neverworld." To our knowledge, there are two versions of Neverworld. The early version of Neverworld (e.g., in Khani et al. (2019)), which contains simplified topography of the Southern Hemisphere. Neverworld-2 in Marques et al. (2022) is a cross-hemisphere domain that mimics the Atlantic Ocean. We believe ISOM has several distinctive aspects from the two.

(1) Though the topography in ISOM is idealized, it is more complex and realistic compared to the early Neverworld and Neverworld-2, focusing on the dynamics in the Southern Ocean.

(2) The early Neverworld uses a layer model (MOM) with only six constant density layers. While it is not entirely fair to directly compare the vertical grids of a layer model and a z-level model, having only six layers is insufficient to capture all the dynamics, particularly the mixed layer. This inadequacy is likely why Neverworld-2 employs 15 layers. The previous version of ISOM had 40 vertical levels. We now increase it to 75 levels in the new version and improve the simulation of the mixed layer.

(3) The highest horizontal resolution of Neverworld in Khani et al. (2019) is 1/16 degree, which is sufficient to produce eddy-rich results for studying mesoscale eddy effects. However, it may not be adequate to explicitly resolve the deformation radius and fully capture the mesoscale-submesoscale interactions, and to obtain MODNS for the design of MOLES. Neverworld-2 in Marques et al. (2022) has the highest resolution of 1/32 degree, demonstrating the capability to capture the entire mesoscale dynamics. Nevertheless, Neverworld-2 oversimplifies the topography in the Southern Ocean and simulates an excessively strong ACC current. We have optimized the topographic settings in the new version of ISOM, significantly enhancing the simulation of ACC transport (increasing from 70 Sv in the previous manuscript to 145 Sv). The new manuscript is currently in preparation, and we hope you will be glad to read it.

Basically, we bridge the gap between model/simulation and "intermediate complexity" (quote from Marques et al. (2022)) for the Southern Ocean (Marques et al. (2022) seems to focus on the Atlantic Ocean). This is the merit of ISOM in terms of model configuration. As a tool, we aim for it to assist the scientific community in gaining a deeper understanding.

Khani et al. (2019), Diagnosing Subgrid Mesoscale Eddy Fluxes With and Without Topography, JAMES

Marques et al. (2022), NeverWorld2: an idealized model hierarchy to investigate ocean mesoscale eddies across resolutions, GMD

Specific comments: Line 31: "as well as generating ..."

Line 35: A reference to Eden and Greatbatch (2008, doi: 10.1016/j.ocemod.2007.09.002) is also appropriate here. 

Section 2.1: Please define all symbols used in the equations, including u, v, θ, etc. You also switch from θ in equation 4 to T in equation 7, and it's unclear why this is named potential temperature rather than simply temperature given the linear equation of state. 

Line 398: Space missing after "tensor". 

Section 5: This should be written in the past tense, "In this paper, we have introduced ...", etc. 

General comment: The paper is full of acronyms which do not help with its readability. Can the authors please work hard to reduce the number of these? 

Thank you very much for the comments regarding the readability of the manuscript. We will address them in the new version.

Data availability: To be truly reproducible, the authors should provide the version number of MITgcm employed in this study.

Thank you for the comment. To be honest, we have forgotten the exact version of MITgcm that was used. We only recall downloading it sometime in 2021 or early 2022. If you know how to determine the version number, please leave a comment. Perhaps we should provide a copy of MITgcm on the website that distributes the datasets. For the previous manuscript, we have provided all the necessary configuration files, scripts, and pickup files, which should be sufficient to conduct the simulation. We are unsure how different the latest version of MITgcm is from the version we used, but it should be compatible (?).

Line 46: Can you explain what you mean by "its mathematical properties are not fully satisfied by the grid discretization of numerical models"? Does this imply something that is needed beyond the mathematical properties derived in Maddison and Marshall (2013, doi: 10.1017/jfm.2013.259), sections 2.3 and 3.3, and which hold for isopycnal thickness-weighted averaging? 

Thank you for the comment. Maddison and Marshall (2013) utilize and generalize the Eulerian-mean framework, but their framework is still Reynolds-averaged. The clue is that the example they provide is an ensemble average. The ensemble average is a classical Reynolds-averaged (RA) method that does not consider the situation of local scale separation. A good example is an ocean model with a horizontal resolution of (let us assume) 30 km. The model grid separates the entire process controlled by the primitive equations into a resolved part and a sub-grid scale part in a local way, and these two parts essentially interact with each other (cross-scale or multi-scale interactions). This is the scenario for the grid discretization of numerical models, especially when the grid scale is located in the eddying regime (not only for the ocean but also for models of other fluid dynamics). The RA-based as well as the Eulerian-mean framework are definitely self-consistent theoretical systems. If one uses them to diagnose something according to their interests, they can always find explanations for phenomena within the framework. However, the RA-based framework does not align with model practice. An extreme example is zonal average and time average (both are classical RA methods). If one uses the zonal average or time average to represent the resolved part of an ocean model, it is impossible to discern the behavior of the interactions between the resolved scale and sub-grid scale of the model. This does not mean RA-based methods are useless or incorrect; it just means RA-based methods are not suitable to guide sub-grid scale parameterization design and analyze interactive behavior near the model grid-scale. I recommend Khani et al. (2019), Buzzicotti et al. (2023), and Xie et al. (2023). All works discuss relevant aspects in a local or coarse-grained sense (or personally, I prefer to call it in the LES sense), which is more aligned with model practice than RA-based methods.

Khani et al. (2019), Diagnosing Subgrid Mesoscale Eddy Fluxes With and Without Topography, JAMES

Buzzicotti et al. (2023), Spatio-Temporal Coarse-Graining Decomposition of the Global Ocean Geostrophic Kinetic Energy, JAMES

Xie et al. (2023), A Multifaceted Isoneutral Eddy Transport Diagnostic Framework and Its Application in the Southern Ocean, JAMES

Section 2.2: Why is the reference ocean depth so shallow at 3km, and why does the model Drake Passage shallows to 1km (by eye, the model Drake Passage looks even shallower in Figure 1a, but this might be an illusion?) The authors subsequently find that the ACC transport is relatively weak at 65Sv, which I don't find surprising given the overly shallow reference ocean depth and Drake Passage depth (despite the authors' comment to the contrary on line 178). A realistic ACC transport should be a prerequisite and I suggest that the authors considering deepening their ocean to at least 4km, and the model Drake Passage to around 3km since there is circumpolar connection at the latter depth. It is also worth noting that with the surface temperature relaxed strongly to 0 degrees at the south, and the northern vertical temperature profile prescribed through the sponge layer, the thermal wind is effectively imposed in the model:  what transport does this give, assuming no flow at 3km? 

Thank you very much for this comment. In the previous manuscript, the reference ocean depth was set at 3 km. This choice was influenced by the use of several idealized channel models focusing on Southern Ocean dynamics that employed the same depth. We shallowed the model depth of the Drake Passage to 1 km based on the shallowness of the realistic Drake Passage outlet. However, as you and another reviewer pointed out, this depth was too shallow and inhibits the development of ACC transport. After testing, we finally set the reference depth to 4 km and adjust the highest bottom topography in the model Drake Passage region from 1 km to 2.5 km in the updated version of ISOM. This adjustment resulted in a time-averaged ACC transport of around 145 Sv, which falls within a reasonable range. Regarding the sponge layer, it is a mature practice in idealized Southern Ocean modeling works, e.g., Abernathey et al. (2011), Abernathey et al. (2013) and Bischoff&Thompson (2014). For ISOM, the setting of the sponge layer is to follow the realistic temperature profile in the north and, together with the surface relaxation boundary condition, to generate slope isopycnal/isothermal surfaces. The relation of zonal velocity and the stratification has been shown in Figure 4 in the previous manuscript. For more information on the sponge layer technique, I highly recommend Abernathey et al. (2011) and Abernathey et al. (2013). In addition, the ACC transport not only depends on the prescribed stratification profile but is also highly sensitive to the surface wind stress. We will add more content about the relevant issues in the new manuscript.

Abernathey et al. (2011), The Dependence of Southern Ocean Meridional Overturning on Wind Stress, JPO

Abernathey et al. (2013), Diagnostics of isopycnal mixing in a circumpolar channel, Ocean Modelling

Bischoff&Thompson (2014), Configuration of a Southern Ocean Storm Track, JPO

Line 147: Does "The topography within the passage has a piecewise linear depth" mean that the depth varies in a linear manner with partial cells, or that you are employing piecewise linear shaved cells? I think the former, but the text might be read as implying the latter. 

Thank you for the comment. It is the former. We will improve the readability in the new manuscripts.

General comment on Section 2 and Figure 1a: I'm not so sure the idealized topography is as novel as the authors imply. Many model studies have employed idealized Drake Paggase and/or ridge, e.g., the choices made here are not so different to those used in the two layer model of Tansley and Marshall (2001, doi: 10.1175/1520-0485(2001)031<3258:OTDOWD>2.0.CO;2). 

Thank you for the comment. As previously explained, our work bridges the gap between model/simulation and "intermediate complexity" for the Southern Ocean. This is the merit of ISOM in terms of model configuration. As for Tansley and Marshall (2001), we disagree with the comment. Our model is significantly more complicated than theirs, particularly in terms of topography and other model configurations. We have refined and optimized the model topography in the latest version to enhance its realism. The manuscript is currently in preparation. We hope you will appreciate the significance of our work as you review the upcoming manuscript.

Line 242: Why do you choose a 7th order advection scheme for the 8km simulation? This is indeed a very accurate scheme with minimal spurious diapycnal mixing, but I would not have expected that to be a major issue for such short integrations, unless you have a problem with undershoots, as described by Hecht (2010, doi:10.1016/j.ocemod.2010.07.005)?  I understand that the 3rd order DST advection scheme is chosen for computational efficiency at higher resolution - does this introduce any issues and, if not, then why not use it for the 8km integrations for consistency? 

Thank you for the comment. The reentrant channel case in MITgcm documentation uses the 7th-order scheme, and we followed it in the previous work for the 8km simulation. When the resolution becomes finer, we found it is too slow to integrate. Thus, we changed to DST-3 scheme. In the simulation of the new ISOM, we use DTS-3 scheme for all the simulations. Another reason to choose DST-3 scheme is because it is also very accurate, as you can see from Figures 2.14-2.16 in MITgcm documentation. As for Hecht (2010), the numerical effect discussed is mainly caused by the second-order centered-in-space advection scheme and is not that relevant to our work.

Line 297: The Rossby number is usually defined as a non-dimensional parameter quantifying the relative importance of the inertial and Coriolis accelerations, or of relative and planetary vorticity. Here, however, it is used as non-dimensionalized relative vorticity, which gives it a fundamentally different meaning. Indeed, on line 300, the authors state that mesoscale processes dominates the flow when |Ro| ≪ O(1), but there are many regions (filaments) in Figure 6(d) where Ro=1 and yet submesoscale processes are prevalent. I therefore strongly advise avoiding this terminology and instead describing what is shown, i.e., ζ/f. 

Thank you for the comment. However, we do not agree with the comment. The definition and usage of Rossby number in our work are correct and can be supported by numerous submesoscale-related studies from various academic journals. We provide a list of several works here:

Su et al. (2018), Ocean submesoscales as a key component of the global heat budget, Nature Communications

Zhang et al. (2023), Submesoscale inverse energy cascade enhances Southern Ocean eddy heat transport, Nature Communications

Zhang et al. (2023), Ocean Modeling with Adaptive REsolution (OMARE; version 1.0) - refactoring the NEMO model (version 4.0.1) with the parallel computing framework of JASMIN-Part 1: Adaptive grid refinement in an idealized double-gyre case, GMD

Hohenegger et al. (2023), ICON-Sapphire: simulating the components of the Earth system and their interactions at kilometer and subkilometer scales, GMD

Schubert et al. (2019), Submesoscale Impacts on Mesoscale Agulhas Dynamics, JAMES

Naveira Garabato et al. (2022), Kinetic Energy Transfers between Mesoscale and Submesoscale Motions in the Open Ocean's Upper Layers, JPO

Schubert et al. (2020), The Submesoscale Kinetic Energy Cascade: Mesoscale Absorption of Submesoscale Mixed Layer Eddies and Frontal Downscale Fluxes, JPO

Luko et al. (2023), Topographically Generated Submesoscale Shear Instabilities Associated with Brazil Current Meanders, JPO

Balwada et al. (2018), Submesoscale Vertical Velocities Enhance Tracer Subduction in an Idealized Antarctic Circumpolar Current, GRL

Thank you very much for your patience and understanding of the purpose of our work. We have optimized the model configuration, which means the simulation needs to be redone. It takes time. Therefore, the new manuscript is still in preparation. At this stage, we will respond to the comments in the following.

Lines 88-89: The authors claim "We emphasize that the focus of the simulations should be on controlling the dynamics of the idealized model rather than on precise comparisons with observations or realistic model results." but without comparisons to direct observations, it would be difficult to assess the impact of advection schemes of momentum on the states realized by MODNS (e.g. Uchida et al., 2022, https://doi.org/10.5194/gmd-15-5829-2022); Thiry et al., 2024, https://doi.org/10.5194/gmd-17-1749-2024). Do the authors intend on testing different advection schemes as part of the ISOM database? It appears that even at 2km resolution, numerical simulations of the ocean are very sensitive to the advection schemes and/or closures used for surface convection.

Thank you for the comment and recommended papers. We agree that the issue you bring up is worth noting.

Even in the most idealized and simplest case, different advection schemes still have varying numerical effects. The collective effects generated by different sub-grid scale parameterizations and/or the same scheme under different parameter settings are often diverse. Therefore, it is understandable that different ocean models (with distinct advection schemes and parameterization settings) yield inconsistent simulation results. What makes the issue even more complicated is that interactions between the advection scheme and the sub-grid scale parameterization exist, such as the relationship between numerical mixing and physical mixing effects mentioned in Holmes et al. (2021). Though it is not directly related to this work, we are pleased to inform you that we are currently using a realistic global ocean model to evaluate the effects of advection and parameterization schemes.

Now back to MODNS. If we conduct systematic testing on the advection and closure schemes, it would be a huge project. Although I personally think this is meaningful, it goes far beyond the current scope of work. This means we must make some choices. Choosing means giving up something, which can sometimes make us and readers feel a bit discouraged, but that is life. We tend to believe that between the advection scheme and sub-grid scale parameterization, the advection scheme has a more fundamental impact. Therefore, when generating MODNS at this stage, we prioritize choosing advection schemes that have excellent simulation performance in idealized examples such as 1-D and 2-D wave transport (such as DST-3 with a flux limiter, as you can see from Figures 2.14-2.16 in MITgcm documentation). As for the closure scheme, since the horizontal resolution itself is sufficient to ensure the natural development of mesoscale processes, the closure scheme of MODNS actually needs to deal with the collective effects of submesoscale and other fine-scale processes that the model cannot capture. The use of parameterization schemes that focus on different processes would obviously lead to differences in the expression of these processes. But they are not the focus of MODNS (if one wants to better study these processes, the most important thing to do is to simulate them at higher resolution as some kind of DNS corresponding to physical processes). Therefore, we tend to use energy dissipation schemes (Smagorinsky and biharmonic scheme) with a clear effect that is easy to understand in MODNS.

In addition, regarding the comparison with observations, we agree with your opinion. We may not have expressed ourselves clearly before. What we mean is that due to the idealization of the model, it is almost impossible to quantitatively align with the observation. However, for such a model like ISOM, it is necessary to obtain reasonable results that are qualitatively consistent with the observations. Based on this principle, we have reflected on our previous work after considering the opinions of all the reviewers. We think that the problem lies in ensuring that ACC transport can be within a reasonable range. Therefore, we optimized the configuration of ISOM and obtained a time average ACC transport of around 145 Sv in the new simulation, which is within a reasonable range qualitatively consistent with the observation. As your subsequent questions are related to this, we will provide relevant responses under the following questions.  The new manuscript is currently in preparation, and we hope you will be glad to read it.

Holmes et al. (2021), The Geography of Numerical Mixing in a Suite of Global Ocean Models, JAMES

Line 135: The vertical resolution is quite low for a model with the horizontal resolution of 2km (e.g. Balwada et al., 2018, https://doi.org/10.1029/2018GL079244; Ajayi et al., 2020, https://doi.org/10.1029/2019JC015827). Does this hinder the development of baroclinic instability and/or flow-bathymetry interactions in ISOM? Furthermore, the depth of 3000m is quite shallower that the actual Southern Ocean. Given this fact, does it make sense to compare the ACC transport to observational estimates as in Fig. 3?

Thank you very much for this comment. The setting of 3 km depth and 40 vertical levels is used in several idealized channel models focusing on Southern Ocean dynamics [e.g., Abernathey et al. (2011), Abernathey et al. (2013), and Bischoff&Thompson (2014)]. We followed them in the previous manuscript. When you pointed out the issue, we realized that we needed to step forward. We have increased the depth to 4 km and used 75 vertical levels (according to Stewart et al. (2017)) in the new ISOM. This helps to improve the simulation of ACC transport.  

Abernathey et al. (2011), The Dependence of Southern Ocean Meridional Overturning on Wind Stress, JPO

Abernathey et al. (2013), Diagnostics of isopycnal mixing in a circumpolar channel, Ocean Modelling

Bischoff&Thompson (2014), Configuration of a Southern Ocean Storm Track, JPO

Stewart et al. (2017), Vertical resolution of baroclinic modes in global ocean models, Ocean Modelling

Line 184: Regarding the intense meridional displacement of the eddies and idealized ACC jet, have the authors played with different meridional wind profiles? In order to acheive a meridionally narrower ACC, some studies employ a quadratic form in sinusoid (e.g. Balwada et al., 2018, https://doi.org/10.1029/2018GL079244).

Thank you very much for the comment. We have used the new ISOM to test the influence of wind profiles. The ACC transport is sensitive to the wind stress forcing. When fixing the magnitude, the meridional wind profile with a quadratic form in sinusoid results in a narrower ACC stream but smaller ACC transport (135 Sv) compared to the normal one (145 Sv). This is because the total energy injected is less, especially in the model Drake Passage. Given the sensitivity, 2-D wind stress would be effective in fine-tuning the ACC transport and flow pattern. However, in such an idealized work like ISOM and Neverworld, it is not that necessary to involve such complexity. We will give more details in the upcoming manuscript.

Section 3.2: The -3 slope is a prediction for the forward enstrophy cascade range but -5/3 is predicted for the inverse energy cascade range. Given that this is a modeling study, spectral cascades of enstrophy and energy should be also shown to convince the reader on the robustness of the spectral slopes.

Thank you for the comment. We will add more spectral content in the upcoming manuscript.

Section 4: Perhaps the authors have in mind to document the eddy transport tensor diagnosed from their passive tracers in another paper but I would like to see how the tensor changes in time as the homogenization of passive tracers takes place. If the tensor were to show strong time dependence, this would be an obstacle towards obtaining robust estimates of the tensor.

Putting the relevant work into another paper is indeed our plan. At present, we have not yet carried out work on the diagnosis of the transport tensor. However, this does not prevent us from discussing this issue here.

The main problem brought about by the homogenization of passive tracers is that the linear problem of solving the transport tensor becomes underdetermined. Therefore, more tracers are needed to figure out "how the tensor changes in time as the homogenization of passive tracers takes place". The combination of four passive tracers in this manuscript aims at maintaining overdetermination or at least well-posedness of the linear system during tracer release intervals (such as one year).

In addition, the time dependence of the transport tensor is inevitable because the eddying flow field varies in time, thus the transport efficiency must vary in time (e.g., Haigh et al. (2020) and Haigh et al. (2021) show snapshots of the transport tensor).

Regarding the robust estimates of the tensor, I guess what you might want to express is that different tracers may produce inconsistent estimation results for transport tensors (Uchida et al. (2023) and Xie et al. (2023)), namely, there may not be a universal/unique transport tensor for all variables. Personally, I think this is a pain point. I hope both our future work and your work will provide an more explicit answer to this.

Haigh et al. (2020), Tracer-based estimates of eddy-induced diffusivities, Deep Sea Research

Haigh et al. (2021), On eddy transport in the ocean. Part I: The diffusion tensor, Ocean Modelling

Uchida et al. (2023), Cautionary tales from the mesoscale eddy transport tensor, Ocean Modelling

Xie et al. (2023), A multifaceted isoneutral eddy transport diagnostic framework and its application in the Southern Ocean, JAMES

Thank you very much for the comments and your patience. The underestimation of ACC transport is worth further optimization. Following your and other reviewers' suggestions, we are optimizing the model configuration and re-conducting all the simulations. However, it is not finished yet, and the new manuscript is still in preparation. At this stage, we will just respond to the comments and discuss what we think below.

1. Figure 3: as it is explained in the manuscript the ACC transport is significantly underestimated in this setup. One possible reason for this performance is that the channel depth in the Drake Passage area is not realistic. I suggest the channel depth to be increased to be around 4 km and that the bottom topography can be set to be around 2.5 km height like the realistic Drake Passage region. I think this changes to the SO configuration setup would improve the poor estimate of ACC transport. I also suggest the author perform an investigation on changes of ACC transport with the surface wind stress and/or initial salinity profile. 

Thank you for this helpful comment. We have optimized the model topography. The reference domain depth has been set to 4 km as recommended. We have tested the setting of the model Drake Passage. We finally adjusted the highest bottom topography in the model Drake Passage region from 1 km (in the previous manuscript) to 2.5 km. The time-averaged ACC transport is now about 145 Sv, which falls within a reasonable range. Another reviewer has also suggested the depth of Drake Passage, and the specific number provided are deeper. Interestingly, we found that a situation between the two suggested values can give a reasonable result. We are very grateful for the suggestions from the reviewers. In addition, we will also add content about changes in ACC transport with the surface wind stress in the new manuscript.

2. Lines 295-300: In this manuscript, Ro ~ 1 is introduced as an indication for submesoscale process. While moderate Ro might be a sign for smaller eddy-like features (in comparison with mesoscale range), but Ro~1 may not be considered as submesoscale "process". I suggest the authors find other metrics when they discuss submesoscale process (see e.g. McWilliams JC. 2016, Submesoscale currents in the ocean. Proc. R. Soc. A 472: 20160117, for more information).

Thank you for the suggestion and paper recommended. We will add more submesoscale-related indicators in the revised manuscript.

3. Equations (1) and (2): in actual LES/DNS setup for atmosphere and ocean, it is important to have a coupled (anisotropic) horizontal and vertical dissipation scheme in horizontal and vertical momentum equations (see e.g. Khani S, ML Waite 2020, An anisotropic subgrid-Scale parameterization for large-eddy simulations of stratified turbulence. Mon. Wea. Rev. 148, 4299-4311). I understand that if the authors might leave implementing this coupled/anisotropic setup for their future work, but I would like they add a short statement about this in their manuscript since their model setup is performed at very resolution 2 km. 

Thank you very much for your understanding and the paper. We agree that a coupled (anisotropic) horizontal and vertical distribution scheme will play a crucial role in high-resolution models. We will add the statement in the manuscript. Personally, I like Khani&Waite (2020). It would be very interesting to conduct an a prior test on the Khani&Waite (2020) scheme using MODNS data generated by ISOM. In addition, I also recommended the work of Khani&Waite (2020) to my collaborators who are planning high-resolution experiments using a realistic multi-resolution ocean model.

4. Figure 5: the horizontal axis in these plots is labelled as 'scale' which might be confusing. I suggest using horizontal wavenumber k_h [1/km] (the authors might consider labeling the top horizontal axis with the 'horizontal scale [km]'.

Minor comments:

Line 365: "submesoscale signals" ---> submesoscale filaments

Lines 397 and 400: put a space before '('.

I suggest a proof-reading before submitting the revision (there are a few grammatical typos in the text). 

Thank you for your patience. We will improve the manuscript and enhance the readability.