## 1.

## Introduction

Short wave radiation (SWR) and its attenuation with depth have a major impact on the vertical distribution of oceanic water temperature, dynamical processes, and ocean–atmosphere interactions. SWR is a sum of the ultraviolet–visible wavelengths [350 to 700 nm, photosynthetically available radiation (PAR)] and infrared (IR) wavelengths (700 to 2500 nm). The IR is usually absorbed in the upper few meters of water. However, PAR can penetrate below 50 m in clear water. In existing physical oceanic models,^{1}2.3.4.5.6.^{–}^{7} the temperature is modeled as a tracer, when its temporal change is described by the advection–diffusion-source equation. The PAR and its attenuation with depth represent a source (forcing) term in the governing equation for the temperature. As a result, in physical ocean models, there are three types of errors related to PAR and its attenuation with the depth: (1) errors in specification of the surface values of PAR, (2) errors in accuracy of the optical model for the PAR attenuation with the depth, and (3) errors in attenuation coefficients of PAR. The remote sensing community has been developing satellite-derived products of PAR and inherent optical properties (IOPs, such as absorption and backscattering, diffuse attenuation coefficient at 490 nm), which are used to estimate PAR attenuation coefficients.^{8}9.^{–}^{10} Also, there are well-established optical models, such as HydroLight and EcoLight,^{11}12.^{–}^{13} to provide accurate attenuation of PAR with the depth. All of these products can certainly improve modeling of PAR and its attenuation with depth in physical oceanic models. However, at present, there is a very limited use of remotely sensed PAR-related products in oceanic models.

In numerical modeling of oceanic processes, the PAR usually comes from the atmospheric model predictions, and in numerous studies it has been reported that atmospheric models show a tendency to overestimate the SWR due to the underestimation of predicted low-level clouds.^{14}15.16.17.18.^{–}^{19} The use of PAR from the atmospheric model is mostly dictated by already in place coupling of the oceanic model to the atmospheric model outputs and also by the fact that the atmospheric model provides forecasts of PAR values that can be used to estimate the corresponding forecast of oceanic conditions, whereas remotely sensed PAR products provide only hindcast or nowcast of PAR values.

The application of optical models based on radiative transfer theory has demonstrated promising results in the idealized channel flow model.^{11} However, optical models based on radiative transfer theory have not found wide applications in the numerical oceanic modeling community, mostly because of a high increase in computational cost when they are used with high resolution, data assimilative physical models. For this reason, most existing oceanic physical models^{1}2.3.4.5.6.^{–}^{7} use a simplified optical model based on an exponential attenuation of PAR with the depth. Specifications of the PAR attenuation coefficient in exponential attenuation schemes are mostly based on: the Jerlov’s water-types classification,^{20} climatological estimates of attenuation coefficients, or from the biological model predictions of light-absorbing and scattering water constituents.^{9}^{,}^{21}22.23.24.^{–}^{25} All of these PAR attenuation coefficient specifications are prone for introducing errors in the attenuation of SWR with the depth.

It is clear that remotely sensed PAR-related products should find wider applications in the oceanic physical models. One way is to assimilate PAR-related products into the oceanic model, when remotely sensed products are merged together with the oceanic model predictions of PAR and its attenuation. In existing data assimilation approaches, this merge is based on estimates of errors in remotely sensed products and corresponding errors in the PAR values and its attenuation in the oceanic model. To achieve this merge, remote sensing and numerical oceanic modeling communities have to better understand the sensitivity of the model heat predictions to the errors in PAR, optical model, and PAR attenuation coefficient.

The sensitivity of the heat modeling to errors introduced by the use of a simplified exponential attenuation scheme was investigated in Ref. 11.

In the present paper, our objective is to study the sensitivity of the modeled upper ocean heat content to other mentioned above two types of errors: those related to the specification of the magnitude of PAR at the surface (from atmospheric model) and those related to the incorrect specification of the PAR attenuation coefficient with depth (prescribed in the oceanic model) in the simplified exponential attenuation model.

Among the questions addressed in this study are: What errors (in surface PAR or in its attenuation), and under what conditions, have a greater impact on the modeled oceanic heat content of the area of interest? Under what conditions do the two considered types of errors combine to reinforce the heat content error? Under what conditions do the PAR errors work in opposition, such that they reduce the total error in the heat content? How do the heat content sensitivities to errors in PAR and its attenuation change in time and space?

To address the questions above, we have developed an approach for estimating the sensitivity of the oceanic model heat content to changes in surface values of PAR and its attenuation with the depth. The approach is based on using the adjoint to the oceanic model governing temperature equation. These sensitivities are valid for any model and modeling area where the temperature change is described as tracer with advection–diffusion-source model. Based on derived sensitivities, we have estimated errors in the upper ocean heat content due to errors in surface PAR values and due to errors in the PAR attenuation coefficient.

The structure of the paper is as follows: Sec. 2 describes the modeling system, observations, and the approach for estimating the model heat content sensitivities. Section 3 is devoted to the results. Discussions and conclusions are presented in Sec. 4.

## 2.

## Methods

## 2.1.

### Modeling System

In the present study, the Monterey Bay model [the domain is shown in Fig. 1(a)] is based on the Navy Coastal Ocean Model (NCOM)^{2}^{,}^{3} and is triply nested inside the global and regional (California Current) NCOM-based models.^{17}^{,}^{18} The model, which is set up on a curvilinear orthogonal grid with resolution ranging from 1 to 4 km, uses a sigma vertical coordinate system with 30 levels. On open boundaries, the model is one-way coupled to a larger scale regional (California Current) NCOM-based model.^{17} The Navy Coupled Ocean Data Assimilation system^{26} is used for assimilation of observations into the Monterey Bay model.^{27} The oceanic model is forced with surface fluxes from the Coupled Ocean-Atmosphere Mesoscale Prediction System (COAMPS) at 3 km horizontal resolution.^{28} The COMAPS shortwave calculation considers absorption by water vapor and ozone. Two types of clouds are considered in the radiative transfer calculation: stratiform and cumulus.

