! !################################################################## !################################################################## !###### ###### !###### Advanced Regional Prediction System (ARPS) ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! PROGRAM arps,25 ! !----------------------------------------------------------------------- ! ! CONTACT: ! ! For further information, contact: ! ! Center for Analysis and Prediction of Storms ! University of Oklahoma ! Sarkeys Energy Center, ! 100 East Boyd ! Norman, OK 73019 USA ! Phone: (405) 325-0453 ! FAX : (405) 325-7614 ! E-mail: kkd@ou.edu or: mxue@ou.edu ! ! User support: arpssupport@ou.edu ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! This program is the driver for the Advanced Regional Prediction ! System (ARPS) model. ! ! The ARPS is a fully three-dimensional, compressible, ! nonhydrostatic model based on the Arakawa C grid. The model ! solves three momentum equations, a thermodynamic energy ! equation, a pressure equation, six moisture equations (water ! vapor, cloud water, rainwater, cloud ice, hail, and snow), ! and also provides for the statistical representation of ! unresolvable processes via a turbulence parameterization scheme ! (user option). The model is written in a terrain-following ! coordinate system, and has provisions for several boundary ! condition options. ! ! For a full description of ARPS, please consult the ARPS 4.0 ! User's Guide. ! !----------------------------------------------------------------------- ! ! MODIFICATION HISTORY: ! ! The modification history can be found in file HISTORY. ! !----------------------------------------------------------------------- ! ! GRID STRUCTURE AND STENCIL NUMBERING CONVENTION ! ! Each of the three velocity components is located at a physical ! boundary in the model. The numbering convention for the ! staggered grid is shown below for each of the coordinate ! directions, where S denotes a scalar variable. ! ! Arakawa C-grid is used by this model. ! ! ! X-Direction (East.......West): ! ! West boundary East boundary ! | | ! <-------- Physical Domain ------> ! +-----+-----+-----+-----+-----+-----+-----+-----+-----+-----+ ! U S U S U S U S U S U ! ! i= 1 1 2 2 . . . NX-2 NX-2 NX-1 NX-1 NX ! ! ! Y-Direction (North......South): ! ! South boundary North boundary ! | | ! <------- Physical Domain -------> ! +-----+-----+-----+-----+-----+-----+-----+-----+-----+-----+ ! V S V S V S V S V S V ! ! j= 1 1 2 2 . . . NY-2 NY-2 NY-1 NY-1 NY ! ! ! Z-Direction: ! ! NZ + W ! | ! | ! NZ-1 + S ! | ! | ! NZ-1 + W -- Physical Top boundary. ! | ! . ! . ! ! | ! 2 + S ! | ! | ! 2 + W -- Physical Ground Level ! | ! | ! 1 + S ! | ! | ! k= 1 + W ! !----------------------------------------------------------------------- ! ! Variable Declarations: ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! Force explicit declarations ! !----------------------------------------------------------------------- ! ! Include files: ! !----------------------------------------------------------------------- INCLUDE 'globcst.inc' ! Global constants that control model ! execution INCLUDE 'bndry.inc' INCLUDE 'timelvls.inc' INCLUDE 'radcst.inc' ! Includes radiation grid size information. INCLUDE 'mp.inc' ! Message passing parameters. INCLUDE 'alloc.inc' ! Memory allocation declaration and interfaces INCLUDE 'nudging.inc' ! Memory allocation declaration and interfaces ! !----------------------------------------------------------------------- ! ! Time dependent variables: ! !----------------------------------------------------------------------- REAL, ALLOCATABLE :: u (:,:,:,:) ! Total u-velocity (m/s) REAL, ALLOCATABLE :: v (:,:,:,:) ! Total v-velocity (m/s) REAL, ALLOCATABLE :: w (:,:,:,:) ! Total w-velocity (m/s) REAL, ALLOCATABLE :: wcont (:,:,:) ! Contravariant vertical velocity (m/s) REAL, ALLOCATABLE :: ptprt (:,:,:,:) ! Perturbation potential temperature ! from that of base state atmosphere (K) REAL, ALLOCATABLE :: pprt (:,:,:,:) ! Perturbation pressure from that ! of base state atmosphere (Pascal) REAL, ALLOCATABLE :: qv (:,:,:,:) ! Water vapor specific humidity (kg/kg) REAL, ALLOCATABLE :: qc (:,:,:,:) ! Cloud water mixing ratio (kg/kg) REAL, ALLOCATABLE :: qr (:,:,:,:) ! Rain water mixing ratio (kg/kg) REAL, ALLOCATABLE :: qi (:,:,:,:) ! Cloud ice mixing ratio (kg/kg) REAL, ALLOCATABLE :: qs (:,:,:,:) ! Snow mixing ratio (kg/kg) REAL, ALLOCATABLE :: qh (:,:,:,:) ! Hail mixing ratio (kg/kg) REAL, ALLOCATABLE :: tke (:,:,:,:) ! Turbulent Kinetic Energy ((m/s)**2) REAL, ALLOCATABLE :: pbldpth(:,:,:) ! Planetary boundary layer depth (m) REAL, ALLOCATABLE :: kmh (:,:,:) ! Horizontal turb. mixing coef. for ! momentum. ( m**2/s ) REAL, ALLOCATABLE :: kmv (:,:,:) ! Vertical turb. mixing coef. for ! momentum. ( m**2/s ) REAL, ALLOCATABLE :: rprntl(:,:,:) ! Reciprocal of Prandtl number ! !----------------------------------------------------------------------- ! ! Time tendencies of certain time dependent variables at the lateral ! boundaries, which are used to set the lateral boundary conditions ! for options 4 and 5. ! !----------------------------------------------------------------------- ! REAL, ALLOCATABLE :: udteb (:,:) ! Time tendency of u at e-boundary (m/s**2) REAL, ALLOCATABLE :: udtwb (:,:) ! Time tendency of u at w-boundary (m/s**2) REAL, ALLOCATABLE :: udtnb (:,:) ! Time tendency of u at n-boundary (m/s**2) REAL, ALLOCATABLE :: udtsb (:,:) ! Time tendency of u at s-boundary (m/s**2) REAL, ALLOCATABLE :: vdteb (:,:) ! Time tendency of v at e-boundary (m/s**2) REAL, ALLOCATABLE :: vdtwb (:,:) ! Time tendency of v at w-boundary (m/s**2) REAL, ALLOCATABLE :: vdtnb (:,:) ! Time tendency of v at n-boundary (m/s**2) REAL, ALLOCATABLE :: vdtsb (:,:) ! Time tendency of v at s-boundary (m/s**2) REAL, ALLOCATABLE :: wdteb (:,:) ! Time tendency of w at e-boundary (m/s**2) REAL, ALLOCATABLE :: wdtwb (:,:) ! Time tendency of w at w-boundary (m/s**2) REAL, ALLOCATABLE :: wdtnb (:,:) ! Time tendency of w at n-boundary (m/s**2) REAL, ALLOCATABLE :: wdtsb (:,:) ! Time tendency of w at s-boundary (m/s**2) REAL, ALLOCATABLE :: pdteb (:,:) ! Time tendency of pprt at e-boundary (Pascal/s) REAL, ALLOCATABLE :: pdtwb (:,:) ! Time tendency of pprt at w-boundary (Pascal/s) REAL, ALLOCATABLE :: pdtnb (:,:) ! Time tendency of pprt at n-boundary (Pascal/s) REAL, ALLOCATABLE :: pdtsb (:,:) ! Time tendency of pprt at s-boundary (Pascal/s) REAL, ALLOCATABLE :: sdteb (:,:) ! Time tendency of a scalar at e-boundary (m/s**2) REAL, ALLOCATABLE :: sdtwb (:,:) ! Time tendency of a scalar at w-boundary (m/s**2) REAL, ALLOCATABLE :: sdtnb (:,:) ! Time tendency of a scalar at n-boundary (m/s**2) REAL, ALLOCATABLE :: sdtsb (:,:) ! Time tendency of a scalar at s-boundary (m/s**2) ! !----------------------------------------------------------------------- ! ! Base state variables: ! !----------------------------------------------------------------------- ! REAL, ALLOCATABLE :: ubar (:,:,:) ! Base state u-velocity (m/s) REAL, ALLOCATABLE :: vbar (:,:,:) ! Base state v-velocity (m/s) REAL, ALLOCATABLE :: ptbar (:,:,:) ! Base state potential temperature (K) REAL, ALLOCATABLE :: pbar (:,:,:) ! Base state pressure (Pascal) REAL, ALLOCATABLE :: ptbari(:,:,:) ! Inverse Base state pot. temperature (K) REAL, ALLOCATABLE :: pbari (:,:,:) ! Inverse Base state pressure (Pascal) REAL, ALLOCATABLE :: rhostr(:,:,:) ! Base state density rhobar times j3. REAL, ALLOCATABLE :: rhostri(:,:,:)! Inv. base state density rhobar times j3. REAL, ALLOCATABLE :: qvbar (:,:,:) ! Base state water vapor specific humidity (kg/kg) REAL, ALLOCATABLE :: ppi (:,:,:) ! Exner function. REAL, ALLOCATABLE :: csndsq(:,:,:) ! Sound wave speed squared. ! !----------------------------------------------------------------------- ! ! Arrays related to model grid definition: ! !----------------------------------------------------------------------- ! REAL, ALLOCATABLE :: x (:) ! The x-coord. of the physical and ! computational grid. Defined at u-point. REAL, ALLOCATABLE :: y (:) ! The y-coord. of the physical and ! computational grid. Defined at v-point. REAL, ALLOCATABLE :: z (:) ! The z-coord. of the computational grid. ! Defined at w-point on the staggered grid. REAL, ALLOCATABLE :: zp (:,:,:) ! The physical height coordinate defined at ! w-point of the staggered grid. REAL, ALLOCATABLE :: hterain(:,:) ! The height of the terrain. REAL, ALLOCATABLE :: mapfct(:,:,:) ! Map factors at scalar, u and v points REAL, ALLOCATABLE :: j1 (:,:,:) ! Coordinate transformation Jacobian defined ! as - d( zp )/d( x ). REAL, ALLOCATABLE :: j2 (:,:,:) ! Coordinate transformation Jacobian defined ! as - d( zp )/d( y ). REAL, ALLOCATABLE :: j3 (:,:,:) ! Coordinate transformation Jacobian defined ! as d( zp )/d( z ). REAL, ALLOCATABLE :: aj3x (:,:,:) ! Coordinate transformation Jacobian defined ! as d( zp )/d( z ) AVERAGED IN THE X-DIR. REAL, ALLOCATABLE :: aj3y (:,:,:) ! Coordinate transformation Jacobian defined ! as d( zp )/d( z ) AVERAGED IN THE Y-DIR. REAL, ALLOCATABLE :: aj3z (:,:,:) ! Coordinate transformation Jacobian defined ! as d( zp )/d( z ) AVERAGED IN THE Z-DIR. REAL, ALLOCATABLE :: j3inv (:,:,:) ! Inverse of J3. REAL, ALLOCATABLE :: trigs1(:) ! Array containing pre-computed trig ! function for fftopt=1. REAL, ALLOCATABLE :: trigs2(:) ! Array containing pre-computed trig ! function for fftopt=1. INTEGER, ALLOCATABLE :: ifax1(:) ! Array containing the factors of nx for ! fftopt=1. INTEGER, ALLOCATABLE :: ifax2(:) ! Array containing the factors of ny for ! fftopt=1. REAL, ALLOCATABLE :: vwork1 (:,:) ! 2-D work array for fftopt=2. REAL, ALLOCATABLE :: vwork2 (:,:) ! 2-D work array for fftopt=2. REAL, ALLOCATABLE :: wsave1 (:) ! Work array for fftopt =2. REAL, ALLOCATABLE :: wsave2 (:) ! Work array for fftopt =2. REAL, ALLOCATABLE :: sinlat(:,:) ! Sin of latitude at each grid point REAL, ALLOCATABLE :: ptcumsrc(:,:,:) ! Source term in pt-equation due ! to cumulus parameterization REAL, ALLOCATABLE :: qcumsrc(:,:,:,:)! Source term in water equations due ! to cumulus parameterization: ! qcumsrc(1,1,1,1) for qv equation ! qcumsrc(1,1,1,2) for qc equation ! qcumsrc(1,1,1,3) for qr equation ! qcumsrc(1,1,1,4) for qi equation ! qcumsrc(1,1,1,5) for qs equation REAL, ALLOCATABLE :: w0avg(:,:,:) ! a closing running average vertical ! velocity in 10min for K-F scheme ! !----------------------------------------------------------------------- ! ! Arrays for soil model. ! !----------------------------------------------------------------------- ! INTEGER, ALLOCATABLE :: soiltyp(:,:,:) ! Soil type at each point REAL, ALLOCATABLE :: stypfrct(:,:,:) ! Fraction of soil type INTEGER, ALLOCATABLE :: vegtyp (:,:) ! Vegetation type at each point REAL, ALLOCATABLE :: lai (:,:) ! Leaf Area Index REAL, ALLOCATABLE :: roufns (:,:) ! Surface roughness REAL, ALLOCATABLE :: veg (:,:) ! Vegetation fraction REAL, ALLOCATABLE :: tsfc (:,:,:) ! Ground temperature at surface (K) REAL, ALLOCATABLE :: tsoil (:,:,:) ! Deep soil temperature (K) REAL, ALLOCATABLE :: wetsfc (:,:,:) ! Surface soil moisture REAL, ALLOCATABLE :: wetdp (:,:,:) ! Deep soil moisture REAL, ALLOCATABLE :: wetcanp(:,:,:) ! Canopy water amount REAL, ALLOCATABLE :: snowdpth(:,:) ! Snow depth (m) REAL, ALLOCATABLE :: qvsfc (:,:,:) ! Effective sfc. qv REAL, ALLOCATABLE :: ptsfc (:,:) ! Ground surface potential temperature (K) REAL, ALLOCATABLE :: raing(:,:) ! Grid supersaturation rain REAL, ALLOCATABLE :: rainc(:,:) ! Cumulus convective rain REAL, ALLOCATABLE :: prcrate(:,:,:) ! precipitation rates (kg/(m**2*s)) ! prcrate(1,1,1) = total precip. rate ! prcrate(1,1,2) = grid scale precip. rate ! prcrate(1,1,3) = cumulative precip. rate ! prcrate(1,1,4) = microphysics precip. rate REAL, ALLOCATABLE :: raincv(:,:) ! K-F convective rainfall (cm) INTEGER, ALLOCATABLE :: nca(:,:) ! K-F counter for CAPE release ! !----------------------------------------------------------------------- ! ! Arrays for radiation ! !----------------------------------------------------------------------- ! REAL, ALLOCATABLE :: radfrc(:,:,:) ! Radiation forcing (K/s) REAL, ALLOCATABLE :: radsw (:,:) ! Solar radiation reaching the surface REAL, ALLOCATABLE :: rnflx (:,:) ! Net radiation flux absorbed by surface REAL, ALLOCATABLE :: rad2d (:,:,:) ! 2-D arrays for radiation (see radcst.inc) REAL, ALLOCATABLE :: radbuf(:) ! Buffer to carry working arrays for ! radiation computing (see radcst.inc) ! !----------------------------------------------------------------------- ! ! Arrays that carry surface fluxes ! !----------------------------------------------------------------------- ! REAL, ALLOCATABLE :: usflx (:,:) ! Surface flux of u-momentum (kg/(m*s**2)) REAL, ALLOCATABLE :: vsflx (:,:) ! Surface flux of v-momentum (kg/(m*s**2)) REAL, ALLOCATABLE :: ptsflx(:,:) ! Surface heat flux (K*kg/(m**2*s)) REAL, ALLOCATABLE :: qvsflx(:,:) ! Surface moisture flux (kg/(m**2*s)) ! !----------------------------------------------------------------------- ! ! Buffer arrays that carry external boundary 3-D arrays ! !----------------------------------------------------------------------- ! REAL, ALLOCATABLE :: exbcbuf( : ) ! EXBC buffer array REAL, ALLOCATABLE :: bcrlx(:,:) ! EXBC relaxation function coefficients ! !----------------------------------------------------------------------- ! ! Arrays for analysis increment updating (a type of nudging) output ! !----------------------------------------------------------------------- ! REAL, ALLOCATABLE :: uincr(:,:,:) REAL, ALLOCATABLE :: vincr(:,:,:) REAL, ALLOCATABLE :: wincr(:,:,:) REAL, ALLOCATABLE :: pincr(:,:,:) REAL, ALLOCATABLE :: ptincr(:,:,:) REAL, ALLOCATABLE :: qvincr(:,:,:) REAL, ALLOCATABLE :: qcincr(:,:,:) REAL, ALLOCATABLE :: qrincr(:,:,:) REAL, ALLOCATABLE :: qiincr(:,:,:) REAL, ALLOCATABLE :: qsincr(:,:,:) REAL, ALLOCATABLE :: qhincr(:,:,:) ! !----------------------------------------------------------------------- ! ! Work arrays that do not carry physical meaning in the code ! !----------------------------------------------------------------------- ! REAL, ALLOCATABLE :: temxy1(:,:) ! Temporary work array, used as phydro REAL, ALLOCATABLE :: tem1 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem2 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem3 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem4 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem5 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem6 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem7 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem8 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem9 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem10 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem11 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem12 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem13 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem14 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem15 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem16 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem17 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem18 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem19 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem20 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem21 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem22 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem23 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem24 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem25 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem26 (:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem1_0(:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem2_0(:,:,:) ! Temporary work array. REAL, ALLOCATABLE :: tem3_0(:,:,:) ! Temporary work array. ! !----------------------------------------------------------------------- ! ! ARPS dimensions: ! ! nx, ny, nz: Dimensions of computatioal grid. When run on MPP ! with PVM or MPI, they represent of the size of the ! subdomains. See below. ! ! Given nx, ny and nz, the physical domain size will be ! xl=(nx-3)*dx by yl=(ny-3)*dy by zh=(nz-3)*dz. The ! variables nx, ny, nz, dx, dy and dz are read in from ! the input file by subroutine INITPARA. ! !----------------------------------------------------------------------- ! INTEGER :: nx ! Number of grid points in the x-direction INTEGER :: ny ! Number of grid points in the y-direction INTEGER :: nz ! Number of grid points in the z-direction INTEGER :: nxndg ! Number of x grid points for nudging (1 or nx) INTEGER :: nyndg ! Number of y grid points for nudging (1 or ny) INTEGER :: nzndg ! Number of z grid points for nudging (1 or nz) ! !----------------------------------------------------------------------- ! ! Soil types. ! !----------------------------------------------------------------------- ! INTEGER :: nstyps ! Number of soil types ! !----------------------------------------------------------------------- ! ! Radiation buffer size. Set equal to ! n2d_radiat*nx*ny + n3d_radiat*nx*ny*nz ! in subroutine inipara if radopt=2, otherwise set to 1. ! !----------------------------------------------------------------------- ! INTEGER :: rbufsz ! !----------------------------------------------------------------------- ! ! External boundary buffer size. Set equal to ! 22*nx*ny*nz ! in subroutine inipara if lbcopt=2, otherwise set to 1. ! !----------------------------------------------------------------------- ! INTEGER :: exbcbufsz INTEGER :: frstep ! Flag for the initial time step INTEGER :: mptr ! Grid identifier. INTEGER :: ierr, i,j,k REAL :: eps ! Small value to compensate the roundoff ! error in checking of the end of time ! integration. DATA eps /0.01/ INTEGER istatus, nxy, nxyz REAL :: twall0, t_walltime ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! !----------------------------------------------------------------------- ! ! First, initialize the model parameters. Most of them are ! global parameters passed among subroutines through common blocks. ! ! These parameters are declared in the include files: ! globcst.inc, bndry.inc, exbc.inc, radcst.inc etc. ! !----------------------------------------------------------------------- ! CALL acct_init CALL set_acct(init_acct) CALL initpara(nx,ny,nz,nstyps) ! Set radiation and external boundary work array buffer sizes: IF (radopt == 2) THEN rbufsz = n2d_radiat*nx*ny + n3d_radiat*nx*ny*nz ELSE rbufsz = 1 END IF IF (lbcopt == 2) THEN exbcbufsz = 22*nx*ny*nz ELSE exbcbufsz = 1 END IF IF ( nudgopt > 0 ) THEN nxndg=nx nyndg=ny nzndg=nz ELSE nxndg=1 nyndg=1 nzndg=1 END IF ! !----------------------------------------------------------------------- ! ! Print a log file in namelist format. ! !----------------------------------------------------------------------- ! CALL acct_interrupt(output_acct) IF (myproc == 0) THEN IF (mp_opt == 0) THEN CALL prtlog(nx,ny,nz,6) CALL prtlog(nx,ny,nz,0) ELSE CALL prtlog(nproc_x*(nx-3)+3,nproc_y*(ny-3)+3,nz,6) CALL prtlog(nproc_x*(nx-3)+3,nproc_y*(ny-3)+3,nz,0) ENDIF END IF CALL acct_stop_inter ! !----------------------------------------------------------------------- ! ! Allocate arrays etc. ! !----------------------------------------------------------------------- ! mptr = 1 ! Set the grid number to 1 for a single grid run. nestgrd = 0 ! Set the grid nesting flag to zero, i.e. no nesting. nxy = nx*ny nxyz = nxy*nz ALLOCATE(u (nx,ny,nz,nt),STAT=istatus) u = 0 ALLOCATE(v (nx,ny,nz,nt),STAT=istatus) v = 0 ALLOCATE(w (nx,ny,nz,nt),STAT=istatus) w = 0 ALLOCATE(wcont(nx,ny,nz ),STAT=istatus) wcont = 0 ALLOCATE(ptprt(nx,ny,nz,nt),STAT=istatus) ptprt = 0 ALLOCATE(pprt (nx,ny,nz,nt),STAT=istatus) pprt = 0 ALLOCATE(qv (nx,ny,nz,nt),STAT=istatus) qv = 0 ALLOCATE(qc (nx,ny,nz,nt),STAT=istatus) qc = 0 ALLOCATE(qr (nx,ny,nz,nt),STAT=istatus) qr = 0 ALLOCATE(qi (nx,ny,nz,nt),STAT=istatus) qi = 0 ALLOCATE(qs (nx,ny,nz,nt),STAT=istatus) qs = 0 ALLOCATE(qh (nx,ny,nz,nt),STAT=istatus) qh = 0 ALLOCATE(tke (nx,ny,nz,nt),STAT=istatus) tke = 0 ALLOCATE(pbldpth(nx,ny,nt),STAT=istatus) pbldpth = 0 ALLOCATE(kmh(nx,ny,nz),STAT=istatus) kmh = 0 ALLOCATE(kmv(nx,ny,nz),STAT=istatus) kmv = 0 ALLOCATE(rprntl(nx,ny,nz),STAT=istatus) rprntl = 0 ALLOCATE(udteb(ny,nz),STAT=istatus) udteb = 0 ALLOCATE(udtwb(ny,nz),STAT=istatus) udtwb = 0 ALLOCATE(udtnb(nx,nz),STAT=istatus) udtnb = 0 ALLOCATE(udtsb(nx,nz),STAT=istatus) udtsb = 0 ALLOCATE(vdteb(ny,nz),STAT=istatus) vdteb = 0 ALLOCATE(vdtwb(ny,nz),STAT=istatus) vdtwb = 0 ALLOCATE(vdtnb(nx,nz),STAT=istatus) vdtnb = 0 ALLOCATE(vdtsb(nx,nz),STAT=istatus) vdtsb = 0 ALLOCATE(wdteb(ny,nz),STAT=istatus) wdteb = 0 ALLOCATE(wdtwb(ny,nz),STAT=istatus) wdtwb = 0 ALLOCATE(wdtnb(nx,nz),STAT=istatus) wdtnb = 0 ALLOCATE(wdtsb(nx,nz),STAT=istatus) wdtsb = 0 ALLOCATE(pdteb(ny,nz),STAT=istatus) pdteb = 0 ALLOCATE(pdtwb(ny,nz),STAT=istatus) pdtwb = 0 ALLOCATE(pdtnb(nx,nz),STAT=istatus) pdtnb = 0 ALLOCATE(pdtsb(nx,nz),STAT=istatus) pdtsb = 0 ALLOCATE(sdteb(ny,nz),STAT=istatus) sdteb = 0 ALLOCATE(sdtwb(ny,nz),STAT=istatus) sdtwb = 0 ALLOCATE(sdtnb(nx,nz),STAT=istatus) sdtnb = 0 ALLOCATE(sdtsb(nx,nz),STAT=istatus) sdtsb = 0 ALLOCATE(ubar (nx,ny,nz),STAT=istatus) ubar = 0 ALLOCATE(vbar (nx,ny,nz),STAT=istatus) vbar = 0 ALLOCATE(ptbar (nx,ny,nz),STAT=istatus) ptbar = 0 ALLOCATE(pbar (nx,ny,nz),STAT=istatus) pbar = 0 ALLOCATE(ptbari(nx,ny,nz),STAT=istatus) ptbari = 0 ALLOCATE(pbari (nx,ny,nz),STAT=istatus) pbari = 0 ALLOCATE(rhostr(nx,ny,nz),STAT=istatus) rhostr = 0 ALLOCATE(rhostri(nx,ny,nz),STAT=istatus) rhostri = 0 ALLOCATE(qvbar (nx,ny,nz),STAT=istatus) qvbar = 0 ALLOCATE(ppi (nx,ny,nz),STAT=istatus) ppi = 0 ALLOCATE(csndsq(nx,ny,nz),STAT=istatus) csndsq = 0 ALLOCATE( x = 0 ALLOCATE( y = 0 ALLOCATE( z = 0 ALLOCATE(zp(nx,ny,nz),STAT=istatus) zp = 0 ALLOCATE(hterain(nx,ny),STAT=istatus) hterain = 0 ALLOCATE(mapfct(nx,ny,8),STAT=istatus) mapfct = 0 ALLOCATE(j1 (nx,ny,nz),STAT=istatus) j1 = 0 ALLOCATE(j2 (nx,ny,nz),STAT=istatus) j2 = 0 ALLOCATE(j3 (nx,ny,nz),STAT=istatus) j3 = 0 ALLOCATE(aj3x (nx,ny,nz),STAT=istatus) aj3x = 0 ALLOCATE(aj3y (nx,ny,nz),STAT=istatus) aj3y = 0 ALLOCATE(aj3z (nx,ny,nz),STAT=istatus) aj3z = 0 ALLOCATE(j3inv(nx,ny,nz),STAT=istatus) j3inv = 0 ALLOCATE(trigs1(3*(nx-1)/2+1),STAT=istatus) trigs1 = 0 ALLOCATE(trigs2(3*(ny-1)/2+1),STAT=istatus) trigs2 = 0 ALLOCATE(ifax1(13),STAT=istatus) ifax1 = 0 ALLOCATE(ifax2(13),STAT=istatus) ifax2 = 0 ALLOCATE(vwork1(nx+1,ny+1),STAT=istatus) vwork1 = 0 ALLOCATE(vwork2(ny,nx+1),STAT=istatus) vwork2 = 0 ALLOCATE(wsave1(3*(ny-1)+15),STAT=istatus) wsave1 = 0 ALLOCATE(wsave2(3*(nx-1)+15),STAT=istatus) wsave2 = 0 ALLOCATE(sinlat(nx,ny),STAT=istatus) sinlat = 0 ALLOCATE(ptcumsrc(nx,ny,nz),STAT=istatus) ptcumsrc = 0 ALLOCATE(qcumsrc(nx,ny,nz,5),STAT=istatus) qcumsrc = 0 ALLOCATE(w0avg(nx,ny,nz),STAT=istatus) w0avg = 0 ALLOCATE(soiltyp (nx,ny,nstyps),STAT=istatus) soiltyp = 0 ALLOCATE(stypfrct(nx,ny,nstyps),STAT=istatus) stypfrct = 0 ALLOCATE(vegtyp(nx,ny),STAT=istatus) vegtyp = 0 ALLOCATE(lai (nx,ny),STAT=istatus) lai = 0 ALLOCATE(roufns(nx,ny),STAT=istatus) roufns = 0 ALLOCATE(veg (nx,ny),STAT=istatus) veg = 0 ALLOCATE(tsfc (nx,ny,0:nstyps),STAT=istatus) tsfc = 0 ALLOCATE(tsoil (nx,ny,0:nstyps),STAT=istatus) tsoil = 0 ALLOCATE(wetsfc (nx,ny,0:nstyps),STAT=istatus) wetsfc = 0 ALLOCATE(wetdp (nx,ny,0:nstyps),STAT=istatus) wetdp = 0 ALLOCATE(wetcanp (nx,ny,0:nstyps),STAT=istatus) wetcanp = 0 ALLOCATE(qvsfc (nx,ny,0:nstyps),STAT=istatus) qvsfc = 0 ALLOCATE(ptsfc (nx,ny),STAT=istatus) ptsfc = 0 ALLOCATE(snowdpth(nx,ny),STAT=istatus) snowdpth = 0 ALLOCATE(raing (nx,ny),STAT=istatus) raing = 0 ALLOCATE(rainc (nx,ny),STAT=istatus) rainc = 0 ALLOCATE(raincv (nx,ny),STAT=istatus) raincv = 0 ALLOCATE(nca (nx,ny),STAT=istatus) nca = 0 ALLOCATE(prcrate(nx,ny,4),STAT=istatus) prcrate = 0 ALLOCATE(radfrc(nx,ny,nz),STAT=istatus) radfrc = 0 ALLOCATE(radsw (nx,ny),STAT=istatus) radsw = 0 ALLOCATE(rnflx (nx,ny),STAT=istatus) rnflx = 0 ALLOCATE(rad2d (nx,ny,nrad2d),STAT=istatus) rad2d = 0 ALLOCATE(radbuf(rbufsz),STAT=istatus) radbuf = 0 ALLOCATE(usflx (nx,ny),STAT=istatus) usflx = 0 ALLOCATE(vsflx (nx,ny),STAT=istatus) vsflx = 0 ALLOCATE(ptsflx(nx,ny),STAT=istatus) ptsflx = 0 ALLOCATE(qvsflx(nx,ny),STAT=istatus) qvsflx = 0 ALLOCATE(exbcbuf(exbcbufsz ),STAT=istatus) exbcbuf = 0 ALLOCATE(bcrlx(nx,ny),STAT=istatus) bcrlx = 0 print *, ' nxndg,nyndg,nzndg = ',nxndg,nyndg,nzndg ALLOCATE(uincr(nxndg,nyndg,nzndg),STAT=istatus) uincr = 0. ALLOCATE(vincr(nxndg,nyndg,nzndg),STAT=istatus) vincr = 0. ALLOCATE(wincr(nxndg,nyndg,nzndg),STAT=istatus) wincr = 0. ALLOCATE(pincr(nxndg,nyndg,nzndg),STAT=istatus) pincr = 0. ALLOCATE(ptincr(nxndg,nyndg,nzndg),STAT=istatus) ptincr = 0. ALLOCATE(qvincr(nxndg,nyndg,nzndg),STAT=istatus) qvincr = 0. ALLOCATE(qcincr(nxndg,nyndg,nzndg),STAT=istatus) qcincr = 0. ALLOCATE(qrincr(nxndg,nyndg,nzndg),STAT=istatus) qrincr = 0. ALLOCATE(qiincr(nxndg,nyndg,nzndg),STAT=istatus) qiincr = 0. ALLOCATE(qsincr(nxndg,nyndg,nzndg),STAT=istatus) qsincr = 0. ALLOCATE(qhincr(nxndg,nyndg,nzndg),STAT=istatus) qhincr = 0. ALLOCATE(temxy1(nx,ny),STAT=istatus) temxy1 = 0 ALLOCATE(tem1 (nx,ny,nz),STAT=istatus) tem1 = 0 ALLOCATE(tem2 (nx,ny,nz),STAT=istatus) tem2 = 0 ALLOCATE(tem3 (nx,ny,nz),STAT=istatus) tem3 = 0 ALLOCATE(tem4 (nx,ny,nz),STAT=istatus) tem4 = 0 ALLOCATE(tem5 (nx,ny,nz),STAT=istatus) tem5 = 0 ALLOCATE(tem6 (nx,ny,nz),STAT=istatus) tem6 = 0 ALLOCATE(tem7 (nx,ny,nz),STAT=istatus) tem7 = 0 ALLOCATE(tem8 (nx,ny,nz),STAT=istatus) tem8 = 0 ALLOCATE(tem9 (nx,ny,nz),STAT=istatus) tem9 = 0 ALLOCATE(tem10(nx,ny,nz),STAT=istatus) tem10 = 0 ALLOCATE(tem11(nx,ny,nz),STAT=istatus) tem11 = 0 ALLOCATE(tem12(nx,ny,nz),STAT=istatus) tem12 = 0 ALLOCATE(tem13(nx,ny,nz),STAT=istatus) tem13 = 0 ALLOCATE(tem14(nx,ny,nz),STAT=istatus) tem14 = 0 ALLOCATE(tem15(nx,ny,nz),STAT=istatus) tem15 = 0 ALLOCATE(tem16(nx,ny,nz),STAT=istatus) tem16 = 0 ALLOCATE(tem17(nx,ny,nz),STAT=istatus) tem17 = 0 ALLOCATE(tem18(nx,ny,nz),STAT=istatus) tem18 = 0 ALLOCATE(tem19(nx,ny,nz),STAT=istatus) tem19 = 0 ALLOCATE(tem20(nx,ny,nz),STAT=istatus) tem20 = 0 ALLOCATE(tem21(nx,ny,nz),STAT=istatus) tem21 = 0 ALLOCATE(tem22(nx,ny,nz),STAT=istatus) tem22 = 0 ALLOCATE(tem23(nx,ny,nz),STAT=istatus) tem23 = 0 ALLOCATE(tem24(nx,ny,nz),STAT=istatus) tem24 = 0 ALLOCATE(tem25(nx,ny,nz),STAT=istatus) tem25 = 0 ALLOCATE(tem26(nx,ny,nz),STAT=istatus) tem26 = 0 ALLOCATE(tem1_0(0:nx,0:ny,0:nz),STAT=istatus) tem1_0 = 0 ALLOCATE(tem2_0(0:nx,0:ny,0:nz),STAT=istatus) tem2_0 = 0 ALLOCATE(tem3_0(0:nx,0:ny,0:nz),STAT=istatus) tem3_0 = 0 ! !----------------------------------------------------------------------- ! ! Set up exbcbuf pointers. ! !----------------------------------------------------------------------- ! IF( lbcopt == 2 .AND. mptr == 1 ) CALL setexbcptr(nx,ny,nz) ! !----------------------------------------------------------------------- ! ! Initialize all arrays and variables, and set initial ! condition for thermal disturbance if desired. ! !----------------------------------------------------------------------- ! CALL initial(mptr,nx,ny,nz,nxndg,nyndg,nzndg,nstyps,exbcbufsz, & u,v,w,wcont,ptprt,pprt,qv,qc,qr,qi,qs,qh,tke, & udteb,udtwb,udtnb,udtsb,vdteb,vdtwb,vdtnb,vdtsb, & wdteb,wdtwb,wdtnb,wdtsb,pdteb,pdtwb,pdtnb,pdtsb, & sdteb,sdtwb,sdtnb,sdtsb, & ubar,vbar,ptbar,pbar,ptbari,pbari,rhostr,rhostri, & qvbar,ppi,csndsq, & x,y,z,zp,hterain,mapfct, & j1,j2,j3,aj3x,aj3y,aj3z,j3inv, & trigs1,trigs2,ifax1,ifax2, & wsave1,wsave2,vwork1,vwork2, & sinlat,kmh,kmv,rprntl, & soiltyp,stypfrct,vegtyp,lai,roufns,veg, & tsfc,tsoil,wetsfc,wetdp,wetcanp,snowdpth,ptsfc,qvsfc, & ptcumsrc,qcumsrc,w0avg,nca, & raincv,raing,rainc,prcrate,exbcbuf,bcrlx,radfrc,radsw, & rnflx,usflx,vsflx,ptsflx,qvsflx, & uincr,vincr,wincr,pincr,ptincr,qvincr, & qcincr,qrincr,qiincr,qsincr,qhincr, & temxy1,tem1,tem2,tem3,tem4,tem5,tem6,tem7, & tem8,tem9,tem10,tem11,tem12,tem13, & tem14,tem15,tem16,tem17,tem18,tem19, & tem20,tem21,tem22,tem23,tem24,tem25,tem26, & tem1_0,tem2_0,tem3_0) ! !----------------------------------------------------------------------- ! ! The current model time (initial time): ! !----------------------------------------------------------------------- ! curtim = tstart IF( myproc == 0 ) WRITE(6,'(1x,a,f13.3,a)') & 'The initial model time is at ', curtim,' seconds.' ! !----------------------------------------------------------------------- ! ! Output the initial fields ! !----------------------------------------------------------------------- ! CALL set_acct(output_acct) CALL initout(mptr,nx,ny,nz,nstyps,exbcbufsz, & u,v,w,wcont,ptprt,pprt,qv,qc,qr,qi,qs,qh,tke, & ubar,vbar,ptbar,pbar,rhostr,qvbar,kmh,kmv, & x,y,z,zp,hterain,mapfct,j1,j2,j3,j3inv, & soiltyp,stypfrct,vegtyp,lai,roufns,veg, & tsfc,tsoil,wetsfc,wetdp,wetcanp,snowdpth, & raing,rainc,prcrate,exbcbuf, & radfrc,radsw,rnflx, & usflx,vsflx,ptsflx,qvsflx, & tem1,tem2,tem3,tem4,tem5,tem6,tem7, & tem8,tem9,tem10,tem11) IF (mp_opt > 0) THEN CALL mpbarrier () END IF CALL set_acct(misc_acct) nstep = nint( curtim/dtbig ) ! !----------------------------------------------------------------------- ! ! Time integration loop begins -----------------------------------> ! !----------------------------------------------------------------------- ! 900 CONTINUE ! Entry point of time step integration. ! !----------------------------------------------------------------------- ! ! Define the time step counter nstep reference to time zero ! (curtim=0.0). This means for tstart.ne.0 case, nstep for the ! first step of integration is not 1. ! !----------------------------------------------------------------------- ! nstep = nstep + 1 IF( (ABS(curtim-tstart) <= 1.0E-10) .AND. (restrt /= 1) ) THEN frstep=1 ! Indicate that this is the initial step of ! model integration. For this step, a forward ! time integration scheme will be used. ELSE ! For non-first step or restart run frstep=0 ! When frstep=0, leapfrog scheme is used. END IF ! !----------------------------------------------------------------------- ! ! Perform one time step integration for all equations: ! On exit of this routine, all time dependent fields are ! advanced by one time step. ! !----------------------------------------------------------------------- ! CALL cordintg(mptr,frstep, nx,ny,nz,nxndg,nyndg,nzndg, & rbufsz,nstyps,exbcbufsz, & u,v,w,wcont,ptprt,pprt,qv,qc,qr,qi,qs,qh, & tke,pbldpth, & udteb,udtwb,udtnb,udtsb,vdteb,vdtwb,vdtnb,vdtsb, & wdteb,wdtwb,wdtnb,wdtsb,pdteb,pdtwb,pdtnb,pdtsb, & sdteb,sdtwb,sdtnb,sdtsb, & ubar,vbar,ptbar,pbar,ptbari,pbari,rhostr,rhostri, & qvbar,ppi,csndsq, & x,y,z,zp, mapfct, & j1,j2,j3,aj3x,aj3y,aj3z,j3inv, & trigs1,trigs2,ifax1,ifax2, & wsave1,wsave2,vwork1,vwork2, & sinlat, kmh,kmv,rprntl, & soiltyp,stypfrct,vegtyp,lai,roufns,veg, & tsfc,tsoil,wetsfc,wetdp,wetcanp, & snowdpth,ptsfc,qvsfc, & ptcumsrc,qcumsrc,raing,rainc,prcrate, & w0avg,nca,raincv, & radfrc,radsw,rnflx,rad2d,radbuf,exbcbuf,bcrlx, & usflx,vsflx,ptsflx,qvsflx, & uincr,vincr,wincr,pincr,ptincr,qvincr, & qcincr,qrincr,qiincr,qsincr,qhincr, & temxy1,tem1,tem2,tem3,tem4,tem5,tem6,tem7, & tem8,tem9,tem10,tem11,tem12,tem13, & tem14,tem15,tem16,tem17,tem18,tem19, & tem20,tem21,tem22,tem23,tem24,tem25,tem26, & tem1_0,tem2_0,tem3_0) CALL set_acct(misc_acct) ! !----------------------------------------------------------------------- ! ! Update physical time and generate output files ! !----------------------------------------------------------------------- ! curtim = curtim + dtbig ! !----------------------------------------------------------------------- ! ! Print timestep marker ! !----------------------------------------------------------------------- ! IF (myproc == 0) THEN WRITE(6,'(a,i8,a,f10.2,a,f8.3,a)') & ' End of integration time step ', nstep, & ', Model time=',curtim,' s (', curtim/3600.0, ' h)' ENDIF CALL set_acct(output_acct) CALL output(mptr,nx,ny,nz,nstyps,exbcbufsz, & u,v,w,wcont,ptprt,pprt,qv,qc,qr,qi,qs,qh,tke, & udteb,udtwb,vdtnb,vdtsb,pdteb,pdtwb,pdtnb,pdtsb, & ubar,vbar,ptbar,pbar,rhostr,qvbar,kmh,kmv, & x,y,z,zp,hterain, mapfct, j1,j2,j3,j3inv, & soiltyp,stypfrct,vegtyp,lai,roufns,veg, & tsfc,tsoil,wetsfc,wetdp,wetcanp,snowdpth,qvsfc, & ptcumsrc,qcumsrc,w0avg,nca,raincv, & raing,rainc,prcrate,exbcbuf, & radfrc,radsw,rnflx, & usflx,vsflx,ptsflx,qvsflx, & tem1, tem2,tem3,tem4,tem5, tem6,tem7, & tem8,tem9,tem10,tem11) CALL set_acct(misc_acct) ! !----------------------------------------------------------------------- ! ! Check the stability of the integration. If the model solution ! appears to be exceeding the linear stability condition, then ! perform a data dump for post-mortem. ! !----------------------------------------------------------------------- ! DO k=1,nz-1 DO j=1,ny-1 DO i=1,nx-1 tem11(i,j,k) = rhostr(i,j,k)*j3inv(i,j,k) END DO END DO END DO CALL chkstab(mptr,nx,ny,nz,nstyps, & u,v,w,wcont,ptprt,pprt,qv,qc,qr,qi,qs,qh,tke, & ubar,vbar,ptbar,pbar,tem11,qvbar,kmh,kmv, & x,y,z,zp,hterain,mapfct, j1,j2,j3, & soiltyp,stypfrct,vegtyp,lai,roufns,veg, & tsfc,tsoil,wetsfc,wetdp,wetcanp,snowdpth, & raing,rainc,prcrate, & radfrc,radsw,rnflx, & usflx,vsflx,ptsflx,qvsflx, & tem1,tem2,tem3,tem4) ! !----------------------------------------------------------------------- ! ! Calculate and apply the adjustment of domain translation ! speed using either cell-tracking or optimal mean speed algorithm. ! !----------------------------------------------------------------------- ! CALL grdtran(nx,ny,nz,ubar,vbar,u,v,w,ptprt,pprt, & qv,qc,qr,qi,qs,qh,qvbar,rhostr,x,y,zp,j3,j3inv, & tem1,tem2,tem3,tem4,tem5,tem6,tem7, & tem8,tem9,tem10,tem11) ! !----------------------------------------------------------------------- ! ! Time integration loop ends, go back to the beginning of the ! time step loop if stop time is not reached. ! !----------------------------------------------------------------------- ! IF( curtim+eps*dtbig < tstop ) GO TO 900 ! !----------------------------------------------------------------------- ! ! End of entire model time integration. The program stops. ! ! Close opened files. ! !----------------------------------------------------------------------- ! IF (hdmpfmt == 5) THEN CALL mclosescheme (gridid, ierr) CALL mclosedataset (dsindex, ierr) ELSE IF (hdmpfmt == 9) THEN CLOSE (nchdmp) END IF CALL acct_finish CALL acct_report_arps IF (myproc == 0) THEN WRITE(6,'(//1x,a,/1x,a,f13.3,a,/1x,a)') & 'ARPS stopped normally in the main program. ', & 'The ending time was ', curtim,' seconds.', & 'Thanks for using ARPS.' WRITE (6,*) "Maxumum memory allocation (in words):", & max_memory_use END IF IF (mp_opt > 0) THEN CALL mpexit(0) END IF STOP END PROGRAM ARPS