subroutine gloopr x (trie_ls,trio_ls, x ls_node,ls_nodes,max_ls_nodes, x lats_nodes_r,global_lats_r, x lonsperlar, x epse,epso,epsedn,epsodn, x snnp1ev,snnp1od, plnev_r,plnod_r, x pddev_r,pddod_r, ! x snnp1ev,snnp1od,ndexev,ndexod, ! x plnev_r,plnod_r,pddev_r,pddod_r,plnew_r,plnow_r, x phour, & xlon,xlat,coszdg,COSZEN, & SLMSK,SHELEG,SNCOVR,SNOALB,ZORL,TSEA,HPRIME, Clu [+1L]: extract snow-free albedo (SFALB) + SFALB, & ALVSF,ALNSF ,ALVWF ,ALNWF,FACSF ,FACWF,CV,CVT , & CVB ,SWH,HLW,SFCNSW,SFCDLW, & FICE ,TISFC, SFCDSW, ! FOR SEA-ICE - XW Nov04 & TSFLW,FLUXR , phy_f3d,slag,sdec,cdec,NBLCK,KDT, & global_times_r) cc #include "f_hpm.h" ! USE MACHINE , ONLY : kind_phys USE FUNCPHYS , ONLY : fpkap USE PHYSCONS, FV => con_fvirt, rerth => con_rerth ! hmhj use module_radiation_driver, only : radinit, grrad use module_radiation_astronomy,only : astronomy ! !! --- for optional spectral band heating outputs !! use module_radsw_parameters, only : NBDSW !! use module_radlw_parameters, only : NBDLW ! use resol_def use layout1 use gg_def use vert_def use date_def use namelist_def use coordinate_def ! hmhj use tracer_const ! hmhj use d3d_def , only : cldcov ! implicit none include 'mpif.h' ! real (kind=kind_phys), parameter :: QMIN =1.0e-10 real (kind=kind_evod), parameter :: Typical_pgr = 95.0 real (kind=kind_evod), parameter :: cons0 = 0.0, cons2 = 2.0 ! ! --- ... inputs: integer, intent(in) :: ls_node(LS_DIM,3), ls_nodes(LS_DIM,NODES), & & max_ls_nodes(NODES), lats_nodes_r(NODES), & & global_lats_r(LATR), lonsperlar(LATR) integer, intent(in) :: NBLCK real(kind=kind_evod), dimension(LEN_TRIE_LS), intent(in) :: & & epse, epsedn, snnp1ev real(kind=kind_evod), dimension(LEN_TRIO_LS), intent(in) :: & & epso, epsodn, snnp1od real(kind=kind_evod), intent(in) :: plnev_r(LEN_TRIE_LS, LATR2) real(kind=kind_evod), intent(in) :: plnod_r(LEN_TRIO_LS, LATR2) real (kind=kind_phys), dimension(LONR,LATS_NODE_R), intent(in) :: & & xlon, xlat, slmsk, sheleg, zorl, tsea, & & alvsf, alnsf, alvwf, alnwf, facsf, facwf, & & cv, cvt, cvb, FICE, tisfc, sncovr, snoalb real (kind=kind_phys), intent(in) :: & & hprime(NMTVR,LONR,LATS_NODE_R), phour, & & phy_f3d(NGPTC,LEVS,NBLCK,LATS_NODE_R,NUM_P3D) ! ! --- ... input and output: real(kind=kind_evod), intent(inout) :: & & trie_ls(LEN_TRIE_LS,2,11*LEVS+3*LEVH+6), & & trio_ls(LEN_TRIO_LS,2,11*LEVS+3*LEVH+6) integer ipt_ls ! hmhj real(kind=kind_evod) reall ! hmhj real(kind=kind_evod) rlcs2(jcap1) ! hmhj real (kind=kind_phys), intent(inout) :: & & fluxr (NFXR,LONR,LATS_NODE_R) ! --- ... inputs but not used anymore: real(kind=kind_evod), intent(in) :: pddev_r(LEN_TRIE_LS,LATR2), & & pddod_r(LEN_TRIO_LS,LATR2) & ! & plnew_r(LEN_TRIE_LS,LATR2), & ! & plnow_r(LEN_TRIO_LS,LATR2) ! & syn_ls_r(4*LS_DIM,LOTS,LATR2) ! integer, intent(in) :: ndexev(LEN_TRIE_LS), ndexod(LEN_TRIO_LS) integer, intent(in) :: KDT ! --- ... outputs: real(kind=kind_evod), intent(out) :: & & global_times_r(LATG,NODES) real(kind=kind_evod) :: & & for_gr_r_1(LONRX*LOTS,LATS_DIM_R), & & dyn_gr_r_1(lonrx*lotd,lats_dim_r), ! hmhj & for_gr_r_2(LONRX*LOTS,LATS_DIM_R), & dyn_gr_r_2(lonrx*lotd,lats_dim_r) ! hmhj real (kind=kind_phys), intent(out) :: & & swh(NGPTC,LEVS,NBLCK,LATS_NODE_R), & & hlw(NGPTC,LEVS,NBLCK,LATS_NODE_R) real (kind=kind_phys),dimension(LONR,LATS_NODE_R), intent(out) :: & & coszdg, coszen, sfcnsw, sfcdlw, tsflw, & & sfcdsw, SFALB real (kind=kind_phys), intent(out) :: slag, sdec, cdec !! --- ... optional spectral band heating rates !! real (kind=kind_phys), optional, intent(out) :: & !! & htrswb(NGPTC,LEVS,NBDSW,NBLCK,LATS_NODE_R), & !! & htrlwb(NGPTC,LEVS,NBDLW,NBLCK,LATS_NODE_R) ! --- ... locals: ! real(kind=kind_phys) :: prsl(NGPTC,LEVS), prdel(NGPTC,LEVS), & real(kind=kind_phys) :: prsl(NGPTC,LEVS), prslk(NGPTC,LEVS), & & prsi(NGPTC,LEVP1), prsik(NGPTC,LEVP1) real (kind=kind_phys) :: si_loc(LEVR+1) ! real (kind=kind_phys) :: si_loc(LEVR+1), prslk(NGPTC,LEVR) real (kind=kind_phys) :: gu(NGPTC,LEVS), gv1(NGPTC,LEVS), & & gt(NGPTC,LEVR), gd (NGPTC,LEVS), & & gr(NGPTC,LEVR), gr1(NGPTC,LEVR,NTRAC-1), & & gphi(NGPTC), glam(NGPTC), gq(NGPTC), & & sumq(NGPTC,LEVR), xcp(NGPTC,LEVR), &! hmhj & gtv(NGPTC,LEVR), gtvx(NGPTC,LEVR), &! hmhj & gtvy(NGPTC,LEVR) ! hmhj real (kind=kind_phys) :: f_ice(NGPTC,LEVS), f_rain(NGPTC,LEVS), & & r_rime(NGPTC,LEVS) real (kind=kind_phys) :: cldcov_v(NGPTC,LEVS), hprime_v(NGPTC), & & fluxr_v(NGPTC,NFXR), vvel(NGPTC,LEVS) real (kind=kind_phys) :: flgmin_l(ngptc), work1, work2 real (kind=kind_phys) :: rinc(5), dtsw, dtlw, solcon real (kind=kind_phys), save :: facoz integer :: njeff, lon, lan, lat, iblk, lon_dim, lons_lat, istrt integer :: idat(8), jdat(8), DAYS(13), iday, imon, midmon, id integer :: lmax integer, save :: icwp, k1oz, k2oz, midm, midp ! --- number of days in a month data DAYS / 31,28,31,30,31,30,31,31,30,31,30,31,30 / ! --- ... control parameters: ! (some of the them may be moved into model namelist) ! --- ICTM=yyyy#, controls time sensitive external data (e.g. CO2, solcon, aerosols, etc) ! integer, parameter :: ICTM = -2 ! same as ICTM=0, but add seasonal cycle from ! ! climatology. no extrapolation. ! integer, parameter :: ICTM = -1 ! use user provided external data set for the ! ! forecast time. no extrapolation. ! integer, parameter :: ICTM = 0 ! use data at initial cond time, if not ! ! available, use latest, no extrapolation. !! integer, parameter :: ICTM = 1 ! use data at the forecast time, if not ! ! available, use latest and extrapolation. ! integer, parameter :: ICTM =yyyy0 ! use yyyy data for the forecast time, ! ! no further data extrapolation. ! integer, parameter :: ICTM =yyyy1 ! use yyyy data for the fcst. if needed, do ! ! extrapolation to match the fcst time. ! --- ISOL controls solar constant data source !! integer, parameter :: ISOL = 0 ! use prescribed solar constant ! integer, parameter :: ISOL = 1 ! use varying solar const with 11-yr cycle ! --- ICO2 controls co2 data source for radiation ! integer, parameter :: ICO2 = 0 ! prescribed global mean value (old opernl) !! integer, parameter :: ICO2 = 1 ! use obs co2 annual mean value only ! integer, parameter :: ICO2 = 2 ! use obs co2 monthly data with 2-d variation ! --- IALB controls surface albedo for sw radiation !! integer, parameter :: IALB = 0 ! use climatology alb, based on sfc type ! integer, parameter :: IALB = 1 ! use modis derived alb (to be developed) ! --- IEMS controls surface emissivity for lw radiation !! integer, parameter :: IEMS = 0 ! use fixed value of 1.0 ! integer, parameter :: IEMS = 1 ! use varying sfc emiss, based on sfc type ! --- IAER controls aerosols scheme selections ! Old definition ! integer, parameter :: IAER = 1 ! opac climatology, without volc forcing ! integer, parameter :: IAER =11 ! opac climatology, with volcanic forcing ! integer, parameter :: IAER = 2 ! gocart prognostic, without volc forcing ! integer, parameter :: IAER =12 ! gocart prognostic, with volcanic forcing ! New definition in this code ! IAER = 0 --> no aerosol effect at all (volc, sw, lw) ! = 1 --> only tropospheric sw aerosols, no trop-lw and volc ! = 10 --> only tropospheric lw aerosols, no trop-sw and volc ! = 11 --> both trop-sw and trop-lw aerosols, no volc ! = 100 --> only strato-volc aeros, no trop-sw and trop-lw ! = 101 --> only sw aeros (trop + volc), no lw aeros ! = 110 --> only lw aeros (trop + volc), no sw aeros ! = 111 --> both sw and lw aeros (trop + volc) ! ! --- IOVR controls cloud overlapping method in radiation: ! integer, parameter :: IOVR_SW = 0 ! sw: random overlap clouds !! integer, parameter :: IOVR_SW = 1 ! sw: max-random overlap clouds ! integer, parameter :: IOVR_LW = 0 ! lw: random overlap clouds !! integer, parameter :: IOVR_LW = 1 ! lw: max-random overlap clouds ! --- iflip indicates model vertical index direction: ! integer, parameter :: IFLIP = 0 ! virtical profile index from top to bottom integer, parameter :: IFLIP = 1 ! virtical profile index from bottom to top ! ! The following parameters are from gbphys ! real (kind=kind_phys), parameter :: dxmax=-16.118095651, & & dxmin=-9.800790154, dxinv=1.0/(dxmax-dxmin) integer :: kr, kt, kd, kq, ku, kv, ierr, dimg, kx, ky integer :: i, j, k, n integer :: kdtphi,kdtlam,ks ! hmhj logical :: lslag, change, lprnt data lslag / .false. /, lprnt / .false. / logical, save :: first, sas_shal data first / .true. / ! --- timers: real*8 :: rtc, timer1, timer2 ! !===> *** ... begin here ! !! cc !cmr ls_node(1,1) ... ls_node(ls_max_node,1) : values of L !cmr ls_node(1,2) ... ls_node(ls_max_node,2) : values of jbasev !cmr ls_node(1,3) ... ls_node(ls_max_node,3) : values of jbasod cc cc c$$$ integer lots,lotd,lota c$$$cc c$$$ parameter ( lots = 5*levs+1*levh+3 ) c$$$ parameter ( lotd = 6*levs+2*levh+0 ) c$$$ parameter ( lota = 3*levs+1*levh+1 ) cc cc integer kap,kar,kat,kau,kav,kdrlam integer ksd,ksplam,kspphi,ksq,ksr,kst integer ksu,ksv,ksz,node cc ! real(kind=kind_evod) spdlat(levs,lats_dim_r) !Moor real(kind=kind_phys) slk(levs) ! real(kind=kind_evod) spdmax_node (levs) ! real(kind=kind_evod) spdmax_nodes(levs,nodes) cc cc ccxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx cc cc cc................................................................ cc syn(1, 0*levs+0*levh+1, lan) ze cc syn(1, 1*levs+0*levh+1, lan) di cc syn(1, 2*levs+0*levh+1, lan) te cc syn(1, 3*levs+0*levh+1, lan) rq cc syn(1, 3*levs+1*levh+1, lan) q cc syn(1, 3*levs+1*levh+2, lan) dpdlam cc syn(1, 3*levs+1*levh+3, lan) dpdphi cc syn(1, 3*levs+1*levh+4, lan) uln cc syn(1, 4*levs+1*levh+4, lan) vln cc................................................................ cc dyn(1, 0*levs+0*levh+1, lan) d(t)/d(phi) cc dyn(1, 1*levs+0*levh+1, lan) d(rq)/d(phi) cc dyn(1, 1*levs+1*levh+1, lan) d(t)/d(lam) cc dyn(1, 2*levs+1*levh+1, lan) d(rq)/d(lam) cc dyn(1, 2*levs+2*levh+1, lan) d(u)/d(lam) cc dyn(1, 3*levs+2*levh+1, lan) d(v)/d(lam) cc dyn(1, 4*levs+2*levh+1, lan) d(u)/d(phi) cc dyn(1, 5*levs+2*levh+1, lan) d(v)/d(phi) cc................................................................ cc anl(1, 0*levs+0*levh+1, lan) w dudt cc anl(1, 1*levs+0*levh+1, lan) x dvdt cc anl(1, 2*levs+0*levh+1, lan) y dtdt cc anl(1, 3*levs+0*levh+1, lan) rt drdt cc anl(1, 3*levs+1*levh+1, lan) z dqdt cc................................................................ cc cc c$$$ parameter(ksz =0*levs+0*levh+1, c$$$ x ksd =1*levs+0*levh+1, c$$$ x kst =2*levs+0*levh+1, c$$$ x ksr =3*levs+0*levh+1, c$$$ x ksq =3*levs+1*levh+1, c$$$ x ksplam =3*levs+1*levh+2, c$$$ x kspphi =3*levs+1*levh+3, c$$$ x ksu =3*levs+1*levh+4, c$$$ x ksv =4*levs+1*levh+4) cc c$$$ parameter(kdtphi =0*levs+0*levh+1, c$$$ x kdrphi =1*levs+0*levh+1, c$$$ x kdtlam =1*levs+1*levh+1, c$$$ x kdrlam =2*levs+1*levh+1, c$$$ x kdulam =2*levs+2*levh+1, c$$$ x kdvlam =3*levs+2*levh+1, c$$$ x kduphi =4*levs+2*levh+1, c$$$ x kdvphi =5*levs+2*levh+1) cc c$$$ parameter(kau =0*levs+0*levh+1, c$$$ x kav =1*levs+0*levh+1, c$$$ x kat =2*levs+0*levh+1, c$$$ x kar =3*levs+0*levh+1, c$$$ x kap =3*levs+1*levh+1) cc cc c$$$ integer P_gz,P_zem,P_dim,P_tem,P_rm,P_qm c$$$ integer P_ze,P_di,P_te,P_rq,P_q,P_dlam,P_dphi,P_uln,P_vln c$$$ integer P_w,P_x,P_y,P_rt,P_zq c$$$cc c$$$cc old common /comfspec/ c$$$ parameter(P_gz = 0*levs+0*levh+1, ! gze/o(lnte/od,2), c$$$ x P_zem = 0*levs+0*levh+2, ! zeme/o(lnte/od,2,levs), c$$$ x P_dim = 1*levs+0*levh+2, ! dime/o(lnte/od,2,levs), c$$$ x P_tem = 2*levs+0*levh+2, ! teme/o(lnte/od,2,levs), c$$$ x P_rm = 3*levs+0*levh+2, ! rme/o(lnte/od,2,levh), c$$$ x P_qm = 3*levs+1*levh+2, ! qme/o(lnte/od,2), c$$$ x P_ze = 3*levs+1*levh+3, ! zee/o(lnte/od,2,levs), c$$$ x P_di = 4*levs+1*levh+3, ! die/o(lnte/od,2,levs), c$$$ x P_te = 5*levs+1*levh+3, ! tee/o(lnte/od,2,levs), c$$$ x P_rq = 6*levs+1*levh+3, ! rqe/o(lnte/od,2,levh), c$$$ x P_q = 6*levs+2*levh+3, ! qe/o(lnte/od,2), c$$$ x P_dlam= 6*levs+2*levh+4, ! dpdlame/o(lnte/od,2), c$$$ x P_dphi= 6*levs+2*levh+5, ! dpdphie/o(lnte/od,2), c$$$ x P_uln = 6*levs+2*levh+6, ! ulne/o(lnte/od,2,levs), c$$$ x P_vln = 7*levs+2*levh+6, ! vlne/o(lnte/od,2,levs), c$$$ x P_w = 8*levs+2*levh+6, ! we/o(lnte/od,2,levs), c$$$ x P_x = 9*levs+2*levh+6, ! xe/o(lnte/od,2,levs), c$$$ x P_y =10*levs+2*levh+6, ! ye/o(lnte/od,2,levs), c$$$ x P_rt =11*levs+2*levh+6, ! rte/o(lnte/od,2,levh), c$$$ x P_zq =11*levs+3*levh+6) ! zqe/o(lnte/od,2) cc cc ! print *,' in gloopr vertcoord_id =',vertcoord_id ! ksz =0*levs+0*levh+1 ksd =1*levs+0*levh+1 kst =2*levs+0*levh+1 ksr =3*levs+0*levh+1 ksq =3*levs+1*levh+1 ksplam =3*levs+1*levh+2 kspphi =3*levs+1*levh+3 ksu =3*levs+1*levh+4 ksv =4*levs+1*levh+4 kdtphi =0*levs+0*levh+1 ! hmhj kdtlam =1*levs+1*levh+1 ! hmhj cc lslag=.false. cc idat = 0 idat(1) = idate(4) idat(2) = idate(2) idat(3) = idate(3) idat(5) = idate(1) rinc = 0. rinc(2) = fhour call w3movdat(rinc, idat, jdat) ! if (ntoz .le. 0) then ! Climatological Ozone! ! ! if(me .eq. 0) WRITE (6,989) jdat(1),jdat(2),jdat(3),jdat(5) ! 989 FORMAT(' UPDATING OZONE FOR ', I4,I3,I3,I3) ! IDAY = jdat(3) IMON = jdat(2) MIDMON = DAYS(IMON)/2 + 1 CHANGE = FIRST .OR. & ( (IDAY .EQ. MIDMON) .AND. (jdat(5).EQ.0) ) ! IF (CHANGE) THEN IF (IDAY .LT. MIDMON) THEN K1OZ = MOD(IMON+10,12) + 1 MIDM = DAYS(K1OZ)/2 + 1 K2OZ = IMON MIDP = DAYS(K1OZ) + MIDMON ELSE K1OZ = IMON MIDM = MIDMON K2OZ = MOD(IMON,12) + 1 MIDP = DAYS(K2OZ)/2 + 1 + DAYS(K1OZ) ENDIF ENDIF ! IF (IDAY .LT. MIDMON) THEN ID = IDAY + DAYS(K1OZ) ELSE ID = IDAY ENDIF FACOZ = real (ID-MIDM) / real (MIDP-MIDM) endif ! if (first) then sas_shal = sashal .and. (.not. ras) ! if( hybrid.or.gen_coord_hybrid ) then ! hmhj if( gen_coord_hybrid ) then ! hmhj si_loc(levr+1) = si(levp1) ! hmhj do k=1,levr ! hmhj si_loc(k) = si(k) ! hmhj enddo ! hmhj else ! hmhj ! --- get some sigma distribution for radiation-cloud initialization !sela si(k)=(ak5(k)+bk5(k)*Typical_pgr)/Typical_pgr !ak(k) bk(k) go top to botto si_loc(levr+1)= ak5(1)/typical_pgr+bk5(1) do k=1,levr si_loc(levr+1-k)= ak5(levp1-levr+k)/typical_pgr & + bk5(levp1-levr+k) enddo endif else do k = 1, levr si_loc(k) = si(k) enddo si_loc(levr+1) = si(levp1) endif ! end_if_hybrid ! --- determin prognostic/diagnostic cloud scheme icwp = 0 if (NTCW > 0) icwp = 1 first = .false. endif ! end_if_first ! !===> *** ... radiation initialization ! dtsw = 3600.0 * fhswr dtlw = 3600.0 * fhlwr call radinit & ! --- input: & ( si_loc, LEVR, IFLIP, NUM_P3D, ICTM, & & ISOL, ICO2, ICWP, IALB, IEMS, IAER, idat, jdat, me ) ! --- output: ( none ) ! !===> *** ... astronomy for sw radiation calculation. ! call astronomy & ! --- inputs: & ( lonsperlar, global_lats_r, sinlat_r, coslat_r, xlon, & ! & fhswr, jdat, deltim, & & fhswr, jdat, & & LONR, LATS_NODE_R, LATR, IPT_LATS_NODE_R, lsswr, me, & ! --- outputs: & solcon, slag, sdec, cdec, coszen, coszdg & & ) ! print *,' returned from astro' ! !===> *** ... spectrum to grid transformation for radiation calculation. ! ----------------------------------- cc call f_hpmstart(61,"gr delnpe") call delnpe(trie_ls(1,1,P_q ), x trio_ls(1,1,P_dphi), x trie_ls(1,1,P_dlam), x epse,epso,ls_node) call f_hpmstop(61) cc call f_hpmstart(62,"gr delnpo") call delnpo(trio_ls(1,1,P_q ), x trie_ls(1,1,P_dphi), x trio_ls(1,1,P_dlam), x epse,epso,ls_node) call f_hpmstop(62) cc cc call f_hpmstart(63,"gr dezouv dozeuv") ! !$omp parallel do shared(trie_ls,trio_ls) !$omp+shared(epsedn,epsodn,snnp1ev,snnp1od,ls_node) !$omp+private(k) do k=1,levs call dezouv(trie_ls(1,1,P_di +k-1), trio_ls(1,1,P_ze +k-1), x trie_ls(1,1,P_uln+k-1), trio_ls(1,1,P_vln+k-1), x epsedn,epsodn,snnp1ev,snnp1od,ls_node) cc call dozeuv(trio_ls(1,1,P_di +k-1), trie_ls(1,1,P_ze +k-1), x trio_ls(1,1,P_uln+k-1), trie_ls(1,1,P_vln+k-1), x epsedn,epsodn,snnp1ev,snnp1od,ls_node) enddo call f_hpmstop(63) cc !sela print*,'completed call to dztouv' cc !cmr call mpi_barrier (mpi_comm_world,ierr) cc CALL countperf(0,5,0.) CALL synctime() CALL countperf(1,5,0.) !! dimg=0 CALL countperf(0,1,0.) cc call f_hpmstart(67,"gr sumfln") cc !sela print*,'begining call to sumfln' call sumflna_r(trie_ls(1,1,P_ze), x trio_ls(1,1,P_ze), x lat1s_r, x plnev_r,plnod_r, x lots,ls_node,latr2, x lslag,lats_dim_a,lots,for_gr_r_1, x ls_nodes,max_ls_nodes, x lats_nodes_r,global_lats_r, x lats_node_r,ipt_lats_node_r,lon_dims_r,dimg, x lonsperlar,lonrx,latr) cc !sela print*,'completed call to sumfln' call f_hpmstop(67) cc CALL countperf(1,1,0.) cc ! ----------------------------------- if( vertcoord_id.eq.3. ) then ! ----------------------------------- CALL countperf(0,1,0.) ! hmhj ! call f_hpmstart(68,"gr sumder2") ! hmhj ! call sumdera_r(trie_ls(1,1,P_te), ! hmhj x trio_ls(1,1,P_te), ! hmhj x lat1s_r, ! hmhj x pddev_r,pddod_r, ! hmhj x levs,ls_node,latr2, ! hmhj x lslag,lats_dim_r,lotd, ! hmhj x dyn_gr_r_1, ! hmhj x ls_nodes,max_ls_nodes, ! hmhj x lats_nodes_r,global_lats_r, ! hmhj x lats_node_r,ipt_lats_node_r,lon_dims_r,dimg, ! hmhj x lonsperlar,lonrx,latr) ! hmhj ! call f_hpmstop(68) ! hmhj ! CALL countperf(1,1,0.) ! hmhj ! -------------------------------- endif ! vertcoord_id=3 ! -------------------------------- ! !cmr call mpi_barrier (mpi_comm_world,ierr) do lan=1,lats_node_r timer1=rtc() cc lat = global_lats_r(ipt_lats_node_r-1+lan) cc lon_dim = lon_dims_r(lan) cc lons_lat = lonsperlar(lat) ! ------------------------------------------------------- if( gen_coord_hybrid .and. vertcoord_id.eq.3. ) then ! ------------------------------------------------------- ! lmax = min(jcap,lons_lat/2) ! hmhj ! ipt_ls=min(lat,latr-lat+1) ! hmhj do i=1,lmax+1 ! hmhj if ( ipt_ls .ge. lat1s_r(i-1) ) then ! hmhj reall=i-1 ! hmhj rlcs2(i)=reall*rcs2_r(ipt_ls)/rerth ! hmhj else ! hmhj rlcs2(i)=cons0 !constant ! hmhj endif ! hmhj enddo ! hmhj ! !$omp parallel do private(k,i) do k=1,levs ! hmhj do i=1,lmax+1 ! hmhj ! ! d(t)/d(lam) ! hmhj dyn_gr_r_1(2*i-1+(kdtlam-2+k)*lon_dim,lan)= ! hmhj x -for_gr_r_1(2*i +(kst -2+k)*lon_dim,lan)*rlcs2(i) ! hmhj dyn_gr_r_1(2*i +(kdtlam-2+k)*lon_dim,lan)= ! hmhj x for_gr_r_1(2*i-1+(kst -2+k)*lon_dim,lan)*rlcs2(i) ! hmhj enddo ! hmhj enddo ! hmhj ! -------------------- endif ! gc and vertcoord_id=3 ! --------------------- ! cc CALL countperf(0,6,0.) !sela print*,' beginning call four2grid',lan CALL FOUR2GRID_thread(for_gr_r_1(1,lan),for_gr_r_2(1,lan), & & lon_dim,lons_lat,lonrx,5*levs+levh+3,lan,me) ! ------------------------------------------------------- if( gen_coord_hybrid.and.vertcoord_id.eq.3. ) then ! hmhj ! ------------------------------------------------------- CALL FOUR2GRID_thread(dyn_gr_r_1(1,lan),dyn_gr_r_2(1,lan), ! hmhj & lon_dim,lons_lat,lonrx,levs,lan,me) ! hmhj CALL FOUR2GRID_thread(dyn_gr_r_1((kdtlam-1)*lon_dim+1,lan), ! hmhj & dyn_gr_r_2((kdtlam-1)*lon_dim+1,lan), ! hmhj & lon_dim,lons_lat,lonrx,levs,lan,me) ! hmhj ! ------------------------- endif ! gc and vertcoord_id=3 ! ------------------------- !sela print*,' completed call four2grid lan=',lan CALL countperf(1,6,0.) !! if( .not. gen_coord_hybrid ) then ! hmhj do k = 1, LEVS kr = (KSR + k - 2) * lon_dim kt = (KST + k - 2) * lon_dim do j = 1, lons_lat if (for_gr_r_2(j+kr,lan) <= 0.0) then for_gr_r_2(j+kr,lan) = QMIN endif for_gr_r_2(j+kt,lan) = for_gr_r_2(j+kt,lan) & & / (1.0 + FV*for_gr_r_2(j+kr,lan)) enddo enddo kq = (KSQ - 1)*lon_dim do j = 1, lons_lat for_gr_r_2(j+kq,lan) = exp( for_gr_r_2(j+kq,lan) ) enddo endif ! hmhj c timer2=rtc() global_times_r(lat,me+1)=timer2-timer1 c$$$ print*,'timeloopr',me,timer1,timer2,global_times_r(lat,me+1) !! enddo !lan ! call f_hpmstart(69,"gr lat_loop2") ! !===> *** ... starting latitude loop ! do lan=1,lats_node_r cc lat = global_lats_r(ipt_lats_node_r-1+lan) cc lons_lat = lonsperlar(lat) !! !$omp parallel do schedule(dynamic,1) private(lon,j,k,lon_dim) !$omp+private(istrt,njeff,iblk,ku,kv,kd,kq,kt,kr,kx,ky,ks,n) !$omp+private(vvel,gu,gv1,gd,gt,gr,gr1,gq,gphi,glam) !$omp+private(gtv,gtvx,gtvy,sumq,xcp) !$omp+private(cldcov_v,hprime_v,fluxr_v,f_ice,f_rain,r_rime) !$omp+private(prslk,prsl,prsik,prsi) DO lon=1,lons_lat,NGPTC !! lon_dim = lon_dims_r(lan) NJEFF = MIN(NGPTC,lons_lat-lon+1) ISTRT = lon if (NGPTC.ne.1) then IBLK = lon/NGPTC+1 else IBLK = lon endif do k = 1, LEVS ku = lon - 1 + (KSU + k - 2)*lon_dim kv = lon - 1 + (KSV + k - 2)*lon_dim kd = lon - 1 + (KSD + k - 2)*lon_dim do j = 1, njeff gu(j,k) = for_gr_r_2(j+ku,lan) gv1(j,k) = for_gr_r_2(j+kv,lan) gd(j,k) = for_gr_r_2(j+kd,lan) enddo enddo if( gen_coord_hybrid ) then ! hmhj do k=1,levr ! hmhj kt = lon - 1 + (KST + k - 2)*lon_dim kr = lon - 1 + (KSR + k - 2)*lon_dim do j=1,njeff ! hmhj gtv(j,k) = for_gr_r_2(j+kt,lan) gr(j,k) = max(qmin, for_gr_r_2(j+kr,lan)) enddo ! hmhj enddo ! hmhj ! -------------------------------------- if( vertcoord_id.eq.3. ) then ! -------------------------------------- do k=1,levr ! hmhj kx = lon - 1 + (kdtlam + k - 2)*lon_dim ky = lon - 1 + (kdtphi + k - 2)*lon_dim do j=1,njeff ! hmhj gtvx(j,k) = dyn_gr_r_2(j+kx,lan) gtvy(j,k) = dyn_gr_r_2(j+ky,lan) enddo ! hmhj enddo ! hmhj ! ----------------------------- endif ! ----------------------------- if( thermodyn_id.eq.3 ) then ! get dry temperature from enthalpy ! hmhj sumq=0.0 ! hmhj xcp=0.0 ! hmhj do i=1,ntrac ! hmhj if( cpi(i).ne.0.0 ) then ! hmhj ks=ksr+(i-1)*levs ! hmhj do k=1,levr ! hmhj kr = lon - 1 + (ks + k - 2)*lon_dim ! hmhj do j=1,njeff ! hmhj sumq(j,k)=sumq(j,k)+for_gr_r_2(j+kr,lan) ! hmhj xcp(j,k)=xcp(j,k)+cpi(i)*for_gr_r_2(j+kr,lan) ! hmhj enddo ! hmhj enddo ! hmhj endif ! hmhj enddo ! hmhj do k=1,levr ! hmhj do j=1,njeff ! hmhj xcp(j,k)=(1.-sumq(j,k))*cpi(0)+xcp(j,k) ! hmhj gt(j,k)=gtv(j,k)/xcp(j,k) ! hmhj enddo ! hmhj enddo ! hmhj else if( thermodyn_id.le.1 ) then ! hmhj ! get dry temperature from virtual temperature ! hmhj do k=1,levr ! hmhj do j=1,njeff ! hmhj gt(j,k) = gtv(j,k) / (1.+fv*gr(j,k)) ! hmhj enddo ! hmhj enddo ! hmhj else ! get dry temperature from dry temperature ! hmhj do k=1,levr ! hmhj do j=1,njeff ! hmhj gt(j,k) = gtv(j,k) ! hmhj enddo ! hmhj enddo endif else ! hmhj ! do k = 1, LEVR kt = lon - 1 + (KST + k - 2)*lon_dim kr = lon - 1 + (KSR + k - 2)*lon_dim do j = 1, njeff gt(j,k) = for_gr_r_2(j+kt,lan) gr(j,k) = for_gr_r_2(j+kr,lan) enddo enddo endif ! ! Remaining tracers ! do n = 1, NTRAC-1 do k = 1, LEVR kr = lon - 1 + (KSR + n*LEVS + k - 2)*lon_dim do j = 1, njeff gr1(j,k,n) = for_gr_r_2(j+kr,lan) enddo enddo enddo kq = lon - 1 + (KSQ - 1)*lon_dim kt = lon - 1 + (KSPPHI - 1)*lon_dim kr = lon - 1 + (KSPLAM - 1)*lon_dim do j = 1, njeff gq (j) = for_gr_r_2(j+kq,lan) gphi(j) = for_gr_r_2(j+kt,lan) glam(j) = for_gr_r_2(j+kr,lan) enddo !! ! --- vertical structure variables: del,si,sl,prslk,prdel ! if( gen_coord_hybrid ) then ! hmhj !Moor call hyb2press_gc(njeff,ngptc,gq,gtv,prsi,prsl,prdel) ! hmhj call hyb2press_gc(njeff,ngptc,gq,gtv,prsi,prsl,prsik,prslk) ! hmhj call omegtes_gc(njeff,ngptc,rcs2_r(min(lat,latr-lat+1)), ! hmhj & gq,gphi,glam,gtv,gtvx,gtvy,gd,gu,gv1,vvel) ! hmhj elseif (hybrid) then !Moor call hyb2press(njeff,ngptc,gq, prsi, prsl,prdel) call hyb2press(njeff,ngptc,gq, prsi, prsl,prsik, prslk) call omegtes(njeff,ngptc,rcs2_r(min(lat,latr-lat+1)), & gq,gphi,glam,gd,gu,gv1,vvel) else !Moor call sig2press(njeff,ngptc,gq,sl,del,si, prsi, prsl,prdel) call sig2press(njeff,ngptc,gq,sl,si,slk,sik, & prsi,prsl,prsik,prslk) CALL countperf(0,12,0.) call omegast3(njeff,ngptc,levs, & gphi,glam,gu,gv1,gd,del, & rcs2_r(min(lat,latr-lat+1)),vvel,gq,sl) endif !..... if (levr .lt. levs) then do j=1,njeff prsi(j,levr+1) = prsi(j,levp1) prsl(j,levr) = (prsi(j,levp1)+prsi(j,levr)) * 0.5 prsik(j,levr+1) = prslk(j,levp1) prslk(j,levr) = fpkap(prsl(j,levr)*1000.0) enddo endif if (ldiag3d) then do k=1,levr do j=1,njeff !Moor prslk(j,k) = fpkap(prsl(j,k)*1000.0) cldcov_v(j,k) = cldcov(k,istrt+j-1,lan) enddo enddo endif ! do j=1,njeff hprime_v(j) = hprime(1,istrt+j-1,lan) enddo ! do k=1,nfxr do j=1,njeff fluxr_v(j,k) = fluxr(k,istrt+j-1,lan) enddo enddo if (NUM_P3D == 3) then do k = 1, LEVR do j = 1, njeff f_ice (j,k) = phy_f3d(j,k,iblk,lan,1) f_rain(j,k) = phy_f3d(j,k,iblk,lan,2) r_rime(j,k) = phy_f3d(j,k,iblk,lan,3) enddo enddo endif work1 = (log(coslat_r(lat) / (lons_lat*latg)) - dxmin) * dxinv work1 = max(0.0, min(1.0,work1)) work2 = flgmin(1)*work1 + flgmin(2)*(1.0-work1) do j=1,njeff flgmin_l(j) = work2 enddo ! *** ... calling radiation driver ! ! lprnt = me .eq. 0 .and. kdt .ge. 120 ! if (lprnt) then ! if (kdt .gt. 85) then ! print *,' calling grrad for me=',me,' lan=',lan,' lat=',lat ! &,' num_p3d=',num_p3d,' snoalb=',snoalb(lon,lan),' lon=',lon ! &,' tsea=',tsea(lon,lan),' sncovr=',sncovr(lon,lan), ! &' sheleg=',sheleg(lon,lan) ! call grrad & ! --- inputs: & ( prsi,prsl,prslk,gt,gr,gr1,vvel,slmsk(lon,lan), & & xlon(lon,lan),xlat(lon,lan),tsea(lon,lan), & & sheleg(lon,lan),sncovr(lon,lan),snoalb(lon,lan), & & zorl(lon,lan),hprime_v, & & alvsf(lon,lan),alnsf(lon,lan),alvwf(lon,lan), & & alnwf(lon,lan),facsf(lon,lan),facwf(lon,lan), & ! fice FOR SEA-ICE XW Nov04 & fice(lon,lan),tisfc(lon,lan), & & solcon,coszen(lon,lan),coszdg(lon,lan),k1oz,k2oz,facoz, & & cv(lon,lan),cvt(lon,lan),cvb(lon,lan), & & IOVR_SW,IOVR_LW,f_ice,f_rain,r_rime,flgmin_l, & & NUM_P3D,NTCW-1,NCLD,NTOZ-1,NTRAC-1,NFXR, & & dtlw,dtsw, lsswr,lslwr,lssav,ldiag3d,sas_shal,norad_precip,& & crick_proof, ccnorm, & ! & dtlw,dtsw, lsswr,lslwr,lssav,ldiag3d, & & NGPTC,njeff,LEVR,IFLIP, me, lprnt, & ! --- outputs: & swh(1,1,iblk,lan),sfcnsw(lon,lan),sfcdsw(lon,lan), & ! sfcdsw FOR SEA-ICE XW Nov04 & sfalb(lon,lan), & ! lu [+1L]: add sfalb & hlw(1,1,iblk,lan),sfcdlw(lon,lan),tsflw(lon,lan), & ! --- input/output: & fluxr_v,cldcov_v & !! --- optional outputs: !! &, HTRSWB=htrswb(1,1,1,iblk,lan), & !! &, HTRLWB=htrlwb(1,1,1,iblk,lan) & & ) ! ! if (lprnt) print *,' returned from grrad for me=',me,' lan=', ! &lan,' lat=',lat,' kdt=',kdt ! ! if (ldiag3d) then do k=1,levr do j=1,njeff cldcov(k,istrt+j-1,lan) = cldcov_v(j,k) enddo enddo endif do k=1,nfxr do j=1,njeff fluxr(k,istrt+j-1,lan) = fluxr_v(j,k) enddo enddo if (levr .lt. levs) then do k=levr+1,levs do j=1,njeff hlw(j,k,iblk,lan) = hlw(j,levr,iblk,lan) swh(j,k,iblk,lan) = swh(j,levr,iblk,lan) enddo enddo endif c$$$ write(2900+lat,*) ' ilon = ',istrt c$$$ write(2900+lat,'("swh",T16,"hlw")') c$$$ do k=1,levs c$$$ write(2900+lat, c$$$ . '(e10.3,T16,e10.3,T31,e10.3)') c$$$ . swh(1,k,iblk,lan),hlw(1,k,iblk,lan) c$$$ enddo !! CALL countperf(1,12,0.) ENDDO ! enddo cc call f_hpmstop(69) !! CALL countperf(0,5,0.) CALL synctime() CALL countperf(1,5,0.) !sela print*,'completed gloopr_v kdt=',kdt !! return end subroutine gloopr