Back to home page

MITgcm
 
 

    

File indexing completed on 2018-03-02 18:40:11 UTC

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