Back to home page

MITgcm
 
 

    

File indexing completed on 2018-03-02 18:45:26 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
4868e56500 Andr*0001 #include "FIZHI_OPTIONS.h"
eec2a90a83 Jean*0002        SUBROUTINE FIZHI_INIT_VARS (myThid)
4868e56500 Andr*0003 c-----------------------------------------------------------------------
                0004 c  Routine to initialise the fizhi state.
eec2a90a83 Jean*0005 c
4868e56500 Andr*0006 c  Input: myThid       - Process number calling this routine
                0007 c
eec2a90a83 Jean*0008 c  Notes:
                0009 c   1) For a Cold Start -
4868e56500 Andr*0010 c      This routine takes the initial condition on the dynamics grid
                0011 c      and interpolates to the physics grid to initialize the state
                0012 c      variables that are on both grids. It initializes the variables
                0013 c      of the turbulence scheme to 0., and the land state from a model
                0014 c      climatology.
                0015 c   2) For a Restart, read the fizhi pickup file
                0016 c   3) The velocity component physics fields are on an A-Grid
                0017 c
                0018 c Calls: dyn2phys (x4)
                0019 c-----------------------------------------------------------------------
eec2a90a83 Jean*0020        IMPLICIT NONE
4868e56500 Andr*0021 #include "SIZE.h"
                0022 #include "fizhi_SIZE.h"
                0023 #include "fizhi_land_SIZE.h"
                0024 #include "GRID.h"
                0025 #include "DYNVARS.h"
                0026 #include "gridalt_mapping.h"
                0027 #include "fizhi_coms.h"
                0028 #include "fizhi_land_coms.h"
                0029 #include "fizhi_earth_coms.h"
                0030 #include "EEPARAMS.h"
                0031 #include "SURFACE.h"
                0032 #include "PARAMS.h"
                0033 #include "chronos.h"
                0034 
eec2a90a83 Jean*0035        INTEGER myThid
4868e56500 Andr*0036 
                0037 c pe on dynamics and physics grid refers to bottom edge
                0038        _RL pephy(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nrphys+1,nSx,nSy)
                0039        _RL pedyn(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1,nSx,nSy)
                0040        _RL windphy(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nrphys,nSx,nSy)
                0041        _RL udyntemp(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
                0042        _RL vdyntemp(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
                0043        _RL tempphy(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nrphys,nSx,nSy)
                0044 
eec2a90a83 Jean*0045        INTEGER i, j, L, bi, bj, Lbotij
                0046        INTEGER im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2
                0047        LOGICAL alarm
                0048        EXTERNAL alarm
4868e56500 Andr*0049 
                0050        im1 = 1-OLx
                0051        im2 = sNx+OLx
                0052        jm1 = 1-OLy
                0053        jm2 = sNy+OLy
                0054        idim1 = 1
                0055        idim2 = sNx
                0056        jdim1 = 1
                0057        jdim2 = sNy
                0058 
                0059 c   First Check to see if we can start a fizhi experiment at current time
                0060 c    All Fizhi alarms must be on for the first time step of a segment
                0061 
                0062       if( .not.alarm('moist') .or. .not.alarm('turb')   .or.
eec2a90a83 Jean*0063      &    .not.alarm('radsw') .or. .not.alarm('radlw') ) then
4868e56500 Andr*0064        print *,' Cant Start Fizhi experiment at ',nymd,' ',nhms
                0065        stop
                0066       endif
                0067 
                0068 C Deal Here with Variables that are on a Fizhi Pickup or need Initialization
                0069 
                0070       IF ( startTime.EQ.baseTime .AND. nIter0.EQ.0 ) THEN
                0071       print *,' In fizhi_init_vars: Beginning of New Experiment '
                0072 
                0073        do bj = myByLo(myThid), myByHi(myThid)
                0074        do bi = myBxLo(myThid), myBxHi(myThid)
                0075 
                0076 C Build pressures on dynamics grid
                0077         do j = 1,sNy
                0078         do i = 1,sNx
                0079          do L = 1,Nr
                0080           pedyn(i,j,L,bi,bj) = 0.
                0081          enddo
                0082         enddo
                0083         enddo
                0084         do j = 1,sNy
                0085         do i = 1,sNx
eec2a90a83 Jean*0086          Lbotij = kSurfC(i,j,bi,bj)
                0087          if(Lbotij.ne.0.)
                0088      &    pedyn(i,j,Lbotij,bi,bj) = Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj)
4868e56500 Andr*0089         enddo
                0090         enddo
                0091         do j = 1,sNy
                0092         do i = 1,sNx
eec2a90a83 Jean*0093          Lbotij = kSurfC(i,j,bi,bj)
4868e56500 Andr*0094          do L = Lbotij+1,Nr+1
                0095           pedyn(i,j,L,bi,bj) = pedyn(i,j,L-1,bi,bj) -
eec2a90a83 Jean*0096      &                        drF(L-1)*hfacC(i,j,L-1,bi,bj)
4868e56500 Andr*0097          enddo
                0098 c Do not use a zero field as the top edge pressure for interpolation
                0099          if(pedyn(i,j,Nr+1,bi,bj).lt.1.e-5)
eec2a90a83 Jean*0100      &                               pedyn(i,j,Nr+1,bi,bj) = 1.e-5
4868e56500 Andr*0101         enddo
                0102         enddo
                0103 C Build pressures on physics grid
                0104         do j = 1,sNy
                0105         do i = 1,sNx
                0106          pephy(i,j,1,bi,bj)=Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj)
                0107          do L = 2,Nrphys+1
                0108           pephy(i,j,L,bi,bj)=pephy(i,j,L-1,bi,bj)-dpphys0(i,j,L-1,bi,bj)
                0109          enddo
                0110 c Do not use a zero field as the top edge pressure for interpolation
                0111          if(pephy(i,j,Nrphys+1,bi,bj).lt.1.e-5)
