idealized_hurricane Module Reference

Detailed Description

Forcing for the idealized hurricane and SCM_idealized_hurricane examples.

Data Types
type	idealized_hurricane_cs
	Container for parameters describing idealized wind structure. More...

Functions/Subroutines
subroutine, public	idealized_hurricane_wind_init (Time, G, US, param_file, CS)
	Initializes wind profile for the SCM idealized hurricane example. More...
subroutine, public	idealized_hurricane_wind_forcing (state, forces, day, G, US, CS)
	Computes the surface wind for the idealized hurricane test cases. More...
subroutine	idealized_hurricane_wind_profile (CS, absf, YY, XX, UOCN, VOCN, Tx, Ty)
	Calculate the wind speed at a location as a function of time. More...
subroutine, public	scm_idealized_hurricane_wind_forcing (state, forces, day, G, US, CS)
	This subroutine is primarily needed as a legacy for reproducing answers. It is included as an additional subroutine rather than padded into the previous routine with flags to ease its eventual removal. Its functionality is replaced with the new routines and it can be deleted when answer changes are acceptable. More...

Variables
character(len=40)	mdl = "idealized_hurricane"
	This module's name. More...

Function/Subroutine Documentation

idealized_hurricane_wind_forcing()

subroutine, public idealized_hurricane::idealized_hurricane_wind_forcing	(	type(surface), intent(in)	state,
		type(mech_forcing), intent(inout)	forces,
		type(time_type), intent(in)	day,
		type(ocean_grid_type), intent(inout)	G,
		type(unit_scale_type), intent(in)	US,
		type(idealized_hurricane_cs), pointer	CS
	)

Computes the surface wind for the idealized hurricane test cases.

Parameters

[in]	state	Surface state structure
[in,out]	forces	A structure with the driving mechanical forces
[in]	day	Time in days
[in,out]	g	Grid structure
[in]	us	A dimensional unit scaling type
	cs	Container for idealized hurricane parameters

Compute storm center location

Computes taux

Computes tauy

Get Ustar