The Monterey Bay model outputs were validated against temperature, salinity, and currents profiles at moorings, high frequency (HF) radar-derived currents, and aircraft sea surface temperatures (SSTs).^{18}^{,}^{27}

## 2.2.

### Approach for Estimation of the Model Heat Content Sensitivities

The governing equation for temperature has the following form in oceanic models:^{1}^{,}^{2}3.4.5.6.^{–}^{7}

## (1)

$$\frac{\partial T}{\partial t}=-\mathbf{U}\nabla T+\mathrm{div}(k\nabla T)-\frac{Qsr}{{\rho}_{0}{c}_{p}}\frac{\partial \gamma}{\partial z},$$Initial condition for Eq. (1) at $t={t}_{0}$ has the form

where $D$ is the modeling domain and $\tau (x,y,z)\in D$.In Eq. (1), the $Qsr$ is a sum of the visible–ultraviolet bands (350 to 700 nm, PAR) and IR bands (700 to 2500 nm)

where $f$ represents a fraction of the visible–ultraviolet bands (PAR) in total solar radiation.The vertical distribution of the PAR and IR is described by two separate exponential functions (see for example Ref. 9)

## (4)

$$Qsr\gamma ={Q}_{\mathrm{PAR}}\text{\hspace{0.17em}}\mathrm{exp}(-{K}_{\mathrm{PAR}}z)+{Q}_{\mathrm{IR}}\text{\hspace{0.17em}}\mathrm{exp}(-{K}_{\mathrm{IR}}z),$$Let us consider the heat content $J$ at verification time ${t}_{f}>{t}_{0}$ in the subdomain $\mathrm{\Omega}$ (area of interest) of the modeling domain $D$ (our interest is in the upper ocean model heat content)

where $d\tau $ is the element of volume in the integral [Eq. (5)].As emphasized in Refs. 9 and 29, ${K}_{\mathrm{PAR}}$ varies with depth and horizontally. In Refs. 21 and 22, ${K}_{\mathrm{PAR}}$ is empirically modeled as functions of chlorophyll while in Ref. 9, ${K}_{\mathrm{PAR}}$ is empirically modeled as functions of IOPs (absorption and backscattering at 490 nm), which are derived from the satellite remotely sensed reflectances. Profiles of absorption and backscattering at 490 nm predicted by the biochemical model are used to estimate ${K}_{\mathrm{PAR}}$ in Ref. 25. It is known that there is high level of uncertainty in estimation of chlorophyll and IOPs from remotely sensed reflectances, as well as, based on IOP derivations from the biochemical model simulations. For example, according to Ref. 30, “the satellite data product accuracy generally accepted by the international missions is $\pm 5\%$ for water-leaving radiances and $\pm 35\%$ for chlorophyll a in the open ocean.” These relatively large uncertainties make the analytical analysis of errors in the heat content due to errors in ${K}_{\mathrm{PAR}}$ challenging, particularly when ${K}_{\mathrm{PAR}}$ varies horizontally as well as vertically. To simplify the analysis of the heat content sensitivity to errors in ${K}_{\mathrm{PAR}}$, we assume that ${K}_{\mathrm{PAR}}$ can vary only horizontally. This includes the case when values of $f$, ${K}_{\mathrm{PAR}}$, and ${K}_{\mathrm{IR}}$ in Eq. (4) are chosen based on the classification of Jerlov’s^{20} (as interpreted by Ref. 31), which is currently employed in many oceanic models.^{1}^{,}^{3}^{,}^{5}^{,}^{6} Oceanic waters are classified based on its clarity, and five values of ${K}_{\mathrm{PAR}}$ are defined ranging from clear to increasingly turbid water (Table 1). The case of horizontally varying ${K}_{\mathrm{PAR}}$ also includes a specification of ${K}_{\mathrm{PAR}}$ based on satellite observations of chlorophyll (see for example Ref. 23).

## Table 1

Parameters for Jerlov’s seawater optical types (as interpreted by Paulson and Simpson31). The ratio of sensitivities SKPARSQPAR [Eqs. (17)–(18)] estimated for Jerlov’s seawater optical types.

Water type | I | IA | IB | II | III |
---|---|---|---|---|---|

$f$ | 0.42 | 0.3785 | 0.3276 | 0.2330 | 0.22 |

${K}_{\mathrm{PAR}}$ | 0.0435 | 0.0491 | 0.0579 | 0.0685 | 0.1265 |

${K}_{\mathrm{IR}}$ | 2.86 | 1.67 | 1. | 0.67 | 0.72 |

$\frac{{S}_{{K}_{\mathrm{PAR}}}}{{S}_{{Q}_{\mathrm{PAR}}}}$ | 0.53 | 0.51 | 0.45 | 0.38 | 0.14 |

