ROMS
Loading...
Searching...
No Matches
lmd_swfrac.F
Go to the documentation of this file.
1#include "cppdefs.h"
2
3#if defined NONLINEAR && (defined LMD_SKPP || defined SOLAR_SOURCE) && \
4 defined solve3d
5
6 SUBROUTINE lmd_swfrac_tile (ng, tile, &
7 & LBi, UBi, LBj, UBj, &
8 & IminS, ImaxS, JminS, JmaxS, &
9 & Zscale, Z, swdk)
10!
11!git $Id$
12!================================================== Hernan G. Arango ===
13! Copyright (c) 2002-2025 The ROMS Group !
14! Licensed under a MIT/X style license !
15! See License_ROMS.md !
16!=======================================================================
17! !
18! This subroutine computes the fraction of solar shortwave flux !
19! penetrating to specified depth (times Zscale) due to exponential !
20! decay in Jerlov water type. !
21! !
22! On Input: !
23! !
24! Zscale scale factor to apply to depth array. !
25! Z vertical height (meters, negative) for !
26! desired solar short-wave fraction. !
27! !
28! On Output: !
29! !
30! swdk shortwave (radiation) fractional decay. !
31! !
32! Reference: !
33! !
34! Paulson, C.A., and J.J. Simpson, 1977: Irradiance meassurements !
35! in the upper ocean, J. Phys. Oceanogr., 7, 952-956. !
36! !
37! This routine was adapted from Bill Large 1995 code. !
38! !
39!=======================================================================
40!
41 USE mod_param
42 USE mod_mixing
43 USE mod_scalars
44!
45 implicit none
46!
47! Imported variable declarations.
48!
49 integer, intent(in) :: ng, tile
50 integer, intent(in) :: LBi, UBi, LBj, UBj
51 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
52
53 real(r8), intent(in) :: Zscale
54
55 real(r8), intent(in) :: Z(IminS:ImaxS,JminS:JmaxS)
56 real(r8), intent(out) :: swdk(IminS:ImaxS,JminS:JmaxS)
57!
58! Local variable declarations.
59!
60 integer :: Jindex, i, j
61
62 real(r8), dimension(IminS:ImaxS) :: fac1, fac2, fac3
63
64# include "set_bounds.h"
65!
66!-----------------------------------------------------------------------
67! Use Paulson and Simpson (1977) two wavelength bands solar
68! absorption model.
69!-----------------------------------------------------------------------
70!
71 DO j=jstr,jend
72 DO i=istr,iend
73 jindex=int(mixing(ng)%Jwtype(i,j))
74 fac1(i)=zscale/lmd_mu1(jindex)
75 fac2(i)=zscale/lmd_mu2(jindex)
76 fac3(i)=lmd_r1(jindex)
77 END DO
78!!DIR$ VECTOR ALWAYS
79 DO i=istr,iend
80 swdk(i,j)=exp(z(i,j)*fac1(i))*fac3(i)+ &
81 & exp(z(i,j)*fac2(i))*(1.0_r8-fac3(i))
82 END DO
83 END DO
84 RETURN
85 END SUBROUTINE lmd_swfrac_tile
86#else
87 SUBROUTINE lmd_swfrac
88 RETURN
89 END SUBROUTINE lmd_swfrac
90#endif
91
subroutine lmd_swfrac_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, zscale, z, swdk)
Definition lmd_swfrac.F:10
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
real(r8), dimension(9) lmd_mu1
real(r8), dimension(9) lmd_mu2
real(r8), dimension(9) lmd_r1