\hypertarget{namespacecoord__adapt}{}\section{coord\+\_\+adapt Module Reference}
\label{namespacecoord__adapt}\index{coord\+\_\+adapt@{coord\+\_\+adapt}}


\subsection{Detailed Description}
Regrid columns for the adaptive coordinate. \subsection*{Data Types}
\begin{DoxyCompactItemize}
\item 
type \mbox{\hyperlink{structcoord__adapt_1_1adapt__cs}{adapt\+\_\+cs}}
\begin{DoxyCompactList}\small\item\em Control structure for adaptive coordinates (\mbox{\hyperlink{namespacecoord__adapt}{coord\+\_\+adapt}}). \end{DoxyCompactList}\end{DoxyCompactItemize}
\subsection*{Functions/\+Subroutines}
\begin{DoxyCompactItemize}
\item 
subroutine, public \mbox{\hyperlink{namespacecoord__adapt_a53eb167f67898f13a7f2258365f1b968}{init\+\_\+coord\+\_\+adapt}} (CS, nk, coordinate\+Resolution, m\+\_\+to\+\_\+H, kg\+\_\+m3\+\_\+to\+\_\+R)
\begin{DoxyCompactList}\small\item\em Initialise an adapt\+\_\+\+CS with parameters. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacecoord__adapt_a92e88dea18cc4b3adef844f58a8fc0ff}{end\+\_\+coord\+\_\+adapt}} (CS)
\begin{DoxyCompactList}\small\item\em Clean up the coordinate control structure. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacecoord__adapt_a35706e0359c3953aac56b160e48b0ef8}{set\+\_\+adapt\+\_\+params}} (CS, adapt\+Time\+Ratio, adapt\+Alpha, adapt\+Zoom, adapt\+Zoom\+Coeff, adapt\+Buoy\+Coeff, adapt\+Drho0, adapt\+Do\+Min)
\begin{DoxyCompactList}\small\item\em This subtroutine can be used to set the parameters for \mbox{\hyperlink{namespacecoord__adapt}{coord\+\_\+adapt}} module. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacecoord__adapt_a7c351df31210e9ecef21fcddb0683003}{build\+\_\+adapt\+\_\+column}} (CS, G, GV, US, tv, i, j, z\+Int, t\+Int, s\+Int, h, z\+Next)
\end{DoxyCompactItemize}


\subsection{Function/\+Subroutine Documentation}
\mbox{\Hypertarget{namespacecoord__adapt_a7c351df31210e9ecef21fcddb0683003}\label{namespacecoord__adapt_a7c351df31210e9ecef21fcddb0683003}} 
\index{coord\+\_\+adapt@{coord\+\_\+adapt}!build\+\_\+adapt\+\_\+column@{build\+\_\+adapt\+\_\+column}}
\index{build\+\_\+adapt\+\_\+column@{build\+\_\+adapt\+\_\+column}!coord\+\_\+adapt@{coord\+\_\+adapt}}
\subsubsection{\texorpdfstring{build\+\_\+adapt\+\_\+column()}{build\_adapt\_column()}}
{\footnotesize\ttfamily subroutine, public coord\+\_\+adapt\+::build\+\_\+adapt\+\_\+column (\begin{DoxyParamCaption}\item[{type(\mbox{\hyperlink{structcoord__adapt_1_1adapt__cs}{adapt\+\_\+cs}}), intent(in)}]{CS,  }\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(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv,  }\item[{integer, intent(in)}]{i,  }\item[{integer, intent(in)}]{j,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, gv \%ke+1), intent(in)}]{z\+Int,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, gv \%ke+1), intent(in)}]{t\+Int,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, gv \%ke+1), intent(in)}]{s\+Int,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, gv \%ke), intent(in)}]{h,  }\item[{real, dimension( gv \%ke+1), intent(inout)}]{z\+Next }\end{DoxyParamCaption})}


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em cs} & The control structure for this module\\
\hline
\mbox{\tt in}  & {\em g} & The ocean\textquotesingle{}s grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em tv} & A structure pointing to various thermodynamic variables\\
\hline
\mbox{\tt in}  & {\em i} & The i-\/index of the column to work on\\
\hline
\mbox{\tt in}  & {\em j} & The j-\/index of the column to work on\\
\hline
\mbox{\tt in}  & {\em zint} & Interface heights \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}.\\
\hline
\mbox{\tt in}  & {\em tint} & Interface temperatures \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em sint} & Interface salinities \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in,out}  & {\em znext} & updated interface positions \\
\hline
\end{DoxyParams}