Appendix presents the derivations of the errors in the model heat content [Eq. (5)] due to errors in surface PAR and its attenuation with the depth. Equations (15) and (16) provide a representation of errors $\delta {J}_{{Q}_{\mathrm{PAR}}}$ and $\delta {J}_{{K}_{\mathrm{PAR}}}$ in the heat content [Eq. (5)] due to normalized errors ${\u03f5}_{{Q}_{\mathrm{PAR}}}$ [Eq. (19)] and ${\u03f5}_{{K}_{\mathrm{PAR}}}$ [Eq. (20)] in surface PAR values (${Q}_{\mathrm{PAR}}$) and in the attenuation coefficient ${K}_{\mathrm{PAR}}$ correspondingly. Equations (15)–(20) provide an estimate of errors in the heat content $J$ [Eq. (5)] at verification time $t={t}_{f}$ due to relative errors in ${Q}_{\mathrm{PAR}}$ and ${K}_{\mathrm{PAR}}$ at any time $t$ between initial time ${t}_{0}$ and verification time ${t}_{f}$. Also, Eqs. (15)–(20) present an estimate of errors in the heat content $J$ [Eq. (5)] due to relative errors in ${Q}_{\mathrm{PAR}}(s,t)$ and ${K}_{\mathrm{PAR}}(s,t)$ at any location on the surface $S$ of the modeling domain D ($s\in S$), even outside of the domain of interest $\mathrm{\Omega}$, where the heat content $J$ [Eq. (5)] is estimated.

## 2.3.

### Observations

## 2.3.1.

#### Satellite ocean color data

The moderate-resolution imaging spectroradiometer (MODIS)-Aqua satellite imagery was processed using the Naval Research Laboratory (NRL) automated processing system (APS). APS is a complete end-to-end system that includes sensor calibration, atmospheric correction (with near-IR correction for coastal waters), and bio-optical inversion. APS incorporates, and is consistent with, the latest NASA MODIS code (SeaDAS).^{32}^{,}^{33} The diffuse attenuation coefficient *Kd* (490) was estimated in accord with Ref. 9 with 1-km resolution for August 15, 2003. Equation (9) from Ref. 10 was used to estimate ${K}_{\mathrm{PAR}}$ from *Kd* (490)

*Kd*(490) with 9-km resolution from Aqua-MODIS imagery for August 2003 were obtained from Ref. 34. The same above Eq. (9) from Ref. 10 was used to estimate monthly values of ${K}_{\mathrm{PAR}}$ from

*Kd*(490).

Both above ${K}_{\mathrm{PAR}}$ distributions (August 15 and monthly mean) were interpolated to the model grid.

## 2.3.2.

#### Mooring data

Observations of winds and PAR from the Monterey Bay Aquarium Research Institute surface moorings M1 (122.02° W, 36.74° N) and M2 (122.40° W, 36.67° N) are used in this study [Fig. 1(a)]. Near-surface 3-m wind speed and direction were measured by a MetSys monitor. The observed PAR was measured by the spectroradiometer mounted on moorings $\sim 3\text{\hspace{0.17em}\hspace{0.17em}}\mathrm{m}$ above the water surface.^{35}

## 2.3.3.

#### HF radar surface currents

Surface current observations used in this study were derived from the California Coastal Ocean Current Mapping Program’s HF radar network. Surface currents were estimated based on inputs from five HF radar sites for August 2003. Vector currents were estimated on a Cartesian grid with a horizontal resolution of 3 km by computing the best-fit vector velocity components using all radial velocity observations within a radius of 3 km for each grid point each hour.^{36}

## 3.

## Results

Comparisons of COAMPS predictions with observations from moorings and aircraft surveys were reported in Refs. 18, 27, 28, 37, and 38. The comparisons during August to September of 2003 demonstrated high complex correlations between the observed and COAMPS-predicted wind velocities, and good agreement in the magnitudes and the extent of each observed upwelling event. We estimated the magnitude of the complex correlation coefficient and angular displacement^{39} between observed and COAMPS model 10 m wind velocities at the location of the mooring M1 [Fig. 1(a)]. The estimated value of complex correlation between observed and model wind velocities is 0.79, which is significant at a 95% confidence level. The angular displacement gives the average counterclockwise angle of the COAMPS wind velocity with respect to the observed wind velocity. The angular displacement between model and observed wind velocities is small and equals $-11\text{\hspace{0.17em}\hspace{0.17em}}\mathrm{deg}$.

At the same time, comparisons of the observed and model-predicted daily-averaged values of PAR [Fig. 1(c)] show an overestimation of the PAR in the COAMPS predictions, especially during the wind relaxation events of August 20 to 22 and September 1 to 3 (overestimation of mean daily PAR values ranges from 20 to $80\text{\hspace{0.17em}\hspace{0.17em}}\mathrm{W}/{\mathrm{m}}^{2}$). The excessive PAR is likely related to the modeling of low-level clouds. Prediction of extensive low-level clouds in the Monterey Bay area during summer, which is believed to be underestimated by COAMPS during this time period,^{28} is a very challenging problem.

To address the research questions raised at the end of Sec. 1, we need to establish “true” values of PAR (${Q}_{\mathrm{PAR}}$) and ${K}_{\mathrm{PAR}}$ against which we estimate errors in the modeled PAR and the model-used values of the PAR attenuation coefficient. The observed PAR values at moorings M1 and M2 can be used to establish magnitudes of errors between the observed and the COAMPS model-predicted values of PAR. In the present study, only satellite-derived values of ${K}_{\mathrm{PAR}}$ are available to us in order to estimate errors in the Monterey Bay model specification of ${K}_{\mathrm{PAR}}$ (from the Jerlov’s water-types classification or climatological estimates of attenuation coefficients). This introduces challenges, because satellite-derived ${K}_{\mathrm{PAR}}$ values are prone to their own errors, and satellite images of ${K}_{\mathrm{PAR}}$ are available mostly during clear days, when errors in the COAMPS model predictions of PAR represent the lower bound of possible errors, whereas during cloudy days, when the model predictions of PAR represent the upper bound of errors, the satellite images of ${K}_{\mathrm{PAR}}$ are not available. To mitigate this, we present results for errors in PAR corresponding to the clear day (specifically August 15) and to the monthly mean of PAR for the month of August. In this case, our analysis of the heat content errors will be for the case of the lower bound error in the COAMPS model-predicted PAR values and for the monthly mean error in the modeled PAR.

