Changeset 3757 for trunk/LMDZ.MARS/libf
- Timestamp:
- May 9, 2025, 9:22:41 AM (2 months ago)
- Location:
- trunk/LMDZ.MARS/libf/phymars
- Files:
-
- 12 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.MARS/libf/phymars/calldrag_noro_mod.F90
r2642 r3757 31 31 use dimradmars_mod, only: ndomainsz 32 32 use drag_noro_mod, only: drag_noro 33 use sugwd_mod, only: sugwd 33 34 USE ioipsl_getin_p_mod, ONLY : getin_p 34 35 -
trunk/LMDZ.MARS/libf/phymars/flusv.F
r1268 r3757 1 MODULE flusv_mod 2 3 IMPLICIT NONE 4 5 CONTAINS 6 1 7 SUBROUTINE flusv(KDLON,nsf,n,omega,g,tau,emis,bh,bsol,fah,fdh) 2 use dimradmars_mod, only: ndlo2, ndlon, nflev3 IMPLICIT NONE4 c....................................................................... 5 c 6 c c alcul des flux ascendant et descendant aux interfaces entre n couches7 c * dans l'infrarouge8 c * B est une fonction lineaire de $\tau$ a l'interieur de chaque couche9 c * le B du sol peut etre different de celui qui correspond au profil10 c de la n-ieme couche11 c * l'hypothese est une hypothese a deux flux isotropes sur chaque12 c hemisphere ("hemispheric constant") + " source function technique"13 c ( voirToon et al. 1988)14 c * le flux descendant en haut de l'atmosphere est nul15 c * l es couches sont numerotees du haut de l'atmosphere vers le sol16 c 17 c in : * KDLON ---> dimension de vectorisation8 use dimradmars_mod, only: ndlo2, ndlon, nflev 9 IMPLICIT NONE 10 c....................................................................... 11 c 12 c compute the upward and downward fluxes at the interface between n layers 13 c * in the infrared 14 c * B is a linear function of $\tau$ in each layer 15 c * B at the surface can be different than what corresponds to the profile 16 c in the n-th layer 17 c * the work hypothes isthat we have two isotropic fluxes for each 18 c hemisphere ("hemispheric constant") + "technical source function" 19 c (see Toon et al. 1988) 20 c * the downward flux at the top of the atmosphere is zero 21 c * layers are numbered from the top of the atmosphere to the ground 22 c 23 c in : * KDLON ---> vectorisation dimension 18 24 c * nsf ---> nsf=0 ==> "hemispheric constant" 19 25 c nsf>0 ==> "hemispheric constant" + "source function" 20 c * n ---> nombre de couches 21 c * omega(i) ---> single scattering albedo pour la i-eme couche 22 c * g(i) ---> asymmetry parameter pour la i-eme couche 23 c * tau(i) ---> epaisseur optique de la i-eme couche 24 c * emis ---> emissivite du sol 25 c * bh(i) ---> luminance du corps noir en haut de la i-eme 26 c couche, bh(n+1) pour la valeur au sol qui 27 c correspond au profil de la n-ieme couche 28 c * bsol ---> luminance du corps noir au sol 29 c 30 c out : * fah(i) ---> flux ascendant en haut de la i-eme couche, 31 c fah(n+1) pour le sol 32 c * fdh(i) ---> flux descendant en haut de la i-eme couche, 33 c fdh(n+1) pour le sol 34 c 35 c....................................................................... 36 c declaration des arguments 37 c 38 INTEGER KDLON,nsf,n 39 REAL omega(NDLO2,n),g(NDLO2,n),tau(NDLO2,n),emis(NDLO2) 40 &,bh(NDLO2,n+1),bsol(NDLO2),fah(NDLO2,n+1),fdh(NDLO2,n+1) 41 c....................................................................... 42 c declaration des variables locales 43 c 44 REAL pi 45 PARAMETER (pi=3.141592653589793E+0) 26 c * n ---> number of layers 27 c * omega(i) ---> single scattering albedo for the i-th layer 28 c * g(i) ---> asymmetry parameter for the i-th layer 29 c * tau(i) ---> optical thickness of the i-th layer 30 c * emis ---> ground emissivity 31 c * bh(i) ---> black body luminance at the top of the i-th layer, 32 c bh(n+1) for the ground value which 33 c corresponds to the profile for the n-th layer 34 c * bsol ---> black body luminance of the ground 35 c 36 c out : * fah(i) ---> upward flux at the top of the i-th layer, 37 c fah(n+1) for the ground 38 c * fdh(i) ---> downward flux at the top of the i-th layer, 39 c fdh(n+1) for the ground 40 c 41 c....................................................................... 42 c arguments 43 c 44 INTEGER,INTENT(IN) :: KDLON,nsf,n 45 REAL,INTENT(IN) :: omega(NDLO2,n),g(NDLO2,n) 46 REAL,INTENT(IN) :: tau(NDLO2,n),emis(NDLO2) 47 REAL,INTENT(IN) :: bh(NDLO2,n+1),bsol(NDLO2) 48 REAL,INTENT(OUT) :: fah(NDLO2,n+1),fdh(NDLO2,n+1) 49 c....................................................................... 50 c local variables 51 c 52 REAL,PARAMETER :: pi=3.141592653589793E+0 46 53 INTEGER iv,i,j 47 54 REAL beta,gama1,gama2,amu1,grgama,b0,b1 … … 58 65 & ,alpha1(NDLON,2*nflev),alpha2(NDLON,2*nflev) 59 66 & ,sigma1(NDLON,2*nflev),sigma2(NDLON,2*nflev) 60 INTEGER nq 61 PARAMETER (nq=8) 62 REAL x(nq),w(nq),gri(NDLON,nq) 63 DATA x/1.9855071751231860E-2 , 0.1016667612931866E+0 67 INTEGER,PARAMETER :: nq=8 68 REAL,PARAMETER :: x(nq) = 69 & (/1.9855071751231860E-2 , 0.1016667612931866E+0 64 70 & , 0.2372337950418355E+0 , 0.4082826787521751E+0 65 71 & , 0.5917173212478250E+0 , 0.7627662049581645E+0 66 & , 0.8983332387068134E+0 , 0.9801449282487682E+0/ 67 DATA w/5.0614268145185310E-2 , 0.1111905172266872E+0 72 & , 0.8983332387068134E+0 , 0.9801449282487682E+0/) 73 REAL,PARAMETER :: w(nq) = 74 & (/5.0614268145185310E-2 , 0.1111905172266872E+0 68 75 & , 0.1568533229389437E+0 , 0.1813418916891810E+0 69 76 & , 0.1813418916891810E+0 , 0.1568533229389437E+0 70 & , 0.1111905172266872E+0 , 5.0614268145185310E-2/ 71 c....................................................................... 72 c 73 c....................................................................... 74 do 10001 i=1,n 75 do 99999 iv=1,KDLON 77 & , 0.1111905172266872E+0 , 5.0614268145185310E-2/) 78 REAL :: gri(NDLON,nq) 79 c....................................................................... 80 c 81 c....................................................................... 82 do i=1,n 83 do iv=1,KDLON 76 84 beta=(1.E+0-g(iv,i))/2.E+0 77 85 gama1=2.E+0*(1.E+0-omega(iv,i)*(1.E+0-beta)) … … 81 89 grgama=(gama1-alambda(iv,i))/gama2 82 90 c 83 c ici une petite bidouille : si l'epaisseur optique d'une couche84 c est trop faible, $dB \over d\tau$ devient tres grand et le schema85 c s'ecroule. dans ces cas, on fait l'hypothese de couche isotherme.91 c small hack here : if the optical depth of a layer is too small, 92 c $dB \over d\tau$ becomes very large and the scheme fails. 93 c In those cases we assume an isothermal layer. 86 94 c 87 95 if (tau(iv,i).gt.1.E-3) then 88 b0=bh(iv,i)89 b1=(bh(iv,i+1)-b0)/tau(iv,i)96 b0=bh(iv,i) 97 b1=(bh(iv,i+1)-b0)/tau(iv,i) 90 98 else 91 b0=(bh(iv,i)+bh(iv,i+1))/2.E+092 b1=0.E+099 b0=(bh(iv,i)+bh(iv,i+1))/2.E+0 100 b1=0.E+0 93 101 endif 94 102 c … … 111 119 sigma2(iv,i)=alpha2(iv,i) 112 120 c 113 99999 continue 114 10001 continue 115 c....................................................................... 116 do 99998iv=1,KDLON117 a(iv,1)=0.E+0118 b(iv,1)=e1(iv,1)119 d(iv,1)=-e2(iv,1)120 e(iv,1)=-cdh(iv,1)121 99998 continue 122 c 123 do 10002i=1,n-1124 j=2*i+1125 do 99997iv=1,KDLON126 a(iv,j)=e2(iv,i)*e3(iv,i)-e4(iv,i)*e1(iv,i)127 b(iv,j)=e1(iv,i)*e1(iv,i+1)-e3(iv,i)*e3(iv,i+1)128 d(iv,j)=e3(iv,i)*e4(iv,i+1)-e1(iv,i)*e2(iv,i+1)129 e(iv,j)=e3(iv,i)*(cah(iv,i+1)-cab(iv,i))130 & +e1(iv,i)*(cdb(iv,i)-cdh(iv,i+1))131 99997 continue 132 10002 continue 133 c 134 do 10003i=1,n-1135 j=2*i136 do 99996iv=1,KDLON137 a(iv,j)=e2(iv,i+1)*e1(iv,i)-e3(iv,i)*e4(iv,i+1)138 b(iv,j)=e2(iv,i)*e2(iv,i+1)-e4(iv,i)*e4(iv,i+1)139 d(iv,j)=e1(iv,i+1)*e4(iv,i+1)-e2(iv,i+1)*e3(iv,i+1)140 e(iv,j)=e2(iv,i+1)*(cah(iv,i+1)-cab(iv,i))141 & +e4(iv,i+1)*(cdb(iv,i)-cdh(iv,i+1))142 99996 continue 143 10003 continue 121 enddo ! of do iv=1,KDLON 122 enddo ! of do i=1,n 123 c....................................................................... 124 do iv=1,KDLON 125 a(iv,1)=0.E+0 126 b(iv,1)=e1(iv,1) 127 d(iv,1)=-e2(iv,1) 128 e(iv,1)=-cdh(iv,1) 129 enddo 130 c 131 do i=1,n-1 132 j=2*i+1 133 do iv=1,KDLON 134 a(iv,j)=e2(iv,i)*e3(iv,i)-e4(iv,i)*e1(iv,i) 135 b(iv,j)=e1(iv,i)*e1(iv,i+1)-e3(iv,i)*e3(iv,i+1) 136 d(iv,j)=e3(iv,i)*e4(iv,i+1)-e1(iv,i)*e2(iv,i+1) 137 e(iv,j)=e3(iv,i)*(cah(iv,i+1)-cab(iv,i)) 138 & +e1(iv,i)*(cdb(iv,i)-cdh(iv,i+1)) 139 enddo 140 enddo ! of do i=1,n-1 141 c 142 do i=1,n-1 143 j=2*i 144 do iv=1,KDLON 145 a(iv,j)=e2(iv,i+1)*e1(iv,i)-e3(iv,i)*e4(iv,i+1) 146 b(iv,j)=e2(iv,i)*e2(iv,i+1)-e4(iv,i)*e4(iv,i+1) 147 d(iv,j)=e1(iv,i+1)*e4(iv,i+1)-e2(iv,i+1)*e3(iv,i+1) 148 e(iv,j)=e2(iv,i+1)*(cah(iv,i+1)-cab(iv,i)) 149 & +e4(iv,i+1)*(cdb(iv,i)-cdh(iv,i+1)) 150 enddo 151 enddo ! of do i=1,n-1 144 152 c 145 153 j=2*n 146 do 99995 iv=1,KDLON 147 a(iv,j)=e1(iv,n)-(1.E+0-emis(iv))*e3(iv,n) 148 b(iv,j)=e2(iv,n)-(1.E+0-emis(iv))*e4(iv,n) 149 d(iv,j)=0.E+0 150 e(iv,j)=emis(iv)*pi*bsol(iv)-cab(iv,n)+(1.E+0-emis(iv))*cdb(iv,n) 151 99995 continue 154 do iv=1,KDLON 155 a(iv,j)=e1(iv,n)-(1.E+0-emis(iv))*e3(iv,n) 156 b(iv,j)=e2(iv,n)-(1.E+0-emis(iv))*e4(iv,n) 157 d(iv,j)=0.E+0 158 e(iv,j)=emis(iv)*pi*bsol(iv)-cab(iv,n) 159 & +(1.E+0-emis(iv))*cdb(iv,n) 160 enddo 152 161 c....................................................................... 153 162 call sys3v(KDLON,2*n,a,b,d,e,y) 154 163 c....................................................................... 155 do 10004i=1,n156 do 99994iv=1,KDLON157 grg(iv,i)=grg(iv,i)*(y(iv,2*i-1)+y(iv,2*i))158 grh(iv,i)=grh(iv,i)*(y(iv,2*i-1)-y(iv,2*i))159 grj(iv,i)=grj(iv,i)*(y(iv,2*i-1)+y(iv,2*i))160 grk(iv,i)=grk(iv,i)*(y(iv,2*i-1)-y(iv,2*i))161 99994 continue 162 10004 continue 163 c....................................................................... 164 c les valeurs de flux "hemispheric constant"164 do i=1,n 165 do iv=1,KDLON 166 grg(iv,i)=grg(iv,i)*(y(iv,2*i-1)+y(iv,2*i)) 167 grh(iv,i)=grh(iv,i)*(y(iv,2*i-1)-y(iv,2*i)) 168 grj(iv,i)=grj(iv,i)*(y(iv,2*i-1)+y(iv,2*i)) 169 grk(iv,i)=grk(iv,i)*(y(iv,2*i-1)-y(iv,2*i)) 170 enddo 171 enddo 172 c....................................................................... 173 c values of "hemispheric constant" fluxes 165 174 c 166 175 IF (nsf.eq.0) THEN 167 do 10005 i=1,n 168 do 99993 iv=1,KDLON 169 fah(iv,i)=e3(iv,i)*y(iv,2*i-1)-e4(iv,i)*y(iv,2*i)+cah(iv,i) 170 fdh(iv,i)=e1(iv,i)*y(iv,2*i-1)-e2(iv,i)*y(iv,2*i)+cdh(iv,i) 171 99993 continue 172 10005 continue 173 do 99992 iv=1,KDLON 174 fah(iv,n+1)=e1(iv,n)*y(iv,2*n-1)+e2(iv,n)*y(iv,2*n)+cab(iv,n) 175 fdh(iv,n+1)=e3(iv,n)*y(iv,2*n-1)+e4(iv,n)*y(iv,2*n)+cdb(iv,n) 176 99992 continue 177 GOTO 10100 178 ENDIF 179 c....................................................................... 180 c passage a "source function" 181 c 182 c on applique une quadrature de dimension nq 183 c x est le vecteur des \mu de la quadrature 184 c w est le vecteur des poids correspondants 185 c nq est fixe en parameter 186 c x et w sont fixes en data 187 c 188 c....................................................................... 189 c on part d'en haut et on descent selon les nq angles pour calculer 190 c tous les flux descendants 191 c 192 do 10006 j=1,nq 193 do 99991 iv=1,KDLON 194 gri(iv,j)=0.E+0 195 99991 continue 196 10006 continue 197 do 99990 iv=1,KDLON 198 fdh(iv,1)=0.E+0 199 99990 continue 200 do 10007 i=1,n 201 do 10008 j=1,nq 202 do 99989 iv=1,KDLON 203 gri(iv,j)=gri(iv,j)*exp(-tau(iv,i)/x(j)) 204 &+grj(iv,i)/(alambda(iv,i)*x(j)+1.E+0) 205 &*(1.E+0-exp(-tau(iv,i)*(alambda(iv,i)+1.E+0/x(j)))) 206 &+grk(iv,i)/(alambda(iv,i)*x(j)-1.E+0) 207 &*(exp(-tau(iv,i)/x(j))-exp(-tau(iv,i)*alambda(iv,i))) 208 &+sigma1(iv,i)*(1.E+0-exp(-tau(iv,i)/x(j))) 209 &+sigma2(iv,i)*(x(j)*exp(-tau(iv,i)/x(j))+tau(iv,i)-x(j)) 210 99989 continue 211 10008 continue 212 do 99988 iv=1,KDLON 213 fdh(iv,i+1)=0.E+0 214 99988 continue 215 do 10009 j=1,nq 216 do 99987 iv=1,KDLON 217 fdh(iv,i+1)=fdh(iv,i+1)+w(j)*x(j)*gri(iv,j) 218 99987 continue 219 10009 continue 220 10007 continue 221 c....................................................................... 222 c on applique la condition de reflexion a sol 223 c 224 do 99986 iv=1,KDLON 225 fah(iv,n+1)=(1.E+0-emis(iv))*fdh(iv,n+1)+pi*emis(iv)*bsol(iv) 226 99986 continue 227 do 10010 j=1,nq 228 do 99985 iv=1,KDLON 229 gri(iv,j)=2.E+0*fah(iv,n+1) 230 99985 continue 231 10010 continue 232 c....................................................................... 233 c on remonte pour calculer tous les flux ascendants 176 do i=1,n 177 do iv=1,KDLON 178 fah(iv,i)=e3(iv,i)*y(iv,2*i-1)-e4(iv,i)*y(iv,2*i)+cah(iv,i) 179 fdh(iv,i)=e1(iv,i)*y(iv,2*i-1)-e2(iv,i)*y(iv,2*i)+cdh(iv,i) 180 enddo 181 enddo 182 do iv=1,KDLON 183 fah(iv,n+1)=e1(iv,n)*y(iv,2*n-1)+e2(iv,n)*y(iv,2*n)+cab(iv,n) 184 fdh(iv,n+1)=e3(iv,n)*y(iv,2*n-1)+e4(iv,n)*y(iv,2*n)+cdb(iv,n) 185 enddo 186 ELSE 187 c....................................................................... 188 c going to the "source function" 189 c 190 c apply a quadrature over nq (fixed parameter) points 191 c x is the vector of the \mu of the quadrature 192 c w is the vector of corresponding weights 193 c x() et w() are fixed parameters 194 c 195 c....................................................................... 196 c start from the top and go down along the nq angles to compute all 197 c downward fluxes 198 c 199 do j=1,nq 200 do iv=1,KDLON 201 gri(iv,j)=0.E+0 202 enddo 203 enddo 204 do iv=1,KDLON 205 fdh(iv,1)=0.E+0 206 enddo 207 do i=1,n 208 do j=1,nq 209 do iv=1,KDLON 210 gri(iv,j)=gri(iv,j)*exp(-tau(iv,i)/x(j)) 211 & +grj(iv,i)/(alambda(iv,i)*x(j)+1.E+0) 212 & *(1.E+0-exp(-tau(iv,i)*(alambda(iv,i)+1.E+0/x(j)))) 213 & +grk(iv,i)/(alambda(iv,i)*x(j)-1.E+0) 214 & *(exp(-tau(iv,i)/x(j))-exp(-tau(iv,i)*alambda(iv,i))) 215 & +sigma1(iv,i)*(1.E+0-exp(-tau(iv,i)/x(j))) 216 & +sigma2(iv,i)*(x(j)*exp(-tau(iv,i)/x(j))+tau(iv,i)-x(j)) 217 enddo ! of do iv=1,KDLON 218 enddo ! of do j=1,nq 219 do iv=1,KDLON 220 fdh(iv,i+1)=0.E+0 221 enddo 222 do j=1,nq 223 do iv=1,KDLON 224 fdh(iv,i+1)=fdh(iv,i+1)+w(j)*x(j)*gri(iv,j) 225 enddo ! of do iv=1,KDLON 226 enddo ! of do j=1,nq 227 enddo ! of do i=1,n 228 c....................................................................... 229 c apply the reflexion condition on the ground 230 c 231 do iv=1,KDLON 232 fah(iv,n+1)=(1.E+0-emis(iv))*fdh(iv,n+1)+pi*emis(iv)*bsol(iv) 233 enddo 234 do j=1,nq 235 do iv=1,KDLON 236 gri(iv,j)=2.E+0*fah(iv,n+1) 237 enddo 238 enddo 239 c....................................................................... 240 c going back up to compute all the upward fluxes 234 241 c 235 do 10011i=n,1,-1236 do 10012j=1,nq237 do 99984iv=1,KDLON238 gri(iv,j)=gri(iv,j)*exp(-tau(iv,i)/x(j))239 & +grg(iv,i)/(alambda(iv,i)*x(j)-1.E+0)240 & *(exp(-tau(iv,i)/x(j))-exp(-tau(iv,i)*alambda(iv,i)))241 & +grh(iv,i)/(alambda(iv,i)*x(j)+1.E+0)242 & *(1.E+0-exp(-tau(iv,i)*(alambda(iv,i)+1.E+0/x(j))))243 & +alpha1(iv,i)*(1.E+0-exp(-tau(iv,i)/x(j)))244 & +alpha2(iv,i)*(x(j)-(tau(iv,i)+x(j))*exp(-tau(iv,i)/x(j)))245 99984 continue 246 10012 continue 247 do 99983iv=1,KDLON248 fah(iv,i)=0.E+0249 99983 continue 250 do 10013j=1,nq251 do 99982iv=1,KDLON252 fah(iv,i)=fah(iv,i)+w(j)*x(j)*gri(iv,j)253 99982 continue 254 10013 continue 255 10011 continue 256 c....................................................................... 257 10100 continue 258 c....................................................................... 259 c 260 c....................................................................... 261 RETURN 262 END 242 do i=n,1,-1 243 do j=1,nq 244 do iv=1,KDLON 245 gri(iv,j)=gri(iv,j)*exp(-tau(iv,i)/x(j)) 246 & +grg(iv,i)/(alambda(iv,i)*x(j)-1.E+0) 247 & *(exp(-tau(iv,i)/x(j))-exp(-tau(iv,i)*alambda(iv,i))) 248 & +grh(iv,i)/(alambda(iv,i)*x(j)+1.E+0) 249 & *(1.E+0-exp(-tau(iv,i)*(alambda(iv,i)+1.E+0/x(j)))) 250 & +alpha1(iv,i)*(1.E+0-exp(-tau(iv,i)/x(j))) 251 & +alpha2(iv,i)*(x(j)-(tau(iv,i)+x(j))*exp(-tau(iv,i)/x(j))) 252 enddo ! of do iv=1,KDLON 253 enddo ! of do j=1,nq 254 do iv=1,KDLON 255 fah(iv,i)=0.E+0 256 enddo 257 do j=1,nq 258 do iv=1,KDLON 259 fah(iv,i)=fah(iv,i)+w(j)*x(j)*gri(iv,j) 260 enddo 261 enddo ! of do j=1,nq 262 enddo ! of do i=n,1,-1 263 c....................................................................... 264 ENDIF ! of IF (nsf.eq.0) 265 c....................................................................... 266 c 267 c....................................................................... 268 269 END SUBROUTINE flusv 263 270 264 271 c *************************************************************** … … 270 277 c....................................................................... 271 278 c 272 c inversion d'un systeme lineaire tridiagonal279 c solve a tridiagonal linear system such that: 273 280 c 274 281 c | b1 d1 | | y1 | | e1 | … … 278 285 c | an bn | | yn | | en | 279 286 c 280 c in : * KDLON --> dimension de vectorisation 281 c * n --> dimension du systeme 282 c * a,b,d,e --> voir dessin 283 c 284 c out : * y --> voir dessin 285 c 286 c....................................................................... 287 c declaration des arguments 288 c 289 INTEGER KDLON,n 290 REAL a(NDLO2,n),b(NDLO2,n),d(NDLO2,n),e(NDLO2,n),y(NDLO2,n) 291 c....................................................................... 292 c declaration des variables locales 293 c 294 INTEGER iv,i 295 REAL as(NDLON,4*nflev),ds(NDLON,4*nflev) 287 c in : * KDLON --> vectorisation dimension 288 c * n --> system size 289 c * a,b,d,e --> coefficients as shown above 290 c 291 c out : * y --> see above 292 c 293 c....................................................................... 294 c arguments 295 c 296 INTEGER,INTENT(IN) :: KDLON,n 297 REAL,INTENT(IN) :: a(NDLO2,n),b(NDLO2,n),d(NDLO2,n),e(NDLO2,n) 298 REAL,INTENT(OUT) :: y(NDLO2,n) 299 c....................................................................... 300 c local variables 301 c 302 INTEGER :: iv,i 303 REAL :: as(NDLON,4*nflev),ds(NDLON,4*nflev) 296 304 & ,x(NDLON,4*nflev) 297 305 c....................................................................... 298 306 c 299 307 c....................................................................... 300 do 99999 iv=1,KDLON 301 as(iv,n)=a(iv,n)/b(iv,n) 302 ds(iv,n)=e(iv,n)/b(iv,n) 303 99999 continue 304 do 10001 i=n-1,1,-1 305 do 99998 iv=1,KDLON 306 x(iv,i)=1.E+0/(b(iv,i)-d(iv,i)*as(iv,i+1)) 307 as(iv,i)=a(iv,i)*x(iv,i) 308 ds(iv,i)=(e(iv,i)-d(iv,i)*ds(iv,i+1))*x(iv,i) 309 99998 continue 310 10001 continue 311 do 99997 iv=1,KDLON 312 y(iv,1)=ds(iv,1) 313 99997 continue 314 do 10002 i=2,n 315 do 99996 iv=1,KDLON 316 y(iv,i)=ds(iv,i)-as(iv,i)*y(iv,i-1) 317 99996 continue 318 10002 continue 319 c....................................................................... 320 c 321 c....................................................................... 322 RETURN 323 END 308 do iv=1,KDLON 309 as(iv,n)=a(iv,n)/b(iv,n) 310 ds(iv,n)=e(iv,n)/b(iv,n) 311 enddo 312 do i=n-1,1,-1 313 do iv=1,KDLON 314 x(iv,i)=1.E+0/(b(iv,i)-d(iv,i)*as(iv,i+1)) 315 as(iv,i)=a(iv,i)*x(iv,i) 316 ds(iv,i)=(e(iv,i)-d(iv,i)*ds(iv,i+1))*x(iv,i) 317 enddo 318 enddo 319 do iv=1,KDLON 320 y(iv,1)=ds(iv,1) 321 enddo 322 do i=2,n 323 do iv=1,KDLON 324 y(iv,i)=ds(iv,i)-as(iv,i)*y(iv,i-1) 325 enddo 326 enddo 327 c....................................................................... 328 c 329 c....................................................................... 330 331 END SUBROUTINE sys3v 332 333 END MODULE flusv_mod -
trunk/LMDZ.MARS/libf/phymars/lwb.F
r1266 r3757 1 module lwb_mod 2 3 implicit none 4 5 contains 6 1 7 subroutine lwb (kdlon,kflev,tlev,tlay,dt0 2 8 . ,bsurf,btop,blay,blev,dblay,dbsublay) … … 132 138 c---------------------------------------------------------------------- 133 139 134 return 135 end 140 end subroutine lwb 141 142 end module lwb_mod -
trunk/LMDZ.MARS/libf/phymars/lwdiff.F
r3726 r3757 1 module lwdiff_mod 2 3 implicit none 4 5 contains 6 1 7 subroutine lwdiff (kdlon,kflev 2 8 . ,pbsur,pbtop,pdbsl … … 8 14 use yomlw_h, only: nlaylte 9 15 use comcstfi_h, only: pi 16 use flusv_mod, only: flusv 10 17 IMPLICIT NONE 11 18 … … 34 41 C --------- 35 42 C 36 integer kdlon,kflev 37 REAL PBSUR(NDLO2,nir), PBTOP(NDLO2,nir) 38 S , PDBSL(NDLO2,nir,KFLEV*2), PEMIS(NDLO2) 39 40 real PFLUC(NDLO2,2,KFLEV+1) 41 real tautotal(ndlon,nflev,nir) 42 real omegtotal(ndlon,nflev,nir), gtotal(ndlon,nflev,nir) 43 integer,intent(in) :: kdlon,kflev 44 real,intent(in) :: PBSUR(NDLO2,nir), PBTOP(NDLO2,nir) 45 real,intent(in) :: PDBSL(NDLO2,nir,KFLEV*2), PEMIS(NDLO2) 46 47 real,intent(out) :: PFLUC(NDLO2,2,KFLEV+1) 48 real,intent(in) :: tautotal(ndlon,nflev,nir) 49 real,intent(in) :: omegtotal(ndlon,nflev,nir) 50 real,intent(in) :: gtotal(ndlon,nflev,nir) 43 51 44 52 C … … 69 77 C -------------- 70 78 C 71 100 CONTINUE 79 72 80 C 73 81 C* 1.1 INITIALIZE LAYER CONTRIBUTIONS 74 82 C ------------------------------ 75 83 C 76 110 CONTINUE 84 77 85 C 78 86 … … 84 92 enddo 85 93 86 DO 112JK = 1 , nlaylte+187 DO 111JL = 1 , KDLON94 DO JK = 1 , nlaylte+1 95 DO JL = 1 , KDLON 88 96 ZADJD(JL,JK) = 0. 89 97 ZADJU(JL,JK) = 0. 90 98 ZDISD(JL,JK) = 0. 91 99 ZDISU(JL,JK) = 0. 92 111 CONTINUE93 112 CONTINUE100 ENDDO 101 ENDDO 94 102 C 95 103 C … … 106 114 C ================================================================== 107 115 C 108 200 CONTINUE 116 109 117 C 110 118 C ------------------------------------------------------------------ … … 133 141 c! boucle sur les bandes hors CO2 134 142 c! 135 DO 10001iir=3,nir143 DO iir=3,nir 136 144 c! 137 145 do jl=1,kdlon … … 202 210 enddo 203 211 ENDDO 204 10001 CONTINUE 212 ENDDO ! of DO iir=3,nir 205 213 206 214 DO J2=1,nlaylte+1 … … 213 221 214 222 215 END 223 END SUBROUTINE lwdiff 224 225 end module lwdiff_mod -
trunk/LMDZ.MARS/libf/phymars/lwflux.F
r3726 r3757 19 19 use dimradmars_mod, only: ndlo2, nir, ndlon, nuco2, nflev 20 20 use yomlw_h, only: nlaylte, xi, xi_ground, gcp 21 use lwdiff_mod, only: lwdiff 21 22 implicit none 22 23 -
trunk/LMDZ.MARS/libf/phymars/lwmain_mod.F
r3726 r3757 22 22 use yomlw_h, only: nlaylte, xi 23 23 use lwi_mod, only: lwi 24 use lwb_mod, only: lwb 25 use lwu_mod, only: lwu 26 use lwxd_mod, only: lwxd 27 use lwxn_mod, only: lwxn 28 use lwxb_mod, only: lwxb 24 29 use lwflux_mod, only: lwflux 25 30 use callkeys_mod, only: ilwb, ilwd, ilwn -
trunk/LMDZ.MARS/libf/phymars/lwtt.F
r1266 r3757 1 module lwtt_mod 2 3 implicit none 4 5 contains 6 1 7 subroutine lwtt (kdlon,u,up,nu,tr) 2 8 … … 59 65 60 66 c---------------------------------------------------------------------- 61 return 62 end 67 68 end subroutine lwtt 69 70 end module lwtt_mod -
trunk/LMDZ.MARS/libf/phymars/lwu.F
r3726 r3757 1 module lwu_mod 2 3 implicit none 4 5 contains 6 1 7 subroutine lwu (kdlon,kflev 2 8 & ,dp,plev,tlay,aerosol … … 208 214 209 215 c---------------------------------------------------------------------- 210 end 216 end subroutine lwu 217 218 end module lwu_mod -
trunk/LMDZ.MARS/libf/phymars/lwxb.F
r1266 r3757 1 module lwxb_mod 2 3 implicit none 4 5 contains 6 1 7 subroutine lwxb (ig0,kdlon,kflev 2 8 . ,emis … … 35 41 use dimradmars_mod, only: ndlo2, nuco2, ndlon, nflev 36 42 use yomlw_h, only: xi, nlaylte 43 use lwtt_mod, only: lwtt 37 44 implicit none 38 45 … … 212 219 213 220 c------------------------------------------------------------------------- 214 return 215 end 221 end subroutine lwxb 222 223 end module lwxb_mod -
trunk/LMDZ.MARS/libf/phymars/lwxd.F
r3726 r3757 1 module lwxd_mod 2 3 implicit none 4 5 contains 6 1 7 subroutine lwxd (ig0,kdlon,kflev,emis 2 8 . ,aer_t,co2_u,co2_up) … … 35 41 use yomlw_h, only: nlaylte, xi, xi_emis 36 42 use callkeys_mod, only: callemis 43 use lwtt_mod, only: lwtt 37 44 implicit none 38 45 … … 247 254 248 255 c---------------------------------------------------------------------- 249 end 256 end subroutine lwxd 257 258 end module lwxd_mod -
trunk/LMDZ.MARS/libf/phymars/lwxn.F
r3726 r3757 1 module lwxn_mod 2 3 implicit none 4 5 contains 6 1 7 subroutine lwxn ( ig0,kdlon,kflev 2 8 . , dp … … 73 79 use yomlw_h, only: nlaylte, xi, xi_ground, xi_emis 74 80 use callkeys_mod, only: linear, alphan, ncouche 81 use lwtt_mod, only: lwtt 75 82 implicit none 76 83 … … 379 386 enddo ! boucle sur jk 380 387 c---------------------------------------------------------------------- 381 return 382 end 383 384 385 388 389 end subroutine lwxn 390 391 end module lwxn_mod -
trunk/LMDZ.MARS/libf/phymars/sugwd.F90
r3754 r3757 1 MODULE sugwd_mod 2 3 IMPLICIT NONE 4 5 CONTAINS 6 1 7 SUBROUTINE SUGWD(nlayer,sigtest) 2 8 ! ============================================================================== … … 70 76 GTSEC=1.E-07 ! Security values for Sub-grid scale anisotropy(OROSETUP,GWSTRESS) 71 77 72 RETURN 73 END 78 END SUBROUTINE SUGWD 79 80 END MODULE sugwd_mod
Note: See TracChangeset
for help on using the changeset viewer.