Changeset 776 for trunk/LMDZ.COMMON/libf/dyn3d
- Timestamp:
- Sep 7, 2012, 2:49:58 PM (12 years ago)
- Location:
- trunk/LMDZ.COMMON/libf/dyn3d
- Files:
-
- 17 edited
- 1 moved
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.COMMON/libf/dyn3d/calfis.F
r108 r776 170 170 PARAMETER(ntetaSTD=3) 171 171 REAL rtetaSTD(ntetaSTD) 172 DATA rtetaSTD/350., 380., 405./ 172 DATA rtetaSTD/350., 380., 405./ ! Earth-specific values, beware !! 173 173 REAL PVteta(ngridmx,ntetaSTD) 174 174 c … … 461 461 if (planet_type=="earth") then 462 462 #ifdef CPP_EARTH 463 ! PVtheta calls tetalevel, which is in the (Earth) physics 463 464 cIM calcul PV a teta=350, 380, 405K 464 465 CALL PVtheta(ngridmx,llm,pucov,pvcov,pteta, … … 482 483 ! ne pose pas de probleme a priori. 483 484 484 #ifdef CPP_PHYS485 486 485 ! write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys 487 486 zdt_split=dtphys/nsplit_phys … … 490 489 zdtfic(:,:)=0. 491 490 zdqfic(:,:,:)=0. 491 492 #ifdef CPP_PHYS 492 493 493 494 do isplit=1,nsplit_phys … … 563 564 zdqfic(:,:,:)=zdqfic(:,:,:)+zdqfi(:,:,:) 564 565 565 enddo 566 enddo ! of do isplit=1,nsplit_phys 567 568 #endif 569 ! #endif of #ifdef CPP_PHYS 570 566 571 zdufi(:,:)=zdufic(:,:)/nsplit_phys 567 572 zdvfi(:,:)=zdvfic(:,:)/nsplit_phys … … 569 574 zdqfi(:,:,:)=zdqfic(:,:,:)/nsplit_phys 570 575 571 #endif572 ! #endif of #ifdef CPP_PHYS573 576 574 577 500 CONTINUE -
trunk/LMDZ.COMMON/libf/dyn3d/ce0l.F90
r492 r776 28 28 IMPLICIT NONE 29 29 #ifndef CPP_EARTH 30 WRITE(*,*)'limit_netcdf: Earth-specific program, needs Earth physics' 30 #include "iniprint.h" 31 WRITE(lunout,*)'limit_netcdf: Earth-specific program, needs Earth physics' 31 32 #else 32 33 !------------------------------------------------------------------------------- -
trunk/LMDZ.COMMON/libf/dyn3d/comvert.h
r127 r776 1 1 ! 2 ! $Id: comvert.h 1 520 2011-05-23 11:37:09Z emillour$2 ! $Id: comvert.h 1625 2012-05-09 13:14:48Z lguez $ 3 3 ! 4 4 !----------------------------------------------------------------------- … … 9 9 & aps(llm),bps(llm),scaleheight 10 10 11 common/comverti/disvert_type 11 common/comverti/disvert_type, pressure_exner 12 12 13 13 real ap ! hybrid pressure contribution at interlayers … … 30 30 ! using 'z2sig.def' (or 'esasig.def) file 31 31 32 logical pressure_exner 33 ! compute pressure inside layers using Exner function, else use mean 34 ! of pressure values at interfaces 35 32 36 !----------------------------------------------------------------------- -
trunk/LMDZ.COMMON/libf/dyn3d/disvert.F90
r128 r776 1 ! $Id: disvert.F90 1 520 2011-05-23 11:37:09Z emillour$1 ! $Id: disvert.F90 1645 2012-07-30 16:01:50Z lguez $ 2 2 3 3 SUBROUTINE disvert(pa,preff,ap,bp,dpres,presnivs,nivsigs,nivsig,scaleheight) 4 4 5 5 ! Auteur : P. Le Van 6 7 use new_unit_m, only: new_unit 8 use ioipsl, only: getin 9 use assert_m, only: assert 6 10 7 11 IMPLICIT NONE … … 18 22 19 23 real,intent(in) :: pa, preff 20 real,intent(out) :: ap(llmp1), bp(llmp1) 24 real,intent(out) :: ap(llmp1) ! in Pa 25 real, intent(out):: bp(llmp1) 21 26 real,intent(out) :: dpres(llm), nivsigs(llm), nivsig(llmp1) 22 27 real,intent(out) :: presnivs(llm) … … 26 31 real zk, zkm1, dzk1, dzk2, k0, k1 27 32 28 INTEGER l 33 INTEGER l, unit 29 34 REAL dsigmin 30 35 REAL alpha, beta, deltaz 31 INTEGER iostat32 36 REAL x 33 37 character(len=*),parameter :: modname="disvert" 34 38 39 character(len=6):: vert_sampling 40 ! (allowed values are "param", "tropo", "strato" and "read") 41 35 42 !----------------------------------------------------------------------- 43 44 print *, "Call sequence information: disvert" 36 45 37 46 ! default scaleheight is 8km for earth 38 47 scaleheight=8. 39 48 40 OPEN(99, file='sigma.def', status='old', form='formatted', iostat=iostat) 49 vert_sampling = merge("strato", "tropo ", ok_strato) ! default value 50 call getin('vert_sampling', vert_sampling) 51 print *, 'vert_sampling = ' // vert_sampling 41 52 42 IF (iostat == 0) THEN 43 ! cas 1 on lit les options dans sigma.def: 53 select case (vert_sampling) 54 case ("param") 55 ! On lit les options dans sigma.def: 56 OPEN(99, file='sigma.def', status='old', form='formatted') 44 57 READ(99, *) scaleheight ! hauteur d'echelle 8. 45 58 READ(99, *) deltaz ! epaiseur de la premiere couche 0.04 … … 69 82 sig(llm+1)=0. 70 83 71 DO l = 1, llm 72 dsig(l) = sig(l)-sig(l+1) 73 end DO 74 ELSE 75 if (ok_strato) then 76 if (llm==39) then 77 dsigmin=0.3 78 else if (llm==50) then 79 dsigmin=1. 80 else 81 write(lunout,*) trim(modname), & 82 ' ATTENTION discretisation z a ajuster' 83 dsigmin=1. 84 endif 85 write(lunout,*) trim(modname), & 86 ' Discretisation verticale DSIGMIN=',dsigmin 87 endif 84 bp(: llm) = EXP(1. - 1. / sig(: llm)**2) 85 bp(llmp1) = 0. 88 86 87 ap = pa * (sig - bp) 88 case("tropo") 89 89 DO l = 1, llm 90 90 x = 2*asin(1.) * (l - 0.5) / (llm + 1) 91 92 IF (ok_strato) THEN 93 dsig(l) =(dsigmin + 7. * SIN(x)**2) & 94 *(0.5*(1.-tanh(1.*(x-asin(1.))/asin(1.))))**2 95 ELSE 96 dsig(l) = 1.0 + 7.0 * SIN(x)**2 97 ENDIF 91 dsig(l) = 1.0 + 7.0 * SIN(x)**2 98 92 ENDDO 99 93 dsig = dsig / sum(dsig) … … 102 96 sig(l) = sig(l+1) + dsig(l) 103 97 ENDDO 104 ENDIF 98 99 bp(1)=1. 100 bp(2: llm) = EXP(1. - 1. / sig(2: llm)**2) 101 bp(llmp1) = 0. 102 103 ap(1)=0. 104 ap(2: llm + 1) = pa * (sig(2: llm + 1) - bp(2: llm + 1)) 105 case("strato") 106 if (llm==39) then 107 dsigmin=0.3 108 else if (llm==50) then 109 dsigmin=1. 110 else 111 write(lunout,*) trim(modname), ' ATTENTION discretisation z a ajuster' 112 dsigmin=1. 113 endif 114 WRITE(LUNOUT,*) trim(modname), 'Discretisation verticale DSIGMIN=',dsigmin 115 116 DO l = 1, llm 117 x = 2*asin(1.) * (l - 0.5) / (llm + 1) 118 dsig(l) =(dsigmin + 7. * SIN(x)**2) & 119 *(0.5*(1.-tanh(1.*(x-asin(1.))/asin(1.))))**2 120 ENDDO 121 dsig = dsig / sum(dsig) 122 sig(llm+1) = 0. 123 DO l = llm, 1, -1 124 sig(l) = sig(l+1) + dsig(l) 125 ENDDO 126 127 bp(1)=1. 128 bp(2: llm) = EXP(1. - 1. / sig(2: llm)**2) 129 bp(llmp1) = 0. 130 131 ap(1)=0. 132 ap(2: llm + 1) = pa * (sig(2: llm + 1) - bp(2: llm + 1)) 133 case("read") 134 ! Read "ap" and "bp". First line is skipped (title line). "ap" 135 ! should be in Pa. First couple of values should correspond to 136 ! the surface, that is : "bp" should be in descending order. 137 call new_unit(unit) 138 open(unit, file="hybrid.txt", status="old", action="read", & 139 position="rewind") 140 read(unit, fmt=*) ! skip title line 141 do l = 1, llm + 1 142 read(unit, fmt=*) ap(l), bp(l) 143 end do 144 close(unit) 145 call assert(ap(1) == 0., ap(llm + 1) == 0., bp(1) == 1., & 146 bp(llm + 1) == 0., "disvert: bad ap or bp values") 147 case default 148 call abort_gcm("disvert", 'Wrong value for "vert_sampling"', 1) 149 END select 105 150 106 151 DO l=1, llm … … 111 156 nivsig(l)= REAL(l) 112 157 ENDDO 113 114 ! .... Calculs de ap(l) et de bp(l) ....115 ! ..... pa et preff sont lus sur les fichiers start par lectba .....116 117 bp(llmp1) = 0.118 119 DO l = 1, llm120 bp(l) = EXP( 1. -1./( sig(l)*sig(l)) )121 ap(l) = pa * ( sig(l) - bp(l) )122 ENDDO123 124 bp(1)=1.125 ap(1)=0.126 127 ap(llmp1) = pa * ( sig(llmp1) - bp(llmp1) )128 158 129 159 write(lunout, *) trim(modname),': BP ' -
trunk/LMDZ.COMMON/libf/dyn3d/dynetat0.F
r492 r776 6 6 7 7 USE infotrac 8 use netcdf, only: nf90_get_var 8 9 IMPLICIT NONE 9 10 … … 28 29 #include "comconst.h" 29 30 #include "comvert.h" 30 #include "comgeom .h"31 #include "comgeom2.h" 31 32 #include "ener.h" 32 33 #include "netcdf.inc" … … 40 41 41 42 CHARACTER*(*) fichnom 42 REAL vcov(i p1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm)43 REAL q(i p1jmp1,llm,nqtot),masse(ip1jmp1,llm)44 REAL ps(i p1jmp1),phis(ip1jmp1)43 REAL vcov(iip1, jjm,llm),ucov(iip1, jjp1,llm),teta(iip1, jjp1,llm) 44 REAL q(iip1,jjp1,llm,nqtot),masse(iip1, jjp1,llm) 45 REAL ps(iip1, jjp1),phis(iip1, jjp1) 45 46 46 47 REAL time … … 70 71 CALL abort 71 72 ENDIF 72 #ifdef NC_DOUBLE 73 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, tab_cntrl) 74 #else 75 ierr = NF_GET_VAR_REAL(nid, nvarid, tab_cntrl) 76 #endif 73 ierr = nf90_get_var(nid, nvarid, tab_cntrl) 77 74 IF (ierr .NE. NF_NOERR) THEN 78 75 write(lunout,*)"dynetat0: Lecture echoue pour <controle>" … … 142 139 CALL abort 143 140 ENDIF 144 #ifdef NC_DOUBLE 145 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlonu) 146 #else 147 ierr = NF_GET_VAR_REAL(nid, nvarid, rlonu) 148 #endif 141 ierr = nf90_get_var(nid, nvarid, rlonu) 149 142 IF (ierr .NE. NF_NOERR) THEN 150 143 write(lunout,*)"dynetat0: Lecture echouee pour <rlonu>" … … 157 150 CALL abort 158 151 ENDIF 159 #ifdef NC_DOUBLE 160 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlatu) 161 #else 162 ierr = NF_GET_VAR_REAL(nid, nvarid, rlatu) 163 #endif 152 ierr = nf90_get_var(nid, nvarid, rlatu) 164 153 IF (ierr .NE. NF_NOERR) THEN 165 154 write(lunout,*)"dynetat0: Lecture echouee pour <rlatu>" … … 172 161 CALL abort 173 162 ENDIF 174 #ifdef NC_DOUBLE 175 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlonv) 176 #else 177 ierr = NF_GET_VAR_REAL(nid, nvarid, rlonv) 178 #endif 163 ierr = nf90_get_var(nid, nvarid, rlonv) 179 164 IF (ierr .NE. NF_NOERR) THEN 180 165 write(lunout,*)"dynetat0: Lecture echouee pour <rlonv>" … … 187 172 CALL abort 188 173 ENDIF 189 #ifdef NC_DOUBLE 190 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlatv) 191 #else 192 ierr = NF_GET_VAR_REAL(nid, nvarid, rlatv) 193 #endif 174 ierr = nf90_get_var(nid, nvarid, rlatv) 194 175 IF (ierr .NE. NF_NOERR) THEN 195 176 write(lunout,*)"dynetat0: Lecture echouee pour rlatv" … … 202 183 CALL abort 203 184 ENDIF 204 #ifdef NC_DOUBLE 205 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, cu) 206 #else 207 ierr = NF_GET_VAR_REAL(nid, nvarid, cu) 208 #endif 185 ierr = nf90_get_var(nid, nvarid, cu) 209 186 IF (ierr .NE. NF_NOERR) THEN 210 187 write(lunout,*)"dynetat0: Lecture echouee pour <cu>" … … 217 194 CALL abort 218 195 ENDIF 219 #ifdef NC_DOUBLE 220 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, cv) 221 #else 222 ierr = NF_GET_VAR_REAL(nid, nvarid, cv) 223 #endif 196 ierr = nf90_get_var(nid, nvarid, cv) 224 197 IF (ierr .NE. NF_NOERR) THEN 225 198 write(lunout,*)"dynetat0: Lecture echouee pour <cv>" … … 232 205 CALL abort 233 206 ENDIF 234 #ifdef NC_DOUBLE 235 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, aire) 236 #else 237 ierr = NF_GET_VAR_REAL(nid, nvarid, aire) 238 #endif 207 ierr = nf90_get_var(nid, nvarid, aire) 239 208 IF (ierr .NE. NF_NOERR) THEN 240 209 write(lunout,*)"dynetat0: Lecture echouee pour <aire>" … … 247 216 CALL abort 248 217 ENDIF 249 #ifdef NC_DOUBLE 250 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, phis) 251 #else 252 ierr = NF_GET_VAR_REAL(nid, nvarid, phis) 253 #endif 218 ierr = nf90_get_var(nid, nvarid, phis) 254 219 IF (ierr .NE. NF_NOERR) THEN 255 220 write(lunout,*)"dynetat0: Lecture echouee pour <phisinit>" … … 262 227 CALL abort 263 228 ENDIF 264 #ifdef NC_DOUBLE 265 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, time) 266 #else 267 ierr = NF_GET_VAR_REAL(nid, nvarid, time) 268 #endif 229 ierr = nf90_get_var(nid, nvarid, time) 269 230 IF (ierr .NE. NF_NOERR) THEN 270 231 write(lunout,*)"dynetat0: Lecture echouee <temps>" … … 277 238 CALL abort 278 239 ENDIF 279 #ifdef NC_DOUBLE 280 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, ucov) 281 #else 282 ierr = NF_GET_VAR_REAL(nid, nvarid, ucov) 283 #endif 240 ierr = nf90_get_var(nid, nvarid, ucov) 284 241 IF (ierr .NE. NF_NOERR) THEN 285 242 write(lunout,*)"dynetat0: Lecture echouee pour <ucov>" … … 292 249 CALL abort 293 250 ENDIF 294 #ifdef NC_DOUBLE 295 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, vcov) 296 #else 297 ierr = NF_GET_VAR_REAL(nid, nvarid, vcov) 298 #endif 251 ierr = nf90_get_var(nid, nvarid, vcov) 299 252 IF (ierr .NE. NF_NOERR) THEN 300 253 write(lunout,*)"dynetat0: Lecture echouee pour <vcov>" … … 307 260 CALL abort 308 261 ENDIF 309 #ifdef NC_DOUBLE 310 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, teta) 311 #else 312 ierr = NF_GET_VAR_REAL(nid, nvarid, teta) 313 #endif 262 ierr = nf90_get_var(nid, nvarid, teta) 314 263 IF (ierr .NE. NF_NOERR) THEN 315 264 write(lunout,*)"dynetat0: Lecture echouee pour <teta>" … … 325 274 & "> est absent" 326 275 write(lunout,*)" Il est donc initialise a zero" 327 q(:,:, iq)=0.276 q(:,:,:,iq)=0. 328 277 ELSE 329 #ifdef NC_DOUBLE 330 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, q(1,1,iq)) 331 #else 332 ierr = NF_GET_VAR_REAL(nid, nvarid, q(1,1,iq)) 333 #endif 278 ierr = NF90_GET_VAR(nid, nvarid, q(:,:,:,iq)) 334 279 IF (ierr .NE. NF_NOERR) THEN 335 280 write(lunout,*)"dynetat0: Lecture echouee pour "//tname(iq) … … 345 290 CALL abort 346 291 ENDIF 347 #ifdef NC_DOUBLE 348 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, masse) 349 #else 350 ierr = NF_GET_VAR_REAL(nid, nvarid, masse) 351 #endif 292 ierr = nf90_get_var(nid, nvarid, masse) 352 293 IF (ierr .NE. NF_NOERR) THEN 353 294 write(lunout,*)"dynetat0: Lecture echouee pour <masse>" … … 360 301 CALL abort 361 302 ENDIF 362 #ifdef NC_DOUBLE 363 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, ps) 364 #else 365 ierr = NF_GET_VAR_REAL(nid, nvarid, ps) 366 #endif 303 ierr = nf90_get_var(nid, nvarid, ps) 367 304 IF (ierr .NE. NF_NOERR) THEN 368 305 write(lunout,*)"dynetat0: Lecture echouee pour <ps>" -
trunk/LMDZ.COMMON/libf/dyn3d/dynredem.F
r492 r776 1 1 ! 2 ! $Id: dynredem.F 1 403 2010-07-01 09:02:53Z fairhead$2 ! $Id: dynredem.F 1635 2012-07-12 11:37:16Z lguez $ 3 3 ! 4 4 c … … 8 8 #endif 9 9 USE infotrac 10 use netcdf95, only: NF95_PUT_VAR 10 11 11 12 IMPLICIT NONE … … 19 20 #include "comconst.h" 20 21 #include "comvert.h" 21 #include "comgeom .h"22 #include "comgeom2.h" 22 23 #include "temps.h" 23 24 #include "ener.h" … … 31 32 c ---------- 32 33 INTEGER iday_end 33 REAL phis(i p1jmp1)34 REAL phis(iip1, jjp1) 34 35 CHARACTER*(*) fichnom 35 36 … … 138 139 c 139 140 ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 27, 140 . "Fichier dem arrage dynamique")141 . "Fichier demmarage dynamique") 141 142 c 142 143 c Definir les dimensions du fichiers: … … 166 167 . "Parametres de controle") 167 168 ierr = NF_ENDDEF(nid) 168 #ifdef NC_DOUBLE 169 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,tab_cntrl) 170 #else 171 ierr = NF_PUT_VAR_REAL (nid,nvarid,tab_cntrl) 172 #endif 169 call NF95_PUT_VAR(nid,nvarid,tab_cntrl) 173 170 c 174 171 ierr = NF_REDEF (nid) … … 183 180 . "Longitudes des points U") 184 181 ierr = NF_ENDDEF(nid) 185 #ifdef NC_DOUBLE 186 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlonu) 187 #else 188 ierr = NF_PUT_VAR_REAL (nid,nvarid,rlonu) 189 #endif 182 call NF95_PUT_VAR(nid,nvarid,rlonu) 190 183 c 191 184 ierr = NF_REDEF (nid) … … 200 193 . "Latitudes des points U") 201 194 ierr = NF_ENDDEF(nid) 202 #ifdef NC_DOUBLE 203 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlatu) 204 #else 205 ierr = NF_PUT_VAR_REAL (nid,nvarid,rlatu) 206 #endif 195 call NF95_PUT_VAR (nid,nvarid,rlatu) 207 196 c 208 197 ierr = NF_REDEF (nid) … … 217 206 . "Longitudes des points V") 218 207 ierr = NF_ENDDEF(nid) 219 #ifdef NC_DOUBLE 220 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlonv) 221 #else 222 ierr = NF_PUT_VAR_REAL (nid,nvarid,rlonv) 223 #endif 208 call NF95_PUT_VAR(nid,nvarid,rlonv) 224 209 c 225 210 ierr = NF_REDEF (nid) … … 234 219 . "Latitudes des points V") 235 220 ierr = NF_ENDDEF(nid) 236 #ifdef NC_DOUBLE 237 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlatv) 238 #else 239 ierr = NF_PUT_VAR_REAL (nid,nvarid,rlatv) 240 #endif 221 call NF95_PUT_VAR(nid,nvarid,rlatv) 241 222 c 242 223 ierr = NF_REDEF (nid) … … 251 232 . "Numero naturel des couches s") 252 233 ierr = NF_ENDDEF(nid) 253 #ifdef NC_DOUBLE 254 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,nivsigs) 255 #else 256 ierr = NF_PUT_VAR_REAL (nid,nvarid,nivsigs) 257 #endif 234 call NF95_PUT_VAR(nid,nvarid,nivsigs) 258 235 c 259 236 ierr = NF_REDEF (nid) … … 268 245 . "Numero naturel des couches sigma") 269 246 ierr = NF_ENDDEF(nid) 270 #ifdef NC_DOUBLE 271 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,nivsig) 272 #else 273 ierr = NF_PUT_VAR_REAL (nid,nvarid,nivsig) 274 #endif 247 call NF95_PUT_VAR(nid,nvarid,nivsig) 275 248 c 276 249 ierr = NF_REDEF (nid) … … 285 258 . "Coefficient A pour hybride") 286 259 ierr = NF_ENDDEF(nid) 287 #ifdef NC_DOUBLE 288 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,ap) 289 #else 290 ierr = NF_PUT_VAR_REAL (nid,nvarid,ap) 291 #endif 260 call NF95_PUT_VAR(nid,nvarid,ap) 292 261 c 293 262 ierr = NF_REDEF (nid) … … 302 271 . "Coefficient B pour hybride") 303 272 ierr = NF_ENDDEF(nid) 304 #ifdef NC_DOUBLE 305 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,bp) 306 #else 307 ierr = NF_PUT_VAR_REAL (nid,nvarid,bp) 308 #endif 273 call NF95_PUT_VAR(nid,nvarid,bp) 309 274 c 310 275 ierr = NF_REDEF (nid) … … 317 282 cIM 220306 END 318 283 ierr = NF_ENDDEF(nid) 319 #ifdef NC_DOUBLE 320 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,presnivs) 321 #else 322 ierr = NF_PUT_VAR_REAL (nid,nvarid,presnivs) 323 #endif 284 call NF95_PUT_VAR(nid,nvarid,presnivs) 324 285 c 325 286 c Coefficients de passage cov. <-> contra. <--> naturel … … 338 299 . "Coefficient de passage pour U") 339 300 ierr = NF_ENDDEF(nid) 340 #ifdef NC_DOUBLE 341 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,cu) 342 #else 343 ierr = NF_PUT_VAR_REAL (nid,nvarid,cu) 344 #endif 301 call NF95_PUT_VAR(nid,nvarid,cu) 345 302 c 346 303 ierr = NF_REDEF (nid) … … 357 314 . "Coefficient de passage pour V") 358 315 ierr = NF_ENDDEF(nid) 359 #ifdef NC_DOUBLE 360 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,cv) 361 #else 362 ierr = NF_PUT_VAR_REAL (nid,nvarid,cv) 363 #endif 316 call NF95_PUT_VAR(nid,nvarid,cv) 364 317 c 365 318 c Aire de chaque maille: … … 378 331 . "Aires de chaque maille") 379 332 ierr = NF_ENDDEF(nid) 380 #ifdef NC_DOUBLE 381 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,aire) 382 #else 383 ierr = NF_PUT_VAR_REAL (nid,nvarid,aire) 384 #endif 333 call NF95_PUT_VAR(nid,nvarid,aire) 385 334 c 386 335 c Geopentiel au sol: … … 399 348 . "Geopotentiel au sol") 400 349 ierr = NF_ENDDEF(nid) 401 #ifdef NC_DOUBLE 402 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,phis) 403 #else 404 ierr = NF_PUT_VAR_REAL (nid,nvarid,phis) 405 #endif 350 call NF95_PUT_VAR(nid,nvarid,phis) 406 351 c 407 352 c Definir les variables pour pouvoir les enregistrer plus tard: … … 524 469 USE infotrac 525 470 USE control_mod 471 use netcdf, only: NF90_get_VAR 472 use netcdf95, only: NF95_PUT_VAR 526 473 527 474 IMPLICIT NONE … … 538 485 #include "iniprint.h" 539 486 487 540 488 INTEGER l 541 REAL vcov(i p1jm,llm),ucov(ip1jmp1,llm)542 REAL teta(i p1jmp1,llm)543 REAL ps(i p1jmp1),masse(ip1jmp1,llm)544 REAL q(i p1jmp1,llm,nqtot)489 REAL vcov(iip1,jjm,llm),ucov(iip1, jjp1,llm) 490 REAL teta(iip1, jjp1,llm) 491 REAL ps(iip1, jjp1),masse(iip1, jjp1,llm) 492 REAL q(iip1, jjp1, llm, nqtot) 545 493 CHARACTER*(*) fichnom 546 494 … … 576 524 CALL abort_gcm(modname,abort_message,ierr) 577 525 ENDIF 578 #ifdef NC_DOUBLE 579 ierr = NF_PUT_VAR1_DOUBLE (nid,nvarid,nb,time) 580 #else 581 ierr = NF_PUT_VAR1_REAL (nid,nvarid,nb,time) 582 #endif 526 call NF95_PUT_VAR(nid,nvarid,time,start=(/nb/)) 583 527 write(lunout,*) "dynredem1: Enregistrement pour ", nb, time 584 528 … … 592 536 CALL abort_gcm(modname,abort_message,ierr) 593 537 ENDIF 594 #ifdef NC_DOUBLE 595 ierr = NF_GET_VAR_DOUBLE(nid, nvarid, tab_cntrl) 596 #else 597 ierr = NF_GET_VAR_REAL(nid, nvarid, tab_cntrl) 598 #endif 538 ierr = NF90_GET_VAR(nid, nvarid, tab_cntrl) 599 539 tab_cntrl(31) = REAL(itau_dyn + itaufin) 600 #ifdef NC_DOUBLE 601 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,tab_cntrl) 602 #else 603 ierr = NF_PUT_VAR_REAL (nid,nvarid,tab_cntrl) 604 #endif 540 call NF95_PUT_VAR(nid,nvarid,tab_cntrl) 605 541 606 542 c Ecriture des champs … … 612 548 CALL abort_gcm(modname,abort_message,ierr) 613 549 ENDIF 614 #ifdef NC_DOUBLE 615 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,ucov) 616 #else 617 ierr = NF_PUT_VAR_REAL (nid,nvarid,ucov) 618 #endif 550 call NF95_PUT_VAR(nid,nvarid,ucov) 619 551 620 552 ierr = NF_INQ_VARID(nid, "vcov", nvarid) … … 624 556 CALL abort_gcm(modname,abort_message,ierr) 625 557 ENDIF 626 #ifdef NC_DOUBLE 627 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,vcov) 628 #else 629 ierr = NF_PUT_VAR_REAL (nid,nvarid,vcov) 630 #endif 558 call NF95_PUT_VAR(nid,nvarid,vcov) 631 559 632 560 ierr = NF_INQ_VARID(nid, "teta", nvarid) … … 636 564 CALL abort_gcm(modname,abort_message,ierr) 637 565 ENDIF 638 #ifdef NC_DOUBLE 639 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,teta) 640 #else 641 ierr = NF_PUT_VAR_REAL (nid,nvarid,teta) 642 #endif 566 call NF95_PUT_VAR(nid,nvarid,teta) 643 567 644 568 IF (type_trac == 'inca') THEN … … 662 586 CALL abort_gcm(modname,abort_message,ierr) 663 587 ENDIF 664 #ifdef NC_DOUBLE 665 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,q(1,1,iq)) 666 #else 667 ierr = NF_PUT_VAR_REAL (nid,nvarid,q(1,1,iq)) 668 #endif 588 call NF95_PUT_VAR(nid,nvarid,q(:,:,:,iq)) 669 589 ELSE ! type_trac = inca 670 590 ! lecture de la valeur du traceur dans start_trac.nc … … 681 601 CALL abort_gcm(modname,abort_message,ierr) 682 602 ENDIF 683 #ifdef NC_DOUBLE 684 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,q(1,1,iq)) 685 #else 686 ierr = NF_PUT_VAR_REAL (nid,nvarid,q(1,1,iq)) 687 #endif 603 call NF95_PUT_VAR(nid,nvarid,q(:,:,:,iq)) 688 604 689 605 ELSE 690 606 write(lunout,*) "dynredem1: ",trim(tname(iq)), 691 607 & " est present dans start_trac.nc" 692 #ifdef NC_DOUBLE 693 ierr = NF_GET_VAR_DOUBLE(nid_trac, nvarid_trac, trac_tmp) 694 #else 695 ierr = NF_GET_VAR_REAL(nid_trac, nvarid_trac, trac_tmp) 696 #endif 608 ierr = NF90_GET_VAR(nid_trac, nvarid_trac, trac_tmp) 697 609 IF (ierr .NE. NF_NOERR) THEN 698 610 abort_message="dynredem1: Lecture echouee pour"// … … 708 620 CALL abort_gcm(modname,abort_message,ierr) 709 621 ENDIF 710 #ifdef NC_DOUBLE 711 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,trac_tmp) 712 #else 713 ierr = NF_PUT_VAR_REAL (nid,nvarid,trac_tmp) 714 #endif 622 call NF95_PUT_VAR(nid, nvarid, trac_tmp) 715 623 716 624 ENDIF ! IF (ierr .NE. NF_NOERR) … … 725 633 CALL abort_gcm(modname,abort_message,ierr) 726 634 ENDIF 727 #ifdef NC_DOUBLE 728 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,q(1,1,iq)) 729 #else 730 ierr = NF_PUT_VAR_REAL (nid,nvarid,q(1,1,iq)) 731 #endif 635 call NF95_PUT_VAR(nid,nvarid,q(:,:,:,iq)) 732 636 ENDIF ! (ierr_file .ne. 2) 733 637 END IF !type_trac … … 742 646 CALL abort_gcm(modname,abort_message,ierr) 743 647 ENDIF 744 #ifdef NC_DOUBLE 745 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,masse) 746 #else 747 ierr = NF_PUT_VAR_REAL (nid,nvarid,masse) 748 #endif 648 call NF95_PUT_VAR(nid,nvarid,masse) 749 649 c 750 650 ierr = NF_INQ_VARID(nid, "ps", nvarid) … … 754 654 CALL abort_gcm(modname,abort_message,ierr) 755 655 ENDIF 756 #ifdef NC_DOUBLE 757 ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,ps) 758 #else 759 ierr = NF_PUT_VAR_REAL (nid,nvarid,ps) 760 #endif 656 call NF95_PUT_VAR(nid,nvarid,ps) 761 657 762 658 ierr = NF_CLOSE(nid) -
trunk/LMDZ.COMMON/libf/dyn3d/etat0_netcdf.F90
r127 r776 1 1 ! 2 ! $Id: etat0_netcdf.F90 1 520 2011-05-23 11:37:09Z emillour$2 ! $Id: etat0_netcdf.F90 1625 2012-05-09 13:14:48Z lguez $ 3 3 ! 4 4 !------------------------------------------------------------------------------- … … 251 251 !******************************************************************************* 252 252 CALL pression(ip1jmp1, ap, bp, psol, p3d) 253 if ( disvert_type.eq.1) then253 if (pressure_exner) then 254 254 CALL exner_hyb(ip1jmp1, psol, p3d, alpha, beta, pks, pk, y) 255 else ! we assume that we are in the disvert_type==2 case255 else 256 256 CALL exner_milieu(ip1jmp1,psol,p3d,beta,pks,pk,y) 257 257 endif -
trunk/LMDZ.COMMON/libf/dyn3d/exner_hyb.F
r127 r776 56 56 ! Sanity check 57 57 if (firstcall) then 58 ! check that vertical discretization is compatible59 ! with this routine60 if (disvert_type.ne.1) then61 call abort_gcm(modname,62 & "this routine should only be called if disvert_type==1",42)63 endif64 65 58 ! sanity checks for Shallow Water case (1 vertical layer) 66 59 if (llm.eq.1) then -
trunk/LMDZ.COMMON/libf/dyn3d/exner_milieu.F
r127 r776 53 53 ! Sanity check 54 54 if (firstcall) then 55 ! check that vertical discretization is compatible56 ! with this routine57 if (disvert_type.ne.2) then58 call abort_gcm(modname,59 & "this routine should only be called if disvert_type==2",42)60 endif61 62 55 ! sanity checks for Shallow Water case (1 vertical layer) 63 56 if (llm.eq.1) then -
trunk/LMDZ.COMMON/libf/dyn3d/gcm.F
r492 r776 21 21 ! A nettoyer. On ne veut qu'une ou deux routines d'interface 22 22 ! dynamique -> physique pour l'initialisation 23 ! Ehouarn: for now these only apply to Earth:23 ! Ehouarn: the following are needed with (parallel) physics: 24 24 #ifdef CPP_PHYS 25 25 USE dimphy 26 26 USE comgeomphy 27 #endif28 #ifdef CPP_EARTH29 27 USE mod_phys_lmdz_para, ONLY : klon_mpi_para_nb 30 28 #endif … … 177 175 ! A nettoyer. On ne veut qu'une ou deux routines d'interface 178 176 ! dynamique -> physique pour l'initialisation 179 ! Ehouarn : temporarily (?) keep this only for Earth180 ! if (planet_type.eq."earth") then181 !#ifdef CPP_EARTH182 177 #ifdef CPP_PHYS 183 178 CALL init_phys_lmdz(iim,jjp1,llm,1,(/(jjm-1)*iim+2/)) 184 179 call initcomgeomphy 185 180 #endif 186 ! endif187 181 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 188 182 c … … 465 459 CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys , 466 460 , latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp ) 467 #endif ! CPP_PHYS461 #endif 468 462 call_iniphys=.false. 469 463 ENDIF ! of IF (call_iniphys.and.(iflag_phys.eq.1)) 470 !#endif471 464 472 465 c numero de stockage pour les fichiers de redemarrage: -
trunk/LMDZ.COMMON/libf/dyn3d/guide_mod.F90
r127 r776 644 644 ! ----------------------------------------------------------------- 645 645 CALL pression( ip1jmp1, ap, bp, psi, p ) 646 if ( disvert_type==1) then646 if (pressure_exner) then 647 647 CALL exner_hyb(ip1jmp1,psi,p,alpha,beta,pks,pk,pkf) 648 else ! we assume that we are in the disvert_type==2 case648 else 649 649 CALL exner_milieu(ip1jmp1,psi,p,beta,pks,pk,pkf) 650 650 endif -
trunk/LMDZ.COMMON/libf/dyn3d/iniacademic.F90
r492 r776 1 1 ! 2 ! $Id: iniacademic.F90 1 529 2011-05-26 15:17:33Z fairhead$2 ! $Id: iniacademic.F90 1625 2012-05-09 13:14:48Z lguez $ 3 3 ! 4 4 SUBROUTINE iniacademic(vcov,ucov,teta,q,masse,ps,phis,time_0) … … 222 222 223 223 CALL pression ( ip1jmp1, ap, bp, ps, p ) 224 if ( disvert_type.eq.1) then224 if (pressure_exner) then 225 225 CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf ) 226 else if (disvert_type.eq.2) then226 else 227 227 call exner_milieu(ip1jmp1,ps,p,beta,pks,pk,pkf) 228 else229 write(abort_message,*) "Wrong value for disvert_type: ", &230 disvert_type231 call abort_gcm(modname,abort_message,0)232 228 endif 233 229 CALL massdair(p,masse) -
trunk/LMDZ.COMMON/libf/dyn3d/iniconst.F90
r775 r776 1 1 ! 2 ! $Id: iniconst.F 1520 2011-05-23 11:37:09Z emillour$2 ! $Id: iniconst.F90 1625 2012-05-09 13:14:48Z lguez $ 3 3 ! 4 4 SUBROUTINE iniconst 5 5 6 6 USE control_mod 7 7 #ifdef CPP_IOIPSL 8 8 use IOIPSL 9 9 #else 10 ! if not using IOIPSL, we still need to use (a local version of) getin11 10 ! if not using IOIPSL, we still need to use (a local version of) getin 11 use ioipsl_getincom 12 12 #endif 13 13 14 IMPLICIT NONE 15 c 16 c P. Le Van 17 c 18 c----------------------------------------------------------------------- 19 c Declarations: 20 c ------------- 21 c 22 #include "dimensions.h" 23 #include "paramet.h" 24 #include "comconst.h" 25 #include "temps.h" 26 #include "comvert.h" 27 #include "iniprint.h" 14 IMPLICIT NONE 15 ! 16 ! P. Le Van 17 ! 18 ! Declarations: 19 ! ------------- 20 ! 21 include "dimensions.h" 22 include "paramet.h" 23 include "comconst.h" 24 include "temps.h" 25 include "comvert.h" 26 include "iniprint.h" 28 27 28 character(len=*),parameter :: modname="iniconst" 29 character(len=80) :: abort_message 30 ! 31 ! 32 ! 33 !----------------------------------------------------------------------- 34 ! dimension des boucles: 35 ! ---------------------- 29 36 30 character(len=*),parameter :: modname="iniconst" 31 character(len=80) :: abort_message 32 c 33 c 34 c 35 c----------------------------------------------------------------------- 36 c dimension des boucles: 37 c ---------------------- 37 im = iim 38 jm = jjm 39 lllm = llm 40 imp1 = iim 41 jmp1 = jjm + 1 42 lllmm1 = llm - 1 43 lllmp1 = llm + 1 38 44 39 im = iim 40 jm = jjm 41 lllm = llm 42 imp1 = iim 43 jmp1 = jjm + 1 44 lllmm1 = llm - 1 45 lllmp1 = llm + 1 45 !----------------------------------------------------------------------- 46 46 47 c----------------------------------------------------------------------- 47 dtphys = iphysiq * dtvr 48 unsim = 1./iim 49 pi = 2.*ASIN( 1. ) 48 50 49 dtphys = iphysiq * dtvr 50 unsim = 1./iim 51 pi = 2.*ASIN( 1. ) 51 !----------------------------------------------------------------------- 52 ! 52 53 53 c----------------------------------------------------------------------- 54 c 54 r = cpp * kappa 55 55 56 r = cpp * kappa 56 write(lunout,*) trim(modname),': R CP Kappa ',r,cpp,kappa 57 ! 58 !----------------------------------------------------------------------- 57 59 58 write(lunout,*) trim(modname),': R CP Kappa ',r,cpp,kappa 59 c 60 c----------------------------------------------------------------------- 60 ! vertical discretization: default behavior depends on planet_type flag 61 if (planet_type=="earth") then 62 disvert_type=1 63 else 64 disvert_type=2 65 endif 66 ! but user can also specify using one or the other in run.def: 67 call getin('disvert_type',disvert_type) 68 write(lunout,*) trim(modname),': disvert_type=',disvert_type 61 69 62 ! vertical discretization: default behavior depends on planet_type flag 63 if (planet_type=="earth") then 64 disvert_type=1 65 else 66 disvert_type=2 67 endif 68 ! but user can also specify using one or the other in run.def: 69 call getin('disvert_type',disvert_type) 70 write(lunout,*) trim(modname),': disvert_type=',disvert_type 71 72 if (disvert_type==1) then 73 ! standard case for Earth (automatic generation of levels) 74 call disvert(pa,preff,ap,bp,dpres,presnivs,nivsigs,nivsig, 75 & scaleheight) 76 else if (disvert_type==2) then 77 ! standard case for planets (levels generated using z2sig.def file) 78 call disvert_noterre 79 else 80 write(abort_message,*) "Wrong value for disvert_type: ", 81 & disvert_type 82 call abort_gcm(modname,abort_message,0) 83 endif 70 pressure_exner = disvert_type == 1 ! default value 71 call getin('pressure_exner', pressure_exner) 84 72 85 END 73 if (disvert_type==1) then 74 ! standard case for Earth (automatic generation of levels) 75 call disvert(pa,preff,ap,bp,dpres,presnivs,nivsigs,nivsig, scaleheight) 76 else if (disvert_type==2) then 77 ! standard case for planets (levels generated using z2sig.def file) 78 call disvert_noterre 79 else 80 write(abort_message,*) "Wrong value for disvert_type: ", disvert_type 81 call abort_gcm(modname,abort_message,0) 82 endif 83 84 END SUBROUTINE iniconst -
trunk/LMDZ.COMMON/libf/dyn3d/inidissip.F90
r270 r776 28 28 ! Local variables: 29 29 REAL fact,zvert(llm),zz 30 REAL zh(ip1jmp1),zu(ip1jmp1),zv(ip1jm),deltap(ip1jmp1,llm) 30 REAL zh(ip1jmp1),zu(ip1jmp1), gx(ip1jmp1), divgra(ip1jmp1) 31 real zv(ip1jm), gy(ip1jm), deltap(ip1jmp1,llm) 31 32 REAL ullm,vllm,umin,vmin,zhmin,zhmax 32 REAL zllm ,z1llm33 REAL zllm 33 34 34 35 INTEGER l,ij,idum,ii … … 78 79 DO l = 1,50 79 80 IF(lstardis) THEN 80 CALL divgrad2(1,zh,deltap,niterh, zh)81 CALL divgrad2(1,zh,deltap,niterh,divgra) 81 82 ELSE 82 CALL divgrad (1,zh,niterh, zh)83 CALL divgrad (1,zh,niterh,divgra) 83 84 ENDIF 84 85 85 CALL minmax(iip1*jjp1,zh,zhmin,zhmax ) 86 87 zllm = ABS( zhmax ) 88 z1llm = 1./zllm 89 DO ij = 1,ip1jmp1 90 zh(ij) = zh(ij)* z1llm 91 ENDDO 86 zllm = ABS(maxval(divgra)) 87 zh = divgra / zllm 92 88 ENDDO 93 89 … … 123 119 !cccc CALL covcont( 1,zu,zv,zu,zv ) 124 120 IF(lstardis) THEN 125 CALL gradiv2( 1,zu,zv,nitergdiv, zu,zv)121 CALL gradiv2( 1,zu,zv,nitergdiv,gx,gy ) 126 122 ELSE 127 CALL gradiv ( 1,zu,zv,nitergdiv, zu,zv)123 CALL gradiv ( 1,zu,zv,nitergdiv,gx,gy ) 128 124 ENDIF 129 125 ELSE 130 126 IF(lstardis) THEN 131 CALL nxgraro2( 1,zu,zv,nitergrot, zu,zv)127 CALL nxgraro2( 1,zu,zv,nitergrot,gx,gy ) 132 128 ELSE 133 CALL nxgrarot( 1,zu,zv,nitergrot, zu,zv)129 CALL nxgrarot( 1,zu,zv,nitergrot,gx,gy ) 134 130 ENDIF 135 131 ENDIF 136 132 137 CALL minmax(iip1*jjp1,zu,umin,ullm ) 138 CALL minmax(iip1*jjm, zv,vmin,vllm ) 139 140 ullm = ABS ( ullm ) 141 vllm = ABS ( vllm ) 142 143 zllm = MAX( ullm,vllm ) 144 z1llm = 1./ zllm 145 DO ij = 1, ip1jmp1 146 zu(ij) = zu(ij)* z1llm 147 ENDDO 148 DO ij = 1, ip1jm 149 zv(ij) = zv(ij)* z1llm 150 ENDDO 133 zllm = max(abs(maxval(gx)), abs(maxval(gy))) 134 zu = gx / zllm 135 zv = gy / zllm 151 136 end DO 152 137 -
trunk/LMDZ.COMMON/libf/dyn3d/inigrads.F
r1 r776 9 9 implicit none 10 10 11 integer if,im,jm,lm,i,j,l ,lnblnk11 integer if,im,jm,lm,i,j,l 12 12 real x(im),y(jm),z(lm),fx,fy,fz,dt 13 13 real xmin,xmax,ymin,ymax … … 40 40 ivar(if)=0 41 41 42 fichier(if)= file(1:lnblnk(file))42 fichier(if)=trim(file) 43 43 44 44 firsttime(if)=.true. … … 70 70 71 71 print*,4*(ifd(if)-iid(if))*(jfd(if)-jid(if)) 72 print*, file(1:lnblnk(file))//'.dat'72 print*,trim(file)//'.dat' 73 73 74 OPEN (unit(if)+1,FILE= file(1:lnblnk(file))//'.dat'74 OPEN (unit(if)+1,FILE=trim(file)//'.dat' 75 75 s ,FORM='unformatted', 76 76 s ACCESS='direct' -
trunk/LMDZ.COMMON/libf/dyn3d/integrd.F
r270 r776 1 1 ! 2 ! $Id: integrd.F 1 550 2011-07-05 09:44:55Z lguez$2 ! $Id: integrd.F 1616 2012-02-17 11:59:00Z emillour $ 3 3 ! 4 4 SUBROUTINE integrd 5 5 $ ( nq,vcovm1,ucovm1,tetam1,psm1,massem1, 6 $ dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis,finvmaold ) 6 $ dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis !,finvmaold 7 & ) 7 8 8 9 use control_mod, only : planet_type … … 34 35 #include "temps.h" 35 36 #include "serre.h" 37 #include "iniprint.h" 36 38 37 39 c Arguments: 38 40 c ---------- 39 41 40 INTEGER nq 41 42 REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm) 43 REAL q(ip1jmp1,llm,nq) 44 REAL ps(ip1jmp1),masse(ip1jmp1,llm),phis(ip1jmp1) 45 46 REAL vcovm1(ip1jm,llm),ucovm1(ip1jmp1,llm) 47 REAL tetam1(ip1jmp1,llm),psm1(ip1jmp1),massem1(ip1jmp1,llm) 48 49 REAL dv(ip1jm,llm),du(ip1jmp1,llm) 50 REAL dteta(ip1jmp1,llm),dp(ip1jmp1) 51 REAL dq(ip1jmp1,llm,nq), finvmaold(ip1jmp1,llm) 42 integer,intent(in) :: nq ! number of tracers to handle in this routine 43 real,intent(inout) :: vcov(ip1jm,llm) ! covariant meridional wind 44 real,intent(inout) :: ucov(ip1jmp1,llm) ! covariant zonal wind 45 real,intent(inout) :: teta(ip1jmp1,llm) ! potential temperature 46 real,intent(inout) :: q(ip1jmp1,llm,nq) ! advected tracers 47 real,intent(inout) :: ps(ip1jmp1) ! surface pressure 48 real,intent(inout) :: masse(ip1jmp1,llm) ! atmospheric mass 49 real,intent(in) :: phis(ip1jmp1) ! ground geopotential !!! unused 50 ! values at previous time step 51 real,intent(inout) :: vcovm1(ip1jm,llm) 52 real,intent(inout) :: ucovm1(ip1jmp1,llm) 53 real,intent(inout) :: tetam1(ip1jmp1,llm) 54 real,intent(inout) :: psm1(ip1jmp1) 55 real,intent(inout) :: massem1(ip1jmp1,llm) 56 ! the tendencies to add 57 real,intent(in) :: dv(ip1jm,llm) 58 real,intent(in) :: du(ip1jmp1,llm) 59 real,intent(in) :: dteta(ip1jmp1,llm) 60 real,intent(in) :: dp(ip1jmp1) 61 real,intent(in) :: dq(ip1jmp1,llm,nq) !!! unused 62 ! real,intent(out) :: finvmaold(ip1jmp1,llm) !!! unused 52 63 53 64 c Local: … … 55 66 56 67 REAL vscr( ip1jm ),uscr( ip1jmp1 ),hscr( ip1jmp1 ),pscr(ip1jmp1) 57 REAL massescr( ip1jmp1,llm ), finvmasse(ip1jmp1,llm) 68 REAL massescr( ip1jmp1,llm ) 69 ! REAL finvmasse(ip1jmp1,llm) 58 70 REAL p(ip1jmp1,llmp1) 59 71 REAL tpn,tps,tppn(iim),tpps(iim) … … 61 73 REAL deltap( ip1jmp1,llm ) 62 74 63 INTEGER l,ij,iq 75 INTEGER l,ij,iq,i,j 64 76 65 77 REAL SSUM … … 88 100 DO ij = 1,ip1jmp1 89 101 IF( ps(ij).LT.0. ) THEN 90 PRINT*,' Au point ij = ',ij, ' , pression sol neg. ', ps(ij) 91 print *, ' dans integrd' 92 stop 1 102 write(lunout,*) "integrd: negative surface pressure ",ps(ij) 103 write(lunout,*) " at node ij =", ij 104 ! since ij=j+(i-1)*jjp1 , we have 105 j=modulo(ij,jjp1) 106 i=1+(ij-j)/jjp1 107 write(lunout,*) " lon = ",rlonv(i)*180./pi, " deg", 108 & " lat = ",rlatu(j)*180./pi, " deg" 109 stop 93 110 ENDIF 94 111 ENDDO … … 110 127 CALL massdair ( p , masse ) 111 128 112 CALL SCOPY( ijp1llm , masse, 1, finvmasse, 1 ) 113 CALL filtreg( finvmasse, jjp1, llm, -2, 2, .TRUE., 1 ) 129 ! Ehouarn : we don't use/need finvmaold and finvmasse, 130 ! so might as well not compute them 131 ! CALL SCOPY( ijp1llm , masse, 1, finvmasse, 1 ) 132 ! CALL filtreg( finvmasse, jjp1, llm, -2, 2, .TRUE., 1 ) 114 133 c 115 134 … … 218 237 ENDDO 219 238 220 221 CALL SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )239 ! Ehouarn: forget about finvmaold 240 ! CALL SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 ) 222 241 223 242 endif ! of if (planet_type.eq."earth") -
trunk/LMDZ.COMMON/libf/dyn3d/leapfrog.F
r500 r776 124 124 125 125 REAL SSUM 126 REAL time_0 , finvmaold(ip1jmp1,llm) 126 REAL time_0 127 ! REAL finvmaold(ip1jmp1,llm) 127 128 128 129 cym LOGICAL lafin … … 243 244 dq(:,:,:)=0. 244 245 CALL pression ( ip1jmp1, ap, bp, ps, p ) 245 if ( disvert_type==1) then246 if (pressure_exner) then 246 247 CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf ) 247 else ! we assume that we are in the disvert_type==2 case248 else 248 249 CALL exner_milieu( ip1jmp1, ps, p, beta, pks, pk, pkf ) 249 250 endif … … 271 272 c ---------------------------------- 272 273 273 1 CONTINUE 274 1 CONTINUE ! Matsuno Forward step begins here 274 275 275 276 jD_cur = jD_ref + day_ini - day_ref + & 276 & i nt (itau * dtvr / daysec)277 & itau/day_step 277 278 jH_cur = jH_ref + start_time + & 278 & (itau * dtvr / daysec - int(itau * dtvr / daysec))279 & mod(itau,day_step)/float(day_step) 279 280 jD_cur = jD_cur + int(jH_cur) 280 281 jH_cur = jH_cur - int(jH_cur) … … 307 308 308 309 c ... P.Le Van .26/04/94 .... 309 310 CALL SCOPY ( ijp1llm, masse, 1, finvmaold, 1 )311 CALL filtreg ( finvmaold ,jjp1, llm, -2,2, .TRUE., 1 )312 313 2 CONTINUE 310 ! Ehouarn: finvmaold is actually not used 311 ! CALL SCOPY ( ijp1llm, masse, 1, finvmaold, 1 ) 312 ! CALL filtreg ( finvmaold ,jjp1, llm, -2,2, .TRUE., 1 ) 313 314 2 CONTINUE ! Matsuno backward or leapfrog step begins here 314 315 315 316 c----------------------------------------------------------------------- … … 357 358 call tpot2t(ijp1llm,teta,temp,pk) 358 359 tsurpk = cpp*temp/pk 360 ! compute geopotential phi() 359 361 CALL geopot ( ip1jmp1, tsurpk , pk , pks, phis , phi ) 360 362 … … 372 374 373 375 ! IF( forward. OR . leapf ) THEN 376 ! Ehouarn: NB: at this point p with ps are not synchronized 377 ! (whereas mass and ps are...) 374 378 IF((.not.forward).OR. leapf ) THEN 375 379 ! Ehouarn: gather mass fluxes during backward Matsuno or LF step … … 398 402 399 403 CALL integrd ( 2,vcovm1,ucovm1,tetam1,psm1,massem1 , 400 $ dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis ,401 $ finvmaold )404 $ dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis ) 405 ! $ finvmaold ) 402 406 403 407 IF ((planet_type.eq."titan").and.(tidal)) then … … 431 435 432 436 CALL pression ( ip1jmp1, ap, bp, ps, p ) 433 if ( disvert_type==1) then437 if (pressure_exner) then 434 438 CALL exner_hyb( ip1jmp1, ps, p,alpha,beta,pks, pk, pkf ) 435 else ! we assume that we are in the disvert_type==2 case439 else 436 440 CALL exner_milieu( ip1jmp1, ps, p, beta, pks, pk, pkf ) 437 441 endif 438 442 439 443 jD_cur = jD_ref + day_ini - day_ref + & 440 & i nt (itau * dtvr / daysec)444 & itau/day_step 441 445 jH_cur = jH_ref + start_time + & 442 & (itau * dtvr / daysec - int(itau * dtvr / daysec))446 & mod(itau,day_step)/float(day_step) 443 447 jD_cur = jD_cur + int(jH_cur) 444 448 jH_cur = jH_cur - int(jH_cur) … … 545 549 546 550 CALL pression ( ip1jmp1, ap, bp, ps, p ) 547 if ( disvert_type==1) then551 if (pressure_exner) then 548 552 CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf ) 549 else ! we assume that we are in the disvert_type==2 case553 else 550 554 CALL exner_milieu( ip1jmp1, ps, p, beta, pks, pk, pkf ) 551 555 endif … … 613 617 ENDDO 614 618 615 DO ij = 1,iim 616 tppn(ij) = aire( ij ) * ps ( ij ) 617 tpps(ij) = aire(ij+ip1jm) * ps (ij+ip1jm) 618 ENDDO 619 tpn = SSUM(iim,tppn,1)/apoln 620 tps = SSUM(iim,tpps,1)/apols 621 622 DO ij = 1, iip1 623 ps( ij ) = tpn 624 ps(ij+ip1jm) = tps 625 ENDDO 626 619 if (1 == 0) then 620 !!! Ehouarn: lines here 1) kill 1+1=2 in the dynamics 621 !!! 2) should probably not be here anyway 622 !!! but are kept for those who would want to revert to previous behaviour 623 DO ij = 1,iim 624 tppn(ij) = aire( ij ) * ps ( ij ) 625 tpps(ij) = aire(ij+ip1jm) * ps (ij+ip1jm) 626 ENDDO 627 tpn = SSUM(iim,tppn,1)/apoln 628 tps = SSUM(iim,tpps,1)/apols 629 630 DO ij = 1, iip1 631 ps( ij ) = tpn 632 ps(ij+ip1jm) = tps 633 ENDDO 634 endif ! of if (1 == 0) 627 635 628 636 END IF ! of IF(apdiss) … … 749 757 750 758 CLOSE(99) 759 !!! Ehouarn: Why not stop here and now? 751 760 ENDIF ! of IF (itau.EQ.itaufin) 752 761 -
trunk/LMDZ.COMMON/libf/dyn3d/wrgrads.F
r1 r776 26 26 c local 27 27 28 integer im,jm,lm,i,j,l, lnblnk,iv,iii,iji,iif,ijf28 integer im,jm,lm,i,j,l,iv,iii,iji,iif,ijf 29 29 30 30 logical writectl … … 59 59 nvar(if)=ivar(if) 60 60 var(ivar(if),if)=name 61 tvar(ivar(if),if)=t itlevar(1:lnblnk(titlevar))61 tvar(ivar(if),if)=trim(titlevar) 62 62 nld(ivar(if),if)=nl 63 63 c print*,'initialisation ecriture de ',var(ivar(if),if) … … 101 101 file=fichier(if) 102 102 c WARNING! on reecrase le fichier .ctl a chaque ecriture 103 open(unit(if),file= file(1:lnblnk(file))//'.ctl'103 open(unit(if),file=trim(file)//'.ctl' 104 104 & ,form='formatted',status='unknown') 105 105 write(unit(if),'(a5,1x,a40)') 106 & 'DSET ','^'// file(1:lnblnk(file))//'.dat'106 & 'DSET ','^'//trim(file)//'.dat' 107 107 108 108 write(unit(if),'(a12)') 'UNDEF 1.0E30'
Note: See TracChangeset
for help on using the changeset viewer.