Table 2 provides comparisons between daily-averaged modeled and observed surface PAR (${Q}_{\mathrm{PAR}}$) values for August 15, 2003, as well as for the August 2003 monthly means at moorings M1 and M2. The errors in the modeled values of ${Q}_{\mathrm{PAR}}$ in comparison to observed values are 7.4% and 17% at moorings M1 and M2 for August 15, 2003 and 27% and 37% at M1 and M2 for the August mean values.

## Table 2

Observed and modeled surface PAR (W/m2) comparisons at moorings M1 and M2.

Daily averaged August 15, 2003 | Mean August 2003 | |||||||
---|---|---|---|---|---|---|---|---|

M1 | Error (%) | M2 | Error (%) | M1 | Error (%) | M2 | Error (%) | |

Observed | 135 | — | 115 | — | 106 | — | 96 | — |

Modeled | 145 | 7.4 | 135 | 17 | 135 | 27 | 132 | 37 |

Two model runs were conducted over the considered time period. Run 1 is the model run forced with SWR from the COAMPS predictions. In run 2, the mean daily values of COAMPS SWR were corrected in accord with observed values at the M1 mooring

where ${Q}_{SR}^{ADJ}(x,y)$ is the adjusted mean daily values of SWR at location $(x,y)$, ${Q}_{SR}^{O}(M1)$ and ${Q}_{SR}^{M}(M1)$ are observed and model-predicted daily mean values of SWR at mooring M1.For both runs, values of $f$, ${K}_{\mathrm{PAR}}$, and ${K}_{\mathrm{IR}}$ in Eq. (4) are chosen based on Jerlov’s water type “IA” (Table 1). In both runs, no observations were assimilated.

Comparisons of the observed and model-predicted temperatures at 1- and 5-m depths are shown in Figs. 1(d) and 1(e). There is a good agreement between model and observations during strong upwelling favorable winds, when the model dynamics are mostly wind-driven and the solar heating was weaker than other days in August. However, during short wind relaxation events of August 20 to 22 and September 1 to 3, there are substantial differences between model and observed temperatures for run 1, when the model is forced with the surface SWR values without corrections. The run 1 temperatures are much warmer than observed temperatures at both considered depths. However, there is a reasonably good agreement between the model and observed temperatures for run 2, forced with corrected SWR in accord with the observations [Eq. (6)]. According to Table 3, the bias of model temperature predictions at 1-m depth is reduced from 1.83 (run 1) to 0.47°C (run 2), and root-mean-square error (RMS) is reduced at half. This demonstrates that the oceanic model deficiencies in predictions of observed surface and subsurface temperature during wind relaxation events are correlated strongly with the fact that the surface SWR values from the atmospheric model are larger than observed values at mooring locations.

## Table 3

Comparisons of model temperature predictions (1-m depth) at M1 during the relaxation event of August 20 to 22, 2003.

Bias (°C) | RMS (°C) | |
---|---|---|

Run 1 | 1.83 | 2.25 |

Run 2 | 0.47 | 0.86 |

As we have mentioned before, in many oceanic models applications, the ${K}_{\mathrm{PAR}}$ values are chosen either based on the Jerlov classification (specifically, using the attenuation coefficient for the Jerlov type IA in Table 1), or ${K}_{\mathrm{PAR}}$ values are derived from climatological (such as monthly mean) data sets. Table 4 illustrates possible errors in such specifications of ${K}_{\mathrm{PAR}}$ for a specific day—August 15, 2003. According to Table 4, monthly mean values of ${K}_{\mathrm{PAR}}$ represent more turbid, murkier water at both moorings locations in comparison to the ${K}_{\mathrm{PAR}}$ for August 15, 2003. At the same time, the Jerlov’s type IA water (Table 4, ${K}_{\mathrm{PAR}}\text{\hspace{0.17em}}=\text{\hspace{0.17em}}0.049$) represents a more clear, bluish water then for the observed ${K}_{\mathrm{PAR}}$ from the Aqua-MODIS data for the August 15, 2003. Spatial comparisons of monthly mean ${K}_{\mathrm{PAR}}$ with the values for the August 15 are shown in Fig. 2. For most locations in the domain, the monthly mean of ${K}_{\mathrm{PAR}}$ represents more turbid waters than the ${K}_{\mathrm{PAR}}$ for the particular day, August 15, and the error is reaching almost 300% inside the bay.

## Table 4

Comparisons of the PAR attenuation coefficients with the KPAR (m−1) derived from the Aqua-MODIS imagery for August 15, 2003).

M1 | Error (%) | M2 | Error (%) | |
---|---|---|---|---|

Aqua-MODIS, August 15 | 0.18 | — | 0.16 | — |

Aqua-MODIS, mean August | 0.53 | 196 | 0.21 | 31 |

Jerlov type IA | 0.049 | $-72$ | 0.049 | $-70$ |