Definition at line 116 of file coord\+\_\+adapt.\+F90.


\begin{DoxyCode}
116   \textcolor{keywordtype}{type}(adapt\_CS),                              \textcolor{keywordtype}{intent(in)}    :: CS\textcolor{comment}{   !< The control structure for this
       module}
117   \textcolor{keywordtype}{type}(ocean\_grid\_type),                       \textcolor{keywordtype}{intent(in)}    :: G\textcolor{comment}{    !< The ocean's grid structure}
118   \textcolor{keywordtype}{type}(verticalGrid\_type),                     \textcolor{keywordtype}{intent(in)}    :: GV\textcolor{comment}{   !< The ocean's vertical grid structure}
119   \textcolor{keywordtype}{type}(unit\_scale\_type),                       \textcolor{keywordtype}{intent(in)}    :: US\textcolor{comment}{   !< A dimensional unit scaling type}
120   \textcolor{keywordtype}{type}(thermo\_var\_ptrs),                       \textcolor{keywordtype}{intent(in)}    :: tv\textcolor{comment}{   !< A structure pointing to various}
121 \textcolor{comment}{                                                                     !! thermodynamic variables}
122   \textcolor{keywordtype}{integer},                                     \textcolor{keywordtype}{intent(in)}    :: i\textcolor{comment}{    !< The i-index of the column to work
       on}
123   \textcolor{keywordtype}{integer},                                     \textcolor{keywordtype}{intent(in)}    :: j\textcolor{comment}{    !< The j-index of the column to work
       on}
124   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV)+1)}, \textcolor{keywordtype}{intent(in)}    :: zInt\textcolor{comment}{ !< Interface heights [H ~> m or kg
       m-2].}
125   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV)+1)}, \textcolor{keywordtype}{intent(in)}    :: tInt\textcolor{comment}{ !< Interface temperatures [degC]}
126   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV)+1)}, \textcolor{keywordtype}{intent(in)}    :: sInt\textcolor{comment}{ !< Interface salinities [ppt]}
127   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV))},   \textcolor{keywordtype}{intent(in)}    :: h\textcolor{comment}{    !< Layer thicknesses [H ~> m or kg
       m-2]}
