\hypertarget{namespacerossby__front__2d__initialization}{}\section{rossby\+\_\+front\+\_\+2d\+\_\+initialization Module Reference}
\label{namespacerossby__front__2d__initialization}\index{rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}}


\subsection{Detailed Description}
Initial conditions for the 2D Rossby front test. 

\hypertarget{namespacerossby__front__2d__initialization_section_Rossby_front_2d}{}\subsection{Description of the 2d Rossby front initial conditions}\label{namespacerossby__front__2d__initialization_section_Rossby_front_2d}
Consistent with a linear equation of state, the system has a constant stratification below the mixed layer, stratified in temperature only. Isotherms are flat below the mixed layer and vertical within. Salinity is constant. The mixed layer has a half sine form so that there are no mixed layer or temperature gradients at the side walls.

Below the mixed layer the potential temperature, $\theta(z)$, is given by \[ \theta(z) = \theta_0 - \Delta \theta \left( z + h_{ML} \right) \] where $ \theta_0 $ and $ \Delta \theta $ are external model parameters.

The depth of the mixed layer, $H_{ML}$ is \[ h_{ML}(y) = h_{min} + \left( h_{max} - h_{min} \right) \cos{\pi y/L} \]. The temperature in mixed layer is given by the reference temperature at $z=h_{ML}$ so that \begin{eqnarray} \theta(y,z) = \theta_0 - \Delta \theta \left( z + h_{ML} \right) & \forall & z < h_{ML}(y) T.B.D. \end{eqnarray} \subsection*{Functions/\+Subroutines}
\begin{DoxyCompactItemize}
\item 
subroutine, public \mbox{\hyperlink{namespacerossby__front__2d__initialization_a1cea6a589ecdbb72113e0922aebcadad}{rossby\+\_\+front\+\_\+initialize\+\_\+thickness}} (h, G, GV, US, param\+\_\+file, just\+\_\+read\+\_\+params)
\begin{DoxyCompactList}\small\item\em Initialization of thicknesses in 2D Rossby front test. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacerossby__front__2d__initialization_a97433b821b82240bbd7237076c06f667}{rossby\+\_\+front\+\_\+initialize\+\_\+temperature\+\_\+salinity}} (T, S, h, G, GV, param\+\_\+file, eqn\+\_\+of\+\_\+state, just\+\_\+read\+\_\+params)
\begin{DoxyCompactList}\small\item\em Initialization of temperature and salinity in the Rossby front test. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacerossby__front__2d__initialization_a4cdf99efb62134cf4ee9b3dac0b72205}{rossby\+\_\+front\+\_\+initialize\+\_\+velocity}} (u, v, h, G, GV, US, param\+\_\+file, just\+\_\+read\+\_\+params)
\begin{DoxyCompactList}\small\item\em Initialization of u and v in the Rossby front test. \end{DoxyCompactList}\item 
real function \mbox{\hyperlink{namespacerossby__front__2d__initialization_a15a0b752df24fbae7deabe844a418239}{ypseudo}} (G, lat)
\begin{DoxyCompactList}\small\item\em Pseudo coordinate across domain used by Hml() and d\+Tdy() returns a coordinate from -\/\+P\+I/2 .. P\+I/2 squashed towards the center of the domain. \end{DoxyCompactList}\item 
real function \mbox{\hyperlink{namespacerossby__front__2d__initialization_aa10adb0378184432ecaa78eb339c6c5a}{hml}} (G, lat)
\begin{DoxyCompactList}\small\item\em Analytic prescription of mixed layer depth in 2d Rossby front test, in the same units as Gmax\+\_\+depth. \end{DoxyCompactList}\item 
real function \mbox{\hyperlink{namespacerossby__front__2d__initialization_a587a5f5c3f4694558d3d5206840ccab2}{dtdy}} (G, dT, lat)
\begin{DoxyCompactList}\small\item\em Analytic prescription of mixed layer temperature gradient in 2d Rossby front test. \end{DoxyCompactList}\end{DoxyCompactItemize}
\subsection*{Variables}
\begin{DoxyCompactItemize}
\item 
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_a09788e02ae4a4b4904ee206ddaa3a2a8}\label{namespacerossby__front__2d__initialization_a09788e02ae4a4b4904ee206ddaa3a2a8}} 
character(len=40) \mbox{\hyperlink{namespacerossby__front__2d__initialization_a09788e02ae4a4b4904ee206ddaa3a2a8}{mdl}} = \char`\"{}Rossby\+\_\+front\+\_\+2d\+\_\+initialization\char`\"{}
\begin{DoxyCompactList}\small\item\em This module\textquotesingle{}s name. \end{DoxyCompactList}\item 
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_a8e0cf37e24942144e2ba69f852d0fad4}\label{namespacerossby__front__2d__initialization_a8e0cf37e24942144e2ba69f852d0fad4}} 
real, parameter \mbox{\hyperlink{namespacerossby__front__2d__initialization_a8e0cf37e24942144e2ba69f852d0fad4}{frontfractionalwidth}} = 0.\+5
\begin{DoxyCompactList}\small\item\em Width of front as fraction of domain \mbox{[}nondim\mbox{]}. \end{DoxyCompactList}\item 
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_a75da0dacb75190f3cf8911be0885c9c7}\label{namespacerossby__front__2d__initialization_a75da0dacb75190f3cf8911be0885c9c7}} 
real, parameter \mbox{\hyperlink{namespacerossby__front__2d__initialization_a75da0dacb75190f3cf8911be0885c9c7}{hmlmin}} = 0.\+25
\begin{DoxyCompactList}\small\item\em Shallowest ML as fractional depth of ocean \mbox{[}nondim\mbox{]}. \end{DoxyCompactList}\item 
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_abaf46ae8dfac538478b49d25ddba4ffb}\label{namespacerossby__front__2d__initialization_abaf46ae8dfac538478b49d25ddba4ffb}} 
real, parameter \mbox{\hyperlink{namespacerossby__front__2d__initialization_abaf46ae8dfac538478b49d25ddba4ffb}{hmlmax}} = 0.\+75
\begin{DoxyCompactList}\small\item\em Deepest ML as fractional depth of ocean \mbox{[}nondim\mbox{]}. \end{DoxyCompactList}\end{DoxyCompactItemize}


