source: dynamico_lmdz/aquaplanet/LMDZ5/libf/phylmd/aeropt.F90 @ 3814

Last change on this file since 3814 was 3814, checked in by ymipsl, 10 years ago

remove all dynamic dependency in LMDZ physics except for the include "dimensions.h"

YM

File size: 4.5 KB
Line 
1
2! $Id: aeropt.F90 1992 2014-03-05 13:19:12Z lguez $
3
4SUBROUTINE aeropt(pplay, paprs, t_seri, msulfate, rhcl, tau_ae, piz_ae, &
5    cg_ae, ai)
6
7  USE dimphy
8  IMPLICIT NONE
9
10
11
12  ! ym#include "dimensions.h"
13  ! ym#include "dimphy.h"
14  include "YOMCST.h"
15
16  ! Arguments:
17
18  REAL, INTENT (IN) :: paprs(klon, klev+1)
19  REAL, INTENT (IN) :: pplay(klon, klev), t_seri(klon, klev)
20  REAL, INTENT (IN) :: msulfate(klon, klev) ! masse sulfate ug SO4/m3  [ug/m^3]
21  REAL, INTENT (IN) :: rhcl(klon, klev) ! humidite relative ciel clair
22  REAL, INTENT (OUT) :: tau_ae(klon, klev, 2) ! epaisseur optique aerosol
23  REAL, INTENT (OUT) :: piz_ae(klon, klev, 2) ! single scattering albedo aerosol
24  REAL, INTENT (OUT) :: cg_ae(klon, klev, 2) ! asymmetry parameter aerosol
25  REAL, INTENT (OUT) :: ai(klon) ! POLDER aerosol index
26
27  ! Local
28
29  INTEGER i, k, inu
30  INTEGER rh_num, nbre_rh
31  PARAMETER (nbre_rh=12)
32  REAL rh_tab(nbre_rh)
33  REAL rh_max, delta, rh
34  PARAMETER (rh_max=95.)
35  DATA rh_tab/0., 10., 20., 30., 40., 50., 60., 70., 80., 85., 90., 95./
36  REAL zrho, zdz
37  REAL taue670(klon) ! epaisseur optique aerosol absorption 550 nm
38  REAL taue865(klon) ! epaisseur optique aerosol extinction 865 nm
39  REAL alpha_aer_sulfate(nbre_rh, 5) !--unit m2/g SO4
40  REAL alphasulfate
41
42  CHARACTER (LEN=20) :: modname = 'aeropt'
43  CHARACTER (LEN=80) :: abort_message
44
45
46  ! Proprietes optiques
47
48  REAL alpha_aer(nbre_rh, 2) !--unit m2/g SO4
49  REAL cg_aer(nbre_rh, 2)
50  DATA alpha_aer/.500130E+01, .500130E+01, .500130E+01, .500130E+01, &
51    .500130E+01, .616710E+01, .826850E+01, .107687E+02, .136976E+02, &
52    .162972E+02, .211690E+02, .354833E+02, .139460E+01, .139460E+01, &
53    .139460E+01, .139460E+01, .139460E+01, .173910E+01, .244380E+01, &
54    .332320E+01, .440120E+01, .539570E+01, .734580E+01, .136038E+02/
55  DATA cg_aer/.619800E+00, .619800E+00, .619800E+00, .619800E+00, &
56    .619800E+00, .662700E+00, .682100E+00, .698500E+00, .712500E+00, &
57    .721800E+00, .734600E+00, .755800E+00, .545600E+00, .545600E+00, &
58    .545600E+00, .545600E+00, .545600E+00, .583700E+00, .607100E+00, &
59    .627700E+00, .645800E+00, .658400E+00, .676500E+00, .708500E+00/
60  DATA alpha_aer_sulfate/4.910, 4.910, 4.910, 4.910, 6.547, 7.373, 8.373, &
61    9.788, 12.167, 14.256, 17.924, 28.433, 1.453, 1.453, 1.453, 1.453, 2.003, &
62    2.321, 2.711, 3.282, 4.287, 5.210, 6.914, 12.305, 4.308, 4.308, 4.308, &
63    4.308, 5.753, 6.521, 7.449, 8.772, 11.014, 12.999, 16.518, 26.772, 3.265, &
64    3.265, 3.265, 3.265, 4.388, 5.016, 5.775, 6.868, 8.745, 10.429, 13.457, &
65    22.538, 2.116, 2.116, 2.116, 2.116, 2.882, 3.330, 3.876, 4.670, 6.059, &
66    7.327, 9.650, 16.883/
67
68  DO i = 1, klon
69    taue670(i) = 0.0
70    taue865(i) = 0.0
71  END DO
72
73  DO k = 1, klev
74    DO i = 1, klon
75      IF (t_seri(i,k)==0) WRITE (*, *) 'aeropt T ', i, k, t_seri(i, k)
76      IF (pplay(i,k)==0) WRITE (*, *) 'aeropt p ', i, k, pplay(i, k)
77      zrho = pplay(i, k)/t_seri(i, k)/rd ! kg/m3
78      zdz = (paprs(i,k)-paprs(i,k+1))/zrho/rg ! m
79      rh = min(rhcl(i,k)*100., rh_max)
80      rh_num = int(rh/10.+1.)
81      IF (rh<0.) THEN
82        abort_message = 'aeropt: RH < 0 not possible'
83        CALL abort_physic(modname, abort_message, 1)
84      END IF
85      IF (rh>85.) rh_num = 10
86      IF (rh>90.) rh_num = 11
87      delta = (rh-rh_tab(rh_num))/(rh_tab(rh_num+1)-rh_tab(rh_num))
88
89      inu = 1
90      tau_ae(i, k, inu) = alpha_aer(rh_num, inu) + delta*(alpha_aer(rh_num+1, &
91        inu)-alpha_aer(rh_num,inu))
92      tau_ae(i, k, inu) = tau_ae(i, k, inu)*msulfate(i, k)*zdz*1.E-6
93      piz_ae(i, k, inu) = 1.0
94      cg_ae(i, k, inu) = cg_aer(rh_num, inu) + delta*(cg_aer(rh_num+1,inu)- &
95        cg_aer(rh_num,inu))
96
97      inu = 2
98      tau_ae(i, k, inu) = alpha_aer(rh_num, inu) + delta*(alpha_aer(rh_num+1, &
99        inu)-alpha_aer(rh_num,inu))
100      tau_ae(i, k, inu) = tau_ae(i, k, inu)*msulfate(i, k)*zdz*1.E-6
101      piz_ae(i, k, inu) = 1.0
102      cg_ae(i, k, inu) = cg_aer(rh_num, inu) + delta*(cg_aer(rh_num+1,inu)- &
103        cg_aer(rh_num,inu))
104      ! jq
105      ! jq for aerosol index
106
107      alphasulfate = alpha_aer_sulfate(rh_num, 4) + &
108        delta*(alpha_aer_sulfate(rh_num+1,4)-alpha_aer_sulfate(rh_num,4)) !--m2/g
109
110      taue670(i) = taue670(i) + alphasulfate*msulfate(i, k)*zdz*1.E-6
111
112      alphasulfate = alpha_aer_sulfate(rh_num, 5) + &
113        delta*(alpha_aer_sulfate(rh_num+1,5)-alpha_aer_sulfate(rh_num,5)) !--m2/g
114
115      taue865(i) = taue865(i) + alphasulfate*msulfate(i, k)*zdz*1.E-6
116
117    END DO
118  END DO
119
120  DO i = 1, klon
121    ai(i) = (-log(max(taue670(i),0.0001)/max(taue865(i), &
122      0.0001))/log(670./865.))*taue865(i)
123  END DO
124
125  RETURN
126END SUBROUTINE aeropt
Note: See TracBrowser for help on using the repository browser.