128   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)},                 \textcolor{keywordtype}{intent(inout)} :: zNext\textcolor{comment}{ !< updated interface positions}
129 
130   \textcolor{comment}{! Local variables}
131   \textcolor{keywordtype}{integer} :: k, nz
132   \textcolor{keywordtype}{real} :: h\_up, b1, b\_denom\_1, d1, depth, nominal\_z, stretching
133   \textcolor{keywordtype}{real} :: drdz  \textcolor{comment}{! The vertical density gradient [R H-1 ~> kg m-4 or m-1]}
134   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: alpha \textcolor{comment}{! drho/dT [R degC-1 ~> kg m-3 degC-1]}
135   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: beta  \textcolor{comment}{! drho/dS [R ppt-1 ~> kg m-3 ppt-1]}
136   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: del2sigma \textcolor{comment}{! Laplacian of in situ density times grid spacing [R ~> kg m-3]}
137   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: dh\_d2s \textcolor{comment}{! Thickness change in response to del2sigma [H ~> m or kg m-2]}
138   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV))} :: kGrid, c1 \textcolor{comment}{! grid diffusivity on layers, and tridiagonal work array}
139 
140   nz = cs%nk
141 
142   \textcolor{comment}{! set bottom and surface of zNext}
143   znext(1) = 0.
144   znext(nz+1) = zint(i,j,nz+1)
145 
146   \textcolor{comment}{! local depth for scaling diffusivity}
147   depth = g%bathyT(i,j) * gv%Z\_to\_H
148 
149   \textcolor{comment}{! initialize del2sigma and the thickness change response to it zero}
150   del2sigma(:) = 0.0 ; dh\_d2s(:) = 0.0
151 
152   \textcolor{comment}{! calculate del-squared of neutral density by a}
153   \textcolor{comment}{! stencilled finite difference}
154   \textcolor{comment}{! TODO: this needs to be adjusted to account for vanished layers near topography}
155 
156   \textcolor{comment}{! up (j-1)}
157   \textcolor{keywordflow}{if} (g%mask2dT(i,j-1) > 0.) \textcolor{keywordflow}{then}
158     \textcolor{keyword}{call }calculate\_density\_derivs( &
159          0.5 * (tint(i,j,2:nz) + tint(i,j-1,2:nz)), &
160          0.5 * (sint(i,j,2:nz) + sint(i,j-1,2:nz)), &
161          0.5 * (zint(i,j,2:nz) + zint(i,j-1,2:nz)) * (gv%H\_to\_RZ * gv%g\_Earth), &
162          alpha, beta, tv%eqn\_of\_state, (/2,nz/) )
163 
164     del2sigma(2:nz) = del2sigma(2:nz) + &
165          (alpha(2:nz) * (tint(i,j-1,2:nz) - tint(i,j,2:nz)) + &
166           beta(2:nz)  * (sint(i,j-1,2:nz) - sint(i,j,2:nz)))
167 \textcolor{keywordflow}{  endif}
168   \textcolor{comment}{! down (j+1)}
169   \textcolor{keywordflow}{if} (g%mask2dT(i,j+1) > 0.) \textcolor{keywordflow}{then}
170     \textcolor{keyword}{call }calculate\_density\_derivs( &
171          0.5 * (tint(i,j,2:nz) + tint(i,j+1,2:nz)), &
172          0.5 * (sint(i,j,2:nz) + sint(i,j+1,2:nz)), &
173          0.5 * (zint(i,j,2:nz) + zint(i,j+1,2:nz)) * (gv%H\_to\_RZ * gv%g\_Earth), &
174          alpha, beta, tv%eqn\_of\_state, (/2,nz/) )
175 
176     del2sigma(2:nz) = del2sigma(2:nz) + &
177          (alpha(2:nz) * (tint(i,j+1,2:nz) - tint(i,j,2:nz)) + &
178           beta(2:nz)  * (sint(i,j+1,2:nz) - sint(i,j,2:nz)))
179 \textcolor{keywordflow}{  endif}
180   \textcolor{comment}{! left (i-1)}
181   \textcolor{keywordflow}{if} (g%mask2dT(i-1,j) > 0.) \textcolor{keywordflow}{then}
182     \textcolor{keyword}{call }calculate\_density\_derivs( &
183          0.5 * (tint(i,j,2:nz) + tint(i-1,j,2:nz)), &
184          0.5 * (sint(i,j,2:nz) + sint(i-1,j,2:nz)), &
185          0.5 * (zint(i,j,2:nz) + zint(i-1,j,2:nz)) * (gv%H\_to\_RZ * gv%g\_Earth), &
186          alpha, beta, tv%eqn\_of\_state, (/2,nz/) )
187 
188     del2sigma(2:nz) = del2sigma(2:nz) + &
189          (alpha(2:nz) * (tint(i-1,j,2:nz) - tint(i,j,2:nz)) + &
190           beta(2:nz)  * (sint(i-1,j,2:nz) - sint(i,j,2:nz)))
191 \textcolor{keywordflow}{  endif}
192   \textcolor{comment}{! right (i+1)}
193   \textcolor{keywordflow}{if} (g%mask2dT(i+1,j) > 0.) \textcolor{keywordflow}{then}
194     \textcolor{keyword}{call }calculate\_density\_derivs( &
195          0.5 * (tint(i,j,2:nz) + tint(i+1,j,2:nz)), &
196          0.5 * (sint(i,j,2:nz) + sint(i+1,j,2:nz)), &
197          0.5 * (zint(i,j,2:nz) + zint(i+1,j,2:nz)) * (gv%H\_to\_RZ * gv%g\_Earth), &
198          alpha, beta, tv%eqn\_of\_state, (/2,nz/) )
199 
200     del2sigma(2:nz) = del2sigma(2:nz) + &
201          (alpha(2:nz) * (tint(i+1,j,2:nz) - tint(i,j,2:nz)) + &
202           beta(2:nz)  * (sint(i+1,j,2:nz) - sint(i,j,2:nz)))
203 \textcolor{keywordflow}{  endif}
204 
205   \textcolor{comment}{! at this point, del2sigma contains the local neutral density curvature at}
206   \textcolor{comment}{! h-points, on interfaces}
207   \textcolor{comment}{! we need to divide by drho/dz to give an interfacial displacement}
208   \textcolor{comment}{!}
209   \textcolor{comment}{! a positive curvature means we're too light relative to adjacent columns,}
210   \textcolor{comment}{! so del2sigma needs to be positive too (push the interface deeper)}
211   \textcolor{keyword}{call }calculate\_density\_derivs(tint(i,j,:), sint(i,j,:), zint(i,j,:) * (gv%H\_to\_RZ * gv%g\_Earth), &
212        alpha, beta, tv%eqn\_of\_state, (/1,nz+1/) )
213   \textcolor{keywordflow}{do} k = 2, nz
214     \textcolor{comment}{! TODO make lower bound here configurable}
215     dh\_d2s(k) = del2sigma(k) * (0.5 * (h(i,j,k-1) + h(i,j,k))) / &
216          max(alpha(k) * (tv%T(i,j,k) - tv%T(i,j,k-1)) + &
217              beta(k)  * (tv%S(i,j,k) - tv%S(i,j,k-1)), 1e-20*us%kg\_m3\_to\_R)
218 
219     \textcolor{comment}{! don't move the interface so far that it would tangle with another}
220     \textcolor{comment}{! interface in the direction we're moving (or exceed a Nyquist limit}
221     \textcolor{comment}{! that could cause oscillations of the interface)}
222     h\_up = merge(h(i,j,k), h(i,j,k-1), dh\_d2s(k) > 0.)
223     dh\_d2s(k) = 0.5 * cs%adaptAlpha * &
224          sign(min(abs(del2sigma(k)), 0.5 * h\_up), dh\_d2s(k))
225 
226     \textcolor{comment}{! update interface positions so we can diffuse them}
227     znext(k) = zint(i,j,k) + dh\_d2s(k)
228 \textcolor{keywordflow}{  enddo}
229 
230   \textcolor{comment}{! solve diffusivity equation to smooth grid}
231   \textcolor{comment}{! upper diagonal coefficients: -kGrid(2:nz)}
232   \textcolor{comment}{! lower diagonal coefficients: -kGrid(1:nz-1)}
233   \textcolor{comment}{! diagonal coefficients:       1 + (kGrid(1:nz-1) + kGrid(2:nz))}
234   \textcolor{comment}{!}
235   \textcolor{comment}{! first, calculate the diffusivities within layers}
236   \textcolor{keywordflow}{do} k = 1, nz
237     \textcolor{comment}{! calculate the dr bit of drdz}
238     drdz = 0.5 * (alpha(k) + alpha(k+1)) * (tint(i,j,k+1) - tint(i,j,k)) + &
239            0.5 * (beta(k)  + beta(k+1))  * (sint(i,j,k+1) - sint(i,j,k))
240     \textcolor{comment}{! divide by dz from the new interface positions}
241     drdz = drdz / (znext(k) - znext(k+1) + gv%H\_subroundoff)
242     \textcolor{comment}{! don't do weird stuff in unstably-stratified regions}
243     drdz = max(drdz, 0.)
244 
245     \textcolor{comment}{! set vertical grid diffusivity}
246     kgrid(k) = (cs%adaptTimeRatio * nz**2 * depth) * &
247          (cs%adaptZoomCoeff / (cs%adaptZoom + 0.5*(znext(k) + znext(k+1))) + &
248          (cs%adaptBuoyCoeff * drdz / cs%adaptDrho0) + &
249          max(1.0 - cs%adaptZoomCoeff - cs%adaptBuoyCoeff, 0.0) / depth)
250 \textcolor{keywordflow}{  enddo}
251 
252   \textcolor{comment}{! initial denominator (first diagonal element)}
253   b1 = 1.0
254   \textcolor{comment}{! initial Q\_1 = 1 - q\_1 = 1 - 0/1}
255   d1 = 1.0
256   \textcolor{comment}{! work on all interior interfaces}
257   \textcolor{keywordflow}{do} k = 2, nz
258     \textcolor{comment}{! calculate numerator of Q\_k}
259     b\_denom\_1 = 1. + d1 * kgrid(k-1)
260     \textcolor{comment}{! update denominator for k}
261     b1 = 1.0 / (b\_denom\_1 + kgrid(k))
262 
263     \textcolor{comment}{! calculate q\_k}
264     c1(k) = kgrid(k) * b1
265     \textcolor{comment}{! update Q\_k = 1 - q\_k}
266     d1 = b\_denom\_1 * b1
267 
268     \textcolor{comment}{! update RHS}
269     znext(k) = b1 * (znext(k) + kgrid(k-1)*znext(k-1))
270 \textcolor{keywordflow}{  enddo}
271   \textcolor{comment}{! final substitution}
272   \textcolor{keywordflow}{do} k = nz, 2, -1
273     znext(k) = znext(k) + c1(k)*znext(k+1)
274 \textcolor{keywordflow}{  enddo}
275 
276   \textcolor{keywordflow}{if} (cs%adaptDoMin) \textcolor{keywordflow}{then}
277     nominal\_z = 0.
278     stretching = zint(i,j,nz+1) / depth
279 
280     \textcolor{keywordflow}{do} k = 2, nz+1
281       nominal\_z = nominal\_z + cs%coordinateResolution(k-1) * stretching
282       \textcolor{comment}{! take the deeper of the calculated and nominal positions}
283       znext(k) = max(znext(k), nominal\_z)
284       \textcolor{comment}{! interface can't go below topography}
285       znext(k) = min(znext(k), zint(i,j,nz+1))
286 \textcolor{keywordflow}{    enddo}
287 \textcolor{keywordflow}{  endif}
\end{DoxyCode}
\mbox{\Hypertarget{namespacecoord__adapt_a92e88dea18cc4b3adef844f58a8fc0ff}\label{namespacecoord__adapt_a92e88dea18cc4b3adef844f58a8fc0ff}} 
\index{coord\+\_\+adapt@{coord\+\_\+adapt}!end\+\_\+coord\+\_\+adapt@{end\+\_\+coord\+\_\+adapt}}
\index{end\+\_\+coord\+\_\+adapt@{end\+\_\+coord\+\_\+adapt}!coord\+\_\+adapt@{coord\+\_\+adapt}}
\subsubsection{\texorpdfstring{end\+\_\+coord\+\_\+adapt()}{end\_coord\_adapt()}}
{\footnotesize\ttfamily subroutine, public coord\+\_\+adapt\+::end\+\_\+coord\+\_\+adapt (\begin{DoxyParamCaption}\item[{type(\mbox{\hyperlink{structcoord__adapt_1_1adapt__cs}{adapt\+\_\+cs}}), pointer}]{CS }\end{DoxyParamCaption})}



Clean up the coordinate control structure. 


\begin{DoxyParams}{Parameters}
{\em cs} & The control structure for this module \\
\hline
\end{DoxyParams}


Definition at line 80 of file coord\+\_\+adapt.\+F90.


\begin{DoxyCode}
80   \textcolor{keywordtype}{type}(adapt\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{  !< The control structure for this module}
81 
82   \textcolor{comment}{! nothing to do}
83   \textcolor{keywordflow}{if} (.not. \textcolor{keyword}{associated}(cs)) \textcolor{keywordflow}{return}
84   \textcolor{keyword}{deallocate}(cs%coordinateResolution)
85   \textcolor{keyword}{deallocate}(cs)
\end{DoxyCode}
\mbox{\Hypertarget{namespacecoord__adapt_a53eb167f67898f13a7f2258365f1b968}\label{namespacecoord__adapt_a53eb167f67898f13a7f2258365f1b968}} 
\index{coord\+\_\+adapt@{coord\+\_\+adapt}!init\+\_\+coord\+\_\+adapt@{init\+\_\+coord\+\_\+adapt}}
\index{init\+\_\+coord\+\_\+adapt@{init\+\_\+coord\+\_\+adapt}!coord\+\_\+adapt@{coord\+\_\+adapt}}
\subsubsection{\texorpdfstring{init\+\_\+coord\+\_\+adapt()}{init\_coord\_adapt()}}
{\footnotesize\ttfamily subroutine, public coord\+\_\+adapt\+::init\+\_\+coord\+\_\+adapt (\begin{DoxyParamCaption}\item[{type(\mbox{\hyperlink{structcoord__adapt_1_1adapt__cs}{adapt\+\_\+cs}}), pointer}]{CS,  }\item[{integer, intent(in)}]{nk,  }\item[{real, dimension(\+:), intent(in)}]{coordinate\+Resolution,  }\item[{real, intent(in)}]{m\+\_\+to\+\_\+H,  }\item[{real, intent(in)}]{kg\+\_\+m3\+\_\+to\+\_\+R }\end{DoxyParamCaption})}



Initialise an adapt\+\_\+\+CS with parameters. 


\begin{DoxyParams}[1]{Parameters}
 & {\em cs} & Unassociated pointer to hold the control structure\\
\hline
\mbox{\tt in}  & {\em nk} & Number of layers in the grid\\
\hline
\mbox{\tt in}  & {\em coordinateresolution} & Nominal near-\/surface resolution \mbox{[}m\mbox{]} or other units specified with m\+\_\+to\+\_\+H\\
\hline
\mbox{\tt in}  & {\em m\+\_\+to\+\_\+h} & A conversion factor from m to the units of thicknesses\\
\hline
\mbox{\tt in}  & {\em kg\+\_\+m3\+\_\+to\+\_\+r} & A conversion factor from kg m-\/3 to the units of density \\
\hline
\end{DoxyParams}


Definition at line 54 of file coord\+\_\+adapt.\+F90.


\begin{DoxyCode}
54   \textcolor{keywordtype}{type}(adapt\_CS),     \textcolor{keywordtype}{pointer}    :: CS\textcolor{comment}{ !< Unassociated pointer to hold the control structure}
55   \textcolor{keywordtype}{integer},            \textcolor{keywordtype}{intent(in)} :: nk\textcolor{comment}{ !< Number of layers in the grid}
56   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:)}, \textcolor{keywordtype}{intent(in)} :: coordinateResolution\textcolor{comment}{ !< Nominal near-surface resolution [m] or}
57 \textcolor{comment}{                                       !! other units specified with m\_to\_H}
58   \textcolor{keywordtype}{real},               \textcolor{keywordtype}{intent(in)} :: m\_to\_H\textcolor{comment}{ !< A conversion factor from m to the units of thicknesses}
59   \textcolor{keywordtype}{real},               \textcolor{keywordtype}{intent(in)} :: kg\_m3\_to\_R\textcolor{comment}{ !< A conversion factor from kg m-3 to the units of density}
60 
61   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{"init\_coord\_adapt: CS already associated"})
62   \textcolor{keyword}{allocate}(cs)
63   \textcolor{keyword}{allocate}(cs%coordinateResolution(nk))
64 
65   cs%nk = nk
66   cs%coordinateResolution(:) = coordinateresolution(:)
67 
68   \textcolor{comment}{! Set real parameter default values}
69   cs%adaptTimeRatio = 1e-1 \textcolor{comment}{! Nondim.}
70   cs%adaptAlpha     = 1.0  \textcolor{comment}{! Nondim.}
71   cs%adaptZoom      = 200.0 * m\_to\_h    \textcolor{comment}{! [H ~> m or kg m-2]}
72   cs%adaptZoomCoeff = 0.0  \textcolor{comment}{! Nondim.}
73   cs%adaptBuoyCoeff = 0.0  \textcolor{comment}{! Nondim.}
74   cs%adaptDrho0     = 0.5 * kg\_m3\_to\_r  \textcolor{comment}{! [R ~> kg m-3]}
75 
\end{DoxyCode}
\mbox{\Hypertarget{namespacecoord__adapt_a35706e0359c3953aac56b160e48b0ef8}\label{namespacecoord__adapt_a35706e0359c3953aac56b160e48b0ef8}} 
\index{coord\+\_\+adapt@{coord\+\_\+adapt}!set\+\_\+adapt\+\_\+params@{set\+\_\+adapt\+\_\+params}}
\index{set\+\_\+adapt\+\_\+params@{set\+\_\+adapt\+\_\+params}!coord\+\_\+adapt@{coord\+\_\+adapt}}
\subsubsection{\texorpdfstring{set\+\_\+adapt\+\_\+params()}{set\_adapt\_params()}}
{\footnotesize\ttfamily subroutine, public coord\+\_\+adapt\+::set\+\_\+adapt\+\_\+params (\begin{DoxyParamCaption}\item[{type(\mbox{\hyperlink{structcoord__adapt_1_1adapt__cs}{adapt\+\_\+cs}}), pointer}]{CS,  }\item[{real, intent(in), optional}]{adapt\+Time\+Ratio,  }\item[{real, intent(in), optional}]{adapt\+Alpha,  }\item[{real, intent(in), optional}]{adapt\+Zoom,  }\item[{real, intent(in), optional}]{adapt\+Zoom\+Coeff,  }\item[{real, intent(in), optional}]{adapt\+Buoy\+Coeff,  }\item[{real, intent(in), optional}]{adapt\+Drho0,  }\item[{logical, intent(in), optional}]{adapt\+Do\+Min }\end{DoxyParamCaption})}



This subtroutine can be used to set the parameters for \mbox{\hyperlink{namespacecoord__adapt}{coord\+\_\+adapt}} module. 


\begin{DoxyParams}[1]{Parameters}
 & {\em cs} & The control structure for this module\\
\hline
\mbox{\tt in}  & {\em adapttimeratio} & Ratio of optimisation and diffusion timescales\\
\hline
\mbox{\tt in}  & {\em adaptalpha} & Nondimensional coefficient determining how much optimisation to apply\\
\hline
\mbox{\tt in}  & {\em adaptzoom} & Near-\/surface zooming depth \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em adaptzoomcoeff} & Near-\/surface zooming coefficient\\
\hline
\mbox{\tt in}  & {\em adaptbuoycoeff} & Stratification-\/dependent diffusion coefficient\\
\hline
\mbox{\tt in}  & {\em adaptdrho0} & Reference density difference for stratification-\/dependent diffusion \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]}\\
\hline
\mbox{\tt in}  & {\em adaptdomin} & If true, form a H\+Y\+C\+O\+M1-\/like mixed layer by preventing interfaces from becoming shallower than the depths set by coordinate\+Resolution \\
\hline
\end{DoxyParams}