eec2a90a83 Jean*0112      &                               pephy(i,j,Nrphys+1,bi,bj) = 1.e-5
4868e56500 Andr*0113         enddo
                0114         enddo
                0115 c
                0116 c Create an initial wind magnitude field on the physics grid -
                0117 c   Use a log wind law with z0=1cm, u*=1 cm/sec,
                0118 c   do units and get u = .025*ln(dP*10), with dP in pa.
                0119         do L = 1,Nrphys
                0120         do j = 1,sNy
                0121         do i = 1,sNx
                0122          windphy(i,j,L,bi,bj) = 0.025 *
eec2a90a83 Jean*0123      &             log((pephy(i,j,1,bi,bj)-pephy(i,j,L+1,bi,bj))*10.)
4868e56500 Andr*0124         enddo
                0125         enddo
                0126         enddo
                0127 
                0128        enddo
                0129        enddo
eec2a90a83 Jean*0130 
4868e56500 Andr*0131 c Create initial fields on phys. grid - Move Dynamics u and v to A-Grid
                0132        call CtoA(myThid,uvel,vvel,maskW,maskS,im1,im2,jm1,jm2,Nr,
eec2a90a83 Jean*0133      &                     nSx,nSy,1,sNx,1,sNy,udyntemp,vdyntemp)
4868e56500 Andr*0134 
                0135        do bj = myByLo(myThid), myByHi(myThid)
                0136        do bi = myBxLo(myThid), myBxHi(myThid)
                0137 
                0138 c Create initial fields on phys. grid - interpolate from dyn. grid
eec2a90a83 Jean*0139         call dyn2phys(udyntemp,pedyn,im1,im2,jm1,jm2,Nr,nSx,nSy,
                0140      & 1,sNx,1,sNy,bi,bj,windphy,pephy,kSurfC,Nrphys,nlperdyn,1,tempphy)
