Skip to content

Commit

Permalink
Add incremental analysis update capability (#2151)
Browse files Browse the repository at this point in the history 
TYPE: new feature

KEYWORDS: data assimilation, incremental analysis update

SOURCE: Min Chen of IUM/CMA and internal

DESCRIPTION OF CHANGES:
Add incremental analysis update capability. In the DA code, code is
added to write out analysis increments for all variables in WRF netCDF
format using auxiliary history output stream #5. In the model, analysis
increments are divided by the number of time steps in a specified time
window and added to the model similar to physics tendencies. The input
stream for the model is 15. The capability is turned on by adding iau =
1 and iau_time_window_sec in &time_control. For example:

```
 auxinput15_inname                   = "wrfiauinp_d01"
 io_form_auxinput15                  = 2
 auxinput15_interval                 = 360,
 iau                                 = 1
 iau_time_window_sec                 = 3600.
```
LIST OF MODIFIED FILES: 
M       Registry/Registry.EM_COMMON.var
M       Registry/Registry.wrfvar
M       Registry/registry.em_shared_collection
A       Registry/registry.iau
M       dyn_em/module_em.F
M       dyn_em/module_first_rk_step_part2.F
M       phys/module_physics_addtendc.F
M       var/build/da.make
M       var/build/depend.txt
M       var/da/da_main/da_solve.inc
M       var/da/da_main/da_wrfvar_io.f90
M       var/da/da_main/da_wrfvar_top.f90
A       var/da/da_main/da_write_anaincrements.inc
M       var/da/da_transfer_model/da_transfer_xatowrf.inc

TESTS CONDUCTED: 
The Jenkins tests have passed.

RELEASE NOTE: This PR adds an incremental analysis update capability. In
the DA code, code is added to write out analysis increments for all
variables in WRF netCDF format using auxiliary history output stream #5.
In the model, analysis increments are divided by the number of time
steps in a specified time window and added to the model similar to
physics tendencies. The input stream for the model is 15. The capability
is turned on by adding iau = 1 and iau_time_window_sec in &time_control.
The way the increments are added to the model is similar to what
described by the paper by Chen et al.
(https://doi-org.cuucar.idm.oclc.org/10.1175/WAF-D-22-0127.1).
  • Loading branch information
@weiwangncar
weiwangncar authored Jan 29, 2025
1 parent 70855a7 commit 6741f01
Show file tree
Hide file tree
Showing 14 changed files with 461 additions and 8 deletions.
14 changes: 14 additions & 0 deletions Registry/Registry.EM_COMMON.var
View file
Original file line number Diff line number Diff line change
Expand Up @@ -347,6 +347,20 @@ state real QFX ij misc 1 - irh "Q
state real REGIME ij misc 1 - irh "REGIME" "FLAGS: 1=Night/Stable, 2=Mechanical Turbulent, 3=Forced Conv, 4=Free Conv" ""
state integer KPBL ij misc 1 - irh "KPBL" "LEVEL OF PBL TOP" ""

# Increment Output
state real u_iau ijk dyn_em 1 X ih5 "U_IAU" "x-wind component inc" "m s-1"
state real v_iau ijk dyn_em 1 Y ih5 "V_IAU" "y-wind component inc" "m s-1"
state real t_iau ijk dyn_em 1 - ih5 "T_IAU" "potential temp inc" "K"
state real w_iau ijk dyn_em 1 - ih5 "W_IAU" "z-wind component inc" "m s-1"
state real qv_iau ijk dyn_em 1 - ih5 "QV_IAU" "water water mixing ratio inc" "kg kg-1"
state real qc_iau ijk dyn_em 1 - ih5 "QC_IAU" "cloud water mixing ratio inc" "kg kg-1"
state real qr_iau ijk dyn_em 1 - ih5 "QR_IAU" "rain water mixing ratio inc" "kg kg-1"
state real qi_iau ijk dyn_em 1 - ih5 "QI_IAU" "ice water mixing ratio inc" "kg kg-1"
state real qs_iau ijk dyn_em 1 - ih5 "QS_IAU" "snow water mixing ratio inc" "kg kg-1"
state real qg_iau ijk dyn_em 1 - ih5 "QG_IAU" "graupel mixing ratio inc" "kg kg-1"
state real ph_iau ijk dyn_em 1 - ih5 "PH_IAU" "perturbation geopotential inc" "m2 s-2"
state real mu_iau ij dyn_em 1 - ih5 "MU_IAU" "dry air mass inc" "pa"

#
#---------------------------------------------------------------------------------------------------------------------------------------
#
Expand Down
2 changes: 1 addition & 1 deletion Registry/Registry.wrfvar
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@ state real - ijkft g_scalar 1 - - -

state real landmask ij misc 1 - i012rhd=(interp_fcnm)u=(copy_fcnm) "LANDMASK" "LAND MASK (1 FOR LAND, 0 FOR WATER)" ""

state real SST ij misc 1 - i01245rh05d=(interp_mask_water_field:lu_index,iswater) "SST" "SEA SURFACE TEMPERATURE" "K"
state real SST ij misc 1 - i01245rd=(interp_mask_water_field:lu_index,iswater) "SST" "SEA SURFACE TEMPERATURE" "K"

# Registry entries that are exclusive to Registry.EM

Expand Down
1 change: 1 addition & 0 deletions Registry/registry.em_shared_collection
View file
Original file line number Diff line number Diff line change
Expand Up @@ -28,4 +28,5 @@ include registry.new3d_wif
include registry.trad_fields
include registry.solar_fields
include registry.diags
include registry.iau
include registry.CMAQ
38 changes: 38 additions & 0 deletions Registry/registry.iau
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,38 @@
# IAU variables

state character iau_time - - - - i{15}r "TIME_IAU" " " " "
state real mu_iau ij misc 1 - i{15}r "MU_IAU" "mu analysis increments array" " "
state real u_iau ikj misc 1 - i{15}r "U_IAU" "u analysis increments array" " "
state real v_iau ikj misc 1 - i{15}r "V_IAU" "v analysis increments array" " "
state real w_iau ikj misc 1 - i{15}r "W_IAU" "w analysis increments array" " "
state real p_iau ikj misc 1 - i{15}r "P_IAU" "p analysis increments array" " "
state real t_iau ikj misc 1 - i{15}r "T_IAU" "t analysis increments array" " "
state real ph_iau ikj misc 1 - i{15}r "PH_IAU" "ph analysis increments array" " "
state real qv_iau ikj misc 1 - i{15}r "QV_IAU" "qv analysis increments array" " "
state real qr_iau ikj misc 1 - i{15}r "QR_IAU" "qr analysis increments array" " "
state real qc_iau ikj misc 1 - i{15}r "QC_IAU" "qc analysis increments array" " "
state real qs_iau ikj misc 1 - i{15}r "QS_IAU" "qs analysis increments array" " "
state real qi_iau ikj misc 1 - i{15}r "QI_IAU" "qice analysis increments array" " "
state real qg_iau ikj misc 1 - i{15}r "QG_IAU" "qgraupel analysis increments array" " "

state real RUIAUTEN ikj misc 1 X r "RUIAUTEN" "X WIND TENDENCY DUE TO IAU" "m s-2"
state real RVIAUTEN ikj misc 1 Y r "RVIAUTEN" "Y WIND TENDENCY DUE TO IAU" "m s-2"
state real RTHIAUTEN ikj misc 1 - r "RTHIAUTEN" "THETA TENDENCY DUE TO IAU" "K s-1"
state real RPHIAUTEN ikj misc 1 - r "RPHIAUTEN" "GEOPOTENTIAL TENDENCY DUE TO IAU" "m2 s-3"
state real RQVIAUTEN ikj misc 1 - r "RQVIAUTEN" "Q_V TENDENCY DUE TO IAU" "kg kg-1 s-1"
state real RQCIAUTEN ikj misc 1 - r "RQCIAUTEN" "Q_C TENDENCY DUE TO IAU" "kg kg-1 s-1"
state real RQRIAUTEN ikj misc 1 - r "RQRIAUTEN" "Q_R TENDENCY DUE TO IAU" "kg kg-1 s-1"
state real RQIIAUTEN ikj misc 1 - r "RQIIAUTEN" "Q_I TENDENCY DUE TO IAU" "kg kg-1 s-1"
state real RQSIAUTEN ikj misc 1 - r "RQSIAUTEN" "Q_S TENDENCY DUE TO IAU" "kg kg-1 s-1"
state real RQGIAUTEN ikj misc 1 - r "RQGIAUTEN" "Q_G TENDENCY DUE TO IAU" "kg kg-1 s-1"
state real RMUIAUTEN ij misc 1 - r "RMUIAUTEN" "MU TENDENCY DUE TO IAU" "Pa s-1"

# IAU namelist options

rconfig integer iau namelist,time_control max_domains 0 irh "analysis increments read" "0/1 ACTIVATE FOR ANALYSIS INCREMENTS UPDATES" ""
rconfig real iau_time_window_sec namelist,time_control max_domains 3600. irh "iau_time_window_sec" "TIME WINDOW OF INCREMENTS ANALYSIS UPDATES" "SECONDS"

# IAU packages

package noiau iau==0 - -
package iau iau==1 - state:u_iau,v_iau,w_iau,p_iau,t_iau,ph_iau,qqv_iau,qqr_iau,qqc_iau,qqs_iau,qqi_iau,qqg_iau,ruiauten,rviauten,rthiauten,rqviauten,rqciauten,rqriauten,rqiiauten,rqsiauten,rqgiauten
144 changes: 144 additions & 0 deletions dyn_em/module_em.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2089,6 +2089,14 @@ SUBROUTINE calculate_phy_tend (config_flags,c1,c2, &
RQISHTEN,RQSSHTEN,RQGSHTEN, &
RUNDGDTEN,RVNDGDTEN,RTHNDGDTEN,RQVNDGDTEN, &
RMUNDGDTEN, &
ruiauten,rviauten,rthiauten, &
rqviauten,rqciauten,rqriauten, &
rqiiauten,rqsiauten,rqgiauten, &
rmuiauten,rphiauten, &
u_iau,v_iau,t_iau,ph_iau, &
mu_iau,qv_iau,qc_iau,qr_iau, &
qi_iau,qs_iau,qg_iau, &
itimestep, &
scalar, scalar_tend, num_scalar, &
tracer, tracer_tend, num_tracer, &
ids,ide, jds,jde, kds,kde, &
Expand Down Expand Up @@ -2163,16 +2171,52 @@ SUBROUTINE calculate_phy_tend (config_flags,c1,c2, &
REAL, DIMENSION( ims:ime, jms:jme ) , &
INTENT(INOUT) :: RMUNDGDTEN

! iau
REAL, DIMENSION( ims:ime , kms:kme , jms:jme ), &
INTENT(INOUT) :: RUIAUTEN, &
RVIAUTEN, &
RTHIAUTEN, &
RPHIAUTEN, &
RQVIAUTEN, &
RQCIAUTEN, &
RQRIAUTEN, &
RQIIAUTEN, &
RQSIAUTEN, &
RQGIAUTEN

REAL, DIMENSION( ims:ime, jms:jme ), &
INTENT(INOUT) :: RMUIAUTEN

REAL, DIMENSION( ims:ime , kms:kme , jms:jme ), &
INTENT(IN ) :: &
u_iau, &
v_iau, &
t_iau, &
ph_iau, &
qv_iau, &
qc_iau, &
qr_iau, &
qs_iau, &
qi_iau, &
qg_iau

REAL, DIMENSION( ims:ime, jms:jme ), &
INTENT(IN ) :: mu_iau

! 4d arrays

REAL ,DIMENSION(ims:ime,kms:kme,jms:jme,num_tracer),INTENT(INOUT) :: tracer
REAL ,DIMENSION(ims:ime,kms:kme,jms:jme,num_tracer),INTENT(INOUT) :: tracer_tend
REAL ,DIMENSION(ims:ime,kms:kme,jms:jme,num_scalar),INTENT(INOUT) :: scalar
REAL ,DIMENSION(ims:ime,kms:kme,jms:jme,num_scalar),INTENT(INOUT) :: scalar_tend

INTEGER :: itimestep
INTEGER :: i,k,j, im
INTEGER :: itf,ktf,jtf,itsu,jtsv

! local
real :: wgt_iau

!-----------------------------------------------------------------------

!<DESCRIPTION>
Expand Down Expand Up @@ -2448,6 +2492,88 @@ SUBROUTINE calculate_phy_tend (config_flags,c1,c2, &
ENDDO
ENDDO

! IAU (Incremental Analysis Updates)

IF (config_flags%iau .gt. 0) THEN

wgt_iau = iau_coef(config_flags%iau_time_window_sec,config_flags%dt,itimestep)

DO J=jts,jtf
DO K=kts,ktf
DO I=its,itf
RUIAUTEN(I,K,J) =(c1(k)*mut(I,J)+c2(k))*U_IAU(I,K,J)*WGT_IAU
RVIAUTEN(I,K,J) =(c1(k)*mut(I,J)+c2(k))*V_IAU(I,K,J)*WGT_IAU
RTHIAUTEN(I,K,J)=(c1(k)*mut(I,J)+c2(k))*T_IAU(I,K,J)*WGT_IAU
RQVIAUTEN(I,K,J)=(c1(k)*mut(I,J)+c2(k))*QV_IAU(I,K,J)*WGT_IAU
ENDDO
ENDDO
ENDDO

DO J=jts,jtf
DO K=kts,ktf+1
DO I=its,itf
RPHIAUTEN(I,K,J)=PH_IAU(I,K,J)*WGT_IAU
ENDDO
ENDDO
ENDDO

DO J=jts,jtf
DO I=its,itf
RMUIAUTEN(I,J) =MU_IAU(I,J)*WGT_IAU
ENDDO
ENDDO

IF (P_QC .ge. PARAM_FIRST_SCALAR)THEN
DO J=jts,jtf
DO K=kts,ktf
DO I=its,itf
RQCIAUTEN(I,K,J)=(c1(k)*mut(I,J)+c2(k))*QC_IAU(I,K,J)*WGT_IAU
ENDDO
ENDDO
ENDDO
ENDIF

IF (P_QR .ge. PARAM_FIRST_SCALAR)THEN
DO J=jts,jtf
DO K=kts,ktf
DO I=its,itf
RQRIAUTEN(I,K,J)=(c1(k)*mut(I,J)+c2(k))*QR_IAU(I,K,J)*WGT_IAU
ENDDO
ENDDO
ENDDO
ENDIF

IF (P_QI .ge. PARAM_FIRST_SCALAR)THEN
DO J=jts,jtf
DO K=kts,ktf
DO I=its,itf
RQIIAUTEN(I,K,J)=(c1(k)*mut(I,J)+c2(k))*QI_IAU(I,K,J)*WGT_IAU
ENDDO
ENDDO
ENDDO
ENDIF

IF(P_QS .ge. PARAM_FIRST_SCALAR)THEN
DO J=jts,jtf
DO K=kts,ktf
DO I=its,itf
RQSIAUTEN(I,K,J)=(c1(k)*mut(I,J)+c2(k))*QS_IAU(I,K,J)*WGT_IAU
ENDDO
ENDDO
ENDDO
ENDIF

IF(P_QG .ge. PARAM_FIRST_SCALAR)THEN
DO J=jts,jtf
DO K=kts,ktf
DO I=its,itf
RQGIAUTEN(I,K,J)=(c1(k)*mut(I,J)+c2(k))*QG_IAU(I,K,J)*WGT_IAU
ENDDO
ENDDO
ENDDO
ENDIF

ENDIF

END SUBROUTINE calculate_phy_tend

Expand Down Expand Up @@ -3018,4 +3144,22 @@ subroutine trajmapproj (grid,config_flags,ts_proj)

end subroutine trajmapproj

FUNCTION iau_coef (iau_time_window_sec, dt, itimestep) result(wgt_iau)

implicit none

! This function returns the coeficient for IAU

INTEGER :: itimestep
INTEGER :: nsteps_iau
REAL :: wgt_iau, iau_time_window_sec, dt

wgt_iau = 0.
nsteps_iau = nint(iau_time_window_sec / dt)
if (itimestep <= nsteps_iau) then
wgt_iau = 1.0/iau_time_window_sec
endif

END FUNCTION iau_coef

END MODULE module_em
12 changes: 12 additions & 0 deletions dyn_em/module_first_rk_step_part2.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -402,6 +402,14 @@ SUBROUTINE first_rk_step_part2 ( grid , config_flags &
grid%rqishten,grid%rqsshten,grid%rqgshten, &
grid%RUNDGDTEN,grid%RVNDGDTEN,grid%RTHNDGDTEN,grid%RQVNDGDTEN, &
grid%RMUNDGDTEN, &
grid%ruiauten,grid%rviauten,grid%rthiauten, &
grid%rqviauten,grid%rqciauten,grid%rqriauten, &
grid%rqiiauten,grid%rqsiauten,grid%rqgiauten, &
grid%rmuiauten,grid%rphiauten, &
grid%u_iau,grid%v_iau,grid%t_iau,grid%ph_iau, &
grid%mu_iau,grid%qv_iau,grid%qc_iau,grid%qr_iau, &
grid%qs_iau,grid%qi_iau,grid%qg_iau, &
grid%itimestep, &
scalar, scalar_tend, num_scalar, &
tracer, tracer_tend, num_tracer, &
ids,ide, jds,jde, kds,kde, &
Expand Down Expand Up @@ -807,6 +815,10 @@ SUBROUTINE first_rk_step_part2 ( grid , config_flags &
grid%RVNDGDTEN,grid%RTHNDGDTEN,grid%RPHNDGDTEN, &
grid%RQVNDGDTEN,grid%RMUNDGDTEN, &
grid%rthfrten,grid%rqvfrten, & ! fire
grid%ruiauten,grid%rviauten,grid%rthiauten, &
grid%rqviauten,grid%rqciauten,grid%rqriauten, &
grid%rqiiauten,grid%rqsiauten,grid%rqgiauten, &
grid%rphiauten,grid%rmuiauten, &
num_moist,num_scalar,config_flags,rk_step, &
grid%adv_moist_cond, &
ids, ide, jds, jde, kds, kde, &
Expand Down
Loading

0 comments on commit 6741f01

Please sign in to comment.