Representation of Bay of Bengal Upper-Ocean Salinity in General Circulation Models BAY OF BENGAL : FROM MONSOONS TO MIXING

Author Posting. © The Oceanography Society, 2016. This article is posted here by permission of The Oceanography Society for personal use, not for redistribution. The definitive version was published in Oceanography 29, no. 2 (2016): 38–49, doi:10.5670/oceanog.2016.37.

general, sea surface temperature (SST) remains over 28°C in the BoB throughout the year, which is favorable for generation of active convection (Gadgil et al., 1984;Graham and Barnett, 1987;Sanilkumar et al., 1994).The shallow halocline and enhanced stratification in the upper ocean are known to maintain the high heat content and SST in the BoB (Shenoi et al., 2002, Mahadevan et al., 2016, in this issue).Further, such salinity structure influences the evolution of mixed layer temperature in regions of near-surface haline stratification (Rao andSivakumar, 1999, 2003;Howden and Murtugudde, 2001;Kara et al., 2003;de Boyer Montégut et al., 2004).
Representation of the BoB salinity structure both in coupled and uncoupled (stand-alone) general circulation models (GCMs) has been challenging due to the complexity arising from the unique physical and dynamical processes in this basin.The East India Coastal Current (EICC) is known to play an important role in maintaining the large-scale mass and salt balances between the BoB and the Arabian Sea (e.g., Durand et al., 2009).Studies suggest that the EICC transports freshwater southward along the east coast of India at the end of summer, and it reaches the southern tip of India by the end of boreal fall (e.g., Durand et al., 2009;Akhil et al., 2014).Apart from this, the winter monsoon current transports low-salinity BoB water to the Arabian Sea (Vinayachandran et al., 2005).Previous ocean model studies (e.g., Akhil et al., 2014) have suggested that erosion of the freshwater tongue along the east coast of India is due to vertical processes rather than horizontal advection.Furthermore, observations have revealed the importance of upward pumping of saltier water to the surface to maintain the salt balance in the southern BoB (Vinayachandran et al., 2013).This suggests the importance of both vertical and lateral processes in controlling the salinity budget in the BoB.
Thus, it is essential to properly represent freshwater flux, current systems, and vertical processes (e.g., Benshila et al., 2014;Akhil et al., 2014;Wilson and Riser, 2016;Parekh et al., 2016;Behara and Vinayachandran, 2016, and references therein) within GCMs so that they simulate realistic surface and subsurface salinity structure in the BoB.Previous modeling studies examined BoB surface salinity variability using sensitivity

INTRODUCTION
The Bay of Bengal (BoB) is a unique tropical basin whose physical properties exhibit extreme variability.This variability is governed by southwesterly monsoon winds from May to September and northeasterly monsoon winds from October to March (e.g., Schott and McCreary, 2001;Shankar et al., 2002).In the BoB, precipitation exceeds evaporation (e.g., Prasad, 1997), and adjoining rivers input a large amount of freshwater (e.g., Subramanian, 1993), making the upper ocean in this region substantially fresher compared to other tropical oceans (e.g., Sengupta et al., 2006).The low-salinity upper layers maintain stable stratification (Narvekar and Prasanna Kumar, 2006;Agarwal et al., 2012), which is largely controlled by surface freshwater fluxes (e.g., Sengupta et al., 2016;Wilson and Riser, 2016).In experiments to understand the relative contributions of various processes (Thompson et al., 2006;Vinayachandran and Nanjundiah, 2009;Akhil et al., 2014;Rahaman et al., 2014).Many GCMs have positive sea surface salinity biases of 2-3 psu over the freshwater-dominated northern BoB.However, factors responsible for such salinity biases in the BoB (both in the surface and in the subsurface) are unclear and need examination.Biases in subsurface salinity and temperature in GCMs can alter the ocean circulation, sea level, vertical mixing, and the coupling between ocean and atmosphere (e.g., Seo et al., 2009;Brown et al., 2013)