To quantify possible errors in the upper ocean model heat content [Eq. (5)] due to errors in surface PAR and its attenuation with the depth, we consider now a specific area of interest $\mathrm{\Omega}$ in Eq. (5) as the area of $3\times 3$ horizontal grid cells (approximately an area of $4\text{\hspace{0.17em}\hspace{0.17em}}\mathrm{km}\times 4\text{\hspace{0.17em}\hspace{0.17em}}\mathrm{km}$) around the M1 mooring down to a depth of 25 m [Fig. 3(a), the area of interest is the box around mooring M1].

The depth of the area of interest is chosen at 25 m because the first optical depth for the clearest water in the Jerlov’s classification is around 23 m, therefore 25 m is deeper than the first optical depth for all water types considered in Table 1. For the verification time ${t}_{f}$, we select August 15, 2003 (verification time ${t}_{f}$ is equal 2003081500). Equation (15) in Appendix estimates errors $\delta {J}_{{Q}_{\mathrm{PAR}}}$ in the model heat content [Eq. (5)] due to normalized errors ${\u03f5}_{{Q}_{\mathrm{PAR}}}$ [Eq. (19)] in surface PAR values (${Q}_{\mathrm{PAR}}$) while Eq. (16) estimates errors $\delta {J}_{{K}_{\mathrm{PAR}}}$ in the model heat content [Eq. (5)] due to normalized errors ${\u03f5}_{{K}_{\mathrm{PAR}}}$ [Eq. (20)] in the attenuation coefficient ${K}_{\mathrm{PAR}}$. Table 2 provides estimates of ${\u03f5}_{{Q}_{\mathrm{PAR}}}$ [Eq. (19)], and Table 4 provides estimates of ${\u03f5}_{{K}_{\mathrm{PAR}}}$ [Eq. (20)].

Based on those errors and Eqs. (15) and (16), Table 4 shows the ratios of errors $\delta {J}_{{K}_{\mathrm{PAR}}}/\delta {J}_{{Q}_{\mathrm{PAR}}}$ in the model heat content. In the first example of Table 5, the error ${\u03f5}_{{Q}_{\mathrm{PAR}}}$ [Eq. (19)] is 7.4% (based on observed and atmospheric model-predicted surface PAR values at mooring M1, Table 2), and the error in attenuation coefficient ${\u03f5}_{{K}_{\mathrm{PAR}}}$ [Eq. (20)] is $-72\%$ (based on difference between Jerlov’s type IA attenuation coefficient and ${K}_{\mathrm{PAR}}$ derived from the Aqua-MODIS for August 15, 2003). The error in ${K}_{\mathrm{PAR}}$ is negative (the attenuation coefficient value in the model is smaller than observed ${K}_{\mathrm{PAR}}$, i.e., shifted toward more clear water), and the ratio of errors in heat content $\delta {J}_{{K}_{\mathrm{PAR}}}/\delta {J}_{{Q}_{\mathrm{PAR}}}$ is also negative ($-0.51$). In example 2 of Table 5, the error ${\u03f5}_{{Q}_{\mathrm{PAR}}}$ is positive but larger than in example 1, while the error ${\u03f5}_{{K}_{\mathrm{PAR}}}$ is negative and is the same as in example 1. As a result, the ratio $\delta {J}_{{K}_{\mathrm{PAR}}}/\delta {J}_{{Q}_{\mathrm{PAR}}}$ is still negative ($-0.14$), and there is also a reduction of the total error in the model heat content. Based on examples 1 and 2, the overestimation of PAR leads to the positive error in the heat content, but underestimation of ${K}_{\mathrm{PAR}}$ (shifted in the model toward more clear water) results in the penetration of this excessive PAR surface heating deeper into the water column and leads to a reduction of the total error in the heat content of the upper 25 m. In other words, the ${Q}_{\mathrm{PAR}}$ and ${K}_{\mathrm{PAR}}$ errors offset each other to some extent.

## Table 5

Estimates of errors in heat content at the mooring M1 (top 25 m).

Examples | QPAR | QPAR | ϵQPAR | KPAR | KPAR | ϵKPAR | δJKPARδJQPAR |
---|---|---|---|---|---|---|---|

Observed W/m2 | Modeled W/m2 | (%) | Observed m−1 | Modeled m−1 | (%) | ||

1 | 135 | 145 | 7.4 | 0.18 | 0.049 | $-72$ | $-0.51$ |

2 | 106 | 135 | 27 | 0.18 | 0.049 | $-72$ | $-0.14$ |

3 | 135 | 145 | 7.4 | 0.18 | 0.53 | 196 | 1.38 |

4 | 106 | 135 | 27 | 0.18 | 0.53 | 196 | 0.38 |

For examples 3 and 4 of Table 5, impacts of errors in PAR and its attenuation coefficient on the model heat content errors are different from examples 1 and 2. The errors in surface PAR (${\u03f5}_{{Q}_{\mathrm{PAR}}}$) are still positive in examples 3 and 4 (the same as in examples 1 and 2, respectively). However, the errors in ${K}_{\mathrm{PAR}}$ are now also positive in contrast to the negative errors in ${K}_{\mathrm{PAR}}$ of examples 1 and 2. The attenuation coefficient in the model is specified based on Aqua-MODIS mean ${K}_{\mathrm{PAR}}$ value for August at mooring M1. In this case, the error in ${K}_{\mathrm{PAR}}$ is positive and around 196% (the model value of the attenuation coefficient is greater than the observed ${K}_{\mathrm{PAR}}$, i.e., shifted toward more turbid water). Corresponding errors in the heat content $\delta {J}_{{K}_{\mathrm{PAR}}}$ and $\delta {J}_{{Q}_{\mathrm{PAR}}}$ are positive, which results in the summation of these two types of errors and an increase in the total error in the model heat content of the upper 25 m.

