Changeset 1795 for LMDZ5/branches/testing/libf/dyn3d
- Timestamp:
- Jul 18, 2013, 10:20:28 AM (11 years ago)
- Location:
- LMDZ5/branches/testing
- Files:
-
- 3 deleted
- 10 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ5/branches/testing
- Property svn:mergeinfo changed
/LMDZ5/trunk merged: 1747-1749,1751,1753-1767,1769,1771-1772,1774-1776,1778-1794
- Property svn:mergeinfo changed
-
LMDZ5/branches/testing/libf/dyn3d/ce0l.F90
r1665 r1795 20 20 USE comgeomphy 21 21 USE infotrac 22 USE indice_sol_mod 22 23 23 24 #ifdef CPP_IOIPSL … … 36 37 #include "dimensions.h" 37 38 #include "paramet.h" 38 #include "indicesol.h"39 !#include "indicesol.h" 39 40 #include "iniprint.h" 40 41 #include "temps.h" -
LMDZ5/branches/testing/libf/dyn3d/comconst.h
r1707 r1795 6 6 7 7 COMMON/comconsti/im,jm,lllm,imp1,jmp1,lllmm1,lllmp1,lcl, & 8 & iflag_top_bound 8 & iflag_top_bound,mode_top_bound 9 9 COMMON/comconstr/dtvr,daysec, & 10 10 & pi,dtphys,dtdiss,rad,r,cpp,kappa,cotot,unsim,g,omeg & … … 30 30 REAL omeg ! (rad/s) rotation rate of the planet 31 31 REAL dissip_factz,dissip_deltaz,dissip_zref 32 INTEGER iflag_top_bound 33 REAL tau_top_bound 32 ! top_bound sponge: 33 INTEGER iflag_top_bound ! sponge type 34 INTEGER mode_top_bound ! sponge mode 35 REAL tau_top_bound ! inverse of sponge characteristic time scale (Hz) 34 36 REAL daylen ! length of solar day, in 'standard' day length 35 37 REAL year_day ! Number of standard days in a year -
LMDZ5/branches/testing/libf/dyn3d/comvert.h
r1669 r1795 23 23 real bps ! hybrid sigma contribution at mid-layers 24 24 real scaleheight ! atmospheric (reference) scale height (km) 25 real pseudoalt ! for planets 25 real pseudoalt ! pseudo-altitude of model levels (km), based on presnivs(), 26 ! preff and scaleheight 26 27 27 28 integer disvert_type ! type of vertical discretization: -
LMDZ5/branches/testing/libf/dyn3d/conf_gcm.F
r1707 r1795 307 307 CALL getin('dissip_zref',dissip_zref ) 308 308 309 ! top_bound sponge: only active if ok_strato=.true. and iflag_top_bound!=0 310 ! iflag_top_bound=0 for no sponge 311 ! iflag_top_bound=1 for sponge over 4 topmost layers 312 ! iflag_top_bound=2 for sponge from top to ~1% of top layer pressure 309 313 iflag_top_bound=1 314 CALL getin('iflag_top_bound',iflag_top_bound) 315 316 ! mode_top_bound : fields towards which sponge relaxation will be done: 317 ! mode_top_bound=0: no relaxation 318 ! mode_top_bound=1: u and v relax towards 0 319 ! mode_top_bound=2: u and v relax towards their zonal mean 320 ! mode_top_bound=3: u,v and pot. temp. relax towards their zonal mean 321 mode_top_bound=3 322 CALL getin('mode_top_bound',mode_top_bound) 323 324 ! top_bound sponge : inverse of charactericstic relaxation time scale for sponge 310 325 tau_top_bound=1.e-5 311 CALL getin('iflag_top_bound',iflag_top_bound)312 326 CALL getin('tau_top_bound',tau_top_bound) 313 327 -
LMDZ5/branches/testing/libf/dyn3d/gcm.F
r1707 r1795 16 16 USE infotrac 17 17 USE control_mod 18 19 #ifdef INCA 20 ! Only INCA needs these informations (from the Earth's physics) 21 USE indice_sol_mod 22 #endif 18 23 19 24 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! … … 77 82 #ifdef INCA 78 83 ! Only INCA needs these informations (from the Earth's physics) 79 #include "indicesol.h"84 !#include "indicesol.h" 80 85 #endif 81 86 INTEGER longcles -
LMDZ5/branches/testing/libf/dyn3d/inigeom.F
r1403 r1795 426 426 radclatm = 0.5* rad * coslatm 427 427 c 428 ai14 = un4rad2 * coslatp * yprp 429 ai23 = un4rad2 * coslatm * yprm 428 430 DO 32 i = 1,iim 429 431 xprp = xprimp025( i ) 430 432 xprm = xprimm025( i ) 431 433 432 ai14 = un4rad2 * coslatp * yprp433 ai23 = un4rad2 * coslatm * yprm434 434 aireij1 ( i,j ) = ai14 * xprp 435 435 aireij2 ( i,j ) = ai23 * xprp -
LMDZ5/branches/testing/libf/dyn3d/leapfrog.F
r1707 r1795 436 436 $ clesphy0, dufi,dvfi,dtetafi,dqfi,dpfi ) 437 437 438 IF (ok_strato) THEN439 CALL top_bound( vcov,ucov,teta,masse,dufi,dvfi,dtetafi)440 ENDIF441 442 438 c ajout des tendances physiques: 443 439 c ------------------------------ … … 445 441 $ ucov, vcov, teta , q ,ps , 446 442 $ dufi, dvfi, dtetafi , dqfi ,dpfi ) 443 444 IF (ok_strato) THEN 445 CALL top_bound( vcov,ucov,teta,masse,dtphys) 446 ENDIF 447 447 448 c 448 449 c Diagnostique de conservation de l'énergie : difference … … 476 477 ! Sponge layer (if any) 477 478 IF (ok_strato) THEN 478 dufi(:,:)=0. 479 dvfi(:,:)=0. 480 dtetafi(:,:)=0. 481 dqfi(:,:,:)=0. 482 dpfi(:)=0. 483 CALL top_bound(vcov,ucov,teta,masse,dufi,dvfi,dtetafi) 484 CALL addfi( dtvr, leapf, forward , 485 $ ucov, vcov, teta , q ,ps , 486 $ dufi, dvfi, dtetafi , dqfi ,dpfi ) 479 ! dufi(:,:)=0. 480 ! dvfi(:,:)=0. 481 ! dtetafi(:,:)=0. 482 ! dqfi(:,:,:)=0. 483 ! dpfi(:)=0. 484 ! CALL top_bound(vcov,ucov,teta,masse,dufi,dvfi,dtetafi) 485 CALL top_bound( vcov,ucov,teta,masse,dtvr) 486 ! CALL addfi( dtvr, leapf, forward , 487 ! $ ucov, vcov, teta , q ,ps , 488 ! $ dufi, dvfi, dtetafi , dqfi ,dpfi ) 487 489 ENDIF ! of IF (ok_strato) 488 490 ENDIF ! of IF (iflag_phys.EQ.2) -
LMDZ5/branches/testing/libf/dyn3d/temps.h
r1665 r1795 13 13 ! INCLUDE 'temps.h' 14 14 15 COMMON/temps/ itaufin, dt, day_ini, day_end, annee_ref, day_ref,&16 & itau_dyn, itau_phy, jD_ref, jH_ref, calend,&17 & start_time15 COMMON/temps/ dt, jD_ref, jH_ref, start_time, & 16 & day_ini, day_end, annee_ref, day_ref, & 17 & itau_dyn, itau_phy, itaufin, calend 18 18 19 19 20 20 INTEGER itaufin 21 INTEGER itau_dyn, itau_phy 21 INTEGER itau_dyn, itau_phy 22 22 INTEGER day_ini, day_end, annee_ref, day_ref 23 23 REAL dt, jD_ref, jH_ref, start_time -
LMDZ5/branches/testing/libf/dyn3d/top_bound.F
r1279 r1795 1 SUBROUTINE top_bound( vcov,ucov,teta,masse, du,dv,dh ) 1 ! 2 ! $Id$ 3 ! 4 SUBROUTINE top_bound(vcov,ucov,teta,masse,dt) 2 5 IMPLICIT NONE 3 6 c … … 24 27 c 25 28 c======================================================================= 26 c-----------------------------------------------------------------------27 c Declarations:28 c -------------29 29 30 ! #include "comgeom.h" 30 ! top_bound sponge layer model: 31 ! Quenching is modeled as: A(t)=Am+A0*exp(-lambda*t) 32 ! where Am is the zonal average of the field (or zero), and lambda the inverse 33 ! of the characteristic quenching/relaxation time scale 34 ! Thus, assuming Am to be time-independent, field at time t+dt is given by: 35 ! A(t+dt)=A(t)-(A(t)-Am)*(1-exp(-lambda*t)) 36 ! Moreover lambda can be a function of model level (see below), and relaxation 37 ! can be toward the average zonal field or just zero (see below). 38 39 ! NB: top_bound sponge is only called from leapfrog if ok_strato=.true. 40 41 ! sponge parameters: (loaded/set in conf_gcm.F ; stored in comconst.h) 42 ! iflag_top_bound=0 for no sponge 43 ! iflag_top_bound=1 for sponge over 4 topmost layers 44 ! iflag_top_bound=2 for sponge from top to ~1% of top layer pressure 45 ! mode_top_bound=0: no relaxation 46 ! mode_top_bound=1: u and v relax towards 0 47 ! mode_top_bound=2: u and v relax towards their zonal mean 48 ! mode_top_bound=3: u,v and pot. temp. relax towards their zonal mean 49 ! tau_top_bound : inverse of charactericstic relaxation time scale at 50 ! the topmost layer (Hz) 51 52 31 53 #include "comdissipn.h" 54 #include "iniprint.h" 32 55 33 56 c Arguments: 34 57 c ---------- 35 58 36 REAL ucov(iip1,jjp1,llm),vcov(iip1,jjm,llm),teta(iip1,jjp1,llm) 37 REAL masse(iip1,jjp1,llm) 38 REAL dv(iip1,jjm,llm),du(iip1,jjp1,llm),dh(iip1,jjp1,llm) 59 real,intent(inout) :: ucov(iip1,jjp1,llm) ! covariant zonal wind 60 real,intent(inout) :: vcov(iip1,jjm,llm) ! covariant meridional wind 61 real,intent(inout) :: teta(iip1,jjp1,llm) ! potential temperature 62 real,intent(in) :: masse(iip1,jjp1,llm) ! mass of atmosphere 63 real,intent(in) :: dt ! time step (s) of sponge model 39 64 40 65 c Local: … … 44 69 REAL uzon(jjp1,llm),vzon(jjm,llm),tzon(jjp1,llm) 45 70 46 INTEGER NDAMP47 PARAMETER (NDAMP=4)48 71 integer i 49 REAL,SAVE :: rdamp(llm) 50 ! & (/(0., i =1,llm-NDAMP),0.125E-5,.25E-5,.5E-5,1.E-5/) 72 REAL,SAVE :: rdamp(llm) ! quenching coefficient 73 real,save :: lambda(llm) ! inverse or quenching time scale (Hz) 51 74 52 75 LOGICAL,SAVE :: first=.true. 53 76 54 77 INTEGER j,l 55 56 57 C CALCUL DES CHAMPS EN MOYENNE ZONALE:58 78 59 79 if (iflag_top_bound.eq.0) return … … 61 81 if (first) then 62 82 if (iflag_top_bound.eq.1) then 63 ! couche eponge dans les 4 dernieres couches du modele64 rdamp(:)=0.65 rdamp(llm)=tau_top_bound66 rdamp(llm-1)=tau_top_bound/2.67 rdamp(llm-2)=tau_top_bound/4.68 rdamp(llm-3)=tau_top_bound/8.83 ! sponge quenching over the topmost 4 atmospheric layers 84 lambda(:)=0. 85 lambda(llm)=tau_top_bound 86 lambda(llm-1)=tau_top_bound/2. 87 lambda(llm-2)=tau_top_bound/4. 88 lambda(llm-3)=tau_top_bound/8. 69 89 else if (iflag_top_bound.eq.2) then 70 ! couce eponge dans toutes les couches de pression plus faible que71 ! 100 fois la pression de la derniere couche72 rdamp(:)=tau_top_bound90 ! sponge quenching over topmost layers down to pressures which are 91 ! higher than 100 times the topmost layer pressure 92 lambda(:)=tau_top_bound 73 93 s *max(presnivs(llm)/presnivs(:)-0.01,0.) 74 94 endif 95 96 ! quenching coefficient rdamp(:) 97 ! rdamp(:)=dt*lambda(:) ! Explicit Euler approx. 98 rdamp(:)=1.-exp(-lambda(:)*dt) 99 100 write(lunout,*)'TOP_BOUND mode',mode_top_bound 101 write(lunout,*)'Sponge layer coefficients' 102 write(lunout,*)'p (Pa) z(km) tau(s) 1./tau (Hz)' 103 do l=1,llm 104 if (rdamp(l).ne.0.) then 105 write(lunout,'(6(1pe12.4,1x))') 106 & presnivs(l),log(preff/presnivs(l))*scaleheight, 107 & 1./lambda(l),lambda(l) 108 endif 109 enddo 75 110 first=.false. 76 print*,'TOP_BOUND rdamp=',rdamp 77 endif 111 endif ! of if (first) 78 112 79 113 CALL massbar(masse,massebx,masseby) 80 114 81 do l=1,llm 115 ! compute zonal average of vcov and u 116 if (mode_top_bound.ge.2) then 117 do l=1,llm 82 118 do j=1,jjm 83 119 vzon(j,l)=0. 84 120 zm=0. 85 121 do i=1,iim 86 ! Rm: on peut travailler directement avec la moyenne zonale de vcov 87 ! plutot qu'avec celle de v car le coefficient cv qui relie les deux 88 ! ne varie qu'en latitude 122 ! NB: we can work using vcov zonal mean rather than v since the 123 ! cv coefficient (which relates the two) only varies with latitudes 89 124 vzon(j,l)=vzon(j,l)+vcov(i,j,l)*masseby(i,j,l) 90 125 zm=zm+masseby(i,j,l) … … 92 127 vzon(j,l)=vzon(j,l)/zm 93 128 enddo 94 enddo129 enddo 95 130 96 do l=1,llm 97 do i=1,iip1 98 do j=1,jjm 99 dv(i,j,l)=dv(i,j,l)-rdamp(l)*(vcov(i,j,l)-vzon(j,l)) 100 enddo 101 enddo 102 enddo 103 104 do l=1,llm 105 do j=2,jjm 131 do l=1,llm 132 do j=2,jjm ! excluding poles 106 133 uzon(j,l)=0. 107 134 zm=0. … … 112 139 uzon(j,l)=uzon(j,l)/zm 113 140 enddo 114 enddo 141 enddo 142 else ! ucov and vcov will relax towards 0 143 vzon(:,:)=0. 144 uzon(:,:)=0. 145 endif ! of if (mode_top_bound.ge.2) 115 146 116 do l=1,llm 117 do j=2,jjm 147 ! compute zonal average of potential temperature, if necessary 148 if (mode_top_bound.ge.3) then 149 do l=1,llm 150 do j=2,jjm ! excluding poles 118 151 zm=0. 119 152 tzon(j,l)=0. … … 124 157 tzon(j,l)=tzon(j,l)/zm 125 158 enddo 126 enddo 159 enddo 160 endif ! of if (mode_top_bound.ge.3) 127 161 128 C AMORTISSEMENTS LINEAIRES: 129 130 do l=1,llm162 if (mode_top_bound.ge.1) then 163 ! Apply sponge quenching on vcov: 164 do l=1,llm 131 165 do i=1,iip1 132 do j=2,jjm 133 du(i,j,l)=du(i,j,l) 134 s -rdamp(l)*(ucov(i,j,l)-cu(i,j)*uzon(j,l)) 135 dh(i,j,l)=dh(i,j,l)-rdamp(l)*(teta(i,j,l)-tzon(j,l)) 166 do j=1,jjm 167 vcov(i,j,l)=vcov(i,j,l) 168 & -rdamp(l)*(vcov(i,j,l)-vzon(j,l)) 136 169 enddo 137 170 enddo 138 enddo 139 171 enddo 140 172 141 RETURN 173 ! Apply sponge quenching on ucov: 174 do l=1,llm 175 do i=1,iip1 176 do j=2,jjm ! excluding poles 177 ucov(i,j,l)=ucov(i,j,l) 178 & -rdamp(l)*(ucov(i,j,l)-cu(i,j)*uzon(j,l)) 179 enddo 180 enddo 181 enddo 182 endif ! of if (mode_top_bound.ge.1) 183 184 if (mode_top_bound.ge.3) then 185 ! Apply sponge quenching on teta: 186 do l=1,llm 187 do i=1,iip1 188 do j=2,jjm ! excluding poles 189 teta(i,j,l)=teta(i,j,l) 190 & -rdamp(l)*(teta(i,j,l)-tzon(j,l)) 191 enddo 192 enddo 193 enddo 194 endif ! of if (mode_top_bound.ge.3) 195 142 196 END
Note: See TracChangeset
for help on using the changeset viewer.