MODELS AND DATA
NCEP-CFSv2 is a fully coupled oceanatmosphere-land model with advanced physics and increased resolution compared to CFSv1 (Saha et al., 2006).Under the National Monsoon Mission of India, CFSv2 has been selected as an operational model for dynamical monsoon prediction over the Indian region (http:// www.tropmet.res.in/monsoon/index.php); this model has also been widely used for global-scale forecasts (Saha et al., 2014).The Global Forecast System's (GFS) atmospheric component has a horizontal resolution of T126 (~100 km) and 64 sigma layers vertically (Saha et al., 2014).The ocean component of CFSv2 is the GFDL-MOM4 model (Griffies et al., 2004).CFSv2 is integrated over a period of 100 years.The last 60 years are utilized for preparing monthly climatology for the present study.More details of the CFSv2

River runoff
Annual mean runoff (single time step; Dai and Trenberth, 2002) Monthly climatology (Dai et al., 2009) Interannual monthly river discharge (derived from data sets provided by Dai et al., 2009, andPapa et al., 2010) Data assimilation technique None None 3D VAR (Derber and Rosati, 1989) Variables assimilated None None Temperature and salinity profiles from Argo, TAO/TRITON, RAMA, PIRATA, XBT, CTD, and NDPB (Ravichandran et al., 2013) Evaporation (E) Computed from bulk formula using model SST and moisture from atmospheric model (GFS) Computed from bulk formula using model SST (COREv2) Computed from bulk formula using model SST

Model domain and time period
Global and last 60 years of 100 years free run is used free run are provided in Roxy (2014).
The stand-alone ocean model MOM5 is forced with the Coordinated Oceanice Reference Experiments version 2 (COREv2) surface forcing data sets for global ocean-ice modeling (Large and Yeager, 2004), developed at the National Center for Atmospheric Research (NCAR).In the model, turbulent heat fluxes are derived from the NCAR bulk formulae using the simulated SST and the prescribed 10 m atmospheric parameters (see Table 1 for details).Monthly climatological river discharge at discrete river mouth locations at 1° × 1° global grid (Dai and Trenberth, 2002) is used.
Finally, ocean reanalysis generated by INC-GODAS is also examined (Ravichandran et al., 2013).INC-GODAS (ocean model is MOM4) is an improved version of NCEP-GODAS (ocean model is MOM3) with a shorter assimilation window (Ravichandran et al., 2013).Temperature and salinity profiles are assimilated at six-hour intervals using all observations from the 10-day assimilation window, and the 3DVAR assimilation technique is used for data assimilation.et al., 1998) moored buoys, expendable bathythermographs (XBTs), and National Data Buoy Program (NDBP) moorings.We note that maintaining such varied observational networks (Venkatesan et al., 2013) is crucial if models are to be substantially constrained by data assimilation.Ravichandran et al. (2013) provide more details of INC-GODAS.The model is forced by momentum flux, heat flux, and freshwater flux from the NCEP-DOE atmospheric reanalysis version 2 (Kanamitsu et al., 2002;hereafter, NCEP R2).INC-GODAS provides three-dimensional ocean structures and is currently being used to provide initial ocean conditions for the operational model CFSv2 for both seasonal and extended range prediction of the Indian summer monsoon (e.g., Sahai et al., 2013).
Vertical mixing in all three models is prescribed using the standard K-profile parameterization (KPP).The KPP scheme (Large et al., 1994) is one of the most widely used parameterization schemes for vertical mixing in ocean GCMs.Unlike bulk mixed layer models, KPP does not assume a well-mixed boundary layer and explicitly predicts an ocean boundary layer depth.Within this boundary layer, the turbulent mixing is parameterized using a nonlocal bulk Richardson number and similarity theory.Below the boundary layer, KPP invokes eddy coefficients to parameterize vertical mixing due to shear instability (when model resolved shear produces small Richardson numbers), breaking of a background internal wave field (assumed to be constant and quite weak), and double diffusion.Table 1 provides details about forcing and specifications of mixing schemes for the three models used in the present study.
The simulated fields from the three models, CFSv2, MOM5, and INC-GODAS, are compared with World Ocean Atlas 2013 (WOA13) temperature and salinity products (Locarnini et al., 2013;Zweng et al., 2013).The latter is a reasonable choice, as it has been found to be consistent with RAMA buoy data (not shown) and Argo data (Fousiya et al., 2015).We also consider the Estimating the Circulation and Climate of the Ocean, Phase II (ECCO2; Menemenlis et al., 2008) reanalysis data, the TropFlux air-sea flux product (Praveen Kumar et al., 2010), and Global Precipitation Climatology Project (GPCP) monthly precipitation data (Adler et al., 2003).The ECCO2 product is used in the present study to estimate circulation bias in the models.It is an automatic differentiation tool that is used to calculate the adjoint code of MITgcm (Marshall et al., 1997;Azaneu et al., 2014).ECCO2 climatology is constructed from the period 1992-2007.
In the present study, the three models (CFSv2, MOM5, and INC-GODAS) are validated against the three data sets (WOA13, ECCO2, and TropFlux).To examine the causes for the salinity biases in these models, we compare the Brunt-Vaisala frequency (N 2 ; Gill, 1982), vertical shear of horizontal currents, and energy required for mixing (ERM).The latter is a metric of the potential energy of the upper water column (Shenoi et al., 2002).ERM is computed as where g is acceleration due to gravity (9.8 m s -2 ), ρ s is surface layer density (kg m -3 ), ρ b is density at any particular depth beneath the surface mixed layer (kg m -3 ), and h is mixed layer depth (m).ERM is a proxy for upper-ocean stratification.