\subsection{Function/\+Subroutine Documentation}
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_a587a5f5c3f4694558d3d5206840ccab2}\label{namespacerossby__front__2d__initialization_a587a5f5c3f4694558d3d5206840ccab2}} 
\index{rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}!dtdy@{dtdy}}
\index{dtdy@{dtdy}!rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}}
\subsubsection{\texorpdfstring{dtdy()}{dtdy()}}
{\footnotesize\ttfamily real function rossby\+\_\+front\+\_\+2d\+\_\+initialization\+::dtdy (\begin{DoxyParamCaption}\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{real, intent(in)}]{dT,  }\item[{real, intent(in)}]{lat }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



Analytic prescription of mixed layer temperature gradient in 2d Rossby front test. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & Grid structure\\
\hline
\mbox{\tt in}  & {\em dt} & Top to bottom temperature difference\\
\hline
\mbox{\tt in}  & {\em lat} & Latitude \\
\hline
\end{DoxyParams}


Definition at line 256 of file Rossby\+\_\+front\+\_\+2d\+\_\+initialization.\+F90.


\begin{DoxyCode}
256   \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: G\textcolor{comment}{     !< Grid structure}
257   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)} :: dT\textcolor{comment}{    !< Top to bottom temperature difference}
258   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)} :: lat\textcolor{comment}{   !< Latitude}
259   \textcolor{comment}{! Local}
260   \textcolor{keywordtype}{real} :: PI, dHML, dHdy
261   \textcolor{keywordtype}{real} :: km = 1.e3 \textcolor{comment}{! AXIS\_UNITS = 'k' (1000 m)}
262 
263   pi = 4.0 * atan(1.0)
264   dhml = 0.5 * ( hmlmax - hmlmin ) * g%max\_depth
265   dhdy = dhml * ( pi / ( frontfractionalwidth * g%len\_lat * km ) ) * cos( ypseudo(g, lat) )
266   dtdy = -( dt / g%max\_depth ) * dhdy
267 
\end{DoxyCode}
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_aa10adb0378184432ecaa78eb339c6c5a}\label{namespacerossby__front__2d__initialization_aa10adb0378184432ecaa78eb339c6c5a}} 
\index{rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}!hml@{hml}}
\index{hml@{hml}!rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}}
\subsubsection{\texorpdfstring{hml()}{hml()}}
{\footnotesize\ttfamily real function rossby\+\_\+front\+\_\+2d\+\_\+initialization\+::hml (\begin{DoxyParamCaption}\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{real, intent(in)}]{lat }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



Analytic prescription of mixed layer depth in 2d Rossby front test, in the same units as Gmax\+\_\+depth. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & Grid structure\\
\hline
\mbox{\tt in}  & {\em lat} & Latitude \\
\hline
\end{DoxyParams}


Definition at line 243 of file Rossby\+\_\+front\+\_\+2d\+\_\+initialization.\+F90.


\begin{DoxyCode}
243   \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: G\textcolor{comment}{   !< Grid structure}
244   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)} :: lat\textcolor{comment}{ !< Latitude}
245   \textcolor{comment}{! Local}
246   \textcolor{keywordtype}{real} :: dHML, HMLmean
247 
248   dhml = 0.5 * ( hmlmax - hmlmin ) * g%max\_depth
249   hmlmean = 0.5 * ( hmlmin + hmlmax ) * g%max\_depth
250   hml = hmlmean + dhml * sin( ypseudo(g, lat) )
\end{DoxyCode}
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_a97433b821b82240bbd7237076c06f667}\label{namespacerossby__front__2d__initialization_a97433b821b82240bbd7237076c06f667}} 
\index{rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}!rossby\+\_\+front\+\_\+initialize\+\_\+temperature\+\_\+salinity@{rossby\+\_\+front\+\_\+initialize\+\_\+temperature\+\_\+salinity}}
\index{rossby\+\_\+front\+\_\+initialize\+\_\+temperature\+\_\+salinity@{rossby\+\_\+front\+\_\+initialize\+\_\+temperature\+\_\+salinity}!rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}}
\subsubsection{\texorpdfstring{rossby\+\_\+front\+\_\+initialize\+\_\+temperature\+\_\+salinity()}{rossby\_front\_initialize\_temperature\_salinity()}}
{\footnotesize\ttfamily subroutine, public rossby\+\_\+front\+\_\+2d\+\_\+initialization\+::rossby\+\_\+front\+\_\+initialize\+\_\+temperature\+\_\+salinity (\begin{DoxyParamCaption}\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed,  g \%ke), intent(out)}]{T,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed,  g \%ke), intent(out)}]{S,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed,  g \%ke), intent(in)}]{h,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(param\+\_\+file\+\_\+type), intent(in)}]{param\+\_\+file,  }\item[{type(eos\+\_\+type), pointer}]{eqn\+\_\+of\+\_\+state,  }\item[{logical, intent(in), optional}]{just\+\_\+read\+\_\+params }\end{DoxyParamCaption})}



