| Description |
This update includes several improvements to the 4D-Var data assimilation algorithms:
- Modified routine allocate_boundary in module mod_boundary.F to allow the allocation of the boundary arrays using both the nonlinear and adjoint ad_LBC and LBC switches. For example,
#if defined ADJOINT || defined TANGENT || defined TL_IOMS
IF ( LBC(iwest,isFsur,ng)%acquire.or. &
& ad_LBC(iwest,isFsur,ng)%acquire) THEN
#else
IF (LBC(iwest,isFsur,ng)%acquire) THEN
#endif
allocate ( BOUNDARY(ng) % zeta_west(LBj:UBj) )
...
END IF
It allows processing of the boundary data for the adjoint, tangent linear or representer models even if it is not needed by the nonlinear model when using different boundary conditions between models. Many thanks to Julia Levin for bringing this to my attention.
- Renumbered the model and message identifiers passed to routine get_state. Added new entries to differentiate the processing of various state tasks:
character (len=48), dimension(17) :: StateMsg = &
& (/'state initial conditions, ', & !01
& 'previous state initial conditions, ', & !02
& 'previous adjoint state solution, ', & !03
& 'latest adjoint state solution, ', & !04
& 'surface forcing and or OBC increments, ', & !05
& 'tangent linear model error forcing, ', & !06
& 'impulse forcing, ', & !07
& 'v-space increments, ', & !08
& 'background state, ', & !09
& 'IC correlation standard deviation, ', & !10
& 'model error correlation standard deviation, ', & !11
& 'OBC correlation standard deviation, ', & !12
& 'surface forcing correlation standard deviation, ', & !13
& 'IC normalization factors, ', & !14
& 'model error normalization factors, ', & !15
& 'OBC normalization factors, ', & !16
& 'surface forcing normalization factors, '/) !17
For Example, during the processing of the correlations standard deviations we now have:
! Initial conditions standard deviation. They are loaded in Tindex=1
! of the e_var(...,Tindex) state variables.
!
STDrec=1
Tindex=1
DO ng=1,Ngrids
CALL get_state (ng, 10, 10, STD(1,ng)%name, STDrec, Tindex)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
END DO
!
! Model error standard deviation. They are loaded in Tindex=2
! of the e_var(...,Tindex) state variables.
!
STDrec=1
Tindex=2
IF (NSA.eq.2) THEN
DO ng=1,Ngrids
CALL get_state (ng, 11, 11, STD(2,ng)%name, STDrec, Tindex)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
END DO
END IF
#ifdef ADJUST_BOUNDARY
!
! Open boundary conditions standard deviation.
!
STDrec=1
Tindex=1
DO ng=1,Ngrids
CALL get_state (ng, 12, 12, STD(3,ng)%name, STDrec, Tindex)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
END DO
#endif
#if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
!
! Surface forcing standard deviation.
!
STDrec=1
Tindex=1
DO ng=1,Ngrids
CALL get_state (ng, 13, 13, STD(4,ng)%name, STDrec, Tindex)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
END DO
#endif
Similarly, to read and process the error covariance normalization factors we have:
NRMrec=1
CALL get_state (ng, 14, 14, NRM(1,ng)%name, NRMrec, 1)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
IF (NSA.eq.2) THEN
CALL get_state (ng, 15, 15, NRM(2,ng)%name, NRMrec, 2)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
END IF
#ifdef ADJUST_BOUNDARY
CALL get_state (ng, 16, 16, NRM(3,ng)%name, NRMrec, 1)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
#endif
#if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
CALL get_state (ng, 17, 17, NRM(4,ng)%name, NRMrec, 1)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
#endif
END IF
WARNING: All the 4D-Var drivers and associated algorithms were modified because of the renumbering of the tasks. The new indices allow to turn off the time messages when processing both standard deviation and normalization factors. It was confusing for some users. In get_state.F, we just need to overwrite the string t_code as follows:
IF (Master) THEN
CALL time_string (INPtime, t_code)
IF ((10.le.model).and.(model.le.17)) THEN
t_code=' ' ! time is meaningless for these fields
END IF
WRITE (Tstring,'(f15.4)') tdays(ng)
IF (ERend.gt.ERstr) THEN
WRITE (stdout,40) string, 'Reading '//TRIM(StateMsg(msg)), &
& t_code, ng, Nrun, TRIM(ADJUSTL(Tstring)), &
& ncname(lstr:lend), InpRec, Tindex
ELSE
WRITE (stdout,50) string, 'Reading '//TRIM(StateMsg(msg)), &
& t_code, ng, TRIM(ADJUSTL(Tstring)), &
& ncname(lstr:lend), InpRec, Tindex
END IF
END IF
- Updated the observation impact and sensitivity algorithms to allow the saving of the analysis data for each inner loop. A new C-preprocessing option IMPACT_INNER is introduced to achieve the saving of analysis quantities into the output MOD NetCDF file. For Example, in def_mod.F we have:
!
! Define observation impact due to initial condition increments.
!
# ifdef IMPACT_INNER
Vinfo( 1)='ObsImpact_IC'
Vinfo( 2)='observation impact due to initial conditions'
vardim(1)=datumDim
vardim(2)=NinnerDim
Vinfo(24)='_FillValue'
Aval(6)=spval
status=def_var(ng, iNLM, DAV(ng)%ncid, varid, nf90_double, &
& 2, vardim, Aval, Vinfo, ncname)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
# else
Vinfo( 1)='ObsImpact_IC'
Vinfo( 2)='observation impact due to initial conditions'
vardim(1)=datumDim
Vinfo(24)='_FillValue'
Aval(6)=spval
status=def_var(ng, iNLM, DAV(ng)%ncid, varid, nf90_double, &
& 1, vardim, Aval, Vinfo, ncname)
IF (FoundError(exit_flag, NoError, __LINE__, &
& __FILE__)) RETURN
# endif
Notice that the output variable ObsImpact_IC can be either one- or two-dimension array.
- Added a new adjoint routine ad_wvelocity.F to allow observation impacts or sensitivities using vertical velocity functionals. Many thanks to Andy Moore for coding the adjoint of ROMS diagnostic routine wvelocity.F. Also, for modifying the algorithms to allow impacts and sensitivities functionals that use vertical velocity.
|
| Description |
Added a new C-processing option WIND_MINUS_CURRENT to compute the effective wind forcing by subtracting surface ocean currents to input wind components in bulk_flux.F:
#ifdef WIND_MINUS_CURRENT
!
! Modify near-surface (2m or 10m) effective winds by subtracting the
! ocean surface current (J Wilkin). See:
!
! Bye, J.A.T. and J.-O. Wolff, 1999: Atmosphere-ocean momentum exchange
! in general circulation models. J. Phys. Oceanogr. 29, 671-691.
!
DO j=Jstr-1,JendR
DO i=Istr-1,IendR
Uwind(i,j)=Uwind(i,j)- &
& 0.5_r8*(u(i,j,N(ng),nrhs)+u(i+1,j,N(ng),nrhs))
Vwind(i,j)=Vwind(i,j)- &
& 0.5_r8*(v(i,j,N(ng),nrhs)+v(i,j+1,N(ng),nrhs))
END DO
END DO
#endif
Many thanks to John Wilkin for coding and evaluating this option.
References:
Bye, J.A.T. and J.-O. Wolff, 1999: Atmosphere-ocean momentum exchange in general circulation models, J. Phys. Oceanogr., 29, 671-691.
Sandery, P.A., G.B. Brassington, A. Craig, and T. Pugh, 2010: Impacts of Ocean-Atmosphere Coupling on the Tropical Cyclone Intensity Change and Ocean Prediction in the Australia Region, Mon. Wea. Rev., 138, 2074-2091, doi: 10.1175/2010MWR3101.1
|
| Description |
Minor update to several files:
- Updated routines calendar and ref_clock of module dateclock.F for the computation of the 360DAY calendar (time_ref=-1). It has a year length of 360 days and every month has 30 days:
In this calendar, the time in days is simply:
TimeInDays = year * 360 + (month - 1) * 30 + (day - 1)
And its inverse:
year = INT(TimeInDays / 360)
yday = INT((TimeInDays - year * 360) + 1)
month = INT(((yday - 1) / 30) + 1)
day = MOD(yday - 1, 30) + 1
In the above equations, the origin (TimeInDays=0) corresponds to 01-Jan-000. However, historically ROMS assumed that TimeInDays=1 corresponded to 01-Jan-0000 instead. So there is a one-day shift when compared with third-party software (like ESMF). The equations can be manipulated to give either origin, but it is confusing. Since the above equations are simpler, effective immediately the origin (TimeInDays=0) occurs on 01-Jan-0000. Therefore, to guarantee compatibility with previous ROMS solutions and datasets with this climatological calendar, the reference date is changed to
'time-units since 0000-12-30 00:00:00'
to fix the one-day shift.
- Updated routines datenum and datevec of module dateclock.F to be consitent with the Matlab version of these functions. For simplicity, the equations coded in ROMS assumed that the origin (datenum=0) occurs on 01-Mar-0000 while in Matlab occurs on 00-00-0000:
Matlab: datenum(0000,00,00)=0 reference date
datenum(0000,01,01)=1
ROMS: datenum(0000,03,01)=0 refecence date: Mar 1, 0000 !
datenum(0000,01,01)=-59
To avoid confusion, an offset of 61 days is added to match Matlab datenum function. The difference between 0000-00-00 00:00:00 and 0000-03-01 00:00:00 is 61 days.
- Added a new option -p ot the build.sh and build.bash compiling scripts to print any macro of the compilation:
# Usage: :::
# :::
# ./build.sh [options] :::
# :::
# Options: :::
# :::
# -j [N] Compile in parallel using N CPUs :::
# omit argument for all available CPUs :::
# :::
# -p macro Prints any Makefile macro value. For example, :::
# :::
# build.sh -p FFLAGS :::
# :::
# -noclean Do not clean already compiled objects :::
- Added an extra argument Mgrids to function load_s2d of file inp_par.F to the number of nested grids in a application to allow coupling to other DATA model in the future:
FUNCTION load_s2d (Nval, Fname, Fdim, line, label, ifile, igrid, &
& Mgrids, Nfiles, Ncount, idim, S)
Notice that in read_phypar.F now we have:
CASE ('FRCNAME')
label='FRC - forcing fields'
Npts=load_s2d(Nval, Cval, Cdim, line, label, ifile, &
& igrid, Ngrids, nFfiles, Ncount, max_Ffiles, &
& FRC)
|