UPPER-OCEAN SALINITY BIASES AND POSSIBLE CAUSES
Figure Figure 2 shows vertical profiles of temperature, salinity, and density during boreal summer and fall, averaged over the BoB for CFSv2, MOM5, and INC-GODAS and compared with WOA13.In general, the temperature is high near the surface and decreases with depth and vice versa for salinity and density.However, upperocean salinity in the BoB is much lower   than in the rest of the tropical Indian Ocean due to the large fresh water influx.CFSv2 displays a large cold bias for the upper 80 m of the water column, a positive bias in salinity from the surface to 100 m depth, and a positive bias in density up to 60 m depth (Figure 2d).In MOM5, the temperature bias is smaller, but salinity and density show a strong negative bias around 70 m depth.In the case of INC-GODAS, biases in temperature and density are high above 60 m (Figure 2f).During fall, the maximum positive bias in salinity is noted around 40 m in MOM5 and CFSv2 and negative bias in INC-GODAS (Figure 2j-l).Biases in temperature, salinity, and density are higher in fall compared to summer.We next examine freshwater fluxes, water column stability, and current profiles in the upper ocean in order to investigate the causes of the biases.Figure 3a shows WOA13 climatology and CFSv2 vertical profiles of the Brunt-Vaisala frequency (N 2 ) averaged over the BoB for summer.CFSv2 displays a too-weak stability for the upper 100 m as compared to WOA13, which is clearly seen in the bias (Figure 3d).The ERM profile is consistent with that for the entire tropical Indian Ocean (Chowdary et al., 2016b) and confirms the weak stability in CFSv2 (Figure 3b,e).Though upperocean salinity is better represented in MOM5 compared to CFSv2, the stability is still too weak in the top 60 m of the water column (Figure 3d,e).In the case of INC-GODAS, the stability is slightly higher than observed and is better represented compared to CFSv2 and MOM5.This is mainly due to the assimilation of both temperature and salinity profiles.Figure 3c shows vertical shear of horizontal currents in CFSv2, MOM5, and INC-GODAS.Shear bias is calculated with respect to the ECCO2 product.It is evident from the bias that vertical shear is higher in all models compared to ECCO2 (Figure 3f). Figure 4 shows evaporation (E) minus precipitation (P) in CFSv2, MOM5, and INC-GODAS and their biases over the BoB during summer and fall.The eastwest gradient in E−P is too strong in CFSv2 compared to observations in both summer and fall.Excess precipitation over evaporation is apparent along the Myanmar coast in CFSv2 during summer, a bias that has also been reported in earlier studies (e.g., Parekh et al., 2016;Chowdary et al., 2016a).However, precipitation in general is underestimated over most of the BoB, which arises from the biased rainfall simulation from the GFS.This strong positive E−P in the model would contribute to the high salinity bias, especially in the surface layer of the BoB in CFSv2 (Figure 2e).Similarly, CFSv2 displays a positive bias in freshwater flux during fall, which is responsible for the positive surface salinity bias (Figure 2k).Consistent with other coupled model studies, precipitation in the atmospheric component of CFSv2 contains strong biases, with overestimation of evaporation resulting in poor representation of surface salinity in the BoB (e.g., Seo et al., 2009).MOM5 displays a weak bias in E−P, which is reflected in a weak positive bias in surface salinity  (Figures 4c and 2e).INC-GODAS also shows a weak positive bias in E−P in most of the area during both summer and fall (Figure 4c).Evaporation changes in both models are partly dependent on simulated SSTs.In response to positive bias in E−P, INC-GODAS shows a weak positive surface salinity bias (Figure 2e).The weak biases in INC-GODAS and MOM5 are partly due to high evaporation, while the strong bias in CFSv2 is due to the bias in the model rainfall.Therefore, accurate representation of the net freshwater flux in the coupled models is critical for maintaining a more realistic upper-ocean salinity structure over the BoB, especially at the surface.Horizontal advection is found to be important for salinity variations in the BoB for most of the year (e.g., Rao and Sivakumar, 2003;Sengupta et al., 2016, and references therein).During summer, ECCO2 shows some regions of northeastward or southeastward surface current averaged for the upper 40 m in the BoB (Figure 5a); these are consistent with Ekman transport by southwesterly monsoonal winds as an important component of average flow.The southward component of the surface current is weak in all models west of 86°E (Figure 5b-d).Westward current bias in the central BoB is particularly strong in CFSv2 and MOM5 and slightly weak in INC-GODAS, suggesting weaker eastward mean flows in all models.This may also reflect a failure to properly capture the surface-confined Ekman transports.During fall, a basinscale cyclonic gyre is seen in ECCO2 surface currents (Figure 5e).The southward current component near the east coast of India (EICC) is strong in ECCO2 but is weaker in CFSv2, MOM5, and INC-GODAS (Figure 5f-h).Coastal current systems such as the EICC are not very well captured in models, mainly due to their coarser resolution.Overall, a tooweak cyclonic gyre during fall is evident in many models compared to ECCO2. Figure 6 illustrates mean meridional salinity advection at the southern boundary of the BoB (10°N) averaged from 85°E-90°E.As compared to ECCO2, during summer, CFSv2 overestimates salinity advection due to monsoon currents at the surface and in the subsurface (Figure 6a,b).This is consistent with the salinity biases in CFSv2 during summer (Figure 2).MOM5 also displays positive biases in salinity advection but bias is weaker (Figure 6c).In INC-GODAS, the strong negative bias in salinity advection in the upper ocean during spring and positive bias during summer are apparent (Figure 6d).Discrepancies in salinity advection or inflow and outflow of saline water from the BoB in GCMs also contribute to the salinity bias.