Definition at line 91 of file coord\+\_\+adapt.\+F90.


\begin{DoxyCode}
91   \textcolor{keywordtype}{type}(adapt\_CS),    \textcolor{keywordtype}{pointer}    :: CS\textcolor{comment}{  !< The control structure for this module}
92   \textcolor{keywordtype}{real},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: adaptTimeRatio\textcolor{comment}{ !< Ratio of optimisation and diffusion timescales}
93   \textcolor{keywordtype}{real},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: adaptAlpha\textcolor{comment}{     !< Nondimensional coefficient determining}
94 \textcolor{comment}{                                                  !! how much optimisation to apply}
95   \textcolor{keywordtype}{real},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: adaptZoom\textcolor{comment}{      !< Near-surface zooming depth [H ~> m or kg m-2]}
96   \textcolor{keywordtype}{real},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: adaptZoomCoeff\textcolor{comment}{ !< Near-surface zooming coefficient}
97   \textcolor{keywordtype}{real},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: adaptBuoyCoeff\textcolor{comment}{ !< Stratification-dependent diffusion coefficient}
98   \textcolor{keywordtype}{real},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: adaptDrho0\textcolor{comment}{  !< Reference density difference for}
99 \textcolor{comment}{                                               !! stratification-dependent diffusion [R ~> kg m-3]}
100   \textcolor{keywordtype}{logical}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: adaptDoMin\textcolor{comment}{  !< If true, form a HYCOM1-like mixed layer by}
101 \textcolor{comment}{                                               !! preventing interfaces from becoming shallower than}
102 \textcolor{comment}{                                               !! the depths set by coordinateResolution}
103 
104   \textcolor{keywordflow}{if} (.not. \textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{"set\_adapt\_params: CS not associated"})
105 
106   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(adapttimeratio)) cs%adaptTimeRatio = adapttimeratio
107   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(adaptalpha)) cs%adaptAlpha = adaptalpha
108   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(adaptzoom)) cs%adaptZoom = adaptzoom
109   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(adaptzoomcoeff)) cs%adaptZoomCoeff = adaptzoomcoeff
110   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(adaptbuoycoeff)) cs%adaptBuoyCoeff = adaptbuoycoeff
111   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(adaptdrho0)) cs%adaptDrho0 = adaptdrho0
112   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(adaptdomin)) cs%adaptDoMin = adaptdomin
\end{DoxyCode}
