From b1ae0eec74309e477f5b75ba7ded653fa8148104 Mon Sep 17 00:00:00 2001 From: jswhit2 Date: Mon, 3 Jan 2022 21:12:00 +0000 Subject: [PATCH 1/5] fix for 4diau with iau_filter_increments=T --- tools/fv_iau_mod.F90 | 16 +++++++--------- 1 file changed, 7 insertions(+), 9 deletions(-) diff --git a/tools/fv_iau_mod.F90 b/tools/fv_iau_mod.F90 index 12b87e70f..c46534894 100644 --- a/tools/fv_iau_mod.F90 +++ b/tools/fv_iau_mod.F90 @@ -305,19 +305,17 @@ subroutine getiauforcing(IPD_Control,IAU_Data) return endif - t1=iau_state%hr1 - IPD_Control%iau_delthrs*0.5 - t2=iau_state%hr1 + IPD_Control%iau_delthrs*0.5 if (IPD_Control%iau_filter_increments) then ! compute increment filter weight - ! t1 beginning of window, t2 end of window + ! IPD_Control%iaufhrs(1) beginning of window, IPD_Control%iaufhrs(nfiles) end of window ! IPD_Control%fhour current time ! in window kstep=-nstep,nstep (2*nstep+1 total) ! time step IPD_control%dtp dtp=IPD_control%dtp nstep = 0.5*IPD_Control%iau_delthrs*3600/dtp ! compute normalized filter weight - kstep = (IPD_Control%fhour-(t1+IPD_Control%iau_delthrs*0.5))*3600./dtp - if (kstep .ge. -nstep .and. kstep .le. nstep) then + kstep = ((IPD_Control%fhour-IPD_Control%iaufhrs(1)) - 0.5*IPD_Control%iau_delthrs)*3600./dtp + if (IPD_Control%fhour >= IPD_Control%iaufhrs(1) .and. IPD_Control%fhour < IPD_Control%iaufhrs(nfiles)) then sx = acos(-1.)*kstep/nstep wx = acos(-1.)*kstep/(nstep+1) if (kstep .ne. 0) then @@ -326,7 +324,7 @@ subroutine getiauforcing(IPD_Control,IAU_Data) wt = 1. endif iau_state%wt = iau_state%wt_normfact*wt - if (is_master()) print *,'filter wt',kstep,IPD_Control%fhour,iau_state%wt + if (is_master()) print *,'filter wt',kstep,IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles),iau_state%wt/iau_state%wt_normfact else iau_state%wt = 0. endif @@ -335,12 +333,11 @@ subroutine getiauforcing(IPD_Control,IAU_Data) if (nfiles.EQ.1) then ! on check to see if we are in the IAU window, no need to update the ! tendencies since they are fixed over the window - if ( IPD_Control%fhour < t1 .or. IPD_Control%fhour >= t2 ) then -! if (is_master()) print *,'no iau forcing',t1,IPD_Control%fhour,t2 + if (IPD_Control%fhour < IPD_Control%iaufhrs(1) .or. IPD_Control%fhour >= IPD_Control%iaufhrs(nfiles)) then IAU_Data%in_interval=.false. else if (IPD_Control%iau_filter_increments) call setiauforcing(IPD_Control,IAU_Data,iau_state%wt) - if (is_master()) print *,'apply iau forcing',t1,IPD_Control%fhour,t2 + if (is_master()) print *,'apply iau forcing',IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles) IAU_Data%in_interval=.true. endif return @@ -352,6 +349,7 @@ subroutine getiauforcing(IPD_Control,IAU_Data) ! if (is_master()) print *,'no iau forcing',IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles) IAU_Data%in_interval=.false. else + if (is_master()) print *,'apply iau forcing',IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles) IAU_Data%in_interval=.true. do k=nfiles,1,-1 if (IPD_Control%iaufhrs(k) > IPD_Control%fhour) then From ddc90d1b51ff3500802ee0a392a6916154446d62 Mon Sep 17 00:00:00 2001 From: Jeffrey S Whitaker Date: Tue, 11 Jan 2022 20:16:43 -0600 Subject: [PATCH 2/5] fix time interval for 3DIAU --- tools/fv_iau_mod.F90 | 25 ++++++++++++++++--------- 1 file changed, 16 insertions(+), 9 deletions(-) diff --git a/tools/fv_iau_mod.F90 b/tools/fv_iau_mod.F90 index c46534894..175782fbe 100644 --- a/tools/fv_iau_mod.F90 +++ b/tools/fv_iau_mod.F90 @@ -253,6 +253,7 @@ subroutine IAU_initialize (IPD_Control, IAU_Data,Init_parm) allocate (iau_state%inc1%tracer_inc(is:ie, js:je, km,ntracers)) iau_state%hr1=IPD_Control%iaufhrs(1) iau_state%wt = 1.0 ! IAU increment filter weights (default 1.0) + iau_state%wt_normfact = 1.0 if (IPD_Control%iau_filter_increments) then ! compute increment filter weights, sum to obtain normalization factor dtp=IPD_control%dtp @@ -305,17 +306,24 @@ subroutine getiauforcing(IPD_Control,IAU_Data) return endif + if (nfiles .eq. 1) then + t1 = IPD_Control%iaufhrs(1)-0.5*IPD_Control%iau_delthrs + t2 = IPD_Control%iaufhrs(1)+0.5*IPD_Control%iau_delthrs + else + t1 = IPD_Control%iaufhrs(1) + t2 = IPD_Control%iaufhrs(nfiles) + endif if (IPD_Control%iau_filter_increments) then ! compute increment filter weight - ! IPD_Control%iaufhrs(1) beginning of window, IPD_Control%iaufhrs(nfiles) end of window + ! t1 isbeginning of window, t2 end of window ! IPD_Control%fhour current time ! in window kstep=-nstep,nstep (2*nstep+1 total) ! time step IPD_control%dtp dtp=IPD_control%dtp nstep = 0.5*IPD_Control%iau_delthrs*3600/dtp ! compute normalized filter weight - kstep = ((IPD_Control%fhour-IPD_Control%iaufhrs(1)) - 0.5*IPD_Control%iau_delthrs)*3600./dtp - if (IPD_Control%fhour >= IPD_Control%iaufhrs(1) .and. IPD_Control%fhour < IPD_Control%iaufhrs(nfiles)) then + kstep = ((IPD_Control%fhour-t1) - 0.5*IPD_Control%iau_delthrs)*3600./dtp + if (IPD_Control%fhour >= t1 .and. IPD_Control%fhour < t2) then sx = acos(-1.)*kstep/nstep wx = acos(-1.)*kstep/(nstep+1) if (kstep .ne. 0) then @@ -324,7 +332,7 @@ subroutine getiauforcing(IPD_Control,IAU_Data) wt = 1. endif iau_state%wt = iau_state%wt_normfact*wt - if (is_master()) print *,'filter wt',kstep,IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles),iau_state%wt/iau_state%wt_normfact + !if (is_master()) print *,'kstep,t1,t,t2,filter wt=',kstep,t1,IPD_Control%fhour,t2,iau_state%wt/iau_state%wt_normfact else iau_state%wt = 0. endif @@ -333,11 +341,10 @@ subroutine getiauforcing(IPD_Control,IAU_Data) if (nfiles.EQ.1) then ! on check to see if we are in the IAU window, no need to update the ! tendencies since they are fixed over the window - if (IPD_Control%fhour < IPD_Control%iaufhrs(1) .or. IPD_Control%fhour >= IPD_Control%iaufhrs(nfiles)) then - IAU_Data%in_interval=.false. + if (IPD_Control%fhour < t1 .or. IPD_Control%fhour >= t2) then else if (IPD_Control%iau_filter_increments) call setiauforcing(IPD_Control,IAU_Data,iau_state%wt) - if (is_master()) print *,'apply iau forcing',IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles) + if (is_master()) print *,'apply iau forcing t1,t,t2,filter wt=',t1,IPD_Control%fhour,t2,iau_state%wt/iau_state%wt_normfact IAU_Data%in_interval=.true. endif return @@ -345,11 +352,11 @@ subroutine getiauforcing(IPD_Control,IAU_Data) if (nfiles > 1) then t2=2 - if (IPD_Control%fhour < IPD_Control%iaufhrs(1) .or. IPD_Control%fhour >= IPD_Control%iaufhrs(nfiles)) then + if (IPD_Control%fhour < t1 .or. IPD_Control%fhour >= t2) then ! if (is_master()) print *,'no iau forcing',IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles) IAU_Data%in_interval=.false. else - if (is_master()) print *,'apply iau forcing',IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles) + if (is_master()) print *,'apply iau forcing t1,t,t2,filter wt=',t1,IPD_Control%fhour,t2,iau_state%wt/iau_state%wt_normfact IAU_Data%in_interval=.true. do k=nfiles,1,-1 if (IPD_Control%iaufhrs(k) > IPD_Control%fhour) then From f93268deff5c78496f5a06d4e9a8f9700c9c23f7 Mon Sep 17 00:00:00 2001 From: Jeffrey S Whitaker Date: Tue, 11 Jan 2022 20:24:39 -0600 Subject: [PATCH 3/5] fix typo in comment --- tools/fv_iau_mod.F90 | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/tools/fv_iau_mod.F90 b/tools/fv_iau_mod.F90 index 175782fbe..eab283100 100644 --- a/tools/fv_iau_mod.F90 +++ b/tools/fv_iau_mod.F90 @@ -315,7 +315,7 @@ subroutine getiauforcing(IPD_Control,IAU_Data) endif if (IPD_Control%iau_filter_increments) then ! compute increment filter weight - ! t1 isbeginning of window, t2 end of window + ! t1 is beginning of window, t2 end of window ! IPD_Control%fhour current time ! in window kstep=-nstep,nstep (2*nstep+1 total) ! time step IPD_control%dtp @@ -341,7 +341,9 @@ subroutine getiauforcing(IPD_Control,IAU_Data) if (nfiles.EQ.1) then ! on check to see if we are in the IAU window, no need to update the ! tendencies since they are fixed over the window - if (IPD_Control%fhour < t1 .or. IPD_Control%fhour >= t2) then + if ( IPD_Control%fhour < t1 .or. IPD_Control%fhour >= t2 ) then +! if (is_master()) print *,'no iau forcing',t1,IPD_Control%fhour,t2 + IAU_Data%in_interval=.false. else if (IPD_Control%iau_filter_increments) call setiauforcing(IPD_Control,IAU_Data,iau_state%wt) if (is_master()) print *,'apply iau forcing t1,t,t2,filter wt=',t1,IPD_Control%fhour,t2,iau_state%wt/iau_state%wt_normfact From dfb0f8065715932f32d931f5188f027af44598bb Mon Sep 17 00:00:00 2001 From: jswhit2 Date: Tue, 1 Feb 2022 03:32:21 +0000 Subject: [PATCH 4/5] fix bug found in review by @MingjingTong-NOAA --- tools/fv_iau_mod.F90 | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/tools/fv_iau_mod.F90 b/tools/fv_iau_mod.F90 index eab283100..9501fa73d 100644 --- a/tools/fv_iau_mod.F90 +++ b/tools/fv_iau_mod.F90 @@ -298,7 +298,7 @@ subroutine getiauforcing(IPD_Control,IAU_Data) implicit none type (IPD_control_type), intent(in) :: IPD_Control type(IAU_external_data_type), intent(inout) :: IAU_Data - real(kind=kind_phys) t1,t2,sx,wx,wt,dtp + real(kind=kind_phys) t1,t2,tnext,sx,wx,wt,dtp integer n,i,j,k,sphum,kstep,nstep IAU_Data%in_interval=.false. @@ -353,7 +353,7 @@ subroutine getiauforcing(IPD_Control,IAU_Data) endif if (nfiles > 1) then - t2=2 + tnext=2 if (IPD_Control%fhour < t1 .or. IPD_Control%fhour >= t2) then ! if (is_master()) print *,'no iau forcing',IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles) IAU_Data%in_interval=.false. @@ -362,16 +362,16 @@ subroutine getiauforcing(IPD_Control,IAU_Data) IAU_Data%in_interval=.true. do k=nfiles,1,-1 if (IPD_Control%iaufhrs(k) > IPD_Control%fhour) then - t2=k + tnext=k endif enddo -! if (is_master()) print *,'t2=',t2 +! if (is_master()) print *,'tnext=',tnext if (IPD_Control%fhour >= iau_state%hr2) then ! need to read in next increment file iau_state%hr1=iau_state%hr2 - iau_state%hr2=IPD_Control%iaufhrs(int(t2)) + iau_state%hr2=IPD_Control%iaufhrs(int(tnext)) iau_state%inc1=iau_state%inc2 - if (is_master()) print *,'reading next increment file',trim(IPD_Control%iau_inc_files(int(t2))) - call read_iau_forcing(IPD_Control,iau_state%inc2,'INPUT/'//trim(IPD_Control%iau_inc_files(int(t2)))) + if (is_master()) print *,'reading next increment file',trim(IPD_Control%iau_inc_files(int(tnext))) + call read_iau_forcing(IPD_Control,iau_state%inc2,'INPUT/'//trim(IPD_Control%iau_inc_files(int(tnext)))) endif call updateiauforcing(IPD_Control,IAU_Data,iau_state%wt) endif From 8ed31c8b2537f37b3d2c38ff5f8b736a092cc9ed Mon Sep 17 00:00:00 2001 From: Jeffrey S Whitaker Date: Tue, 1 Feb 2022 14:27:12 -0600 Subject: [PATCH 5/5] change tnext to integer variable itnext --- tools/fv_iau_mod.F90 | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/tools/fv_iau_mod.F90 b/tools/fv_iau_mod.F90 index 9501fa73d..703640a4f 100644 --- a/tools/fv_iau_mod.F90 +++ b/tools/fv_iau_mod.F90 @@ -298,8 +298,8 @@ subroutine getiauforcing(IPD_Control,IAU_Data) implicit none type (IPD_control_type), intent(in) :: IPD_Control type(IAU_external_data_type), intent(inout) :: IAU_Data - real(kind=kind_phys) t1,t2,tnext,sx,wx,wt,dtp - integer n,i,j,k,sphum,kstep,nstep + real(kind=kind_phys) t1,t2,sx,wx,wt,dtp + integer n,i,j,k,sphum,kstep,nstep,itnext IAU_Data%in_interval=.false. if (nfiles.LE.0) then @@ -353,7 +353,7 @@ subroutine getiauforcing(IPD_Control,IAU_Data) endif if (nfiles > 1) then - tnext=2 + itnext=2 if (IPD_Control%fhour < t1 .or. IPD_Control%fhour >= t2) then ! if (is_master()) print *,'no iau forcing',IPD_Control%iaufhrs(1),IPD_Control%fhour,IPD_Control%iaufhrs(nfiles) IAU_Data%in_interval=.false. @@ -362,16 +362,16 @@ subroutine getiauforcing(IPD_Control,IAU_Data) IAU_Data%in_interval=.true. do k=nfiles,1,-1 if (IPD_Control%iaufhrs(k) > IPD_Control%fhour) then - tnext=k + itnext=k endif enddo -! if (is_master()) print *,'tnext=',tnext +! if (is_master()) print *,'itnext=',itnext if (IPD_Control%fhour >= iau_state%hr2) then ! need to read in next increment file iau_state%hr1=iau_state%hr2 - iau_state%hr2=IPD_Control%iaufhrs(int(tnext)) + iau_state%hr2=IPD_Control%iaufhrs(itnext) iau_state%inc1=iau_state%inc2 - if (is_master()) print *,'reading next increment file',trim(IPD_Control%iau_inc_files(int(tnext))) - call read_iau_forcing(IPD_Control,iau_state%inc2,'INPUT/'//trim(IPD_Control%iau_inc_files(int(tnext)))) + if (is_master()) print *,'reading next increment file',trim(IPD_Control%iau_inc_files(itnext)) + call read_iau_forcing(IPD_Control,iau_state%inc2,'INPUT/'//trim(IPD_Control%iau_inc_files(itnext))) endif call updateiauforcing(IPD_Control,IAU_Data,iau_state%wt) endif