Definition at line 197 of file Idealized_Hurricane.F90.

  type(surface),                intent(in)    :: state  !< Surface state structure
  type(mech_forcing),           intent(inout) :: forces !< A structure with the driving mechanical forces
  type(time_type),              intent(in)    :: day    !< Time in days
  type(ocean_grid_type),        intent(inout) :: G      !< Grid structure
  type(unit_scale_type),        intent(in)    :: US     !< A dimensional unit scaling type
  type(idealized_hurricane_CS), pointer       :: CS     !< Container for idealized hurricane parameters
 
  ! Local variables
  integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq
  integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
 
  real :: TX,TY       !< wind stress
  real :: Uocn, Vocn  !< Surface ocean velocity components
  real :: LAT, LON    !< Grid location
  real :: YY, XX      !< storm relative position
  real :: XC, YC      !< Storm center location
  real :: f           !< Coriolis
  real :: fbench      !< The benchmark 'f' value
  real :: fbench_fac  !< A factor that is set to 0 to use the
                      !!  benchmark 'f' value
  real :: rel_tau_fac !< A factor that is set to 0 to disable
                      !!  current relative stress calculation
 
  ! Bounds for loops and memory allocation
  is = g%isc    ; ie = g%iec    ; js = g%jsc    ; je = g%jec
  isq = g%IscB  ; ieq = g%IecB  ; jsq = g%JscB  ; jeq = g%JecB
  isd = g%isd   ; ied = g%ied   ; jsd = g%jsd   ; jed = g%jed
  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
 
  ! Allocate the forcing arrays, if necessary.
  call allocate_mech_forcing(g, forces, stress=.true., ustar=.true.)
 
  if (cs%relative_tau) then
     rel_tau_fac = 1.
  else
     rel_tau_fac = 0. !Multiplied to 0 surface current
  endif
 
  !> Compute storm center location
  xc = cs%Hurr_cen_X0 + (time_type_to_real(day)*cs%hurr_translation_spd*&
       cos(cs%hurr_translation_dir))
  yc = cs%Hurr_cen_Y0 + (time_type_to_real(day)*cs%hurr_translation_spd*&
       sin(cs%hurr_translation_dir))
 
 
  if (cs%BR_Bench) then
     ! f reset to value used in generated wind for benchmark test
     fbench = 5.5659e-05
     fbench_fac = 0.0
  else
     fbench = 0.0
     fbench_fac = 1.0
  endif
 
  !> Computes taux
  do j=js,je
    do i=is-1,ieq
      uocn = state%u(i,j)*rel_tau_fac
      vocn = 0.25*(state%v(i,j)+state%v(i+1,j-1)&
             +state%v(i+1,j)+state%v(i,j-1))*rel_tau_fac
      f = abs(0.5*us%s_to_T*(g%CoriolisBu(i,j)+g%CoriolisBu(i,j-1)))*fbench_fac + fbench
      ! Calculate position as a function of time.
      if (cs%SCM_mode) then
        yy = yc + cs%dy_from_center
        xx = xc
      else
        lat = g%geoLatCu(i,j)*1000. ! Convert Lat from km to m.
        lon = g%geoLonCu(i,j)*1000. ! Convert Lon from km to m.
        yy = lat - yc
        xx = lon - xc
      endif
      call idealized_hurricane_wind_profile(cs,f,yy,xx,uocn,vocn,tx,ty)
      forces%taux(i,j) = g%mask2dCu(i,j) * us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z * tx
    enddo
  enddo
  !> Computes tauy
  do j=js-1,jeq
    do i=is,ie
      uocn = 0.25*(state%u(i,j)+state%u(i-1,j+1)&
            +state%u(i-1,j)+state%u(i,j+1))*rel_tau_fac
      vocn = state%v(i,j)*rel_tau_fac
      f = abs(0.5*us%s_to_T*(g%CoriolisBu(i-1,j)+g%CoriolisBu(i,j)))*fbench_fac + fbench
      ! Calculate position as a function of time.
      if (cs%SCM_mode) then
        yy = yc + cs%dy_from_center
        xx = xc
      else
        lat = g%geoLatCv(i,j)*1000. ! Convert Lat from km to m.
        lon = g%geoLonCv(i,j)*1000. ! Convert Lon from km to m.
        yy = lat - yc
        xx = lon - xc
      endif
      call idealized_hurricane_wind_profile(cs, f, yy, xx, uocn, vocn, tx, ty)
      forces%tauy(i,j) = g%mask2dCv(i,j) * us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z * ty
    enddo
  enddo
 
  !> Get Ustar
  do j=js,je
    do i=is,ie
      !  This expression can be changed if desired, but need not be.
      forces%ustar(i,j) = g%mask2dT(i,j) * sqrt(us%L_to_Z * (cs%gustiness/cs%Rho0 + &
              sqrt(0.5*(forces%taux(i-1,j)**2 + forces%taux(i,j)**2) + &
                   0.5*(forces%tauy(i,j-1)**2 + forces%tauy(i,j)**2))/cs%Rho0))
    enddo
  enddo
 
  return

References mom_forcing_type::allocate_mech_forcing(), and idealized_hurricane_wind_profile().

Here is the call graph for this function:

idealized_hurricane_wind_init()

subroutine, public idealized_hurricane::idealized_hurricane_wind_init	(	type(time_type), intent(in)	Time,
		type(ocean_grid_type), intent(in)	G,
		type(unit_scale_type), intent(in)	US,
		type(param_file_type), intent(in)	param_file,
		type(idealized_hurricane_cs), pointer	CS
	)

Initializes wind profile for the SCM idealized hurricane example.

Parameters

[in]	time	Model time
[in]	g	Grid structure
[in]	us	A dimensional unit scaling type
[in]	param_file	Input parameter structure
	cs	Parameter container for this module

Definition at line 94 of file Idealized_Hurricane.F90.

  type(time_type),               intent(in) :: Time   !< Model time
  type(ocean_grid_type),         intent(in) :: G      !< Grid structure
  type(unit_scale_type),         intent(in) :: US     !< A dimensional unit scaling type
  type(param_file_type),         intent(in) :: param_file !< Input parameter structure
  type(idealized_hurricane_CS),  pointer    :: CS     !< Parameter container for this module
 
  real :: DP, C
 
