| Directory: | ./ |
|---|---|
| File: | dyn/leapfrog.f |
| Date: | 2022-01-11 19:19:34 |
| Exec | Total | Coverage | |
|---|---|---|---|
| Lines: | 138 | 243 | 56.8% |
| Branches: | 111 | 218 | 50.9% |
| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | ! | ||
| 2 | ! $Id: leapfrog.F 3947 2021-06-21 15:15:02Z emillour $ | ||
| 3 | ! | ||
| 4 | c | ||
| 5 | c | ||
| 6 | 37441 | SUBROUTINE leapfrog(ucov,vcov,teta,ps,masse,phis,q,time_0) | |
| 7 | |||
| 8 | |||
| 9 | cIM : pour sortir les param. du modele dans un fis. netcdf 110106 | ||
| 10 | use IOIPSL | ||
| 11 | USE infotrac, ONLY: nqtot,ok_iso_verif | ||
| 12 | USE guide_mod, ONLY : guide_main | ||
| 13 | USE write_field, ONLY: writefield | ||
| 14 | USE control_mod, ONLY: nday, day_step, planet_type, offline, | ||
| 15 | & iconser, iphysiq, iperiod, dissip_period, | ||
| 16 | & iecri, ip_ebil_dyn, ok_dynzon, ok_dyn_ins, | ||
| 17 | & periodav, ok_dyn_ave, output_grads_dyn | ||
| 18 | use exner_hyb_m, only: exner_hyb | ||
| 19 | use exner_milieu_m, only: exner_milieu | ||
| 20 | USE comvert_mod, ONLY: ap,bp,pressure_exner,presnivs | ||
| 21 | USE comconst_mod, ONLY: cpp, dtphys, dtvr, pi, ihf | ||
| 22 | USE logic_mod, ONLY: iflag_phys,ok_guide,forward,leapf,apphys, | ||
| 23 | & statcl,conser,apdiss,purmats,ok_strato | ||
| 24 | USE temps_mod, ONLY: jD_ref,jH_ref,itaufin,day_ini,day_ref, | ||
| 25 | & start_time,dt | ||
| 26 | |||
| 27 | IMPLICIT NONE | ||
| 28 | |||
| 29 | c ...... Version du 10/01/98 .......... | ||
| 30 | |||
| 31 | c avec coordonnees verticales hybrides | ||
| 32 | c avec nouveaux operat. dissipation * ( gradiv2,divgrad2,nxgraro2 ) | ||
| 33 | |||
| 34 | c======================================================================= | ||
| 35 | c | ||
| 36 | c Auteur: P. Le Van /L. Fairhead/F.Hourdin | ||
| 37 | c ------- | ||
| 38 | c | ||
| 39 | c Objet: | ||
| 40 | c ------ | ||
| 41 | c | ||
| 42 | c GCM LMD nouvelle grille | ||
| 43 | c | ||
| 44 | c======================================================================= | ||
| 45 | c | ||
| 46 | c ... Dans inigeom , nouveaux calculs pour les elongations cu , cv | ||
| 47 | c et possibilite d'appeler une fonction f(y) a derivee tangente | ||
| 48 | c hyperbolique a la place de la fonction a derivee sinusoidale. | ||
| 49 | |||
| 50 | c ... Possibilite de choisir le shema pour l'advection de | ||
| 51 | c q , en modifiant iadv dans traceur.def (10/02) . | ||
| 52 | c | ||
| 53 | c Pour Van-Leer + Vapeur d'eau saturee, iadv(1)=4. (F.Codron,10/99) | ||
| 54 | c Pour Van-Leer iadv=10 | ||
| 55 | c | ||
| 56 | c----------------------------------------------------------------------- | ||
| 57 | c Declarations: | ||
| 58 | c ------------- | ||
| 59 | |||
| 60 | include "dimensions.h" | ||
| 61 | include "paramet.h" | ||
| 62 | include "comdissnew.h" | ||
| 63 | include "comgeom.h" | ||
| 64 | include "description.h" | ||
| 65 | include "iniprint.h" | ||
| 66 | include "academic.h" | ||
| 67 | |||
| 68 | REAL,INTENT(IN) :: time_0 ! not used | ||
| 69 | |||
| 70 | c dynamical variables: | ||
| 71 | REAL,INTENT(INOUT) :: ucov(ip1jmp1,llm) ! zonal covariant wind | ||
| 72 | REAL,INTENT(INOUT) :: vcov(ip1jm,llm) ! meridional covariant wind | ||
| 73 | REAL,INTENT(INOUT) :: teta(ip1jmp1,llm) ! potential temperature | ||
| 74 | REAL,INTENT(INOUT) :: ps(ip1jmp1) ! surface pressure (Pa) | ||
| 75 | REAL,INTENT(INOUT) :: masse(ip1jmp1,llm) ! air mass | ||
| 76 | REAL,INTENT(INOUT) :: phis(ip1jmp1) ! geopotentiat at the surface | ||
| 77 | REAL,INTENT(INOUT) :: q(ip1jmp1,llm,nqtot) ! advected tracers | ||
| 78 | |||
| 79 | REAL p (ip1jmp1,llmp1 ) ! interlayer pressure | ||
| 80 | REAL pks(ip1jmp1) ! exner at the surface | ||
| 81 | REAL pk(ip1jmp1,llm) ! exner at mid-layer | ||
| 82 | REAL pkf(ip1jmp1,llm) ! filtered exner at mid-layer | ||
| 83 | REAL phi(ip1jmp1,llm) ! geopotential | ||
| 84 | REAL w(ip1jmp1,llm) ! vertical velocity | ||
| 85 | |||
| 86 | real zqmin,zqmax | ||
| 87 | |||
| 88 | c variables dynamiques intermediaire pour le transport | ||
| 89 | REAL pbaru(ip1jmp1,llm),pbarv(ip1jm,llm) !flux de masse | ||
| 90 | |||
| 91 | c variables dynamiques au pas -1 | ||
| 92 | REAL vcovm1(ip1jm,llm),ucovm1(ip1jmp1,llm) | ||
| 93 | REAL tetam1(ip1jmp1,llm),psm1(ip1jmp1) | ||
| 94 | REAL massem1(ip1jmp1,llm) | ||
| 95 | |||
| 96 | c tendances dynamiques | ||
| 97 | REAL dv(ip1jm,llm),du(ip1jmp1,llm) | ||
| 98 | 1 | REAL dteta(ip1jmp1,llm),dq(ip1jmp1,llm,nqtot),dp(ip1jmp1) | |
| 99 | |||
| 100 | c tendances de la dissipation | ||
| 101 | REAL dvdis(ip1jm,llm),dudis(ip1jmp1,llm) | ||
| 102 | REAL dtetadis(ip1jmp1,llm) | ||
| 103 | |||
| 104 | c tendances physiques | ||
| 105 | REAL dvfi(ip1jm,llm),dufi(ip1jmp1,llm) | ||
| 106 | 1 | REAL dtetafi(ip1jmp1,llm),dqfi(ip1jmp1,llm,nqtot),dpfi(ip1jmp1) | |
| 107 | |||
| 108 | c variables pour le fichier histoire | ||
| 109 | REAL dtav ! intervalle de temps elementaire | ||
| 110 | |||
| 111 | REAL tppn(iim),tpps(iim),tpn,tps | ||
| 112 | c | ||
| 113 | INTEGER itau,itaufinp1,iav | ||
| 114 | ! INTEGER iday ! jour julien | ||
| 115 | REAL time | ||
| 116 | |||
| 117 | REAL SSUM | ||
| 118 | ! REAL finvmaold(ip1jmp1,llm) | ||
| 119 | |||
| 120 | cym LOGICAL lafin | ||
| 121 | LOGICAL :: lafin=.false. | ||
| 122 | INTEGER ij,iq,l | ||
| 123 | INTEGER ik | ||
| 124 | |||
| 125 | real time_step, t_wrt, t_ops | ||
| 126 | |||
| 127 | ! REAL rdayvrai,rdaym_ini | ||
| 128 | ! jD_cur: jour julien courant | ||
| 129 | ! jH_cur: heure julienne courante | ||
| 130 | REAL :: jD_cur, jH_cur | ||
| 131 | INTEGER :: an, mois, jour | ||
| 132 | REAL :: secondes | ||
| 133 | |||
| 134 | LOGICAL first,callinigrads | ||
| 135 | cIM : pour sortir les param. du modele dans un fis. netcdf 110106 | ||
| 136 | save first | ||
| 137 | data first/.true./ | ||
| 138 | real dt_cum | ||
| 139 | character*10 infile | ||
| 140 | integer zan, tau0, thoriid | ||
| 141 | integer nid_ctesGCM | ||
| 142 | save nid_ctesGCM | ||
| 143 | real degres | ||
| 144 | real rlong(iip1), rlatg(jjp1) | ||
| 145 | real zx_tmp_2d(iip1,jjp1) | ||
| 146 | integer ndex2d(iip1*jjp1) | ||
| 147 | logical ok_sync | ||
| 148 | parameter (ok_sync = .true.) | ||
| 149 | logical physic | ||
| 150 | |||
| 151 | data callinigrads/.true./ | ||
| 152 | character*10 string10 | ||
| 153 | |||
| 154 | REAL :: flxw(ip1jmp1,llm) ! flux de masse verticale | ||
| 155 | |||
| 156 | c+jld variables test conservation energie | ||
| 157 | REAL ecin(ip1jmp1,llm),ecin0(ip1jmp1,llm) | ||
| 158 | C Tendance de la temp. potentiel d (theta)/ d t due a la | ||
| 159 | C tansformation d'energie cinetique en energie thermique | ||
| 160 | C cree par la dissipation | ||
| 161 | REAL dtetaecdt(ip1jmp1,llm) | ||
| 162 | REAL vcont(ip1jm,llm),ucont(ip1jmp1,llm) | ||
| 163 | REAL vnat(ip1jm,llm),unat(ip1jmp1,llm) | ||
| 164 | REAL d_h_vcol, d_qt, d_qw, d_ql, d_ec | ||
| 165 | CHARACTER*15 ztit | ||
| 166 | !IM INTEGER ip_ebil_dyn ! PRINT level for energy conserv. diag. | ||
| 167 | !IM SAVE ip_ebil_dyn | ||
| 168 | !IM DATA ip_ebil_dyn/0/ | ||
| 169 | c-jld | ||
| 170 | |||
| 171 | character*80 dynhist_file, dynhistave_file | ||
| 172 | character(len=*),parameter :: modname="leapfrog" | ||
| 173 | character*80 abort_message | ||
| 174 | |||
| 175 | logical dissip_conservative | ||
| 176 | save dissip_conservative | ||
| 177 | data dissip_conservative/.true./ | ||
| 178 | |||
| 179 | LOGICAL prem | ||
| 180 | save prem | ||
| 181 | DATA prem/.true./ | ||
| 182 | INTEGER testita | ||
| 183 | PARAMETER (testita = 9) | ||
| 184 | |||
| 185 | logical , parameter :: flag_verif = .false. | ||
| 186 | |||
| 187 | |||
| 188 | integer itau_w ! pas de temps ecriture = itap + itau_phy | ||
| 189 | |||
| 190 | |||
| 191 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | if (nday>=0) then |
| 192 | 1 | itaufin = nday*day_step | |
| 193 | else | ||
| 194 | ✗ | itaufin = -nday | |
| 195 | endif | ||
| 196 | 1 | itaufinp1 = itaufin +1 | |
| 197 | 1 | itau = 0 | |
| 198 | physic=.true. | ||
| 199 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
|
1 | if (iflag_phys==0.or.iflag_phys==2) physic=.false. |
| 200 | |||
| 201 | c iday = day_ini+itau/day_step | ||
| 202 | c time = REAL(itau-(iday-day_ini)*day_step)/day_step+time_0 | ||
| 203 | c IF(time.GT.1.) THEN | ||
| 204 | c time = time-1. | ||
| 205 | c iday = iday+1 | ||
| 206 | c ENDIF | ||
| 207 | |||
| 208 | |||
| 209 | c----------------------------------------------------------------------- | ||
| 210 | c On initialise la pression et la fonction d'Exner : | ||
| 211 | c -------------------------------------------------- | ||
| 212 | |||
| 213 |
6/6✓ Branch 0 taken 5 times.
✓ Branch 1 taken 1 times.
✓ Branch 2 taken 195 times.
✓ Branch 3 taken 5 times.
✓ Branch 4 taken 212355 times.
✓ Branch 5 taken 195 times.
|
212556 | dq(:,:,:)=0. |
| 214 | 1 | CALL pression ( ip1jmp1, ap, bp, ps, p ) | |
| 215 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | if (pressure_exner) then |
| 216 | 1 | CALL exner_hyb( ip1jmp1, ps, p, pks, pk, pkf ) | |
| 217 | else | ||
| 218 | ✗ | CALL exner_milieu( ip1jmp1, ps, p, pks, pk, pkf ) | |
| 219 | endif | ||
| 220 | |||
| 221 | c----------------------------------------------------------------------- | ||
| 222 | c Debut de l'integration temporelle: | ||
| 223 | c ---------------------------------- | ||
| 224 | |||
| 225 | 1 CONTINUE ! Matsuno Forward step begins here | ||
| 226 | |||
| 227 | c date: (NB: date remains unchanged for Backward step) | ||
| 228 | c ----- | ||
| 229 | |||
| 230 | jD_cur = jD_ref + day_ini - day_ref + & | ||
| 231 | 481 | & (itau+1)/day_step | |
| 232 | jH_cur = jH_ref + start_time + & | ||
| 233 | 481 | & mod(itau+1,day_step)/float(day_step) | |
| 234 | 481 | jD_cur = jD_cur + int(jH_cur) | |
| 235 | 481 | jH_cur = jH_cur - int(jH_cur) | |
| 236 | |||
| 237 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 481 times.
|
481 | if (ok_iso_verif) then |
| 238 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 321') | |
| 239 | endif !if (ok_iso_verif) then | ||
| 240 | |||
| 241 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 481 times.
|
481 | if (ok_guide) then |
| 242 | ✗ | call guide_main(itau,ucov,vcov,teta,q,masse,ps) | |
| 243 | endif | ||
| 244 | |||
| 245 | |||
| 246 | c | ||
| 247 | c IF( MOD( itau, 10* day_step ).EQ.0 ) THEN | ||
| 248 | c CALL test_period ( ucov,vcov,teta,q,p,phis ) | ||
| 249 | c PRINT *,' ---- Test_period apres continue OK ! -----', itau | ||
| 250 | c ENDIF | ||
| 251 | c | ||
| 252 | |||
| 253 | ! Save fields obtained at previous time step as '...m1' | ||
| 254 | 481 | CALL SCOPY( ijmllm ,vcov , 1, vcovm1 , 1 ) | |
| 255 | 481 | CALL SCOPY( ijp1llm,ucov , 1, ucovm1 , 1 ) | |
| 256 | 481 | CALL SCOPY( ijp1llm,teta , 1, tetam1 , 1 ) | |
| 257 | 481 | CALL SCOPY( ijp1llm,masse, 1, massem1, 1 ) | |
| 258 | 481 | CALL SCOPY( ip1jmp1, ps , 1, psm1 , 1 ) | |
| 259 | |||
| 260 | 481 | forward = .TRUE. | |
| 261 | 481 | leapf = .FALSE. | |
| 262 | 481 | dt = dtvr | |
| 263 | |||
| 264 | c ... P.Le Van .26/04/94 .... | ||
| 265 | ! Ehouarn: finvmaold is actually not used | ||
| 266 | ! CALL SCOPY ( ijp1llm, masse, 1, finvmaold, 1 ) | ||
| 267 | ! CALL filtreg ( finvmaold ,jjp1, llm, -2,2, .TRUE., 1 ) | ||
| 268 | |||
| 269 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 481 times.
|
481 | if (ok_iso_verif) then |
| 270 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 400') | |
| 271 | endif !if (ok_iso_verif) then | ||
| 272 | |||
| 273 | 2 CONTINUE ! Matsuno backward or leapfrog step begins here | ||
| 274 | |||
| 275 | c----------------------------------------------------------------------- | ||
| 276 | |||
| 277 | c date: (NB: only leapfrog step requires recomputing date) | ||
| 278 | c ----- | ||
| 279 | |||
| 280 |
2/2✓ Branch 0 taken 1920 times.
✓ Branch 1 taken 961 times.
|
2881 | IF (leapf) THEN |
| 281 | jD_cur = jD_ref + day_ini - day_ref + | ||
| 282 | 1920 | & (itau+1)/day_step | |
| 283 | jH_cur = jH_ref + start_time + | ||
| 284 | 1920 | & mod(itau+1,day_step)/float(day_step) | |
| 285 | 1920 | jD_cur = jD_cur + int(jH_cur) | |
| 286 | 1920 | jH_cur = jH_cur - int(jH_cur) | |
| 287 | ENDIF | ||
| 288 | |||
| 289 | |||
| 290 | c gestion des appels de la physique et des dissipations: | ||
| 291 | c ------------------------------------------------------ | ||
| 292 | c | ||
| 293 | c ... P.Le Van ( 6/02/95 ) .... | ||
| 294 | |||
| 295 | 2881 | apphys = .FALSE. | |
| 296 | 2881 | statcl = .FALSE. | |
| 297 | 2881 | conser = .FALSE. | |
| 298 | 2881 | apdiss = .FALSE. | |
| 299 | |||
| 300 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2881 times.
|
2881 | IF( purmats ) THEN |
| 301 | ! Purely Matsuno time stepping | ||
| 302 | ✗ | IF( MOD(itau,iconser) .EQ.0.AND. forward ) conser = .TRUE. | |
| 303 | ✗ | IF( MOD(itau,dissip_period ).EQ.0.AND..NOT.forward ) | |
| 304 | ✗ | s apdiss = .TRUE. | |
| 305 | IF( MOD(itau,iphysiq ).EQ.0.AND..NOT.forward | ||
| 306 | ✗ | s .and. physic ) apphys = .TRUE. | |
| 307 | ELSE | ||
| 308 | ! Leapfrog/Matsuno time stepping | ||
| 309 |
2/2✓ Branch 0 taken 81 times.
✓ Branch 1 taken 2800 times.
|
2881 | IF( MOD(itau ,iconser) .EQ. 0 ) conser = .TRUE. |
| 310 |
3/4✓ Branch 0 taken 480 times.
✓ Branch 1 taken 2401 times.
✓ Branch 2 taken 480 times.
✗ Branch 3 not taken.
|
2881 | IF( MOD(itau+1,dissip_period).EQ.0 .AND. .NOT. forward ) |
| 311 | 480 | s apdiss = .TRUE. | |
| 312 |
3/4✓ Branch 0 taken 480 times.
✓ Branch 1 taken 2401 times.
✓ Branch 2 taken 480 times.
✗ Branch 3 not taken.
|
2881 | IF( MOD(itau+1,iphysiq).EQ.0.AND.physic ) apphys=.TRUE. |
| 313 | END IF | ||
| 314 | |||
| 315 | ! Ehouarn: for Shallow Water case (ie: 1 vertical layer), | ||
| 316 | ! supress dissipation step | ||
| 317 | if (llm.eq.1) then | ||
| 318 | apdiss=.false. | ||
| 319 | endif | ||
| 320 | |||
| 321 | |||
| 322 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2881 times.
|
2881 | if (ok_iso_verif) then |
| 323 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 589') | |
| 324 | endif !if (ok_iso_verif) then | ||
| 325 | |||
| 326 | c----------------------------------------------------------------------- | ||
| 327 | c calcul des tendances dynamiques: | ||
| 328 | c -------------------------------- | ||
| 329 | |||
| 330 | ! compute geopotential phi() | ||
| 331 | 2881 | CALL geopot ( ip1jmp1, teta , pk , pks, phis , phi ) | |
| 332 | |||
| 333 | 2881 | time = jD_cur + jH_cur | |
| 334 | CALL caldyn | ||
| 335 | $ ( itau,ucov,vcov,teta,ps,masse,pk,pkf,phis , | ||
| 336 | 2881 | $ phi,conser,du,dv,dteta,dp,w, pbaru,pbarv, time ) | |
| 337 | |||
| 338 | |||
| 339 | c----------------------------------------------------------------------- | ||
| 340 | c calcul des tendances advection des traceurs (dont l'humidite) | ||
| 341 | c ------------------------------------------------------------- | ||
| 342 | |||
| 343 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2881 times.
|
2881 | if (ok_iso_verif) then |
| 344 | call check_isotopes_seq(q,ip1jmp1, | ||
| 345 | ✗ | & 'leapfrog 686: avant caladvtrac') | |
| 346 | endif !if (ok_iso_verif) then | ||
| 347 | |||
| 348 |
4/4✓ Branch 0 taken 2400 times.
✓ Branch 1 taken 481 times.
✓ Branch 2 taken 1920 times.
✓ Branch 3 taken 480 times.
|
2881 | IF( forward. OR . leapf ) THEN |
| 349 | ! Ehouarn: NB: fields sent to advtrac are those at the beginning of the time step | ||
| 350 | CALL caladvtrac(q,pbaru,pbarv, | ||
| 351 | * p, masse, dq, teta, | ||
| 352 | 2401 | . flxw, pk) | |
| 353 | !write(*,*) 'caladvtrac 346' | ||
| 354 | |||
| 355 | |||
| 356 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2401 times.
|
2401 | IF (offline) THEN |
| 357 | Cmaf stokage du flux de masse pour traceurs OFF-LINE | ||
| 358 | |||
| 359 | CALL fluxstokenc(pbaru,pbarv,masse,teta,phi,phis, | ||
| 360 | ✗ | . dtvr, itau) | |
| 361 | |||
| 362 | |||
| 363 | ENDIF ! of IF (offline) | ||
| 364 | c | ||
| 365 | ENDIF ! of IF( forward. OR . leapf ) | ||
| 366 | |||
| 367 | |||
| 368 | c----------------------------------------------------------------------- | ||
| 369 | c integrations dynamique et traceurs: | ||
| 370 | c ---------------------------------- | ||
| 371 | |||
| 372 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2881 times.
|
2881 | if (ok_iso_verif) then |
| 373 | ✗ | write(*,*) 'leapfrog 720' | |
| 374 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 756') | |
| 375 | endif !if (ok_iso_verif) then | ||
| 376 | |||
| 377 | CALL integrd ( nqtot,vcovm1,ucovm1,tetam1,psm1,massem1 , | ||
| 378 | 2881 | $ dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis ) | |
| 379 | ! $ finvmaold ) | ||
| 380 | |||
| 381 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2881 times.
|
2881 | if (ok_iso_verif) then |
| 382 | ✗ | write(*,*) 'leapfrog 724' | |
| 383 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 762') | |
| 384 | endif !if (ok_iso_verif) then | ||
| 385 | |||
| 386 | c .P.Le Van (26/04/94 ajout de finvpold dans l'appel d'integrd) | ||
| 387 | c | ||
| 388 | c----------------------------------------------------------------------- | ||
| 389 | c calcul des tendances physiques: | ||
| 390 | c ------------------------------- | ||
| 391 | c ######## P.Le Van ( Modif le 6/02/95 ) ########### | ||
| 392 | c | ||
| 393 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2881 times.
|
2881 | IF( purmats ) THEN |
| 394 | ✗ | IF( itau.EQ.itaufin.AND..NOT.forward ) lafin = .TRUE. | |
| 395 | ELSE | ||
| 396 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 2880 times.
|
2881 | IF( itau+1. EQ. itaufin ) lafin = .TRUE. |
| 397 | ENDIF | ||
| 398 | c | ||
| 399 | c | ||
| 400 |
2/2✓ Branch 0 taken 480 times.
✓ Branch 1 taken 2401 times.
|
2881 | IF( apphys ) THEN |
| 401 | c | ||
| 402 | c ....... Ajout P.Le Van ( 17/04/96 ) ........... | ||
| 403 | c | ||
| 404 | |||
| 405 | 480 | CALL pression ( ip1jmp1, ap, bp, ps, p ) | |
| 406 |
1/2✓ Branch 0 taken 480 times.
✗ Branch 1 not taken.
|
480 | if (pressure_exner) then |
| 407 | 480 | CALL exner_hyb( ip1jmp1, ps, p,pks, pk, pkf ) | |
| 408 | else | ||
| 409 | ✗ | CALL exner_milieu( ip1jmp1, ps, p, pks, pk, pkf ) | |
| 410 | endif | ||
| 411 | |||
| 412 | ! Appel a geopot ajoute le 2014/05/08 pour garantir la convergence numerique | ||
| 413 | ! avec dyn3dmem | ||
| 414 | 480 | CALL geopot ( ip1jmp1, teta , pk , pks, phis , phi ) | |
| 415 | |||
| 416 | ! rdaym_ini = itau * dtvr / daysec | ||
| 417 | ! rdayvrai = rdaym_ini + day_ini | ||
| 418 | ! jD_cur = jD_ref + day_ini - day_ref | ||
| 419 | ! $ + int (itau * dtvr / daysec) | ||
| 420 | ! jH_cur = jH_ref + & | ||
| 421 | ! & (itau * dtvr / daysec - int(itau * dtvr / daysec)) | ||
| 422 | jD_cur = jD_ref + day_ini - day_ref + & | ||
| 423 | 480 | & (itau+1)/day_step | |
| 424 | |||
| 425 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
|
480 | IF (planet_type .eq."generic") THEN |
| 426 | ! AS: we make jD_cur to be pday | ||
| 427 | ✗ | jD_cur = int(day_ini + itau/day_step) | |
| 428 | ENDIF | ||
| 429 | |||
| 430 | jH_cur = jH_ref + start_time + & | ||
| 431 | 480 | & mod(itau+1,day_step)/float(day_step) | |
| 432 | 480 | jD_cur = jD_cur + int(jH_cur) | |
| 433 | 480 | jH_cur = jH_cur - int(jH_cur) | |
| 434 | ! write(lunout,*)'itau, jD_cur = ', itau, jD_cur, jH_cur | ||
| 435 | ! call ju2ymds(jD_cur+jH_cur, an, mois, jour, secondes) | ||
| 436 | ! write(lunout,*)'current date = ',an, mois, jour, secondes | ||
| 437 | |||
| 438 | c rajout debug | ||
| 439 | c lafin = .true. | ||
| 440 | |||
| 441 | |||
| 442 | c Inbterface avec les routines de phylmd (phymars ... ) | ||
| 443 | c ----------------------------------------------------- | ||
| 444 | |||
| 445 | c+jld | ||
| 446 | |||
| 447 | c Diagnostique de conservation de l'�nergie : initialisation | ||
| 448 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
|
480 | IF (ip_ebil_dyn.ge.1 ) THEN |
| 449 | ✗ | ztit='bil dyn' | |
| 450 | ! Ehouarn: be careful, diagedyn is Earth-specific! | ||
| 451 | ✗ | IF (planet_type.eq."earth") THEN | |
| 452 | CALL diagedyn(ztit,2,1,1,dtphys | ||
| 453 | ✗ | & , ucov , vcov , ps, p ,pk , teta , q(:,:,1), q(:,:,2)) | |
| 454 | ENDIF | ||
| 455 | ENDIF ! of IF (ip_ebil_dyn.ge.1 ) | ||
| 456 | c-jld | ||
| 457 | cIM decommenter les 6 lignes suivantes pour sortir quelques parametres dynamiques de LMDZ | ||
| 458 | cIM uncomment next 6 lines to get some parameters for LMDZ dynamics | ||
| 459 | c IF (first) THEN | ||
| 460 | c first=.false. | ||
| 461 | c#include "ini_paramLMDZ_dyn.h" | ||
| 462 | c ENDIF | ||
| 463 | c | ||
| 464 | c#include "write_paramLMDZ_dyn.h" | ||
| 465 | c | ||
| 466 | ! #endif of #ifdef 1 | ||
| 467 | CALL calfis( lafin , jD_cur, jH_cur, | ||
| 468 | $ ucov,vcov,teta,q,masse,ps,p,pk,phis,phi , | ||
| 469 | $ du,dv,dteta,dq, | ||
| 470 | 480 | $ flxw,dufi,dvfi,dtetafi,dqfi,dpfi ) | |
| 471 | c ajout des tendances physiques: | ||
| 472 | c ------------------------------ | ||
| 473 | CALL addfi( dtphys, leapf, forward , | ||
| 474 | $ ucov, vcov, teta , q ,ps , | ||
| 475 | 480 | $ dufi, dvfi, dtetafi , dqfi ,dpfi ) | |
| 476 | ! since addfi updates ps(), also update p(), masse() and pk() | ||
| 477 | 480 | CALL pression (ip1jmp1,ap,bp,ps,p) | |
| 478 | 480 | CALL massdair(p,masse) | |
| 479 |
1/2✓ Branch 0 taken 480 times.
✗ Branch 1 not taken.
|
480 | if (pressure_exner) then |
| 480 | 480 | CALL exner_hyb(ip1jmp1,ps,p,pks,pk,pkf) | |
| 481 | else | ||
| 482 | ✗ | CALL exner_milieu(ip1jmp1,ps,p,pks,pk,pkf) | |
| 483 | endif | ||
| 484 | |||
| 485 |
1/2✓ Branch 0 taken 480 times.
✗ Branch 1 not taken.
|
480 | IF (ok_strato) THEN |
| 486 | 480 | CALL top_bound( vcov,ucov,teta,masse,dtphys) | |
| 487 | ENDIF | ||
| 488 | |||
| 489 | c | ||
| 490 | c Diagnostique de conservation de l'�nergie : difference | ||
| 491 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
|
480 | IF (ip_ebil_dyn.ge.1 ) THEN |
| 492 | ✗ | ztit='bil phys' | |
| 493 | ✗ | IF (planet_type.eq."earth") THEN | |
| 494 | CALL diagedyn(ztit,2,1,1,dtphys | ||
| 495 | ✗ | & , ucov , vcov , ps, p ,pk , teta , q(:,:,1), q(:,:,2)) | |
| 496 | ENDIF | ||
| 497 | ENDIF ! of IF (ip_ebil_dyn.ge.1 ) | ||
| 498 | |||
| 499 | ENDIF ! of IF( apphys ) | ||
| 500 | |||
| 501 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2881 times.
|
2881 | IF(iflag_phys.EQ.2) THEN ! "Newtonian" case |
| 502 | ! Academic case : Simple friction and Newtonan relaxation | ||
| 503 | ! ------------------------------------------------------- | ||
| 504 | ✗ | DO l=1,llm | |
| 505 | ✗ | DO ij=1,ip1jmp1 | |
| 506 | teta(ij,l)=teta(ij,l)-dtvr* | ||
| 507 | ✗ | & (teta(ij,l)-tetarappel(ij,l))*(knewt_g+knewt_t(l)*clat4(ij)) | |
| 508 | ENDDO | ||
| 509 | ENDDO ! of DO l=1,llm | ||
| 510 | |||
| 511 | ✗ | if (planet_type.eq."giant") then | |
| 512 | ! add an intrinsic heat flux at the base of the atmosphere | ||
| 513 | ✗ | teta(:,1)=teta(:,1)+dtvr*aire(:)*ihf/cpp/masse(:,1) | |
| 514 | endif | ||
| 515 | |||
| 516 | ✗ | call friction(ucov,vcov,dtvr) | |
| 517 | |||
| 518 | ! Sponge layer (if any) | ||
| 519 | ✗ | IF (ok_strato) THEN | |
| 520 | ! dufi(:,:)=0. | ||
| 521 | ! dvfi(:,:)=0. | ||
| 522 | ! dtetafi(:,:)=0. | ||
| 523 | ! dqfi(:,:,:)=0. | ||
| 524 | ! dpfi(:)=0. | ||
| 525 | ! CALL top_bound(vcov,ucov,teta,masse,dufi,dvfi,dtetafi) | ||
| 526 | ✗ | CALL top_bound( vcov,ucov,teta,masse,dtvr) | |
| 527 | ! CALL addfi( dtvr, leapf, forward , | ||
| 528 | ! $ ucov, vcov, teta , q ,ps , | ||
| 529 | ! $ dufi, dvfi, dtetafi , dqfi ,dpfi ) | ||
| 530 | ENDIF ! of IF (ok_strato) | ||
| 531 | ENDIF ! of IF (iflag_phys.EQ.2) | ||
| 532 | |||
| 533 | |||
| 534 | c-jld | ||
| 535 | |||
| 536 | 2881 | CALL pression ( ip1jmp1, ap, bp, ps, p ) | |
| 537 |
1/2✓ Branch 0 taken 2881 times.
✗ Branch 1 not taken.
|
2881 | if (pressure_exner) then |
| 538 | 2881 | CALL exner_hyb( ip1jmp1, ps, p, pks, pk, pkf ) | |
| 539 | else | ||
| 540 | ✗ | CALL exner_milieu( ip1jmp1, ps, p, pks, pk, pkf ) | |
| 541 | endif | ||
| 542 | 2881 | CALL massdair(p,masse) | |
| 543 | |||
| 544 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2881 times.
|
2881 | if (ok_iso_verif) then |
| 545 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 1196') | |
| 546 | endif !if (ok_iso_verif) then | ||
| 547 | |||
| 548 | c----------------------------------------------------------------------- | ||
| 549 | c dissipation horizontale et verticale des petites echelles: | ||
| 550 | c ---------------------------------------------------------- | ||
| 551 | |||
| 552 |
2/2✓ Branch 0 taken 480 times.
✓ Branch 1 taken 2401 times.
|
2881 | IF(apdiss) THEN |
| 553 | |||
| 554 | |||
| 555 | c calcul de l'energie cinetique avant dissipation | ||
| 556 | 480 | call covcont(llm,ucov,vcov,ucont,vcont) | |
| 557 | 480 | call enercin(vcov,ucov,vcont,ucont,ecin0) | |
| 558 | |||
| 559 | c dissipation | ||
| 560 | 480 | CALL dissip(vcov,ucov,teta,p,dvdis,dudis,dtetadis) | |
| 561 |
4/4✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 20386080 times.
✓ Branch 3 taken 18720 times.
|
20405280 | ucov=ucov+dudis |
| 562 |
4/4✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 19768320 times.
✓ Branch 3 taken 18720 times.
|
19787520 | vcov=vcov+dvdis |
| 563 | c teta=teta+dtetadis | ||
| 564 | |||
| 565 | |||
| 566 | c------------------------------------------------------------------------ | ||
| 567 |
1/2✓ Branch 0 taken 480 times.
✗ Branch 1 not taken.
|
480 | if (dissip_conservative) then |
| 568 | C On rajoute la tendance due a la transform. Ec -> E therm. cree | ||
| 569 | C lors de la dissipation | ||
| 570 | 480 | call covcont(llm,ucov,vcov,ucont,vcont) | |
| 571 | 480 | call enercin(vcov,ucov,vcont,ucont,ecin) | |
| 572 |
4/4✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 20386080 times.
✓ Branch 3 taken 18720 times.
|
20405280 | dtetaecdt= (ecin0-ecin)/ pk |
| 573 | c teta=teta+dtetaecdt | ||
| 574 |
4/4✓ Branch 0 taken 480 times.
✓ Branch 1 taken 18720 times.
✓ Branch 2 taken 20386080 times.
✓ Branch 3 taken 18720 times.
|
20405280 | dtetadis=dtetadis+dtetaecdt |
| 575 | endif | ||
| 576 |
4/4✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 20386080 times.
✓ Branch 3 taken 18720 times.
|
20405280 | teta=teta+dtetadis |
| 577 | c------------------------------------------------------------------------ | ||
| 578 | |||
| 579 | |||
| 580 | c ....... P. Le Van ( ajout le 17/04/96 ) ........... | ||
| 581 | c ... Calcul de la valeur moyenne, unique de h aux poles ..... | ||
| 582 | c | ||
| 583 | |||
| 584 |
2/2✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
|
19200 | DO l = 1, llm |
| 585 |
2/2✓ Branch 0 taken 599040 times.
✓ Branch 1 taken 18720 times.
|
617760 | DO ij = 1,iim |
| 586 | 599040 | tppn(ij) = aire( ij ) * teta( ij ,l) | |
| 587 | 617760 | tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l) | |
| 588 | ENDDO | ||
| 589 | 18720 | tpn = SSUM(iim,tppn,1)/apoln | |
| 590 | 18720 | tps = SSUM(iim,tpps,1)/apols | |
| 591 | |||
| 592 |
2/2✓ Branch 0 taken 617760 times.
✓ Branch 1 taken 18720 times.
|
636960 | DO ij = 1, iip1 |
| 593 | 617760 | teta( ij ,l) = tpn | |
| 594 | 636480 | teta(ij+ip1jm,l) = tps | |
| 595 | ENDDO | ||
| 596 | ENDDO | ||
| 597 | |||
| 598 | if (1 == 0) then | ||
| 599 | !!! Ehouarn: lines here 1) kill 1+1=2 in the dynamics | ||
| 600 | !!! 2) should probably not be here anyway | ||
| 601 | !!! but are kept for those who would want to revert to previous behaviour | ||
| 602 | DO ij = 1,iim | ||
| 603 | tppn(ij) = aire( ij ) * ps ( ij ) | ||
| 604 | tpps(ij) = aire(ij+ip1jm) * ps (ij+ip1jm) | ||
| 605 | ENDDO | ||
| 606 | tpn = SSUM(iim,tppn,1)/apoln | ||
| 607 | tps = SSUM(iim,tpps,1)/apols | ||
| 608 | |||
| 609 | DO ij = 1, iip1 | ||
| 610 | ps( ij ) = tpn | ||
| 611 | ps(ij+ip1jm) = tps | ||
| 612 | ENDDO | ||
| 613 | endif ! of if (1 == 0) | ||
| 614 | |||
| 615 | END IF ! of IF(apdiss) | ||
| 616 | |||
| 617 | c ajout debug | ||
| 618 | c IF( lafin ) then | ||
| 619 | c abort_message = 'Simulation finished' | ||
| 620 | c call abort_gcm(modname,abort_message,0) | ||
| 621 | c ENDIF | ||
| 622 | |||
| 623 | c ******************************************************************** | ||
| 624 | c ******************************************************************** | ||
| 625 | c .... fin de l'integration dynamique et physique pour le pas itau .. | ||
| 626 | c ******************************************************************** | ||
| 627 | c ******************************************************************** | ||
| 628 | |||
| 629 | c preparation du pas d'integration suivant ...... | ||
| 630 | |||
| 631 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2881 times.
|
2881 | if (ok_iso_verif) then |
| 632 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 1509') | |
| 633 | endif !if (ok_iso_verif) then | ||
| 634 | |||
| 635 |
1/2✓ Branch 0 taken 2881 times.
✗ Branch 1 not taken.
|
2881 | IF ( .NOT.purmats ) THEN |
| 636 | c ........................................................ | ||
| 637 | c .............. schema matsuno + leapfrog .............. | ||
| 638 | c ........................................................ | ||
| 639 | |||
| 640 |
4/4✓ Branch 0 taken 2400 times.
✓ Branch 1 taken 481 times.
✓ Branch 2 taken 1920 times.
✓ Branch 3 taken 480 times.
|
2881 | IF(forward. OR. leapf) THEN |
| 641 | 2401 | itau= itau + 1 | |
| 642 | c iday= day_ini+itau/day_step | ||
| 643 | c time= REAL(itau-(iday-day_ini)*day_step)/day_step+time_0 | ||
| 644 | c IF(time.GT.1.) THEN | ||
| 645 | c time = time-1. | ||
| 646 | c iday = iday+1 | ||
| 647 | c ENDIF | ||
| 648 | ENDIF | ||
| 649 | |||
| 650 | |||
| 651 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 2880 times.
|
2881 | IF( itau. EQ. itaufinp1 ) then |
| 652 | if (flag_verif) then | ||
| 653 | write(79,*) 'ucov',ucov | ||
| 654 | write(80,*) 'vcov',vcov | ||
| 655 | write(81,*) 'teta',teta | ||
| 656 | write(82,*) 'ps',ps | ||
| 657 | write(83,*) 'q',q | ||
| 658 | WRITE(85,*) 'q1 = ',q(:,:,1) | ||
| 659 | WRITE(86,*) 'q3 = ',q(:,:,3) | ||
| 660 | endif | ||
| 661 | |||
| 662 | 1 | abort_message = 'Simulation finished' | |
| 663 | |||
| 664 | 1 | call abort_gcm(modname,abort_message,0) | |
| 665 | ENDIF | ||
| 666 | c----------------------------------------------------------------------- | ||
| 667 | c ecriture du fichier histoire moyenne: | ||
| 668 | c ------------------------------------- | ||
| 669 | |||
| 670 |
3/4✓ Branch 0 taken 2400 times.
✓ Branch 1 taken 480 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2400 times.
|
2880 | IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN |
| 671 | IF(itau.EQ.itaufin) THEN | ||
| 672 | iav=1 | ||
| 673 | ELSE | ||
| 674 | iav=0 | ||
| 675 | ENDIF | ||
| 676 | |||
| 677 | ! ! Ehouarn: re-compute geopotential for outputs | ||
| 678 | 480 | CALL geopot(ip1jmp1,teta,pk,pks,phis,phi) | |
| 679 | |||
| 680 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
|
480 | IF (ok_dynzon) THEN |
| 681 | CALL bilan_dyn(2,dtvr*iperiod,dtvr*day_step*periodav, | ||
| 682 | ✗ | & ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q) | |
| 683 | END IF | ||
| 684 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
|
480 | IF (ok_dyn_ave) THEN |
| 685 | CALL writedynav(itau,vcov, | ||
| 686 | ✗ | & ucov,teta,pk,phi,q,masse,ps,phis) | |
| 687 | ENDIF | ||
| 688 | |||
| 689 | ENDIF ! of IF((MOD(itau,iperiod).EQ.0).OR.(itau.EQ.itaufin)) | ||
| 690 | |||
| 691 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2880 times.
|
2880 | if (ok_iso_verif) then |
| 692 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 1584') | |
| 693 | endif !if (ok_iso_verif) then | ||
| 694 | |||
| 695 | c----------------------------------------------------------------------- | ||
| 696 | c ecriture de la bande histoire: | ||
| 697 | c ------------------------------ | ||
| 698 | |||
| 699 |
1/2✓ Branch 0 taken 2880 times.
✗ Branch 1 not taken.
|
2880 | IF( MOD(itau,iecri).EQ.0) THEN |
| 700 | ! Ehouarn: output only during LF or Backward Matsuno | ||
| 701 |
4/4✓ Branch 0 taken 960 times.
✓ Branch 1 taken 1920 times.
✓ Branch 2 taken 480 times.
✓ Branch 3 taken 480 times.
|
2880 | if (leapf.or.(.not.leapf.and.(.not.forward))) then |
| 702 | 2400 | CALL geopot(ip1jmp1,teta,pk,pks,phis,phi) | |
| 703 | 2400 | unat=0. | |
| 704 |
2/2✓ Branch 0 taken 93600 times.
✓ Branch 1 taken 2400 times.
|
96000 | do l=1,llm |
| 705 |
2/2✓ Branch 0 taken 95752800 times.
✓ Branch 1 taken 93600 times.
|
95846400 | unat(iip2:ip1jm,l)=ucov(iip2:ip1jm,l)/cu(iip2:ip1jm) |
| 706 |
2/2✓ Branch 0 taken 98841600 times.
✓ Branch 1 taken 93600 times.
|
98937600 | vnat(:,l)=vcov(:,l)/cv(:) |
| 707 | enddo | ||
| 708 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2400 times.
|
2400 | if (ok_dyn_ins) then |
| 709 | ! write(lunout,*) "leapfrog: call writehist, itau=",itau | ||
| 710 | ✗ | CALL writehist(itau,vcov,ucov,teta,phi,q,masse,ps,phis) | |
| 711 | ! call WriteField('ucov',reshape(ucov,(/iip1,jmp1,llm/))) | ||
| 712 | ! call WriteField('vcov',reshape(vcov,(/iip1,jjm,llm/))) | ||
| 713 | ! call WriteField('teta',reshape(teta,(/iip1,jmp1,llm/))) | ||
| 714 | ! call WriteField('ps',reshape(ps,(/iip1,jmp1/))) | ||
| 715 | ! call WriteField('masse',reshape(masse,(/iip1,jmp1,llm/))) | ||
| 716 | endif ! of if (ok_dyn_ins) | ||
| 717 | ! For some Grads outputs of fields | ||
| 718 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2400 times.
|
2400 | if (output_grads_dyn) then |
| 719 | ! | ||
| 720 | ! $Header$ | ||
| 721 | ! | ||
| 722 | ✗ | if (callinigrads) then | |
| 723 | |||
| 724 | ✗ | string10='dyn' | |
| 725 | call inigrads(1,iip1 | ||
| 726 | s ,rlonv,180./pi,-180.,180.,jjp1,rlatu,-90.,90.,180./pi | ||
| 727 | s ,llm,presnivs,1. | ||
| 728 | ✗ | s ,dtvr*iperiod,string10,'dyn_zon ') | |
| 729 | |||
| 730 | ✗ | callinigrads=.false. | |
| 731 | |||
| 732 | |||
| 733 | endif | ||
| 734 | |||
| 735 | ✗ | string10='ps' | |
| 736 | ✗ | CALL wrgrads(1,1,ps,string10,string10) | |
| 737 | |||
| 738 | ✗ | string10='u' | |
| 739 | ✗ | CALL wrgrads(1,llm,unat,string10,string10) | |
| 740 | ✗ | string10='v' | |
| 741 | ✗ | CALL wrgrads(1,llm,vnat,string10,string10) | |
| 742 | ✗ | string10='teta' | |
| 743 | ✗ | CALL wrgrads(1,llm,teta,string10,string10) | |
| 744 | ✗ | do iq=1,nqtot | |
| 745 | ✗ | string10='q' | |
| 746 | ✗ | write(string10(2:2),'(i1)') iq | |
| 747 | ✗ | CALL wrgrads(1,llm,q(:,:,iq),string10,string10) | |
| 748 | enddo | ||
| 749 | |||
| 750 | endif | ||
| 751 | endif ! of if (leapf.or.(.not.leapf.and.(.not.forward))) | ||
| 752 | ENDIF ! of IF(MOD(itau,iecri).EQ.0) | ||
| 753 | |||
| 754 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 2879 times.
|
2880 | IF(itau.EQ.itaufin) THEN |
| 755 | |||
| 756 | |||
| 757 | ! if (planet_type.eq."earth") then | ||
| 758 | ! Write an Earth-format restart file | ||
| 759 | CALL dynredem1("restart.nc",start_time, | ||
| 760 | 1 | & vcov,ucov,teta,q,masse,ps) | |
| 761 | ! endif ! of if (planet_type.eq."earth") | ||
| 762 | |||
| 763 | 1 | CLOSE(99) | |
| 764 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
|
1 | if (ok_guide) then |
| 765 | ! set ok_guide to false to avoid extra output | ||
| 766 | ! in following forward step | ||
| 767 | ✗ | ok_guide=.false. | |
| 768 | endif | ||
| 769 | !!! Ehouarn: Why not stop here and now? | ||
| 770 | ENDIF ! of IF (itau.EQ.itaufin) | ||
| 771 | |||
| 772 | c----------------------------------------------------------------------- | ||
| 773 | c gestion de l'integration temporelle: | ||
| 774 | c ------------------------------------ | ||
| 775 | |||
| 776 |
2/2✓ Branch 0 taken 480 times.
✓ Branch 1 taken 2400 times.
|
2880 | IF( MOD(itau,iperiod).EQ.0 ) THEN |
| 777 | GO TO 1 | ||
| 778 |
2/2✓ Branch 0 taken 960 times.
✓ Branch 1 taken 1440 times.
|
2400 | ELSE IF ( MOD(itau-1,iperiod). EQ. 0 ) THEN |
| 779 | |||
| 780 |
2/2✓ Branch 0 taken 480 times.
✓ Branch 1 taken 480 times.
|
960 | IF( forward ) THEN |
| 781 | c fin du pas forward et debut du pas backward | ||
| 782 | |||
| 783 | 480 | forward = .FALSE. | |
| 784 | 480 | leapf = .FALSE. | |
| 785 | 480 | GO TO 2 | |
| 786 | |||
| 787 | ELSE | ||
| 788 | c fin du pas backward et debut du premier pas leapfrog | ||
| 789 | |||
| 790 | 480 | leapf = .TRUE. | |
| 791 | 480 | dt = 2.*dtvr | |
| 792 | 480 | GO TO 2 | |
| 793 | END IF ! of IF (forward) | ||
| 794 | ELSE | ||
| 795 | |||
| 796 | c ...... pas leapfrog ..... | ||
| 797 | |||
| 798 | 1440 | leapf = .TRUE. | |
| 799 | 1440 | dt = 2.*dtvr | |
| 800 | 1440 | GO TO 2 | |
| 801 | END IF ! of IF (MOD(itau,iperiod).EQ.0) | ||
| 802 | ! ELSEIF (MOD(itau-1,iperiod).EQ.0) | ||
| 803 | |||
| 804 | ELSE ! of IF (.not.purmats) | ||
| 805 | |||
| 806 | ✗ | if (ok_iso_verif) then | |
| 807 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 1664') | |
| 808 | endif !if (ok_iso_verif) then | ||
| 809 | |||
| 810 | c ........................................................ | ||
| 811 | c .............. schema matsuno ............... | ||
| 812 | c ........................................................ | ||
| 813 | ✗ | IF( forward ) THEN | |
| 814 | |||
| 815 | ✗ | itau = itau + 1 | |
| 816 | c iday = day_ini+itau/day_step | ||
| 817 | c time = REAL(itau-(iday-day_ini)*day_step)/day_step+time_0 | ||
| 818 | c | ||
| 819 | c IF(time.GT.1.) THEN | ||
| 820 | c time = time-1. | ||
| 821 | c iday = iday+1 | ||
| 822 | c ENDIF | ||
| 823 | |||
| 824 | ✗ | forward = .FALSE. | |
| 825 | ✗ | IF( itau. EQ. itaufinp1 ) then | |
| 826 | ✗ | abort_message = 'Simulation finished' | |
| 827 | ✗ | call abort_gcm(modname,abort_message,0) | |
| 828 | ENDIF | ||
| 829 | GO TO 2 | ||
| 830 | |||
| 831 | ELSE ! of IF(forward) i.e. backward step | ||
| 832 | |||
| 833 | ✗ | if (ok_iso_verif) then | |
| 834 | ✗ | call check_isotopes_seq(q,ip1jmp1,'leapfrog 1698') | |
| 835 | endif !if (ok_iso_verif) then | ||
| 836 | |||
| 837 | ✗ | IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN | |
| 838 | IF(itau.EQ.itaufin) THEN | ||
| 839 | iav=1 | ||
| 840 | ELSE | ||
| 841 | iav=0 | ||
| 842 | ENDIF | ||
| 843 | |||
| 844 | ! ! Ehouarn: re-compute geopotential for outputs | ||
| 845 | ✗ | CALL geopot(ip1jmp1,teta,pk,pks,phis,phi) | |
| 846 | |||
| 847 | ✗ | IF (ok_dynzon) THEN | |
| 848 | CALL bilan_dyn(2,dtvr*iperiod,dtvr*day_step*periodav, | ||
| 849 | ✗ | & ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q) | |
| 850 | ENDIF | ||
| 851 | ✗ | IF (ok_dyn_ave) THEN | |
| 852 | CALL writedynav(itau,vcov, | ||
| 853 | ✗ | & ucov,teta,pk,phi,q,masse,ps,phis) | |
| 854 | ENDIF | ||
| 855 | |||
| 856 | ENDIF ! of IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) | ||
| 857 | |||
| 858 | ✗ | IF(MOD(itau,iecri ).EQ.0) THEN | |
| 859 | c IF(MOD(itau,iecri*day_step).EQ.0) THEN | ||
| 860 | ✗ | CALL geopot(ip1jmp1,teta,pk,pks,phis,phi) | |
| 861 | ✗ | unat=0. | |
| 862 | ✗ | do l=1,llm | |
| 863 | ✗ | unat(iip2:ip1jm,l)=ucov(iip2:ip1jm,l)/cu(iip2:ip1jm) | |
| 864 | ✗ | vnat(:,l)=vcov(:,l)/cv(:) | |
| 865 | enddo | ||
| 866 | ✗ | if (ok_dyn_ins) then | |
| 867 | ! write(lunout,*) "leapfrog: call writehist (b)", | ||
| 868 | ! & itau,iecri | ||
| 869 | ✗ | CALL writehist(itau,vcov,ucov,teta,phi,q,masse,ps,phis) | |
| 870 | endif ! of if (ok_dyn_ins) | ||
| 871 | ! For some Grads outputs | ||
| 872 | ✗ | if (output_grads_dyn) then | |
| 873 | ! | ||
| 874 | ! $Header$ | ||
| 875 | ! | ||
| 876 | ✗ | if (callinigrads) then | |
| 877 | |||
| 878 | ✗ | string10='dyn' | |
| 879 | call inigrads(1,iip1 | ||
| 880 | s ,rlonv,180./pi,-180.,180.,jjp1,rlatu,-90.,90.,180./pi | ||
| 881 | s ,llm,presnivs,1. | ||
| 882 | ✗ | s ,dtvr*iperiod,string10,'dyn_zon ') | |
| 883 | |||
| 884 | ✗ | callinigrads=.false. | |
| 885 | |||
| 886 | |||
| 887 | endif | ||
| 888 | |||
| 889 | ✗ | string10='ps' | |
| 890 | ✗ | CALL wrgrads(1,1,ps,string10,string10) | |
| 891 | |||
| 892 | ✗ | string10='u' | |
| 893 | ✗ | CALL wrgrads(1,llm,unat,string10,string10) | |
| 894 | ✗ | string10='v' | |
| 895 | ✗ | CALL wrgrads(1,llm,vnat,string10,string10) | |
| 896 | ✗ | string10='teta' | |
| 897 | ✗ | CALL wrgrads(1,llm,teta,string10,string10) | |
| 898 | ✗ | do iq=1,nqtot | |
| 899 | ✗ | string10='q' | |
| 900 | ✗ | write(string10(2:2),'(i1)') iq | |
| 901 | ✗ | CALL wrgrads(1,llm,q(:,:,iq),string10,string10) | |
| 902 | enddo | ||
| 903 | |||
| 904 | endif | ||
| 905 | |||
| 906 | ENDIF ! of IF(MOD(itau,iecri ).EQ.0) | ||
| 907 | |||
| 908 | ✗ | IF(itau.EQ.itaufin) THEN | |
| 909 | ! if (planet_type.eq."earth") then | ||
| 910 | CALL dynredem1("restart.nc",start_time, | ||
| 911 | ✗ | & vcov,ucov,teta,q,masse,ps) | |
| 912 | ! endif ! of if (planet_type.eq."earth") | ||
| 913 | ✗ | if (ok_guide) then | |
| 914 | ! set ok_guide to false to avoid extra output | ||
| 915 | ! in following forward step | ||
| 916 | ✗ | ok_guide=.false. | |
| 917 | endif | ||
| 918 | ENDIF ! of IF(itau.EQ.itaufin) | ||
| 919 | |||
| 920 | ✗ | forward = .TRUE. | |
| 921 | ✗ | GO TO 1 | |
| 922 | |||
| 923 | ENDIF ! of IF (forward) | ||
| 924 | |||
| 925 | END IF ! of IF(.not.purmats) | ||
| 926 | |||
| 927 | END | ||
| 928 |