The ratios of sensitivities $\frac{{S}_{{K}_{\mathrm{PAR}}}}{{S}_{{Q}_{\mathrm{PAR}}}}$ [Eqs. (17) and (18)] estimated for the five Jerlov’s water types are presented in Table 1. The ratio $\frac{{S}_{{K}_{\mathrm{PAR}}}}{{S}_{{Q}_{\mathrm{PAR}}}}$ decreases as the water type is changing from clear to increasingly turbid water. This indicates that with the decrease in water clarity, the relative error in surface PAR (${\u03f5}_{{Q}_{\mathrm{PAR}}}$) introduces a larger error in the model heat content than the relative error in the attenuation coefficient (${\u03f5}_{{K}_{\mathrm{PAR}}}$) of the same magnitude. For the Jerlov’s type 2 water (notation IA in Table 1, ${K}_{\mathrm{PAR}}=0.049$), the ratio $\frac{{S}_{{K}_{\mathrm{PAR}}}}{{S}_{{Q}_{\mathrm{PAR}}}}=0.51$, which means that the relative error in the surface PAR introduces error twice as large into the model heat content than the same magnitude relative error in the attenuation coefficient ${K}_{\mathrm{PAR}}$. For more turbid water Jerlov’s type 5 (notation “III” in Table 1) (${K}_{\mathrm{PAR}}=0.127$, $\frac{{S}_{{K}_{\mathrm{PAR}}}}{{S}_{{Q}_{\mathrm{PAR}}}}=0.14$), the relative error in surface PAR introduces error seven times as large into the model heat content than the same magnitude relative error in the attenuation coefficient ${K}_{\mathrm{PAR}}$.

The spatial distribution of the HF radar surface currents [Fig. 3(a)] and MODIS-Aqua SSTs [Fig. 3(b)] shows a development of a cold, strong, southward flowing jet along the entrance to the bay during the upwelling event around August 15, 2003. This cold jet separates warmer water masses of cyclonic circulation inside the bay and warmer water masses of the anticyclonic circulation outside the bay [Fig. 3(b)].

Figures 3(c) and 3(d) show sensitivities ${S}_{{Q}_{\mathrm{PAR}}}$ [Eq. (17)] and ${S}_{{K}_{\mathrm{PAR}}}$ [Eq. (18)] for the time of 24 h prior to the verification time ${t}_{f}=2003081500$. These maps present the sensitivities of the model heat content [Eq. (5)] to errors in the surface PAR and ${K}_{\mathrm{PAR}}$ prior to the verification time. Sensitivities maps are shown for Jerlov types IA (${K}_{\mathrm{PAR}}=0.049$). Back in time, Fig. 3 shows the spread of high sensitivities in the following directions: (1) northward, along the path of the southward flowing cold upwelling jet, (2) offshore (cyclonically) along the edge between the cold jet and warm anticyclonic offshore eddy, and (3) into the northern part of the bay (upwelling shadow area). This is in accord with circulation during upwelling favorable conditions of August 15, 2003, described above, when the water masses at mooring M1 are mostly influenced by the advection of the upwelled cold water from the upwelling center to the north of the Monterey Bay, interaction between upwelled water of the cold jet and offshore warm anticyclonic circulation, and by the advection of the warmer water from the northern part of the bay along the entrance to the bay.

The solution of the adjoint problem formulated in Appendix is nonnegative, and in this case, the sensitivity ${S}_{{Q}_{\mathrm{PAR}}}$ [Fig. 3(c)] of the heat content [Eq. (5)] to the error in the surface PAR (${Q}_{\mathrm{PAR}}$) is positive in accord with [Eq. (17)]. In contrast to sensitivities to surface PAR errors [Fig. 3(c)], sensitivities to the errors in ${K}_{\mathrm{PAR}}$ [${S}_{{K}_{\mathrm{PAR}}}$, Figure 3(d)] have positive and negative values. They are positive in shallower water for locations on the shelf in the northern part of the bay and negative in deeper offshore waters to the west of the M1 mooring. Areas of positive values of ${S}_{{Q}_{\mathrm{PAR}}}$ [Eq. (17)] and ${S}_{{K}_{\mathrm{PAR}}}$ [Eq. (18)] [Figs. 3(c) and 3(d)] indicate that when errors in ${Q}_{\mathrm{PAR}}$ (${\u03f5}_{{Q}_{\mathrm{PAR}}}$) and in ${K}_{\mathrm{PAR}}$ (${\u03f5}_{{K}_{\mathrm{PAR}}}$) are both positive (or both negative), their impacts on the model heat content error will add up and will result in an increase of the total model heat content error in comparison to the case when errors ${\u03f5}_{{Q}_{\mathrm{PAR}}}$ and ${\u03f5}_{{K}_{\mathrm{PAR}}}$ have opposite signs.

## 4.

## Discussions and Conclusions