! This include declares and sets the variable "version".
#include "version_variable.h"
 
  if (associated(cs)) then
    call mom_error(fatal, "idealized_hurricane_wind_init called "// &
                          "with an associated control structure.")
    return
  endif
 
  allocate(cs)
 
  cs%pi = 4.0*atan(1.0)
  cs%Deg2Rad = cs%pi/180.
 
  ! Read all relevant parameters and write them to the model log.
  call log_version(param_file, mdl, version, "")
 
  ! Parameters for computing a wind profile
  call get_param(param_file, mdl, "IDL_HURR_RHO_AIR", cs%rho_a, &
                 "Air density used to compute the idealized hurricane "//&
                 "wind profile.", units='kg/m3', default=1.2)
  call get_param(param_file, mdl, "IDL_HURR_AMBIENT_PRESSURE", &
                 cs%pressure_ambient, "Ambient pressure used in the "//&
                 "idealized hurricane wind profile.", units='Pa', &
                 default=101200.)
  call get_param(param_file, mdl, "IDL_HURR_CENTRAL_PRESSURE", &
                 cs%pressure_central, "Central pressure used in the "//&
                 "idealized hurricane wind profile.", units='Pa', &
                 default=96800.)
  call get_param(param_file, mdl, "IDL_HURR_RAD_MAX_WIND", &
                 cs%rad_max_wind, "Radius of maximum winds used in the "//&
                 "idealized hurricane wind profile.", units='m', &
                 default=50.e3)
  call get_param(param_file, mdl, "IDL_HURR_MAX_WIND", cs%max_windspeed, &
                 "Maximum wind speed used in the idealized hurricane"// &
                 "wind profile.", units='m/s', default=65.)
  call get_param(param_file, mdl, "IDL_HURR_TRAN_SPEED", cs%hurr_translation_spd, &
                 "Translation speed of hurricane used in the idealized "//&
                 "hurricane wind profile.", units='m/s', default=5.0)
  call get_param(param_file, mdl, "IDL_HURR_TRAN_DIR", cs%hurr_translation_dir, &
                 "Translation direction (towards) of hurricane used in the "//&
                 "idealized hurricane wind profile.", units='degrees', &
                 default=180.0)
  cs%hurr_translation_dir = cs%hurr_translation_dir * cs%Deg2Rad
  call get_param(param_file, mdl, "IDL_HURR_X0", cs%Hurr_cen_X0, &
                 "Idealized Hurricane initial X position", &
                 units='m', default=0.)
  call get_param(param_file, mdl, "IDL_HURR_Y0", cs%Hurr_cen_Y0, &
                 "Idealized Hurricane initial Y position", &
                 units='m', default=0.)
  call get_param(param_file, mdl, "IDL_HURR_TAU_CURR_REL", cs%relative_tau, &
                 "Current relative stress switch "//&
                 "used in the idealized hurricane wind profile.", &
                 units='', default=.false.)
 
  ! Parameters for SCM mode
  call get_param(param_file, mdl, "IDL_HURR_SCM_BR_BENCH", cs%BR_BENCH, &
                 "Single column mode benchmark case switch, which is "// &
                 "invoking a modification (bug) in the wind profile meant to "//&
                 "reproduce a previous implementation.", units='', default=.false.)
  call get_param(param_file, mdl, "IDL_HURR_SCM", cs%SCM_MODE, &
                 "Single Column mode switch "//&
                 "used in the SCM idealized hurricane wind profile.", &
                 units='', default=.false.)
  call get_param(param_file, mdl, "IDL_HURR_SCM_LOCY", cs%DY_from_center, &
                 "Y distance of station used in the SCM idealized hurricane "//&
                 "wind profile.", units='m', default=50.e3)
 
  ! The following parameters are model run-time parameters which are used
  ! and logged elsewhere and so should not be logged here. The default
  ! value should be consistent with the rest of the model.
  call get_param(param_file, mdl, "RHO_0", cs%Rho0, &
                 "The mean ocean density used with BOUSSINESQ true to "//&
                 "calculate accelerations and the mass for conservation "//&
                 "properties, or with BOUSSINSEQ false to convert some "//&
                 "parameters from vertical units of m to kg m-2.", &
                 units="kg m-3", default=1035.0, scale=us%kg_m3_to_R, do_not_log=.true.)
  call get_param(param_file, mdl, "GUST_CONST", cs%gustiness, &
                 "The background gustiness in the winds.", units="Pa", &
                 default=0.0, scale=us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z, do_not_log=.true.)
 
 
  if (cs%BR_BENCH) then
    cs%rho_a = 1.2
  endif
  dp = cs%pressure_ambient - cs%pressure_central
  c = cs%max_windspeed / sqrt( dp )
  cs%Holland_B = c**2 * cs%rho_a * exp(1.0)
  cs%Holland_A = (cs%rad_max_wind)**cs%Holland_B
  cs%Holland_AxBxDP = cs%Holland_A*cs%Holland_B*dp
 