Initialization of temperature and salinity in the Rossby front test. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & Grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\
\hline
\mbox{\tt out}  & {\em t} & Potential temperature \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt out}  & {\em s} & Salinity \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h} & Thickness \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em param\+\_\+file} & Parameter file handle\\
\hline
 & {\em eqn\+\_\+of\+\_\+state} & Equation of state structure\\
\hline
\mbox{\tt in}  & {\em just\+\_\+read\+\_\+params} & If present and true, this call will only read parameters without changing T \& S. \\
\hline
\end{DoxyParams}


Definition at line 114 of file Rossby\+\_\+front\+\_\+2d\+\_\+initialization.\+F90.


\begin{DoxyCode}
114   \textcolor{keywordtype}{type}(ocean\_grid\_type),                     \textcolor{keywordtype}{intent(in)}  :: G\textcolor{comment}{  !< Grid structure}
115   \textcolor{keywordtype}{type}(verticalGrid\_type),                   \textcolor{keywordtype}{intent(in)}  :: GV\textcolor{comment}{ !< The ocean's vertical grid structure.}
116   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G), SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: T\textcolor{comment}{  !< Potential temperature [degC]}
117   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G), SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: S\textcolor{comment}{  !< Salinity [ppt]}
118   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G), SZK\_(G))}, \textcolor{keywordtype}{intent(in)}  :: h\textcolor{comment}{  !< Thickness [H ~> m or kg m-2]}
119   \textcolor{keywordtype}{type}(param\_file\_type),                     \textcolor{keywordtype}{intent(in)}  :: param\_file\textcolor{comment}{   !< Parameter file handle}
120   \textcolor{keywordtype}{type}(EOS\_type),                            \textcolor{keywordtype}{pointer}     :: eqn\_of\_state\textcolor{comment}{ !< Equation of state structure}
121   \textcolor{keywordtype}{logical},       \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: just\_read\_params\textcolor{comment}{ !< If present and true, this call will}
122 \textcolor{comment}{                                                      !! only read parameters without changing T & S.}
123 
124   \textcolor{keywordtype}{integer}   :: i, j, k, is, ie, js, je, nz
125   \textcolor{keywordtype}{real}      :: T\_ref, S\_ref \textcolor{comment}{! Reference salinity and temerature within surface layer}
126   \textcolor{keywordtype}{real}      :: T\_range      \textcolor{comment}{! Range of salinities and temperatures over the vertical}
127   \textcolor{keywordtype}{real}      :: y, zc, zi, dTdz
128   \textcolor{keywordtype}{logical} :: just\_read    \textcolor{comment}{! If true, just read parameters but set nothing.}
129   \textcolor{keywordtype}{character(len=40)} :: verticalCoordinate
130   \textcolor{keywordtype}{real}      :: PI                   \textcolor{comment}{! 3.1415926... calculated as 4*atan(1)}
131 
132   is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
133 
134   just\_read = .false. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(just\_read\_params)) just\_read = just\_read\_params
135 
136   \textcolor{keyword}{call }get\_param(param\_file, mdl,\textcolor{stringliteral}{"REGRIDDING\_COORDINATE\_MODE"}, verticalcoordinate, &
137             default=default\_coordinate\_mode, do\_not\_log=just\_read)
138   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"S\_REF"}, s\_ref, \textcolor{stringliteral}{'Reference salinity'}, &
139                  default=35.0, units=\textcolor{stringliteral}{'1e-3'}, do\_not\_log=just\_read)
140   \textcolor{keyword}{call }get\_param(param\_file, mdl,\textcolor{stringliteral}{"T\_REF"},t\_ref,\textcolor{stringliteral}{'Reference temperature'},units=\textcolor{stringliteral}{'C'},&
141                  fail\_if\_missing=.not.just\_read, do\_not\_log=just\_read)
142   \textcolor{keyword}{call }get\_param(param\_file, mdl,\textcolor{stringliteral}{"T\_RANGE"},t\_range,\textcolor{stringliteral}{'Initial temperature range'},&
143                  units=\textcolor{stringliteral}{'C'}, default=0.0, do\_not\_log=just\_read)
144 
145   \textcolor{keywordflow}{if} (just\_read) \textcolor{keywordflow}{return} \textcolor{comment}{! All run-time parameters have been read, so return.}
146 
147   t(:,:,:) = 0.0
148   s(:,:,:) = s\_ref
149   dtdz = t\_range / g%max\_depth
150 
151   \textcolor{keywordflow}{do} j = g%jsc,g%jec ; \textcolor{keywordflow}{do} i = g%isc,g%iec
152     zi = 0.
153     \textcolor{keywordflow}{do} k = 1, nz
154       zi = zi - h(i,j,k)              \textcolor{comment}{! Bottom interface position}
155       zc = gv%H\_to\_Z * (zi - 0.5*h(i,j,k))    \textcolor{comment}{! Position of middle of cell}
156       zc = min( zc, -hml(g, g%geoLatT(i,j)) ) \textcolor{comment}{! Bound by depth of mixed layer}
157       t(i,j,k) = t\_ref + dtdz * zc \textcolor{comment}{! Linear temperature profile}
158 \textcolor{keywordflow}{    enddo}
159 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
160 
\end{DoxyCode}
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_a1cea6a589ecdbb72113e0922aebcadad}\label{namespacerossby__front__2d__initialization_a1cea6a589ecdbb72113e0922aebcadad}} 
\index{rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}!rossby\+\_\+front\+\_\+initialize\+\_\+thickness@{rossby\+\_\+front\+\_\+initialize\+\_\+thickness}}
\index{rossby\+\_\+front\+\_\+initialize\+\_\+thickness@{rossby\+\_\+front\+\_\+initialize\+\_\+thickness}!rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}}
\subsubsection{\texorpdfstring{rossby\+\_\+front\+\_\+initialize\+\_\+thickness()}{rossby\_front\_initialize\_thickness()}}
{\footnotesize\ttfamily subroutine, public rossby\+\_\+front\+\_\+2d\+\_\+initialization\+::rossby\+\_\+front\+\_\+initialize\+\_\+thickness (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)), intent(out)}]{h,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(param\+\_\+file\+\_\+type), intent(in)}]{param\+\_\+file,  }\item[{logical, intent(in), optional}]{just\+\_\+read\+\_\+params }\end{DoxyParamCaption})}



