Changeset 5246 for LMDZ6/trunk/libf/dynphy_lonlat
- Timestamp:
- Oct 21, 2024, 2:58:45 PM (7 weeks ago)
- Location:
- LMDZ6/trunk/libf/dynphy_lonlat
- Files:
-
- 6 moved
-
calfis.F90 (moved) (moved from LMDZ6/trunk/libf/dynphy_lonlat/calfis.F) (1 diff)
-
calfis_loc.F90 (moved) (moved from LMDZ6/trunk/libf/dynphy_lonlat/calfis_loc.F) (1 diff)
-
gr_dyn_fi.f90 (moved) (moved from LMDZ6/trunk/libf/dynphy_lonlat/gr_dyn_fi.F) (1 diff)
-
gr_dyn_fi_p.F90 (moved) (moved from LMDZ6/trunk/libf/dynphy_lonlat/gr_dyn_fi_p.F) (1 diff)
-
gr_fi_dyn.f90 (moved) (moved from LMDZ6/trunk/libf/dynphy_lonlat/gr_fi_dyn.F) (1 diff)
-
gr_fi_dyn_p.F90 (moved) (moved from LMDZ6/trunk/libf/dynphy_lonlat/gr_fi_dyn_p.F) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/trunk/libf/dynphy_lonlat/calfis.F90
r5245 r5246 2 2 ! $Id$ 3 3 ! 4 C 5 C 6 SUBROUTINE calfis(lafin, 7 $ jD_cur, jH_cur,8 $ pucov,9 $ pvcov,10 $ pteta,11 $ pq,12 $ pmasse,13 $ pps,14 $ pp,15 $ ppk,16 $ pphis,17 $ pphi,18 $ pducov,19 $ pdvcov,20 $ pdteta,21 $ pdq,22 $ flxw,23 $ pdufi,24 $ pdvfi,25 $ pdhfi,26 $ pdqfi,27 $pdpsfi)28 c 29 c Auteur : P. Le Van, F. Hourdin 30 c.........31 32 4 ! 5 ! 6 SUBROUTINE calfis(lafin, & 7 jD_cur, jH_cur, & 8 pucov, & 9 pvcov, & 10 pteta, & 11 pq, & 12 pmasse, & 13 pps, & 14 pp, & 15 ppk, & 16 pphis, & 17 pphi, & 18 pducov, & 19 pdvcov, & 20 pdteta, & 21 pdq, & 22 flxw, & 23 pdufi, & 24 pdvfi, & 25 pdhfi, & 26 pdqfi, & 27 pdpsfi) 28 ! 29 ! Auteur : P. Le Van, F. Hourdin 30 ! ......... 31 USE infotrac, ONLY: nqtot, tracers 32 USE control_mod, ONLY: planet_type, nsplit_phys 33 33 #ifdef CPP_PHYS 34 USE callphysiq_mod, ONLY: call_physiq 35 #endif 36 USE comconst_mod, ONLY: cpp, daysec, dtphys, dtvr, kappa, pi 37 USE comvert_mod, ONLY: preff, presnivs 38 39 IMPLICIT NONE 40 c======================================================================= 41 c 42 c 1. rearrangement des tableaux et transformation 43 c variables dynamiques > variables physiques 44 c 2. calcul des termes physiques 45 c 3. retransformation des tendances physiques en tendances dynamiques 46 c 47 c remarques: 48 c ---------- 49 c 50 c - les vents sont donnes dans la physique par leurs composantes 51 c naturelles. 52 c - la variable thermodynamique de la physique est une variable 53 c intensive : T 54 c pour la dynamique on prend T * ( preff / p(l) ) **kappa 55 c - les deux seules variables dependant de la geometrie necessaires 56 c pour la physique sont la latitude pour le rayonnement et 57 c l'aire de la maille quand on veut integrer une grandeur 58 c horizontalement. 59 c - les points de la physique sont les points scalaires de la 60 c la dynamique; numerotation: 61 c 1 pour le pole nord 62 c (jjm-1)*iim pour l'interieur du domaine 63 c ngridmx pour le pole sud 64 c ---> ngridmx=2+(jjm-1)*iim 65 c 66 c Input : 67 c ------- 68 c pucov covariant zonal velocity 69 c pvcov covariant meridional velocity 70 c pteta potential temperature 71 c pps surface pressure 72 c pmasse masse d'air dans chaque maille 73 c pts surface temperature (K) 74 c callrad clef d'appel au rayonnement 75 c 76 c Output : 77 c -------- 78 c pdufi tendency for the natural zonal velocity (ms-1) 79 c pdvfi tendency for the natural meridional velocity 80 c pdhfi tendency for the potential temperature 81 c pdtsfi tendency for the surface temperature 82 c 83 c pdtrad radiative tendencies \ both input 84 c pfluxrad radiative fluxes / and output 85 c 86 c======================================================================= 87 c 88 c----------------------------------------------------------------------- 89 c 90 c 0. Declarations : 91 c ------------------ 92 93 include "dimensions.h" 94 include "paramet.h" 95 96 INTEGER ngridmx 97 PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm ) 98 99 include "comgeom2.h" 100 include "iniprint.h" 101 102 c Arguments : 103 c ----------- 104 LOGICAL,INTENT(IN) :: lafin ! .true. for the very last call to physics 105 REAL,INTENT(IN):: jD_cur, jH_cur 106 REAL,INTENT(IN) :: pvcov(iip1,jjm,llm) ! covariant meridional velocity 107 REAL,INTENT(IN) :: pucov(iip1,jjp1,llm) ! covariant zonal velocity 108 REAL,INTENT(IN) :: pteta(iip1,jjp1,llm) ! potential temperature 109 REAL,INTENT(IN) :: pmasse(iip1,jjp1,llm) ! mass in each cell ! not used 110 REAL,INTENT(IN) :: pq(iip1,jjp1,llm,nqtot) ! tracers 111 REAL,INTENT(IN) :: pphis(iip1,jjp1) ! surface geopotential 112 REAL,INTENT(IN) :: pphi(iip1,jjp1,llm) ! geopotential 113 114 REAL,INTENT(IN) :: pdvcov(iip1,jjm,llm) ! dynamical tendency on vcov 115 REAL,INTENT(IN) :: pducov(iip1,jjp1,llm) ! dynamical tendency on ucov 116 REAL,INTENT(IN) :: pdteta(iip1,jjp1,llm) ! dynamical tendency on teta 117 ! NB: pdteta is used only to compute pcvgt which is in fact not used... 118 REAL,INTENT(IN) :: pdq(iip1,jjp1,llm,nqtot) ! dynamical tendency on tracers 119 ! NB: pdq is only used to compute pcvgq which is in fact not used... 120 121 REAL,INTENT(IN) :: pps(iip1,jjp1) ! surface pressure (Pa) 122 REAL,INTENT(IN) :: pp(iip1,jjp1,llmp1) ! pressure at mesh interfaces (Pa) 123 REAL,INTENT(IN) :: ppk(iip1,jjp1,llm) ! Exner at mid-layer 124 REAL,INTENT(IN) :: flxw(iip1,jjp1,llm) ! Vertical mass flux on lower mesh interfaces (kg/s) (on llm because flxw(:,:,llm+1)=0) 125 126 ! tendencies (in */s) from the physics 127 REAL,INTENT(OUT) :: pdvfi(iip1,jjm,llm) ! tendency on covariant meridional wind 128 REAL,INTENT(OUT) :: pdufi(iip1,jjp1,llm) ! tendency on covariant zonal wind 129 REAL,INTENT(OUT) :: pdhfi(iip1,jjp1,llm) ! tendency on potential temperature (K/s) 130 REAL,INTENT(OUT) :: pdqfi(iip1,jjp1,llm,nqtot) ! tendency on tracers 131 REAL,INTENT(OUT) :: pdpsfi(iip1,jjp1) ! tendency on surface pressure (Pa/s) 132 133 134 c Local variables : 135 c ----------------- 136 137 INTEGER i,j,l,ig0,ig,iq,itr 138 REAL zpsrf(ngridmx) 139 REAL zplev(ngridmx,llm+1),zplay(ngridmx,llm) 140 REAL zphi(ngridmx,llm),zphis(ngridmx) 141 c 142 REAL zrot(iip1,jjm,llm) ! AdlC May 2014 143 REAL zufi(ngridmx,llm), zvfi(ngridmx,llm) 144 REAL zrfi(ngridmx,llm) ! relative wind vorticity 145 REAL ztfi(ngridmx,llm),zqfi(ngridmx,llm,nqtot) 146 REAL zpk(ngridmx,llm) 147 c 148 REAL pcvgu(ngridmx,llm), pcvgv(ngridmx,llm) 149 REAL pcvgt(ngridmx,llm), pcvgq(ngridmx,llm,2) 150 c 151 REAL zdufi(ngridmx,llm),zdvfi(ngridmx,llm) 152 REAL zdtfi(ngridmx,llm),zdqfi(ngridmx,llm,nqtot) 153 REAL zdpsrf(ngridmx) 154 c 155 REAL zdufic(ngridmx,llm),zdvfic(ngridmx,llm) 156 REAL zdtfic(ngridmx,llm),zdqfic(ngridmx,llm,nqtot) 157 REAL jH_cur_split,zdt_split 158 LOGICAL debut_split,lafin_split 159 INTEGER isplit 160 161 REAL zsin(iim),zcos(iim),z1(iim) 162 REAL zsinbis(iim),zcosbis(iim),z1bis(iim) 163 REAL unskap, pksurcp 164 c 165 REAL flxwfi(ngridmx,llm) ! Flux de masse verticale sur la grille physiq 166 c 167 168 REAL SSUM 169 170 LOGICAL,SAVE :: firstcal=.true., debut=.true. 171 ! REAL rdayvrai 172 173 c 174 c----------------------------------------------------------------------- 175 c 176 c 1. Initialisations : 177 c -------------------- 178 c 179 c 180 IF ( firstcal ) THEN 181 debut = .TRUE. 182 IF (ngridmx.NE.2+(jjm-1)*iim) THEN 183 write(lunout,*) 'STOP dans calfis' 184 write(lunout,*) 185 & 'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim' 186 write(lunout,*) ' ngridmx jjm iim ' 187 write(lunout,*) ngridmx,jjm,iim 188 call abort_gcm("calfis", "", 1) 189 ENDIF 190 ELSE 191 debut = .FALSE. 192 ENDIF ! of IF (firstcal) 193 194 c 195 c 196 c----------------------------------------------------------------------- 197 c 40. transformation des variables dynamiques en variables physiques: 198 c --------------------------------------------------------------- 199 200 c 41. pressions au sol (en Pascals) 201 c ---------------------------------- 202 203 204 zpsrf(1) = pps(1,1) 205 206 ig0 = 2 207 DO j = 2,jjm 208 CALL SCOPY( iim,pps(1,j),1,zpsrf(ig0), 1 ) 209 ig0 = ig0+iim 210 ENDDO 211 212 zpsrf(ngridmx) = pps(1,jjp1) 213 214 215 c 42. pression intercouches et fonction d'Exner: 216 c 217 c ----------------------------------------------------------------- 218 c .... zplev definis aux (llm +1) interfaces des couches .... 219 c .... zplay definis aux ( llm ) milieux des couches .... 220 c ----------------------------------------------------------------- 221 222 c ... Exner = cp * ( p(l) / preff ) ** kappa .... 223 c 224 unskap = 1./ kappa 225 c 226 DO l = 1, llm 227 zpk( 1,l ) = ppk(1,1,l) 228 zplev( 1,l ) = pp(1,1,l) 229 ig0 = 2 230 DO j = 2, jjm 231 DO i =1, iim 232 zpk( ig0,l ) = ppk(i,j,l) 233 zplev( ig0,l ) = pp(i,j,l) 234 ig0 = ig0 +1 235 ENDDO 236 ENDDO 237 zpk( ngridmx,l ) = ppk(1,jjp1,l) 238 zplev( ngridmx,l ) = pp(1,jjp1,l) 239 ENDDO 240 zplev( 1,llmp1 ) = pp(1,1,llmp1) 241 ig0 = 2 242 DO j = 2, jjm 243 DO i =1, iim 244 zplev( ig0,llmp1 ) = pp(i,j,llmp1) 245 ig0 = ig0 +1 246 ENDDO 247 ENDDO 248 zplev( ngridmx,llmp1 ) = pp(1,jjp1,llmp1) 249 c 250 c 251 252 c 43. temperature naturelle (en K) et pressions milieux couches . 253 c --------------------------------------------------------------- 254 255 DO l=1,llm 256 257 pksurcp = ppk(1,1,l) / cpp 258 zplay(1,l) = preff * pksurcp ** unskap 259 ztfi(1,l) = pteta(1,1,l) * pksurcp 260 pcvgt(1,l) = pdteta(1,1,l) * pksurcp / pmasse(1,1,l) 261 ig0 = 2 262 263 DO j = 2, jjm 264 DO i = 1, iim 265 pksurcp = ppk(i,j,l) / cpp 266 zplay(ig0,l) = preff * pksurcp ** unskap 267 ztfi(ig0,l) = pteta(i,j,l) * pksurcp 268 pcvgt(ig0,l) = pdteta(i,j,l) * pksurcp / pmasse(i,j,l) 269 ig0 = ig0 + 1 270 ENDDO 271 ENDDO 272 273 pksurcp = ppk(1,jjp1,l) / cpp 274 zplay(ig0,l) = preff * pksurcp ** unskap 275 ztfi (ig0,l) = pteta(1,jjp1,l) * pksurcp 276 pcvgt(ig0,l) = pdteta(1,jjp1,l) * pksurcp/ pmasse(1,jjp1,l) 277 278 ENDDO 279 280 c 43.bis traceurs 281 c --------------- 282 c 283 itr=0 284 DO iq=1,nqtot 285 IF(.NOT.tracers(iq)%isAdvected) CYCLE 286 itr = itr + 1 287 DO l=1,llm 288 zqfi(1,l,itr) = pq(1,1,l,iq) 289 ig0 = 2 290 DO j=2,jjm 291 DO i = 1, iim 292 zqfi(ig0,l,itr) = pq(i,j,l,iq) 293 ig0 = ig0 + 1 294 ENDDO 295 ENDDO 296 zqfi(ig0,l,itr) = pq(1,jjp1,l,iq) 34 USE callphysiq_mod, ONLY: call_physiq 35 #endif 36 USE comconst_mod, ONLY: cpp, daysec, dtphys, dtvr, kappa, pi 37 USE comvert_mod, ONLY: preff, presnivs 38 39 IMPLICIT NONE 40 !======================================================================= 41 ! 42 ! 1. rearrangement des tableaux et transformation 43 ! variables dynamiques > variables physiques 44 ! 2. calcul des termes physiques 45 ! 3. retransformation des tendances physiques en tendances dynamiques 46 ! 47 ! remarques: 48 ! ---------- 49 ! 50 ! - les vents sont donnes dans la physique par leurs composantes 51 ! naturelles. 52 ! - la variable thermodynamique de la physique est une variable 53 ! intensive : T 54 ! pour la dynamique on prend T * ( preff / p(l) ) **kappa 55 ! - les deux seules variables dependant de la geometrie necessaires 56 ! pour la physique sont la latitude pour le rayonnement et 57 ! l'aire de la maille quand on veut integrer une grandeur 58 ! horizontalement. 59 ! - les points de la physique sont les points scalaires de la 60 ! la dynamique; numerotation: 61 ! 1 pour le pole nord 62 ! (jjm-1)*iim pour l'interieur du domaine 63 ! ngridmx pour le pole sud 64 ! ---> ngridmx=2+(jjm-1)*iim 65 ! 66 ! Input : 67 ! ------- 68 ! pucov covariant zonal velocity 69 ! pvcov covariant meridional velocity 70 ! pteta potential temperature 71 ! pps surface pressure 72 ! pmasse masse d'air dans chaque maille 73 ! pts surface temperature (K) 74 ! callrad clef d'appel au rayonnement 75 ! 76 ! Output : 77 ! -------- 78 ! pdufi tendency for the natural zonal velocity (ms-1) 79 ! pdvfi tendency for the natural meridional velocity 80 ! pdhfi tendency for the potential temperature 81 ! pdtsfi tendency for the surface temperature 82 ! 83 ! pdtrad radiative tendencies \ both input 84 ! pfluxrad radiative fluxes / and output 85 ! 86 !======================================================================= 87 ! 88 !----------------------------------------------------------------------- 89 ! 90 ! 0. Declarations : 91 ! ------------------ 92 93 include "dimensions.h" 94 include "paramet.h" 95 96 INTEGER :: ngridmx 97 PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm ) 98 99 include "comgeom2.h" 100 include "iniprint.h" 101 102 ! Arguments : 103 ! ----------- 104 LOGICAL,INTENT(IN) :: lafin ! .true. for the very last call to physics 105 REAL,INTENT(IN):: jD_cur, jH_cur 106 REAL,INTENT(IN) :: pvcov(iip1,jjm,llm) ! covariant meridional velocity 107 REAL,INTENT(IN) :: pucov(iip1,jjp1,llm) ! covariant zonal velocity 108 REAL,INTENT(IN) :: pteta(iip1,jjp1,llm) ! potential temperature 109 REAL,INTENT(IN) :: pmasse(iip1,jjp1,llm) ! mass in each cell ! not used 110 REAL,INTENT(IN) :: pq(iip1,jjp1,llm,nqtot) ! tracers 111 REAL,INTENT(IN) :: pphis(iip1,jjp1) ! surface geopotential 112 REAL,INTENT(IN) :: pphi(iip1,jjp1,llm) ! geopotential 113 114 REAL,INTENT(IN) :: pdvcov(iip1,jjm,llm) ! dynamical tendency on vcov 115 REAL,INTENT(IN) :: pducov(iip1,jjp1,llm) ! dynamical tendency on ucov 116 REAL,INTENT(IN) :: pdteta(iip1,jjp1,llm) ! dynamical tendency on teta 117 ! ! NB: pdteta is used only to compute pcvgt which is in fact not used... 118 REAL,INTENT(IN) :: pdq(iip1,jjp1,llm,nqtot) ! dynamical tendency on tracers 119 ! ! NB: pdq is only used to compute pcvgq which is in fact not used... 120 121 REAL,INTENT(IN) :: pps(iip1,jjp1) ! surface pressure (Pa) 122 REAL,INTENT(IN) :: pp(iip1,jjp1,llmp1) ! pressure at mesh interfaces (Pa) 123 REAL,INTENT(IN) :: ppk(iip1,jjp1,llm) ! Exner at mid-layer 124 REAL,INTENT(IN) :: flxw(iip1,jjp1,llm) ! Vertical mass flux on lower mesh interfaces (kg/s) (on llm because flxw(:,:,llm+1)=0) 125 126 ! ! tendencies (in */s) from the physics 127 REAL,INTENT(OUT) :: pdvfi(iip1,jjm,llm) ! tendency on covariant meridional wind 128 REAL,INTENT(OUT) :: pdufi(iip1,jjp1,llm) ! tendency on covariant zonal wind 129 REAL,INTENT(OUT) :: pdhfi(iip1,jjp1,llm) ! tendency on potential temperature (K/s) 130 REAL,INTENT(OUT) :: pdqfi(iip1,jjp1,llm,nqtot) ! tendency on tracers 131 REAL,INTENT(OUT) :: pdpsfi(iip1,jjp1) ! tendency on surface pressure (Pa/s) 132 133 134 ! Local variables : 135 ! ----------------- 136 137 INTEGER :: i,j,l,ig0,ig,iq,itr 138 REAL :: zpsrf(ngridmx) 139 REAL :: zplev(ngridmx,llm+1),zplay(ngridmx,llm) 140 REAL :: zphi(ngridmx,llm),zphis(ngridmx) 141 ! 142 REAL :: zrot(iip1,jjm,llm) ! AdlC May 2014 143 REAL :: zufi(ngridmx,llm), zvfi(ngridmx,llm) 144 REAL :: zrfi(ngridmx,llm) ! relative wind vorticity 145 REAL :: ztfi(ngridmx,llm),zqfi(ngridmx,llm,nqtot) 146 REAL :: zpk(ngridmx,llm) 147 ! 148 REAL :: pcvgu(ngridmx,llm), pcvgv(ngridmx,llm) 149 REAL :: pcvgt(ngridmx,llm), pcvgq(ngridmx,llm,2) 150 ! 151 REAL :: zdufi(ngridmx,llm),zdvfi(ngridmx,llm) 152 REAL :: zdtfi(ngridmx,llm),zdqfi(ngridmx,llm,nqtot) 153 REAL :: zdpsrf(ngridmx) 154 ! 155 REAL :: zdufic(ngridmx,llm),zdvfic(ngridmx,llm) 156 REAL :: zdtfic(ngridmx,llm),zdqfic(ngridmx,llm,nqtot) 157 REAL :: jH_cur_split,zdt_split 158 LOGICAL :: debut_split,lafin_split 159 INTEGER :: isplit 160 161 REAL :: zsin(iim),zcos(iim),z1(iim) 162 REAL :: zsinbis(iim),zcosbis(iim),z1bis(iim) 163 REAL :: unskap, pksurcp 164 ! 165 REAL :: flxwfi(ngridmx,llm) ! Flux de masse verticale sur la grille physiq 166 ! 167 168 REAL :: SSUM 169 170 LOGICAL,SAVE :: firstcal=.true., debut=.true. 171 ! REAL rdayvrai 172 173 ! 174 !----------------------------------------------------------------------- 175 ! 176 ! 1. Initialisations : 177 ! -------------------- 178 ! 179 ! 180 IF ( firstcal ) THEN 181 debut = .TRUE. 182 IF (ngridmx.NE.2+(jjm-1)*iim) THEN 183 write(lunout,*) 'STOP dans calfis' 184 write(lunout,*) & 185 'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim' 186 write(lunout,*) ' ngridmx jjm iim ' 187 write(lunout,*) ngridmx,jjm,iim 188 call abort_gcm("calfis", "", 1) 189 ENDIF 190 ELSE 191 debut = .FALSE. 192 ENDIF ! of IF (firstcal) 193 194 ! 195 ! 196 !----------------------------------------------------------------------- 197 ! 40. transformation des variables dynamiques en variables physiques: 198 ! --------------------------------------------------------------- 199 200 ! 41. pressions au sol (en Pascals) 201 ! ---------------------------------- 202 203 204 zpsrf(1) = pps(1,1) 205 206 ig0 = 2 207 DO j = 2,jjm 208 CALL SCOPY( iim,pps(1,j),1,zpsrf(ig0), 1 ) 209 ig0 = ig0+iim 210 ENDDO 211 212 zpsrf(ngridmx) = pps(1,jjp1) 213 214 215 ! 42. pression intercouches et fonction d'Exner: 216 ! 217 ! ----------------------------------------------------------------- 218 ! .... zplev definis aux (llm +1) interfaces des couches .... 219 ! .... zplay definis aux ( llm ) milieux des couches .... 220 ! ----------------------------------------------------------------- 221 222 ! ... Exner = cp * ( p(l) / preff ) ** kappa .... 223 ! 224 unskap = 1./ kappa 225 ! 226 DO l = 1, llm 227 zpk( 1,l ) = ppk(1,1,l) 228 zplev( 1,l ) = pp(1,1,l) 229 ig0 = 2 230 DO j = 2, jjm 231 DO i =1, iim 232 zpk( ig0,l ) = ppk(i,j,l) 233 zplev( ig0,l ) = pp(i,j,l) 234 ig0 = ig0 +1 297 235 ENDDO 298 236 ENDDO 299 300 c convergence dynamique pour les traceurs "EAU" 301 ! Earth-specific treatment of first 2 tracers (water) 302 if (planet_type=="earth") then 303 DO iq=1,2 304 DO l=1,llm 305 pcvgq(1,l,iq)= pdq(1,1,l,iq) / pmasse(1,1,l) 306 ig0 = 2 307 DO j=2,jjm 308 DO i = 1, iim 309 pcvgq(ig0,l,iq) = pdq(i,j,l,iq) / pmasse(i,j,l) 310 ig0 = ig0 + 1 311 ENDDO 312 ENDDO 313 pcvgq(ig0,l,iq)= pdq(1,jjp1,l,iq) / pmasse(1,jjp1,l) 314 ENDDO 315 ENDDO 316 endif ! of if (planet_type=="earth") 317 318 319 c Geopotentiel calcule par rapport a la surface locale: 320 c ----------------------------------------------------- 321 322 CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,pphi,zphi) 323 CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,pphis,zphis) 324 DO l=1,llm 325 DO ig=1,ngridmx 326 zphi(ig,l)=zphi(ig,l)-zphis(ig) 237 zpk( ngridmx,l ) = ppk(1,jjp1,l) 238 zplev( ngridmx,l ) = pp(1,jjp1,l) 239 ENDDO 240 zplev( 1,llmp1 ) = pp(1,1,llmp1) 241 ig0 = 2 242 DO j = 2, jjm 243 DO i =1, iim 244 zplev( ig0,llmp1 ) = pp(i,j,llmp1) 245 ig0 = ig0 +1 327 246 ENDDO 328 247 ENDDO 329 330 c .... Calcul de la vitesse verticale ( en Pa*m*s ou Kg/s ) .... 331 c JG : ancien calcule de omega utilise dans physiq.F. Maintenant le flux 332 c de masse est calclue dans advtrac.F 333 c DO l=1,llm 334 c pvervel(1,l)=pw(1,1,l) * g /apoln 335 c ig0=2 336 c DO j=2,jjm 337 c DO i = 1, iim 338 c pvervel(ig0,l) = pw(i,j,l) * g * unsaire(i,j) 339 c ig0 = ig0 + 1 340 c ENDDO 341 c ENDDO 342 c pvervel(ig0,l)=pw(1,jjp1,l) * g /apols 343 c ENDDO 344 345 c 346 c 45. champ u: 347 c ------------ 348 349 DO 50 l=1,llm 350 351 DO 25 j=2,jjm 352 ig0 = 1+(j-2)*iim 353 zufi(ig0+1,l)= 0.5 * 354 $ ( pucov(iim,j,l)/cu(iim,j) + pucov(1,j,l)/cu(1,j) ) 355 pcvgu(ig0+1,l)= 0.5 * 356 $ ( pducov(iim,j,l)/cu(iim,j) + pducov(1,j,l)/cu(1,j) ) 357 DO 10 i=2,iim 358 zufi(ig0+i,l)= 0.5 * 359 $ ( pucov(i-1,j,l)/cu(i-1,j) + pucov(i,j,l)/cu(i,j) ) 360 pcvgu(ig0+i,l)= 0.5 * 361 $ ( pducov(i-1,j,l)/cu(i-1,j) + pducov(i,j,l)/cu(i,j) ) 362 10 CONTINUE 363 25 CONTINUE 364 365 50 CONTINUE 366 367 368 C Alvaro de la Camara (May 2014) 369 C 46.1 Calcul de la vorticite et passage sur la grille physique 370 C -------------------------------------------------------------- 371 DO l=1,llm 372 do i=1,iim 373 do j=1,jjm 374 zrot(i,j,l) = (pvcov(i+1,j,l) - pvcov(i,j,l) 375 $ + pucov(i,j+1,l) - pucov(i,j,l)) 376 $ / (cu(i,j)+cu(i,j+1)) 377 $ / (cv(i+1,j)+cv(i,j)) *4 378 enddo 379 enddo 248 zplev( ngridmx,llmp1 ) = pp(1,jjp1,llmp1) 249 ! 250 ! 251 252 ! 43. temperature naturelle (en K) et pressions milieux couches . 253 ! --------------------------------------------------------------- 254 255 DO l=1,llm 256 257 pksurcp = ppk(1,1,l) / cpp 258 zplay(1,l) = preff * pksurcp ** unskap 259 ztfi(1,l) = pteta(1,1,l) * pksurcp 260 pcvgt(1,l) = pdteta(1,1,l) * pksurcp / pmasse(1,1,l) 261 ig0 = 2 262 263 DO j = 2, jjm 264 DO i = 1, iim 265 pksurcp = ppk(i,j,l) / cpp 266 zplay(ig0,l) = preff * pksurcp ** unskap 267 ztfi(ig0,l) = pteta(i,j,l) * pksurcp 268 pcvgt(ig0,l) = pdteta(i,j,l) * pksurcp / pmasse(i,j,l) 269 ig0 = ig0 + 1 270 ENDDO 271 ENDDO 272 273 pksurcp = ppk(1,jjp1,l) / cpp 274 zplay(ig0,l) = preff * pksurcp ** unskap 275 ztfi (ig0,l) = pteta(1,jjp1,l) * pksurcp 276 pcvgt(ig0,l) = pdteta(1,jjp1,l) * pksurcp/ pmasse(1,jjp1,l) 277 278 ENDDO 279 280 ! 43.bis traceurs 281 ! --------------- 282 ! 283 itr=0 284 DO iq=1,nqtot 285 IF(.NOT.tracers(iq)%isAdvected) CYCLE 286 itr = itr + 1 287 DO l=1,llm 288 zqfi(1,l,itr) = pq(1,1,l,iq) 289 ig0 = 2 290 DO j=2,jjm 291 DO i = 1, iim 292 zqfi(ig0,l,itr) = pq(i,j,l,iq) 293 ig0 = ig0 + 1 294 ENDDO 295 ENDDO 296 zqfi(ig0,l,itr) = pq(1,jjp1,l,iq) 297 ENDDO 298 ENDDO 299 300 ! convergence dynamique pour les traceurs "EAU" 301 ! Earth-specific treatment of first 2 tracers (water) 302 if (planet_type=="earth") then 303 DO iq=1,2 304 DO l=1,llm 305 pcvgq(1,l,iq)= pdq(1,1,l,iq) / pmasse(1,1,l) 306 ig0 = 2 307 DO j=2,jjm 308 DO i = 1, iim 309 pcvgq(ig0,l,iq) = pdq(i,j,l,iq) / pmasse(i,j,l) 310 ig0 = ig0 + 1 311 ENDDO 312 ENDDO 313 pcvgq(ig0,l,iq)= pdq(1,jjp1,l,iq) / pmasse(1,jjp1,l) 314 ENDDO 315 ENDDO 316 endif ! of if (planet_type=="earth") 317 318 319 ! Geopotentiel calcule par rapport a la surface locale: 320 ! ----------------------------------------------------- 321 322 CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,pphi,zphi) 323 CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,pphis,zphis) 324 DO l=1,llm 325 DO ig=1,ngridmx 326 zphi(ig,l)=zphi(ig,l)-zphis(ig) 327 ENDDO 328 ENDDO 329 330 ! .... Calcul de la vitesse verticale ( en Pa*m*s ou Kg/s ) .... 331 ! JG : ancien calcule de omega utilise dans physiq.F. Maintenant le flux 332 ! de masse est calclue dans advtrac.F 333 ! DO l=1,llm 334 ! pvervel(1,l)=pw(1,1,l) * g /apoln 335 ! ig0=2 336 ! DO j=2,jjm 337 ! DO i = 1, iim 338 ! pvervel(ig0,l) = pw(i,j,l) * g * unsaire(i,j) 339 ! ig0 = ig0 + 1 340 ! ENDDO 341 ! ENDDO 342 ! pvervel(ig0,l)=pw(1,jjp1,l) * g /apols 343 ! ENDDO 344 345 ! 346 ! 45. champ u: 347 ! ------------ 348 349 DO l=1,llm 350 351 DO j=2,jjm 352 ig0 = 1+(j-2)*iim 353 zufi(ig0+1,l)= 0.5 * & 354 ( pucov(iim,j,l)/cu(iim,j) + pucov(1,j,l)/cu(1,j) ) 355 pcvgu(ig0+1,l)= 0.5 * & 356 ( pducov(iim,j,l)/cu(iim,j) + pducov(1,j,l)/cu(1,j) ) 357 DO i=2,iim 358 zufi(ig0+i,l)= 0.5 * & 359 ( pucov(i-1,j,l)/cu(i-1,j) + pucov(i,j,l)/cu(i,j) ) 360 pcvgu(ig0+i,l)= 0.5 * & 361 ( pducov(i-1,j,l)/cu(i-1,j) + pducov(i,j,l)/cu(i,j) ) 362 END DO 363 END DO 364 365 END DO 366 367 368 ! Alvaro de la Camara (May 2014) 369 ! 46.1 Calcul de la vorticite et passage sur la grille physique 370 ! -------------------------------------------------------------- 371 DO l=1,llm 372 do i=1,iim 373 do j=1,jjm 374 zrot(i,j,l) = (pvcov(i+1,j,l) - pvcov(i,j,l) & 375 + pucov(i,j+1,l) - pucov(i,j,l)) & 376 / (cu(i,j)+cu(i,j+1)) & 377 / (cv(i+1,j)+cv(i,j)) *4 378 enddo 379 enddo 380 ENDDO 381 382 ! 46.champ v: 383 ! ----------- 384 385 DO l=1,llm 386 DO j=2,jjm 387 ig0=1+(j-2)*iim 388 DO i=1,iim 389 zvfi(ig0+i,l)= 0.5 * & 390 ( pvcov(i,j-1,l)/cv(i,j-1) + pvcov(i,j,l)/cv(i,j) ) 391 pcvgv(ig0+i,l)= 0.5 * & 392 ( pdvcov(i,j-1,l)/cv(i,j-1) + pdvcov(i,j,l)/cv(i,j) ) 393 ENDDO 394 zrfi(ig0 + 1,l)= 0.25 *(zrot(iim,j-1,l)+zrot(iim,j,l) & 395 +zrot(1,j-1,l)+zrot(1,j,l)) 396 DO i=2,iim 397 zrfi(ig0 + i,l)= 0.25 *(zrot(i-1,j-1,l)+zrot(i-1,j,l) & 398 +zrot(i,j-1,l)+zrot(i,j,l)) ! AdlC MAY 2014 399 ENDDO 400 ENDDO 401 ENDDO 402 403 404 ! 47. champs de vents aux pole nord 405 ! ------------------------------ 406 ! U = 1 / pi * integrale [ v * cos(long) * d long ] 407 ! V = 1 / pi * integrale [ v * sin(long) * d long ] 408 409 DO l=1,llm 410 411 z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,1,l)/cv(1,1) 412 z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,1,l)/cv(1,1) 413 DO i=2,iim 414 z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1) 415 z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,1,l)/cv(i,1) 416 ENDDO 417 418 DO i=1,iim 419 zcos(i) = COS(rlonv(i))*z1(i) 420 zcosbis(i)= COS(rlonv(i))*z1bis(i) 421 zsin(i) = SIN(rlonv(i))*z1(i) 422 zsinbis(i)= SIN(rlonv(i))*z1bis(i) 423 ENDDO 424 425 zufi(1,l) = SSUM(iim,zcos,1)/pi 426 pcvgu(1,l) = SSUM(iim,zcosbis,1)/pi 427 zvfi(1,l) = SSUM(iim,zsin,1)/pi 428 pcvgv(1,l) = SSUM(iim,zsinbis,1)/pi 429 zrfi(1, l) = 0. 430 ENDDO 431 432 433 ! 48. champs de vents aux pole sud: 434 ! --------------------------------- 435 ! U = 1 / pi * integrale [ v * cos(long) * d long ] 436 ! V = 1 / pi * integrale [ v * sin(long) * d long ] 437 438 DO l=1,llm 439 440 z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,jjm,l)/cv(1,jjm) 441 z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,jjm,l)/cv(1,jjm) 442 DO i=2,iim 443 z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm) 444 z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,jjm,l)/cv(i,jjm) 445 ENDDO 446 447 DO i=1,iim 448 zcos(i) = COS(rlonv(i))*z1(i) 449 zcosbis(i) = COS(rlonv(i))*z1bis(i) 450 zsin(i) = SIN(rlonv(i))*z1(i) 451 zsinbis(i) = SIN(rlonv(i))*z1bis(i) 452 ENDDO 453 454 zufi(ngridmx,l) = SSUM(iim,zcos,1)/pi 455 pcvgu(ngridmx,l) = SSUM(iim,zcosbis,1)/pi 456 zvfi(ngridmx,l) = SSUM(iim,zsin,1)/pi 457 pcvgv(ngridmx,l) = SSUM(iim,zsinbis,1)/pi 458 zrfi(ngridmx, l) = 0. 459 ENDDO 460 ! 461 ! On change de grille, dynamique vers physiq, pour le flux de masse verticale 462 CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,flxw,flxwfi) 463 464 !----------------------------------------------------------------------- 465 ! Appel de la physique: 466 ! --------------------- 467 468 469 470 ! write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys 471 zdt_split=dtphys/nsplit_phys 472 zdufic(:,:)=0. 473 zdvfic(:,:)=0. 474 zdtfic(:,:)=0. 475 zdqfic(:,:,:)=0. 476 477 #ifdef CPP_PHYS 478 479 do isplit=1,nsplit_phys 480 481 jH_cur_split=jH_cur+(isplit-1) * dtvr / (daysec *nsplit_phys) 482 debut_split=debut.and.isplit==1 483 lafin_split=lafin.and.isplit==nsplit_phys 484 485 ! if (planet_type=="earth") then 486 CALL call_physiq(ngridmx,llm,nqtot,tracers(:)%name, & 487 debut_split,lafin_split, & 488 jD_cur,jH_cur_split,zdt_split, & 489 zplev,zplay, & 490 zpk,zphi,zphis, & 491 presnivs, & 492 zufi,zvfi,zrfi,ztfi,zqfi, & 493 flxwfi,pducov, & 494 zdufi,zdvfi,zdtfi,zdqfi,zdpsrf) 495 ! 496 ! else if ( planet_type=="generic" ) then 497 ! 498 ! CALL physiq (ngridmx, !! ngrid 499 ! . llm, !! nlayer 500 ! . nqtot, !! nq 501 ! . tracers(:)%name,!! tracer names from dynamical core (given in infotrac) 502 ! . debut_split, !! firstcall 503 ! . lafin_split, !! lastcall 504 ! . jD_cur, !! pday. see leapfrog 505 ! . jH_cur_split, !! ptime "fraction of day" 506 ! . zdt_split, !! ptimestep 507 ! . zplev, !! pplev 508 ! . zplay, !! pplay 509 ! . zphi, !! pphi 510 ! . zufi, !! pu 511 ! . zvfi, !! pv 512 ! . ztfi, !! pt 513 ! . zqfi, !! pq 514 ! . flxwfi, !! pw !! or 0. anyway this is for diagnostic. not used in physiq. 515 ! . zdufi, !! pdu 516 ! . zdvfi, !! pdv 517 ! . zdtfi, !! pdt 518 ! . zdqfi, !! pdq 519 ! . zdpsrf, !! pdpsrf 520 ! . tracerdyn) !! tracerdyn <-- utilite ??? 521 ! 522 ! endif ! of if (planet_type=="earth") 523 524 zufi(:,:)=zufi(:,:)+zdufi(:,:)*zdt_split 525 zvfi(:,:)=zvfi(:,:)+zdvfi(:,:)*zdt_split 526 ztfi(:,:)=ztfi(:,:)+zdtfi(:,:)*zdt_split 527 zqfi(:,:,:)=zqfi(:,:,:)+zdqfi(:,:,:)*zdt_split 528 529 zdufic(:,:)=zdufic(:,:)+zdufi(:,:) 530 zdvfic(:,:)=zdvfic(:,:)+zdvfi(:,:) 531 zdtfic(:,:)=zdtfic(:,:)+zdtfi(:,:) 532 zdqfic(:,:,:)=zdqfic(:,:,:)+zdqfi(:,:,:) 533 534 enddo ! of do isplit=1,nsplit_phys 535 536 #endif 537 ! of #ifdef CPP_PHYS 538 539 zdufi(:,:)=zdufic(:,:)/nsplit_phys 540 zdvfi(:,:)=zdvfic(:,:)/nsplit_phys 541 zdtfi(:,:)=zdtfic(:,:)/nsplit_phys 542 zdqfi(:,:,:)=zdqfic(:,:,:)/nsplit_phys 543 544 545 500 CONTINUE 546 547 !----------------------------------------------------------------------- 548 ! transformation des tendances physiques en tendances dynamiques: 549 ! --------------------------------------------------------------- 550 551 ! tendance sur la pression : 552 ! ----------------------------------- 553 554 CALL gr_fi_dyn(1,ngridmx,iip1,jjp1,zdpsrf,pdpsfi) 555 ! 556 ! 62. enthalpie potentielle 557 ! --------------------- 558 559 DO l=1,llm 560 561 DO i=1,iip1 562 pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l) 563 pdhfi(i,jjp1,l) = cpp * zdtfi(ngridmx,l)/ ppk(i,jjp1,l) 564 ENDDO 565 566 DO j=2,jjm 567 ig0=1+(j-2)*iim 568 DO i=1,iim 569 pdhfi(i,j,l) = cpp * zdtfi(ig0+i,l) / ppk(i,j,l) 570 ENDDO 571 pdhfi(iip1,j,l) = pdhfi(1,j,l) 572 ENDDO 573 574 ENDDO 575 576 577 ! 62. humidite specifique 578 ! --------------------- 579 ! Ehouarn: removed this useless bit: was overwritten at step 63 anyways 580 ! DO iq=1,nqtot 581 ! DO l=1,llm 582 ! DO i=1,iip1 583 ! pdqfi(i,1,l,iq) = zdqfi(1,l,iq) 584 ! pdqfi(i,jjp1,l,iq) = zdqfi(ngridmx,l,iq) 585 ! ENDDO 586 ! DO j=2,jjm 587 ! ig0=1+(j-2)*iim 588 ! DO i=1,iim 589 ! pdqfi(i,j,l,iq) = zdqfi(ig0+i,l,iq) 590 ! ENDDO 591 ! pdqfi(iip1,j,l,iq) = pdqfi(1,j,l,iq) 592 ! ENDDO 593 ! ENDDO 594 ! ENDDO 595 596 ! 63. traceurs 597 ! ------------ 598 ! initialisation des tendances 599 pdqfi(:,:,:,:)=0. 600 ! 601 itr = 0 602 DO iq=1,nqtot 603 IF(.NOT.tracers(iq)%isAdvected) CYCLE 604 itr = itr + 1 605 DO l=1,llm 606 DO i=1,iip1 607 pdqfi(i,1,l,iq) = zdqfi(1,l,itr) 608 pdqfi(i,jjp1,l,iq) = zdqfi(ngridmx,l,itr) 609 ENDDO 610 DO j=2,jjm 611 ig0=1+(j-2)*iim 612 DO i=1,iim 613 pdqfi(i,j,l,iq) = zdqfi(ig0+i,l,itr) 614 ENDDO 615 pdqfi(iip1,j,l,iq) = pdqfi(1,j,l,itr) 616 ENDDO 617 ENDDO 618 ENDDO 619 620 ! 65. champ u: 621 ! ------------ 622 623 DO l=1,llm 624 625 DO i=1,iip1 626 pdufi(i,1,l) = 0. 627 pdufi(i,jjp1,l) = 0. 628 ENDDO 629 630 DO j=2,jjm 631 ig0=1+(j-2)*iim 632 DO i=1,iim-1 633 pdufi(i,j,l)= & 634 0.5*(zdufi(ig0+i,l)+zdufi(ig0+i+1,l))*cu(i,j) 635 ENDDO 636 pdufi(iim,j,l)= & 637 0.5*(zdufi(ig0+1,l)+zdufi(ig0+iim,l))*cu(iim,j) 638 pdufi(iip1,j,l)=pdufi(1,j,l) 639 ENDDO 640 641 ENDDO 642 643 644 ! 67. champ v: 645 ! ------------ 646 647 DO l=1,llm 648 649 DO j=2,jjm-1 650 ig0=1+(j-2)*iim 651 DO i=1,iim 652 pdvfi(i,j,l)= & 653 0.5*(zdvfi(ig0+i,l)+zdvfi(ig0+i+iim,l))*cv(i,j) 654 ENDDO 655 pdvfi(iip1,j,l) = pdvfi(1,j,l) 656 ENDDO 657 ENDDO 658 659 660 ! 68. champ v pres des poles: 661 ! --------------------------- 662 ! v = U * cos(long) + V * SIN(long) 663 664 DO l=1,llm 665 666 DO i=1,iim 667 pdvfi(i,1,l)= & 668 zdufi(1,l)*COS(rlonv(i))+zdvfi(1,l)*SIN(rlonv(i)) 669 pdvfi(i,jjm,l)=zdufi(ngridmx,l)*COS(rlonv(i)) & 670 +zdvfi(ngridmx,l)*SIN(rlonv(i)) 671 pdvfi(i,1,l)= & 672 0.5*(pdvfi(i,1,l)+zdvfi(i+1,l))*cv(i,1) 673 pdvfi(i,jjm,l)= & 674 0.5*(pdvfi(i,jjm,l)+zdvfi(ngridmx-iip1+i,l))*cv(i,jjm) 380 675 ENDDO 381 676 382 c 46.champ v: 383 c ----------- 384 385 DO l=1,llm 386 DO j=2,jjm 387 ig0=1+(j-2)*iim 388 DO i=1,iim 389 zvfi(ig0+i,l)= 0.5 * 390 $ ( pvcov(i,j-1,l)/cv(i,j-1) + pvcov(i,j,l)/cv(i,j) ) 391 pcvgv(ig0+i,l)= 0.5 * 392 $ ( pdvcov(i,j-1,l)/cv(i,j-1) + pdvcov(i,j,l)/cv(i,j) ) 393 ENDDO 394 zrfi(ig0 + 1,l)= 0.25 *(zrot(iim,j-1,l)+zrot(iim,j,l) 395 & +zrot(1,j-1,l)+zrot(1,j,l)) 396 DO i=2,iim 397 zrfi(ig0 + i,l)= 0.25 *(zrot(i-1,j-1,l)+zrot(i-1,j,l) 398 $ +zrot(i,j-1,l)+zrot(i,j,l)) ! AdlC MAY 2014 399 ENDDO 400 ENDDO 401 ENDDO 402 403 404 c 47. champs de vents aux pole nord 405 c ------------------------------ 406 c U = 1 / pi * integrale [ v * cos(long) * d long ] 407 c V = 1 / pi * integrale [ v * sin(long) * d long ] 408 409 DO l=1,llm 410 411 z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,1,l)/cv(1,1) 412 z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,1,l)/cv(1,1) 413 DO i=2,iim 414 z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1) 415 z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,1,l)/cv(i,1) 416 ENDDO 417 418 DO i=1,iim 419 zcos(i) = COS(rlonv(i))*z1(i) 420 zcosbis(i)= COS(rlonv(i))*z1bis(i) 421 zsin(i) = SIN(rlonv(i))*z1(i) 422 zsinbis(i)= SIN(rlonv(i))*z1bis(i) 423 ENDDO 424 425 zufi(1,l) = SSUM(iim,zcos,1)/pi 426 pcvgu(1,l) = SSUM(iim,zcosbis,1)/pi 427 zvfi(1,l) = SSUM(iim,zsin,1)/pi 428 pcvgv(1,l) = SSUM(iim,zsinbis,1)/pi 429 zrfi(1, l) = 0. 430 ENDDO 431 432 433 c 48. champs de vents aux pole sud: 434 c --------------------------------- 435 c U = 1 / pi * integrale [ v * cos(long) * d long ] 436 c V = 1 / pi * integrale [ v * sin(long) * d long ] 437 438 DO l=1,llm 439 440 z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,jjm,l)/cv(1,jjm) 441 z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,jjm,l)/cv(1,jjm) 442 DO i=2,iim 443 z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm) 444 z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,jjm,l)/cv(i,jjm) 445 ENDDO 446 447 DO i=1,iim 448 zcos(i) = COS(rlonv(i))*z1(i) 449 zcosbis(i) = COS(rlonv(i))*z1bis(i) 450 zsin(i) = SIN(rlonv(i))*z1(i) 451 zsinbis(i) = SIN(rlonv(i))*z1bis(i) 452 ENDDO 453 454 zufi(ngridmx,l) = SSUM(iim,zcos,1)/pi 455 pcvgu(ngridmx,l) = SSUM(iim,zcosbis,1)/pi 456 zvfi(ngridmx,l) = SSUM(iim,zsin,1)/pi 457 pcvgv(ngridmx,l) = SSUM(iim,zsinbis,1)/pi 458 zrfi(ngridmx, l) = 0. 459 ENDDO 460 c 461 c On change de grille, dynamique vers physiq, pour le flux de masse verticale 462 CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,flxw,flxwfi) 463 464 c----------------------------------------------------------------------- 465 c Appel de la physique: 466 c --------------------- 467 468 469 470 ! write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys 471 zdt_split=dtphys/nsplit_phys 472 zdufic(:,:)=0. 473 zdvfic(:,:)=0. 474 zdtfic(:,:)=0. 475 zdqfic(:,:,:)=0. 476 477 #ifdef CPP_PHYS 478 479 do isplit=1,nsplit_phys 480 481 jH_cur_split=jH_cur+(isplit-1) * dtvr / (daysec *nsplit_phys) 482 debut_split=debut.and.isplit==1 483 lafin_split=lafin.and.isplit==nsplit_phys 484 485 ! if (planet_type=="earth") then 486 CALL call_physiq(ngridmx,llm,nqtot,tracers(:)%name, 487 & debut_split,lafin_split, 488 & jD_cur,jH_cur_split,zdt_split, 489 & zplev,zplay, 490 & zpk,zphi,zphis, 491 & presnivs, 492 & zufi,zvfi,zrfi,ztfi,zqfi, 493 & flxwfi,pducov, 494 & zdufi,zdvfi,zdtfi,zdqfi,zdpsrf) 495 ! 496 ! else if ( planet_type=="generic" ) then 497 ! 498 ! CALL physiq (ngridmx, !! ngrid 499 ! . llm, !! nlayer 500 ! . nqtot, !! nq 501 ! . tracers(:)%name,!! tracer names from dynamical core (given in infotrac) 502 ! . debut_split, !! firstcall 503 ! . lafin_split, !! lastcall 504 ! . jD_cur, !! pday. see leapfrog 505 ! . jH_cur_split, !! ptime "fraction of day" 506 ! . zdt_split, !! ptimestep 507 ! . zplev, !! pplev 508 ! . zplay, !! pplay 509 ! . zphi, !! pphi 510 ! . zufi, !! pu 511 ! . zvfi, !! pv 512 ! . ztfi, !! pt 513 ! . zqfi, !! pq 514 ! . flxwfi, !! pw !! or 0. anyway this is for diagnostic. not used in physiq. 515 ! . zdufi, !! pdu 516 ! . zdvfi, !! pdv 517 ! . zdtfi, !! pdt 518 ! . zdqfi, !! pdq 519 ! . zdpsrf, !! pdpsrf 520 ! . tracerdyn) !! tracerdyn <-- utilite ??? 521 ! 522 ! endif ! of if (planet_type=="earth") 523 524 zufi(:,:)=zufi(:,:)+zdufi(:,:)*zdt_split 525 zvfi(:,:)=zvfi(:,:)+zdvfi(:,:)*zdt_split 526 ztfi(:,:)=ztfi(:,:)+zdtfi(:,:)*zdt_split 527 zqfi(:,:,:)=zqfi(:,:,:)+zdqfi(:,:,:)*zdt_split 528 529 zdufic(:,:)=zdufic(:,:)+zdufi(:,:) 530 zdvfic(:,:)=zdvfic(:,:)+zdvfi(:,:) 531 zdtfic(:,:)=zdtfic(:,:)+zdtfi(:,:) 532 zdqfic(:,:,:)=zdqfic(:,:,:)+zdqfi(:,:,:) 533 534 enddo ! of do isplit=1,nsplit_phys 535 536 #endif 537 ! of #ifdef CPP_PHYS 538 539 zdufi(:,:)=zdufic(:,:)/nsplit_phys 540 zdvfi(:,:)=zdvfic(:,:)/nsplit_phys 541 zdtfi(:,:)=zdtfic(:,:)/nsplit_phys 542 zdqfi(:,:,:)=zdqfic(:,:,:)/nsplit_phys 543 544 545 500 CONTINUE 546 547 c----------------------------------------------------------------------- 548 c transformation des tendances physiques en tendances dynamiques: 549 c --------------------------------------------------------------- 550 551 c tendance sur la pression : 552 c ----------------------------------- 553 554 CALL gr_fi_dyn(1,ngridmx,iip1,jjp1,zdpsrf,pdpsfi) 555 c 556 c 62. enthalpie potentielle 557 c --------------------- 558 559 DO l=1,llm 560 561 DO i=1,iip1 562 pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l) 563 pdhfi(i,jjp1,l) = cpp * zdtfi(ngridmx,l)/ ppk(i,jjp1,l) 564 ENDDO 565 566 DO j=2,jjm 567 ig0=1+(j-2)*iim 568 DO i=1,iim 569 pdhfi(i,j,l) = cpp * zdtfi(ig0+i,l) / ppk(i,j,l) 570 ENDDO 571 pdhfi(iip1,j,l) = pdhfi(1,j,l) 572 ENDDO 573 574 ENDDO 575 576 577 c 62. humidite specifique 578 c --------------------- 579 ! Ehouarn: removed this useless bit: was overwritten at step 63 anyways 580 ! DO iq=1,nqtot 581 ! DO l=1,llm 582 ! DO i=1,iip1 583 ! pdqfi(i,1,l,iq) = zdqfi(1,l,iq) 584 ! pdqfi(i,jjp1,l,iq) = zdqfi(ngridmx,l,iq) 585 ! ENDDO 586 ! DO j=2,jjm 587 ! ig0=1+(j-2)*iim 588 ! DO i=1,iim 589 ! pdqfi(i,j,l,iq) = zdqfi(ig0+i,l,iq) 590 ! ENDDO 591 ! pdqfi(iip1,j,l,iq) = pdqfi(1,j,l,iq) 592 ! ENDDO 593 ! ENDDO 594 ! ENDDO 595 596 c 63. traceurs 597 c ------------ 598 C initialisation des tendances 599 pdqfi(:,:,:,:)=0. 600 C 601 itr = 0 602 DO iq=1,nqtot 603 IF(.NOT.tracers(iq)%isAdvected) CYCLE 604 itr = itr + 1 605 DO l=1,llm 606 DO i=1,iip1 607 pdqfi(i,1,l,iq) = zdqfi(1,l,itr) 608 pdqfi(i,jjp1,l,iq) = zdqfi(ngridmx,l,itr) 609 ENDDO 610 DO j=2,jjm 611 ig0=1+(j-2)*iim 612 DO i=1,iim 613 pdqfi(i,j,l,iq) = zdqfi(ig0+i,l,itr) 614 ENDDO 615 pdqfi(iip1,j,l,iq) = pdqfi(1,j,l,itr) 616 ENDDO 617 ENDDO 618 ENDDO 619 620 c 65. champ u: 621 c ------------ 622 623 DO l=1,llm 624 625 DO i=1,iip1 626 pdufi(i,1,l) = 0. 627 pdufi(i,jjp1,l) = 0. 628 ENDDO 629 630 DO j=2,jjm 631 ig0=1+(j-2)*iim 632 DO i=1,iim-1 633 pdufi(i,j,l)= 634 $ 0.5*(zdufi(ig0+i,l)+zdufi(ig0+i+1,l))*cu(i,j) 635 ENDDO 636 pdufi(iim,j,l)= 637 $ 0.5*(zdufi(ig0+1,l)+zdufi(ig0+iim,l))*cu(iim,j) 638 pdufi(iip1,j,l)=pdufi(1,j,l) 639 ENDDO 640 641 ENDDO 642 643 644 c 67. champ v: 645 c ------------ 646 647 DO l=1,llm 648 649 DO j=2,jjm-1 650 ig0=1+(j-2)*iim 651 DO i=1,iim 652 pdvfi(i,j,l)= 653 $ 0.5*(zdvfi(ig0+i,l)+zdvfi(ig0+i+iim,l))*cv(i,j) 654 ENDDO 655 pdvfi(iip1,j,l) = pdvfi(1,j,l) 656 ENDDO 657 ENDDO 658 659 660 c 68. champ v pres des poles: 661 c --------------------------- 662 c v = U * cos(long) + V * SIN(long) 663 664 DO l=1,llm 665 666 DO i=1,iim 667 pdvfi(i,1,l)= 668 $ zdufi(1,l)*COS(rlonv(i))+zdvfi(1,l)*SIN(rlonv(i)) 669 pdvfi(i,jjm,l)=zdufi(ngridmx,l)*COS(rlonv(i)) 670 $ +zdvfi(ngridmx,l)*SIN(rlonv(i)) 671 pdvfi(i,1,l)= 672 $ 0.5*(pdvfi(i,1,l)+zdvfi(i+1,l))*cv(i,1) 673 pdvfi(i,jjm,l)= 674 $ 0.5*(pdvfi(i,jjm,l)+zdvfi(ngridmx-iip1+i,l))*cv(i,jjm) 675 ENDDO 676 677 pdvfi(iip1,1,l) = pdvfi(1,1,l) 678 pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l) 679 680 ENDDO 681 682 c----------------------------------------------------------------------- 677 pdvfi(iip1,1,l) = pdvfi(1,1,l) 678 pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l) 679 680 ENDDO 681 682 !----------------------------------------------------------------------- 683 683 684 684 700 CONTINUE 685 686 687 688 689 END 685 686 firstcal = .FALSE. 687 688 RETURN 689 END SUBROUTINE calfis -
LMDZ6/trunk/libf/dynphy_lonlat/calfis_loc.F90
r5245 r5246 2 2 ! $Id$ 3 3 ! 4 C 5 C 6 SUBROUTINE calfis_loc(lafin, 7 $ jD_cur, jH_cur,8 $ pucov,9 $ pvcov,10 $ pteta,11 $ pq,12 $ pmasse,13 $ pps,14 $ pp,15 $ ppk,16 $ pphis,17 $ pphi,18 $ pducov,19 $ pdvcov,20 $ pdteta,21 $ pdq,22 $ flxw,23 $ pdufi,24 $ pdvfi,25 $ pdhfi,26 $ pdqfi,27 $pdpsfi)4 ! 5 ! 6 SUBROUTINE calfis_loc(lafin, & 7 jD_cur, jH_cur, & 8 pucov, & 9 pvcov, & 10 pteta, & 11 pq, & 12 pmasse, & 13 pps, & 14 pp, & 15 ppk, & 16 pphis, & 17 pphi, & 18 pducov, & 19 pdvcov, & 20 pdteta, & 21 pdq, & 22 flxw, & 23 pdufi, & 24 pdvfi, & 25 pdhfi, & 26 pdqfi, & 27 pdpsfi) 28 28 #ifdef CPP_PHYS 29 ! If using physics30 c 31 c Auteur : P. Le Van, F. Hourdin 32 c.........33 34 35 36 37 38 29 ! If using physics 30 ! 31 ! Auteur : P. Le Van, F. Hourdin 32 ! ......... 33 USE dimphy 34 USE mod_phys_lmdz_mpi_data, mpi_root_xx=>mpi_master 35 USE mod_phys_lmdz_omp_data, ONLY: klon_omp, klon_omp_begin 36 USE mod_const_mpi, ONLY: COMM_LMDZ 37 USE mod_interface_dyn_phys 38 USE IOPHY 39 39 #endif 40 40 USE lmdz_mpi 41 41 42 42 #ifdef CPP_PARA 43 USE parallel_lmdz,ONLY:omp_chunk,using_mpi,jjb_u,jje_u,jjb_v,jje_v44 $,jj_begin_dyn=>jj_begin,jj_end_dyn=>jj_end45 46 47 43 USE parallel_lmdz,ONLY:omp_chunk,using_mpi,jjb_u,jje_u,jjb_v,jje_v & 44 ,jj_begin_dyn=>jj_begin,jj_end_dyn=>jj_end 45 USE Write_Field 46 Use Write_field_p 47 USE Times 48 48 #endif 49 50 49 USE infotrac, ONLY: nqtot, tracers 50 USE control_mod, ONLY: planet_type, nsplit_phys 51 51 #ifdef CPP_PHYS 52 53 #endif 54 55 52 USE callphysiq_mod, ONLY: call_physiq 53 #endif 54 USE comvert_mod, ONLY: preff, presnivs 55 USE comconst_mod, ONLY: cpp, daysec, dtphys, dtvr, kappa, pi 56 56 57 57 #ifdef CPP_PARA 58 59 c=======================================================================60 c 61 c1. rearrangement des tableaux et transformation62 cvariables dynamiques > variables physiques63 c2. calcul des termes physiques64 c3. retransformation des tendances physiques en tendances dynamiques65 c 66 cremarques:67 c----------68 c 69 c - les vents sont donnes dans la physique par leurs composantes 70 cnaturelles.71 c- la variable thermodynamique de la physique est une variable72 c intensive : T 73 cpour la dynamique on prend T * ( preff / p(l) ) **kappa74 c- les deux seules variables dependant de la geometrie necessaires75 c pour la physique sont la latitude pour le rayonnement et 76 c l'aire de la maille quand on veut integrer une grandeur 77 chorizontalement.78 c - les points de la physique sont les points scalaires de la 79 cla dynamique; numerotation:80 c1 pour le pole nord81 c(jjm-1)*iim pour l'interieur du domaine82 cngridmx pour le pole sud83 c---> ngridmx=2+(jjm-1)*iim84 c 85 cInput :86 c-------87 cecritphy frequence d'ecriture (en jours)de histphy88 cpucov covariant zonal velocity89 c pvcov covariant meridional velocity 90 cpteta potential temperature91 cpps surface pressure92 cpmasse masse d'air dans chaque maille93 cpts surface temperature (K)94 ccallrad clef d'appel au rayonnement95 c 96 cOutput :97 c--------98 cpdufi tendency for the natural zonal velocity (ms-1)99 c pdvfi tendency for the natural meridional velocity 100 cpdhfi tendency for the potential temperature101 cpdtsfi tendency for the surface temperature102 c 103 cpdtrad radiative tendencies \ both input104 cpfluxrad radiative fluxes / and output105 c 106 c=======================================================================107 c 108 c-----------------------------------------------------------------------109 c 110 c0. Declarations :111 c------------------112 113 114 115 116 INTEGERngridmx117 118 119 120 121 cArguments :122 c-----------123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 58 IMPLICIT NONE 59 !======================================================================= 60 ! 61 ! 1. rearrangement des tableaux et transformation 62 ! variables dynamiques > variables physiques 63 ! 2. calcul des termes physiques 64 ! 3. retransformation des tendances physiques en tendances dynamiques 65 ! 66 ! remarques: 67 ! ---------- 68 ! 69 ! - les vents sont donnes dans la physique par leurs composantes 70 ! naturelles. 71 ! - la variable thermodynamique de la physique est une variable 72 ! intensive : T 73 ! pour la dynamique on prend T * ( preff / p(l) ) **kappa 74 ! - les deux seules variables dependant de la geometrie necessaires 75 ! pour la physique sont la latitude pour le rayonnement et 76 ! l'aire de la maille quand on veut integrer une grandeur 77 ! horizontalement. 78 ! - les points de la physique sont les points scalaires de la 79 ! la dynamique; numerotation: 80 ! 1 pour le pole nord 81 ! (jjm-1)*iim pour l'interieur du domaine 82 ! ngridmx pour le pole sud 83 ! ---> ngridmx=2+(jjm-1)*iim 84 ! 85 ! Input : 86 ! ------- 87 ! ecritphy frequence d'ecriture (en jours)de histphy 88 ! pucov covariant zonal velocity 89 ! pvcov covariant meridional velocity 90 ! pteta potential temperature 91 ! pps surface pressure 92 ! pmasse masse d'air dans chaque maille 93 ! pts surface temperature (K) 94 ! callrad clef d'appel au rayonnement 95 ! 96 ! Output : 97 ! -------- 98 ! pdufi tendency for the natural zonal velocity (ms-1) 99 ! pdvfi tendency for the natural meridional velocity 100 ! pdhfi tendency for the potential temperature 101 ! pdtsfi tendency for the surface temperature 102 ! 103 ! pdtrad radiative tendencies \ both input 104 ! pfluxrad radiative fluxes / and output 105 ! 106 !======================================================================= 107 ! 108 !----------------------------------------------------------------------- 109 ! 110 ! 0. Declarations : 111 ! ------------------ 112 113 include "dimensions.h" 114 include "paramet.h" 115 116 INTEGER :: ngridmx 117 PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm ) 118 119 include "comgeom2.h" 120 include "iniprint.h" 121 ! Arguments : 122 ! ----------- 123 LOGICAL,INTENT(IN) :: lafin ! .true. for the very last call to physics 124 REAL,INTENT(IN):: jD_cur, jH_cur 125 REAL,INTENT(IN):: pvcov(iip1,jjb_v:jje_v,llm) ! covariant meridional velocity 126 REAL,INTENT(IN):: pucov(iip1,jjb_u:jje_u,llm) ! covariant zonal velocity 127 REAL,INTENT(IN):: pteta(iip1,jjb_u:jje_u,llm) ! potential temperature 128 REAL,INTENT(IN):: pmasse(iip1,jjb_u:jje_u,llm) ! mass in each cell ! not used 129 REAL,INTENT(IN):: pq(iip1,jjb_u:jje_u,llm,nqtot) ! tracers 130 REAL,INTENT(IN):: pphis(iip1,jjb_u:jje_u) ! surface geopotential 131 REAL,INTENT(IN):: pphi(iip1,jjb_u:jje_u,llm) ! geopotential 132 133 REAL,INTENT(IN) :: pdvcov(iip1,jjb_v:jje_v,llm) ! dynamical tendency on vcov ! not used 134 REAL,INTENT(IN) :: pducov(iip1,jjb_u:jje_u,llm) ! dynamical tendency on ucov 135 REAL,INTENT(IN) :: pdteta(iip1,jjb_u:jje_u,llm) ! dynamical tendency on teta ! not used 136 REAL,INTENT(IN) :: pdq(iip1,jjb_u:jje_u,llm,nqtot) ! dynamical tendency on tracers ! not used 137 138 REAL,INTENT(IN) :: pps(iip1,jjb_u:jje_u) ! surface pressure (Pa) 139 REAL,INTENT(IN) :: pp(iip1,jjb_u:jje_u,llmp1) ! pressure at mesh interfaces (Pa) 140 REAL,INTENT(IN) :: ppk(iip1,jjb_u:jje_u,llm) ! Exner at mid-layer 141 REAL,INTENT(IN) :: flxw(iip1,jjb_u:jje_u,llm) ! Vertical mass flux on lower mesh interfaces (kg/s) (on llm because flxw(:,:,llm+1)=0) 142 143 ! ! tendencies (in */s) from the physics 144 REAL,INTENT(OUT) :: pdvfi(iip1,jjb_v:jje_v,llm) ! tendency on covariant meridional wind 145 REAL,INTENT(OUT) :: pdufi(iip1,jjb_u:jje_u,llm) ! tendency on covariant zonal wind 146 REAL,INTENT(OUT) :: pdhfi(iip1,jjb_u:jje_u,llm) ! tendency on potential temperature (K/s) 147 REAL,INTENT(OUT) :: pdqfi(iip1,jjb_u:jje_u,llm,nqtot) ! tendency on tracers 148 REAL,INTENT(OUT) :: pdpsfi(iip1,jjb_u:jje_u) ! tendency on surface pressure (Pa/s) 149 149 150 150 #ifdef CPP_PHYS 151 ! Ehouarn: for now calfis_p needs some informations from physics to compile 152 c Local variables : 153 c ----------------- 154 155 INTEGER i,j,l,ig0,ig,iq,itr 156 REAL,ALLOCATABLE,SAVE :: zpsrf(:) 157 REAL,ALLOCATABLE,SAVE :: zplev(:,:),zplay(:,:) 158 REAL,ALLOCATABLE,SAVE :: zphi(:,:),zphis(:) 159 c 160 REAL zrot(iip1,jjb_v:jje_v,llm) ! AdlC May 2014 161 REAL,ALLOCATABLE,SAVE :: zufi(:,:), zvfi(:,:), zrfi(:,:) 162 REAL,ALLOCATABLE,SAVE :: ztfi(:,:),zqfi(:,:,:) 163 REAL,ALLOCATABLE,SAVE :: zpk(:,:) 164 c 165 REAL,ALLOCATABLE,SAVE :: pcvgu(:,:), pcvgv(:,:) 166 REAL,ALLOCATABLE,SAVE :: pcvgt(:,:), pcvgq(:,:,:) 167 c 168 REAL,ALLOCATABLE,SAVE :: zdufi(:,:),zdvfi(:,:) 169 REAL,ALLOCATABLE,SAVE :: zdtfi(:,:),zdqfi(:,:,:) 170 REAL,ALLOCATABLE,SAVE :: zdpsrf(:) 171 REAL,SAVE,ALLOCATABLE :: flxwfi(:,:) ! Flux de masse verticale sur la grille physiq 172 173 c 174 REAL,ALLOCATABLE,SAVE :: zplev_omp(:,:) 175 REAL,ALLOCATABLE,SAVE :: zplay_omp(:,:) 176 REAL,ALLOCATABLE,SAVE :: zpk_omp(:,:) 177 REAL,ALLOCATABLE,SAVE :: zphi_omp(:,:) 178 REAL,ALLOCATABLE,SAVE :: zphis_omp(:) 179 REAL,ALLOCATABLE,SAVE :: presnivs_omp(:) 180 REAL,ALLOCATABLE,SAVE :: zufi_omp(:,:) 181 REAL,ALLOCATABLE,SAVE :: zvfi_omp(:,:) 182 REAL,ALLOCATABLE,SAVE :: zrfi_omp(:,:) 183 REAL,ALLOCATABLE,SAVE :: ztfi_omp(:,:) 184 REAL,ALLOCATABLE,SAVE :: zqfi_omp(:,:,:) 185 REAL,ALLOCATABLE,SAVE :: zdufi_omp(:,:) 186 REAL,ALLOCATABLE,SAVE :: zdvfi_omp(:,:) 187 REAL,ALLOCATABLE,SAVE :: zdtfi_omp(:,:) 188 REAL,ALLOCATABLE,SAVE :: zdqfi_omp(:,:,:) 189 REAL,ALLOCATABLE,SAVE :: zdpsrf_omp(:) 190 REAL,SAVE,ALLOCATABLE :: flxwfi_omp(:,:) ! Flux de masse verticale sur la grille physiq 191 192 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 193 ! Introduction du splitting (FH) 194 ! Question pour Yann : 195 ! J'ai été surpris au début que les tableaux zufi_omp, zdufi_omp n'co soitent 196 ! en SAVE. Je crois comprendre que c'est parce que tu voulais qu'il 197 ! soit allocatable (plutot par exemple que de passer une dimension 198 ! dépendant du process en argument des routines) et que, du coup, 199 ! le SAVE évite d'avoir à refaire l'allocation à chaque appel. 200 ! Tu confirmes ? 201 ! J'ai suivi le même principe pour les zdufic_omp 202 ! Mais c'est surement bien que tu controles. 203 ! 204 205 REAL,ALLOCATABLE,SAVE :: zdufic_omp(:,:) 206 REAL,ALLOCATABLE,SAVE :: zdvfic_omp(:,:) 207 REAL,ALLOCATABLE,SAVE :: zdtfic_omp(:,:) 208 REAL,ALLOCATABLE,SAVE :: zdqfic_omp(:,:,:) 209 REAL jH_cur_split,zdt_split 210 LOGICAL debut_split,lafin_split 211 INTEGER isplit 212 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 213 214 c$OMP THREADPRIVATE(zplev_omp,zplay_omp,zpk_omp,zphi_omp,zphis_omp, 215 c$OMP+ presnivs_omp,zufi_omp,zvfi_omp,ztfi_omp, 216 c$OMP+ zrfi_omp,zqfi_omp,zdufi_omp,zdvfi_omp, 217 c$OMP+ zdtfi_omp,zdqfi_omp,zdpsrf_omp,flxwfi_omp, 218 c$OMP+ zdufic_omp,zdvfic_omp,zdtfic_omp,zdqfic_omp) 219 220 LOGICAL,SAVE :: first_omp=.true. 221 c$OMP THREADPRIVATE(first_omp) 222 223 REAL zsin(iim),zcos(iim),z1(iim) 224 REAL zsinbis(iim),zcosbis(iim),z1bis(iim) 225 REAL unskap, pksurcp 226 c 227 REAL SSUM 228 229 LOGICAL,SAVE :: firstcal=.true., debut=.true. 230 c$OMP THREADPRIVATE(firstcal,debut) 231 232 REAL,SAVE,dimension(1:iim,1:llm):: du_send,du_recv,dv_send,dv_recv 233 INTEGER :: ierr 234 INTEGER,dimension(MPI_STATUS_SIZE,4) :: Status 235 INTEGER, dimension(4) :: Req 236 REAL,ALLOCATABLE,SAVE:: zdufi2(:,:),zdvfi2(:,:) 237 integer :: k,kstart,kend 238 INTEGER :: offset 239 INTEGER :: jjb,jje 240 241 c 242 c----------------------------------------------------------------------- 243 c 244 c 1. Initialisations : 245 c -------------------- 246 c 247 248 klon=klon_mpi 249 250 c 251 IF ( firstcal ) THEN 252 debut = .TRUE. 253 IF (ngridmx.NE.2+(jjm-1)*iim) THEN 254 write(lunout,*) 'STOP dans calfis' 255 write(lunout,*) 256 & 'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim' 257 write(lunout,*) ' ngridmx jjm iim ' 258 write(lunout,*) ngridmx,jjm,iim 259 call abort_gcm("calfis_loc", "", 1) 260 ENDIF 261 c$OMP MASTER 262 ALLOCATE(zpsrf(klon)) 263 ALLOCATE(zplev(klon,llm+1),zplay(klon,llm)) 264 ALLOCATE(zphi(klon,llm),zphis(klon)) 265 ALLOCATE(zufi(klon,llm), zvfi(klon,llm),zrfi(klon,llm)) 266 ALLOCATE(ztfi(klon,llm),zqfi(klon,llm,nqtot)) 267 ALLOCATE(pcvgu(klon,llm), pcvgv(klon,llm)) 268 ALLOCATE(pcvgt(klon,llm), pcvgq(klon,llm,2)) 269 ALLOCATE(zdufi(klon,llm),zdvfi(klon,llm)) 270 ALLOCATE(zdtfi(klon,llm),zdqfi(klon,llm,nqtot)) 271 ALLOCATE(zdpsrf(klon)) 272 ALLOCATE(zdufi2(klon+iim,llm),zdvfi2(klon+iim,llm)) 273 ALLOCATE(flxwfi(klon,llm)) 274 ALLOCATE(zpk(klon,llm)) 275 c$OMP END MASTER 276 c$OMP BARRIER 277 ELSE 278 debut = .FALSE. 279 ENDIF 280 281 c 282 c 283 c----------------------------------------------------------------------- 284 c 40. transformation des variables dynamiques en variables physiques: 285 c --------------------------------------------------------------- 286 287 c 41. pressions au sol (en Pascals) 288 c ---------------------------------- 289 290 c$OMP MASTER 291 call start_timer(timer_physic) 292 c$OMP END MASTER 293 294 c$OMP MASTER 295 !CDIR ON_ADB(index_i) 296 !CDIR ON_ADB(index_j) 297 do ig0=1,klon 298 i=index_i(ig0) 299 j=index_j(ig0) 300 zpsrf(ig0)=pps(i,j) 151 ! Ehouarn: for now calfis_p needs some informations from physics to compile 152 ! Local variables : 153 ! ----------------- 154 155 INTEGER :: i,j,l,ig0,ig,iq,itr 156 REAL,ALLOCATABLE,SAVE :: zpsrf(:) 157 REAL,ALLOCATABLE,SAVE :: zplev(:,:),zplay(:,:) 158 REAL,ALLOCATABLE,SAVE :: zphi(:,:),zphis(:) 159 ! 160 REAL :: zrot(iip1,jjb_v:jje_v,llm) ! AdlC May 2014 161 REAL,ALLOCATABLE,SAVE :: zufi(:,:), zvfi(:,:), zrfi(:,:) 162 REAL,ALLOCATABLE,SAVE :: ztfi(:,:),zqfi(:,:,:) 163 REAL,ALLOCATABLE,SAVE :: zpk(:,:) 164 ! 165 REAL,ALLOCATABLE,SAVE :: pcvgu(:,:), pcvgv(:,:) 166 REAL,ALLOCATABLE,SAVE :: pcvgt(:,:), pcvgq(:,:,:) 167 ! 168 REAL,ALLOCATABLE,SAVE :: zdufi(:,:),zdvfi(:,:) 169 REAL,ALLOCATABLE,SAVE :: zdtfi(:,:),zdqfi(:,:,:) 170 REAL,ALLOCATABLE,SAVE :: zdpsrf(:) 171 REAL,SAVE,ALLOCATABLE :: flxwfi(:,:) ! Flux de masse verticale sur la grille physiq 172 173 ! 174 REAL,ALLOCATABLE,SAVE :: zplev_omp(:,:) 175 REAL,ALLOCATABLE,SAVE :: zplay_omp(:,:) 176 REAL,ALLOCATABLE,SAVE :: zpk_omp(:,:) 177 REAL,ALLOCATABLE,SAVE :: zphi_omp(:,:) 178 REAL,ALLOCATABLE,SAVE :: zphis_omp(:) 179 REAL,ALLOCATABLE,SAVE :: presnivs_omp(:) 180 REAL,ALLOCATABLE,SAVE :: zufi_omp(:,:) 181 REAL,ALLOCATABLE,SAVE :: zvfi_omp(:,:) 182 REAL,ALLOCATABLE,SAVE :: zrfi_omp(:,:) 183 REAL,ALLOCATABLE,SAVE :: ztfi_omp(:,:) 184 REAL,ALLOCATABLE,SAVE :: zqfi_omp(:,:,:) 185 REAL,ALLOCATABLE,SAVE :: zdufi_omp(:,:) 186 REAL,ALLOCATABLE,SAVE :: zdvfi_omp(:,:) 187 REAL,ALLOCATABLE,SAVE :: zdtfi_omp(:,:) 188 REAL,ALLOCATABLE,SAVE :: zdqfi_omp(:,:,:) 189 REAL,ALLOCATABLE,SAVE :: zdpsrf_omp(:) 190 REAL,SAVE,ALLOCATABLE :: flxwfi_omp(:,:) ! Flux de masse verticale sur la grille physiq 191 192 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 193 ! Introduction du splitting (FH) 194 ! Question pour Yann : 195 ! J'ai été surpris au début que les tableaux zufi_omp, zdufi_omp n'co soitent 196 ! en SAVE. Je crois comprendre que c'est parce que tu voulais qu'il 197 ! soit allocatable (plutot par exemple que de passer une dimension 198 ! dépendant du process en argument des routines) et que, du coup, 199 ! le SAVE évite d'avoir à refaire l'allocation à chaque appel. 200 ! Tu confirmes ? 201 ! J'ai suivi le même principe pour les zdufic_omp 202 ! Mais c'est surement bien que tu controles. 203 ! 204 205 REAL,ALLOCATABLE,SAVE :: zdufic_omp(:,:) 206 REAL,ALLOCATABLE,SAVE :: zdvfic_omp(:,:) 207 REAL,ALLOCATABLE,SAVE :: zdtfic_omp(:,:) 208 REAL,ALLOCATABLE,SAVE :: zdqfic_omp(:,:,:) 209 REAL :: jH_cur_split,zdt_split 210 LOGICAL :: debut_split,lafin_split 211 INTEGER :: isplit 212 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 213 214 !$OMP THREADPRIVATE(zplev_omp,zplay_omp,zpk_omp,zphi_omp,zphis_omp, & 215 !$OMP presnivs_omp,zufi_omp,zvfi_omp,ztfi_omp, & 216 !$OMP zrfi_omp,zqfi_omp,zdufi_omp,zdvfi_omp, & 217 !$OMP zdtfi_omp,zdqfi_omp,zdpsrf_omp,flxwfi_omp, & 218 !$OMP zdufic_omp,zdvfic_omp,zdtfic_omp,zdqfic_omp) 219 220 LOGICAL,SAVE :: first_omp=.true. 221 !$OMP THREADPRIVATE(first_omp) 222 223 REAL :: zsin(iim),zcos(iim),z1(iim) 224 REAL :: zsinbis(iim),zcosbis(iim),z1bis(iim) 225 REAL :: unskap, pksurcp 226 ! 227 REAL :: SSUM 228 229 LOGICAL,SAVE :: firstcal=.true., debut=.true. 230 !$OMP THREADPRIVATE(firstcal,debut) 231 232 REAL,SAVE,dimension(1:iim,1:llm):: du_send,du_recv,dv_send,dv_recv 233 INTEGER :: ierr 234 INTEGER,dimension(MPI_STATUS_SIZE,4) :: Status 235 INTEGER, dimension(4) :: Req 236 REAL,ALLOCATABLE,SAVE:: zdufi2(:,:),zdvfi2(:,:) 237 integer :: k,kstart,kend 238 INTEGER :: offset 239 INTEGER :: jjb,jje 240 241 ! 242 !----------------------------------------------------------------------- 243 ! 244 ! 1. Initialisations : 245 ! -------------------- 246 ! 247 248 klon=klon_mpi 249 250 ! 251 IF ( firstcal ) THEN 252 debut = .TRUE. 253 IF (ngridmx.NE.2+(jjm-1)*iim) THEN 254 write(lunout,*) 'STOP dans calfis' 255 write(lunout,*) & 256 'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim' 257 write(lunout,*) ' ngridmx jjm iim ' 258 write(lunout,*) ngridmx,jjm,iim 259 call abort_gcm("calfis_loc", "", 1) 260 ENDIF 261 !$OMP MASTER 262 ALLOCATE(zpsrf(klon)) 263 ALLOCATE(zplev(klon,llm+1),zplay(klon,llm)) 264 ALLOCATE(zphi(klon,llm),zphis(klon)) 265 ALLOCATE(zufi(klon,llm), zvfi(klon,llm),zrfi(klon,llm)) 266 ALLOCATE(ztfi(klon,llm),zqfi(klon,llm,nqtot)) 267 ALLOCATE(pcvgu(klon,llm), pcvgv(klon,llm)) 268 ALLOCATE(pcvgt(klon,llm), pcvgq(klon,llm,2)) 269 ALLOCATE(zdufi(klon,llm),zdvfi(klon,llm)) 270 ALLOCATE(zdtfi(klon,llm),zdqfi(klon,llm,nqtot)) 271 ALLOCATE(zdpsrf(klon)) 272 ALLOCATE(zdufi2(klon+iim,llm),zdvfi2(klon+iim,llm)) 273 ALLOCATE(flxwfi(klon,llm)) 274 ALLOCATE(zpk(klon,llm)) 275 !$OMP END MASTER 276 !$OMP BARRIER 277 ELSE 278 debut = .FALSE. 279 ENDIF 280 281 ! 282 ! 283 !----------------------------------------------------------------------- 284 ! 40. transformation des variables dynamiques en variables physiques: 285 ! --------------------------------------------------------------- 286 287 ! 41. pressions au sol (en Pascals) 288 ! ---------------------------------- 289 290 !$OMP MASTER 291 call start_timer(timer_physic) 292 !$OMP END MASTER 293 294 !$OMP MASTER 295 !CDIR ON_ADB(index_i) 296 !CDIR ON_ADB(index_j) 297 do ig0=1,klon 298 i=index_i(ig0) 299 j=index_j(ig0) 300 zpsrf(ig0)=pps(i,j) 301 enddo 302 !$OMP END MASTER 303 304 305 ! 42. pression intercouches : 306 ! 307 ! ----------------------------------------------------------------- 308 ! .... zplev definis aux (llm +1) interfaces des couches .... 309 ! .... zplay definis aux ( llm ) milieux des couches .... 310 ! ----------------------------------------------------------------- 311 312 ! ... Exner = cp * ( p(l) / preff ) ** kappa .... 313 ! 314 unskap = 1./ kappa 315 ! 316 ! print *,omp_rank,'klon--->',klon 317 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 318 DO l = 1, llmp1 319 !CDIR ON_ADB(index_i) 320 !CDIR ON_ADB(index_j) 321 do ig0=1,klon 322 i=index_i(ig0) 323 j=index_j(ig0) 324 zplev( ig0,l ) = pp(i,j,l) 325 enddo 326 ENDDO 327 !$OMP END DO NOWAIT 328 329 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 330 DO l=1,llm 331 do ig0=1,klon 332 i=index_i(ig0) 333 j=index_j(ig0) 334 zpk(ig0,l)=ppk(i,j,l) 335 enddo 336 ENDDO 337 !$OMP END DO NOWAIT 338 339 ! 340 ! 341 342 ! 43. temperature naturelle (en K) et pressions milieux couches . 343 ! --------------------------------------------------------------- 344 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 345 DO l=1,llm 346 !CDIR ON_ADB(index_i) 347 !CDIR ON_ADB(index_j) 348 do ig0=1,klon 349 i=index_i(ig0) 350 j=index_j(ig0) 351 pksurcp = ppk(i,j,l) / cpp 352 zplay(ig0,l) = preff * pksurcp ** unskap 353 ztfi(ig0,l) = pteta(i,j,l) * pksurcp 354 enddo 355 356 ENDDO 357 !$OMP END DO NOWAIT 358 359 ! 43.bis traceurs 360 ! --------------- 361 ! 362 363 itr = 0 364 DO iq=1,nqtot 365 IF(.NOT.tracers(iq)%isAdvected) CYCLE 366 itr = itr + 1 367 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 368 DO l=1,llm 369 !CDIR ON_ADB(index_i) 370 !CDIR ON_ADB(index_j) 371 do ig0=1,klon 372 i=index_i(ig0) 373 j=index_j(ig0) 374 zqfi(ig0,l,itr) = pq(i,j,l,iq) 375 enddo 376 ENDDO 377 !$OMP END DO NOWAIT 378 ENDDO 379 380 381 ! Geopotentiel calcule par rapport a la surface locale: 382 ! ----------------------------------------------------- 383 384 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 385 DO l=1,llm 386 !CDIR ON_ADB(index_i) 387 !CDIR ON_ADB(index_j) 388 do ig0=1,klon 389 i=index_i(ig0) 390 j=index_j(ig0) 391 zphi(ig0,l) = pphi(i,j,l) 392 enddo 393 ENDDO 394 !$OMP END DO NOWAIT 395 396 ! CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,pphi,zphi) 397 398 !$OMP MASTER 399 !CDIR ON_ADB(index_i) 400 !CDIR ON_ADB(index_j) 401 do ig0=1,klon 402 i=index_i(ig0) 403 j=index_j(ig0) 404 zphis(ig0) = pphis(i,j) 405 enddo 406 !$OMP END MASTER 407 408 409 ! CALL gr_dyn_fi_p(1,iip1,jjp1,klon,pphis,zphis) 410 411 !$OMP BARRIER 412 413 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 414 DO l=1,llm 415 DO ig=1,klon 416 zphi(ig,l)=zphi(ig,l)-zphis(ig) 417 ENDDO 418 ENDDO 419 !$OMP END DO NOWAIT 420 421 422 ! 423 ! 45. champ u: 424 ! ------------ 425 426 kstart=1 427 kend=klon 428 429 if (is_north_pole_dyn) kstart=2 430 if (is_south_pole_dyn) kend=klon-1 431 432 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 433 DO l=1,llm 434 !CDIR ON_ADB(index_i) 435 !CDIR ON_ADB(index_j) 436 !CDIR SPARSE 437 do ig0=kstart,kend 438 i=index_i(ig0) 439 j=index_j(ig0) 440 if (i==1) then 441 zufi(ig0,l)= 0.5 *( pucov(iim,j,l)/cu(iim,j) & 442 + pucov(1,j,l)/cu(1,j) ) 443 else 444 zufi(ig0,l)= 0.5*( pucov(i-1,j,l)/cu(i-1,j) & 445 + pucov(i,j,l)/cu(i,j) ) 446 endif 447 enddo 448 ENDDO 449 !$OMP END DO NOWAIT 450 451 ! 452 ! Alvaro de la Camara (May 2014) 453 ! 46.1 Calcul de la vorticite et passage sur la grille physique 454 ! -------------------------------------------------------------- 455 456 jjb=jj_begin_dyn-1 457 jje=jj_end_dyn+1 458 if (is_north_pole_dyn) jjb=1 459 if (is_south_pole_dyn) jje=jjm 460 461 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 462 463 DO l=1,llm 464 do i=1,iim 465 do j=jjb,jje 466 zrot(i,j,l) = (pvcov(i+1,j,l) - pvcov(i,j,l) & 467 + pucov(i,j+1,l) - pucov(i,j,l)) & 468 / (cu(i,j)+cu(i,j+1)) & 469 / (cv(i+1,j)+cv(i,j)) *4 301 470 enddo 302 c$OMP END MASTER 303 304 305 c 42. pression intercouches : 306 c 307 c ----------------------------------------------------------------- 308 c .... zplev definis aux (llm +1) interfaces des couches .... 309 c .... zplay definis aux ( llm ) milieux des couches .... 310 c ----------------------------------------------------------------- 311 312 c ... Exner = cp * ( p(l) / preff ) ** kappa .... 313 c 314 unskap = 1./ kappa 315 c 316 c print *,omp_rank,'klon--->',klon 317 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 318 DO l = 1, llmp1 319 !CDIR ON_ADB(index_i) 320 !CDIR ON_ADB(index_j) 321 do ig0=1,klon 471 enddo 472 ENDDO 473 474 475 ! 46.2champ v: 476 ! ----------- 477 478 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 479 DO l=1,llm 480 !CDIR ON_ADB(index_i) 481 !CDIR ON_ADB(index_j) 482 DO ig0=kstart,kend 483 i=index_i(ig0) 484 j=index_j(ig0) 485 zvfi(ig0,l)= 0.5 *( pvcov(i,j-1,l)/cv(i,j-1) & 486 + pvcov(i,j,l)/cv(i,j) ) 487 if (j==1 .OR. j==jjp1) then ! AdlC MAY 2014 488 zrfi(ig0,l) = 0 ! AdlC MAY 2014 489 else 490 if(i==1)then 491 zrfi(ig0,l)= 0.25 *(zrot(iim,j-1,l)+zrot(iim,j,l) & 492 +zrot(1,j-1,l)+zrot(1,j,l)) ! AdlC MAY 2014 493 else 494 zrfi(ig0,l)= 0.25 *(zrot(i-1,j-1,l)+zrot(i-1,j,l) & 495 +zrot(i,j-1,l)+zrot(i,j,l)) ! AdlC MAY 2014 496 endif 497 endif 498 499 500 ENDDO 501 ENDDO 502 !$OMP END DO NOWAIT 503 504 ! 47. champs de vents aux pole nord 505 ! ------------------------------ 506 ! U = 1 / pi * integrale [ v * cos(long) * d long ] 507 ! V = 1 / pi * integrale [ v * sin(long) * d long ] 508 509 if (is_north_pole_dyn) then 510 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 511 DO l=1,llm 512 513 z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,1,l)/cv(1,1) 514 DO i=2,iim 515 z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1) 516 ENDDO 517 518 DO i=1,iim 519 zcos(i) = COS(rlonv(i))*z1(i) 520 zsin(i) = SIN(rlonv(i))*z1(i) 521 ENDDO 522 523 zufi(1,l) = SSUM(iim,zcos,1)/pi 524 zvfi(1,l) = SSUM(iim,zsin,1)/pi 525 zrfi(1,l) = 0. 526 527 ENDDO 528 !$OMP END DO NOWAIT 529 endif 530 531 532 ! 48. champs de vents aux pole sud: 533 ! --------------------------------- 534 ! U = 1 / pi * integrale [ v * cos(long) * d long ] 535 ! V = 1 / pi * integrale [ v * sin(long) * d long ] 536 537 if (is_south_pole_dyn) then 538 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 539 DO l=1,llm 540 541 z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,jjm,l)/cv(1,jjm) 542 DO i=2,iim 543 z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm) 544 ENDDO 545 546 DO i=1,iim 547 zcos(i) = COS(rlonv(i))*z1(i) 548 zsin(i) = SIN(rlonv(i))*z1(i) 549 ENDDO 550 551 zufi(klon,l) = SSUM(iim,zcos,1)/pi 552 zvfi(klon,l) = SSUM(iim,zsin,1)/pi 553 zrfi(klon,l) = 0. 554 ENDDO 555 !$OMP END DO NOWAIT 556 endif 557 558 ! On change de grille, dynamique vers physiq, pour le flux de masse verticale 559 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 560 DO l=1,llm 561 !CDIR ON_ADB(index_i) 562 !CDIR ON_ADB(index_j) 563 do ig0=1,klon 564 i=index_i(ig0) 565 j=index_j(ig0) 566 flxwfi(ig0,l) = flxw(i,j,l) 567 enddo 568 ENDDO 569 !$OMP END DO NOWAIT 570 571 ! CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,flxw,flxwfi) 572 573 !----------------------------------------------------------------------- 574 ! Appel de la physique: 575 ! --------------------- 576 577 578 !$OMP BARRIER 579 if (first_omp) then 580 klon=klon_omp 581 582 allocate(zplev_omp(klon,llm+1)) 583 allocate(zplay_omp(klon,llm)) 584 allocate(zpk_omp(klon,llm)) 585 allocate(zphi_omp(klon,llm)) 586 allocate(zphis_omp(klon)) 587 allocate(presnivs_omp(llm)) 588 allocate(zufi_omp(klon,llm)) 589 allocate(zvfi_omp(klon,llm)) 590 allocate(zrfi_omp(klon,llm)) ! LG Ari 2014 591 allocate(ztfi_omp(klon,llm)) 592 allocate(zqfi_omp(klon,llm,nqtot)) 593 allocate(zdufi_omp(klon,llm)) 594 allocate(zdvfi_omp(klon,llm)) 595 allocate(zdtfi_omp(klon,llm)) 596 allocate(zdqfi_omp(klon,llm,nqtot)) 597 allocate(zdufic_omp(klon,llm)) 598 allocate(zdvfic_omp(klon,llm)) 599 allocate(zdtfic_omp(klon,llm)) 600 allocate(zdqfic_omp(klon,llm,nqtot)) 601 allocate(zdpsrf_omp(klon)) 602 allocate(flxwfi_omp(klon,llm)) 603 first_omp=.false. 604 endif 605 606 607 klon=klon_omp 608 offset=klon_omp_begin-1 609 610 do l=1,llm+1 611 do i=1,klon 612 zplev_omp(i,l)=zplev(offset+i,l) 613 enddo 614 enddo 615 616 do l=1,llm 617 do i=1,klon 618 zplay_omp(i,l)=zplay(offset+i,l) 619 enddo 620 enddo 621 622 do l=1,llm 623 do i=1,klon 624 zpk_omp(i,l)=zpk(offset+i,l) 625 enddo 626 enddo 627 628 do l=1,llm 629 do i=1,klon 630 zphi_omp(i,l)=zphi(offset+i,l) 631 enddo 632 enddo 633 634 do i=1,klon 635 zphis_omp(i)=zphis(offset+i) 636 enddo 637 638 639 do l=1,llm 640 presnivs_omp(l)=presnivs(l) 641 enddo 642 643 do l=1,llm 644 do i=1,klon 645 zufi_omp(i,l)=zufi(offset+i,l) 646 enddo 647 enddo 648 649 do l=1,llm 650 do i=1,klon 651 zvfi_omp(i,l)=zvfi(offset+i,l) 652 enddo 653 enddo 654 655 do l=1,llm 656 do i=1,klon 657 zrfi_omp(i,l)=zrfi(offset+i,l) 658 enddo 659 enddo 660 661 do l=1,llm 662 do i=1,klon 663 ztfi_omp(i,l)=ztfi(offset+i,l) 664 enddo 665 enddo 666 667 do iq=1,nqtot 668 do l=1,llm 669 do i=1,klon 670 zqfi_omp(i,l,iq)=zqfi(offset+i,l,iq) 671 enddo 672 enddo 673 enddo 674 675 do l=1,llm 676 do i=1,klon 677 zdufi_omp(i,l)=zdufi(offset+i,l) 678 enddo 679 enddo 680 681 do l=1,llm 682 do i=1,klon 683 zdvfi_omp(i,l)=zdvfi(offset+i,l) 684 enddo 685 enddo 686 687 do l=1,llm 688 do i=1,klon 689 zdtfi_omp(i,l)=zdtfi(offset+i,l) 690 enddo 691 enddo 692 693 do iq=1,nqtot 694 do l=1,llm 695 do i=1,klon 696 zdqfi_omp(i,l,iq)=zdqfi(offset+i,l,iq) 697 enddo 698 enddo 699 enddo 700 701 do i=1,klon 702 zdpsrf_omp(i)=zdpsrf(offset+i) 703 enddo 704 705 do l=1,llm 706 do i=1,klon 707 flxwfi_omp(i,l)=flxwfi(offset+i,l) 708 enddo 709 enddo 710 711 !$OMP BARRIER 712 713 714 !$OMP MASTER 715 ! write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys 716 !$OMP END MASTER 717 zdt_split=dtphys/nsplit_phys 718 zdufic_omp(:,:)=0. 719 zdvfic_omp(:,:)=0. 720 zdtfic_omp(:,:)=0. 721 zdqfic_omp(:,:,:)=0. 722 723 #ifdef CPP_PHYS 724 do isplit=1,nsplit_phys 725 726 jH_cur_split=jH_cur+(isplit-1) * dtvr / (daysec *nsplit_phys) 727 debut_split=debut.and.isplit==1 728 lafin_split=lafin.and.isplit==nsplit_phys 729 730 CALL call_physiq(klon,llm,nqtot,tracers(:)%name, & 731 debut_split,lafin_split, & 732 jD_cur,jH_cur_split,zdt_split, & 733 zplev_omp,zplay_omp, & 734 zpk_omp,zphi_omp,zphis_omp, & 735 presnivs_omp, & 736 zufi_omp,zvfi_omp,zrfi_omp,ztfi_omp,zqfi_omp, & 737 flxwfi_omp,pducov, & 738 zdufi_omp,zdvfi_omp,zdtfi_omp,zdqfi_omp, & 739 zdpsrf_omp) 740 741 742 zufi_omp(:,:)=zufi_omp(:,:)+zdufi_omp(:,:)*zdt_split 743 zvfi_omp(:,:)=zvfi_omp(:,:)+zdvfi_omp(:,:)*zdt_split 744 ztfi_omp(:,:)=ztfi_omp(:,:)+zdtfi_omp(:,:)*zdt_split 745 zqfi_omp(:,:,:)=zqfi_omp(:,:,:)+zdqfi_omp(:,:,:)*zdt_split 746 747 zdufic_omp(:,:)=zdufic_omp(:,:)+zdufi_omp(:,:) 748 zdvfic_omp(:,:)=zdvfic_omp(:,:)+zdvfi_omp(:,:) 749 zdtfic_omp(:,:)=zdtfic_omp(:,:)+zdtfi_omp(:,:) 750 zdqfic_omp(:,:,:)=zdqfic_omp(:,:,:)+zdqfi_omp(:,:,:) 751 752 enddo 753 754 #endif 755 ! of #ifdef CPP_PHYS 756 757 758 zdufi_omp(:,:)=zdufic_omp(:,:)/nsplit_phys 759 zdvfi_omp(:,:)=zdvfic_omp(:,:)/nsplit_phys 760 zdtfi_omp(:,:)=zdtfic_omp(:,:)/nsplit_phys 761 zdqfi_omp(:,:,:)=zdqfic_omp(:,:,:)/nsplit_phys 762 763 !$OMP BARRIER 764 765 do l=1,llm+1 766 do i=1,klon 767 zplev(offset+i,l)=zplev_omp(i,l) 768 enddo 769 enddo 770 771 do l=1,llm 772 do i=1,klon 773 zplay(offset+i,l)=zplay_omp(i,l) 774 enddo 775 enddo 776 777 do l=1,llm 778 do i=1,klon 779 zphi(offset+i,l)=zphi_omp(i,l) 780 enddo 781 enddo 782 783 784 do i=1,klon 785 zphis(offset+i)=zphis_omp(i) 786 enddo 787 788 789 do l=1,llm 790 presnivs(l)=presnivs_omp(l) 791 enddo 792 793 do l=1,llm 794 do i=1,klon 795 zufi(offset+i,l)=zufi_omp(i,l) 796 enddo 797 enddo 798 799 do l=1,llm 800 do i=1,klon 801 zvfi(offset+i,l)=zvfi_omp(i,l) 802 enddo 803 enddo 804 805 do l=1,llm 806 do i=1,klon 807 ztfi(offset+i,l)=ztfi_omp(i,l) 808 enddo 809 enddo 810 811 do iq=1,nqtot 812 do l=1,llm 813 do i=1,klon 814 zqfi(offset+i,l,iq)=zqfi_omp(i,l,iq) 815 enddo 816 enddo 817 enddo 818 819 do l=1,llm 820 do i=1,klon 821 zdufi(offset+i,l)=zdufi_omp(i,l) 822 enddo 823 enddo 824 825 do l=1,llm 826 do i=1,klon 827 zdvfi(offset+i,l)=zdvfi_omp(i,l) 828 enddo 829 enddo 830 831 do l=1,llm 832 do i=1,klon 833 zdtfi(offset+i,l)=zdtfi_omp(i,l) 834 enddo 835 enddo 836 837 do iq=1,nqtot 838 do l=1,llm 839 do i=1,klon 840 zdqfi(offset+i,l,iq)=zdqfi_omp(i,l,iq) 841 enddo 842 enddo 843 enddo 844 845 do i=1,klon 846 zdpsrf(offset+i)=zdpsrf_omp(i) 847 enddo 848 849 850 klon=klon_mpi 851 500 CONTINUE 852 !$OMP BARRIER 853 854 !$OMP MASTER 855 call stop_timer(timer_physic) 856 !$OMP END MASTER 857 858 IF (using_mpi) THEN 859 860 if (MPI_rank>0) then 861 862 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 863 DO l=1,llm 864 du_send(1:iim,l)=zdufi(1:iim,l) 865 dv_send(1:iim,l)=zdvfi(1:iim,l) 866 ENDDO 867 !$OMP END DO NOWAIT 868 869 !$OMP BARRIER 870 871 !$OMP MASTER 872 !$OMP CRITICAL (MPI) 873 call MPI_ISSEND(du_send,iim*llm,MPI_REAL8,MPI_Rank-1,401, & 874 COMM_LMDZ,Req(1),ierr) 875 call MPI_ISSEND(dv_send,iim*llm,MPI_REAL8,MPI_Rank-1,402, & 876 COMM_LMDZ,Req(2),ierr) 877 !$OMP END CRITICAL (MPI) 878 !$OMP END MASTER 879 880 !$OMP BARRIER 881 882 endif 883 884 if (MPI_rank<MPI_Size-1) then 885 !$OMP BARRIER 886 887 !$OMP MASTER 888 !$OMP CRITICAL (MPI) 889 call MPI_IRECV(du_recv,iim*llm,MPI_REAL8,MPI_Rank+1,401, & 890 COMM_LMDZ,Req(3),ierr) 891 call MPI_IRECV(dv_recv,iim*llm,MPI_REAL8,MPI_Rank+1,402, & 892 COMM_LMDZ,Req(4),ierr) 893 !$OMP END CRITICAL (MPI) 894 !$OMP END MASTER 895 896 endif 897 898 !$OMP BARRIER 899 900 901 !$OMP MASTER 902 !$OMP CRITICAL (MPI) 903 if (MPI_rank>0 .and. MPI_rank< MPI_Size-1) then 904 call MPI_WAITALL(4,Req(1),Status,ierr) 905 else if (MPI_rank>0) then 906 call MPI_WAITALL(2,Req(1),Status,ierr) 907 else if (MPI_rank <MPI_Size-1) then 908 call MPI_WAITALL(2,Req(3),Status,ierr) 909 endif 910 !$OMP END CRITICAL (MPI) 911 !$OMP END MASTER 912 913 !$OMP BARRIER 914 915 ENDIF ! using_mpi 916 917 918 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 919 DO l=1,llm 920 921 zdufi2(1:klon,l)=zdufi(1:klon,l) 922 zdufi2(klon+1:klon+iim,l)=du_recv(1:iim,l) 923 924 zdvfi2(1:klon,l)=zdvfi(1:klon,l) 925 zdvfi2(klon+1:klon+iim,l)=dv_recv(1:iim,l) 926 927 pdhfi(:,jj_begin,l)=0 928 pdqfi(:,jj_begin,l,:)=0 929 pdufi(:,jj_begin,l)=0 930 pdvfi(:,jj_begin,l)=0 931 932 if (.not. is_south_pole_dyn) then 933 pdhfi(:,jj_end:jj_end+1,l)=0 934 pdqfi(:,jj_end:jj_end+1,l,:)=0 935 pdufi(:,jj_end:jj_end+1,l)=0 936 pdvfi(:,jj_end:jj_end+1,l)=0 937 endif 938 939 ENDDO 940 !$OMP END DO NOWAIT 941 942 !$OMP MASTER 943 pdpsfi(:,jj_begin)=0 944 945 if (.not. is_south_pole_dyn) then 946 pdpsfi(:,jj_end:jj_end+1)=0 947 endif 948 !$OMP END MASTER 949 !----------------------------------------------------------------------- 950 ! transformation des tendances physiques en tendances dynamiques: 951 ! --------------------------------------------------------------- 952 953 ! tendance sur la pression : 954 ! ----------------------------------- 955 ! CALL gr_fi_dyn_p(1,klon,iip1,jjp1,zdpsrf,pdpsfi) 956 957 !$OMP MASTER 958 kstart=1 959 kend=klon 960 961 if (is_north_pole_dyn) kstart=2 962 if (is_south_pole_dyn) kend=klon-1 963 964 !CDIR ON_ADB(index_i) 965 !CDIR ON_ADB(index_j) 966 !cdir NODEP 967 do ig0=kstart,kend 968 i=index_i(ig0) 969 j=index_j(ig0) 970 pdpsfi(i,j) = zdpsrf(ig0) 971 if (i==1) pdpsfi(iip1,j) = zdpsrf(ig0) 972 enddo 973 974 if (is_north_pole_dyn) then 975 DO i=1,iip1 976 pdpsfi(i,1) = zdpsrf(1) 977 enddo 978 endif 979 980 if (is_south_pole_dyn) then 981 DO i=1,iip1 982 pdpsfi(i,jjp1) = zdpsrf(klon) 983 ENDDO 984 endif 985 !$OMP END MASTER 986 !c$OMP BARRIER 987 988 ! 989 ! 62. enthalpie potentielle 990 ! --------------------- 991 992 kstart=1 993 kend=klon 994 995 if (is_north_pole_dyn) kstart=2 996 if (is_south_pole_dyn) kend=klon-1 997 998 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 999 DO l=1,llm 1000 1001 !CDIR ON_ADB(index_i) 1002 !CDIR ON_ADB(index_j) 1003 !cdir NODEP 1004 do ig0=kstart,kend 1005 i=index_i(ig0) 1006 j=index_j(ig0) 1007 pdhfi(i,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l) 1008 if (i==1) pdhfi(iip1,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l) 1009 enddo 1010 1011 if (is_north_pole_dyn) then 1012 DO i=1,iip1 1013 pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l) 1014 enddo 1015 endif 1016 1017 if (is_south_pole_dyn) then 1018 DO i=1,iip1 1019 pdhfi(i,jjp1,l) = cpp * zdtfi(klon,l)/ ppk(i,jjp1,l) 1020 ENDDO 1021 endif 1022 ENDDO 1023 !$OMP END DO NOWAIT 1024 1025 ! 62. humidite specifique 1026 ! --------------------- 1027 ! Ehouarn: removed this useless bit: was overwritten at step 63 anyways 1028 ! DO iq=1,nqtot 1029 !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1030 ! DO l=1,llm 1031 !!!cdir NODEP 1032 ! do ig0=kstart,kend 1033 ! i=index_i(ig0) 1034 ! j=index_j(ig0) 1035 ! pdqfi(i,j,l,iq) = zdqfi(ig0,l,iq) 1036 ! if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,iq) 1037 ! enddo 1038 ! 1039 ! if (is_north_pole_dyn) then 1040 ! do i=1,iip1 1041 ! pdqfi(i,1,l,iq) = zdqfi(1,l,iq) 1042 ! enddo 1043 ! endif 1044 ! 1045 ! if (is_south_pole_dyn) then 1046 ! do i=1,iip1 1047 ! pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,iq) 1048 ! enddo 1049 ! endif 1050 ! ENDDO 1051 !c$OMP END DO NOWAIT 1052 ! ENDDO 1053 1054 ! 63. traceurs 1055 ! ------------ 1056 ! initialisation des tendances 1057 1058 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1059 DO l=1,llm 1060 pdqfi(:,jj_begin:jj_end,l,:)=0. 1061 ENDDO 1062 !$OMP END DO NOWAIT 1063 1064 ! 1065 !cdir NODEP 1066 itr = 0 1067 DO iq=1,nqtot 1068 IF(.NOT.tracers(iq)%isAdvected) CYCLE 1069 itr = itr + 1 1070 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1071 DO l=1,llm 1072 !CDIR ON_ADB(index_i) 1073 !CDIR ON_ADB(index_j) 1074 !cdir NODEP 1075 DO ig0=kstart,kend 322 1076 i=index_i(ig0) 323 1077 j=index_j(ig0) 324 zplev( ig0,l ) = pp(i,j,l) 325 enddo 1078 pdqfi(i,j,l,iq) = zdqfi(ig0,l,itr) 1079 if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,itr) 1080 ENDDO 1081 1082 IF (is_north_pole_dyn) then 1083 DO i=1,iip1 1084 pdqfi(i,1,l,iq) = zdqfi(1,l,itr) 1085 ENDDO 1086 ENDIF 1087 1088 IF (is_south_pole_dyn) then 1089 DO i=1,iip1 1090 pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,itr) 1091 ENDDO 1092 ENDIF 1093 1094 ENDDO 1095 !$OMP END DO NOWAIT 1096 ENDDO 1097 1098 ! 65. champ u: 1099 ! ------------ 1100 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1101 DO l=1,llm 1102 !CDIR ON_ADB(index_i) 1103 !CDIR ON_ADB(index_j) 1104 !cdir NODEP 1105 do ig0=kstart,kend 1106 i=index_i(ig0) 1107 j=index_j(ig0) 1108 1109 if (i/=iim) then 1110 pdufi(i,j,l)=0.5*(zdufi2(ig0,l)+zdufi2(ig0+1,l))*cu(i,j) 1111 endif 1112 1113 if (i==1) then 1114 pdufi(iim,j,l)=0.5*( zdufi2(ig0,l) & 1115 + zdufi2(ig0+iim-1,l))*cu(iim,j) 1116 pdufi(iip1,j,l)=0.5*(zdufi2(ig0,l)+zdufi2(ig0+1,l))*cu(i,j) 1117 endif 1118 1119 enddo 1120 1121 if (is_north_pole_dyn) then 1122 DO i=1,iip1 1123 pdufi(i,1,l) = 0. 1124 ENDDO 1125 endif 1126 1127 if (is_south_pole_dyn) then 1128 DO i=1,iip1 1129 pdufi(i,jjp1,l) = 0. 1130 ENDDO 1131 endif 1132 1133 ENDDO 1134 !$OMP END DO NOWAIT 1135 1136 ! 67. champ v: 1137 ! ------------ 1138 1139 kstart=1 1140 kend=klon 1141 1142 if (is_north_pole_dyn) kstart=2 1143 if (is_south_pole_dyn) kend=klon-1-iim 1144 1145 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1146 DO l=1,llm 1147 !CDIR ON_ADB(index_i) 1148 !CDIR ON_ADB(index_j) 1149 !cdir NODEP 1150 do ig0=kstart,kend 1151 i=index_i(ig0) 1152 j=index_j(ig0) 1153 pdvfi(i,j,l)=0.5*(zdvfi2(ig0,l)+zdvfi2(ig0+iim,l))*cv(i,j) 1154 if (i==1) pdvfi(iip1,j,l) = 0.5*(zdvfi2(ig0,l)+ & 1155 zdvfi2(ig0+iim,l)) & 1156 *cv(i,j) 1157 enddo 1158 1159 ENDDO 1160 !$OMP END DO NOWAIT 1161 1162 1163 ! 68. champ v pres des poles: 1164 ! --------------------------- 1165 ! v = U * cos(long) + V * SIN(long) 1166 1167 if (is_north_pole_dyn) then 1168 1169 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1170 DO l=1,llm 1171 1172 DO i=1,iim 1173 pdvfi(i,1,l)= & 1174 zdufi(1,l)*COS(rlonv(i))+zdvfi(1,l)*SIN(rlonv(i)) 1175 1176 pdvfi(i,1,l)= & 1177 0.5*(pdvfi(i,1,l)+zdvfi(i+1,l))*cv(i,1) 326 1178 ENDDO 327 c$OMP END DO NOWAIT 328 329 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 330 DO l=1,llm 331 do ig0=1,klon 332 i=index_i(ig0) 333 j=index_j(ig0) 334 zpk(ig0,l)=ppk(i,j,l) 335 enddo 336 ENDDO 337 c$OMP END DO NOWAIT 338 339 c 340 c 341 342 c 43. temperature naturelle (en K) et pressions milieux couches . 343 c --------------------------------------------------------------- 344 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 345 DO l=1,llm 346 !CDIR ON_ADB(index_i) 347 !CDIR ON_ADB(index_j) 348 do ig0=1,klon 349 i=index_i(ig0) 350 j=index_j(ig0) 351 pksurcp = ppk(i,j,l) / cpp 352 zplay(ig0,l) = preff * pksurcp ** unskap 353 ztfi(ig0,l) = pteta(i,j,l) * pksurcp 354 enddo 355 356 ENDDO 357 c$OMP END DO NOWAIT 358 359 c 43.bis traceurs 360 c --------------- 361 c 362 363 itr = 0 364 DO iq=1,nqtot 365 IF(.NOT.tracers(iq)%isAdvected) CYCLE 366 itr = itr + 1 367 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 368 DO l=1,llm 369 !CDIR ON_ADB(index_i) 370 !CDIR ON_ADB(index_j) 371 do ig0=1,klon 372 i=index_i(ig0) 373 j=index_j(ig0) 374 zqfi(ig0,l,itr) = pq(i,j,l,iq) 375 enddo 376 ENDDO 377 c$OMP END DO NOWAIT 378 ENDDO 379 380 381 c Geopotentiel calcule par rapport a la surface locale: 382 c ----------------------------------------------------- 383 384 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 385 DO l=1,llm 386 !CDIR ON_ADB(index_i) 387 !CDIR ON_ADB(index_j) 388 do ig0=1,klon 389 i=index_i(ig0) 390 j=index_j(ig0) 391 zphi(ig0,l) = pphi(i,j,l) 392 enddo 393 ENDDO 394 c$OMP END DO NOWAIT 395 396 c CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,pphi,zphi) 397 398 c$OMP MASTER 399 !CDIR ON_ADB(index_i) 400 !CDIR ON_ADB(index_j) 401 do ig0=1,klon 402 i=index_i(ig0) 403 j=index_j(ig0) 404 zphis(ig0) = pphis(i,j) 405 enddo 406 c$OMP END MASTER 407 408 409 c CALL gr_dyn_fi_p(1,iip1,jjp1,klon,pphis,zphis) 410 411 c$OMP BARRIER 412 413 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 414 DO l=1,llm 415 DO ig=1,klon 416 zphi(ig,l)=zphi(ig,l)-zphis(ig) 417 ENDDO 418 ENDDO 419 c$OMP END DO NOWAIT 420 421 422 c 423 c 45. champ u: 424 c ------------ 425 426 kstart=1 427 kend=klon 428 429 if (is_north_pole_dyn) kstart=2 430 if (is_south_pole_dyn) kend=klon-1 431 432 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 433 DO l=1,llm 434 !CDIR ON_ADB(index_i) 435 !CDIR ON_ADB(index_j) 436 !CDIR SPARSE 437 do ig0=kstart,kend 438 i=index_i(ig0) 439 j=index_j(ig0) 440 if (i==1) then 441 zufi(ig0,l)= 0.5 *( pucov(iim,j,l)/cu(iim,j) 442 $ + pucov(1,j,l)/cu(1,j) ) 443 else 444 zufi(ig0,l)= 0.5*( pucov(i-1,j,l)/cu(i-1,j) 445 $ + pucov(i,j,l)/cu(i,j) ) 446 endif 447 enddo 448 ENDDO 449 c$OMP END DO NOWAIT 450 451 c 452 C Alvaro de la Camara (May 2014) 453 C 46.1 Calcul de la vorticite et passage sur la grille physique 454 C -------------------------------------------------------------- 455 456 jjb=jj_begin_dyn-1 457 jje=jj_end_dyn+1 458 if (is_north_pole_dyn) jjb=1 459 if (is_south_pole_dyn) jje=jjm 460 461 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 462 463 DO l=1,llm 464 do i=1,iim 465 do j=jjb,jje 466 zrot(i,j,l) = (pvcov(i+1,j,l) - pvcov(i,j,l) 467 $ + pucov(i,j+1,l) - pucov(i,j,l)) 468 $ / (cu(i,j)+cu(i,j+1)) 469 $ / (cv(i+1,j)+cv(i,j)) *4 470 enddo 471 enddo 472 ENDDO 473 474 475 c 46.2champ v: 476 c ----------- 477 478 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 479 DO l=1,llm 480 !CDIR ON_ADB(index_i) 481 !CDIR ON_ADB(index_j) 482 DO ig0=kstart,kend 483 i=index_i(ig0) 484 j=index_j(ig0) 485 zvfi(ig0,l)= 0.5 *( pvcov(i,j-1,l)/cv(i,j-1) 486 $ + pvcov(i,j,l)/cv(i,j) ) 487 if (j==1 .OR. j==jjp1) then ! AdlC MAY 2014 488 zrfi(ig0,l) = 0 ! AdlC MAY 2014 489 else 490 if(i==1)then 491 zrfi(ig0,l)= 0.25 *(zrot(iim,j-1,l)+zrot(iim,j,l) 492 $ +zrot(1,j-1,l)+zrot(1,j,l)) ! AdlC MAY 2014 493 else 494 zrfi(ig0,l)= 0.25 *(zrot(i-1,j-1,l)+zrot(i-1,j,l) 495 $ +zrot(i,j-1,l)+zrot(i,j,l)) ! AdlC MAY 2014 496 endif 497 endif 498 499 500 ENDDO 501 ENDDO 502 c$OMP END DO NOWAIT 503 504 c 47. champs de vents aux pole nord 505 c ------------------------------ 506 c U = 1 / pi * integrale [ v * cos(long) * d long ] 507 c V = 1 / pi * integrale [ v * sin(long) * d long ] 508 509 if (is_north_pole_dyn) then 510 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 511 DO l=1,llm 512 513 z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,1,l)/cv(1,1) 514 DO i=2,iim 515 z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1) 516 ENDDO 517 518 DO i=1,iim 519 zcos(i) = COS(rlonv(i))*z1(i) 520 zsin(i) = SIN(rlonv(i))*z1(i) 521 ENDDO 522 523 zufi(1,l) = SSUM(iim,zcos,1)/pi 524 zvfi(1,l) = SSUM(iim,zsin,1)/pi 525 zrfi(1,l) = 0. 526 527 ENDDO 528 c$OMP END DO NOWAIT 529 endif 530 531 532 c 48. champs de vents aux pole sud: 533 c --------------------------------- 534 c U = 1 / pi * integrale [ v * cos(long) * d long ] 535 c V = 1 / pi * integrale [ v * sin(long) * d long ] 536 537 if (is_south_pole_dyn) then 538 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 539 DO l=1,llm 540 541 z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,jjm,l)/cv(1,jjm) 542 DO i=2,iim 543 z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm) 544 ENDDO 545 546 DO i=1,iim 547 zcos(i) = COS(rlonv(i))*z1(i) 548 zsin(i) = SIN(rlonv(i))*z1(i) 549 ENDDO 550 551 zufi(klon,l) = SSUM(iim,zcos,1)/pi 552 zvfi(klon,l) = SSUM(iim,zsin,1)/pi 553 zrfi(klon,l) = 0. 554 ENDDO 555 c$OMP END DO NOWAIT 556 endif 557 558 c On change de grille, dynamique vers physiq, pour le flux de masse verticale 559 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 560 DO l=1,llm 561 !CDIR ON_ADB(index_i) 562 !CDIR ON_ADB(index_j) 563 do ig0=1,klon 564 i=index_i(ig0) 565 j=index_j(ig0) 566 flxwfi(ig0,l) = flxw(i,j,l) 567 enddo 568 ENDDO 569 c$OMP END DO NOWAIT 570 571 c CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,flxw,flxwfi) 572 573 c----------------------------------------------------------------------- 574 c Appel de la physique: 575 c --------------------- 576 577 578 c$OMP BARRIER 579 if (first_omp) then 580 klon=klon_omp 581 582 allocate(zplev_omp(klon,llm+1)) 583 allocate(zplay_omp(klon,llm)) 584 allocate(zpk_omp(klon,llm)) 585 allocate(zphi_omp(klon,llm)) 586 allocate(zphis_omp(klon)) 587 allocate(presnivs_omp(llm)) 588 allocate(zufi_omp(klon,llm)) 589 allocate(zvfi_omp(klon,llm)) 590 allocate(zrfi_omp(klon,llm)) ! LG Ari 2014 591 allocate(ztfi_omp(klon,llm)) 592 allocate(zqfi_omp(klon,llm,nqtot)) 593 allocate(zdufi_omp(klon,llm)) 594 allocate(zdvfi_omp(klon,llm)) 595 allocate(zdtfi_omp(klon,llm)) 596 allocate(zdqfi_omp(klon,llm,nqtot)) 597 allocate(zdufic_omp(klon,llm)) 598 allocate(zdvfic_omp(klon,llm)) 599 allocate(zdtfic_omp(klon,llm)) 600 allocate(zdqfic_omp(klon,llm,nqtot)) 601 allocate(zdpsrf_omp(klon)) 602 allocate(flxwfi_omp(klon,llm)) 603 first_omp=.false. 604 endif 605 606 607 klon=klon_omp 608 offset=klon_omp_begin-1 609 610 do l=1,llm+1 611 do i=1,klon 612 zplev_omp(i,l)=zplev(offset+i,l) 613 enddo 614 enddo 615 616 do l=1,llm 617 do i=1,klon 618 zplay_omp(i,l)=zplay(offset+i,l) 619 enddo 620 enddo 621 622 do l=1,llm 623 do i=1,klon 624 zpk_omp(i,l)=zpk(offset+i,l) 625 enddo 626 enddo 627 628 do l=1,llm 629 do i=1,klon 630 zphi_omp(i,l)=zphi(offset+i,l) 631 enddo 632 enddo 633 634 do i=1,klon 635 zphis_omp(i)=zphis(offset+i) 636 enddo 637 638 639 do l=1,llm 640 presnivs_omp(l)=presnivs(l) 641 enddo 642 643 do l=1,llm 644 do i=1,klon 645 zufi_omp(i,l)=zufi(offset+i,l) 646 enddo 647 enddo 648 649 do l=1,llm 650 do i=1,klon 651 zvfi_omp(i,l)=zvfi(offset+i,l) 652 enddo 653 enddo 654 655 do l=1,llm 656 do i=1,klon 657 zrfi_omp(i,l)=zrfi(offset+i,l) 658 enddo 659 enddo 660 661 do l=1,llm 662 do i=1,klon 663 ztfi_omp(i,l)=ztfi(offset+i,l) 664 enddo 665 enddo 666 667 do iq=1,nqtot 668 do l=1,llm 669 do i=1,klon 670 zqfi_omp(i,l,iq)=zqfi(offset+i,l,iq) 671 enddo 672 enddo 673 enddo 674 675 do l=1,llm 676 do i=1,klon 677 zdufi_omp(i,l)=zdufi(offset+i,l) 678 enddo 679 enddo 680 681 do l=1,llm 682 do i=1,klon 683 zdvfi_omp(i,l)=zdvfi(offset+i,l) 684 enddo 685 enddo 686 687 do l=1,llm 688 do i=1,klon 689 zdtfi_omp(i,l)=zdtfi(offset+i,l) 690 enddo 691 enddo 692 693 do iq=1,nqtot 694 do l=1,llm 695 do i=1,klon 696 zdqfi_omp(i,l,iq)=zdqfi(offset+i,l,iq) 697 enddo 698 enddo 699 enddo 700 701 do i=1,klon 702 zdpsrf_omp(i)=zdpsrf(offset+i) 703 enddo 704 705 do l=1,llm 706 do i=1,klon 707 flxwfi_omp(i,l)=flxwfi(offset+i,l) 708 enddo 709 enddo 710 711 c$OMP BARRIER 712 713 714 !$OMP MASTER 715 ! write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys 716 !$OMP END MASTER 717 zdt_split=dtphys/nsplit_phys 718 zdufic_omp(:,:)=0. 719 zdvfic_omp(:,:)=0. 720 zdtfic_omp(:,:)=0. 721 zdqfic_omp(:,:,:)=0. 722 723 #ifdef CPP_PHYS 724 do isplit=1,nsplit_phys 725 726 jH_cur_split=jH_cur+(isplit-1) * dtvr / (daysec *nsplit_phys) 727 debut_split=debut.and.isplit==1 728 lafin_split=lafin.and.isplit==nsplit_phys 729 730 CALL call_physiq(klon,llm,nqtot,tracers(:)%name, 731 & debut_split,lafin_split, 732 & jD_cur,jH_cur_split,zdt_split, 733 & zplev_omp,zplay_omp, 734 & zpk_omp,zphi_omp,zphis_omp, 735 & presnivs_omp, 736 & zufi_omp,zvfi_omp,zrfi_omp,ztfi_omp,zqfi_omp, 737 & flxwfi_omp,pducov, 738 & zdufi_omp,zdvfi_omp,zdtfi_omp,zdqfi_omp, 739 & zdpsrf_omp) 740 741 742 zufi_omp(:,:)=zufi_omp(:,:)+zdufi_omp(:,:)*zdt_split 743 zvfi_omp(:,:)=zvfi_omp(:,:)+zdvfi_omp(:,:)*zdt_split 744 ztfi_omp(:,:)=ztfi_omp(:,:)+zdtfi_omp(:,:)*zdt_split 745 zqfi_omp(:,:,:)=zqfi_omp(:,:,:)+zdqfi_omp(:,:,:)*zdt_split 746 747 zdufic_omp(:,:)=zdufic_omp(:,:)+zdufi_omp(:,:) 748 zdvfic_omp(:,:)=zdvfic_omp(:,:)+zdvfi_omp(:,:) 749 zdtfic_omp(:,:)=zdtfic_omp(:,:)+zdtfi_omp(:,:) 750 zdqfic_omp(:,:,:)=zdqfic_omp(:,:,:)+zdqfi_omp(:,:,:) 751 752 enddo 753 1179 1180 pdvfi(iip1,1,l) = pdvfi(1,1,l) 1181 1182 ENDDO 1183 !$OMP END DO NOWAIT 1184 1185 endif 1186 1187 if (is_south_pole_dyn) then 1188 1189 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1190 DO l=1,llm 1191 1192 DO i=1,iim 1193 pdvfi(i,jjm,l)=zdufi(klon,l)*COS(rlonv(i)) & 1194 +zdvfi(klon,l)*SIN(rlonv(i)) 1195 1196 pdvfi(i,jjm,l)= & 1197 0.5*(pdvfi(i,jjm,l)+zdvfi(klon-iip1+i,l))*cv(i,jjm) 1198 ENDDO 1199 1200 pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l) 1201 1202 ENDDO 1203 !$OMP END DO NOWAIT 1204 1205 endif 1206 !----------------------------------------------------------------------- 1207 1208 700 CONTINUE 1209 1210 firstcal = .FALSE. 1211 1212 #else 1213 call abort_gcm("calfis_loc", & 1214 "calfis_p: for now can only work with parallel physics", 1) 754 1215 #endif 755 ! of #ifdef CPP_PHYS 756 757 758 zdufi_omp(:,:)=zdufic_omp(:,:)/nsplit_phys 759 zdvfi_omp(:,:)=zdvfic_omp(:,:)/nsplit_phys 760 zdtfi_omp(:,:)=zdtfic_omp(:,:)/nsplit_phys 761 zdqfi_omp(:,:,:)=zdqfic_omp(:,:,:)/nsplit_phys 762 763 c$OMP BARRIER 764 765 do l=1,llm+1 766 do i=1,klon 767 zplev(offset+i,l)=zplev_omp(i,l) 768 enddo 769 enddo 770 771 do l=1,llm 772 do i=1,klon 773 zplay(offset+i,l)=zplay_omp(i,l) 774 enddo 775 enddo 776 777 do l=1,llm 778 do i=1,klon 779 zphi(offset+i,l)=zphi_omp(i,l) 780 enddo 781 enddo 782 783 784 do i=1,klon 785 zphis(offset+i)=zphis_omp(i) 786 enddo 787 788 789 do l=1,llm 790 presnivs(l)=presnivs_omp(l) 791 enddo 792 793 do l=1,llm 794 do i=1,klon 795 zufi(offset+i,l)=zufi_omp(i,l) 796 enddo 797 enddo 798 799 do l=1,llm 800 do i=1,klon 801 zvfi(offset+i,l)=zvfi_omp(i,l) 802 enddo 803 enddo 804 805 do l=1,llm 806 do i=1,klon 807 ztfi(offset+i,l)=ztfi_omp(i,l) 808 enddo 809 enddo 810 811 do iq=1,nqtot 812 do l=1,llm 813 do i=1,klon 814 zqfi(offset+i,l,iq)=zqfi_omp(i,l,iq) 815 enddo 816 enddo 817 enddo 818 819 do l=1,llm 820 do i=1,klon 821 zdufi(offset+i,l)=zdufi_omp(i,l) 822 enddo 823 enddo 824 825 do l=1,llm 826 do i=1,klon 827 zdvfi(offset+i,l)=zdvfi_omp(i,l) 828 enddo 829 enddo 830 831 do l=1,llm 832 do i=1,klon 833 zdtfi(offset+i,l)=zdtfi_omp(i,l) 834 enddo 835 enddo 836 837 do iq=1,nqtot 838 do l=1,llm 839 do i=1,klon 840 zdqfi(offset+i,l,iq)=zdqfi_omp(i,l,iq) 841 enddo 842 enddo 843 enddo 844 845 do i=1,klon 846 zdpsrf(offset+i)=zdpsrf_omp(i) 847 enddo 848 849 850 klon=klon_mpi 851 500 CONTINUE 852 c$OMP BARRIER 853 854 c$OMP MASTER 855 call stop_timer(timer_physic) 856 c$OMP END MASTER 857 858 IF (using_mpi) THEN 859 860 if (MPI_rank>0) then 861 862 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 863 DO l=1,llm 864 du_send(1:iim,l)=zdufi(1:iim,l) 865 dv_send(1:iim,l)=zdvfi(1:iim,l) 866 ENDDO 867 c$OMP END DO NOWAIT 868 869 c$OMP BARRIER 870 871 c$OMP MASTER 872 !$OMP CRITICAL (MPI) 873 call MPI_ISSEND(du_send,iim*llm,MPI_REAL8,MPI_Rank-1,401, 874 & COMM_LMDZ,Req(1),ierr) 875 call MPI_ISSEND(dv_send,iim*llm,MPI_REAL8,MPI_Rank-1,402, 876 & COMM_LMDZ,Req(2),ierr) 877 !$OMP END CRITICAL (MPI) 878 c$OMP END MASTER 879 880 c$OMP BARRIER 881 882 endif 883 884 if (MPI_rank<MPI_Size-1) then 885 c$OMP BARRIER 886 887 c$OMP MASTER 888 !$OMP CRITICAL (MPI) 889 call MPI_IRECV(du_recv,iim*llm,MPI_REAL8,MPI_Rank+1,401, 890 & COMM_LMDZ,Req(3),ierr) 891 call MPI_IRECV(dv_recv,iim*llm,MPI_REAL8,MPI_Rank+1,402, 892 & COMM_LMDZ,Req(4),ierr) 893 !$OMP END CRITICAL (MPI) 894 c$OMP END MASTER 895 896 endif 897 898 c$OMP BARRIER 899 900 901 c$OMP MASTER 902 !$OMP CRITICAL (MPI) 903 if (MPI_rank>0 .and. MPI_rank< MPI_Size-1) then 904 call MPI_WAITALL(4,Req(1),Status,ierr) 905 else if (MPI_rank>0) then 906 call MPI_WAITALL(2,Req(1),Status,ierr) 907 else if (MPI_rank <MPI_Size-1) then 908 call MPI_WAITALL(2,Req(3),Status,ierr) 909 endif 910 !$OMP END CRITICAL (MPI) 911 c$OMP END MASTER 912 913 c$OMP BARRIER 914 915 ENDIF ! using_mpi 916 917 918 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 919 DO l=1,llm 920 921 zdufi2(1:klon,l)=zdufi(1:klon,l) 922 zdufi2(klon+1:klon+iim,l)=du_recv(1:iim,l) 923 924 zdvfi2(1:klon,l)=zdvfi(1:klon,l) 925 zdvfi2(klon+1:klon+iim,l)=dv_recv(1:iim,l) 926 927 pdhfi(:,jj_begin,l)=0 928 pdqfi(:,jj_begin,l,:)=0 929 pdufi(:,jj_begin,l)=0 930 pdvfi(:,jj_begin,l)=0 931 932 if (.not. is_south_pole_dyn) then 933 pdhfi(:,jj_end:jj_end+1,l)=0 934 pdqfi(:,jj_end:jj_end+1,l,:)=0 935 pdufi(:,jj_end:jj_end+1,l)=0 936 pdvfi(:,jj_end:jj_end+1,l)=0 937 endif 938 939 ENDDO 940 c$OMP END DO NOWAIT 941 942 c$OMP MASTER 943 pdpsfi(:,jj_begin)=0 944 945 if (.not. is_south_pole_dyn) then 946 pdpsfi(:,jj_end:jj_end+1)=0 947 endif 948 c$OMP END MASTER 949 c----------------------------------------------------------------------- 950 c transformation des tendances physiques en tendances dynamiques: 951 c --------------------------------------------------------------- 952 953 c tendance sur la pression : 954 c ----------------------------------- 955 c CALL gr_fi_dyn_p(1,klon,iip1,jjp1,zdpsrf,pdpsfi) 956 957 c$OMP MASTER 958 kstart=1 959 kend=klon 960 961 if (is_north_pole_dyn) kstart=2 962 if (is_south_pole_dyn) kend=klon-1 963 964 !CDIR ON_ADB(index_i) 965 !CDIR ON_ADB(index_j) 966 !cdir NODEP 967 do ig0=kstart,kend 968 i=index_i(ig0) 969 j=index_j(ig0) 970 pdpsfi(i,j) = zdpsrf(ig0) 971 if (i==1) pdpsfi(iip1,j) = zdpsrf(ig0) 972 enddo 973 974 if (is_north_pole_dyn) then 975 DO i=1,iip1 976 pdpsfi(i,1) = zdpsrf(1) 977 enddo 978 endif 979 980 if (is_south_pole_dyn) then 981 DO i=1,iip1 982 pdpsfi(i,jjp1) = zdpsrf(klon) 983 ENDDO 984 endif 985 c$OMP END MASTER 986 cc$OMP BARRIER 987 988 c 989 c 62. enthalpie potentielle 990 c --------------------- 991 992 kstart=1 993 kend=klon 994 995 if (is_north_pole_dyn) kstart=2 996 if (is_south_pole_dyn) kend=klon-1 997 998 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 999 DO l=1,llm 1000 1001 !CDIR ON_ADB(index_i) 1002 !CDIR ON_ADB(index_j) 1003 !cdir NODEP 1004 do ig0=kstart,kend 1005 i=index_i(ig0) 1006 j=index_j(ig0) 1007 pdhfi(i,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l) 1008 if (i==1) pdhfi(iip1,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l) 1009 enddo 1010 1011 if (is_north_pole_dyn) then 1012 DO i=1,iip1 1013 pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l) 1014 enddo 1015 endif 1016 1017 if (is_south_pole_dyn) then 1018 DO i=1,iip1 1019 pdhfi(i,jjp1,l) = cpp * zdtfi(klon,l)/ ppk(i,jjp1,l) 1020 ENDDO 1021 endif 1022 ENDDO 1023 c$OMP END DO NOWAIT 1024 1025 c 62. humidite specifique 1026 c --------------------- 1027 ! Ehouarn: removed this useless bit: was overwritten at step 63 anyways 1028 ! DO iq=1,nqtot 1029 !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1030 ! DO l=1,llm 1031 !!!cdir NODEP 1032 ! do ig0=kstart,kend 1033 ! i=index_i(ig0) 1034 ! j=index_j(ig0) 1035 ! pdqfi(i,j,l,iq) = zdqfi(ig0,l,iq) 1036 ! if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,iq) 1037 ! enddo 1038 ! 1039 ! if (is_north_pole_dyn) then 1040 ! do i=1,iip1 1041 ! pdqfi(i,1,l,iq) = zdqfi(1,l,iq) 1042 ! enddo 1043 ! endif 1044 ! 1045 ! if (is_south_pole_dyn) then 1046 ! do i=1,iip1 1047 ! pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,iq) 1048 ! enddo 1049 ! endif 1050 ! ENDDO 1051 !c$OMP END DO NOWAIT 1052 ! ENDDO 1053 1054 c 63. traceurs 1055 c ------------ 1056 C initialisation des tendances 1057 1058 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1059 DO l=1,llm 1060 pdqfi(:,jj_begin:jj_end,l,:)=0. 1061 ENDDO 1062 c$OMP END DO NOWAIT 1063 1064 C 1065 !cdir NODEP 1066 itr = 0 1067 DO iq=1,nqtot 1068 IF(.NOT.tracers(iq)%isAdvected) CYCLE 1069 itr = itr + 1 1070 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1071 DO l=1,llm 1072 !CDIR ON_ADB(index_i) 1073 !CDIR ON_ADB(index_j) 1074 !cdir NODEP 1075 DO ig0=kstart,kend 1076 i=index_i(ig0) 1077 j=index_j(ig0) 1078 pdqfi(i,j,l,iq) = zdqfi(ig0,l,itr) 1079 if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,itr) 1080 ENDDO 1081 1082 IF (is_north_pole_dyn) then 1083 DO i=1,iip1 1084 pdqfi(i,1,l,iq) = zdqfi(1,l,itr) 1085 ENDDO 1086 ENDIF 1087 1088 IF (is_south_pole_dyn) then 1089 DO i=1,iip1 1090 pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,itr) 1091 ENDDO 1092 ENDIF 1093 1094 ENDDO 1095 c$OMP END DO NOWAIT 1096 ENDDO 1097 1098 c 65. champ u: 1099 c ------------ 1100 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1101 DO l=1,llm 1102 !CDIR ON_ADB(index_i) 1103 !CDIR ON_ADB(index_j) 1104 !cdir NODEP 1105 do ig0=kstart,kend 1106 i=index_i(ig0) 1107 j=index_j(ig0) 1108 1109 if (i/=iim) then 1110 pdufi(i,j,l)=0.5*(zdufi2(ig0,l)+zdufi2(ig0+1,l))*cu(i,j) 1111 endif 1112 1113 if (i==1) then 1114 pdufi(iim,j,l)=0.5*( zdufi2(ig0,l) 1115 $ + zdufi2(ig0+iim-1,l))*cu(iim,j) 1116 pdufi(iip1,j,l)=0.5*(zdufi2(ig0,l)+zdufi2(ig0+1,l))*cu(i,j) 1117 endif 1118 1119 enddo 1120 1121 if (is_north_pole_dyn) then 1122 DO i=1,iip1 1123 pdufi(i,1,l) = 0. 1124 ENDDO 1125 endif 1126 1127 if (is_south_pole_dyn) then 1128 DO i=1,iip1 1129 pdufi(i,jjp1,l) = 0. 1130 ENDDO 1131 endif 1132 1133 ENDDO 1134 c$OMP END DO NOWAIT 1135 1136 c 67. champ v: 1137 c ------------ 1138 1139 kstart=1 1140 kend=klon 1141 1142 if (is_north_pole_dyn) kstart=2 1143 if (is_south_pole_dyn) kend=klon-1-iim 1144 1145 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1146 DO l=1,llm 1147 !CDIR ON_ADB(index_i) 1148 !CDIR ON_ADB(index_j) 1149 !cdir NODEP 1150 do ig0=kstart,kend 1151 i=index_i(ig0) 1152 j=index_j(ig0) 1153 pdvfi(i,j,l)=0.5*(zdvfi2(ig0,l)+zdvfi2(ig0+iim,l))*cv(i,j) 1154 if (i==1) pdvfi(iip1,j,l) = 0.5*(zdvfi2(ig0,l)+ 1155 $ zdvfi2(ig0+iim,l)) 1156 $ *cv(i,j) 1157 enddo 1158 1159 ENDDO 1160 c$OMP END DO NOWAIT 1161 1162 1163 c 68. champ v pres des poles: 1164 c --------------------------- 1165 c v = U * cos(long) + V * SIN(long) 1166 1167 if (is_north_pole_dyn) then 1168 1169 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1170 DO l=1,llm 1171 1172 DO i=1,iim 1173 pdvfi(i,1,l)= 1174 $ zdufi(1,l)*COS(rlonv(i))+zdvfi(1,l)*SIN(rlonv(i)) 1175 1176 pdvfi(i,1,l)= 1177 $ 0.5*(pdvfi(i,1,l)+zdvfi(i+1,l))*cv(i,1) 1178 ENDDO 1179 1180 pdvfi(iip1,1,l) = pdvfi(1,1,l) 1181 1182 ENDDO 1183 c$OMP END DO NOWAIT 1184 1185 endif 1186 1187 if (is_south_pole_dyn) then 1188 1189 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 1190 DO l=1,llm 1191 1192 DO i=1,iim 1193 pdvfi(i,jjm,l)=zdufi(klon,l)*COS(rlonv(i)) 1194 $ +zdvfi(klon,l)*SIN(rlonv(i)) 1195 1196 pdvfi(i,jjm,l)= 1197 $ 0.5*(pdvfi(i,jjm,l)+zdvfi(klon-iip1+i,l))*cv(i,jjm) 1198 ENDDO 1199 1200 pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l) 1201 1202 ENDDO 1203 c$OMP END DO NOWAIT 1204 1205 endif 1206 c----------------------------------------------------------------------- 1207 1208 700 CONTINUE 1209 1210 firstcal = .FALSE. 1211 1212 #else 1213 call abort_gcm("calfis_loc", 1214 & "calfis_p: for now can only work with parallel physics", 1) 1215 #endif 1216 ! of #ifdef CPP_PHYS 1216 ! of #ifdef CPP_PHYS 1217 1217 #endif 1218 ! of #ifdef CPP_PARA1219 END 1218 ! of #ifdef CPP_PARA 1219 END SUBROUTINE calfis_loc -
LMDZ6/trunk/libf/dynphy_lonlat/gr_dyn_fi.f90
r5245 r5246 2 2 ! $Header$ 3 3 ! 4 5 6 c=======================================================================7 cpassage d'un champ de la grille scalaire a la grille physique8 c=======================================================================4 SUBROUTINE gr_dyn_fi(nfield,im,jm,ngrid,pdyn,pfi) 5 IMPLICIT NONE 6 !======================================================================= 7 ! passage d'un champ de la grille scalaire a la grille physique 8 !======================================================================= 9 9 10 c-----------------------------------------------------------------------11 cdeclarations:12 c-------------10 !----------------------------------------------------------------------- 11 ! declarations: 12 ! ------------- 13 13 14 INTEGERim,jm,ngrid,nfield15 REALpdyn(im,jm,nfield)16 REALpfi(ngrid,nfield)14 INTEGER :: im,jm,ngrid,nfield 15 REAL :: pdyn(im,jm,nfield) 16 REAL :: pfi(ngrid,nfield) 17 17 18 INTEGERj,ifield,ig18 INTEGER :: j,ifield,ig 19 19 20 c-----------------------------------------------------------------------21 ccalcul:22 c-------20 !----------------------------------------------------------------------- 21 ! calcul: 22 ! ------- 23 23 24 25 26 27 ctraitement des poles28 29 24 IF (ngrid.NE.2+(jm-2)*(im-1)) then 25 call abort_gcm("gr_dyn_fi", 'probleme de dim', 1) 26 end if 27 ! traitement des poles 28 CALL SCOPY(nfield,pdyn,im*jm,pfi,ngrid) 29 CALL SCOPY(nfield,pdyn(1,jm,1),im*jm,pfi(ngrid,1),ngrid) 30 30 31 ctraitement des point normaux32 33 34 35 36 37 31 ! traitement des point normaux 32 DO ifield=1,nfield 33 DO j=2,jm-1 34 ig=2+(j-2)*(im-1) 35 CALL SCOPY(im-1,pdyn(1,j,ifield),1,pfi(ig,ifield),1) 36 ENDDO 37 ENDDO 38 38 39 40 END 39 RETURN 40 END SUBROUTINE gr_dyn_fi -
LMDZ6/trunk/libf/dynphy_lonlat/gr_dyn_fi_p.F90
r5245 r5246 2 2 ! $Id$ 3 3 ! 4 4 SUBROUTINE gr_dyn_fi_p(nfield,im,jm,ngrid,pdyn,pfi) 5 5 #ifdef CPP_PARA 6 ! Interface with parallel physics,7 8 9 10 11 c=======================================================================12 cpassage d'un champ de la grille scalaire a la grille physique13 c=======================================================================6 ! Interface with parallel physics, 7 USE mod_interface_dyn_phys 8 USE dimphy 9 USE parallel_lmdz 10 IMPLICIT NONE 11 !======================================================================= 12 ! passage d'un champ de la grille scalaire a la grille physique 13 !======================================================================= 14 14 15 c-----------------------------------------------------------------------16 cdeclarations:17 c-------------15 !----------------------------------------------------------------------- 16 ! declarations: 17 ! ------------- 18 18 19 INTEGERim,jm,ngrid,nfield20 REALpdyn(im,jm,nfield)21 REALpfi(ngrid,nfield)19 INTEGER :: im,jm,ngrid,nfield 20 REAL :: pdyn(im,jm,nfield) 21 REAL :: pfi(ngrid,nfield) 22 22 23 INTEGERi,j,ig,l23 INTEGER :: i,j,ig,l 24 24 25 c-----------------------------------------------------------------------26 ccalcul:27 c-------25 !----------------------------------------------------------------------- 26 ! calcul: 27 ! ------- 28 28 29 cIF(ngrid.NE.2+(jm-2)*(im-1)) STOP 'probleme de dim'30 ctraitement des poles31 ctraitement des point normaux32 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)33 DO l=1,nfield34 35 36 37 38 39 40 c$OMP END DO NOWAIT29 ! IF(ngrid.NE.2+(jm-2)*(im-1)) STOP 'probleme de dim' 30 ! traitement des poles 31 ! traitement des point normaux 32 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 33 DO l=1,nfield 34 DO ig=1,klon 35 i=index_i(ig) 36 j=index_j(ig) 37 pfi(ig,l)=pdyn(i,j,l) 38 ENDDO 39 ENDDO 40 !$OMP END DO NOWAIT 41 41 #endif 42 ! of #ifdef CPP_PARA43 44 END 42 ! of #ifdef CPP_PARA 43 RETURN 44 END SUBROUTINE gr_dyn_fi_p -
LMDZ6/trunk/libf/dynphy_lonlat/gr_fi_dyn.f90
r5245 r5246 2 2 ! $Header$ 3 3 ! 4 5 6 c=======================================================================7 cpassage d'un champ de la grille scalaire a la grille physique8 c=======================================================================4 SUBROUTINE gr_fi_dyn(nfield,ngrid,im,jm,pfi,pdyn) 5 IMPLICIT NONE 6 !======================================================================= 7 ! passage d'un champ de la grille scalaire a la grille physique 8 !======================================================================= 9 9 10 c-----------------------------------------------------------------------11 cdeclarations:12 c-------------10 !----------------------------------------------------------------------- 11 ! declarations: 12 ! ------------- 13 13 14 INTEGERim,jm,ngrid,nfield15 REALpdyn(im,jm,nfield)16 REALpfi(ngrid,nfield)14 INTEGER :: im,jm,ngrid,nfield 15 REAL :: pdyn(im,jm,nfield) 16 REAL :: pfi(ngrid,nfield) 17 17 18 INTEGERi,j,ifield,ig18 INTEGER :: i,j,ifield,ig 19 19 20 c-----------------------------------------------------------------------21 ccalcul:22 c-------20 !----------------------------------------------------------------------- 21 ! calcul: 22 ! ------- 23 23 24 25 ctraitement des poles26 27 28 29 24 DO ifield=1,nfield 25 ! traitement des poles 26 DO i=1,im 27 pdyn(i,1,ifield)=pfi(1,ifield) 28 pdyn(i,jm,ifield)=pfi(ngrid,ifield) 29 ENDDO 30 30 31 ctraitement des point normaux32 33 34 35 36 37 31 ! traitement des point normaux 32 DO j=2,jm-1 33 ig=2+(j-2)*(im-1) 34 CALL SCOPY(im-1,pfi(ig,ifield),1,pdyn(1,j,ifield),1) 35 pdyn(im,j,ifield)=pdyn(1,j,ifield) 36 ENDDO 37 ENDDO 38 38 39 40 END 39 RETURN 40 END SUBROUTINE gr_fi_dyn -
LMDZ6/trunk/libf/dynphy_lonlat/gr_fi_dyn_p.F90
r5245 r5246 2 2 ! $Id$ 3 3 ! 4 4 SUBROUTINE gr_fi_dyn_p(nfield,ngrid,im,jm,pfi,pdyn) 5 5 #ifdef CPP_PARA 6 ! Interface with parallel physics,7 8 9 10 11 c=======================================================================12 cpassage d'un champ de la grille scalaire a la grille physique13 c=======================================================================6 ! Interface with parallel physics, 7 USE mod_interface_dyn_phys 8 USE dimphy 9 USE parallel_lmdz 10 IMPLICIT NONE 11 !======================================================================= 12 ! passage d'un champ de la grille scalaire a la grille physique 13 !======================================================================= 14 14 15 c-----------------------------------------------------------------------16 cdeclarations:17 c-------------15 !----------------------------------------------------------------------- 16 ! declarations: 17 ! ------------- 18 18 19 INTEGERim,jm,ngrid,nfield20 REALpdyn(im,jm,nfield)21 REALpfi(ngrid,nfield)19 INTEGER :: im,jm,ngrid,nfield 20 REAL :: pdyn(im,jm,nfield) 21 REAL :: pfi(ngrid,nfield) 22 22 23 INTEGERi,j,ifield,ig23 INTEGER :: i,j,ifield,ig 24 24 25 c-----------------------------------------------------------------------26 ccalcul:27 c-------28 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)29 25 !----------------------------------------------------------------------- 26 ! calcul: 27 ! ------- 28 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 29 DO ifield=1,nfield 30 30 31 32 33 34 35 36 31 do ig=1,klon 32 i=index_i(ig) 33 j=index_j(ig) 34 pdyn(i,j,ifield)=pfi(ig,ifield) 35 if (i==1) pdyn(im,j,ifield)=pdyn(i,j,ifield) 36 enddo 37 37 38 ctraitement des poles39 40 41 42 43 44 45 46 47 48 49 50 51 52 c$OMP END DO NOWAIT38 ! traitement des poles 39 if (pole_nord) then 40 do i=1,im 41 pdyn(i,1,ifield)=pdyn(1,1,ifield) 42 enddo 43 endif 44 45 if (pole_sud) then 46 do i=1,im 47 pdyn(i,jm,ifield)=pdyn(1,jm,ifield) 48 enddo 49 endif 50 51 ENDDO 52 !$OMP END DO NOWAIT 53 53 #endif 54 ! of #ifdef CPP_PARA55 56 END 54 ! of #ifdef CPP_PARA 55 RETURN 56 END SUBROUTINE gr_fi_dyn_p
Note: See TracChangeset
for help on using the changeset viewer.
