MOM6
user_surface_forcing.F90
Go to the documentation of this file.
1 !> Template for user to code up surface forcing.
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
7 use mom_diag_mediator, only : register_diag_field, diag_ctrl, safe_alloc_ptr
8 use mom_domains, only : pass_var, pass_vector, agrid
9 use mom_error_handler, only : mom_error, fatal, warning, is_root_pe
13 use mom_grid, only : ocean_grid_type
14 use mom_io, only : file_exists, read_data
15 use mom_time_manager, only : time_type, operator(+), operator(/)
19 use mom_variables, only : surface
20 
21 implicit none ; private
22 
24 
25 !> This control structure should be used to store any run-time variables
26 !! associated with the user-specified forcing.
27 !!
28 !! It can be readily modified for a specific case, and because it is private there
29 !! will be no changes needed in other code (although they will have to be recompiled).
30 type, public :: user_surface_forcing_cs ; private
31  ! The variables in the cannonical example are used for some common
32  ! cases, but do not need to be used.
33 
34  logical :: use_temperature !< If true, temperature and salinity are used as state variables.
35  logical :: restorebuoy !< If true, use restoring surface buoyancy forcing.
36  real :: rho0 !< The density used in the Boussinesq approximation [R ~> kg m-3].
37  real :: g_earth !< The gravitational acceleration [L2 Z-1 s-2 ~> m s-2].
38  real :: flux_const !< The restoring rate at the surface [Z T-1 ~> m s-1].
39  real :: gust_const !< A constant unresolved background gustiness
40  !! that contributes to ustar [R L Z T-1 ~> Pa].
41 
42  type(diag_ctrl), pointer :: diag !< A structure that is used to regulate the
43  !! timing of diagnostic output.
45 
46 contains
47 
48 !> This subroutine sets the surface wind stresses, forces%taux and forces%tauy, in [R Z L T-2 ~> Pa].
49 !! These are the stresses in the direction of the model grid (i.e. the same
50 !! direction as the u- and v- velocities).
51 subroutine user_wind_forcing(sfc_state, forces, day, G, US, CS)
52  type(surface), intent(inout) :: sfc_state !< A structure containing fields that
53  !! describe the surface state of the ocean.
54  type(mech_forcing), intent(inout) :: forces !< A structure with the driving mechanical forces
55  type(time_type), intent(in) :: day !< The time of the fluxes
56  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure
57  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
58  type(user_surface_forcing_cs), pointer :: cs !< A pointer to the control structure returned
59  !! by a previous call to user_surface_forcing_init
60 
61  ! Local variables
62  integer :: i, j, is, ie, js, je, isq, ieq, jsq, jeq
63 
64  ! When modifying the code, comment out this error message. It is here
65  ! so that the original (unmodified) version is not accidentally used.
66  call mom_error(fatal, "User_wind_surface_forcing: " // &
67  "User forcing routine called without modification." )
68 
69  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
70  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
71 
72  ! Allocate the forcing arrays, if necessary.
73  call allocate_mech_forcing(g, forces, stress=.true., ustar=.true.)
74 
75  ! Set the surface wind stresses, in units of [R L Z T-1 ~> Pa]. A positive taux
76  ! accelerates the ocean to the (pseudo-)east.
77 
78  ! The i-loop extends to is-1 so that taux can be used later in the
79  ! calculation of ustar - otherwise the lower bound would be Isq.
80  do j=js,je ; do i=is-1,ieq
81  ! Change this to the desired expression.
82  forces%taux(i,j) = g%mask2dCu(i,j) * 0.0*us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z
83  enddo ; enddo
84  do j=js-1,jeq ; do i=is,ie
85  forces%tauy(i,j) = g%mask2dCv(i,j) * 0.0 ! Change this to the desired expression.
86  enddo ; enddo
87 
88  ! Set the surface friction velocity, in units of [Z T-1 ~> m s-1]. ustar
89  ! is always positive.
90  if (associated(forces%ustar)) then ; do j=js,je ; do i=is,ie
91  ! This expression can be changed if desired, but need not be.
92  forces%ustar(i,j) = g%mask2dT(i,j) * sqrt((cs%gust_const + &
93  sqrt(0.5*(forces%taux(i-1,j)**2 + forces%taux(i,j)**2) + &
94  0.5*(forces%tauy(i,j-1)**2 + forces%tauy(i,j)**2))) * (us%L_to_Z/cs%Rho0))
95  enddo ; enddo ; endif
96 
97 end subroutine user_wind_forcing
98 
99 !> This subroutine specifies the current surface fluxes of buoyancy or
100 !! temperature and fresh water. It may also be modified to add
101 !! surface fluxes of user provided tracers.
102 subroutine user_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS)
103  type(surface), intent(inout) :: sfc_state !< A structure containing fields that
104  !! describe the surface state of the ocean.
105  type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields
106  type(time_type), intent(in) :: day !< The time of the fluxes
107  real, intent(in) :: dt !< The amount of time over which
108  !! the fluxes apply [s]
109  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
110  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
111  type(user_surface_forcing_cs), pointer :: cs !< A pointer to the control structure returned
112  !! by a previous call to user_surface_forcing_init
113 
114 ! This subroutine specifies the current surface fluxes of buoyancy or
115 ! temperature and fresh water. It may also be modified to add
116 ! surface fluxes of user provided tracers.
117 
118 ! When temperature is used, there are long list of fluxes that need to be
119 ! set - essentially the same as for a full coupled model, but most of these
120 ! can be simply set to zero. The net fresh water flux should probably be
121 ! set in fluxes%evap and fluxes%lprec, with any salinity restoring
122 ! appearing in fluxes%vprec, and the other water flux components
123 ! (fprec, lrunoff and frunoff) left as arrays full of zeros.
124 ! Evap is usually negative and precip is usually positive. All heat fluxes
125 ! are in W m-2 and positive for heat going into the ocean. All fresh water
126 ! fluxes are in [R Z T-1 ~> kg m-2 s-1] and positive for water moving into the ocean.
127 
128  ! Local variables
129  real :: temp_restore ! The temperature that is being restored toward [degC].
130  real :: salin_restore ! The salinity that is being restored toward [ppt]
131  real :: density_restore ! The potential density that is being restored
132  ! toward [kg m-3].
133  real :: rhoxcp ! The mean density times the heat capacity [J m-3 degC-1].
134  real :: rho0_mks ! The mean density in MKS units [kg m-3]
135  real :: buoy_rest_const ! A constant relating density anomalies to the
136  ! restoring buoyancy flux [L2 m3 T-3 kg-1 ~> m5 s-3 kg-1].
137 
138  integer :: i, j, is, ie, js, je
139  integer :: isd, ied, jsd, jed
140 
141  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
142  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
143  rho0_mks = cs%Rho0 * us%R_to_kg_m3
144 
145  ! When modifying the code, comment out this error message. It is here
146  ! so that the original (unmodified) version is not accidentally used.
147  call mom_error(fatal, "User_buoyancy_surface_forcing: " // &
148  "User forcing routine called without modification." )
149 
150  ! Allocate and zero out the forcing arrays, as necessary. This portion is
151  ! usually not changed.
152  if (cs%use_temperature) then
153  call safe_alloc_ptr(fluxes%evap, isd, ied, jsd, jed)
154  call safe_alloc_ptr(fluxes%lprec, isd, ied, jsd, jed)
155  call safe_alloc_ptr(fluxes%fprec, isd, ied, jsd, jed)
156  call safe_alloc_ptr(fluxes%lrunoff, isd, ied, jsd, jed)
157  call safe_alloc_ptr(fluxes%frunoff, isd, ied, jsd, jed)
158  call safe_alloc_ptr(fluxes%vprec, isd, ied, jsd, jed)
159 
160  call safe_alloc_ptr(fluxes%sw, isd, ied, jsd, jed)
161  call safe_alloc_ptr(fluxes%lw, isd, ied, jsd, jed)
162  call safe_alloc_ptr(fluxes%latent, isd, ied, jsd, jed)
163  call safe_alloc_ptr(fluxes%sens, isd, ied, jsd, jed)
164  else ! This is the buoyancy only mode.
165  call safe_alloc_ptr(fluxes%buoy, isd, ied, jsd, jed)
166  endif
167 
168 
169  ! MODIFY THE CODE IN THE FOLLOWING LOOPS TO SET THE BUOYANCY FORCING TERMS.
170 
171  if ( cs%use_temperature ) then
172  ! Set whichever fluxes are to be used here. Any fluxes that
173  ! are always zero do not need to be changed here.
174  do j=js,je ; do i=is,ie
175  ! Fluxes of fresh water through the surface are in units of [R Z T-1 ~> kg m-2 s-1]
176  ! and are positive downward - i.e. evaporation should be negative.
177  fluxes%evap(i,j) = -0.0 * g%mask2dT(i,j)
178  fluxes%lprec(i,j) = 0.0 * g%mask2dT(i,j)
179 
180  ! vprec will be set later, if it is needed for salinity restoring.
181  fluxes%vprec(i,j) = 0.0
182 
183  ! Heat fluxes are in units of W m-2 and are positive into the ocean.
184  fluxes%lw(i,j) = 0.0 * g%mask2dT(i,j)
185  fluxes%latent(i,j) = 0.0 * g%mask2dT(i,j)
186  fluxes%sens(i,j) = 0.0 * g%mask2dT(i,j)
187  fluxes%sw(i,j) = 0.0 * g%mask2dT(i,j)
188  enddo ; enddo
189  else ! This is the buoyancy only mode.
190  do j=js,je ; do i=is,ie
191  ! fluxes%buoy is the buoyancy flux into the ocean [L2 T-3 ~> m2 s-3]. A positive
192  ! buoyancy flux is of the same sign as heating the ocean.
193  fluxes%buoy(i,j) = 0.0 * g%mask2dT(i,j)
194  enddo ; enddo
195  endif
196 
197  if (cs%restorebuoy) then
198  if (cs%use_temperature) then
199  call safe_alloc_ptr(fluxes%heat_added, isd, ied, jsd, jed)
200  ! When modifying the code, comment out this error message. It is here
201  ! so that the original (unmodified) version is not accidentally used.
202  call mom_error(fatal, "User_buoyancy_surface_forcing: " // &
203  "Temperature and salinity restoring used without modification." )
204 
205  rhoxcp = rho0_mks * fluxes%C_p
206  do j=js,je ; do i=is,ie
207  ! Set Temp_restore and Salin_restore to the temperature (in degC) and
208  ! salinity (in PSU or ppt) that are being restored toward.
209  temp_restore = 0.0
210  salin_restore = 0.0
211 
212  fluxes%heat_added(i,j) = (g%mask2dT(i,j) * (rhoxcp * us%Z_to_m*us%s_to_T*cs%Flux_const)) * &
213  (temp_restore - sfc_state%SST(i,j))
214  fluxes%vprec(i,j) = - (g%mask2dT(i,j) * (cs%Rho0*cs%Flux_const)) * &
215  ((salin_restore - sfc_state%SSS(i,j)) / (0.5 * (salin_restore + sfc_state%SSS(i,j))))
216  enddo ; enddo
217  else
218  ! When modifying the code, comment out this error message. It is here
219  ! so that the original (unmodified) version is not accidentally used.
220  call mom_error(fatal, "User_buoyancy_surface_forcing: " // &
221  "Buoyancy restoring used without modification." )
222 
223  ! The -1 is because density has the opposite sign to buoyancy.
224  buoy_rest_const = -1.0 * (cs%G_Earth * cs%Flux_const) / rho0_mks
225  do j=js,je ; do i=is,ie
226  ! Set density_restore to an expression for the surface potential
227  ! density [kg m-3] that is being restored toward.
228  density_restore = 1030.0
229 
230  fluxes%buoy(i,j) = g%mask2dT(i,j) * buoy_rest_const * &
231  (density_restore - sfc_state%sfc_density(i,j))
232  enddo ; enddo
233  endif
234  endif ! end RESTOREBUOY
235 
236 end subroutine user_buoyancy_forcing
237 
238 !> This subroutine initializes the USER_surface_forcing module
239 subroutine user_surface_forcing_init(Time, G, US, param_file, diag, CS)
240  type(time_type), intent(in) :: time !< The current model time
241  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
242  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
243  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
244  type(diag_ctrl), target, intent(in) :: diag !< A structure that is used to regulate diagnostic output.
245  type(user_surface_forcing_cs), pointer :: cs !< A pointer that is set to point to
246  !! the control structure for this module
247 
248 ! This include declares and sets the variable "version".
249 #include "version_variable.h"
250  character(len=40) :: mdl = "user_surface_forcing" ! This module's name.
251 
252  if (associated(cs)) then
253  call mom_error(warning, "USER_surface_forcing_init called with an associated "// &
254  "control structure.")
255  return
256  endif
257  allocate(cs)
258  cs%diag => diag
259 
260  ! Read all relevant parameters and write them to the model log.
261  call log_version(param_file, mdl, version, "")
262  call get_param(param_file, mdl, "ENABLE_THERMODYNAMICS", cs%use_temperature, &
263  "If true, Temperature and salinity are used as state "//&
264  "variables.", default=.true.)
265 
266  call get_param(param_file, mdl, "G_EARTH", cs%G_Earth, &
267  "The gravitational acceleration of the Earth.", &
268  units="m s-2", default = 9.80, scale=us%m_to_L**2*us%Z_to_m*us%T_to_s**2)
269  call get_param(param_file, mdl, "RHO_0", cs%Rho0, &
270  "The mean ocean density used with BOUSSINESQ true to "//&
271  "calculate accelerations and the mass for conservation "//&
272  "properties, or with BOUSSINSEQ false to convert some "//&
273  "parameters from vertical units of m to kg m-2.", &
274  units="kg m-3", default=1035.0, scale=us%kg_m3_to_R)
275  call get_param(param_file, mdl, "GUST_CONST", cs%gust_const, &
276  "The background gustiness in the winds.", &
277  units="Pa", default=0.02, scale=us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z)
278 
279  call get_param(param_file, mdl, "RESTOREBUOY", cs%restorebuoy, &
280  "If true, the buoyancy fluxes drive the model back "//&
281  "toward some specified surface state with a rate "//&
282  "given by FLUXCONST.", default= .false.)
283  if (cs%restorebuoy) then
284  call get_param(param_file, mdl, "FLUXCONST", cs%Flux_const, &
285  "The constant that relates the restoring surface fluxes "//&
286  "to the relative surface anomalies (akin to a piston "//&
287  "velocity). Note the non-MKS units.", &
288  units="m day-1", scale=us%m_to_Z*us%T_to_s, fail_if_missing=.true.)
289  ! Convert CS%Flux_const from m day-1 to m s-1.
290  cs%Flux_const = cs%Flux_const / 86400.0
291  endif
292 
293 end subroutine user_surface_forcing_init
294 
295 !! \namespace user_surface_forcing
296 !!
297 !! This file contains the subroutines that a user should modify to
298 !! to set the surface wind stresses and fluxes of buoyancy or
299 !! temperature and fresh water. They are called when the run-time
300 !! parameters WIND_CONFIG or BUOY_CONFIG are set to "USER". The
301 !! standard version has simple examples, along with run-time error
302 !! messages that will cause the model to abort if this code has no
303 !! been modified. This code is intended for use with relatively
304 !! simple specifications of the forcing. For more complicated forms,
305 !! it is probably a good idea to read the forcing from input files
306 !! using "file" for WIND_CONFIG and BUOY_CONFIG.
307 !!
308 !! USER_wind_forcing() should set the surface wind stresses (taux and
309 !! tauy) perhaps along with the surface friction velocity (ustar).
310 !!
311 !! USER_buoyancy() forcing is used to set the surface buoyancy
312 !! forcing, which may include a number of fresh water flux fields
313 !! (evap, lprec, fprec, lrunoff, frunoff, and
314 !! vprec) and the surface heat fluxes (sw, lw, latent and sens)
315 !! if temperature and salinity are state variables, or it may simply
316 !! be the buoyancy flux if it is not. This routine also has coded a
317 !! restoring to surface values of temperature and salinity.
318 
319 end module user_surface_forcing
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_forcing_type::allocate_forcing_type
subroutine, public allocate_forcing_type(G, fluxes, water, heat, ustar, press, shelf, iceberg, salt)
Conditionally allocate fields within the forcing type.
Definition: MOM_forcing_type.F90:2811
mom_forcing_type::mech_forcing
Structure that contains pointers to the mechanical forcing at the surface used to drive the liquid oc...
Definition: MOM_forcing_type.F90:187
mom_variables::surface
Pointers to various fields which may be used describe the surface state of MOM, and which will be ret...
Definition: MOM_variables.F90:38
user_surface_forcing::user_buoyancy_forcing
subroutine, public user_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS)
This subroutine specifies the current surface fluxes of buoyancy or temperature and fresh water....
Definition: user_surface_forcing.F90:103
mom_diag_mediator::query_averaging_enabled
logical function, public query_averaging_enabled(diag_cs, time_int, time_end)
Call this subroutine to determine whether the averaging is currently enabled. .true....
Definition: MOM_diag_mediator.F90:1850
user_surface_forcing::user_surface_forcing_cs
This control structure should be used to store any run-time variables associated with the user-specif...
Definition: user_surface_forcing.F90:30
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_diag_mediator
The subroutines here provide convenient wrappers to the fms diag_manager interfaces with additional d...
Definition: MOM_diag_mediator.F90:3
mom_forcing_type::allocate_mech_forcing
subroutine, public allocate_mech_forcing(G, forces, stress, ustar, shelf, press, iceberg)
Conditionally allocate fields within the mechanical forcing type.
Definition: MOM_forcing_type.F90:2879
mom_diag_mediator::register_diag_field
integer function, public register_diag_field(module_name, field_name, axes_in, init_time, long_name, units, missing_value, range, mask_variant, standard_name, verbose, do_not_log, err_msg, interp_method, tile_count, cmor_field_name, cmor_long_name, cmor_units, cmor_standard_name, cell_methods, x_cell_method, y_cell_method, v_cell_method, conversion, v_extensive)
Returns the "diag_mediator" handle for a group (native, CMOR, z-coord, ...) of diagnostics derived fr...
Definition: MOM_diag_mediator.F90:1878
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_diag_mediator::post_data
Make a diagnostic available for averaging or output.
Definition: MOM_diag_mediator.F90:70
mom_domains::pass_vector
Do a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:54
mom_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
user_surface_forcing::user_surface_forcing_init
subroutine, public user_surface_forcing_init(Time, G, US, param_file, diag, CS)
This subroutine initializes the USER_surface_forcing module.
Definition: user_surface_forcing.F90:240
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_tracer_flow_control
Orchestrates the registration and calling of tracer packages.
Definition: MOM_tracer_flow_control.F90:2
user_surface_forcing
Template for user to code up surface forcing.
Definition: user_surface_forcing.F90:2
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_tracer_flow_control::tracer_flow_control_cs
The control structure for orchestrating the calling of tracer packages.
Definition: MOM_tracer_flow_control.F90:75
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_forcing_type::forcing
Structure that contains pointers to the boundary forcing used to drive the liquid ocean simulated by ...
Definition: MOM_forcing_type.F90:50
mom_io::file_exists
Indicate whether a file exists, perhaps with domain decomposition.
Definition: MOM_io.F90:68
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
user_surface_forcing::user_wind_forcing
subroutine, public user_wind_forcing(sfc_state, forces, day, G, US, CS)
This subroutine sets the surface wind stresses, forcestaux and forcestauy, in [R Z L T-2 ~> Pa]....
Definition: user_surface_forcing.F90:52
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_tracer_flow_control::call_tracer_set_forcing
subroutine, public call_tracer_set_forcing(state, fluxes, day_start, day_interval, G, CS)
This subroutine calls the individual tracer modules' subroutines to specify or read quantities relate...
Definition: MOM_tracer_flow_control.F90:382
mom_diag_mediator::diag_ctrl
The following data type a list of diagnostic fields an their variants, as well as variables that cont...
Definition: MOM_diag_mediator.F90:239