4868e56500 Andr*0141 c   Note: Interpolation gives bottom-up arrays (level 1 is bottom),
                0142 c         Physics works top-down. so -> need to flip arrays
                0143         do L = 1,Nrphys
                0144         do j = 1,sNy
                0145         do i = 1,sNx
                0146          uphy(i,j,Nrphys+1-L,bi,bj) = tempphy(i,j,L,bi,bj)
                0147         enddo
                0148         enddo
                0149         enddo
eec2a90a83 Jean*0150         call dyn2phys(vdyntemp,pedyn,im1,im2,jm1,jm2,Nr,nSx,nSy,
                0151      & 1,sNx,1,sNy,bi,bj,windphy,pephy,kSurfC,Nrphys,nlperdyn,1,tempphy)
4868e56500 Andr*0152         do L = 1,Nrphys
                0153         do j = 1,sNy
                0154         do i = 1,sNx
                0155          vphy(i,j,Nrphys+1-L,bi,bj) = tempphy(i,j,L,bi,bj)
                0156         enddo
                0157         enddo
                0158         enddo
eec2a90a83 Jean*0159         call dyn2phys(theta,pedyn,im1,im2,jm1,jm2,Nr,nSx,nSy,
                0160      & 1,sNx,1,sNy,bi,bj,windphy,pephy,kSurfC,Nrphys,nlperdyn,0,tempphy)
4868e56500 Andr*0161         do L = 1,Nrphys
                0162         do j = 1,sNy
                0163         do i = 1,sNx
                0164          thphy(i,j,Nrphys+1-L,bi,bj) = tempphy(i,j,L,bi,bj)
                0165         enddo
                0166         enddo
                0167         enddo
                0168 
eec2a90a83 Jean*0169         call dyn2phys(salt,pedyn,im1,im2,jm1,jm2,Nr,nSx,nSy,
                0170      & 1,sNx,1,sNy,bi,bj,windphy,pephy,kSurfC,Nrphys,nlperdyn,0,tempphy)
4868e56500 Andr*0171         do L = 1,Nrphys
                0172         do j = 1,sNy
                0173         do i = 1,sNx
                0174          sphy(i,j,Nrphys+1-L,bi,bj) = tempphy(i,j,L,bi,bj)
                0175         enddo
                0176         enddo
                0177         enddo
                0178 
                0179 c Zero out fizhi tendency arrays on the fizhi grid
                0180         do L = 1,Nrphys
                0181         do j = 1,sNy
                0182         do i = 1,sNx
                0183          duphy(i,j,L,bi,bj) = 0.
                0184          dvphy(i,j,L,bi,bj) = 0.
                0185          dthphy(i,j,L,bi,bj) = 0.
                0186          dsphy(i,j,L,bi,bj) = 0.
                0187         enddo
                0188         enddo
                0189         enddo
                0190 
                0191 c Zero out fizhi tendency arrays on the dynamics grid
                0192         do L = 1,Nr
                0193         do j = jm1,jm2
                0194         do i = im1,im2
                0195          guphy(i,j,L,bi,bj) = 0.
                0196          gvphy(i,j,L,bi,bj) = 0.
                0197          gthphy(i,j,L,bi,bj) = 0.
                0198          gsphy(i,j,L,bi,bj) = 0.
                0199         enddo
                0200         enddo
                0201         enddo
                0202 
                0203 c Initialize vegetation tile tke, xlmt, khmt, xxmt, yymt, ctmt, zetamt,
                0204         if( (nhms.eq.nhms0) .and. (nymd.eq.nymd0) ) then
                0205          print *,' Cold Start: Zero out Turb second moments '
                0206          do i = 1,nchp
                0207           ctmt(i,bi,bj) = 0.
                0208           xxmt(i,bi,bj) = 0.
                0209           yymt(i,bi,bj) = 0.
                0210           zetamt(i,bi,bj) = 0.
                0211          enddo
                0212          do L = 1,Nrphys
                0213          do i = 1,nchp
                0214           tke(i,L,bi,bj) = 0.
                0215           xlmt(i,L,bi,bj) = 0.
                0216           khmt(i,L,bi,bj) = 0.
                0217          enddo
                0218          enddo
                0219         else
                0220          print *,' Need initial Values for TKE - dont have them! '
                0221          stop
                0222         endif