DISCUSSION
The model salinity bias in the BoB could be due to a combination of inaccurate river forcing or freshwater forcing and the incorrect physical representation of the mixing.Our analysis reveals that the salinity bias, especially in the subsurface,  June-July-August-September October-November does not arise primarily from the misrepresentation of river discharge.For example, we considered models with mean annual (CFSv2) and mean seasonal (MOM5, INC-GODAS) varying river discharges, but all three models show similar biases in subsurface structure, suggesting improper vertical mixing in these models is a likely culprit.
Although the stand-alone model (MOM5) was forced with seasonally varying river discharge data, there is no significant reduction in salinity biases, as anticipated by previous modeling studies (e.g., Durand et al., 2011).This issue is important and warrants further investigation.Because river discharges are a result of complex estuarine dynamics that are absent in ocean models, a dedicated observational campaign and sensitivity studies are needed to improve understanding of freshwater forcing due to river discharge (e.g., Mahadevan et al., 2016, in this issue).For example, river insertion thickness in the ocean models might be a major factor, as it represents the thickness of the water column into which tracers from rivers are inserted.Tracer insertion is usually done in the mixed layer, so any misrepresentation of mixed layer thickness in the models would result in unrealistic insertion of river water.Another potentially important problem is inability to treat river water properties such as temperature, salinity, and turbidity.Of course, it is of utmost importance to develop better river discharge data at higher temporal resolution.When combined with in situ observations, satellite altimetry products can be useful for estimating river discharges at higher temporal frequency.Overall, a lack of accurate estimates of river runoff (amount, time scale, and properties) has resulted in a lack of consensus in the modeling community on the specification of river runoff; this remains the key challenge for numerical modeling of BoB circulation.
Ocean models with vertical shear parameterizations can successfully reproduce many general features of the tropical ocean (Packnowski and Philander, 1981).However, the BoB has several unusual features that may influence both the balance of dynamics that set turbulent mixing rates and upper-ocean stratification and the most appropriate parameterization strategies.There are two categories of concern here: ways in which KPP does not accurately represent the in situ vertical mixing processes, and processes that are simply not included or represented in existing models.In the former category of concern, one of the primary issues is resolution.In the BoB, many of the observed fresh water layers are on the order of 10 m thick or even thinner (Lucas et al., 2016, in this issue); models that have a 10 m vertical resolution will never be able to represent such layers accurately, and will always mix "too much." High priority should be thus given to increasing vertical resolution in the upper ocean for operational forecast models.Additionally, even with sufficient resolution, it is not clear that KPP, tuned for weaker stratification, will accurately represent entrainment and erosion at the extremely sharp haloclines found in the upper BoB in certain seasons (Lucas et al. and Jinadasa et al., 2016, both in this issue;Sengupta et al., 2016).Finally, KPP does not include a variety of vertical mixing processes, such as Langmuir turbulence, and this is likely important in mixing both scalers and momentum, though efforts are underway to remedy that omission.
In the second category of concern, which is accurate representation of in situ vertical mixing, a recent observational campaign and numerical process study have highlighted the importance of a variety of fundamentally three-dimensional phenomena that may play important roles in setting stratification but are not generally represented in the current generation of GCMs.Much of the surface layer in the BoB is freshwater stratified, so subject to the "interior ocean mixing" portion of KPP.Internal wave breaking in this layer is likely much more related to locally wind-generated near-inertial internal waves than the background levels assumed in KPP (Johnston et al., 2016, in this issue), and needs to be parameterized differently.A recent approach developed by Jochum et al. (2013) to explicitly include the effects of near-inertial motions for upperocean mixing in the tropics has shown great promise.In general, the shear that produces subcritical Richardson numbers in the ocean is almost never directly related to the types of mesoscale features resolved by the models presented here; it is instead related to internal waves, internal tides, and frontal processes, each of which must be properly parameterized (MacKinnon et al., 2013).
Many of the observations strongly suggest submesoscale instabilities (see Lucas et al. and Sarkar et al., 2016, both in this issue;Wijesekera et al., in press).
At times, such instabilities induce significant subduction of one water mass beneath another, increasing local stratification and inhibiting vertical mixing; lack of inclusion of such re-stratifying processes may be one reason for the low stratification/ high mixing model biases reported here.Many submesoscale instabilities preferentially occur at sharp lateral density fronts, of the sort not explicitly resolved by the models considered here; hence, their effects must be parameterized in a different manner.Fox-Kemper et al. (2008a,b, 2011) parameterized one element in this family of restratifying processes; however, it is not yet clear whether the particular baroclinic instability represented in that scheme is prevalent in the BoB.Ongoing close collaboration between observationalists, process modelers, and the GCM community, of the sort facilitated by the Air-Sea Interactions Regional Initiative-Ocean Mixing and Monsoon (ASIRI-OMM) experiments, is imperative for development of the next generation of successful forecast models for the Bay of Bengal.

