[3816] | 1 | module bulk_flux_m |
---|
| 2 | |
---|
| 3 | implicit none |
---|
| 4 | |
---|
| 5 | contains |
---|
| 6 | |
---|
| 7 | subroutine bulk_flux(tkt, tks, taur, dter, dser, t_int, s_int, ds_ns, dt_ns, & |
---|
| 8 | u, t_ocean_1, s1, rain, hf, hlb, rnl, tau, rhoa, xlv, rf, dtime, rns) |
---|
| 9 | |
---|
| 10 | use config_ocean_skin_m, only: jwarm, jcool, rain_effect |
---|
| 11 | use Microlayer_m, only: Microlayer |
---|
| 12 | use mom_flux_rain_m, only: mom_flux_rain |
---|
| 13 | use Near_Surface_m, only: Near_Surface, depth |
---|
| 14 | use therm_expans_m, only: therm_expans |
---|
| 15 | |
---|
| 16 | real, intent(out):: tkt(:) |
---|
| 17 | ! thickness of cool skin (microlayer), in m |
---|
| 18 | |
---|
| 19 | real, intent(out):: tks(:) |
---|
| 20 | ! thickness of mass diffusion layer (microlayer), in m |
---|
| 21 | |
---|
| 22 | real, intent(out):: taur(:) ! momentum flux due to rain, in Pa |
---|
| 23 | |
---|
| 24 | real, intent(out):: dter(:) |
---|
| 25 | ! Temperature variation in the diffusive microlayer, that is |
---|
| 26 | ! ocean-air interface temperature minus subskin temperature. In K. |
---|
| 27 | |
---|
| 28 | real, intent(out):: dser(:) |
---|
| 29 | ! Salinity variation in the diffusive microlayer, that is ocean-air |
---|
| 30 | ! interface salinity minus subskin salinity. In ppt. |
---|
| 31 | |
---|
| 32 | real, intent(out):: t_int(:) ! interface temperature, in K |
---|
| 33 | real, intent(out):: s_int(:) ! interface salinity, in ppt |
---|
| 34 | |
---|
| 35 | real, intent(inout):: ds_ns(:) |
---|
| 36 | ! "delta salinity near surface". Salinity variation in the |
---|
| 37 | ! near-surface turbulent layer. That is subskin salinity minus |
---|
| 38 | ! foundation salinity. In ppt. |
---|
| 39 | |
---|
| 40 | real, intent(inout):: dt_ns(:) |
---|
| 41 | ! "delta temperature near surface". Temperature variation in the |
---|
| 42 | ! near-surface turbulent layer. That is subskin temperature minus |
---|
| 43 | ! foundation temperature. (Can be negative.) In K. |
---|
| 44 | |
---|
| 45 | real, intent(in):: u(:) |
---|
| 46 | ! Wind speed relative to the sea surface, i. e. taking current |
---|
| 47 | ! vector into account. In m s-1. |
---|
| 48 | |
---|
| 49 | real, intent(in):: t_ocean_1(:) ! input sea temperature, at depth_1, in K |
---|
| 50 | real, intent(in):: S1(:) ! salinity at depth_1, in ppt |
---|
| 51 | |
---|
| 52 | real, intent(in):: rain(:) |
---|
| 53 | ! rain mass flux, averaged on a timestep, in kg m-2 s-1 |
---|
| 54 | |
---|
| 55 | real, intent(in):: hf(:) |
---|
| 56 | ! turbulent part of sensible heat flux, positive upward, in W m-2 |
---|
| 57 | |
---|
| 58 | real, intent(in):: hlb(:) |
---|
| 59 | ! latent heat flux at the surface, positive upward (W m-2) |
---|
| 60 | |
---|
| 61 | real, intent(in):: rnl(:) |
---|
| 62 | ! net longwave radiation, positive upward, in W m-2 |
---|
| 63 | |
---|
| 64 | real, intent(in):: tau(:) |
---|
| 65 | ! wind stress at the surface, turbulent part only, in Pa |
---|
| 66 | |
---|
| 67 | real, intent(in):: rhoa(:) ! density of moist air (kg / m3) |
---|
| 68 | real, intent(in):: xlv(:) ! latent heat of evaporation (J / kg) |
---|
| 69 | |
---|
| 70 | real, intent(in):: rf(:) |
---|
| 71 | ! sensible heat flux at the surface due to rainfall, in W m-2, |
---|
| 72 | ! positive upward |
---|
| 73 | |
---|
| 74 | real, intent(in):: dtime ! time step, in s |
---|
| 75 | real, intent(in):: rns(:) ! net downward shortwave radiation, in W m-2 |
---|
| 76 | |
---|
| 77 | ! Local: |
---|
| 78 | |
---|
| 79 | real al(size(t_ocean_1)) ! water thermal expansion coefficient (in K-1) |
---|
| 80 | real dels(size(t_ocean_1)), null_array(size(t_ocean_1)) |
---|
| 81 | integer iter |
---|
| 82 | real t_subskin(size(t_ocean_1)) ! subskin temperature, in K |
---|
| 83 | real s_subskin(size(t_ocean_1)) ! subskin salinity, in ppt |
---|
| 84 | |
---|
| 85 | real, parameter:: fxp = 1. - (0.28 * 0.014 & |
---|
| 86 | + 0.27 * 0.357 * (1. - exp(- depth / 0.357)) & |
---|
| 87 | + .45 * 12.82 * (1.- exp(- depth / 12.82))) / depth |
---|
| 88 | ! underflow ! fxp = 1. - (0.28 * 0.014 * (1. - exp(- depth / 0.014)) & |
---|
| 89 | ! Soloviev solar absorption profile |
---|
| 90 | ! H. Bellenger 2016 |
---|
| 91 | |
---|
| 92 | real tau_with_min(size(t_ocean_1)) |
---|
| 93 | ! modified wind stress, avoiding very low values |
---|
| 94 | |
---|
| 95 | real, parameter:: tau_0 = 1e-3 ! in N m-2 |
---|
| 96 | |
---|
| 97 | !------------------------------------------------------------------- |
---|
| 98 | |
---|
| 99 | if (rain_effect) then |
---|
| 100 | taur = mom_flux_rain(u, rain) |
---|
| 101 | else |
---|
| 102 | if (jwarm .or. jcool) null_array = 0. |
---|
| 103 | taur = 0. |
---|
| 104 | end if |
---|
| 105 | |
---|
| 106 | if (jwarm .or. jcool) tau_with_min = tau + tau_0 * (1. - exp(- tau_0 / tau)) |
---|
| 107 | |
---|
| 108 | if (Jwarm) then |
---|
| 109 | if (rain_effect) then |
---|
| 110 | call Near_Surface(al, t_subskin, s_subskin, ds_ns, dt_ns, & |
---|
| 111 | tau_with_min, taur, hlb, rhoa, xlv, dtime, t_ocean_1, s1, rain, & |
---|
| 112 | q_pwp = fxp * rns - (hf + hlb + rnl + rf)) |
---|
| 113 | else |
---|
| 114 | call Near_Surface(al, t_subskin, s_subskin, ds_ns, dt_ns, & |
---|
| 115 | tau_with_min, taur, hlb, rhoa, xlv, dtime, t_ocean_1, s1, & |
---|
| 116 | rain = null_array, q_pwp = fxp * rns - (hf + hlb + rnl)) |
---|
| 117 | end if |
---|
| 118 | else |
---|
| 119 | if (Jcool) al = therm_expans(t_ocean_1) |
---|
| 120 | t_subskin = t_ocean_1 |
---|
| 121 | s_subskin = s1 |
---|
| 122 | end if |
---|
| 123 | |
---|
| 124 | if (Jcool) then |
---|
| 125 | ! First guess: |
---|
| 126 | tkt = 0.001 |
---|
| 127 | tks = 5e-4 |
---|
| 128 | |
---|
| 129 | do iter = 1, 3 |
---|
| 130 | ! Cool skin |
---|
| 131 | dels = rns * (0.065 + 11. * tkt - 6.6e-5 / tkt & |
---|
| 132 | * (1. - exp(- tkt / 8e-4))) ! equation 16 Ohlmann |
---|
| 133 | if (rain_effect) then |
---|
| 134 | call Microlayer(dter, dser, tkt, tks, hlb, tau_with_min, & |
---|
| 135 | s_subskin, al, xlv, taur, rf, rain, & |
---|
| 136 | qcol = rnl + hf + hlb - dels) |
---|
| 137 | else |
---|
| 138 | call Microlayer(dter, dser, tkt, tks, hlb, tau_with_min, & |
---|
| 139 | s_subskin, al, xlv, taur, rf = null_array, & |
---|
| 140 | rain = null_array, qcol = rnl + hf + hlb - dels) |
---|
| 141 | end if |
---|
| 142 | end do |
---|
| 143 | else |
---|
| 144 | tkt = 0. |
---|
| 145 | tks = 0. |
---|
| 146 | dter = 0. |
---|
| 147 | dser = 0. |
---|
| 148 | end if |
---|
| 149 | |
---|
| 150 | t_int = t_subskin + dter |
---|
| 151 | s_int = s_subskin + dser |
---|
| 152 | |
---|
| 153 | end subroutine bulk_flux |
---|
| 154 | |
---|
| 155 | end module bulk_flux_m |
---|