Changeset 125


Ignore:
Timestamp:
02/01/13 01:52:22 (11 years ago)
Author:
dubos
Message:

Ready for energy-conserving variant

File:
1 edited

Legend:

Unmodified
Added
Removed

  • codes/icosagcm/trunk/src/caldyn_gcm.f90

    r123 r125  
    22  USE icosa 
    33 
    4   TYPE(t_field),POINTER :: f_out_u(:) 
    5   REAL(rstd),POINTER :: out_u(:,:) 
     4  TYPE(t_field),POINTER :: f_out_u(:), f_p(:), f_rhodz(:), f_qu(:) 
     5  REAL(rstd),POINTER :: out_u(:,:), p(:,:), rhodz(:,:), qu(:,:) 
    66 
    77  TYPE(t_field),POINTER :: f_buf_i(:), f_buf_ulon(:), f_buf_ulat(:), f_buf_u3d(:) 
    … …  
    1010  PUBLIC init_caldyn, caldyn, write_output_fields 
    1111 
    12   INTEGER :: caldyn_hydrostat 
     12  INTEGER :: caldyn_hydrostat, caldyn_conserv 
    1313 
    1414CONTAINS 
    … …  
    2020    CHARACTER(len=255) :: def 
    2121   
     22    def='direct' 
     23    CALL getin('caldyn_conserv',def) 
     24    SELECT CASE(TRIM(def)) 
     25    CASE('energy') 
     26       caldyn_conserv=1 
     27    CASE('direct') 
     28       caldyn_conserv=2 
     29    CASE DEFAULT 
     30       PRINT*,'Bad selector for variable caldyn_conserv : <', TRIM(def),'> options are <energy>, <enstrophy>' 
     31       STOP 
     32    END SELECT 
     33 
    2234    def='direct' 
    2335    CALL getin('caldyn_exner',def) 
    … …  
    5365 
    5466  CALL allocate_field(f_out_u,field_u,type_real,llm)  
     67  CALL allocate_field(f_p,field_t,type_real,llm+1) 
     68  CALL allocate_field(f_rhodz,field_t,type_real,llm)  
     69  CALL allocate_field(f_qu,field_u,type_real,llm)  
    5570   
    5671  CALL allocate_field(f_buf_i,field_t,type_real,llm) 
    … …  
    118133  TYPE(t_field),POINTER :: f_du(:) 
    119134 
    120   REAL(rstd),POINTER :: phis(:) 
    121   REAL(rstd),POINTER :: ps(:) 
    122   REAL(rstd),POINTER :: theta_rhodz(:,:) 
    123   REAL(rstd),POINTER :: u(:,:) 
    124   REAL(rstd),POINTER :: dps(:) 
    125   REAL(rstd),POINTER :: dtheta_rhodz(:,:) 
    126   REAL(rstd),POINTER :: du(:,:) 
     135  REAL(rstd),POINTER :: phis(:), ps(:), dps(:) 
     136  REAL(rstd),POINTER :: theta_rhodz(:,:), dtheta_rhodz(:,:) 
     137  REAL(rstd),POINTER :: u(:,:), du(:,:) 
     138  REAL(rstd),POINTER :: p(:,:), rhodz(:,:), qu(:,:) 
    127139  INTEGER :: ind,ij 
    128140 
    … …  
    137149      CALL swap_geometry(ind) 
    138150       
    139       out_u=f_out_u(ind)  
    140151      phis=f_phis(ind) 
    141152      ps=f_ps(ind) 
     153      dps=f_dps(ind) 
    142154      theta_rhodz=f_theta_rhodz(ind) 
     155      dtheta_rhodz=f_dtheta_rhodz(ind) 
     156      rhodz=f_rhodz(ind) 
     157      p=f_p(ind) 
     158      qu=f_qu(ind) 
    143159      u=f_u(ind) 
    144       dps=f_dps(ind) 
    145       dtheta_rhodz=f_dtheta_rhodz(ind) 
    146160      du=f_du(ind) 
     161      out_u=f_out_u(ind)  
    147162!$OMP PARALLEL DEFAULT(SHARED) 
    148       CALL compute_caldyn(phis, ps, theta_rhodz, u, dps, dtheta_rhodz, du) 
     163      CALL compute_caldyn(ps, u, p,rhodz,qu, phis, theta_rhodz, dps, dtheta_rhodz, du) 
    149164!$OMP END PARALLEL 
    150165    ENDDO 
    … …  
    161176   
    162177   
    163   SUBROUTINE compute_caldyn(phis, ps, theta_rhodz, u, dps, dtheta_rhodz, du) 
     178  SUBROUTINE compute_caldyn(ps, u, p,rhodz,qu, phis, theta_rhodz, dps, dtheta_rhodz, du) 
    164179  USE icosa 
    165180  USE disvert_mod 
    166181  USE exner_mod 
    167182  IMPLICIT NONE 
    168     REAL(rstd),INTENT(IN) :: phis(iim*jjm) 
    169     REAL(rstd),INTENT(IN) :: u(iim*3*jjm,llm) 
    170     REAL(rstd),INTENT(IN) :: theta_rhodz(iim*jjm,llm) 
    171     REAL(rstd),INTENT(IN) :: ps(iim*jjm) 
     183    REAL(rstd),INTENT(IN)  :: phis(iim*jjm) 
     184    REAL(rstd),INTENT(IN)  :: u(iim*3*jjm,llm) 
     185    REAL(rstd),INTENT(IN)  :: theta_rhodz(iim*jjm,llm) 
     186    REAL(rstd),INTENT(IN)  :: ps(iim*jjm) 
     187 
     188    REAL(rstd),INTENT(OUT) :: p(iim*jjm,llm+1) 
     189    REAL(rstd),INTENT(OUT) :: rhodz(iim*jjm,llm) 
     190    REAL(rstd),INTENT(OUT) :: qu(iim*3*jjm,llm) 
     191 
    172192    REAL(rstd),INTENT(OUT) :: du(iim*3*jjm,llm) 
    173     REAL(rstd),INTENT(OUT):: dtheta_rhodz(iim*jjm,llm) 
    174     REAL(rstd),INTENT(OUT):: dps(iim*jjm) 
     193    REAL(rstd),INTENT(OUT) :: dtheta_rhodz(iim*jjm,llm) 
     194    REAL(rstd),INTENT(OUT) :: dps(iim*jjm) 
    175195     
    176196    INTEGER :: i,j,ij,l 
    … …  
    180200 
    181201    REAL(rstd),ALLOCATABLE,SAVE :: theta(:,:)  ! potential temperature 
    182     REAL(rstd),ALLOCATABLE,SAVE :: p(:,:)  ! pression 
    183202    REAL(rstd),ALLOCATABLE,SAVE :: pk(:,:)   ! Exner function 
    184203    REAL(rstd),ALLOCATABLE,SAVE :: pks(:) 
    185  ! Intermediate variable to compute exner function 
     204! Intermediate variable to compute exner function 
    186205    REAL(rstd),ALLOCATABLE,SAVE :: alpha(:,:) 
    187206    REAL(rstd),ALLOCATABLE,SAVE :: beta(:,:) 
    188  !    
     207    
    189208    REAL(rstd),ALLOCATABLE,SAVE :: phi(:,:)    ! geopotential 
    190     REAL(rstd),ALLOCATABLE,SAVE :: rhodz(:,:)   ! mass density    
    191209    REAL(rstd),ALLOCATABLE,SAVE :: Fe(:,:)   ! mass flux    
    192210    REAL(rstd),ALLOCATABLE,SAVE :: Ftheta(:,:) ! theta flux    
    … …  
    202220!$OMP MASTER   
    203221!    IF (first) THEN 
    204       ALLOCATE(theta(iim*jjm,llm))  ! potential temperature 
    205       ALLOCATE(p(iim*jjm,llm+1))  ! pression 
    206       ALLOCATE(pk(iim*jjm,llm))   ! Exner function 
     222      ALLOCATE(theta(iim*jjm,llm))    ! potential temperature 
     223      ALLOCATE(pk(iim*jjm,llm))       ! Exner function 
    207224      ALLOCATE(pks(iim*jjm)) 
    208225      ALLOCATE(alpha(iim*jjm,llm)) 
    209226      ALLOCATE(beta(iim*jjm,llm)) 
    210       ALLOCATE(phi(iim*jjm,llm))    ! geopotential 
    211       ALLOCATE(rhodz(iim*jjm,llm))   ! mass density    
    212       ALLOCATE(Fe(3*iim*jjm,llm))   ! mass flux    
     227      ALLOCATE(phi(iim*jjm,llm))      ! geopotential 
     228      ALLOCATE(Fe(3*iim*jjm,llm))     ! mass flux    
    213229      ALLOCATE(Ftheta(3*iim*jjm,llm)) ! theta flux    
    214230      ALLOCATE(convm(iim*jjm,llm))    ! mass flux convergence 
    215       ALLOCATE(w(iim*jjm,llm))       ! vertical velocity       
    216       ALLOCATE(qv(2*iim*jjm,llm))       ! potential velocity   
    217       ALLOCATE(berni(iim*jjm,llm))   ! bernouilli term 
     231      ALLOCATE(w(iim*jjm,llm))        ! vertical velocity       
     232      ALLOCATE(qv(2*iim*jjm,llm))     ! potential velocity   
     233      ALLOCATE(berni(iim*jjm,llm))    ! bernouilli term 
    218234!      first=.FALSE. 
    219235!    ENDIF 
    … …  
    292308        ENDDO 
    293309      ENDDO 
    294     ENDDO 
    295  
     310 
     311    DO j=jj_begin,jj_end 
     312      DO i=ii_begin,ii_end 
     313        ij=(j-1)*iim+i 
     314        qu(ij+u_right,l) = 0.5*(qv(ij+z_rdown,l)+qv(ij+z_rup,l))   
     315        qu(ij+u_lup,l) = 0.5*(qv(ij+z_up,l)+qv(ij+z_lup,l))   
     316        qu(ij+u_ldown,l) = 0.5*(qv(ij+z_ldown,l)+qv(ij+z_down,l))   
     317     END DO 
     318  END DO 
     319   
     320ENDDO 
     321 
     322!----------- needed next : rhodz, p --------------! 
    296323 
    297324  DO l = 1, llm 
    … …  
    378405        ij=(j-1)*iim+i 
    379406 
    380        du(ij+u_right,l) = 0.5*(qv(ij+z_rdown,l)+qv(ij+z_rup,l))/de(ij+u_right) *      &   
     407       du(ij+u_right,l) = qu(ij+u_right,l)/de(ij+u_right) *      &   
    381408                      ( wee(ij+u_right,1,1)*Fe(ij+u_rup,l)+                           & 
    382409                        wee(ij+u_right,2,1)*Fe(ij+u_lup,l)+                           & 
    … …  
    391418 
    392419 
    393        du(ij+u_lup,l) = 0.5*(qv(ij+z_up,l)+qv(ij+z_lup,l))/de(ij+u_lup) *         &   
     420       du(ij+u_lup,l) = qu(ij+u_lup,l)/de(ij+u_lup) *         &   
    394421                      ( wee(ij+u_lup,1,1)*Fe(ij+u_left,l)+                        & 
    395422                        wee(ij+u_lup,2,1)*Fe(ij+u_ldown,l)+                       & 
    … …  
    404431 
    405432 
    406        du(ij+u_ldown,l) = 0.5*(qv(ij+z_ldown,l)+qv(ij+z_down,l))/de(ij+u_ldown) *   &   
     433       du(ij+u_ldown,l) = qu(ij+u_ldown,l)/de(ij+u_ldown) *   &   
    407434                      ( wee(ij+u_ldown,1,1)*Fe(ij+u_rdown,l)+                      & 
    408435                        wee(ij+u_ldown,2,1)*Fe(ij+u_right,l)+                      & 
    … …  
    536563!!$OMP MASTER   
    537564    DEALLOCATE(theta)  ! potential temperature 
    538     DEALLOCATE(p)  ! pression 
    539565    DEALLOCATE(pk)   ! Exner function 
    540566    DEALLOCATE(pks) 
    … …  
    542568    DEALLOCATE(beta) 
    543569    DEALLOCATE(phi)    ! geopotential 
    544 !    DEALLOCATE(mass)   ! mass     
    545     DEALLOCATE(rhodz)   ! mass density    
    546570    DEALLOCATE(Fe)   ! mass flux    
    547571    DEALLOCATE(Ftheta) ! theta flux    
Note: See TracChangeset for help on using the changeset viewer.