References mdl, and mom_error_handler::mom_error().

Here is the call graph for this function:

idealized_hurricane_wind_profile()

subroutine idealized_hurricane::idealized_hurricane_wind_profile ( type(idealized_hurricane_cs), pointer  CS, real, intent(in)  absf, real, intent(in)  YY, real, intent(in)  XX, real, intent(in)  UOCN, real, intent(in)  VOCN, real, intent(out)  Tx, real, intent(out)  Ty  )

private

Calculate the wind speed at a location as a function of time.

Parameters

	cs	Container for SCM parameters
[in]	absf	Input coriolis magnitude
[in]	yy	Location in m relative to center y
[in]	xx	Location in m relative to center x
[in]	uocn	X surface current
[in]	vocn	Y surface current
[out]	tx	X stress
[out]	ty	Y stress

Definition at line 309 of file Idealized_Hurricane.F90.

  ! Author: Brandon Reichl
  ! Date: Nov-20-2014
  !       Aug-14-2018 Generalized for non-SCM configuration
 
  ! Input parameters
  type(idealized_hurricane_CS), &
        pointer     :: CS   !< Container for SCM parameters
  real, intent(in)  :: absf !<Input coriolis magnitude
  real, intent(in)  :: YY   !< Location in m relative to center y
  real, intent(in)  :: XX   !< Location in m relative to center x
  real, intent(in)  :: UOCN !< X surface current
  real, intent(in)  :: VOCN !< Y surface current
  real, intent(out) :: Tx   !< X stress
  real, intent(out) :: Ty   !< Y stress
 
  ! Local variables
 
  ! Wind profile terms
  real :: U10
  real :: radius
  real :: radius10
  real :: radius_km
  real :: radiusB
  real :: fcor
  real :: du10
  real :: du
  real :: dv
  real :: CD
 
  !Wind angle variables
  real :: Alph !< The resulting inflow angle (positive outward)
  real :: Rstr
  real :: A0
  real :: A1
  real :: P1
  real :: Adir
  real :: V_TS
  real :: U_TS
 
  ! Implementing Holland (1980) parameteric wind profile
 
  radius = sqrt(xx**2 + yy**2)
 
  !/ BGR
  ! rkm - r converted to km for Holland prof.
  !       used in km due to error, correct implementation should
  !       not need rkm, but to match winds w/ experiment this must
  !       be maintained.  Causes winds far from storm center to be a
  !       couple of m/s higher than the correct Holland prof.
  if (cs%BR_Bench) then
    radius_km = radius/1000.
  else
    ! if not comparing to benchmark, then use correct Holland prof.
    radius_km = radius
  endif
  radiusb = (radius)**cs%Holland_B
 
  !/
  ! Calculate U10 in the interior (inside of 10x radius of maximum wind),
  ! while adjusting U10 to 0 outside of 12x radius of maximum wind.
  if ( (radius/cs%rad_max_wind .gt. 0.001) .and. &
       (radius/cs%rad_max_wind .lt. 10.) ) then
    u10 = sqrt(cs%Holland_AxBxDP*exp(-cs%Holland_A/radiusb)/(cs%rho_A*radiusb)&
                +0.25*(radius_km*absf)**2) - 0.5*radius_km*absf
  elseif ( (radius/cs%rad_max_wind .gt. 10.) .and. &
           (radius/cs%rad_max_wind .lt. 15.) ) then
 
    radius10 = cs%rad_max_wind*10.
 
    if (cs%BR_Bench) then
      radius_km = radius10/1000.
    else
      radius_km = radius10
    endif
    radiusb=radius10**cs%Holland_B
 
    u10 = (sqrt(cs%Holland_AxBxDp*exp(-cs%Holland_A/radiusb)/(cs%rho_A*radiusb)&
                  +0.25*(radius_km*absf)**2)-0.5*radius_km*absf) &
           * (15.-radius/cs%rad_max_wind)/5.
  else
    u10 = 0.
  endif
  adir = atan2(yy,xx)
  !\
 
  ! Wind angle model following Zhang and Ulhorn (2012)
  ! ALPH is inflow angle positive outward.
  rstr = min(10.,radius / cs%rad_max_wind)
  a0 = -0.9*rstr - 0.09*cs%max_windspeed - 14.33
  a1 = -a0*(0.04*rstr + 0.05*cs%Hurr_translation_spd + 0.14)
  p1 = (6.88*rstr - 9.60*cs%Hurr_translation_spd + 85.31) * cs%Deg2Rad
  alph = a0 - a1*cos(cs%hurr_translation_dir-adir-p1)
  if ( (radius/cs%rad_max_wind.gt.10.) .and.&
       (radius/cs%rad_max_wind.lt.15.) ) then
     alph = alph*(15.0-radius/cs%rad_max_wind)/5.
  elseif (radius/cs%rad_max_wind.gt.15.) then
     alph = 0.0
  endif
  alph = alph * cs%Deg2Rad
 
  ! Calculate translation speed components
  u_ts = cs%hurr_translation_spd/2.*cos(cs%hurr_translation_dir)
  v_ts = cs%hurr_translation_spd/2.*sin(cs%hurr_translation_dir)
 
  ! Set output (relative) winds
  du = u10*sin(adir-cs%Pi-alph) - uocn + u_ts
  dv = u10*cos(adir-alph) - vocn + v_ts
 
  !  Use a simple drag coefficient as a function of U10 (from Sullivan et al., 2010)
  du10 = sqrt(du**2+dv**2)
  if (du10.lt.11.) then
     cd = 1.2e-3
  elseif (du10.lt.20.0) then
     cd = (0.49 + 0.065*u10)*1.e-3
  else
     cd = 1.8e-3
  endif
 
  ! Compute stress vector
  tx = cs%rho_A * cd * sqrt(du**2 + dv**2) * du
  ty = cs%rho_A * cd * sqrt(du**2 + dv**2) * dv
 

