| 1 |
jmc |
1.1 |
module radiation_mod |
| 2 |
|
|
|
| 3 |
|
|
! ================================================================================== |
| 4 |
|
|
! ================================================================================== |
| 5 |
|
|
|
| 6 |
|
|
! use fms_mod, only: open_file, check_nml_error, & |
| 7 |
|
|
! mpp_pe, close_file |
| 8 |
|
|
|
| 9 |
|
|
use gcm_params_mod, only: gcm_LEN_MBUF, gcm_SQZ_R, gcm_stdMsgUnit |
| 10 |
|
|
use constants_mod, only: stefan, cp_air, grav, pstd_mks |
| 11 |
|
|
|
| 12 |
|
|
! use diag_manager_mod, only: register_diag_field, send_data |
| 13 |
|
|
|
| 14 |
|
|
! use time_manager_mod, only: time_type, & |
| 15 |
|
|
! operator(+), operator(-), operator(/=) |
| 16 |
|
|
|
| 17 |
|
|
!================================================================================== |
| 18 |
|
|
implicit none |
| 19 |
|
|
private |
| 20 |
|
|
!================================================================================== |
| 21 |
|
|
|
| 22 |
|
|
! version information |
| 23 |
|
|
|
| 24 |
jmc |
1.2 |
character(len=128) :: version='$Id: radiation_mod.F90,v 1.1 2012/09/11 03:53:05 jmc Exp $' |
| 25 |
jmc |
1.1 |
character(len=128) :: tag='homemade' |
| 26 |
|
|
|
| 27 |
|
|
!================================================================================== |
| 28 |
|
|
|
| 29 |
|
|
! public interfaces |
| 30 |
|
|
|
| 31 |
|
|
public :: radiation_init, radiation_down, radiation_up, radiation_end |
| 32 |
|
|
!================================================================================== |
| 33 |
|
|
|
| 34 |
|
|
! module variables |
| 35 |
|
|
logical :: initialized =.false. |
| 36 |
|
|
|
| 37 |
|
|
real :: solar_constant = 1360.0 |
| 38 |
|
|
real :: del_sol = 1.4 |
| 39 |
|
|
! modif omp: winter/summer hemisphere |
| 40 |
|
|
real :: del_sw = 0.0 |
| 41 |
|
|
real :: ir_tau_eq = 6.0 |
| 42 |
|
|
real :: ir_tau_pole = 1.5 |
| 43 |
|
|
real :: atm_abs = 0.0 |
| 44 |
|
|
real :: sw_diff = 0.0 |
| 45 |
|
|
real :: linear_tau = 0.1 |
| 46 |
|
|
real :: albedo_value = 0.06 |
| 47 |
|
|
real :: window = 0.0 ! spectral window transparent to LW |
| 48 |
|
|
real :: wv_exponent = 4.0 |
| 49 |
|
|
real :: solar_exponent = 4.0 |
| 50 |
|
|
|
| 51 |
|
|
!real, allocatable, dimension(:) :: ss, solar, tau_0, solar_tau_0, p2 |
| 52 |
|
|
!real, allocatable, dimension(:,:) :: b_surf |
| 53 |
|
|
!real, allocatable, dimension(:,:,:) :: b, tdt_rad, entrop_rad, tdt_sw |
| 54 |
|
|
!!real, allocatable, dimension(:,:,:) :: up, down, net, solar_down, flux_rad, flux_sw |
| 55 |
|
|
!!real, allocatable, dimension(:,:) :: dtrans |
| 56 |
|
|
!real, allocatable, dimension(:,:,:) :: up, net, flux_rad, flux_sw |
| 57 |
|
|
!real, allocatable, dimension(:,:) :: tau, solar_tau |
| 58 |
|
|
!real, allocatable, dimension(:,:) :: olr, swin |
| 59 |
|
|
!real, allocatable, dimension(:,:) :: albedo |
| 60 |
|
|
|
| 61 |
|
|
real, save :: pi, deg_to_rad , rad_to_deg |
| 62 |
|
|
|
| 63 |
|
|
namelist/radiation_nml/ solar_constant, del_sol, & |
| 64 |
|
|
ir_tau_eq, ir_tau_pole, atm_abs, sw_diff, & |
| 65 |
|
|
linear_tau, del_sw, albedo_value, window, wv_exponent, & |
| 66 |
|
|
solar_exponent |
| 67 |
|
|
|
| 68 |
|
|
!================================================================================== |
| 69 |
|
|
!-------------------- diagnostics fields ------------------------------- |
| 70 |
|
|
|
| 71 |
|
|
!integer :: id_olr, id_swdn_sfc, id_swdn_toa, id_lwdn_sfc, id_lwup_sfc, & |
| 72 |
|
|
! id_tdt_rad, id_flux_rad, id_flux_lw, id_flux_sw, id_entrop_rad, id_tdt_sw |
| 73 |
|
|
|
| 74 |
|
|
character(len=10), parameter :: mod_name = 'two_stream' |
| 75 |
|
|
|
| 76 |
|
|
real :: missing_value = -999. |
| 77 |
|
|
|
| 78 |
|
|
contains |
| 79 |
|
|
|
| 80 |
|
|
! ================================================================================== |
| 81 |
|
|
! ================================================================================== |
| 82 |
|
|
|
| 83 |
|
|
!subroutine radiation_init(is, ie, js, je, num_levels, axes, Time) |
| 84 |
|
|
subroutine radiation_init(is, ie, js, je, num_levels, nSx,nSy, axes, Time, myThid) |
| 85 |
|
|
|
| 86 |
|
|
!------------------------------------------------------------------------------------- |
| 87 |
|
|
integer, intent(in), dimension(4) :: axes |
| 88 |
|
|
!type(time_type), intent(in) :: Time |
| 89 |
|
|
real, intent(in) :: Time |
| 90 |
|
|
integer, intent(in) :: is, ie, js, je, num_levels |
| 91 |
|
|
integer, intent(in) :: nSx, nSy |
| 92 |
|
|
integer, intent(in) :: myThid |
| 93 |
|
|
!------------------------------------------------------------------------------------- |
| 94 |
|
|
integer, dimension(3) :: half = (/1,2,4/) |
| 95 |
|
|
!integer :: ierr, io |
| 96 |
|
|
integer :: iUnit |
| 97 |
|
|
CHARACTER*(gcm_LEN_MBUF) :: msgBuf |
| 98 |
|
|
!----------------------------------------------------------------------------------------- |
| 99 |
|
|
! read namelist and copy to logfile |
| 100 |
|
|
|
| 101 |
|
|
! _BARRIER |
| 102 |
|
|
! _BEGIN_MASTER(myThid) |
| 103 |
|
|
CALL BARRIER(myThid) |
| 104 |
|
|
IF ( myThid.EQ.1 ) THEN |
| 105 |
|
|
|
| 106 |
|
|
WRITE(msgBuf,'(A)') 'RADIATION_INIT: opening data.atm_gray' |
| 107 |
|
|
CALL PRINT_MESSAGE( msgBuf, gcm_stdMsgUnit, gcm_SQZ_R, myThid ) |
| 108 |
|
|
CALL OPEN_COPY_DATA_FILE( & |
| 109 |
|
|
'data.atm_gray', 'RADIATION_INIT', & |
| 110 |
|
|
iUnit, & |
| 111 |
|
|
myThid ) |
| 112 |
|
|
! Read parameters from open data file |
| 113 |
|
|
READ(UNIT=iUnit,NML=radiation_nml) |
| 114 |
|
|
WRITE(msgBuf,'(A)') & |
| 115 |
|
|
'RADIATION_INIT: finished reading data.atm_gray' |
| 116 |
|
|
CALL PRINT_MESSAGE( msgBuf, gcm_stdMsgUnit, gcm_SQZ_R, myThid ) |
| 117 |
|
|
! Close the open data file |
| 118 |
|
|
CLOSE(iUnit) |
| 119 |
|
|
|
| 120 |
|
|
!unit = open_file ('input.nml', action='read') |
| 121 |
|
|
!ierr=1 |
| 122 |
|
|
!do while (ierr /= 0) |
| 123 |
|
|
! read (unit, nml=radiation_nml, iostat=io, end=10) |
| 124 |
|
|
! ierr = check_nml_error (io, 'radiation_nml') |
| 125 |
|
|
!enddo |
| 126 |
|
|
!10 call close_file (unit) |
| 127 |
|
|
|
| 128 |
|
|
!unit = open_file ('logfile.out', action='append') |
| 129 |
|
|
!if ( mpp_pe() == 0 ) then |
| 130 |
|
|
! write (unit,'(/,80("="),/(a))') trim(version), trim(tag) |
| 131 |
|
|
! write (unit, nml=radiation_nml) |
| 132 |
|
|
!endif |
| 133 |
|
|
!call close_file (unit) |
| 134 |
|
|
|
| 135 |
|
|
pi = 4.0*atan(1.) |
| 136 |
|
|
deg_to_rad = 2.*pi/360. |
| 137 |
|
|
rad_to_deg = 360.0/2./pi |
| 138 |
|
|
|
| 139 |
|
|
initialized = .true. |
| 140 |
|
|
|
| 141 |
|
|
!allocate (b (ie-is+1, je-js+1, num_levels)) |
| 142 |
|
|
!allocate (tdt_rad (ie-is+1, je-js+1, num_levels)) |
| 143 |
|
|
!allocate (tdt_sw (ie-is+1, je-js+1, num_levels)) |
| 144 |
|
|
!allocate (entrop_rad (ie-is+1, je-js+1, num_levels)) |
| 145 |
|
|
|
| 146 |
|
|
!allocate (up (ie-is+1, je-js+1, num_levels+1)) |
| 147 |
|
|
!!allocate (down (ie-is+1, je-js+1, num_levels+1)) |
| 148 |
|
|
!allocate (net (ie-is+1, je-js+1, num_levels+1)) |
| 149 |
|
|
!!allocate (solar_down (ie-is+1, je-js+1, num_levels+1)) |
| 150 |
|
|
!allocate (flux_rad (ie-is+1, je-js+1, num_levels+1)) |
| 151 |
|
|
!allocate (flux_sw (ie-is+1, je-js+1, num_levels+1)) |
| 152 |
|
|
|
| 153 |
|
|
!allocate (b_surf (ie-is+1, je-js+1)) |
| 154 |
|
|
!allocate (olr (ie-is+1, je-js+1)) |
| 155 |
|
|
!allocate (swin (ie-is+1, je-js+1)) |
| 156 |
|
|
!allocate (albedo (ie-is+1, je-js+1)) |
| 157 |
|
|
|
| 158 |
|
|
!!allocate (dtrans (je-js+1, num_levels)) |
| 159 |
|
|
|
| 160 |
|
|
!allocate (tau (je-js+1, num_levels+1)) |
| 161 |
|
|
!allocate (solar_tau (je-js+1, num_levels+1)) |
| 162 |
|
|
|
| 163 |
|
|
!allocate (ss (je-js+1)) |
| 164 |
|
|
!allocate (solar (je-js+1)) |
| 165 |
|
|
!allocate (tau_0 (je-js+1)) |
| 166 |
|
|
!allocate (solar_tau_0 (je-js+1)) |
| 167 |
|
|
!allocate (p2 (je-js+1)) |
| 168 |
|
|
|
| 169 |
|
|
ENDIF |
| 170 |
|
|
CALL BARRIER(myThid) |
| 171 |
|
|
|
| 172 |
|
|
!----------------------------------------------------------------------- |
| 173 |
|
|
!------------ initialize diagnostic fields --------------- |
| 174 |
|
|
|
| 175 |
|
|
! id_olr = & |
| 176 |
|
|
! register_diag_field ( mod_name, 'olr', axes(1:2), Time, & |
| 177 |
|
|
! 'outgoing longwave radiation', & |
| 178 |
|
|
! 'watts/m2', missing_value=missing_value ) |
| 179 |
|
|
! id_swdn_sfc = & |
| 180 |
|
|
! register_diag_field ( mod_name, 'swdn_sfc', axes(1:2), Time, & |
| 181 |
|
|
! 'SW flux down at surface', & |
| 182 |
|
|
! 'watts/m2', missing_value=missing_value ) |
| 183 |
|
|
! id_swdn_toa = & |
| 184 |
|
|
! register_diag_field ( mod_name, 'swdn_toa', axes(1:2), Time, & |
| 185 |
|
|
! 'SW flux down at TOA', & |
| 186 |
|
|
! 'watts/m2', missing_value=missing_value ) |
| 187 |
|
|
! id_lwup_sfc = & |
| 188 |
|
|
! register_diag_field ( mod_name, 'lwup_sfc', axes(1:2), Time, & |
| 189 |
|
|
! 'LW flux up at surface', & |
| 190 |
|
|
! 'watts/m2', missing_value=missing_value ) |
| 191 |
|
|
|
| 192 |
|
|
! id_lwdn_sfc = & |
| 193 |
|
|
! register_diag_field ( mod_name, 'lwdn_sfc', axes(1:2), Time, & |
| 194 |
|
|
! 'LW flux down at surface', & |
| 195 |
|
|
! 'watts/m2', missing_value=missing_value ) |
| 196 |
|
|
|
| 197 |
|
|
! id_tdt_rad = & |
| 198 |
|
|
! register_diag_field ( mod_name, 'tdt_rad', axes(1:3), Time, & |
| 199 |
|
|
! 'Temperature tendency due to radiation', & |
| 200 |
|
|
! 'K/s', missing_value=missing_value ) |
| 201 |
|
|
|
| 202 |
|
|
! id_tdt_sw = & |
| 203 |
|
|
! register_diag_field ( mod_name, 'tdt_sw', axes(1:3), Time, & |
| 204 |
|
|
! 'Temperature tendency due to SW radiation', & |
| 205 |
|
|
! 'K/s', missing_value=missing_value ) |
| 206 |
|
|
|
| 207 |
|
|
! id_flux_rad = & |
| 208 |
|
|
! register_diag_field ( mod_name, 'flux_rad', axes(half), Time, & |
| 209 |
|
|
! 'Total radiative flux (positive up)', & |
| 210 |
|
|
! 'W/m^2', missing_value=missing_value ) |
| 211 |
|
|
! id_flux_lw = & |
| 212 |
|
|
! register_diag_field ( mod_name, 'flux_lw', axes(half), Time, & |
| 213 |
|
|
! 'Net longwave radiative flux (positive up)', & |
| 214 |
|
|
! 'W/m^2', missing_value=missing_value ) |
| 215 |
|
|
! id_flux_sw = & |
| 216 |
|
|
! register_diag_field ( mod_name, 'flux_sw', axes(half), Time, & |
| 217 |
|
|
! 'Net shortwave radiative flux (positive up)', & |
| 218 |
|
|
! 'W/m^2', missing_value=missing_value ) |
| 219 |
|
|
! id_entrop_rad = & |
| 220 |
|
|
! register_diag_field ( mod_name, 'entrop_rad', axes(1:3), Time, & |
| 221 |
|
|
! 'Entropy production by radiation', & |
| 222 |
|
|
! '1/s', missing_value=missing_value ) |
| 223 |
|
|
|
| 224 |
|
|
return |
| 225 |
|
|
end subroutine radiation_init |
| 226 |
|
|
|
| 227 |
|
|
! ================================================================================== |
| 228 |
|
|
|
| 229 |
|
|
subroutine radiation_down (is, js, Time_diag, lat, p_half, t, & |
| 230 |
|
|
! net_surf_sw_down, surf_lw_down) |
| 231 |
|
|
net_surf_sw_down, surf_lw_down, & |
| 232 |
|
|
albedo, dtrans, b, down, solar_down, & |
| 233 |
|
|
myThid ) |
| 234 |
|
|
|
| 235 |
|
|
! Begin the radiation calculation by computing downward fluxes. |
| 236 |
|
|
! This part of the calculation does not depend on the surface temperature. |
| 237 |
|
|
|
| 238 |
|
|
integer, intent(in) :: is, js |
| 239 |
|
|
!type(time_type), intent(in) :: Time_diag |
| 240 |
|
|
real, intent(in) :: Time_diag |
| 241 |
|
|
real, intent(in) , dimension(:,:) :: lat |
| 242 |
|
|
real, intent(out) , dimension(:,:) :: net_surf_sw_down |
| 243 |
|
|
real, intent(out) , dimension(:,:) :: surf_lw_down |
| 244 |
|
|
real, intent(in) , dimension(:,:,:) :: t, p_half |
| 245 |
|
|
real, intent(out), dimension(:,:) :: albedo |
| 246 |
|
|
!real, intent(out), dimension(:,:) :: dtrans |
| 247 |
|
|
real, intent(out), dimension(:,:,:) :: dtrans |
| 248 |
|
|
real, intent(out), dimension(:,:,:) :: b |
| 249 |
|
|
real, intent(out), dimension(:,:,:) :: down |
| 250 |
|
|
real, intent(out), dimension(:,:,:) :: solar_down |
| 251 |
|
|
integer, intent(in) :: myThid |
| 252 |
|
|
|
| 253 |
|
|
integer :: i, j, k, n |
| 254 |
|
|
integer :: im, jm |
| 255 |
|
|
|
| 256 |
|
|
logical :: used |
| 257 |
|
|
|
| 258 |
|
|
! ------------------------------------------------------------------------- |
| 259 |
|
|
real, allocatable, dimension(:,:) :: swin |
| 260 |
|
|
!real, allocatable, dimension(:) :: ss, solar, tau_0, solar_tau_0, p2 |
| 261 |
|
|
!real, allocatable, dimension(:,:) :: tau, solar_tau |
| 262 |
|
|
real, allocatable, dimension(:,:) :: ss, solar, tau_0, solar_tau_0, p2 |
| 263 |
|
|
real, allocatable, dimension(:,:,:) :: tau, solar_tau |
| 264 |
|
|
! ------------------------------------------------------------------------- |
| 265 |
|
|
|
| 266 |
|
|
n = size(t,3) |
| 267 |
|
|
im = size(t,1) |
| 268 |
|
|
jm = size(t,2) |
| 269 |
|
|
|
| 270 |
|
|
! ------------------------------------------------------------------------- |
| 271 |
|
|
allocate (swin (im, jm)) |
| 272 |
|
|
!allocate (tau (jm, n+1)) |
| 273 |
|
|
!allocate (solar_tau (jm, n+1)) |
| 274 |
|
|
!allocate (ss (jm)) |
| 275 |
|
|
!allocate (solar (jm)) |
| 276 |
|
|
!allocate (tau_0 (jm)) |
| 277 |
|
|
!allocate (solar_tau_0 (jm)) |
| 278 |
|
|
!allocate (p2 (jm)) |
| 279 |
|
|
allocate (ss (im, jm)) |
| 280 |
|
|
allocate (solar (im, jm)) |
| 281 |
|
|
allocate (tau_0 (im, jm)) |
| 282 |
|
|
allocate (solar_tau_0 (im, jm)) |
| 283 |
|
|
allocate (p2 (im, jm)) |
| 284 |
|
|
allocate (tau (im, jm, n+1)) |
| 285 |
|
|
allocate (solar_tau (im, jm, n+1)) |
| 286 |
|
|
! ------------------------------------------------------------------------- |
| 287 |
|
|
|
| 288 |
|
|
!ss = sin(lat(1,:)) |
| 289 |
|
|
ss = sin(lat) |
| 290 |
|
|
p2 = (1. - 3.*ss*ss)/4. |
| 291 |
|
|
|
| 292 |
|
|
solar = 0.25*solar_constant*(1.0 + del_sol*p2 + del_sw * ss) |
| 293 |
|
|
|
| 294 |
|
|
tau_0 = ir_tau_eq + (ir_tau_pole - ir_tau_eq)*ss*ss |
| 295 |
|
|
|
| 296 |
|
|
solar_tau_0 = (1.0 - sw_diff*ss*ss)*atm_abs |
| 297 |
|
|
|
| 298 |
|
|
! set a constant albedo for testing |
| 299 |
|
|
albedo(:,:) = albedo_value |
| 300 |
|
|
|
| 301 |
|
|
do k = 1, n+1 |
| 302 |
|
|
|
| 303 |
|
|
! tau(:,k) = tau_0(:) * (linear_tau * p_half(1,:,k)/p_half(1,:,n+1) + (1.0 - linear_tau) & |
| 304 |
|
|
! * (p_half(1,:,k)/p_half(1,:,n+1))**wv_exponent) |
| 305 |
|
|
! solar_tau(:,k) = solar_tau_0(:)*(p_half(1,:,k)/p_half(1,:,n+1))**solar_exponent |
| 306 |
|
|
|
| 307 |
|
|
tau(:,:,k) = tau_0(:,:) * ( & |
| 308 |
|
|
linear_tau * p_half(:,:,k)/p_half(:,:,n+1) & |
| 309 |
|
|
+ (1.0 - linear_tau)*(p_half(:,:,k)/p_half(:,:,n+1))**wv_exponent & |
| 310 |
|
|
) |
| 311 |
|
|
solar_tau(:,:,k) = solar_tau_0(:,:) & |
| 312 |
|
|
*(p_half(:,:,k)/p_half(:,:,n+1))**solar_exponent |
| 313 |
|
|
|
| 314 |
|
|
end do |
| 315 |
|
|
|
| 316 |
|
|
! no radiation from spectral window |
| 317 |
|
|
b = (1.0-window)*stefan*t*t*t*t |
| 318 |
|
|
|
| 319 |
|
|
do k = 1, n |
| 320 |
|
|
dtrans(:,:,k) = exp(-(tau(:,:,k+1)-tau(:,:,k))) |
| 321 |
|
|
end do |
| 322 |
|
|
|
| 323 |
|
|
down(:,:,1) = 0.0 |
| 324 |
|
|
do k = 1,n |
| 325 |
|
|
! do j =1, size(t,2) |
| 326 |
|
|
! down(:,j,k+1) = down(:,j,k)*dtrans(j,k) + b(:,j,k)*(1.0 - dtrans(j,k)) |
| 327 |
|
|
! end do |
| 328 |
|
|
down(:,:,k+1) = down(:,:,k)*dtrans(:,:,k) + b(:,:,k)*(1.0 - dtrans(:,:,k)) |
| 329 |
|
|
end do |
| 330 |
|
|
|
| 331 |
|
|
!do i = 1, size(t,1) |
| 332 |
|
|
!do j = 1, size(t,2) |
| 333 |
|
|
! do k = 1,n+1 |
| 334 |
|
|
! solar_down(i,j,k) = solar(j)*exp(-solar_tau(j,k)) |
| 335 |
|
|
! end do |
| 336 |
|
|
!end do |
| 337 |
|
|
!end do |
| 338 |
|
|
do k = 1,n+1 |
| 339 |
|
|
solar_down(:,:,k) = solar(:,:)*exp(-solar_tau(:,:,k)) |
| 340 |
|
|
end do |
| 341 |
|
|
|
| 342 |
|
|
surf_lw_down = down(:,:,n+1) |
| 343 |
|
|
net_surf_sw_down = solar_down(:,:,n+1)*(1. - albedo(:,:)) |
| 344 |
|
|
swin = solar_down(:,:,1) |
| 345 |
|
|
|
| 346 |
|
|
!------- downward sw flux surface ------- |
| 347 |
|
|
! if ( id_swdn_sfc > 0 ) then |
| 348 |
|
|
! used = send_data ( id_swdn_sfc, net_surf_sw_down, Time_diag) |
| 349 |
|
|
! endif |
| 350 |
|
|
!------- incoming sw flux toa ------- |
| 351 |
|
|
! if ( id_swdn_toa > 0 ) then |
| 352 |
|
|
! used = send_data ( id_swdn_toa, swin, Time_diag) |
| 353 |
|
|
! endif |
| 354 |
|
|
!------- downward lw flux surface ------- |
| 355 |
|
|
! if ( id_lwdn_sfc > 0 ) then |
| 356 |
|
|
! used = send_data ( id_lwdn_sfc, surf_lw_down, Time_diag) |
| 357 |
|
|
! endif |
| 358 |
|
|
|
| 359 |
|
|
! ------------------------------------------------------------------------- |
| 360 |
|
|
deallocate (swin) |
| 361 |
|
|
deallocate (tau, solar_tau) |
| 362 |
|
|
deallocate (ss, solar, tau_0, solar_tau_0, p2) |
| 363 |
|
|
! ------------------------------------------------------------------------- |
| 364 |
|
|
|
| 365 |
|
|
return |
| 366 |
|
|
end subroutine radiation_down |
| 367 |
|
|
|
| 368 |
|
|
! ================================================================================== |
| 369 |
|
|
|
| 370 |
|
|
!subroutine radiation_up (is, js, Time_diag, lat, p_half, t_surf, t, tdt) |
| 371 |
jmc |
1.2 |
subroutine radiation_up ( is, js, Time_diag, lat, p_half, t_surf, t, tdt, olr, & |
| 372 |
|
|
albedo, dtrans, b, down, solar_down, & |
| 373 |
jmc |
1.1 |
myThid ) |
| 374 |
|
|
|
| 375 |
|
|
! Now complete the radiation calculation by computing the upward and net fluxes. |
| 376 |
|
|
|
| 377 |
|
|
integer, intent(in) :: is, js |
| 378 |
|
|
!type(time_type), intent(in) :: Time_diag |
| 379 |
|
|
real, intent(in) :: Time_diag |
| 380 |
|
|
real, intent(in) , dimension(:,:) :: lat |
| 381 |
|
|
real, intent(in) , dimension(:,:) :: t_surf |
| 382 |
|
|
real, intent(in) , dimension(:,:,:) :: t, p_half |
| 383 |
|
|
real, intent(inout), dimension(:,:,:) :: tdt |
| 384 |
jmc |
1.2 |
real, intent(out), dimension(:,:) :: olr |
| 385 |
jmc |
1.1 |
real, intent(in), dimension(:,:) :: albedo |
| 386 |
|
|
!real, intent(in), dimension(:,:) :: dtrans |
| 387 |
|
|
real, intent(in), dimension(:,:,:) :: dtrans |
| 388 |
|
|
real, intent(in), dimension(:,:,:) :: b |
| 389 |
|
|
real, intent(in), dimension(:,:,:) :: down |
| 390 |
|
|
real, intent(in), dimension(:,:,:) :: solar_down |
| 391 |
|
|
integer, intent(in) :: myThid |
| 392 |
|
|
|
| 393 |
|
|
integer :: i, j, k, n |
| 394 |
|
|
integer :: im, jm |
| 395 |
|
|
|
| 396 |
|
|
logical :: used |
| 397 |
|
|
|
| 398 |
|
|
! ------------------------------------------------------------------------- |
| 399 |
|
|
real, allocatable, dimension(:,:) :: b_surf |
| 400 |
|
|
real, allocatable, dimension(:,:,:) :: tdt_rad, entrop_rad, tdt_sw |
| 401 |
|
|
real, allocatable, dimension(:,:,:) :: up, net, flux_rad, flux_sw |
| 402 |
|
|
! ------------------------------------------------------------------------- |
| 403 |
|
|
|
| 404 |
|
|
n = size(t,3) |
| 405 |
|
|
im = size(t,1) |
| 406 |
|
|
jm = size(t,2) |
| 407 |
|
|
|
| 408 |
|
|
! ------------------------------------------------------------------------- |
| 409 |
|
|
allocate (tdt_rad (im, jm, n)) |
| 410 |
|
|
allocate (tdt_sw (im, jm, n)) |
| 411 |
|
|
allocate (entrop_rad (im, jm, n)) |
| 412 |
|
|
allocate (up (im, jm, n+1)) |
| 413 |
|
|
allocate (net (im, jm, n+1)) |
| 414 |
|
|
allocate (flux_rad (im, jm, n+1)) |
| 415 |
|
|
allocate (flux_sw (im, jm, n+1)) |
| 416 |
|
|
allocate (b_surf (im, jm)) |
| 417 |
|
|
! ------------------------------------------------------------------------- |
| 418 |
|
|
|
| 419 |
|
|
! total flux from surface |
| 420 |
|
|
b_surf = stefan*t_surf*t_surf*t_surf*t_surf |
| 421 |
|
|
|
| 422 |
|
|
! first deal with non-window upward flux |
| 423 |
|
|
up(:,:,n+1) = b_surf*(1.0-window) |
| 424 |
|
|
!do k = n,1,-1 |
| 425 |
|
|
! do j = 1, size(t,2) |
| 426 |
|
|
! up(:,j,k) = up(:,j,k+1)*dtrans(j,k) + b(:,j,k)*(1.0 - dtrans(j,k)) |
| 427 |
|
|
! end do |
| 428 |
|
|
!end do |
| 429 |
|
|
do k = n,1,-1 |
| 430 |
|
|
up(:,:,k) = up(:,:,k+1)*dtrans(:,:,k) + b(:,:,k)*(1.0 - dtrans(:,:,k)) |
| 431 |
|
|
end do |
| 432 |
|
|
|
| 433 |
|
|
! add upward flux in spectral window |
| 434 |
|
|
do k = 1,n+1 |
| 435 |
|
|
up(:,:,k) = up(:,:,k) + b_surf(:,:)*window |
| 436 |
|
|
end do |
| 437 |
|
|
|
| 438 |
|
|
do k = 1,n+1 |
| 439 |
|
|
net(:,:,k) = up(:,:,k)-down(:,:,k) |
| 440 |
|
|
flux_sw(:,:,k) = albedo(:,:)*solar_down(:,:,n+1) - solar_down(:,:,k) |
| 441 |
|
|
flux_rad(:,:,k) = net(:,:,k) + flux_sw(:,:,k) |
| 442 |
|
|
end do |
| 443 |
|
|
|
| 444 |
|
|
do k = 1,n |
| 445 |
|
|
tdt_rad(:,:,k) = (net(:,:,k+1) - net(:,:,k) - solar_down(:,:,k+1) + solar_down(:,:,k)) & |
| 446 |
|
|
*grav/(cp_air*(p_half(:,:,k+1)-p_half(:,:,k))) |
| 447 |
|
|
tdt_sw(:,:,k) = (- solar_down(:,:,k+1) + solar_down(:,:,k)) & |
| 448 |
|
|
*grav/(cp_air*(p_half(:,:,k+1)-p_half(:,:,k))) |
| 449 |
|
|
tdt(:,:,k) = tdt(:,:,k) + tdt_rad(:,:,k) |
| 450 |
|
|
end do |
| 451 |
|
|
|
| 452 |
|
|
olr = up(:,:,1) |
| 453 |
|
|
|
| 454 |
|
|
!------- outgoing lw flux toa (olr) ------- |
| 455 |
|
|
! if ( id_olr > 0 ) then |
| 456 |
|
|
! used = send_data ( id_olr, olr, Time_diag) |
| 457 |
|
|
! endif |
| 458 |
|
|
!------- upward lw flux surface ------- |
| 459 |
|
|
! if ( id_lwup_sfc > 0 ) then |
| 460 |
|
|
! used = send_data ( id_lwup_sfc, b_surf, Time_diag) |
| 461 |
|
|
! endif |
| 462 |
|
|
!------- temperature tendency due to radiation ------------ |
| 463 |
|
|
! if ( id_tdt_rad > 0 ) then |
| 464 |
|
|
! used = send_data ( id_tdt_rad, tdt_rad, Time_diag) |
| 465 |
|
|
! endif |
| 466 |
|
|
! if ( id_tdt_sw > 0 ) then |
| 467 |
|
|
! used = send_data ( id_tdt_sw, tdt_sw, Time_diag) |
| 468 |
|
|
! endif |
| 469 |
|
|
!------- total radiative flux (at half levels) ----------- |
| 470 |
|
|
! if ( id_flux_rad > 0 ) then |
| 471 |
|
|
! used = send_data ( id_flux_rad, flux_rad, Time_diag) |
| 472 |
|
|
! endif |
| 473 |
|
|
!------- longwave radiative flux (at half levels) -------- |
| 474 |
|
|
! if ( id_flux_lw > 0 ) then |
| 475 |
|
|
! used = send_data ( id_flux_lw, net, Time_diag) |
| 476 |
|
|
! endif |
| 477 |
|
|
! if ( id_flux_sw > 0 ) then |
| 478 |
|
|
! used = send_data ( id_flux_sw, flux_sw, Time_diag) |
| 479 |
|
|
! endif |
| 480 |
|
|
! if ( id_entrop_rad > 0 ) then |
| 481 |
|
|
! do k=1,n |
| 482 |
|
|
! entrop_rad(:,:,k) =tdt_rad(:,:,k)/t(:,:,k)*p_half(:,:,n+1)/1.e5 |
| 483 |
|
|
! end do |
| 484 |
|
|
! used = send_data ( id_entrop_rad, entrop_rad, Time_diag) |
| 485 |
|
|
! endif |
| 486 |
|
|
|
| 487 |
|
|
! ------------------------------------------------------------------------- |
| 488 |
|
|
deallocate (tdt_rad, tdt_sw, entrop_rad) |
| 489 |
jmc |
1.2 |
deallocate (up, net, flux_rad, flux_sw) |
| 490 |
jmc |
1.1 |
deallocate (b_surf) |
| 491 |
|
|
! ------------------------------------------------------------------------- |
| 492 |
|
|
|
| 493 |
|
|
return |
| 494 |
|
|
end subroutine radiation_up |
| 495 |
|
|
|
| 496 |
|
|
! ================================================================================== |
| 497 |
|
|
|
| 498 |
|
|
subroutine radiation_end |
| 499 |
|
|
|
| 500 |
|
|
!deallocate (b, tdt_rad, tdt_sw, entrop_rad) |
| 501 |
|
|
!!deallocate (up, down, net, solar_down, flux_rad, flux_sw) |
| 502 |
|
|
!deallocate (up, net, flux_rad, flux_sw) |
| 503 |
|
|
!deallocate (b_surf, olr, swin, albedo) |
| 504 |
|
|
!!deallocate (dtrans) |
| 505 |
|
|
!deallocate (tau, solar_tau) |
| 506 |
|
|
!deallocate (ss, solar, tau_0, solar_tau_0, p2) |
| 507 |
|
|
|
| 508 |
|
|
end subroutine radiation_end |
| 509 |
|
|
|
| 510 |
|
|
! ================================================================================== |
| 511 |
|
|
|
| 512 |
|
|
end module radiation_mod |
Parent Directory
| 
Revision Log
| 
Revision Graph