Initialization of thicknesses in 2D Rossby front test. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & Grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & Vertical grid structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt out}  & {\em h} & The thickness that is being initialized \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em param\+\_\+file} & A structure indicating the open file to parse for model parameter values.\\
\hline
\mbox{\tt in}  & {\em just\+\_\+read\+\_\+params} & If present and true, this call will only read parameters without changing h. \\
\hline
\end{DoxyParams}


Definition at line 40 of file Rossby\+\_\+front\+\_\+2d\+\_\+initialization.\+F90.


\begin{DoxyCode}
40   \textcolor{keywordtype}{type}(ocean\_grid\_type),   \textcolor{keywordtype}{intent(in)}  :: G\textcolor{comment}{           !< Grid structure}
41   \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)}  :: GV\textcolor{comment}{          !< Vertical grid structure}
42   \textcolor{keywordtype}{type}(unit\_scale\_type),   \textcolor{keywordtype}{intent(in)}  :: US\textcolor{comment}{          !< A dimensional unit scaling type}
43   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV))}, &
44                            \textcolor{keywordtype}{intent(out)} :: h\textcolor{comment}{           !< The thickness that is being initialized [H ~> m or
       kg m-2]}
45   \textcolor{keywordtype}{type}(param\_file\_type),   \textcolor{keywordtype}{intent(in)}  :: param\_file\textcolor{comment}{  !< A structure indicating the open file}
46 \textcolor{comment}{                                                      !! to parse for model parameter values.}
47   \textcolor{keywordtype}{logical},       \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: just\_read\_params\textcolor{comment}{ !< If present and true, this call will}
48 \textcolor{comment}{                                                      !! only read parameters without changing h.}
49 
50   \textcolor{keywordtype}{integer} :: i, j, k, is, ie, js, je, nz
51   \textcolor{keywordtype}{real}    :: Tz, Dml, eta, stretch, h0
52   \textcolor{keywordtype}{real}    :: min\_thickness, T\_range
53   \textcolor{keywordtype}{real}    :: dRho\_dT      \textcolor{comment}{! The partial derivative of density with temperature [R degC-1 ~> kg m-3 degC-1]}
54   \textcolor{keywordtype}{logical} :: just\_read    \textcolor{comment}{! If true, just read parameters but set nothing.}
55   \textcolor{keywordtype}{character(len=40)} :: verticalCoordinate
56 
57   is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
58 
59   just\_read = .false. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(just\_read\_params)) just\_read = just\_read\_params
60 
61   \textcolor{keywordflow}{if} (.not.just\_read) &
62     \textcolor{keyword}{call }mom\_mesg(\textcolor{stringliteral}{"Rossby\_front\_2d\_initialization.F90, Rossby\_front\_initialize\_thickness: setting
       thickness"})
63 
64   \textcolor{keywordflow}{if} (.not.just\_read) \textcolor{keyword}{call }log\_version(param\_file, mdl, version, \textcolor{stringliteral}{""})
65   \textcolor{comment}{! Read parameters needed to set thickness}
66   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"MIN\_THICKNESS"}, min\_thickness, &
67                  \textcolor{stringliteral}{'Minimum layer thickness'},units=\textcolor{stringliteral}{'m'},default=1.e-3, do\_not\_log=just\_read)
68   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"REGRIDDING\_COORDINATE\_MODE"}, verticalcoordinate, &
69                  default=default\_coordinate\_mode, do\_not\_log=just\_read)
70   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"T\_RANGE"}, t\_range, \textcolor{stringliteral}{'Initial temperature range'}, &
71                  units=\textcolor{stringliteral}{'C'}, default=0.0, do\_not\_log=just\_read)
72   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DRHO\_DT"}, drho\_dt, default=-0.2, scale=us%kg\_m3\_to\_R, do\_not\_log=.true.)
73 
74   \textcolor{keywordflow}{if} (just\_read) \textcolor{keywordflow}{return} \textcolor{comment}{! All run-time parameters have been read, so return.}
75 
76   tz = t\_range / g%max\_depth
77 
78   \textcolor{keywordflow}{select case} ( coordinatemode(verticalcoordinate) )
79 
80     \textcolor{keywordflow}{case} (regridding\_layer, regridding\_rho)
81       \textcolor{keywordflow}{do} j = g%jsc,g%jec ; \textcolor{keywordflow}{do} i = g%isc,g%iec
82         dml = hml( g, g%geoLatT(i,j) )
83         eta = -( -drho\_dt / gv%Rho0 ) * tz * 0.5 * ( dml * dml )
84         stretch = ( ( g%max\_depth + eta ) / g%max\_depth )
85         h0 = ( g%max\_depth / \textcolor{keywordtype}{real(nz)} ) * stretch
86         \textcolor{keywordflow}{do} k = 1, nz
87           h(i,j,k) = h0 * gv%Z\_to\_H
88 \textcolor{keywordflow}{        enddo}
89 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo}
90 
91     \textcolor{keywordflow}{case} (regridding\_zstar, regridding\_sigma)
92       \textcolor{keywordflow}{do} j = g%jsc,g%jec ; \textcolor{keywordflow}{do} i = g%isc,g%iec
93         dml = hml( g, g%geoLatT(i,j) )
94         eta = -( -drho\_dt / gv%Rho0 ) * tz * 0.5 * ( dml * dml )
95         stretch = ( ( g%max\_depth + eta ) / g%max\_depth )
96         h0 = ( g%max\_depth / \textcolor{keywordtype}{real(nz)} ) * stretch
97         \textcolor{keywordflow}{do} k = 1, nz
98           h(i,j,k) = h0 * gv%Z\_to\_H
99 \textcolor{keywordflow}{        enddo}
100 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo}
101 
102 \textcolor{keywordflow}{    case default}
103       \textcolor{keyword}{call }mom\_error(fatal,\textcolor{stringliteral}{"Rossby\_front\_initialize: "}// &
104       \textcolor{stringliteral}{"Unrecognized i.c. setup - set REGRIDDING\_COORDINATE\_MODE"})
105 
106 \textcolor{keywordflow}{  end select}
107 
\end{DoxyCode}
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_a4cdf99efb62134cf4ee9b3dac0b72205}\label{namespacerossby__front__2d__initialization_a4cdf99efb62134cf4ee9b3dac0b72205}} 
\index{rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}!rossby\+\_\+front\+\_\+initialize\+\_\+velocity@{rossby\+\_\+front\+\_\+initialize\+\_\+velocity}}
\index{rossby\+\_\+front\+\_\+initialize\+\_\+velocity@{rossby\+\_\+front\+\_\+initialize\+\_\+velocity}!rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}}
\subsubsection{\texorpdfstring{rossby\+\_\+front\+\_\+initialize\+\_\+velocity()}{rossby\_front\_initialize\_velocity()}}
{\footnotesize\ttfamily subroutine, public rossby\+\_\+front\+\_\+2d\+\_\+initialization\+::rossby\+\_\+front\+\_\+initialize\+\_\+velocity (\begin{DoxyParamCaption}\item[{real, dimension( g \%isdb\+: g \%iedb, g \%jsd\+: g \%jed, g \%ke), intent(out)}]{u,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsdb\+: g \%jedb, g \%ke), intent(out)}]{v,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed,  g \%ke), intent(in)}]{h,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(param\+\_\+file\+\_\+type), intent(in)}]{param\+\_\+file,  }\item[{logical, intent(in), optional}]{just\+\_\+read\+\_\+params }\end{DoxyParamCaption})}