eec2a90a83 Jean*0223         turbStart(bi,bj) = .TRUE.
4868e56500 Andr*0224 
eec2a90a83 Jean*0225 c Now initialize vegetation tile land state too - tcanopy, etc...
4868e56500 Andr*0226 c       call fizhi_init_vegsurftiles( nymd,nhms, 'D', myThid )
                0227 c Now initialize land state too - tcanopy, etc... SET FOR NOW,
                0228 c                                              READ CLIM FOR REAL
                0229         do i = 1,nchp
                0230         tcanopy(i,bi,bj) = 283.
                0231         tdeep(i,bi,bj) = 282.5
                0232         ecanopy(i,bi,bj) = 2.e-2
                0233         swetshal(i,bi,bj) = 0.6
                0234         swetroot(i,bi,bj) = 0.5
                0235         swetdeep(i,bi,bj) = 0.5
                0236         capac(i,bi,bj) = 0.
                0237         snodep(i,bi,bj) = 0.
                0238         enddo
                0239 
                0240 c Now initialize fizhi arrays that will be on a pickup
                0241         print *,' Initialize fizhi arrays that will be on pickup '
                0242         imstturblw(bi,bj) = 0
                0243         imstturbsw(bi,bj) = 0
                0244         iras(bi,bj) = 0
                0245         nlwcld(bi,bj) = 0
                0246         nlwlz(bi,bj) = 0
                0247         nswcld(bi,bj) = 0
                0248         nswlz(bi,bj) = 0
                0249         do L = 1,Nrphys
                0250         do j = 1,sNy
                0251         do i = 1,sNx
                0252          swlz(i,j,L,bi,bj) = 0.
                0253          lwlz(i,j,L,bi,bj) = 0.
                0254          qliqavesw(i,j,L,bi,bj) = 0.
                0255          qliqavelw(i,j,L,bi,bj) = 0.
                0256          fccavesw(i,j,L,bi,bj) = 0.
                0257          fccavelw(i,j,L,bi,bj) = 0.
                0258          cldtot_sw(i,j,L,bi,bj) = 0.
                0259          cldras_sw(i,j,L,bi,bj) = 0.
                0260          cldlsp_sw(i,j,L,bi,bj) = 0.
                0261          cldtot_lw(i,j,L,bi,bj) = 0.
                0262          cldras_lw(i,j,L,bi,bj) = 0.
                0263          cldlsp_lw(i,j,L,bi,bj) = 0.
                0264         enddo
                0265         enddo
                0266         enddo
                0267         do j = 1,sNy
                0268         do i = 1,sNx
                0269          rainlsp(i,j,bi,bj) = 0.
                0270          raincon(i,j,bi,bj) = 0.
                0271          snowfall(i,j,bi,bj) = 0.
                0272         enddo
                0273         enddo
                0274 
                0275        enddo
                0276        enddo
                0277 
                0278       ELSE
                0279       print *,' In fizhi_init_vars: Read from restart '
eec2a90a83 Jean*0280 
4868e56500 Andr*0281 C--   Read fizhi package state variables from pickup file
                0282 
                0283        call fizhi_read_pickup( nIter0, myThid )
                0284        CALL FIZHI_READ_VEGTILES( nIter0, 'D', myThid )
eec2a90a83 Jean*0285        do bj = myByLo(myThid), myByHi(myThid)
                0286        do bi = myBxLo(myThid), myBxHi(myThid)
                0287          turbStart(bi,bj) = .FALSE.
                0288        enddo
                0289        enddo
4868e56500 Andr*0290 
                0291       ENDIF
                0292 
eec2a90a83 Jean*0293       RETURN
                0294       END
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated from https://github.com/MITgcm/MITgcm by the 2.2.1-MITgcm-0.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!