We have presented sensitivity analyses to investigate the impact that errors in model surface PAR and vertical attenuation of PAR have on the upper ocean model heat content. In the Monterey Bay area, we have shown that with a decrease in water clarity, the relative error in surface PAR introduces a larger error in the modeled upper ocean heat content than the same magnitude relative error in the attenuation coefficient. For Jerlov’s type IA water (attenuation coefficient is $0.049\text{\hspace{0.17em}\hspace{0.17em}}{\mathrm{m}}^{-1}$), the relative error in surface PAR introduces error twice as large into the model heat content as the same magnitude relative error in the attenuation coefficient. For the more turbid water Jerlov’s type III (attenuation coefficient is $0.127\text{\hspace{0.17em}\hspace{0.17em}}{\mathrm{m}}^{-1}$), the relative error in surface PAR introduces error seven times as large into the model heat content than the same magnitude relative error in the attenuation coefficient. As a result, it is critical to assimilate the available remotely sensed data products of PAR into the atmospheric model (to improve atmospheric model predictions of PAR used by the oceanic model) or to directly assimilate PAR products into the oceanic model by correcting surface heat flux conditions of the oceanic model. In data assimilation, one merges model predictions with observations based on their respective errors. To achieve this merge, remote sensing and numerical oceanic modeling communities have to better understand the sensitivity of the model heat predictions to the errors in PAR and ${K}_{\mathrm{PAR}}$, which is the paper’s objective. The remote sensing community is familiar with errors in PAR and ${K}_{\mathrm{PAR}}$, however, they are not familiar with how these errors interact in the modeling of the upper ocean heat budget, and how they impact the modeling of the upper ocean heat budget in most existing numerical physical models, and the paper provides this information.

We have presented sensitivities maps of the upper ocean model heat content to the errors in surface PAR and ${K}_{\mathrm{PAR}}$. While the sensitivities to the errors in surface PAR are all positive, sensitivities to the errors in ${K}_{\mathrm{PAR}}$ have positive and negative values. They are positive in shallower water for locations on the shelf in the northern part of the bay and negative in deeper offshore waters to the west of the M1 mooring. Areas of positive values of both sensitivities indicate that when errors in surface PAR and in ${K}_{\mathrm{PAR}}$ are both positive (or both negative), their impacts on the model heat content error will add up and will result in an increase of the total model heat content error in comparison to the case when errors in PAR and ${K}_{\mathrm{PAR}}$ have opposite signs.

Derived equations for sensitivities are valid for any model where the temperature change is described as a tracer with the advection–diffusion-source model. Based on derived sensitivities, the impact and interaction of errors in PAR and ${K}_{\mathrm{PAR}}$ can be estimated for any modeling domain of interest.

The future of oceanic modeling is in the development of the coupled atmosphere-ocean-biochemical modeling systems with the assimilation of physical, bio-optical observations to improve the initialization of coupled systems and their parameters. The assimilation of available satellite and *in situ* physical, bio-optical observations should include the assimilation of PAR, as well as ${K}_{\mathrm{PAR}}$ observations. At present, the fidelity of these coupled modeling systems and existing data assimilation approaches are not sufficient for reliable prediction of surface PAR, as well as, of bio-optical properties that impact the modeling of ${K}_{\mathrm{PAR}}$ in the coupled systems (such as estimation of backscattering and partitioned absorption coefficients). With this in mind, it is important to understand how errors in surface PAR (from the atmospheric component of the coupled system), and ${K}_{\mathrm{PAR}}$ (from the biochemical component of the system) will affect the oceanic model predictions, as for example, the upper ocean model heat content. Sensitivities such as those presented here provide a capability to study and understand the impact of different system components on oceanic model predictions.

## Appendices

## Appendix

Consider the following variable $\lambda $ satisfying the adjoint to Eq. (1):

## (7)

$$-\frac{\partial \lambda}{\partial t}=\mathbf{U}\nabla \lambda +\mathrm{div}(k\nabla \lambda ),$$## (8)

$${\lambda |}_{t={t}_{f}}=1\text{\hspace{0.17em}}\text{in}\text{\hspace{0.17em}}\mathrm{\Omega}\phantom{\rule[-0.0ex]{1em}{0.0ex}}\text{and}\phantom{\rule[-0.0ex]{1em}{0.0ex}}{\lambda |}_{t={t}_{f}}=0\text{\hspace{0.17em}}\text{outside of}\text{\hspace{0.17em}}\mathrm{\Omega}\text{\hspace{0.17em}}\text{in}\text{\hspace{0.17em}}D,$$## (9)

$$\lambda (s,t)=0\phantom{\rule[-0.0ex]{1em}{0.0ex}}\text{and}\phantom{\rule[-0.0ex]{1em}{0.0ex}}\frac{\partial \lambda}{\partial n}(s,t)=0\phantom{\rule[-0.0ex]{1em}{0.0ex}}\text{for}\text{\hspace{0.17em}\hspace{0.17em}}s\in {S}_{d}-\text{boundary of the}\text{\hspace{0.17em}}D\text{\hspace{0.17em}}\text{domain}.$$We multiply Eq. (1) by $\lambda $ and subtract Eq. (7) multiplied by $T$ and integrate over domain $D$ and time $[{t}_{0},{t}_{f}]$, and after integrations by parts we get

## (10)

$${\int}_{D}[T({t}_{f})\lambda ({t}_{f})-{T}_{0}\lambda ({t}_{0})]\mathrm{d}\tau =-{\int}_{D}{\int}_{0}^{{t}_{f}}\lambda \frac{Qsr}{{\rho}_{0}{c}_{p}}\frac{\partial \gamma}{\partial z}\mathrm{d}t\text{\hspace{0.17em}}\mathrm{d}\tau .$$## (11)

$$J={\rho}_{0}{c}_{p}{\int}_{\mathrm{\Omega}}T({t}_{f},\tau )\mathrm{d}\tau ={\rho}_{0}{c}_{p}{\int}_{D}{T}_{0}(\tau )\lambda ({t}_{0},\tau )\mathrm{d}\tau -{\int}_{D}{\int}_{{t}_{0}}^{{t}_{f}}\lambda (t,\tau )Qsr\frac{\partial \gamma}{\partial z}\mathrm{d}t\text{\hspace{0.17em}}\mathrm{d}\tau .$$As stated in Sec. 2.2, the $Qsr$ is a sum of the visible–ultraviolet bands (350 to 700 nm, PAR) and IR bands (700 to 2500 nm), and the vertical distribution of the PAR and IR is described by two exponential functions [see Eqs. (3) and (4)].