Initialization of u and v in the Rossby front test. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & Grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & Vertical grid structure\\
\hline
\mbox{\tt out}  & {\em u} & i-\/component of velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}\\
\hline
\mbox{\tt out}  & {\em v} & j-\/component of velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h} & Thickness \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em param\+\_\+file} & A structure indicating the open file to parse for model parameter values.\\
\hline
\mbox{\tt in}  & {\em just\+\_\+read\+\_\+params} & If present and true, this call will only read parameters without setting u \& v. \\
\hline
\end{DoxyParams}


Definition at line 166 of file Rossby\+\_\+front\+\_\+2d\+\_\+initialization.\+F90.


\begin{DoxyCode}
166   \textcolor{keywordtype}{type}(ocean\_grid\_type),      \textcolor{keywordtype}{intent(in)}  :: G\textcolor{comment}{  !< Grid structure}
167   \textcolor{keywordtype}{type}(verticalGrid\_type),    \textcolor{keywordtype}{intent(in)}  :: GV\textcolor{comment}{ !< Vertical grid structure}
168   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, &
169                               \textcolor{keywordtype}{intent(out)} :: u\textcolor{comment}{  !< i-component of velocity [L T-1 ~> m s-1]}
170   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, &
171                               \textcolor{keywordtype}{intent(out)} :: v\textcolor{comment}{  !< j-component of velocity [L T-1 ~> m s-1]}
172   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G), SZK\_(G))}, &
173                               \textcolor{keywordtype}{intent(in)}  :: h\textcolor{comment}{  !< Thickness [H ~> m or kg m-2]}
174   \textcolor{keywordtype}{type}(unit\_scale\_type),      \textcolor{keywordtype}{intent(in)}  :: US\textcolor{comment}{ !< A dimensional unit scaling type}
175   \textcolor{keywordtype}{type}(param\_file\_type),      \textcolor{keywordtype}{intent(in)}  :: param\_file\textcolor{comment}{ !< A structure indicating the open file}
176 \textcolor{comment}{                                                !! to parse for model parameter values.}
177   \textcolor{keywordtype}{logical},          \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: just\_read\_params\textcolor{comment}{ !< If present and true, this call}
178 \textcolor{comment}{                                                !! will only read parameters without setting u & v.}
179 
180   \textcolor{keywordtype}{real}    :: y            \textcolor{comment}{! Non-dimensional coordinate across channel, 0..pi}
181   \textcolor{keywordtype}{real}    :: T\_range      \textcolor{comment}{! Range of salinities and temperatures over the vertical}
182   \textcolor{keywordtype}{real}    :: dUdT         \textcolor{comment}{! Factor to convert dT/dy into dU/dz, g*alpha/f [L2 Z-1 T-1 degC-1 ~> m s-1
       degC-1]}
