+Add particle code option to advect with uhtr (#492)

* +Add particle code option to advect with uhtr

The particle code has so far used the same velocity has was used in the dynamics step. I would like to add the option for the particle code to use uhtr/h and vhtr/h, so that the velocities used to advect particles may include the effects of parameterized eddies.

To make this work, I have added a flag that controls which velocity to use and moved the particles_run step to take place after uhtr and vhtr are defined. The interfaces in the code and in config_src/external are updated to pass this information to the drifters package.

config_src/external/drifters/MOM_particles.F90

@@ -28,14 +28,19 @@ subroutine particles_init(parts, Grid, Time, dt, u, v, h)
 end subroutine particles_init
 
 !> The main driver the steps updates particles
-subroutine particles_run(parts, time, uo, vo, ho, tv, stagger)
+subroutine particles_run(parts, time, uo, vo, ho, tv, use_uh, stagger)
   ! Arguments
   type(particles), pointer :: parts !< Container for all types and memory
   type(time_type), intent(in) :: time !< Model time
-  real, dimension(:,:,:), intent(in) :: uo !< Ocean zonal velocity [L T-1 ~>m s-1]
-  real, dimension(:,:,:), intent(in) :: vo !< Ocean meridional velocity [L T-1~> m s-1]
+  real, dimension(:,:,:), intent(in) :: uo !< If use_uh is false, ocean zonal velocity [L T-1 ~>m s-1].
+                                           !! If use_uh is true, accumulated zonal thickness fluxes
+                                           !! that are used to advect tracers [H L2 ~> m3 or kg]
+  real, dimension(:,:,:), intent(in) :: vo !< If use_uh is false, ocean meridional velocity [L T-1 ~>m s-1].
+                                           !! If use_uh is true, accumulated meridional thickness fluxes
+                                           !! that are used to advect tracers [H L2 ~> m3 or kg]
   real, dimension(:,:,:), intent(in) :: ho !< Ocean layer thickness [H ~> m or kg m-2]
   type(thermo_var_ptrs),  intent(in) :: tv !< structure containing pointers to available thermodynamic fields
+  logical :: use_uh !< Flag for whether u and v are weighted by thickness
   integer, optional, intent(in) :: stagger !< Flag for whether velocities are staggered
 
 end subroutine particles_run

src/core/MOM.F90

@@ -360,6 +360,7 @@ module MOM
                                 !! higher values use more appropriate expressions that differ at
                                 !! roundoff for non-Boussinesq cases.
   logical :: use_particles      !< Turns on the particles package
+  logical :: use_uh_particles   !< particles are advected by uh/h
   logical :: use_dbclient       !< Turns on the database client used for ML inference/analysis
   character(len=10) :: particle_type !< Particle types include: surface(default), profiling and sail drone.
 
@@ -1266,10 +1267,6 @@ subroutine step_MOM_dynamics(forces, p_surf_begin, p_surf_end, dt, dt_thermo, &
     enddo; enddo
   endif
 
-  if (CS%use_particles .and. CS%do_dynamics) then ! Run particles whether or not stepping is split
-    call particles_run(CS%particles, Time_local, CS%u, CS%v, CS%h, CS%tv) ! Run the particles model
-  endif
-
 
   if ((CS%thickness_diffuse .or. CS%interface_filter) .and. &
       .not.CS%thickness_diffuse_first) then
@@ -1331,6 +1328,17 @@ subroutine step_MOM_dynamics(forces, p_surf_begin, p_surf_end, dt, dt_thermo, &
   endif
   call disable_averaging(CS%diag)
 
+  if (CS%use_particles .and. CS%do_dynamics .and. CS%use_uh_particles) then
+    !Run particles using thickness-weighted velocity
+    call particles_run(CS%particles, Time_local, CS%uhtr, CS%vhtr, CS%h, &
+        CS%tv, CS%use_uh_particles)
+  elseif (CS%use_particles .and. CS%do_dynamics) then
+    !Run particles using unweighted velocity
+    call particles_run(CS%particles, Time_local, CS%u, CS%v, CS%h, &
+                       CS%tv, CS%use_uh_particles)
+  endif
+
+
   ! Advance the dynamics time by dt.
   CS%t_dyn_rel_adv = CS%t_dyn_rel_adv + dt
   CS%n_dyn_steps_in_adv = CS%n_dyn_steps_in_adv + 1
@@ -2440,7 +2448,8 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, &
 
   call get_param(param_file, "MOM", "USE_PARTICLES", CS%use_particles, &
                  "If true, use the particles package.", default=.false.)
-
+  call get_param(param_file, "MOM", "USE_UH_PARTICLES", CS%use_uh_particles, &
+                 "If true, use the uh velocity in the particles package.",default=.false.)
   CS%ensemble_ocean=.false.
   call get_param(param_file, "MOM", "ENSEMBLE_OCEAN", CS%ensemble_ocean, &
                  "If False, The model is being run in serial mode as a single realization. "//&