diff --git a/src/2d/bdy_2d.F90 b/src/2d/bdy_2d.F90 index b60eca30df6285d1db7bd1ccc61af16b3a27000a..4ca6213af7f10fb546015fd03aaa80d7a2705e04 100644 --- a/src/2d/bdy_2d.F90 +++ b/src/2d/bdy_2d.F90 @@ -390,16 +390,18 @@ select case (bdy_2d_type(l)) case (ZERO_GRADIENT,CLAMPED_VEL,FLATHER_VEL) do j = wfj(n),wlj(n) - z(i,j) = z(i+1,j) + if ( az(i+1,j) .ne. 0 ) z(i,j) = z(i+1,j) end do case (SOMMERFELD) do j = wfj(n),wlj(n) + if ( az(i+1,j) .ne. 0 ) then cfl = sqrt(g*_HALF_*(D(i,j)+D(i+1,j)))*dtm/DXU z(i,j) = ( & (_ONE_ - _TWO_*cfl*(_ONE_-theta))*z (i ,j) & +(_ONE_ + _TWO_*cfl*(_ONE_-theta))*zo(i+1,j) & -(_ONE_ - _TWO_*cfl*theta )*z (i+1,j) & )/(_ONE_ + _TWO_*cfl*theta ) + end if end do case (CLAMPED_ELEV,CLAMPED) do j = wfj(n),wlj(n) @@ -409,13 +411,16 @@ end do case (FLATHER_ELEV) do j = wfj(n),wlj(n) + if ( az(i+1,j) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i,j)+D(i+1,j)) ! Note (KK): note approximation of sse at vel-time stage - a = ramp*bdy_data(kl) & - - _TWO_/sqrt(g*depth)*(U(i,j)-ramp*bdy_data_u(kl)*depth) - z(i,j) = max(a,-H(i,j)+min_depth) + a = _TWO_/sqrt(g*depth)*(U(i,j)-ramp*bdy_data_u(kl)*depth) + else + a = _ZERO_ + end if + z(i,j) = max(ramp*bdy_data(kl)-a,-H(i,j)+min_depth) k = k+1 kl = kl + 1 end do @@ -429,16 +434,18 @@ select case (bdy_2d_type(l)) case (ZERO_GRADIENT,CLAMPED_VEL,FLATHER_VEL) do i = nfi(n),nli(n) - z(i,j) = z(i,j-1) + if ( az(i,j-1) .ne. 0 ) z(i,j) = z(i,j-1) end do case (SOMMERFELD) do i = nfi(n),nli(n) + if ( az(i,j-1) .ne. 0 ) then cfl = sqrt(g*_HALF_*(D(i,j-1)+D(i,j)))*dtm/DYVJM1 z(i,j) = ( & (_ONE_ - _TWO_*cfl*(_ONE_-theta))*z (i,j ) & +(_ONE_ + _TWO_*cfl*(_ONE_-theta))*zo(i,j-1) & -(_ONE_ - _TWO_*cfl*theta )*z (i,j-1) & )/(_ONE_ + _TWO_*cfl*theta ) + end if end do case (CLAMPED_ELEV,CLAMPED) do i = nfi(n),nli(n) @@ -448,13 +455,16 @@ end do case (FLATHER_ELEV) do i = nfi(n),nli(n) + if ( az(i,j-1) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i,j-1)+D(i,j)) ! Note (KK): note approximation of sse at vel-time stage - a = ramp*bdy_data(kl) & - + _TWO_/sqrt(g*depth)*(V(i,j-1)-ramp*bdy_data_v(kl)*depth) - z(i,j) = max(a,-H(i,j)+min_depth) + a = _TWO_/sqrt(g*depth)*(V(i,j-1)-ramp*bdy_data_v(kl)*depth) + else + a = _ZERO_ + end if + z(i,j) = max(ramp*bdy_data(kl)+a,-H(i,j)+min_depth) k = k+1 kl = kl + 1 end do @@ -468,16 +478,18 @@ select case (bdy_2d_type(l)) case (ZERO_GRADIENT,CLAMPED_VEL,FLATHER_VEL) do j = efj(n),elj(n) - z(i,j) = z(i-1,j) + if ( az(i-1,j) .ne. 0 ) z(i,j) = z(i-1,j) end do case (SOMMERFELD) do j = efj(n),elj(n) + if ( az(i-1,j) .ne. 0 ) then cfl = sqrt(g*_HALF_*(D(i-1,j)+D(i,j)))*dtm/DXUIM1 z(i,j) = ( & (_ONE_ - _TWO_*cfl*(_ONE_-theta))*z (i ,j) & +(_ONE_ + _TWO_*cfl*(_ONE_-theta))*zo(i-1,j) & -(_ONE_ - _TWO_*cfl*theta )*z (i-1,j) & )/(_ONE_ + _TWO_*cfl*theta ) + end if end do case (CLAMPED_ELEV,CLAMPED) do j = efj(n),elj(n) @@ -487,13 +499,16 @@ end do case (FLATHER_ELEV) do j = efj(n),elj(n) + if ( az(i-1,j) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i-1,j)+D(i,j)) ! Note (KK): note approximation of sse at vel-time stage - a = ramp*bdy_data(kl) & - + _TWO_/sqrt(g*depth)*(U(i-1,j)-ramp*bdy_data_u(kl)*depth) - z(i,j) = max(a,-H(i,j)+min_depth) + a = _TWO_/sqrt(g*depth)*(U(i-1,j)-ramp*bdy_data_u(kl)*depth) + else + a = _ZERO_ + end if + z(i,j) = max(ramp*bdy_data(kl)+a,-H(i,j)+min_depth) k = k+1 kl = kl + 1 end do @@ -507,16 +522,18 @@ select case (bdy_2d_type(l)) case (ZERO_GRADIENT,CLAMPED_VEL,FLATHER_VEL) do i = sfi(n),sli(n) - z(i,j) = z(i,j+1) + if ( az(i,j+1) .ne. 0 ) z(i,j) = z(i,j+1) end do case (SOMMERFELD) do i = sfi(n),sli(n) + if ( az(i,j+1) .ne. 0 ) then cfl = sqrt(g*_HALF_*(D(i,j)+D(i,j+1)))*dtm/DYV z(i,j) = ( & (_ONE_ - _TWO_*cfl*(_ONE_-theta))*z (i,j ) & +(_ONE_ + _TWO_*cfl*(_ONE_-theta))*zo(i,j+1) & -(_ONE_ - _TWO_*cfl*theta )*z (i,j+1) & )/(_ONE_ + _TWO_*cfl*theta ) + end if end do case (CLAMPED_ELEV,CLAMPED) do i = sfi(n),sli(n) @@ -526,13 +543,16 @@ end do case (FLATHER_ELEV) do i = sfi(n),sli(n) + if ( az(i,j+1) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i,j)+D(i,j+1)) ! Note (KK): note approximation of sse at vel-time stage - a = ramp*bdy_data(kl) & - - _TWO_/sqrt(g*depth)*(V(i,j)-ramp*bdy_data_v(kl)*depth) - z(i,j) = max(a,-H(i,j)+min_depth) + a = _TWO_/sqrt(g*depth)*(V(i,j)-ramp*bdy_data_v(kl)*depth) + else + a = _ZERO_ + end if + z(i,j) = max(ramp*bdy_data(kl)-a,-H(i,j)+min_depth) k = k+1 kl = kl + 1 end do @@ -627,6 +647,7 @@ select case (bdy_2d_type(l)) case (FLATHER_VEL) do j = wfj(n),wlj(n) + if ( az(i+1,j) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i,j)+D(i+1,j)) @@ -636,11 +657,13 @@ if ( waveforcing_method .ne. NO_WAVES ) then UEuler(i,j) = U(i,j) - UStokes(i,j) end if + end if k = k+1 kl = kl + 1 end do case (CLAMPED_VEL,CLAMPED) do j = wfj(n),wlj(n) + if ( az(i+1,j) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i,j)+D(i+1,j)) @@ -648,6 +671,7 @@ if ( waveforcing_method .ne. NO_WAVES ) then UEuler(i,j) = U(i,j) - UStokes(i,j) end if + end if k = k+1 kl = kl + 1 end do @@ -662,6 +686,7 @@ select case (bdy_2d_type(l)) case (FLATHER_VEL) do j = efj(n),elj(n) + if ( az(i-1,j) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i-1,j)+D(i,j)) @@ -671,11 +696,13 @@ if ( waveforcing_method .ne. NO_WAVES ) then UEuler(i-1,j) = U(i-1,j) - UStokes(i-1,j) end if + end if k = k+1 kl = kl + 1 end do case (CLAMPED_VEL,CLAMPED) do j = efj(n),elj(n) + if ( az(i-1,j) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i-1,j)+D(i,j)) @@ -683,6 +710,7 @@ if ( waveforcing_method .ne. NO_WAVES ) then UEuler(i-1,j) = U(i-1,j) - UStokes(i-1,j) end if + end if k = k+1 kl = kl + 1 end do @@ -700,6 +728,7 @@ select case (bdy_2d_type(l)) case (FLATHER_VEL) do i = nfi(n),nli(n) + if ( az(i,j-1) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i,j-1)+D(i,j)) @@ -709,11 +738,13 @@ if ( waveforcing_method .ne. NO_WAVES ) then VEuler(i,j-1) = V(i,j-1) - VStokes(i,j-1) end if + end if k = k+1 kl = kl + 1 end do case (CLAMPED_VEL,CLAMPED) do i = nfi(n),nli(n) + if ( az(i,j-1) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i,j-1)+D(i,j)) @@ -721,6 +752,7 @@ if ( waveforcing_method .ne. NO_WAVES ) then VEuler(i,j-1) = V(i,j-1) - VStokes(i,j-1) end if + end if k = k+1 kl = kl + 1 end do @@ -735,6 +767,7 @@ select case (bdy_2d_type(l)) case (FLATHER_VEL) do i = sfi(n),sli(n) + if ( az(i,j+1) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i,j)+D(i,j+1)) @@ -744,11 +777,13 @@ if ( waveforcing_method .ne. NO_WAVES ) then VEuler(i,j) = V(i,j) - VStokes(i,j) end if + end if k = k+1 kl = kl + 1 end do case (CLAMPED_VEL,CLAMPED) do i = sfi(n),sli(n) + if ( az(i,j+1) .ne. 0 ) then ! Note (KK): approximate interface depths at vel-time stage ! by spatial mean at last sse-time stage depth = _HALF_*(D(i,j)+D(i,j+1)) @@ -756,6 +791,7 @@ if ( waveforcing_method .ne. NO_WAVES ) then VEuler(i,j) = V(i,j) - VStokes(i,j) end if + end if k = k+1 kl = kl + 1 end do diff --git a/src/2d/m2d.F90 b/src/2d/m2d.F90 index aae261db7012e05014e1b9697df232d677ef82cc..8931fee23bc61ef3bbd6103e208f651167569b73 100644 --- a/src/2d/m2d.F90 +++ b/src/2d/m2d.F90 @@ -167,7 +167,6 @@ integer :: i,j integer :: elev_method=1 REALTYPE :: elev_const=_ZERO_ - integer,parameter :: rk = kind(_ONE_) character(LEN = PATH_MAX) :: elev_file='elev.nc' namelist /m2d/ & elev_method,elev_const,elev_file, & @@ -434,8 +433,8 @@ ! This is only needed for proper flexible output where (az .eq. 0) - z = -9999.0d0 - zo = -9999.0d0 + z = -9999._rk + zo = -9999._rk end where return diff --git a/src/2d/variables_2d.F90 b/src/2d/variables_2d.F90 index a07fd328ba70c3563dd9706972ecd3bdc07aec37..c0596e3770eb6da85139ce7b2fe19c6413a4236e 100644 --- a/src/2d/variables_2d.F90 +++ b/src/2d/variables_2d.F90 @@ -23,6 +23,7 @@ IMPLICIT NONE ! ! !PUBLIC DATA MEMBERS: + integer, parameter :: rk = kind(_ONE_) REALTYPE :: dtm REALTYPE,dimension(:,:),pointer :: zo,z logical :: deformC=.false. @@ -121,7 +122,6 @@ ! ! !LOCAL VARIABLES: integer :: rc - integer,parameter :: rk = kind(_ONE_) !EOP !------------------------------------------------------------------------- !BOC @@ -174,12 +174,12 @@ end if - z = -9999*_ONE_; zo =_ZERO_ + z = -9999._rk ; zo =_ZERO_ U = _ZERO_; DU = _ZERO_; Uint = _ZERO_; UEx = _ZERO_ V = _ZERO_; DV = _ZERO_; Vint = _ZERO_; VEx = _ZERO_ - velx = -9999.0 ; vely = -9999.0 + velx = -9999._rk ; vely = -9999._rk where (az .gt. 0) velx = _ZERO_ ; vely = _ZERO_ end where @@ -336,7 +336,6 @@ ! ! !LOCAL VARIABLES: logical :: used - integer,parameter :: rk = kind(_ONE_) !EOP !----------------------------------------------------------------------- !BOC diff --git a/src/3d/m3d.F90 b/src/3d/m3d.F90 index 12dd39514327c1f8522b8a71c0b9fff15ff59658..3c62de9fc2f910867c91f8eb8164a3b7a6d9013e 100644 --- a/src/3d/m3d.F90 +++ b/src/3d/m3d.F90 @@ -442,8 +442,8 @@ num(i,j,:) = 1.e-15 nuh(i,j,:) = 1.e-15 #ifndef NO_BAROCLINIC - S(i,j,:) = -9999.0 - T(i,j,:) = -9999.0 + S(i,j,:) = -9999._rk + T(i,j,:) = -9999._rk #endif end if end do diff --git a/src/3d/salinity.F90 b/src/3d/salinity.F90 index f197db256bb3a12bd873693fd3f90733713ce28b..df19d4c2e0143d78d2677e6679b2a78f49993fed 100644 --- a/src/3d/salinity.F90 +++ b/src/3d/salinity.F90 @@ -19,7 +19,7 @@ use domain, only: imin,jmin,imax,jmax,kmax,H,az,dry_z !KB use get_field, only: get_3d_field use variables_2d, only: fwf_int - use variables_3d, only: S,hn,kmin + use variables_3d, only: rk,S,hn,kmin use meteo, only: metforcing,met_method,nudge_sss,sss use meteo, only: METEO_CONST,METEO_FROMFILE,METEO_FROMEXT use halo_zones, only: update_3d_halo,wait_halo,D_TAG,H_TAG @@ -305,8 +305,8 @@ stop 'init_salinity' end select - S(:,:,0) = -9999*_ONE_ - forall(i=imin:imax,j=jmin:jmax, az(i,j).eq.0) S(i,j,:) = -9999*_ONE_ + S(:,:,0) = -9999._rk + forall(i=imin:imax,j=jmin:jmax, az(i,j).eq.0) S(i,j,:) = -9999._rk call update_3d_halo(S,S,az,imin,jmin,imax,jmax,kmax,D_TAG) call wait_halo(D_TAG) diff --git a/src/3d/temperature.F90 b/src/3d/temperature.F90 index 5ff525ea48cc2b155917d35d181f0e659c3e5a07..5d665f1772be7f7d9a6146c711750ea86ed413c1 100644 --- a/src/3d/temperature.F90 +++ b/src/3d/temperature.F90 @@ -20,7 +20,7 @@ use domain, only: imin,jmin,imax,kmax,jmax,H,az,dry_z use domain, only: ill,ihl,jll,jhl use domain, only: ilg,ihg,jlg,jhg - use variables_3d, only: T,rad,hn,kmin,A,g1,g2,heatflux_net + use variables_3d, only: rk,T,rad,hn,kmin,A,g1,g2,heatflux_net use meteo, only: metforcing,met_method,nudge_sst,sst use meteo, only: METEO_CONST,METEO_FROMFILE,METEO_FROMEXT !KB use get_field, only: get_3d_field @@ -375,8 +375,8 @@ end interface stop 'init_temperature' end select - T(:,:,0) = -9999*_ONE_ - forall(i=imin:imax,j=jmin:jmax, az(i,j).eq.0) T(i,j,:) = -9999*_ONE_ + T(:,:,0) = -9999._rk + forall(i=imin:imax,j=jmin:jmax, az(i,j).eq.0) T(i,j,:) = -9999._rk call update_3d_halo(T,T,az,imin,jmin,imax,jmax,kmax,D_TAG) call wait_halo(D_TAG) diff --git a/src/3d/variables_3d.F90 b/src/3d/variables_3d.F90 index 05e626c7fd98b4f53d3b9493ce6bbc34819e839c..cbe9613016d8755f677bf77085a9903ebfe4524e 100644 --- a/src/3d/variables_3d.F90 +++ b/src/3d/variables_3d.F90 @@ -117,6 +117,7 @@ IMPLICIT NONE ! ! !PUBLIC DATA MEMBERS: + integer, parameter :: rk = kind(_ONE_) REALTYPE :: dt,cnpar=0.9 REALTYPE :: avmback=_ZERO_,avhback=_ZERO_ logical :: do_numerical_analyses_3d=.false. @@ -213,7 +214,6 @@ ! ! !LOCAL VARIABLES: integer :: i,j, rc - integer,parameter :: rk = kind(_ONE_) !EOP !------------------------------------------------------------------------- !BOC @@ -315,8 +315,8 @@ #endif ! must be nonzero for gotm_fabm in case of calc_temp=F - g1 = -9999*_ONE_ - g2 = -9999*_ONE_ + g1 = -9999._rk + g2 = -9999._rk #ifdef STRUCTURE_FRICTION sf = _ZERO_ @@ -513,7 +513,7 @@ ! Original author(s): Karsten Bolding & Jorn Bruggeman ! ! !LOCAL VARIABLES: - integer,parameter :: rk = kind(_ONE_) +! !EOP !----------------------------------------------------------------------- !BOC diff --git a/src/3d/vertical_coordinates.F90 b/src/3d/vertical_coordinates.F90 index a1b56245bad30976164fd03339a887fee88df21e..fddb5ed6d1dfceee1eb4028d0d92a87eee1b456e 100644 --- a/src/3d/vertical_coordinates.F90 +++ b/src/3d/vertical_coordinates.F90 @@ -183,9 +183,9 @@ stop ! calculate the z-coordinate of the cell centers ! references to mean sea level zc(:,:,0)=-H(:,:) - zc(:,:,1)=-H(:,:) + 0.5*hn(:,:,1) + zc(:,:,1)=-H(:,:) + _HALF_*hn(:,:,1) do k=2,kmax - zc(:,:,k)=zc(:,:,k-1)+0.5*(hn(:,:,k-1)+hn(:,:,k)) + zc(:,:,k)=zc(:,:,k-1)+_HALF_*(hn(:,:,k-1)+hn(:,:,k)) end do #ifdef SLICE_MODEL diff --git a/src/coupling/getm_esmf.F90 b/src/coupling/getm_esmf.F90 index c0208d439d521a95e8f680b6306136b33058348d..d55e021986a285464f9e6fa5a8b8191fb1878e61 100644 --- a/src/coupling/getm_esmf.F90 +++ b/src/coupling/getm_esmf.F90 @@ -1,5 +1,5 @@ #include "cppdefs.h" -!#define ALLEXPORT +#define _GETM_NUOPC_ !----------------------------------------------------------------------- !BOP ! @@ -97,6 +97,16 @@ "wind_x_velocity_at_10m" character(len=*),parameter :: name_windV = & "wind_y_velocity_at_10m" + character(len=*),parameter :: name_airT2 = & + "air_temperature_at_2m" + character(len=*),parameter :: name_humr = & + "relative_humidity" + character(len=*),parameter :: name_hums = & + "specific_humidity" + character(len=*),parameter :: name_dev2 = & + "dew_point_temperature" + character(len=*),parameter :: name_tcc = & + "total_cloud_cover" ! ! !REVISION HISTORY: ! Original author(s): Knut Klingbeil @@ -461,26 +471,28 @@ InitializePhaseMap(1) = "IPDv00p1=1" InitializePhaseMap(2) = "IPDv00p2=2" -! Note (KK): NUOPC attributes are purpose="Instance" -#if 1 +#ifdef _GETM_NUOPC_ call NUOPC_CompAttributeAdd(getmComp,(/"InitializePhaseMap"/),rc=rc) -#else - call ESMF_AttributeAdd(getmComp,convention="NUOPC",purpose="General", & - attrList=(/"InitializePhaseMap"/),rc=rc) -#endif abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) -#if 1 call NUOPC_CompAttributeSet(getmComp,"InitializePhaseMap", & InitializePhaseMap,rc=rc) + abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) + if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) #else +! Note (KK): NUOPC attributes are purpose="Instance" + call ESMF_AttributeAdd(getmComp,convention="NUOPC",purpose="General", & + attrList=(/"InitializePhaseMap"/),rc=rc) + abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) + if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) + call ESMF_AttributeSet(getmComp,name="InitializePhaseMap", & valueList=InitializePhaseMap, & convention="NUOPC",purpose="General",rc=rc) -#endif abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) +#endif call toc(TIM_ESMF) @@ -531,7 +543,6 @@ ! ! !LOCAL VARIABLES type(ESMF_Clock) :: getmClock - type(ESMF_Grid) :: getmGrid2D type(ESMF_Time) :: getmRefTime,getmStartTime,getmStopTime type(ESMF_TimeInterval) :: getmTimeStep type(ESMF_VM) :: vm @@ -708,9 +719,10 @@ end if - call getmComp_init_grid(getmComp,getmGrid2D) - call init_importStateP1(getmComp,getmGrid2D,importState) - call init_exportStateP1(getmComp,getmGrid2D,exportState) + call getmComp_init_grid(getmComp) + call NUOPC_FieldDictionarySetAutoAdd(.true.) + call init_importStateP1(getmComp,importState) + call init_exportStateP1(getmComp,exportState) call toc(TIM_ESMF) @@ -775,8 +787,11 @@ call ESMF_LogWrite(trim(name)//"::InitializeP2...",ESMF_LOGMSG_TRACE) + call init_importStateP2(getmComp,importState) call init_exportStateP2(getmComp,exportState) + call read_importState(getmComp,importState) + ! If the initial Export state needs to be filled, do it here. call getmComp_update_grid(getmComp) call update_exportState(getmComp,exportState) @@ -1048,7 +1063,7 @@ ! !ROUTINE: getmComp_init_grid - Creates Grid ! ! !INTERFACE: - subroutine getmComp_init_grid(getmComp,getmGrid2D) + subroutine getmComp_init_grid(getmComp) ! ! !DESCRIPTION: ! @@ -1067,9 +1082,6 @@ ! !INPUT/OUTPUT PARAMETERS: type(ESMF_GridComp) :: getmComp ! -! !OUTPUT PARAMETERS: - type(ESMF_Grid),intent(out) :: getmGrid2D -! ! !REVISION HISTORY: ! Original Author(s): Knut Klingbeil ! @@ -1089,7 +1101,7 @@ type(ESMF_Array) :: array type(ESMF_CoordSys_Flag) :: coordSys type(ESMF_DistGrid) :: getmDistGrid2D,getmDistGrid3D,distgrid - type(ESMF_Grid) :: getmGrid3D + type(ESMF_Grid) :: getmGrid3D,getmGrid2D type(ESMF_StaggerLoc) :: staggerloc type(ESMF_VM) :: vm ! Note (KK): ESMF_ARRAY's are deep classes, that persist after return. @@ -1282,12 +1294,7 @@ ! staggerAlign and staggerEdgeWidth's. ! internal call to ESMF_GridCreateFrmDistGrid() getmGrid2D = ESMF_GridCreate(getmDistGrid2D,name=trim(name)//"Grid2D", & -#if 0 -! bug in ESMF gridAlign=(/1,1/), & -#else - gridEdgeLWidth=(/1,1/), & -#endif coordSys=coordSys, & coordDimCount=int(coordDimCount(1:2)), & coordDimMap=int(coordDimMap(1:2,1:2)), & @@ -1296,13 +1303,7 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) getmGrid3D = ESMF_GridCreate(getmDistGrid3D,name=trim(name)//"Grid3D", & -#if 0 -! bug in ESMF gridAlign=(/1,1,1/), & -#else - gridEdgeLWidth=(/1,1,1/), & -#endif - coordSys=coordSys, & coordDimCount=coordDimCount, & coordDimMap=coordDimMap, & @@ -1858,22 +1859,20 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) ! !ROUTINE: init_importStateP1 - ! ! !INTERFACE: - subroutine init_importStateP1(getmComp,getmGrid2D,importState) + subroutine init_importStateP1(getmComp,importState) ! ! !DESCRIPTION: ! ! !USES: + use initialise ,only: runtype use domain ,only: grid_type use meteo ,only: met_method,calc_met,METEO_FROMEXT - use meteo ,only: airp,u10,v10,constant_cd + use meteo ,only: airp,u10,v10,t2,tcc use waves ,only: waveforcing_method,WAVES_FROMEXT use waves ,only: waves_ramp use variables_waves,only: waveH,waveK,waveT IMPLICIT NONE ! -! !INPUT PARAMETERS: - type(ESMF_Grid),intent(in) :: getmGrid2D -! ! !INPUT/OUTPUT PARAMETERS: type(ESMF_GridComp) :: getmComp type(ESMF_State) :: importState @@ -1882,7 +1881,8 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) ! Original Author(s): Knut Klingbeil ! ! !LOCAL VARIABLES - type(ESMF_Grid) :: getmGrid3D + type(ESMF_Grid) :: getmGrid3D,getmGrid2D + type(ESMF_Grid), allocatable :: gridList(:) integer :: rc logical :: abort,frc ! @@ -1895,24 +1895,36 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) write(debug,*) 'init_importStateP1() # ',Ncall #endif - call ESMF_GridCompGet(getmComp,grid=getmGrid3D,rc=rc) + call ESMF_GridCompGet(getmComp, gridList=gridList, rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) - call StateAddField(importState,"gridSetField2D",getmGrid2D) + getmGrid3D = gridList(1) + getmGrid2D = gridList(2) if (met_method .eq. METEO_FROMEXT) then !call StateAddField(importState,trim(name_swr ),getmGrid2D,units="W m-2") call StateAddField(importState,trim(name_slp ),getmGrid2D, & farray2D=airp,units="Pa") if (calc_met) then - constant_cd = .true. ! force allocation of new memory if grid rotation needs to be removed frc = (grid_type .ne. 2) call StateAddField(importState,trim(name_windU ),getmGrid2D, & farray2D=u10,units="m s-1",frc=frc) call StateAddField(importState,trim(name_windV ),getmGrid2D, & farray2D=v10,units="m s-1",frc=frc) + if (runtype .gt. 2) then + call StateAddField(importState,trim(name_airT2 ),getmGrid2D, & + farray2D=t2,units="K") + call StateAddField(importState,trim(name_hums ),getmGrid2D, & + units="kg/kg") + call StateAddField(importState,trim(name_humr ),getmGrid2D, & + units="%") + call StateAddField(importState,trim(name_dev2 ),getmGrid2D, & + units="K") + call StateAddField(importState,trim(name_tcc ),getmGrid2D, & + farray2D=tcc,units="") + end if end if ! calc_met end if ! meteo @@ -1929,6 +1941,8 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) farray2D=waveT,units="s") end if + !call ESMF_StatePrint(importState) + #ifdef DEBUG write(debug,*) 'Leaving init_importStateP1()' write(debug,*) @@ -1940,6 +1954,85 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) !----------------------------------------------------------------------- !BOP ! +! !ROUTINE: init_importStateP2 - +! +! !INTERFACE: + subroutine init_importStateP2(getmComp,importState) +! +! !DESCRIPTION: +! +! !USES: + use initialise ,only: runtype + use meteo ,only: met_method,calc_met,METEO_FROMEXT + use meteo ,only: hum + use meteo ,only: hum_method,RELATIVE_HUM,WET_BULB,DEW_POINT,SPECIFIC_HUM + IMPLICIT NONE +! +! !INPUT/OUTPUT PARAMETERS: + type(ESMF_GridComp) :: getmComp + type(ESMF_State) :: importState +! +! !REVISION HISTORY: +! Original Author(s): Knut Klingbeil +! +! !LOCAL VARIABLES + type(ESMF_Grid) :: getmGrid3D,getmGrid2D + type(ESMF_Grid), allocatable :: gridList(:) + real(ESMF_KIND_R8), pointer :: p2dr(:,:) + integer :: rc + logical :: abort,frc +! +!EOP +!----------------------------------------------------------------------- +!BOC +#ifdef DEBUG + integer, save :: Ncall = 0 + Ncall = Ncall+1 + write(debug,*) 'init_importStateP2() # ',Ncall +#endif + + if (met_method .eq. METEO_FROMEXT) then + if (calc_met) then + if (runtype .gt. 2) then + p2dr => NULL() + call StateCompleteConnectedField(importState,trim(name_hums ), & + farray2d=hum,p2dr=p2dr) + if (associated(p2dr)) then + hum_method = SPECIFIC_HUM + else + call StateCompleteConnectedField(importState,trim(name_humr ), & + farray2d=hum,p2dr=p2dr) + if (associated(p2dr)) then + hum_method = RELATIVE_HUM + else + call StateCompleteConnectedField(importState,trim(name_dev2 ), & + farray2d=hum,p2dr=p2dr) + if (associated(p2dr)) then + hum_method = DEW_POINT + else + hum_method = -1 + call ESMF_LogWrite('hum_method=-1', & + ESMF_LOGMSG_WARNING,line=__LINE__,file=FILENAME) + end if + end if + end if + end if + end if ! calc_met + end if ! meteo + + !call ESMF_StatePrint(importState) + +#ifdef DEBUG + write(debug,*) 'Leaving init_importStateP2()' + write(debug,*) +#endif + return + + end subroutine init_importStateP2 +!EOC +!----------------------------------------------------------------------- +!BOP +! ! !ROUTINE: read_importState - ! ! !INTERFACE: @@ -1948,9 +2041,11 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) ! !DESCRIPTION: ! ! !USES: + use initialise ,only: runtype use domain ,only: grid_type,convc,cosconv,sinconv use meteo ,only: met_method,calc_met,METEO_FROMEXT,new_meteo - use meteo ,only: airp,u10,v10 + use meteo ,only: hum_method,RELATIVE_HUM,WET_BULB,DEW_POINT,SPECIFIC_HUM + use meteo ,only: airp,u10,v10,t2,hum,tcc use waves ,only: waveforcing_method,WAVES_FROMEXT,new_waves use waves ,only: waves_ramp use variables_waves,only: coswavedir,sinwavedir,waveH,waveK,waveT @@ -1997,6 +2092,25 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) u10 = cosconv*windU - sinconv*windV v10 = sinconv*windU + cosconv*windV end if + if (runtype .gt. 2) then + call StateReadCompleteField(importState,trim(name_airT2 ),& + farray2d=t2) + + select case(hum_method) + case (SPECIFIC_HUM) + call StateReadCompleteField(importState,trim(name_hums ),& + farray2d=hum) + case (RELATIVE_HUM) + call StateReadCompleteField(importState,trim(name_humr ),& + farray2d=hum) + case (DEW_POINT) + call StateReadCompleteField(importState,trim(name_dev2 ),& + farray2d=hum) + end select + + call StateReadCompleteField(importState,trim(name_tcc ),& + farray2d=tcc) + end if end if ! calc_met end if ! meteo @@ -2035,7 +2149,7 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) ! !ROUTINE: init_exportStateP1 - ! ! !INTERFACE: - subroutine init_exportStateP1(getmComp,getmGrid2D,exportState) + subroutine init_exportStateP1(getmComp,exportState) ! ! !DESCRIPTION: ! @@ -2049,9 +2163,6 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) use waves ,only: WAVES_FROMWIND,WAVES_FROMFILE IMPLICIT NONE ! -! !INPUT PARAMETERS: - type(ESMF_Grid),intent(in) :: getmGrid2D -! ! !INPUT/OUTPUT PARAMETERS: type(ESMF_GridComp) :: getmComp type(ESMF_State) :: exportState @@ -2060,7 +2171,8 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) ! Original Author(s): Knut Klingbeil ! ! !LOCAL VARIABLES - type(ESMF_Grid) :: getmGrid3D + type(ESMF_Grid) :: getmGrid3D,getmGrid2D + type(ESMF_Grid), allocatable :: gridList(:) integer :: rc logical :: abort ! @@ -2073,60 +2185,61 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) write(debug,*) 'init_exportStateP1() # ',Ncall #endif - call ESMF_GridCompGet(getmComp,grid=getmGrid3D,rc=rc) + call ESMF_GridCompGet(getmComp, gridList=gridList, rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) - call StateAddField(exportState,"gridSetField2D",getmGrid2D) + getmGrid3D = gridList(1) + getmGrid2D = gridList(2) -#ifdef ALLEXPORT - call StateAddField(exportState,trim(name_depth ),getmGrid2D,units="m") - call StateAddField(exportState,trim(name_h3D ),getmGrid3D,units="m") - call StateAddField(exportState,trim(name_hbot ),getmGrid2D,units="m") - call StateAddField(exportState,trim(name_U2D ),getmGrid2D,units="m s-1") - call StateAddField(exportState,trim(name_U3D ),getmGrid3D,units="m s-1") - call StateAddField(exportState,trim(name_Ubot ),getmGrid2D,units="m s-1") - call StateAddField(exportState,trim(name_V2D ),getmGrid2D,units="m s-1") - call StateAddField(exportState,trim(name_V3D ),getmGrid3D,units="m s-1") - call StateAddField(exportState,trim(name_Vbot ),getmGrid2D,units="m s-1") + !call StateAddField(exportState,trim(name_depth ),getmGrid2D,units="m") + !call StateAddField(exportState,trim(name_h3D ),getmGrid3D,units="m") + !call StateAddField(exportState,trim(name_hbot ),getmGrid2D,units="m") + !call StateAddField(exportState,trim(name_U2D ),getmGrid2D,units="m s-1") + !call StateAddField(exportState,trim(name_U3D ),getmGrid3D,units="m s-1") + !call StateAddField(exportState,trim(name_Ubot ),getmGrid2D,units="m s-1") + !call StateAddField(exportState,trim(name_V2D ),getmGrid2D,units="m s-1") + !call StateAddField(exportState,trim(name_V3D ),getmGrid3D,units="m s-1") + !call StateAddField(exportState,trim(name_Vbot ),getmGrid2D,units="m s-1") if (runtype .ge. 2) then #ifndef NO_3D - call StateAddField(exportState,trim(name_eps3D ),getmGrid3D,units="m2 s-3",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) - call StateAddField(exportState,trim(name_epsbot ),getmGrid2D,units="m2 s-3") - call StateAddField(exportState,trim(name_num3D ),getmGrid3D,units="m2 s-1",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) - call StateAddField(exportState,trim(name_numbot ),getmGrid2D,units="m2 s-1") - call StateAddField(exportState,trim(name_SS3D ),getmGrid3D,units="s-2",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) - call StateAddField(exportState,trim(name_taubmax),getmGrid2D,units="Pa") - call StateAddField(exportState,trim(name_tke3D ),getmGrid3D,units="m2 s-2",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) - call StateAddField(exportState,trim(name_tkebot ),getmGrid2D,units="m2 s-2") + !call StateAddField(exportState,trim(name_eps3D ),getmGrid3D,units="m2 s-3",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) + !call StateAddField(exportState,trim(name_epsbot ),getmGrid2D,units="m2 s-3") + !call StateAddField(exportState,trim(name_num3D ),getmGrid3D,units="m2 s-1",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) + !call StateAddField(exportState,trim(name_numbot ),getmGrid2D,units="m2 s-1") + !call StateAddField(exportState,trim(name_SS3D ),getmGrid3D,units="s-2",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) + !call StateAddField(exportState,trim(name_taubmax),getmGrid2D,units="Pa") + !call StateAddField(exportState,trim(name_tke3D ),getmGrid3D,units="m2 s-2",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) + !call StateAddField(exportState,trim(name_tkebot ),getmGrid2D,units="m2 s-2") if (runtype .ge. 3) then #ifndef NO_BAROCLINIC - call StateAddField(exportState,trim(name_NN3D ),getmGrid3D,units="s-2",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) - call StateAddField(exportState,trim(name_nuh3D ),getmGrid3D,units="m",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) + !call StateAddField(exportState,trim(name_NN3D ),getmGrid3D,units="s-2",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) + !call StateAddField(exportState,trim(name_nuh3D ),getmGrid3D,units="m",staggerloc=ESMF_STAGGERLOC_CENTER_VFACE) if (calc_salt) then - call StateAddField(exportState,trim(name_S3D ),getmGrid3D,units="") + !call StateAddField(exportState,trim(name_S3D ),getmGrid3D,units="") end if if (calc_temp) then - call StateAddField(exportState,trim(name_T3D ),getmGrid3D,units="degC") - call StateAddField(exportState,trim(name_Tbot ),getmGrid2D,units="degC") + !call StateAddField(exportState,trim(name_T3D ),getmGrid3D,units="degC") + !call StateAddField(exportState,trim(name_Tbot ),getmGrid2D,units="degC") end if #endif end if #endif end if if (met_method.eq.METEO_CONST .or. met_method.eq.METEO_FROMFILE) then - call StateAddField(exportState,trim(name_swr ),getmGrid2D,units="W m-2") - call StateAddField(exportState,trim(name_windU ),getmGrid2D,units="m s-1") - call StateAddField(exportState,trim(name_windV ),getmGrid2D,units="m s-1") + !call StateAddField(exportState,trim(name_swr ),getmGrid2D,units="W m-2") + !call StateAddField(exportState,trim(name_windU ),getmGrid2D,units="m s-1") + !call StateAddField(exportState,trim(name_windV ),getmGrid2D,units="m s-1") end if if ( waveforcing_method.eq.WAVES_FROMWIND & .or. waveforcing_method.eq.WAVES_FROMFILE) then - call StateAddField(exportState,trim(name_waveDir),getmGrid2D,units="rad") - call StateAddField(exportState,trim(name_waveH ),getmGrid2D,units="m") - call StateAddField(exportState,trim(name_waveK ),getmGrid2D,units="m-1") - call StateAddField(exportState,trim(name_waveT ),getmGrid2D,units="s") + !call StateAddField(exportState,trim(name_waveDir),getmGrid2D,units="rad") + !call StateAddField(exportState,trim(name_waveH ),getmGrid2D,units="m") + !call StateAddField(exportState,trim(name_waveK ),getmGrid2D,units="m-1") + !call StateAddField(exportState,trim(name_waveT ),getmGrid2D,units="s") end if -#endif + + !call ESMF_StatePrint(exportState) #ifdef DEBUG write(debug,*) 'Leaving init_exportStateP1()' @@ -2186,7 +2299,6 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) write(debug,*) 'init_exportStateP2() # ',Ncall #endif -#ifdef ALLEXPORT ! force allocation of new memory if Stokes drift needs to be removed frc = (waveforcing_method.ne.NO_WAVES .and. waves_method.ne.WAVES_NOSTOKES) @@ -2363,7 +2475,8 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) call StateCompleteConnectedField(exportState,trim(name_waveT ), & farray2D=waveT) end if ! waves -#endif + + !call ESMF_StatePrint(exportState) #ifdef DEBUG write(debug,*) 'Leaving init_exportStateP2()' @@ -2433,7 +2546,6 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) write(debug,*) 'update_exportState() # ',Ncall #endif -#ifdef ALLEXPORT ! force allocation of new memory if Stokes drift needs to be removed frc = (waveforcing_method.ne.NO_WAVES .and. waves_method.ne.WAVES_NOSTOKES) @@ -2705,7 +2817,6 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) call StateSetCompleteField(exportState,trim(name_waveT ), & farray2D=waveT) end if ! waves -#endif #ifdef DEBUG write(debug,*) 'Leaving update_exportState()' @@ -2857,8 +2968,6 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) xxArray2D = ESMF_ArrayCreate(distgrid,xx1D, & indexflag=ESMF_INDEX_DELOCAL, & -! totalLWidth=(/HALO+1/), & - totalUWidth=(/HALO /), & rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) @@ -2866,8 +2975,6 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) yxArray2D = ESMF_ArrayCreate(distgrid,yx1D, & indexflag=ESMF_INDEX_DELOCAL, & distgridToArrayMap=(/0,1/), & -! totalLWidth=(/HALO+1/), & - totalUWidth=(/HALO /), & rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) @@ -2900,8 +3007,6 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) xxArray3D = ESMF_ArrayCreate(distgrid,xx1D, & indexflag=ESMF_INDEX_DELOCAL, & -! totalLWidth=(/HALO+1/), & - totalUWidth=(/HALO /), & rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) @@ -2909,8 +3014,6 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) yxArray3D = ESMF_ArrayCreate(distgrid,yx1D, & indexflag=ESMF_INDEX_DELOCAL, & distgridToArrayMap=(/0,1,0/), & -! totalLWidth=(/HALO+1/), & - totalUWidth=(/HALO /), & rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) @@ -3002,16 +3105,12 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) xxArray2D = ESMF_ArrayCreate(distgrid,xx2D, & indexflag=ESMF_INDEX_DELOCAL, & -! totalLWidth=(/HALO+1,HALO+1/), & - totalUWidth=(/HALO ,HALO /), & rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) yxArray2D = ESMF_ArrayCreate(distgrid,yx2D, & indexflag=ESMF_INDEX_DELOCAL, & -! totalLWidth=(/HALO+1,HALO+1/), & - totalUWidth=(/HALO ,HALO /), & rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) @@ -3043,16 +3142,12 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) xxArray3D = ESMF_ArrayCreate(distgrid,xx2D, & indexflag=ESMF_INDEX_DELOCAL, & -! totalLWidth=(/HALO+1,HALO+1/), & - totalUWidth=(/HALO ,HALO /), & rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) yxArray3D = ESMF_ArrayCreate(distgrid,yx2D, & indexflag=ESMF_INDEX_DELOCAL, & -! totalLWidth=(/HALO+1,HALO+1/), & - totalUWidth=(/HALO ,HALO /), & rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) @@ -3203,6 +3298,31 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) return end if +#ifdef _GETM_NUOPC_ + if (NUOPC_FieldDictionaryHasEntry(trim(name))) then + call NUOPC_Advertise(state, trim(name), Units=units, rc=rc) + abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) + if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) + + call ESMF_StateGet(state, trim(name), field=field, rc=rc) + abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) + if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) + else + call NUOPC_Advertise(state, trim(name), rc=rc) + abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) + if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) + + call ESMF_StateGet(state, trim(name), field=field, rc=rc) + abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) + if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) + + if (present(units)) then + call ESMF_AttributeSet(field,'units',trim(units),rc=rc) + abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) + if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) + end if + end if +#else field = ESMF_FieldEmptyCreate(name=trim(name),rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) @@ -3221,6 +3341,7 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) call NUOPC_SetAttribute(field,"Connected","false",rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) +#endif call ESMF_FieldEmptySet(field,grid,staggerloc=staggerloc,rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) @@ -3301,9 +3422,11 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) end if +#ifndef _GETM_NUOPC_ call ESMF_StateAdd(state,(/field/),rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) +#endif #ifdef DEBUG write(debug,*) 'Leaving StateAddField()' @@ -3355,6 +3478,7 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) type(ESMF_FieldStatus_Flag) :: status type(ESMF_Grid) :: grid type(ESMF_StaggerLoc) :: staggerloc + type(ESMF_StateItem_Flag) :: itemType real(ESMF_KIND_R8),pointer :: p2dr_(:,:),p3dr_(:,:,:) REALTYPE :: getmreal integer :: rc,dimCount,klen @@ -3370,6 +3494,13 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) write(debug,*) 'StateCompleteConnectedField() # ',Ncall #endif + call ESMF_StateGet(state,trim(name),itemType,rc=rc) + abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) + if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) + +! skip non-existing fields + if (itemType .eq. ESMF_STATEITEM_NOTFOUND) return + call ESMF_StateGet(state,trim(name),field,rc=rc) abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) @@ -3505,6 +3636,10 @@ if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) end if + call NUOPC_Realize(state, field, rc) + abort = ESMF_LogFoundError(rc,line=__LINE__,file=FILENAME) + if (abort) call ESMF_Finalize(endflag=ESMF_END_ABORT) + #ifdef DEBUG write(debug,*) 'Leaving StateCompleteConnectedField()' write(debug,*) diff --git a/src/getm/initialise.F90 b/src/getm/initialise.F90 index cc29f56852df97408320619c56eac10ef91955c2..261e51f7dde4f3a74503c9dd1281bcd4a3fa7fb0 100644 --- a/src/getm/initialise.F90 +++ b/src/getm/initialise.F90 @@ -360,10 +360,11 @@ #endif end if #endif - call init_output_processing() call init_les(runtype) + call init_output_processing() + call init_register_all_variables(runtype) allocate(type_getm_host::output_manager_host) diff --git a/src/getm/register_all_variables.F90 b/src/getm/register_all_variables.F90 index 58d646bb8eba4a8ff364510206c52a9dcfe7037e..1b902b130f2edfe1b9fc9fb3d0cd6e506f66c7be 100644 --- a/src/getm/register_all_variables.F90 +++ b/src/getm/register_all_variables.F90 @@ -11,6 +11,12 @@ ! ! !USES: use field_manager + use variables_2d, only: register_2d_variables + use variables_3d, only: register_3d_variables +#ifdef _FABM_ + use getm_fabm, only: register_fabm_variables +#endif + use output_processing, only: register_processed_variables, finalize_register_processed_variables IMPLICIT NONE ! ! default: all is private. @@ -229,6 +235,7 @@ #ifdef _FABM_ call finalize_register_fabm_variables(fm) #endif + call finalize_register_processed_variables(fm) return end subroutine finalize_register_all_variables diff --git a/src/les/les_smagorinsky_sym.F90 b/src/les/les_smagorinsky_sym.F90 new file mode 100644 index 0000000000000000000000000000000000000000..068bc8ed3b6b2beaa42e18931a4211097d83db0c --- /dev/null +++ b/src/les/les_smagorinsky_sym.F90 @@ -0,0 +1,225 @@ +#include "cppdefs.h" +!----------------------------------------------------------------------- +!BOP +! +! !ROUTINE: les_smagorinsky_sym - \label{les_smagorinsky_sym} +! +! !INTERFACE: + subroutine les_smagorinsky_sym(dudxC,dudxV, & +#ifndef SLICE_MODEL + dvdyC,dvdyU, & +#endif + shearX,shearU, & + AmC,AmX,AmU,AmV) + +! Note (KK): keep in sync with interface in les.F90 +! +! !DESCRIPTION: +! +! +! !USES: + use variables_les, only: SmagC2_2d,SmagX2_2d,SmagU2_2d,SmagV2_2d + use domain, only: imin,imax,jmin,jmax,az,ax,au,av + use domain, only: dxc,dyc,dxx,dyx,dxu,dyu,dxv,dyv + use getm_timers, only: tic,toc,TIM_SMAG2D +!$ use omp_lib + IMPLICIT NONE +! +! !INPUT PARAMETERS: + REALTYPE,dimension(E2DFIELD),intent(in) :: dudxC,dudxV +#ifndef SLICE_MODEL + REALTYPE,dimension(E2DFIELD),intent(in) :: dvdyC,dvdyU +#endif + REALTYPE,dimension(E2DFIELD),intent(in) :: shearX,shearU +! +! !OUTPUT PARAMETERS: + REALTYPE,dimension(E2DFIELD),intent(out),optional :: AmC,AmX,AmU,AmV +! +! !REVISION HISTORY: +! Original author(s): Knut Klingbeil +! +! !LOCAL VARIABLES: + REALTYPE :: dudx,dvdy + integer :: i,j + +!EOP +!----------------------------------------------------------------------- +!BOC +#ifdef DEBUG + integer, save :: Ncall = 0 + Ncall = Ncall+1 + write(debug,*) 'les_smagorinsky_sym() # ',Ncall +#endif +#ifdef SLICE_MODEL + j = jmax/2 ! this MUST NOT be changed!!! +#endif + call tic(TIM_SMAG2D) + +!$OMP PARALLEL DEFAULT(SHARED) & +!$OMP FIRSTPRIVATE(j) & +!$OMP PRIVATE(i,dudx,dvdy) + + if (present(AmC)) then +!$OMP DO SCHEDULE(RUNTIME) +#ifndef SLICE_MODEL + do j=jmin-1,jmax+1 +#endif + do i=imin-1,imax+1 +! Note (KK): AmC(az=1) needed in uv_diffusion + if (az(i,j) .eq. 1) then +! interpolate shearC +! Note (KK): in W/E open boundary cells shearC(az=2) would +! require shearU outside open boundary +! in N/S open boundary cells shearC(az=2) would +! require shearU(au=3) +! however shearC(az=2) not needed + AmC(i,j) = ( & + dudxC(i,j) & +#ifndef SLICE_MODEL + - dvdyC(i,j) & +#endif + )**2 & + + (_HALF_*(shearU(i-1,j) + shearU(i,j)))**2 + AmC(i,j) = SmagC2_2d(i,j)*DXC*DYC*sqrt(AmC(i,j)) + end if + end do +#ifndef SLICE_MODEL + end do +#endif +!$OMP END DO NOWAIT + +#ifdef SLICE_MODEL +!$OMP BARRIER +!$OMP SINGLE + AmC(imin-1:imax+1,j+1) = AmC(imin-1:imax+1,j) +!$OMP END SINGLE +#endif + + end if + + if (present(AmX)) then +!$OMP DO SCHEDULE(RUNTIME) +#ifndef SLICE_MODEL + do j=jmin-1,jmax +#endif + do i=imin-1,imax + if (ax(i,j) .eq. 1) then +! interpolate dudxX and dvdyX + if (av(i,j).eq.3 .or. av(i+1,j).eq.3) then +! Note (KK): western/eastern open bdy + dudx = _ZERO_ + else + dudx = _HALF_*(dudxV(i,j) + dudxV(i+1,j )) + end if +#ifndef SLICE_MODEL + if (au(i,j).eq.3 .or. au(i,j+1).eq.3) then +! Note (KK): northern/southern open bdy + dvdy = _ZERO_ + else + dvdy = _HALF_*(dvdyU(i,j) + dvdyU(i ,j+1)) + end if +#endif + AmX(i,j) = ( & + dudx & +#ifndef SLICE_MODEL + - dvdy & +#endif + )**2 & + + shearX(i,j)**2 + AmX(i,j) = SmagX2_2d(i,j)*DXX*DYX*sqrt(AmX(i,j)) + end if + end do +#ifndef SLICE_MODEL + end do +#endif +!$OMP END DO NOWAIT + +#ifdef SLICE_MODEL +!$OMP BARRIER +!$OMP SINGLE + AmX(imin-1:imax,j-1) = AmX(imin-1:imax,j) + AmX(imin-1:imax,j+1) = AmX(imin-1:imax,j) +!$OMP END SINGLE +#endif + + end if + + if (present(AmU)) then +!$OMP DO SCHEDULE(RUNTIME) +#ifndef SLICE_MODEL + do j=jmin-1,jmax+1 +#endif + do i=imin-1,imax +! Note (KK): shearU(au=3) not available +! (however we only need AmU(au=[1|2]) in tracer_diffusion) + if(au(i,j).eq.1 .or. au(i,j).eq.2) then +! interpolate dudxU (see deformation_rates) + if (au(i,j) .eq. 1) then + dudx = _HALF_*(dudxC(i,j) + dudxC(i+1,j)) + else + dudx = _ZERO_ + end if + AmU(i,j) = ( & + dudx & +#ifndef SLICE_MODEL + - dvdyU(i,j) & +#endif + )**2 & + + shearU(i,j)**2 + AmU(i,j) = SmagU2_2d(i,j)*DXU*DYU*sqrt(AmU(i,j)) + end if + end do +#ifndef SLICE_MODEL + end do +#endif +!$OMP END DO NOWAIT + +#ifdef SLICE_MODEL +!$OMP BARRIER +!$OMP SINGLE + AmU(imin-1:imax,j+1) = AmU(imin-1:imax,j) +!$OMP END SINGLE +#endif + + end if + +#ifndef SLICE_MODEL + if (present(AmV)) then +!$OMP DO SCHEDULE(RUNTIME) + do j=jmin-1,jmax + do i=imin-1,imax+1 +! Note (KK): we only need AmV(av=[1|2]) in tracer_diffusion +! (shearV(av=3) cannot be calculated anyway) + if(av(i,j).eq.1 .or. av(i,j).eq.2) then +! interpolate dvdyV and shearV (see deformation_rates) + if (av(i,j) .eq. 1) then + dvdy = _HALF_*(dvdyC(i,j) + dvdyC(i,j+1)) + else + dvdy = _ZERO_ + end if + AmV(i,j) = ( & + dudxV(i,j) & + - dvdy & + )**2 & + + (_HALF_*(shearX(i-1,j) + shearX(i,j)))**2 + AmV(i,j) = SmagV2_2d(i,j)*DXV*DYV*sqrt(AmV(i,j)) + end if + end do + end do +!$OMP END DO NOWAIT + end if +#endif + +!$OMP END PARALLEL + + call toc(TIM_SMAG2D) +#ifdef DEBUG + write(debug,*) 'Leaving les_smagorinsky_sym()' + write(debug,*) +#endif + return + end subroutine les_smagorinsky_sym +!EOC +!----------------------------------------------------------------------- +! Copyright (C) 2011 - Hans Burchard and Karsten Bolding ! +!----------------------------------------------------------------------- diff --git a/src/ncdf/ncdf_meteo.F90 b/src/ncdf/ncdf_meteo.F90 index c9726f3c6f79aafe678e708a2b8a35944e8e9813..fcb023e987b6faf15588e5894fae97ecd9769ccb 100644 --- a/src/ncdf/ncdf_meteo.F90 +++ b/src/ncdf/ncdf_meteo.F90 @@ -148,14 +148,8 @@ "dimensions do not match") end if il = ilg ; jl = jlg ; ih = ihg ; jh = jhg - else - il = 1 ; jl = 1 ; ih = iextr ; jh = jextr end if - start(1) = il; start(2) = jl; - edges(1) = ih-il+1; edges(2) = jh-jl+1; - edges(3) = 1 - allocate(beta(E2DFIELD),stat=err) if (err /= 0) & stop 'init_meteo_input_ncdf: Error allocating memory (beta)' @@ -174,6 +168,7 @@ call to_rotated_lat_lon(southpole,olon,olat,rlon,rlat,x) beta = x end if + il = 1 ; jl = 1 ; ih = iextr ; jh = jextr else if (met_lat(1) .gt. met_lat(2)) then LEVEL3 'Reverting lat-axis and setting grid_scan to 0' @@ -228,8 +223,21 @@ call getm_error("init_meteo_input_ncdf()", & "Some interpolation coefficients are not valid") end if + + il = minval(gridmap(:,:,1),mask=(gridmap(:,:,1).gt.0)) + jl = minval(gridmap(:,:,2),mask=(gridmap(:,:,2).gt.0)) + ih = min( maxval(gridmap(:,:,1))+1 , iextr ) + jh = min( maxval(gridmap(:,:,2))+1 , jextr ) + + where( gridmap(:,:,1).gt.0 ) gridmap(:,:,1) = gridmap(:,:,1) - il + 1 + where( gridmap(:,:,2).gt.0 ) gridmap(:,:,2) = gridmap(:,:,2) - jl + 1 + end if + start(1) = il; start(2) = jl; + edges(1) = ih-il+1; edges(2) = jh-jl+1; + edges(3) = 1 + allocate(d_airp(E2DFIELD),stat=rc) if (rc /= 0) stop 'init_meteo_input_ncdf: Error allocating memory (d_airp)' d_airp = -9999*_ONE_ @@ -846,7 +854,7 @@ STDERR 'grid_north_pole_longitude ',southpole(2) ! ! !LOCAL VARIABLES: integer :: i,j,err - REALTYPE,dimension(edges(1),edges(2)) :: wrk,wrk_dp + REALTYPE,dimension(edges(1),edges(2)) :: wrk!,wrk_dp REALTYPE :: angle,uu,vv,sinconv,cosconv !EOP !----------------------------------------------------------------------- @@ -866,8 +874,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) end if else !KBKwrk_dp = _ZERO_ - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,airp_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,airp_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,airp_input) end if end if @@ -881,8 +891,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) evap_input(ill:ihl,jll:jhl) = wrk end if else - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,evap_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,evap_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,evap_input) end if if (evap_factor .ne. _ONE_) then evap_input = evap_input * evap_factor @@ -899,8 +911,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) precip_input(ill:ihl,jll:jhl) = wrk end if else - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,precip_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,precip_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,precip_input) end if if (precip_factor .ne. _ONE_) then precip_input = precip_input * precip_factor @@ -920,8 +934,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) end if else !KBKwrk_dp = _ZERO_ - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,u10_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,u10_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,u10_input) end if end if @@ -936,8 +952,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) end if else !KBKwrk_dp = _ZERO_ - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,v10_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,v10_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,v10_input) end if end if @@ -976,8 +994,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) end if else !KBKwrk_dp = _ZERO_ - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,t2_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,t2_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,t2_input) end if end if @@ -992,8 +1012,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) end if else !KBKwrk_dp = _ZERO_ - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,hum_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,hum_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,hum_input) end if end if @@ -1008,8 +1030,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) end if else !KBKwrk_dp = _ZERO_ - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,tcc_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,tcc_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,tcc_input) end if end if @@ -1064,8 +1088,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) sst_input(ill:ihl,jll:jhl) = wrk end if else if (calc_met) then - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,sst_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,sst_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,sst_input) end if end if @@ -1079,8 +1105,10 @@ STDERR 'grid_north_pole_longitude ',southpole(2) sss_input(ill:ihl,jll:jhl) = wrk end if else if (calc_met) then - call copy_var(grid_scan,wrk,wrk_dp) - call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,sss_input) + !call copy_var(grid_scan,wrk,wrk_dp) + !call do_grid_interpol(az,wrk_dp,gridmap,ti,ui,sss_input) + call flip_var(wrk) + call do_grid_interpol(az,wrk,gridmap,ti,ui,sss_input) end if end if @@ -1094,6 +1122,35 @@ STDERR 'grid_north_pole_longitude ',southpole(2) end subroutine read_data !EOC +!----------------------------------------------------------------------- +!BOP +! +! !IROUTINE: flip_var - +! +! !INTERFACE: + subroutine flip_var(var) + IMPLICIT NONE +! +! !DESCRIPTION: +! +! !INPUT/OUTPUT PARAMETERS: + REALTYPE, intent(inout) :: var(edges(1),edges(2)) +! +! !REVISION HISTORY: +! +! !LOCAL VARIABLES: +! +!EOP +!----------------------------------------------------------------------- + + select case (grid_scan) + case (0) + var = var(:,edges(2):1:-1) + end select + return + end subroutine flip_var +!EOC + !----------------------------------------------------------------------- !BOP ! @@ -1125,9 +1182,12 @@ STDERR 'grid_north_pole_longitude ',southpole(2) ! ! !LOCAL VARIABLES: integer :: i,j + integer :: iextr,jextr !EOP !----------------------------------------------------------------------- + iextr = edges(1) ; jextr = edges(2) + select case (grid_scan) case (0) do j=1,jextr diff --git a/src/output/output_processing.F90 b/src/output/output_processing.F90 index 76fdb744080c317aeff70089f5d84f8bf6d53cad..1029d7f445375e19843b380374f9cc5944eed232 100644 --- a/src/output/output_processing.F90 +++ b/src/output/output_processing.F90 @@ -18,18 +18,22 @@ private ! ! !PUBLIC DATA FUNCTIONS: - public init_output_processing, register_processed_variables, do_output_processing + public init_output_processing, do_output_processing + public register_processed_variables, finalize_register_processed_variables ! ! !PUBLIC DATA MEMBERS: +! +! !PRIVATE DATA MEMBERS: REALTYPE, dimension(:,:), allocatable, target :: u2d, v2d - REALTYPE, dimension(:,:), allocatable, target :: u2d_destag, v2d_destag REALTYPE, dimension(:,:,:), allocatable, target :: u3d, v3d + REALTYPE, dimension(:,:), allocatable, target :: u2d_destag, v2d_destag REALTYPE, dimension(:,:,:), allocatable, target :: u3d_destag, v3d_destag -! -! !PRIVATE DATA MEMBERS: - logical, target:: u2d_use, v2d_use + + logical :: u2d_used, v2d_used + logical :: u3d_used, v3d_used + logical, target :: u2d_now, v2d_now + logical, target :: u3d_now, v3d_now logical, target:: u2d_destag_use, v2d_destag_use - logical, target:: u3d_use, v3d_use logical, target:: u3d_destag_use, v3d_destag_use integer, parameter :: rk = kind(_ONE_) ! @@ -65,12 +69,7 @@ !EOP !------------------------------------------------------------------------- !BOC - allocate(u2d(E2DFIELD),stat=rc) - if (rc /= 0) stop 'init_output_processing: Error allocating memory (u2d)' - u2d = 0._rk - allocate(v2d(E2DFIELD),stat=rc) - if (rc /= 0) stop 'init_output_processing: Error allocating memory (v2d)' - v2d = 0._rk + allocate(u2d_destag(E2DFIELD),stat=rc) if (rc /= 0) stop 'init_output_processing: Error allocating memory (u2d_destag)' u2d_destag = 0._rk @@ -79,10 +78,6 @@ v2d_destag = 0._rk #if 0 - allocate(u3d(I3DFIELD),stat=rc) - if (rc /= 0) stop 'init_output_processing: Error allocating memory (u3d)' - allocate(v3d(I3DFIELD),stat=rc) - if (rc /= 0) stop 'init_output_processing: Error allocating memory (v3d)' allocate(u3d_destag(I3DFIELD),stat=rc) if (rc /= 0) stop 'init_output_processing: Error allocating memory (u3d_destag)' allocate(v3d_destag(I3DFIELD),stat=rc) @@ -117,8 +112,15 @@ !BOC LEVEL2 'register_processed_variables()' - call fm%register('u2d', 'm/s', 'velocity in local x-direction', standard_name='', data2d=u2d(_2D_W_), fill_value=-9999._rk, category='velocities', used_now=u2d_use) - call fm%register('v2d', 'm/s', 'velocity in local y-direction', standard_name='', data2d=v2d(_2D_W_), fill_value=-9999._rk, category='velocities', used_now=v2d_use) + call fm%register('u2d', 'm/s', 'velocity in local x-direction', standard_name='', fill_value=-9999._rk, category='velocities', output_level=output_level_debug, used=u2d_used, used_now=u2d_now) + call fm%register('v2d', 'm/s', 'velocity in local y-direction', standard_name='', fill_value=-9999._rk, category='velocities', output_level=output_level_debug, used=v2d_used, used_now=v2d_now) + +#ifndef NO_3D + call fm%register('u3d', 'm/s', 'velocity in local x-direction (3D)', standard_name='', dimensions=(/id_dim_z/), fill_value=-9999._rk, category='velocities', output_level=output_level_debug, used=u3d_used, used_now=u3d_now) + call fm%register('v3d', 'm/s', 'velocity in local y-direction (3D)', standard_name='', dimensions=(/id_dim_z/), fill_value=-9999._rk, category='velocities', output_level=output_level_debug, used=v3d_used, used_now=v3d_now) +#endif + + call fm%register('u2d-destag', 'm/s', 'velocity in local x-direction(destag)', standard_name='', data2d=u2d_destag(_2D_W_), fill_value=-9999._rk, category='velocities',output_level=output_level_debug, used_now=u2d_destag_use) call fm%register('v2d-destag', 'm/s', 'velocity in local y-direction(destag)', standard_name='', data2d=v2d_destag(_2D_W_), fill_value=-9999._rk, category='velocities',output_level=output_level_debug, used_now=v2d_destag_use) @@ -126,6 +128,67 @@ end subroutine register_processed_variables !EOC +!----------------------------------------------------------------------- +!BOP +! +! !ROUTINE: finalize_register_processed_variables() - send optional variables. +! +! !INTERFACE: + subroutine finalize_register_processed_variables(fm) +! +! !DESCRIPTION: +! +! !USES: + use field_manager + IMPLICIT NONE +! +! !INPUT PARAMETERS: + type (type_field_manager) :: fm +! +! !REVISION HISTORY: +! Original author(s): Knut Klingbeil +! +! !LOCAL VARIABLES: + integer :: rc +!EOP +!----------------------------------------------------------------------- +!BOC + LEVEL1 'finalize_register_processed_variables()' + + if (u2d_used) then + allocate(u2d(E2DFIELD),stat=rc) + if (rc /= 0) stop 'finalize_register_processed_variables: Error allocating memory (u2d)' + u2d = 0._rk + call fm%send_data('u2d', u2d(_2D_W_)) + end if + + if (v2d_used) then + allocate(v2d(E2DFIELD),stat=rc) + if (rc /= 0) stop 'finalize_register_processed_variables: Error allocating memory (v2d)' + v2d = 0._rk + call fm%send_data('v2d', v2d(_2D_W_)) + end if + +#ifndef NO_3D + if (u3d_used) then + allocate(u3d(I3DFIELD),stat=rc) + if (rc /= 0) stop 'finalize_register_processed_variables: Error allocating memory (u3d)' + u3d = 0._rk + call fm%send_data('u3d', u3d(_3D_W_)) + end if + + if (v3d_used) then + allocate(v3d(I3DFIELD),stat=rc) + if (rc /= 0) stop 'finalize_register_processed_variables: Error allocating memory (v3d)' + v3d = 0._rk + call fm%send_data('v3d', v3d(_3D_W_)) + end if +#endif + + return + end subroutine finalize_register_processed_variables +!EOC + !----------------------------------------------------------------------- !BOP ! !IROUTINE: do_output_processing - read required variables @@ -137,7 +200,9 @@ use domain, only: az, au, av use variables_2d, only: z,D use variables_2d, only: U,V,DU,DV +#ifndef NO_3D use variables_3d, only: kmin,hn,uu,hun,vv,hvn +#endif IMPLICIT NONE ! ! !DESCRIPTION: @@ -154,29 +219,40 @@ !BOC ! 2D - velocities - if (u2d_use .and. v2d_use) then + + if (u2d_now) then call to_2d_vel(imin,jmin,imax,jmax,au,U,DU,vel_missing, & imin,jmin,imax,jmax,u2d) + end if + + if (v2d_now) then call to_2d_vel(imin,jmin,imax,jmax,av,V,DV,vel_missing, & imin,jmin,imax,jmax,v2d) end if - if (u2d_destag_use .and. v2d_destag_use) then - call to_2d_u(imin,jmin,imax,jmax,az,U,DU,vel_missing, & - imin,jmin,imax,jmax,u2d_destag) - call to_2d_v(imin,jmin,imax,jmax,az,V,DV,vel_missing, & - imin,jmin,imax,jmax,v2d_destag) - end if -#if 0 ! 3D - velocities #ifndef NO_3D - if (allocated(u3d) .and. allocated(v3d)) then + if (u3d_now) then call to_3d_uu(imin,jmin,imax,jmax,kmin,kmax,az, & hun,uu,vel_missing,u3d) + end if + + if (v3d_now) then call to_3d_vv (imin,jmin,imax,jmax,kmin,kmax,az, & hvn,vv,vel_missing,v3d) end if +#endif + + + if (u2d_destag_use .and. v2d_destag_use) then + call to_2d_u(imin,jmin,imax,jmax,az,U,DU,vel_missing, & + imin,jmin,imax,jmax,u2d_destag) + call to_2d_v(imin,jmin,imax,jmax,az,V,DV,vel_missing, & + imin,jmin,imax,jmax,v2d_destag) + end if +#if 0 +#ifndef NO_3D if (allocated(u3d_destag) .and. allocated(v3d_destag)) then call to_3d_vel(imin,jmin,imax,jmax,kmin,kmax,au, & hun,uu,vel_missing,u3d_destag)