- Timestamp:
- Jul 29, 2014, 12:19:12 PM (10 years ago)
- Location:
- trunk/WRFV3/lmdz
- Files:
-
- 8 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/WRFV3/lmdz/calltherm.F90
r1 r142 126 126 endif 127 127 128 ! L. Fita, LMD July 2014. Initializing tendencies and variables 129 d_t_the = 0. 130 d_q_the = 0. 131 d_u_the = 0. 132 d_v_the = 0. 133 Alp = 0. 134 Ale = 0. 135 128 136 ! Incrementer le compteur de la physique 129 137 itap = itap + 1 … … 132 140 ! =================== 133 141 ! print*,'thermiques: WARNING on passe t au lieu de t_seri' 134 135 142 136 143 ! On prend comme valeur initiale des thermiques la valeur du pas … … 237 244 238 245 ! print*,'THERM iflag_thermas_ed=',iflag_thermals_ed 246 239 247 CALL thermcell_main(itap,klon,klev,zdt & 240 248 & ,pplay,paprs,pphi,debut & … … 270 278 ! Ce serait bien de changer, mai en prenant le temps de vérifier que ca 271 279 ! fait bien ce qu'on croit. 272 273 280 flag_bidouille_stratocu=iflag_thermals<=12.or.iflag_thermals==14.or.iflag_thermals==16.or.iflag_thermals==18 274 281 … … 300 307 ENDDO 301 308 fm_therm(:,klev+1)=0. 302 303 304 309 305 310 ! accumulation de la tendance -
trunk/WRFV3/lmdz/diagphy_mod.F90
r31 r142 79 79 !C 80 80 integer pas 81 ! L. Fita, LMD July 201482 CHARACTER(LEN=50) :: errmsg, fname83 84 fname = 'diagphy'85 errmsg = 'ERROR -- error -- ERROR -- error'86 81 87 82 save pas … … 89 84 !$OMP THREADPRIVATE(pas) 90 85 !C 86 87 ! L. Fita, LMD July 2014 88 CHARACTER(LEN=50) :: errmsg, fname, varname 89 LOGICAL :: found 90 REAL :: largest 91 92 fname = 'diagphy' 93 errmsg = 'ERROR -- error -- ERROR -- error' 94 largest = 10.e4 95 91 96 pas=pas+1 92 97 stops=0. … … 153 158 154 159 ! L. Fita, LMD July 2014. Checking for consistency 155 IF (fs_bound .NE. fs_bound .OR. ABS(fs_bound) > 100000.) THEN156 PRINT *,TRIM(errmsg) 157 PRINT *,' ' + TRIM(fname) +': Wrong fs_bound= ',fs_bound,' !!!'160 IF (fs_bound .NE. fs_bound .OR. ABS(fs_bound) > largest) THEN 161 PRINT *,TRIM(errmsg) 162 PRINT *,' ' // TRIM(fname) // ': Wrong fs_bound= ',fs_bound,' !!!' 158 163 PRINT *,' fs_bound: Heat flux at atm. boundaries' 159 164 PRINT *,' fs_bound = stops-stopl - (ssols+ssoll)+ssens+sfront' 160 PRINT *,' stops-stopl= ',stops-stopl 165 PRINT *,' airetot= ',airetot 166 IF (airetot .NE. airetot .OR. ABS(airetot) > largest*10.e6) THEN 167 varname = 'airephy' 168 CALL check_var(fname, varname, airephy, klon, largest*10.e6, .FALSE.) 169 END IF 170 PRINT *,' stops= ',stops 171 IF (stops .NE. stops .OR. ABS(stops) > largest) THEN 172 varname = 'tops' 173 CALL check_var(fname, varname, tops, klon, largest, .FALSE.) 174 END IF 175 PRINT *,' stopl= ',stopl 176 IF (stopl .NE. stopl .OR. ABS(stopl) > largest) THEN 177 varname = 'topl' 178 CALL check_var(fname, varname, topl, klon, largest, .FALSE.) 179 END IF 161 180 PRINT *,' ssols= ',ssols 181 IF (ssols .NE. ssols .OR. ABS(ssols) > largest) THEN 182 varname = 'sols' 183 CALL check_var(fname, varname, sols, klon, largest, .FALSE.) 184 END IF 162 185 PRINT *,' ssoll= ',ssoll 186 IF (ssoll .NE. ssoll .OR. ABS(ssoll) > largest) THEN 187 varname = 'soll' 188 CALL check_var(fname, varname, soll, klon, largest, .FALSE.) 189 END IF 163 190 PRINT *,' ssens= ',ssens 191 IF (ssens .NE. ssens .OR. ABS(ssens) > largest) THEN 192 varname = 'sens' 193 CALL check_var(fname, varname, sens, klon, largest, .FALSE.) 194 END IF 164 195 PRINT *,' sfront= ',sfront 165 STOP 166 END IF 167 168 IF (fq_bound .NE. fq_bound .OR. ABS(fq_bound) > 100000.) THEN 169 PRINT *,TRIM(errmsg) 170 PRINT *,' ' + TRIM(fname) + ': Wrong fs_bound= ',fs_bound,' !!!' 196 IF (sfront .NE. sfront .OR. ABS(sfront) > largest) THEN 197 varname = 'evap' 198 CALL check_var(fname, varname, evap, klon, largest, .FALSE.) 199 varname = 'rain_fall' 200 CALL check_var(fname, varname, rain_fall, klon, largest, .FALSE.) 201 varname = 'snow_fall' 202 CALL check_var(fname, varname, snow_fall, klon, largest, .FALSE.) 203 varname = 'ts' 204 CALL check_var(fname, varname, ts, klon, largest, .FALSE.) 205 END IF 206 STOP 207 END IF 208 209 IF (d_etp_tot .NE. d_etp_tot .OR. ABS(d_etp_tot) > largest) THEN 210 PRINT *,TRIM(errmsg) 211 PRINT *,' ' // TRIM(fname) // ': Wrong d_etp_tot= ',d_etp_tot,' !!!' 212 PRINT *,' d_etp_tot: Heat flux equivalent to atmospheric enthalpy' // & 213 ' change (W/m2)' 214 PRINT *,' d_etp_tot = input value !' 215 STOP 216 END IF 217 218 IF (fq_bound .NE. fq_bound .OR. ABS(fq_bound) > largest) THEN 219 PRINT *,TRIM(errmsg) 220 PRINT *,' ' // TRIM(fname) // ': Wrong fq_bound= ',fs_bound,' !!!' 171 221 PRINT *,' fq_bound: Watter flux at atm. boundaries' 172 222 PRINT *,' fq_bound = evap_tot - rain_fall_tot -snow_fall_tot' 173 223 PRINT *,' evap_tot= ',evap_tot 224 PRINT *,' airetot= ',airetot 225 IF (airetot .NE. airetot .OR. ABS(airetot) > largest*10.e5) THEN 226 varname = 'airephy' 227 CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.) 228 END IF 229 IF (evap_tot .NE. evap_tot .OR. ABS(evap_tot) > largest) THEN 230 varname = 'evap' 231 CALL check_var(fname, varname, evap, klon, largest, .FALSE.) 232 END IF 174 233 PRINT *,' rain_fall_tot= ',rain_fall_tot 234 IF (rain_fall_tot .NE. rain_fall_tot .OR. ABS(rain_fall_tot) > largest) THEN 235 varname = 'rain_fall' 236 CALL check_var(fname, varname, rain_fall, klon, largest, .FALSE.) 237 END IF 175 238 PRINT *,' snow_fall_tot= ',snow_fall_tot 239 IF (snow_fall_tot .NE. snow_fall_tot .OR. ABS(snow_fall_tot) > largest) THEN 240 varname = 'snow_fall' 241 CALL check_var(fname, varname, snow_fall, klon, largest, .FALSE.) 242 END IF 243 STOP 244 END IF 245 246 IF (d_qt_tot .NE. d_qt_tot .OR. ABS(d_qt_tot) > largest) THEN 247 PRINT *,TRIM(errmsg) 248 PRINT *,' ' // TRIM(fname) // ': Wrong d_qt_tot= ',d_qt_tot,' !!!' 249 PRINT *,' d_qt_tot: Mass flux equivalent to atmospheric watter mass' // & 250 ' change (kg/m2/s)' 251 PRINT *,' d_qt_tot = input value !' 176 252 STOP 177 253 END IF … … 184 260 185 261 end SUBROUTINE diagphy 262 263 !L. Fita, LMD 2004. Check variable 264 SUBROUTINE check_var(funcn, varn, var, sizev, bigvalue, stoprun) 265 ! Subroutine to check the consistency of a variable 266 ! * NaN value: by definition is variable /= variable 267 ! * bigvalue: threshold for the variable 268 269 IMPLICIT NONE 270 271 #include "dimensions.h" 272 273 INTEGER, INTENT(IN) :: sizev 274 CHARACTER(LEN=50), INTENT(IN) :: funcn, varn 275 REAL, DIMENSION(sizev), INTENT(IN) :: var 276 REAL, INTENT(IN) :: bigvalue 277 LOGICAL, INTENT(IN) :: stoprun 278 279 ! Local 280 INTEGER :: i, wrongi, xpt, ypt 281 CHARACTER(LEN=50) :: errmsg 282 LOGICAL :: found 283 REAL, DIMENSION(sizev) :: wrongvalues 284 INTEGER, DIMENSION(sizev) :: wronggridpt 285 286 !!!!!!! Variables 287 ! funcn: at which functino of part of the program variable is checked 288 ! varn: name of the variable 289 ! var: variable to check 290 ! sizev: size of the variable 291 ! bigvalue: biggest attenaible value for the variable 292 ! stoprun: Should the run stop if it founds a problem? 293 294 errmsg = 'ERROR -- error -- ERROR -- error' 295 296 found = .FALSE. 297 wrongi = 0 298 DO i=1,sizev 299 IF (var(i) /= var(i) .OR. ABS(var(i)) > bigvalue ) THEN 300 IF (wrongi == 0) found = .TRUE. 301 wrongi = wrongi + 1 302 wrongvalues(wrongi) = var(i) 303 wronggridpt(wrongi) = i 304 END IF 305 END DO 306 307 IF (found) THEN 308 PRINT *,TRIM(errmsg) 309 PRINT *," at '" // TRIM(funcn) // "' variable '" //TRIM(varn)// & 310 "' is wrong in Nvalues= ",wrongi,' at i (x, y) value___' 311 DO i=1,wrongi 312 ypt = INT(wronggridpt(i)/wiim) + 1 313 xpt = wronggridpt(i) - (ypt-1)*wiim 314 PRINT *,wronggridpt(i), '(',xpt,', ',ypt,')', wrongvalues(i) 315 END DO 316 IF (stoprun) THEN 317 STOP 318 END IF 319 END IF 320 321 RETURN 322 323 END SUBROUTINE check_var 324 325 SUBROUTINE check_var3D(funcn, varn, var, sizev, zsize, bigvalue, stoprun) 326 ! Subroutine to check the consistency of a 3D LMDSZ - variable (klon, klev) ! 327 ! * NaN value: by definition is variable /= variable 328 ! * bigvalue: threshold for the variable 329 330 IMPLICIT NONE 331 332 #include "dimensions.h" 333 334 INTEGER, INTENT(IN) :: sizev, zsize 335 CHARACTER(LEN=50), INTENT(IN) :: funcn, varn 336 REAL, DIMENSION(sizev,zsize), INTENT(IN) :: var 337 REAL, INTENT(IN) :: bigvalue 338 LOGICAL, INTENT(IN) :: stoprun 339 340 ! Local 341 INTEGER :: i, k, wrongi, xpt, ypt 342 CHARACTER(LEN=50) :: errmsg 343 LOGICAL :: found 344 REAL, DIMENSION(sizev*zsize) :: wrongvalues 345 INTEGER, DIMENSION(sizev*zsize,2) :: wronggridpt 346 347 !!!!!!! Variables 348 ! funcn: at which functino of part of the program variable is checked 349 ! varn: name of the variable 350 ! var: variable to check 351 ! sizev: size of the variable 352 ! zsize: vertical size of the variable 353 ! bigvalue: biggest attenaible value for the variable 354 ! stoprun: Should the run stop if it founds a problem? 355 356 errmsg = 'ERROR -- error -- ERROR -- error' 357 358 found = .FALSE. 359 wrongi = 0 360 DO i=1,sizev 361 DO k=1,zsize 362 IF (var(i,k) /= var(i,k) .OR. ABS(var(i,k)) > bigvalue ) THEN 363 IF (wrongi == 0) found = .TRUE. 364 wrongi = wrongi + 1 365 wrongvalues(wrongi) = var(i,k) 366 wronggridpt(wrongi,1) = i 367 wronggridpt(wrongi,2) = k 368 END IF 369 END DO 370 END DO 371 372 IF (found) THEN 373 PRINT *,TRIM(errmsg) 374 PRINT *," at '" // TRIM(funcn) // "' variable '" //TRIM(varn)// & 375 "' is wrong in Nvalues= ",wrongi,' at i (x,y) k value___' 376 DO i=1,wrongi 377 ypt = INT(wronggridpt(i,1)/wiim) + 1 378 xpt = wronggridpt(i,1) - (ypt-1)*wiim 379 PRINT *,wronggridpt(i,1), '(',xpt,', ',ypt,')', wronggridpt(i,2), wrongvalues(i) 380 END DO 381 IF (stoprun) THEN 382 STOP 383 END IF 384 END IF 385 386 RETURN 387 388 END SUBROUTINE check_var3D 389 186 390 187 391 !C====================================================================== … … 307 511 !$OMP THREADPRIVATE(h_vcol_pre, h_dair_pre, h_qw_pre, h_ql_pre) 308 512 !$OMP THREADPRIVATE(h_qs_pre, qw_pre, ql_pre, qs_pre , ec_pre) 513 514 ! L. Fita, LMD July 2014 515 CHARACTER(LEN=50) :: errmsg, fname, varname 516 LOGICAL :: found 517 REAL :: largest 518 519 fname = 'diagetpq' 520 errmsg = 'ERROR -- error -- ERROR -- error' 521 largest = 10.e4 309 522 !c====================================================================== 310 523 !C … … 449 662 ec_pre (idiag) = ec_tot 450 663 !C 664 IF (d_h_vcol .NE. d_h_vcol .OR. ABS(d_h_vcol) > largest) THEN 665 PRINT *,TRIM(errmsg) 666 PRINT *,' ' // TRIM(fname) // ': Wrong d_h_vcol= ',d_h_vcol,' !!!' 667 PRINT *,' d_h_vcol: Heat flux (W/m2) define as the Enthalpy change' // & 668 ' (J/m2) during one time step (dtime) for the whole atmosphere (air,' // & 669 ' watter vapour, liquid and solid)' 670 PRINT *,' d_h_vcol = (h_vcol_tot - h_vcol_pre(idiag2) )' 671 PRINT *,' h_vcol_tot= ',h_vcol_tot 672 IF (h_vcol_tot .NE. h_vcol_tot .OR. ABS(h_vcol_tot) > largest) THEN 673 PRINT *,' h_vcol_tot = h_dair_tot+h_qw_tot+h_ql_tot+h_qs_tot' 674 PRINT *,' airetot= ',airetot 675 IF (airetot .NE. airetot .OR. ABS(airetot) > largest*10.e5) THEN 676 varname = 'airephy' 677 CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.) 678 END IF 679 PRINT *,' h_dair_tot= ',h_dair_tot 680 IF (h_dair_tot .NE. h_dair_tot .OR. ABS(h_dair_tot) > largest*10.e6) THEN 681 varname = 'zh_dair_col' 682 CALL check_var(fname, varname, zh_dair_col, klon, largest*10.e6, .FALSE.) 683 PRINT *,' zh_dair_col = func(q,ql,qs,paprs,t)' 684 varname = 'q' 685 CALL check_var3D(fname, varname, q, klon, klev, largest, .FALSE.) 686 varname = 'ql' 687 CALL check_var3D(fname, varname, ql, klon, klev, largest, .FALSE.) 688 varname = 'qs' 689 CALL check_var3D(fname, varname, qs, klon, klev, largest, .FALSE.) 690 varname = 'paprs' 691 CALL check_var3D(fname, varname, paprs, klon, klev, largest*10.e6, .FALSE.) 692 varname = 't' 693 CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.) 694 END IF 695 PRINT *,' h_qw_tot= ',h_qw_tot 696 IF (h_qw_tot .NE. h_qw_tot .OR. ABS(h_qw_tot) > largest*10.e4) THEN 697 varname = 'zh_qw_col' 698 CALL check_var(fname, varname, zh_qw_col, klon, largest*10.e4, .FALSE.) 699 END IF 700 PRINT *,' h_ql_tot= ',h_ql_tot 701 IF (h_ql_tot .NE. h_ql_tot .OR. ABS(h_ql_tot) > largest*10.e2) THEN 702 varname = 'zh_ql_col' 703 CALL check_var(fname, varname, zh_ql_col, klon, largest*10.e2, .FALSE.) 704 END IF 705 PRINT *,' h_qs_tot= ',h_qs_tot 706 IF (h_qs_tot .NE. h_qs_tot .OR. ABS(h_qs_tot) > largest) THEN 707 varname = 'zh_qs_col' 708 CALL check_var(fname, varname, zh_qs_col, klon, largest, .FALSE.) 709 END IF 710 END IF 711 STOP 712 END IF 713 714 IF (qw_tot .NE. qw_tot .OR. ABS(qw_tot) > largest) THEN 715 PRINT *,TRIM(errmsg) 716 PRINT *,' ' // TRIM(fname) // ': Wrong qw_tot= ',qw_tot,' !!!' 717 PRINT *,' qw_tot: total mass of watter vapour (kg/m2)' 718 PRINT *,' qw_tot = f(q,ql,qs,airephy,t,paprs)' 719 varname = 'airephy' 720 CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.) 721 varname = 'q' 722 CALL check_var3D(fname, varname, q, klon, klev, largest, .FALSE.) 723 varname = 'ql' 724 CALL check_var3D(fname, varname, ql, klon, klev, largest, .FALSE.) 725 varname = 'qs' 726 CALL check_var3D(fname, varname, qs, klon, klev, largest, .FALSE.) 727 varname = 't' 728 CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.) 729 varname = 'paprs' 730 CALL check_var3D(fname, varname, paprs, klon, klev, largest, .FALSE.) 731 STOP 732 END IF 733 734 IF (ql_tot .NE. ql_tot .OR. ABS(ql_tot) > largest) THEN 735 PRINT *,TRIM(errmsg) 736 PRINT *,' ' // TRIM(fname) // ': Wrong ql_tot= ',ql_tot,' !!!' 737 PRINT *,' ql_tot: total mass of liquid watter (kg/m2)' 738 PRINT *,' ql_tot = f(ql,airephy,t,paprs)' 739 varname = 'airephy' 740 CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.) 741 varname = 'ql' 742 CALL check_var3D(fname, varname, ql, klon, klev, largest, .FALSE.) 743 varname = 't' 744 CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.) 745 varname = 'paprs' 746 CALL check_var3D(fname, varname, paprs, klon, klev, largest, .FALSE.) 747 STOP 748 END IF 749 750 IF (qs_tot .NE. qs_tot .OR. ABS(qs_tot) > largest) THEN 751 PRINT *,TRIM(errmsg) 752 PRINT *,' ' // TRIM(fname) // ': Wrong qs_tot= ',qs_tot,' !!!' 753 PRINT *,' qs_tot: total mass of solid watter (kg/m2)' 754 PRINT *,' qs_tot = f(qs,airephy,t,paprs)' 755 varname = 'airephy' 756 CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.) 757 varname = 'qs' 758 CALL check_var3D(fname, varname, qs, klon, klev, largest, .FALSE.) 759 varname = 't' 760 CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.) 761 varname = 'paprs' 762 CALL check_var3D(fname, varname, paprs, klon, klev, largest, .FALSE.) 763 STOP 764 END IF 765 766 IF (ec_tot .NE. ec_tot .OR. ABS(ec_tot) > largest*10.e4) THEN 767 PRINT *,TRIM(errmsg) 768 PRINT *,' ' // TRIM(fname) // ': Wrong ec_tot= ',ec_tot,' !!!' 769 PRINT *,' ec_tot: total cinetic energy (kg/m2)' 770 PRINT *,' ec_tot = f(u,v,airephy,t,paprs)' 771 varname = 'airephy' 772 CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.) 773 varname = 'u' 774 CALL check_var3D(fname, varname, u, klon, klev, largest, .FALSE.) 775 varname = 'v' 776 CALL check_var3D(fname, varname, v, klon, klev, largest, .FALSE.) 777 varname = 't' 778 CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.) 779 varname = 'paprs' 780 CALL check_var3D(fname, varname, paprs, klon, klev, largest, .FALSE.) 781 STOP 782 END IF 783 451 784 RETURN 785 452 786 END SUBROUTINE diagetpq 453 787 -
trunk/WRFV3/lmdz/physiq.F90
r6 r142 1299 1299 integer iostat 1300 1300 1301 ! L luis1301 ! L. Fita, LMD. Point for checkings... 1302 1302 INTEGER :: llp 1303 1303 llp = 644 … … 1307 1307 ! 1308 1308 CALL phys_cal_update(jD_cur,jH_cur) 1309 PRINT *,' Lluis physiq: jD_cur: ',jD_cur, ' jH_cur: ',jH_cur, &1310 ' days_elapsed: ',days_elapsed1311 1309 1312 1310 !c====================================================================== … … 1962 1960 & , fs_bound, fq_bound ) 1963 1961 END IF 1962 PRINT *,' Lluis 1 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 1964 1963 1965 1964 !c Diagnostiquer la tendance dynamique … … 2122 2121 !C 2123 2122 END IF 2124 2123 PRINT *,' Lluis 2 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 2125 2124 !c 2126 2125 !c========================================================================= … … 2168 2167 zdtime=dtime*REAL(radpas) ! pas de temps du rayonnement (s) 2169 2168 CALL zenang(zlongi,jH_cur,zdtime,rlat,rlon,rmu0,fract) 2170 PRINT *,' Lluis: radpas: ',radpas,' fract: ',fract(llp), &2171 ' zdtime: ',zdtime,' jH_cur: ',jH_cur,' rmu0: ',rmu0(llp),&2172 ' dtime: ',dtime,' nbapp_rad: ',nbapp_rad2173 PRINT *,' Lluis @rlat@: ',rlat2174 PRINT *,' Lluis @rlon@: ',rlon2175 PRINT *,' Lluis @rmu0@: ',rmu02176 PRINT *,' Lluis @fract@: ',fract2177 2169 ELSE 2178 2170 CALL angle(zlongi, rlat, fract, rmu0) 2179 2171 ENDIF 2180 2172 ENDIF 2181 PRINT *,' Lluis zenang_an rlat: ',rlat(llp), ' rlon:',rlon(llp), &2182 ' rmu0: ',rmu0(llp),' fract: ',fract(llp),' jH_cur: ',jH_cur, &2183 ' cycle_diurne: ',cycle_diurne2184 2173 2185 2174 if (mydebug) then … … 2275 2264 2276 2265 ENDIF 2266 PRINT *,' Lluis 3 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 2277 2267 !c =================================================================== c 2278 2268 !c Calcul de Qsat … … 2635 2625 END IF 2636 2626 !C 2627 PRINT *,' Lluis 4 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 2637 2628 IF (check) THEN 2638 2629 za = qcheck(klon,klev,paprs,q_seri,ql_seri,airephy) … … 2747 2738 & , fs_bound, fq_bound ) 2748 2739 END IF 2749 2740 PRINT *,' Lluis 5 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 2750 2741 !c print*,'apres callwake iflag_cldcon=', iflag_cldcon 2751 2742 !c … … 2771 2762 !c detr_therm(:,:)=0. 2772 2763 !c 2764 2773 2765 IF(prt_level>9)WRITE(lunout,*) & 2774 2766 & 'AVANT LA CONVECTION SECHE , iflag_thermals=' & … … 2958 2950 2959 2951 endif 2952 2960 2953 !c 2961 2954 !c=================================================================== … … 2972 2965 & , fs_bound, fq_bound ) 2973 2966 END IF 2974 2967 PRINT *,' Lluis 6 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 2975 2968 2976 2969 !c------------------------------------------------------------------------- … … 3039 3032 & , fs_bound, fq_bound ) 3040 3033 END IF 3034 PRINT *,' Lluis 7 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 3041 3035 if (mydebug) then 3042 3036 call writefield_phy('u_seri',u_seri,llm) … … 3045 3039 call writefield_phy('q_seri',q_seri,llm) 3046 3040 endif 3047 3048 3041 !c 3049 3042 !c------------------------------------------------------------------- … … 3271 3264 & , fs_bound, fq_bound ) 3272 3265 END IF 3266 PRINT *,' Lluis 8 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 3273 3267 !c 3274 3268 !c Calculer l'humidite relative pour diagnostique … … 3643 3637 & , fs_bound, fq_bound ) 3644 3638 END IF 3639 PRINT *,' Lluis 9 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 3645 3640 !c 3646 3641 !c … … 3817 3812 & , fs_bound, fq_bound ) 3818 3813 END IF 3814 PRINT *,' Lluis 10 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 3819 3815 !c 3820 3816 !c … … 3953 3949 !C Donc la somme de ces 2 variations devrait etre nulle. 3954 3950 3951 PRINT *,' LLuis sending: d_h_vcol: ',d_h_vcol,' d_qt ',d_qt,' d_ec ',d_ec 3952 3955 3953 call diagphy(airephy,ztit,ip_ebil_phy & 3956 3954 & , topsw, toplw, solsw, sollw, sens & … … 3962 3960 !C 3963 3961 END IF 3962 PRINT *,' Lluis 11 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy 3964 3963 PRINT *,'Lluis Reaching the SORTIES point' 3965 3964 -
trunk/WRFV3/lmdz/thermcell_dq.F90
r1 r142 37 37 CHARACTER (LEN=20) :: modname='thermcell_dq' 38 38 CHARACTER (LEN=80) :: abort_message 39 40 39 41 40 ! Old explicite scheme … … 115 114 fqa(:,1)=0. ; fqa(:,nlay)=0. 116 115 117 118 116 ! Trace species evolution 119 117 if (impl==0) then -
trunk/WRFV3/lmdz/thermcell_flux2.F90
r1 r142 150 150 enddo 151 151 152 153 154 152 ! Test provisoire : pour comprendre pourquoi on corrige plein de fois 155 153 ! le cas fm6, on commence par regarder une premiere fois avant les … … 192 190 endif 193 191 enddo 194 195 192 196 193 !------------------------------------------------------------------------- … … 461 458 enddo 462 459 enddo 460 463 461 if (labort_gcm) then 464 462 ig=igout -
trunk/WRFV3/lmdz/thermcell_height.F90
r1 r142 64 64 enddo 65 65 66 ! L. Fita, LMD July 2014. Fixing the If when zw2 == 0 66 67 do l=1,nlay 67 68 do ig=1,ngrid 68 69 if (l.le.lmax(ig)) then 69 if (zw2(ig,l).lt.0.)then 70 print*,'pb2 zw2<0' 71 endif 70 ! Previous version 71 ! if (zw2(ig,l).lt.0.)then 72 ! print*,'pb2 zw2<0' 73 ! endif 74 ! zw2(ig,l)=sqrt(zw2(ig,l)) 75 ! wmax(ig)=max(wmax(ig),zw2(ig,l)) 76 ! else 77 ! zw2(ig,l)=0. 78 ! endif 79 ! New version 80 if (zw2(ig,l).lt.0.)then 81 print*,'pb2 zw2<0' 82 zw2(ig,l)=0. 83 wmax(ig)=max(wmax(ig),zw2(ig,l)) 84 else 72 85 zw2(ig,l)=sqrt(zw2(ig,l)) 73 86 wmax(ig)=max(wmax(ig),zw2(ig,l)) 87 endif 74 88 else 75 89 zw2(ig,l)=0. 76 90 endif 91 77 92 enddo 78 93 enddo -
trunk/WRFV3/lmdz/thermcell_main.F90
r1 r142 105 105 106 106 INTEGER ig,k,l,ll,ierr 107 real zsortie1d( klon)108 INTEGER lmax( klon),lmin(klon),lalim(klon)109 INTEGER lmix( klon)110 INTEGER lmix_bis( klon)111 real linter( klon)112 real zmix( klon)113 real zmax( klon),zw2(klon,klev+1),ztva(klon,klev),zw_est(klon,klev+1),ztva_est(klon,klev)114 ! real fraca( klon,klev)115 116 real zmax_sec( klon)107 real zsortie1d(ngrid) 108 INTEGER lmax(ngrid),lmin(ngrid),lalim(ngrid) 109 INTEGER lmix(ngrid) 110 INTEGER lmix_bis(ngrid) 111 real linter(ngrid) 112 real zmix(ngrid) 113 real zmax(ngrid),zw2(ngrid,nlay+1),ztva(ngrid,nlay),zw_est(ngrid,nlay+1),ztva_est(ngrid,nlay) 114 ! real fraca(ngrid,nlay) 115 116 real zmax_sec(ngrid) 117 117 !on garde le zmax du pas de temps precedent 118 real zmax0( klon)118 real zmax0(ngrid) 119 119 !FH/IM save zmax0 120 120 121 121 real lambda 122 122 123 real zlev( klon,klev+1),zlay(klon,klev)124 real deltaz( klon,klev)125 REAL zh( klon,klev)126 real zthl( klon,klev),zdthladj(klon,klev)127 REAL ztv( klon,klev)128 real zu( klon,klev),zv(klon,klev),zo(klon,klev)129 real zl( klon,klev)130 real zsortie( klon,klev)131 real zva( klon,klev)132 real zua( klon,klev)133 real zoa( klon,klev)134 135 real zta( klon,klev)136 real zha( klon,klev)137 real fraca( klon,klev+1)123 real zlev(ngrid,nlay+1),zlay(ngrid,nlay) 124 real deltaz(ngrid,nlay) 125 REAL zh(ngrid,nlay) 126 real zthl(ngrid,nlay),zdthladj(ngrid,nlay) 127 REAL ztv(ngrid,nlay) 128 real zu(ngrid,nlay),zv(ngrid,nlay),zo(ngrid,nlay) 129 real zl(ngrid,nlay) 130 real zsortie(ngrid,nlay) 131 real zva(ngrid,nlay) 132 real zua(ngrid,nlay) 133 real zoa(ngrid,nlay) 134 135 real zta(ngrid,nlay) 136 real zha(ngrid,nlay) 137 real fraca(ngrid,nlay+1) 138 138 real zf,zf2 139 real thetath2( klon,klev),wth2(klon,klev),wth3(klon,klev)140 real q2( klon,klev)139 real thetath2(ngrid,nlay),wth2(ngrid,nlay),wth3(ngrid,nlay) 140 real q2(ngrid,nlay) 141 141 ! FH probleme de dimensionnement avec l'allocation dynamique 142 142 ! common/comtherm/thetath2,wth2 143 real wq( klon,klev)144 real wthl( klon,klev)145 real wthv( klon,klev)143 real wq(ngrid,nlay) 144 real wthl(ngrid,nlay) 145 real wthv(ngrid,nlay) 146 146 147 real ratqscth( klon,klev)147 real ratqscth(ngrid,nlay) 148 148 real var 149 149 real vardiff 150 real ratqsdiff( klon,klev)150 real ratqsdiff(ngrid,nlay) 151 151 152 152 logical sorties 153 real rho( klon,klev),rhobarz(klon,klev),masse(klon,klev)154 real zpspsk( klon,klev)155 156 real wmax( klon)157 real wmax_tmp( klon)158 real wmax_sec( klon)159 real fm0( klon,klev+1),entr0(klon,klev),detr0(klon,klev)160 real fm( klon,klev+1),entr(klon,klev),detr(klon,klev)161 162 real ztla( klon,klev),zqla(klon,klev),zqta(klon,klev)153 real rho(ngrid,nlay),rhobarz(ngrid,nlay),masse(ngrid,nlay) 154 real zpspsk(ngrid,nlay) 155 156 real wmax(ngrid) 157 real wmax_tmp(ngrid) 158 real wmax_sec(ngrid) 159 real fm0(ngrid,nlay+1),entr0(ngrid,nlay),detr0(ngrid,nlay) 160 real fm(ngrid,nlay+1),entr(ngrid,nlay),detr(ngrid,nlay) 161 162 real ztla(ngrid,nlay),zqla(ngrid,nlay),zqta(ngrid,nlay) 163 163 !niveau de condensation 164 integer nivcon( klon)165 real zcon( klon)164 integer nivcon(ngrid) 165 real zcon(ngrid) 166 166 REAL CHI 167 real zcon2( klon)168 real pcon( klon)169 real zqsat( klon,klev)170 real zqsatth( klon,klev)171 172 real f_star( klon,klev+1),entr_star(klon,klev)173 real detr_star( klon,klev)174 real alim_star_tot( klon)175 real alim_star( klon,klev)176 real alim_star_clos( klon,klev)177 real f( klon), f0(klon)167 real zcon2(ngrid) 168 real pcon(ngrid) 169 real zqsat(ngrid,nlay) 170 real zqsatth(ngrid,nlay) 171 172 real f_star(ngrid,nlay+1),entr_star(ngrid,nlay) 173 real detr_star(ngrid,nlay) 174 real alim_star_tot(ngrid) 175 real alim_star(ngrid,nlay) 176 real alim_star_clos(ngrid,nlay) 177 real f(ngrid), f0(ngrid) 178 178 !FH/IM save f0 179 real zlevinter( klon)179 real zlevinter(ngrid) 180 180 logical debut 181 181 real seuil 182 real csc( klon,klev)182 real csc(ngrid,nlay) 183 183 184 184 !!! nrlmd le 10/04/2012 185 185 186 186 !------Entrées 187 real pbl_tke( klon,klev+1,nbsrf)188 real pctsrf( klon,nbsrf)189 real omega( klon,klev)190 real airephy( klon)187 real pbl_tke(ngrid,nlay+1,nbsrf) 188 real pctsrf(ngrid,nbsrf) 189 real omega(ngrid,nlay) 190 real airephy(ngrid) 191 191 !------Sorties 192 real zlcl( klon),fraca0(klon),w0(klon),w_conv(klon)193 real therm_tke_max0( klon),env_tke_max0(klon)194 real n2( klon),s2(klon)195 real ale_bl_stat( klon)196 real therm_tke_max( klon,klev),env_tke_max(klon,klev)197 real alp_bl_det( klon),alp_bl_fluct_m(klon),alp_bl_fluct_tke(klon),alp_bl_conv(klon),alp_bl_stat(klon)192 real zlcl(ngrid),fraca0(ngrid),w0(ngrid),w_conv(ngrid) 193 real therm_tke_max0(ngrid),env_tke_max0(ngrid) 194 real n2(ngrid),s2(ngrid) 195 real ale_bl_stat(ngrid) 196 real therm_tke_max(ngrid,nlay),env_tke_max(ngrid,nlay) 197 real alp_bl_det(ngrid),alp_bl_fluct_m(ngrid),alp_bl_fluct_tke(ngrid),alp_bl_conv(ngrid),alp_bl_stat(ngrid) 198 198 !------Local 199 199 integer nsrf 200 real rhobarz0( klon) ! Densité au LCL201 logical ok_lcl( klon) ! Existence du LCL des thermiques202 integer klcl( klon) ! Niveau du LCL203 real interp( klon) ! Coef d'interpolation pour le LCL200 real rhobarz0(ngrid) ! Densité au LCL 201 logical ok_lcl(ngrid) ! Existence du LCL des thermiques 202 integer klcl(ngrid) ! Niveau du LCL 203 real interp(ngrid) ! Coef d'interpolation pour le LCL 204 204 !--Triggering 205 205 real Su ! Surface unité: celle d'un updraft élémentaire … … 215 215 real zmax_moy_coef 216 216 parameter(zmax_moy_coef=0.33) 217 real depth( klon) ! Epaisseur moyenne du cumulus218 real w_max( klon) ! Vitesse max statistique219 real s_max( klon)217 real depth(ngrid) ! Epaisseur moyenne du cumulus 218 real w_max(ngrid) ! Vitesse max statistique 219 real s_max(ngrid) 220 220 !--Closure 221 real pbl_tke_max( klon,klev) ! Profil de TKE moyenne222 real pbl_tke_max0( klon) ! TKE moyenne au LCL223 real w_ls( klon,klev) ! Vitesse verticale grande échelle (m/s)221 real pbl_tke_max(ngrid,nlay) ! Profil de TKE moyenne 222 real pbl_tke_max0(ngrid) ! TKE moyenne au LCL 223 real w_ls(ngrid,nlay) ! Vitesse verticale grande échelle (m/s) 224 224 real coef_m ! On considère un rendement pour alp_bl_fluct_m 225 225 parameter(coef_m=1.) … … 231 231 ! 232 232 !nouvelles variables pour la convection 233 real Ale_bl( klon)234 real Alp_bl( klon)235 real alp_int( klon),dp_int(klon),zdp236 real ale_int( klon)237 integer n_int( klon)238 real fm_tot( klon)239 real wght_th( klon,klev)240 integer lalim_conv( klon)233 real Ale_bl(ngrid) 234 real Alp_bl(ngrid) 235 real alp_int(ngrid),dp_int(ngrid),zdp 236 real ale_int(ngrid) 237 integer n_int(ngrid) 238 real fm_tot(ngrid) 239 real wght_th(ngrid,nlay) 240 integer lalim_conv(ngrid) 241 241 !v1d logical therm 242 242 !v1d save therm … … 250 250 EXTERNAL SCOPY 251 251 ! 252 253 ! L. Fita, LMD July 2014. Initializing variables. 254 ! Some not initializated according to values: iflag_trig_bl, iflag_clos_bl 255 zdthladj = 0. 256 pbl_tke_max0 = 0. 257 fraca0 = 0. 258 w_conv = 0. 259 w0 = 0. 260 therm_tke_max0 = 0. 261 env_tke_max0 = 0. 262 alp_bl_det = 0. 263 alp_bl_fluct_m = 0. 264 alp_bl_fluct_tke = 0. 265 alp_bl_conv = 0. 266 alp_bl_stat = 0. 267 interp = 0. 268 klcl = 0 252 269 253 270 !----------------------------------------------------------------------- … … 306 323 307 324 sorties=.true. 308 IF(ngrid.NE. klon) THEN325 IF(ngrid.NE.ngrid) THEN 309 326 PRINT* 310 327 PRINT*,'STOP dans convadj' 311 328 PRINT*,'ngrid =',ngrid 312 PRINT*,' klon =',klon329 PRINT*,'ngrid =',ngrid 313 330 ENDIF 314 331 ! 315 332 ! write(lunout,*)'WARNING thermcell_main f0=max(f0,1.e-2)' 316 do ig=1, klon333 do ig=1,ngrid 317 334 f0(ig)=max(f0(ig),1.e-2) 318 335 zmax0(ig)=max(zmax0(ig),40.) … … 331 348 CALL thermcell_env(ngrid,nlay,po,pt,pu,pv,pplay, & 332 349 & pplev,zo,zh,zl,ztv,zthl,zu,zv,zpspsk,zqsat,lev_out) 333 350 334 351 if (prt_level.ge.1) print*,'thermcell_main apres thermcell_env' 335 352 … … 360 377 enddo 361 378 zlev(:,1)=0. 362 zlev(:,nlay+1)=(2.*pphi(:, klev)-pphi(:,klev-1))/RG379 zlev(:,nlay+1)=(2.*pphi(:,nlay)-pphi(:,nlay-1))/RG 363 380 do l=1,nlay 364 381 zlay(:,l)=pphi(:,l)/RG … … 516 533 ! print*,"ZMAX ",lalim,lmin,linter,lmix,lmax,zmax,zmax0,zmix,wmax 517 534 518 519 520 535 call test_ltherm(ngrid,nlay,pplev,pplay,lalim,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lalim ') 521 536 call test_ltherm(ngrid,nlay,pplev,pplay,lmin ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lmin ') … … 547 562 endif 548 563 549 550 551 552 564 ! Choix de la fonction d'alimentation utilisee pour la fermeture. 553 565 ! Apparemment sans importance … … 584 596 !deduction des flux 585 597 !------------------------------------------------------------------------------- 586 587 598 CALL thermcell_flux2(ngrid,nlay,ptimestep,masse, & 588 599 & lalim,lmax,alim_star, & … … 620 631 call thermcell_dq(ngrid,nlay,dqimpl,ptimestep,fm0,entr0,masse, & 621 632 & po,pdoadj,zoa,lev_out) 622 623 633 !------------------------------------------------------------------ 624 634 ! Calcul de la fraction de l'ascendance 625 635 !------------------------------------------------------------------ 626 do ig=1, klon636 do ig=1,ngrid 627 637 fraca(ig,1)=0. 628 638 fraca(ig,nlay+1)=0. 629 639 enddo 630 640 do l=2,nlay 631 do ig=1, klon641 do ig=1,ngrid 632 642 if (zw2(ig,l).gt.1.e-10) then 633 643 fraca(ig,l)=fm(ig,l)/(rhobarz(ig,l)*zw2(ig,l)) … … 753 763 enddo 754 764 enddo 765 755 766 !calcul des flux: q, thetal et thetav 756 767 do l=1,nlay … … 768 779 do ig=1,ngrid 769 780 ok_lcl(ig)=.false. 770 if ( (pcon(ig) .gt. pplay(ig, klev-1)) .and. (pcon(ig) .lt. pplay(ig,1)) ) ok_lcl(ig)=.true.781 if ( (pcon(ig) .gt. pplay(ig,nlay-1)) .and. (pcon(ig) .lt. pplay(ig,1)) ) ok_lcl(ig)=.true. 771 782 enddo 772 783 … … 775 786 do ig=1,ngrid 776 787 if (ok_lcl(ig)) then 777 if ((pplay(ig,l) .ge. pcon(ig)) .and. (pplay(ig,l+1) .le. pcon(ig))) then 778 klcl(ig)=l 779 interp(ig)=(pcon(ig)-pplay(ig,klcl(ig)))/(pplay(ig,klcl(ig)+1)-pplay(ig,klcl(ig))) 788 ! if ((pplay(ig,l) .ge. pcon(ig)) .and. (pplay(ig,l+1) .le. pcon(ig))) then 789 ! L. Fita, LMD July 2014. Adding avoiding divisons by zero... 790 if ((pplay(ig,l).ge.pcon(ig)) .and. (pplay(ig,l+1).le.pcon(ig)) .and. & 791 (klcl(ig).gt.0).and.(klcl(ig)+1.le.nlay-1) .and. (klcl(ig)+1.gt.0) ) then 792 klcl(ig)=l 793 interp(ig)=(pcon(ig)-pplay(ig,klcl(ig)))/(pplay(ig,klcl(ig)+1)-pplay(ig,klcl(ig))) 780 794 endif 781 795 endif … … 794 808 do ig =1,ngrid 795 809 zmax(ig)=pphi(ig,lmax(ig))/rg 796 if (ok_lcl(ig)) then 797 rhobarz0(ig)=rhobarz(ig,klcl(ig))+(rhobarz(ig,klcl(ig)+1) & 810 ! if (ok_lcl(ig)) then 811 ! L. Fita, LMD July 2014. Adding avoiding divisons by zero... 812 if ( ok_lcl(ig) .and. (klcl(ig).gt.0) .and. (klcl(ig)+1.le.nlay-1) .and. & 813 (klcl(ig)+1.gt.0) ) then 814 rhobarz0(ig)=rhobarz(ig,klcl(ig))+(rhobarz(ig,klcl(ig)+1) & 798 815 & -rhobarz(ig,klcl(ig)))*interp(ig) 799 zlcl(ig)=(pplev(ig,1)-pcon(ig))/(rhobarz0(ig)*RG)800 zlcl(ig)=min(zlcl(ig),zmax(ig)) ! Si zlcl > zmax alors on pose zlcl = zmax816 zlcl(ig)=(pplev(ig,1)-pcon(ig))/(rhobarz0(ig)*RG) 817 zlcl(ig)=min(zlcl(ig),zmax(ig)) ! Si zlcl > zmax alors on pose zlcl = zmax 801 818 else 802 819 rhobarz0(ig)=0. … … 868 885 fraca0(ig)=0. 869 886 w0(ig)=0. 887 ! L. Fita, LMD July 2014. Adding zero value for stability issues 888 w_conv(ig) = 0. 870 889 !!jyg le 27/04/2012 871 890 !! zlcl(ig)=0. … … 933 952 ENDIF ! iflag_trig_bl 934 953 ! print *,'ENDIF iflag_trig_bl' !!jyg 935 936 954 !------------Closure------------------ 937 955 … … 991 1009 lalim_conv(:)=lalim(:) 992 1010 993 do k=1, klev1011 do k=1,nlay 994 1012 do ig=1,ngrid 995 1013 if (k<=lalim_conv(ig)) fm_tot(ig)=fm_tot(ig)+fm(ig,k) … … 999 1017 ! assez bizarre car, si on est dans la couche d'alim et que alim_star et 1000 1018 ! plus petit que 1.e-10, on prend wght_th=1. 1001 do k=1, klev1019 do k=1,nlay 1002 1020 do ig=1,ngrid 1003 1021 if (k<=lalim_conv(ig).and.alim_star(ig,k)>1.e-10) then … … 1059 1077 ratqsdiff(:,:)=0. 1060 1078 1061 do l=1, klev1079 do l=1,nlay 1062 1080 do ig=1,ngrid 1063 1081 if (l<=lalim(ig)) then … … 1069 1087 if (prt_level.ge.1) print*,'14f OK convect8' 1070 1088 1071 do l=1, klev1089 do l=1,nlay 1072 1090 do ig=1,ngrid 1073 1091 if (l<=lalim(ig)) then … … 1105 1123 !----------------------------------------------------------------------------- 1106 1124 1107 subroutine test_ltherm( klon,klev,pplev,pplay,long,seuil,ztv,po,ztva,zqla,f_star,zw2,comment)1125 subroutine test_ltherm(ngrid,nlay,pplev,pplay,long,seuil,ztv,po,ztva,zqla,f_star,zw2,comment) 1108 1126 IMPLICIT NONE 1109 1127 #include "iniprint.h" 1110 1128 1111 integer i, k, klon,klev1112 real pplev( klon,klev+1),pplay(klon,klev)1113 real ztv( klon,klev)1114 real po( klon,klev)1115 real ztva( klon,klev)1116 real zqla( klon,klev)1117 real f_star( klon,klev)1118 real zw2( klon,klev)1119 integer long( klon)1129 integer i, k, ngrid,nlay 1130 real pplev(ngrid,nlay+1),pplay(ngrid,nlay) 1131 real ztv(ngrid,nlay) 1132 real po(ngrid,nlay) 1133 real ztva(ngrid,nlay) 1134 real zqla(ngrid,nlay) 1135 real f_star(ngrid,nlay) 1136 real zw2(ngrid,nlay) 1137 integer long(ngrid) 1120 1138 real seuil 1121 1139 character*21 comment … … 1127 1145 1128 1146 ! test sur la hauteur des thermiques ... 1129 do i=1, klon1147 do i=1,ngrid 1130 1148 !IMtemp if (pplay(i,long(i)).lt.seuil*pplev(i,1)) then 1131 1149 if (prt_level.ge.10) then 1132 1150 print*,'WARNING ',comment,' au point ',i,' K= ',long(i) 1133 1151 print*,' K P(MB) THV(K) Qenv(g/kg)THVA QLA(g/kg) F* W2' 1134 do k=1, klev1152 do k=1,nlay 1135 1153 write(6,'(i3,7f10.3)') k,pplay(i,k),ztv(i,k),1000*po(i,k),ztva(i,k),1000*zqla(i,k),f_star(i,k),zw2(i,k) 1136 1154 enddo -
trunk/WRFV3/lmdz/thermcell_plume.F90
r1 r142 90 90 REAL zw2fact,expa 91 91 Zsat=.false. 92 92 93 ! Initialisation 93 94 RLvCp = RLVTT/RCPD … … 167 168 enddo 168 169 alim_star_tot(:)=1. 169 170 170 171 171 !------------------------------------------------------------------------------ … … 192 192 enddo 193 193 ! 194 195 194 !============================================================================== 196 195 !boucle de calcul de la vitesse verticale dans le thermique … … 539 538 REAL fact_gamma,fact_epsilon,fact_gamma2,fact_epsilon2 540 539 REAL c2(ngrid,klev) 540 541 541 Zsat=.false. 542 542 ! Initialisation … … 603 603 do ig=1,ngrid 604 604 if (ztv(ig,l)> ztv(ig,l+1) .and. ztv(ig,1)>=ztv(ig,l) ) then 605 alim_star(ig,l)=MAX((ztv(ig,l)-ztv(ig,l+1)),0.) & 605 ! L. Fita, LMD July 2014. Adding IF for zlev == 0. 606 IF ( zlev(ig,l+1) < 0.) THEN 607 alim_star(ig,l) = 0. 608 ELSE 609 alim_star(ig,l)=MAX((ztv(ig,l)-ztv(ig,l+1)),0.) & 606 610 & *sqrt(zlev(ig,l+1)) 607 lalim(ig)=l+1 611 lalim(ig)=l+1 612 END IF 608 613 alim_star_tot(ig)=alim_star_tot(ig)+alim_star(ig,l) 609 614 endif … … 618 623 enddo 619 624 alim_star_tot(:)=1. 620 621 622 623 625 !------------------------------------------------------------------------------ 624 626 ! Calcul dans la premiere couche … … 686 688 687 689 !------------------------------------------------ 688 !AJAM:nouveau calcul de w ²690 !AJAM:nouveau calcul de w\B2 689 691 !------------------------------------------------ 690 692 zdz=zlev(ig,l+1)-zlev(ig,l) … … 741 743 endif 742 744 enddo 743 744 745 745 746 !----------------------------------------------------------------------------
Note: See TracChangeset
for help on using the changeset viewer.