SUMMARY
Upper-ocean stratification in the Bay of Bengal is highly dependent upon salinity, which is substantially lower to the north of 14°N than to the south.Representation of BoB salinity in GCMs is challenging because of the complex ocean processes that need to be accurately captured by the models.In this study, we examined the seasonal mean biases in BoB upperocean structure, with a focus on salinity fields, in a coupled model ( CFSv2), an ocean model (MOM5), and a reanalysis product (INC-GODAS) during summer and fall.A positive salinity bias of about 1 psu is noted north of 14°N in the BoB in CFSv2 over the upper 20 m of the water column.Area-averaged salinity over the BoB is high up to 100 m depth, and there is a maximum bias at the surface.At the same time, the temperature bias is negative.In MOM5, a negative (weakly positive) bias in salinity and density around 70 m (at the surface) is noted, and temperature bias is low.The INC-GODAS product displays relatively low bias in temperature and salinity compared to CFSv2 and MOM5.
The processes responsible for upperocean salinity biases such as fresh water flux, vertical mixing, and horizontal advection were examined in all models.Excessive evaporation over precipitation and prescribed annual mean river discharge in CFSv2 are mainly responsible for the positive salinity bias.MOM5 showed low-salinity bias in the upper 30 m, and this could be due to seasonally varying river discharge and prescribed E−P forcing.However, the salinity bias is higher in the subsurface at around 80 m in summer and 40 m in fall.Salinity is better represented in INC-GODAS in the upper 100 m over the BoB, due to assimilation of temperature and salinity profiles from Argo and other observations.CFSv2 and MOM5 show weaker stability with positive bias in density within the top few meters.Furthermore, in both models, vertical shear of the horizontal current is higher compared to ECCO2.This suggests excess mixing of lowsalinity water from the surface with subsurface water in MOM5.In the model, a low Richardson number triggers excess vertical mixing.INC-GODAS salinity is better represented over BoB, but vertical shear of horizontal currents is not well captured.Furthermore, advection of high-salinity water in CFSv2 at the southern boundary of the BoB is consistent with the high-salinity bias.In MOM5, the salinity advection bias is weaker, and the converse is true for INC-GODAS.These results suggest that treatment of the river forcing (vertical mixing, freshwater flux, and lateral advection) is an important factor for improving representation of the BoB salinity fields in coupled and standalone ocean models.To estimate biases accurately, high-frequency spatiotemporal observations over the BoB from field campaigns such as those described in this volume would be very useful.