By substituting Eqs. (3) and (4) into Eq. (11), we get

## (12)

$$J={\rho}_{0}{c}_{p}{\int}_{D}{T}_{0}(\tau )\lambda ({t}_{0},\tau )\mathrm{d}\tau -{\int}_{{t}_{0}}^{{t}_{f}}{\int}_{S}{\int}_{0}^{H(s)}\lambda (s,z,t)[{Q}_{\mathrm{PAR}}\frac{\partial \text{\hspace{0.17em}}\mathrm{exp}(-{K}_{\mathrm{PAR}}z)}{\partial z}+{Q}_{\mathrm{IR}}\frac{\partial \text{\hspace{0.17em}}\mathrm{exp}(-{K}_{\mathrm{IR}}z)}{\partial z}]\mathrm{d}z\text{\hspace{0.17em}}\mathrm{d}s\text{\hspace{0.17em}}\mathrm{d}t,$$Let us introduce notations $\delta {Q}_{\mathrm{PAR}}(s,t)$ and $\delta {K}_{\mathrm{PAR}}(s,t)$, which are errors in ${Q}_{\mathrm{PAR}}$ and ${K}_{\mathrm{PAR}}$ at time $t$ (${t}_{0}\le t\le {t}_{f}$) and at particular location $s\in S$. In this case, the corresponding errors in the heat content $J$ [Eq. (5)] can be written in the following form:

## (13)

$$\delta {J}_{{Q}_{\mathrm{PAR}}}=\frac{\partial J}{\partial {Q}_{\mathrm{PAR}}(s,t)}\xb7\delta {Q}_{\mathrm{PAR}}(s,t),$$## (14)

$$\delta {J}_{{K}_{\mathrm{PAR}}}=\frac{\partial J}{\partial {K}_{\mathrm{PAR}}(s,t)}\xb7\delta {K}_{\mathrm{PAR}}(s,t),$$Based on Eq. (12), and after simple transformations, we have

where## (17)

$${S}_{{Q}_{\mathrm{PAR}}}={Q}_{\mathrm{PAR}}(s,t){K}_{\mathrm{PAR}}(s,t){\int}_{0}^{H(s)}\lambda (s,z,t){e}^{-Kpar(s,t)z}\mathrm{d}\mathrm{z},$$## (18)

$${S}_{{K}_{\mathrm{PAR}}}={Q}_{\mathrm{PAR}}(s,t){K}_{\mathrm{PAR}}(s,t){\int}_{0}^{H(s)}[1-{K}_{\mathrm{PAR}}(s,t)z]\lambda (s,z,t){e}^{-Kpar(s,t)z}\mathrm{d}z,$$## (19)

$${\u03f5}_{{Q}_{\mathrm{PAR}}}=\frac{\delta {Q}_{\mathrm{PAR}}(s,t)}{{Q}_{\mathrm{PAR}}(s,t)},$$## (20)

$${\u03f5}_{{K}_{\mathrm{PAR}}}=\frac{\delta {K}_{\mathrm{PAR}}(s,t)}{{K}_{\mathrm{PAR}}(s,t)},$$In accord with Eqs. (13)–(18), units of errors $\delta {J}_{{Q}_{\mathrm{PAR}}}$ and $\delta {J}_{{K}_{\mathrm{PAR}}}$, as well as sensitivities ${S}_{{Q}_{\mathrm{PAR}}}$ and ${S}_{{K}_{\mathrm{PAR}}}$ are in joules.

## Acknowledgments

This research was funded through the NRL under program element 61153N. Computer time for the numerical simulations was provided through a grant from the Department of Defense High Performance Computing Initiative. Request for access to the data presented in this paper can be sent to igor.shulman@nrlssc.navy.mil. This paper is an NRL contribution: 7330-15-2894.

## References

## Biography

**Igor Shulman** received his PhD in mechanics of fluid, gas, and plasma from the Russian University of Oil and Gas Technology, Russia. He has been an oceanographer with the Naval Research Laboratory since 2003. His research interests include data assimilation and modeling, development of coupled physical and bio-optical models, and bioluminescence. He has been PI on numerous interdisciplinary projects funded by ONR, NOPP, and MURI. He has published more than 70 articles in refereed journals.

**Richard W. Gould Jr.** received his PhD in oceanography from Texas A&M University and has 30 years of oceanographic and remote sensing experience and 54 publications. He has developed multispectral/hyperspectral ocean color algorithms, and has transitioned satellite products to the Navy and NASA. Research interests include satellite algorithm development, uncertainty analyses, optical water mass classification, coastal hypoxia, and physical/bio-optical coupling. He is a head of the bio-optical/physical processes and remote sensing section at the Naval Research Laboratory.

**Stephanie Anderson** received her BS degree in resource development from the University of Rhode Island. She has worked at the Naval Research Laboratory for over 20 years and before that with the Environmental Protection Agency's Narragansett Marine Research Laboratory for eight years. Her current research interests are the analysis and interpretation of oceanographic remote sensing and *in situ* observations and ocean model-observation comparisons and validation. She has published around 16 articles in refereed journals.

**Peter Sakalaukus** received his MS degree in computer science from the University of Southern Mississippi (USM) in 2005. He has been with NRL for six years, and before with USM for 19 years. His research interests include physical, bio-optical modeling, data analysis, and interpretation. He has published 25 articles in refereed journals.