Changeset 108 for trunk/libf/dyn3dpar


Ignore:
Timestamp:
Apr 12, 2011, 11:16:02 AM (14 years ago)
Author:
slebonnois
Message:

Sebastien Lebonnois: sponge layer et dissip horizontale.

Location:
trunk/libf/dyn3dpar
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • trunk/libf/dyn3dpar/addfi_p.F

    r7 r108  
    2828c      pdufi(ip1jmp1,llm)     |
    2929c      pdvfi(ip1jm,llm)       |   respective
    30 c      pdhfi(ip1jmp1)         |      tendencies
     30c      pdhfi(ip1jmp1)         |      tendencies  (unit/s)
    3131c      pdtsfi(ip1jmp1)        |
    3232c

  • trunk/libf/dyn3dpar/calfis_p.F

    r97 r108  
    952952     
    953953! ADAPTATION GCM POUR CP(T)
    954 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    955       DO l=1,llm
    956         ztfi(1:klon,l)=ztfi(1:klon,l)+zdtfi(1:klon,l)*dtphys
    957       ENDDO
    958 !$OMP END DO
    959954      call t2tpot_p(ngridmx,llm,ztfi,zteta,zpk)
    960955

  • trunk/libf/dyn3dpar/comconst.h

    r52 r108  
    66
    77      COMMON/comconsti/im,jm,lllm,imp1,jmp1,lllmm1,lllmp1,lcl,          &
    8      &                 iflag_top_bound
     8     &                 iflag_top_bound,mode_top_bound
    99      COMMON/comconstr/dtvr,daysec,                                     &
    1010     & pi,dtphys,dtdiss,rad,r,cpp,kappa,cotot,unsim,g,omeg              &
    11      &                   ,dissip_factz,dissip_deltaz,dissip_zref        &
    12      &                   ,tau_top_bound,                                &
     11     & ,dissip_fac_mid,dissip_fac_up,dissip_deltaz,dissip_hdelta        &
     12     & ,dissip_pupstart  ,tau_top_bound,                                &
    1313     & daylen,year_day,molmass, ihf
    1414      COMMON/cpdetvenus/nu_venus,t0_venus
    … …  
    2828      REAL g ! (m/s2) gravity
    2929      REAL omeg ! (rad/s) rotation rate of the planet
    30       REAL dissip_factz,dissip_deltaz,dissip_zref
    31       INTEGER iflag_top_bound
     30      REAL dissip_fac_mid,dissip_fac_up,dissip_deltaz,dissip_hdelta
     31      REAL dissip_pupstart
     32      INTEGER iflag_top_bound,mode_top_bound
    3233      REAL tau_top_bound
    3334      REAL daylen ! length of solar day, in 'standard' day length

  • trunk/libf/dyn3dpar/conf_gcm.F

    r97 r108  
    314314! Parametres controlant la variation sur la verticale des constantes de
    315315! dissipation.
    316 ! Pour le moment actifs uniquement dans la version a 39 niveaux
    317 ! avec ok_strato=y
    318 
    319        dissip_factz=4.
    320        dissip_deltaz=10.
    321        dissip_zref=30.
    322        CALL getin('dissip_factz',dissip_factz )
     316! Actifs uniquement avec ok_strato=y
     317
     318       dissip_fac_mid=2.
     319       dissip_fac_up=10.
     320       dissip_deltaz=10.! Intervalle (km) pour le changement mid / up
     321       dissip_hdelta=5. ! scale height (km) dans la zone de la transition(m)
     322       dissip_pupstart=1.e3  ! pression (Pa) au bas la transition mid / up
     323       CALL getin('dissip_fac_mid',dissip_fac_mid )
     324       CALL getin('dissip_fac_up',dissip_fac_up )
    323325       CALL getin('dissip_deltaz',dissip_deltaz )
    324        CALL getin('dissip_zref',dissip_zref )
    325 
     326       CALL getin('dissip_hdelta',dissip_hdelta )
     327       CALL getin('dissip_pupstart',dissip_pupstart )
     328
     329! Parametres controlant la couche eponge
     330! Actifs uniquement avec ok_strato=y
     331!   mode = 0 : pas de sponge
     332!   mode = 1 : u et v -> 0
     333!   mode = 2 : u et v -> moyenne zonale
     334!   mode = 3 : u, v et h -> moyenne zonale
    326335       iflag_top_bound=1
     336       mode_top_bound=3
    327337       tau_top_bound=1.e-5
    328338       CALL getin('iflag_top_bound',iflag_top_bound)
     339       CALL getin('mode_top_bound',mode_top_bound)
    329340       CALL getin('tau_top_bound',tau_top_bound)
    330341

  • trunk/libf/dyn3dpar/leapfrog_p.F

    r101 r108  
    113113      REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: dqfi
    114114
    115        !! Aymeric -- cp(T) comme dans leapfrog.F, SAVE OK ???
    116       REAL,SAVE :: duspg(ip1jmp1,llm) ! for bilan_dyn
    117 
     115c   tendances top_bound (sponge layer)
     116      REAL,SAVE :: dvtop(ip1jm,llm),dutop(ip1jmp1,llm)
     117      REAL,SAVE :: dtetatop(ip1jmp1,llm)
     118      REAL,SAVE :: dptop(ip1jmp1)
     119      REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: dqtop
    118120
    119121c   variables pour le fichier histoire
    … …  
    242244            ALLOCATE(dqfi(ip1jmp1,llm,nqtot))
    243245            ALLOCATE(dqfi_tmp(iip1,llm,nqtot))
     246            ALLOCATE(dqtop(ip1jmp1,llm,nqtot))
    244247         END IF
    245248c$OMP END MASTER     
    … …  
    251254
    252255c$OMP MASTER
     256c INITIALISATIONS
     257        dudis(:,:)   =0.
     258        dvdis(:,:)   =0.
     259        dtetadis(:,:)=0.
     260        dutop(:,:)   =0.
     261        dvtop(:,:)   =0.
     262        dtetatop(:,:)=0.
     263        dqtop(:,:,:) =0.
     264        dptop(:)     =0.
     265        dufi(:,:)   =0.
     266        dvfi(:,:)   =0.
     267        dtetafi(:,:)=0.
     268        dqfi(:,:,:) =0.
     269        dpfi(:)     =0.
    253270      dq(:,:,:)=0.
     271
    254272      CALL pression ( ip1jmp1, ap, bp, ps, p       )
    255273      CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
    … …  
    914932c      ajout des tendances physiques:
    915933c      ------------------------------
    916          IF (ok_strato) THEN
    917            CALL top_bound_p( vcov,ucov,teta,masse,dufi,dvfi,dtetafi)
    918          ENDIF
    919        
    920934          CALL addfi_p( dtphys, leapf, forward   ,
    921935     $                  ucov, vcov, teta , q   ,ps ,
    922936     $                 dufi, dvfi, dtetafi , dqfi ,dpfi  )
    923937
     938c      Couche superieure :
     939c      -------------------
     940         IF (ok_strato) THEN
     941           CALL top_bound_p( vcov,ucov,teta,masse,dutop,dvtop,dtetatop)
     942           CALL addfi_p( dtphys, leapf, forward   ,
     943     $                  ucov, vcov, teta , q   ,ps ,
     944     $                 dutop, dvtop, dtetatop , dqtop ,dptop  )
     945
     946         ENDIF
     947       
    924948c$OMP BARRIER
    925949c$OMP MASTER
    … …  
    14441468                 CALL bilan_dyn_p(dtvr*iperiod,dtvr*day_step*periodav,
    14451469     &                 ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,
    1446      &                 du,dudis,duspg,dufi)
     1470     &                 du,dudis,dutop,dufi)
    14471471c$OMP END MASTER
    14481472              ENDIF !ok_dynzon
    … …  
    16701694                 CALL bilan_dyn_p(dtvr*iperiod,dtvr*day_step*periodav,
    16711695     &                 ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,
    1672      &                 du,dudis,duspg,dufi)
     1696     &                 du,dudis,dutop,dufi)
    16731697
    16741698c$OMP END MASTER

  • trunk/libf/dyn3dpar/top_bound_p.F

    r1 r108  
    6464             rdamp(llm-3)=tau_top_bound/8.
    6565         else if (iflag_top_bound == 2) then
    66 ! couce eponge dans toutes les couches de pression plus faible que
     66! couche eponge dans toutes les couches de pression plus faible que
    6767! 100 fois la pression de la derniere couche
    6868             rdamp(:)=tau_top_bound
    … …  
    7070         endif
    7171         first=.false.
    72          print*,'TOP_BOUND rdamp=',rdamp
     72         print*,'TOP_BOUND mode',mode_top_bound
     73         print*,'Coeffs pour la couche eponge a l equateur'
     74         print*,'p (Pa)  z(km)  tau (s)   dt*rdamp'
     75         do l=1,llm
     76           if (rdamp(l).ne.0.) then
     77            zkm        = phi(iip1/2,jjp1/2,l)/(1000*g)
     78          print*,presnivs(l),zkm,
     79     .          1./rdamp(l),
     80     .          dt*rdamp(l)
     81           endif
     82         enddo
    7383c$OMP END MASTER
    7484c$OMP BARRIER
    … …  
    7787
    7888      CALL massbar_p(masse,massebx,masseby)
     89
     90c   mode = 0 : pas de sponge
     91c   mode = 1 : u et v -> 0
     92c   mode = 2 : u et v -> moyenne zonale
     93c   mode = 3 : u, v et h -> moyenne zonale
     94
     95C POUR V
     96
    7997C  CALCUL DES CHAMPS EN MOYENNE ZONALE:
    8098
    … …  
    84102
    85103c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
    86       do l=1,llm
     104      if (mode_top_bound.ge.2) then
     105       do l=1,llm
    87106        do j=jjb,jje
    88107          zm=0.
    … …  
    97116          vzon(j,l)=vzon(j,l)/zm
    98117        enddo
    99       enddo
     118       enddo
     119      else
     120       do l=1,llm
     121        do j=jjb,jje
     122          vzon(j,l)=0.
     123        enddo
     124       enddo
     125      endif
    100126c$OMP END DO NOWAIT   
    101127
    102 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
    103       do l=1,llm
     128C   AMORTISSEMENTS LINEAIRES:
     129
     130c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
     131      if (mode_top_bound.ge.1) then
     132       do l=1,llm
    104133        do j=jjb,jje
    105134          do i=1,iip1
    106             dv(i,j,l)=dv(i,j,l)-rdamp(l)*(vcov(i,j,l)-vzon(j,l))
    107           enddo
    108         enddo
    109       enddo
    110 c$OMP END DO NOWAIT
     135            dv(i,j,l)= -rdamp(l)*(vcov(i,j,l)-vzon(j,l))
     136          enddo
     137        enddo
     138       enddo
     139      endif
     140c$OMP END DO NOWAIT
     141
     142C POUR U ET H
     143
     144C  CALCUL DES CHAMPS EN MOYENNE ZONALE:
    111145
    112146      jjb=jj_begin
    … …  
    116150
    117151c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
    118       do l=1,llm
     152      if (mode_top_bound.ge.2) then
     153       do l=1,llm
    119154        do j=jjb,jje
    120155          uzon(j,l)=0.
    … …  
    126161          uzon(j,l)=uzon(j,l)/zm
    127162        enddo
    128       enddo
     163       enddo
     164      else
     165       do l=1,llm
     166        do j=jjb,jje
     167          uzon(j,l)=0.
     168        enddo
     169       enddo
     170      endif
    129171c$OMP END DO NOWAIT
    130172
    131173c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)   
    132       do l=1,llm
     174      if (mode_top_bound.ge.3) then
     175       do l=1,llm
    133176        do j=jjb,jje
    134177          zm=0.
    … …  
    140183          tzon(j,l)=tzon(j,l)/zm
    141184        enddo
    142       enddo
     185       enddo
     186      endif
    143187c$OMP END DO NOWAIT
    144188
    … …  
    146190
    147191c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
    148       do l=1,llm
     192      if (mode_top_bound.ge.1) then
     193       do l=1,llm
    149194        do j=jjb,jje
    150195          do i=1,iip1
    151             du(i,j,l)=du(i,j,l)
    152      s               -rdamp(l)*(ucov(i,j,l)-cu(i,j)*uzon(j,l))
    153             dh(i,j,l)=dh(i,j,l)-rdamp(l)*(teta(i,j,l)-tzon(j,l))
    154           enddo
    155        enddo
    156       enddo
     196            du(i,j,l)= -rdamp(l)*(ucov(i,j,l)-cu(i,j)*uzon(j,l))
     197          enddo
     198        enddo
     199       enddo
     200      endif
     201c$OMP END DO NOWAIT
     202     
     203c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
     204      if (mode_top_bound.ge.3) then
     205       do l=1,llm
     206        do j=jjb,jje
     207          do i=1,iip1
     208            dh(i,j,l)= -rdamp(l)*(teta(i,j,l)-tzon(j,l))
     209          enddo
     210        enddo
     211       enddo
     212      endif
    157213c$OMP END DO NOWAIT
    158214     
Note: See TracChangeset for help on using the changeset viewer.