Referenced by idealized_hurricane_wind_forcing().

Here is the caller graph for this function:

scm_idealized_hurricane_wind_forcing()

subroutine, public idealized_hurricane::scm_idealized_hurricane_wind_forcing	(	type(surface), intent(in)	state,
		type(mech_forcing), intent(inout)	forces,
		type(time_type), intent(in)	day,
		type(ocean_grid_type), intent(inout)	G,
		type(unit_scale_type), intent(in)	US,
		type(idealized_hurricane_cs), pointer	CS
	)

This subroutine is primarily needed as a legacy for reproducing answers. It is included as an additional subroutine rather than padded into the previous routine with flags to ease its eventual removal. Its functionality is replaced with the new routines and it can be deleted when answer changes are acceptable.

Parameters

[in]	state	Surface state structure
[in,out]	forces	A structure with the driving mechanical forces
[in]	day	Time in days
[in,out]	g	Grid structure
[in]	us	A dimensional unit scaling type
	cs	Container for SCM parameters

Definition at line 438 of file Idealized_Hurricane.F90.

  type(surface),                intent(in)    :: state  !< Surface state structure
  type(mech_forcing),           intent(inout) :: forces !< A structure with the driving mechanical forces
  type(time_type),              intent(in)    :: day    !< Time in days
  type(ocean_grid_type),        intent(inout) :: G      !< Grid structure
  type(unit_scale_type),        intent(in)    :: US     !< A dimensional unit scaling type
  type(idealized_hurricane_CS), pointer       :: CS     !< Container for SCM parameters
  ! Local variables
  integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq
  integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
  real :: pie, Deg2Rad
  real :: U10, A, B, C, r, f, du10, rkm ! For wind profile expression
  real :: xx, t0 !for location
  real :: dp, rB
  real :: Cd ! Air-sea drag coefficient
  real :: Uocn, Vocn ! Surface ocean velocity components
  real :: dU, dV ! Air-sea differential motion
  !Wind angle variables
  real :: Alph,Rstr, A0, A1, P1, Adir, transdir, V_TS, U_TS
  logical :: BR_Bench
  ! Bounds for loops and memory allocation
  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
 
  ! Allocate the forcing arrays, if necessary.
 
  call allocate_mech_forcing(g, forces, stress=.true., ustar=.true.)
  pie = 4.0*atan(1.0) ; deg2rad = pie/180.
  !/ BR
  ! Implementing Holland (1980) parameteric wind profile
  !------------------------------------------------------|
  br_bench = .true.   !true if comparing to LES runs     |
  t0 = 129600.        !TC 'eye' crosses (0,0) at 36 hours|
  transdir = pie      !translation direction (-x)        |
  !------------------------------------------------------|
  dp = cs%pressure_ambient - cs%pressure_central
  c = cs%max_windspeed / sqrt( dp )
  b = c**2 * cs%rho_a * exp(1.0)
  if (br_bench) then
     ! rho_a reset to value used in generated wind for benchmark test
     b = c**2 * 1.2 * exp(1.0)
  endif
  a = (cs%rad_max_wind/1000.)**b
  f = us%s_to_T*g%CoriolisBu(is,js) ! f=f(x,y) but in the SCM is constant
  if (br_bench) then
     ! f reset to value used in generated wind for benchmark test
     f = 5.5659e-05  !### A constant value in s-1.
  endif
  !/ BR
  ! Calculate x position as a function of time.
  xx = ( t0 - time_type_to_real(day)) * cs%hurr_translation_spd * cos(transdir)
  r = sqrt(xx**2 + cs%DY_from_center**2)
  !/ BR
  ! rkm - r converted to km for Holland prof.
  !       used in km due to error, correct implementation should
  !       not need rkm, but to match winds w/ experiment this must
  !       be maintained.  Causes winds far from storm center to be a
  !       couple of m/s higher than the correct Holland prof.
  if (br_bench) then
     rkm = r/1000.
     rb = (rkm)**b
  else
     ! if not comparing to benchmark, then use correct Holland prof.
     rkm = r
     rb = r**b
  endif
  !/ BR
  ! Calculate U10 in the interior (inside of 10x radius of maximum wind),
  ! while adjusting U10 to 0 outside of 12x radius of maximum wind.
  ! Note that rho_a is set to 1.2 following generated wind for experiment
  if (r/cs%rad_max_wind > 0.001 .AND. r/cs%rad_max_wind < 10.) then
     u10 = sqrt( a*b*dp*exp(-a/rb)/(1.2*rb) + 0.25*(rkm*f)**2 ) - 0.5*rkm*f
  elseif (r/cs%rad_max_wind > 10. .AND. r/cs%rad_max_wind < 12.) then
     r=cs%rad_max_wind*10.
     if (br_bench) then
        rkm = r/1000.
        rb=rkm**b
     else
        rkm = r
        rb = r**b
     endif
     u10 = ( sqrt( a*b*dp*exp(-a/rb)/(1.2*rb) + 0.25*(rkm*f)**2 ) - 0.5*rkm*f) &
           * (12. - r/cs%rad_max_wind)/2.
  else
     u10 = 0.
  endif
  adir = atan2(cs%DY_from_center,xx)
 
  !/ BR
  ! Wind angle model following Zhang and Ulhorn (2012)
  ! ALPH is inflow angle positive outward.
  rstr = min(10.,r / cs%rad_max_wind)
  a0 = -0.9*rstr -0.09*cs%max_windspeed - 14.33
  a1 = -a0 *(0.04*rstr +0.05*cs%hurr_translation_spd+0.14)
  p1 = (6.88*rstr -9.60*cs%hurr_translation_spd+85.31)*pie/180.
  alph = a0 - a1*cos( (transdir - adir ) - p1)
  if (r/cs%rad_max_wind > 10. .AND. r/cs%rad_max_wind < 12.) then
     alph = alph* (12. - r/cs%rad_max_wind)/2.
  elseif (r/cs%rad_max_wind > 12.) then
     alph = 0.0
  endif
  alph = alph * deg2rad
 !/BR
  ! Prepare for wind calculation
  ! X_TS is component of translation speed added to wind vector
  ! due to background steering wind.
  u_ts = cs%hurr_translation_spd/2.*cos(transdir)
  v_ts = cs%hurr_translation_spd/2.*sin(transdir)
 
  ! Set the surface wind stresses, in [Pa]. A positive taux
  ! accelerates the ocean to the (pseudo-)east.
  !   The i-loop extends to is-1 so that taux can be used later in the
  ! calculation of ustar - otherwise the lower bound would be Isq.
  do j=js,je ; do i=is-1,ieq
    !/BR
    ! Turn off surface current for stress calculation to be
    ! consistent with test case.
    uocn = 0.!state%u(I,j)
    vocn = 0.!0.25*( (state%v(i,J) + state%v(i+1,J-1)) &
             !    +(state%v(i+1,J) + state%v(i,J-1)) )
    !/BR
    ! Wind vector calculated from location/direction (sin/cos flipped b/c
    ! cyclonic wind is 90 deg. phase shifted from position angle).
    du = u10*sin(adir-pie-alph) - uocn + u_ts
    dv = u10*cos(adir-alph) - vocn + v_ts
    !/----------------------------------------------------|
    !BR
    !  Add a simple drag coefficient as a function of U10 |
    !/----------------------------------------------------|
    du10=sqrt(du**2+dv**2)
    if (du10 < 11.) then
       cd = 1.2e-3
    elseif (du10 < 20.) then
       cd = (0.49 + 0.065 * u10 )*0.001
    else
       cd = 0.0018
    endif
    forces%taux(i,j) = cs%rho_a * us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z * &
                       g%mask2dCu(i,j) * cd*sqrt(du**2+dv**2)*du
  enddo ; enddo
  !/BR
  ! See notes above
  do j=js-1,jeq ; do i=is,ie
    uocn = 0.!0.25*( (state%u(I,j) + state%u(I-1,j+1)) &
             !    +(state%u(I-1,j) + state%u(I,j+1)) )
    vocn = 0.!state%v(i,J)
    du = u10*sin(adir-pie-alph) - uocn + u_ts
    dv = u10*cos(adir-alph) - vocn + v_ts
    du10=sqrt(du**2+dv**2)
    if (du10 < 11.) then
       cd = 1.2e-3
    elseif (du10 < 20.) then
       cd = (0.49 + 0.065 * u10 )*0.001
    else
       cd = 0.0018
    endif
    forces%tauy(i,j) = cs%rho_a * us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z * &
                       g%mask2dCv(i,j) * cd*du10*dv
  enddo ; enddo
  ! Set the surface friction velocity [m s-1]. ustar is always positive.
  do j=js,je ; do i=is,ie
    !  This expression can be changed if desired, but need not be.
    forces%ustar(i,j) = g%mask2dT(i,j) * sqrt(us%L_to_Z * (cs%gustiness/cs%Rho0 + &
            sqrt(0.5*(forces%taux(i-1,j)**2 + forces%taux(i,j)**2) + &
                 0.5*(forces%tauy(i,j-1)**2 + forces%tauy(i,j)**2))/cs%Rho0))
  enddo ; enddo
 

References mom_forcing_type::allocate_mech_forcing().

Here is the call graph for this function:

Variable Documentation

mdl

character(len=40) idealized_hurricane::mdl = "idealized_hurricane"

This module's name.

Definition at line 88 of file Idealized_Hurricane.F90.

character(len=40)  :: mdl = "idealized_hurricane" !< This module's name.

Referenced by idealized_hurricane_wind_init().