183   \textcolor{keywordtype}{real}    :: dRho\_dT      \textcolor{comment}{! The partial derivative of density with temperature [R degC-1 ~> kg m-3 degC-1]}
184   \textcolor{keywordtype}{real}    :: Dml, zi, zc, zm \textcolor{comment}{! Depths [Z ~> m].}
185   \textcolor{keywordtype}{real}    :: f            \textcolor{comment}{! The local Coriolis parameter [T-1 ~> s-1]}
186   \textcolor{keywordtype}{real}    :: Ty           \textcolor{comment}{! The meridional temperature gradient [degC L-1 ~> degC m-1]}
187   \textcolor{keywordtype}{real}    :: hAtU         \textcolor{comment}{! Interpolated layer thickness [Z ~> m].}
188   \textcolor{keywordtype}{integer} :: i, j, k, is, ie, js, je, nz
189   \textcolor{keywordtype}{logical} :: just\_read    \textcolor{comment}{! If true, just read parameters but set nothing.}
190   \textcolor{keywordtype}{character(len=40)} :: verticalCoordinate
191 
192   is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
193 
194   just\_read = .false. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(just\_read\_params)) just\_read = just\_read\_params
195 
196   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"REGRIDDING\_COORDINATE\_MODE"}, verticalcoordinate, &
197                  default=default\_coordinate\_mode, do\_not\_log=just\_read)
198   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"T\_RANGE"}, t\_range, \textcolor{stringliteral}{'Initial temperature range'}, &
199                  units=\textcolor{stringliteral}{'C'}, default=0.0, do\_not\_log=just\_read)
200   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DRHO\_DT"}, drho\_dt, default=-0.2, scale=us%kg\_m3\_to\_R, do\_not\_log=.true.)
201 
202   \textcolor{keywordflow}{if} (just\_read) \textcolor{keywordflow}{return} \textcolor{comment}{! All run-time parameters have been read, so return.}
203 
204   v(:,:,:) = 0.0
205   u(:,:,:) = 0.0
206 
207   \textcolor{keywordflow}{do} j = g%jsc,g%jec ; \textcolor{keywordflow}{do} i = g%isc-1,g%iec+1
208     f = 0.5* (g%CoriolisBu(i,j) + g%CoriolisBu(i,j-1) )
209     dudt = 0.0 ; \textcolor{keywordflow}{if} (abs(f) > 0.0) &
210       dudt = ( gv%g\_Earth*drho\_dt ) / ( f * gv%Rho0 )
211     dml = hml( g, g%geoLatT(i,j) )
212     ty = us%L\_to\_m*dtdy( g, t\_range, g%geoLatT(i,j) )
213     zi = 0.
214     \textcolor{keywordflow}{do} k = 1, nz
215       hatu = 0.5*(h(i,j,k)+h(i+1,j,k)) * gv%H\_to\_Z
216       zi = zi - hatu              \textcolor{comment}{! Bottom interface position}
217       zc = zi - 0.5*hatu          \textcolor{comment}{! Position of middle of cell}
218       zm = max( zc + dml, 0. )    \textcolor{comment}{! Height above bottom of mixed layer}
219       u(i,j,k) = dudt * ty * zm   \textcolor{comment}{! Thermal wind starting at base of ML}
220 \textcolor{keywordflow}{    enddo}
221 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
222 
\end{DoxyCode}
\mbox{\Hypertarget{namespacerossby__front__2d__initialization_a15a0b752df24fbae7deabe844a418239}\label{namespacerossby__front__2d__initialization_a15a0b752df24fbae7deabe844a418239}} 
\index{rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}!ypseudo@{ypseudo}}
\index{ypseudo@{ypseudo}!rossby\+\_\+front\+\_\+2d\+\_\+initialization@{rossby\+\_\+front\+\_\+2d\+\_\+initialization}}
\subsubsection{\texorpdfstring{ypseudo()}{ypseudo()}}
{\footnotesize\ttfamily real function rossby\+\_\+front\+\_\+2d\+\_\+initialization\+::ypseudo (\begin{DoxyParamCaption}\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{real, intent(in)}]{lat }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



Pseudo coordinate across domain used by Hml() and d\+Tdy() returns a coordinate from -\/\+P\+I/2 .. P\+I/2 squashed towards the center of the domain. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & Grid structure\\
\hline
\mbox{\tt in}  & {\em lat} & Latitude \\
\hline
\end{DoxyParams}


Definition at line 229 of file Rossby\+\_\+front\+\_\+2d\+\_\+initialization.\+F90.


\begin{DoxyCode}
229   \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: G\textcolor{comment}{   !< Grid structure}
230   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)} :: lat\textcolor{comment}{ !< Latitude}
231   \textcolor{comment}{! Local}
232   \textcolor{keywordtype}{real} :: y, PI
233 
234   pi = 4.0 * atan(1.0)
235   ypseudo = ( ( lat - g%south\_lat ) / g%len\_lat ) - 0.5 \textcolor{comment}{! -1/2 .. 1/.2}
236   ypseudo = pi * max(-0.5, min(0.5, ypseudo / frontfractionalwidth))
\end{DoxyCode}