Indo
5° × 0.5° and vertically 10 m in the upper 220 m Horizontally 0.25° × 0.25° and vertically 10 m in the upper 220 m. 0.5° degree in zonal and meridional; meridional resolution is 0.25° within 10°S to 10°N and vertical resolution is 10 m in the upper 220 m Radiation fluxes RRTM (Rapid Radiative Transfer Mode) radiation scheme Daily fluxes from COREv2 Daily NCEP-R2 fluxes

FIGURE 1 .
FIGURE 1. Depth-latitude plots of salinity (mean; shaded, psu) and bias (black contours, psu) averaged between 85°E and 90°E over the Bay of Bengal region in National Centers for Environmental Prediction-Climate Forecast System version 2 (CFSv2; left column), Modular Ocean Model version-5 (MOM5; middle column), and Indian National Centre for Ocean Information Services Global Ocean Data Assimilation System (INC-GODAS; right column) for annual mean (a, f, k), March-April-May (b, g, l), June-July-August-September (c, h, m),October-November (d, i, n), and December-January-February (e, j, o).Biases are calculated with respect to World Ocean Atlas 2013 (WOA13).White contours represent the WOA13 salinity.

FIGURE 2 .FIGURE 3 .
FIGURE 2. Area averaged Bay of Bengal vertical profiles of (a) temperature (°C), (b) salinity (psu), and (c) density (kg m -3 ) for CFSv2, MOM5, INC-GODAS, and WOA13 during summer (June-July-August-September).(d-f) Plots show respective biases.Similar to June through September, October-November vertical profiles of temperature, salinity, and density and their biases are displayed in (g) to (l).Vertical profiles are displayed for the upper 100 m.Bias is calculated with respect to WOA13.
. Thus, accurate representation of the vertical structure of salinity in GCMs is important.