ROMS
Loading...
Searching...
No Matches
mod_fourdvar Module Reference

Data Types

type  t_fourdvar
 

Functions/Subroutines

subroutine, public allocate_fourdvar
 
subroutine, public deallocate_fourdvar
 
subroutine, public initialize_fourdvar
 

Variables

type(t_fourdvar), dimension(:), allocatable fourdvar
 
integer, dimension(:), allocatable obstype
 
integer, dimension(:), allocatable obsstate2type
 
integer, dimension(:), allocatable obstype2state
 
integer, dimension(:), allocatable obsprov
 
real(r8), dimension(:), allocatable obsangler
 
real(r8), dimension(:), allocatable obserr
 
real(r8), dimension(:), allocatable obsmeta
 
real(r8), dimension(:), allocatable obsscale
 
real(r8), dimension(:), allocatable obsvetting
 
real(r8), dimension(:), allocatable obsval
 
real(dp), dimension(:), allocatable tobs
 
real(r8), dimension(:), allocatable xobs
 
real(r8), dimension(:), allocatable yobs
 
real(r8), dimension(:), allocatable zobs
 
real(r8), dimension(:), allocatable admodval
 
real(r8), dimension(:), allocatable nlmodval
 
real(r8), dimension(:), allocatable misfit
 
real(r8), dimension(:), allocatable unvetted
 
real(r8), dimension(:), allocatable uradial
 
real(r8), dimension(:), allocatable vradial
 
real(r8), dimension(:,:), allocatable admodval_s
 
real(r8), dimension(:,:,:), allocatable harnoldi
 
real(r8), dimension(:,:), allocatable gmze
 
real(r8), dimension(:), allocatable cg_beta0
 
real(r8), dimension(:,:), allocatable jobs
 
real(r8), dimension(:,:,:), allocatable zcglwk
 
real(r8), dimension(:,:,:), allocatable vcglwk
 
real(r8), dimension(:), allocatable jb0
 
real(r8), dimension(:,:,:), allocatable vcglev
 
real(r8), dimension(:,:), allocatable zgrad0
 
real(r8), dimension(:), allocatable vgrad0
 
real(r8), dimension(:), allocatable vgrad0s
 
real(r8), dimension(:), allocatable gdgw
 
real(r8), dimension(:), allocatable vgnorm
 
real(r8), dimension(:), allocatable cg_innov
 
real(r8), dimension(:,:), allocatable cg_dla
 
real(r8), dimension(:,:), allocatable ad_zcglwk
 
real(r8), dimension(:), allocatable ad_zgrad0
 
real(r8), dimension(:,:), allocatable ad_vcglwk
 
real(r8), dimension(:), allocatable ad_vgrad0
 
real(r8), dimension(:), allocatable ad_vgrad0s
 
real(r8), dimension(:), allocatable ad_jb0
 
real(r8), dimension(:), allocatable ad_obserr
 
real(r8), dimension(:), allocatable ad_obsscale
 
real(r8), dimension(:), allocatable ad_cg_innov
 
real(r8), dimension(:), allocatable ad_cg_pxsave
 
real(r8), dimension(:), allocatable ad_cg_pxout
 
real(r8), dimension(:,:), allocatable ad_obsval
 
real(r8), dimension(:,:), allocatable tl_zcglwk
 
real(r8), dimension(:), allocatable tl_zgrad0
 
real(r8), dimension(:,:), allocatable tl_vcglwk
 
real(r8), dimension(:), allocatable tl_vgrad0
 
real(r8), dimension(:), allocatable tl_vgrad0s
 
real(r8), dimension(:), allocatable tl_jb0
 
real(r8), dimension(:), allocatable tl_obserr
 
real(r8), dimension(:), allocatable tl_obsscale
 
real(r8), dimension(:), allocatable tl_cg_innov
 
real(r8), dimension(:), allocatable tl_cg_pxsave
 
real(r8), dimension(:), allocatable tl_cg_pxout
 
real(r8), dimension(:), allocatable tl_obsval
 
logical, dimension(:), allocatable load_zobs
 
logical, dimension(:), allocatable wrote_zobs
 
integer mobs
 
integer, dimension(:), allocatable nstatevar
 
integer, dimension(:), allocatable nobsvar
 
character(len=40), dimension(:), allocatable obsname
 
integer nextraobs = 0
 
integer, dimension(:), allocatable extraindex
 
character(len=40), dimension(:), allocatable extraname
 
integer, parameter nweights = 8
 
integer, dimension(nweights), parameter ioffset = (/ 0, 1, 0, 1, 0, 1, 0, 1 /)
 
integer, dimension(nweights), parameter joffset = (/ 0, 0, 1, 1, 0, 0, 1, 1 /)
 
integer, dimension(nweights), parameter koffset = (/ 0, 0, 0, 0, 1, 1, 1, 1 /)
 
integer, dimension(:), allocatable ndatum
 
integer, dimension(:), allocatable nsurvey
 
integer, dimension(:), allocatable obssurvey
 
integer, dimension(:), allocatable nobs
 
integer, dimension(:), allocatable nstrobs
 
integer, dimension(:), allocatable nendobs
 
integer, dimension(:), allocatable nmethod
 
integer, dimension(:), allocatable rscheme
 
integer nrandom = 1000
 
integer nposti
 
integer nimpact
 
integer nvct
 
real(r8), dimension(:), allocatable optimality
 
logical, dimension(:), allocatable processobs
 
logical, dimension(:), allocatable lsadd
 
logical, dimension(:), allocatable wrtnlmod
 
logical, dimension(:), allocatable wrtobsscale
 
logical, dimension(:), allocatable wrtrpmod
 
logical, dimension(:), allocatable wrttlmod
 
logical, dimension(:), allocatable wrtmisfit
 
logical, dimension(:), allocatable wrtimpact_tot
 
logical, dimension(:), allocatable wrtimpact_ic
 
logical, dimension(:), allocatable wrtimpact_fc
 
logical, dimension(:), allocatable wrtimpact_bc
 
logical, dimension(:), allocatable haveadmod
 
logical, dimension(:), allocatable havenlmod
 
logical, dimension(:), allocatable havetlmod
 
logical, dimension(:), allocatable haveobsmeta
 
real(r8), dimension(:), allocatable bgerr
 
real(r8), dimension(:), allocatable nlincrement
 
real(r8), dimension(:), allocatable innovation
 
real(r8), dimension(:), allocatable residual
 
real(r8), dimension(:), allocatable bgthresh
 
real(r8), dimension(:), allocatable jact_s
 
integer, dimension(:,:), allocatable nhsteps
 
integer, dimension(:,:), allocatable nvsteps
 
integer, dimension(:,:), allocatable nhstepsb
 
integer, dimension(:,:), allocatable nvstepsb
 
real(r8), dimension(:,:), allocatable dtsizeh
 
real(r8), dimension(:,:), allocatable dtsizev
 
real(r8), dimension(:,:), allocatable dtsizehb
 
real(r8), dimension(:,:), allocatable dtsizevb
 
real(r8), dimension(:), allocatable khmin
 
real(r8), dimension(:), allocatable khmax
 
real(r8), dimension(:), allocatable kvmin
 
real(r8), dimension(:), allocatable kvmax
 
real(dp), dimension(:), allocatable ad_cg_beta
 
real(dp), dimension(:), allocatable tl_cg_beta
 
real(dp), dimension(:,:), allocatable cg_alpha
 
real(dp), dimension(:,:), allocatable cg_beta
 
real(dp), dimension(:,:), allocatable cg_tau
 
real(dpritzmaxerr
 
real(dp), dimension(:), allocatable ritz
 
real(dp), dimension(:,:,:), allocatable cg_zv
 
real(dp), dimension(:), allocatable ad_cg_delta
 
real(dp), dimension(:), allocatable ad_cg_qg
 
real(dp), dimension(:), allocatable ad_cg_gnorm
 
real(dp), dimension(:), allocatable tl_cg_delta
 
real(dp), dimension(:), allocatable tl_cg_qg
 
real(dp), dimension(:), allocatable tl_cg_gnorm_v
 
real(dp), dimension(:), allocatable ad_cg_gnorm_v
 
real(dp), dimension(:,:), allocatable cg_delta
 
real(dp), dimension(:,:), allocatable cg_gamma
 
real(dp), dimension(:), allocatable cg_gnorm
 
real(dp), dimension(:), allocatable cg_gnorm_v
 
real(dp), dimension(:), allocatable cg_gnorm_y
 
real(dp), dimension(:,:), allocatable cg_greduc
 
real(dp), dimension(:,:), allocatable cg_qg
 
real(dp), dimension(:,:), allocatable cg_tmatrix
 
real(dp), dimension(:,:), allocatable cg_zu
 
real(dp), dimension(:,:), allocatable cg_ritz
 
real(dp), dimension(:,:), allocatable cg_ritzerr
 
logical, dimension(:), allocatable wrtzetaref
 
real(r8), dimension(:,:), allocatable ae_delta
 
real(r8), dimension(:,:), allocatable ae_beta
 
real(r8), dimension(:,:), allocatable ae_zv
 
real(r8), dimension(:), allocatable ae_trace
 
real(r8), dimension(:), allocatable ae_gnorm
 
real(r8), dimension(:,:), allocatable ae_tmatrix
 
real(r8), dimension(:,:), allocatable ae_zu
 
real(r8), dimension(:,:), allocatable ae_ritz
 
real(r8), dimension(:,:), allocatable ae_ritzerr
 
real(r8), dimension(:,:), allocatable zlanczos_coef
 
real(r8), dimension(:,:), allocatable zlanczos_inv
 
real(r8), dimension(:,:), allocatable zlanczos_err
 
real(r8), dimension(:), allocatable zlanczos_diag
 
real(r8), dimension(:), allocatable zlanczos_offdiag
 
real(r8), dimension(:,:), allocatable gsmatrix
 
real(r8), dimension(:,:), allocatable gsmatinv
 
logical lhessianev
 
logical lhotstart
 
logical lprecond
 
logical lritz
 
logical laugweak =.FALSE.
 
integer nritzev
 
integer nconvritz = 0
 
real(r8cg_gammam1 = 0.0_r8
 
real(r8cg_sigmam1 = 0.0_r8
 
real(r8cg_rnorm = 0.0_r8
 
real(r8graderr
 
real(r8hevecerr
 
real(r8dotproduct
 
real(r8addotproduct
 
integer ig1count
 
integer ig2count
 
real(r8), dimension(1000) g1
 
real(r8), dimension(1000) g2
 
logical, dimension(:), allocatable wrtforce
 
logical, dimension(:), allocatable lweakrelax
 
logical, dimension(:), allocatable lsenfct
 
logical, dimension(:), allocatable lsen4dvar
 
logical, dimension(:), allocatable lobspace
 
integer, dimension(:), allocatable ntimes_ana
 
integer, dimension(:), allocatable ntimes_fct
 
real(r8), dimension(:), allocatable forcetime
 

Function/Subroutine Documentation

◆ allocate_fourdvar()

subroutine, public mod_fourdvar::allocate_fourdvar

Definition at line 923 of file mod_fourdvar.F.

924!
925!=======================================================================
926! !
927! This routine allocates several variables in module "mod_fourdvar" !
928! for all nested grids. !
929! !
930!=======================================================================
931!
932!-----------------------------------------------------------------------
933! Allocate module variables due to nested grids.
934!-----------------------------------------------------------------------
935!
936 IF (.not.allocated(load_zobs)) THEN
937 allocate ( load_zobs(ngrids) )
938 dmem(1)=dmem(1)+real(ngrids,r8)
939 END IF
940 IF (.not.allocated(wrote_zobs)) THEN
941 allocate ( wrote_zobs(ngrids) )
942 dmem(1)=dmem(1)+real(ngrids,r8)
943 END IF
944
945 IF (.not.allocated(nstatevar)) THEN
946 allocate ( nstatevar(ngrids) )
947 dmem(1)=dmem(1)+real(ngrids,r8)
948 END IF
949 IF (.not.allocated(nobsvar)) THEN
950 allocate ( nobsvar(ngrids) )
951 dmem(1)=dmem(1)+real(ngrids,r8)
952 END IF
953 IF (.not.allocated(ndatum)) THEN
954 allocate ( ndatum(ngrids) )
955 dmem(1)=dmem(1)+real(ngrids,r8)
956 END IF
957 IF (.not.allocated(nsurvey)) THEN
958 allocate ( nsurvey(ngrids) )
959 dmem(1)=dmem(1)+real(ngrids,r8)
960 END IF
961 IF (.not.allocated(obssurvey)) THEN
962 allocate ( obssurvey(ngrids) )
963 dmem(1)=dmem(1)+real(ngrids,r8)
964 END IF
965 IF (.not.allocated(nobs)) THEN
966 allocate ( nobs(ngrids) )
967 dmem(1)=dmem(1)+real(ngrids,r8)
968 END IF
969
970 IF (.not.allocated(nstrobs)) THEN
971 allocate ( nstrobs(ngrids) )
972 dmem(1)=dmem(1)+real(ngrids,r8)
973 END IF
974 IF (.not.allocated(nendobs)) THEN
975 allocate ( nendobs(ngrids) )
976 dmem(1)=dmem(1)+real(ngrids,r8)
977 END IF
978
979 IF (.not.allocated(nmethod)) THEN
980 allocate ( nmethod(ngrids) )
981 dmem(1)=dmem(1)+real(ngrids,r8)
982 END IF
983 IF (.not.allocated(rscheme)) THEN
984 allocate ( rscheme(ngrids) )
985 dmem(1)=dmem(1)+real(ngrids,r8)
986 END IF
987
988 IF (.not.allocated(optimality)) THEN
989 allocate ( optimality(ngrids) )
990 dmem(1)=dmem(1)+real(ngrids,r8)
991 END IF
992 IF (.not.allocated(processobs)) THEN
993 allocate ( processobs(ngrids) )
994 dmem(1)=dmem(1)+real(ngrids,r8)
995 END IF
996
997# ifdef SP4DVAR
998 IF (.not.allocated(lsadd)) THEN
999 allocate ( lsadd(ngrids) )
1000 dmem(1)=dmem(1)+real(ngrids,r8)
1001 END IF
1002# endif
1003
1004 IF (.not.allocated(wrtnlmod)) THEN
1005 allocate ( wrtnlmod(ngrids) )
1006 dmem(1)=dmem(1)+real(ngrids,r8)
1007 END IF
1008 IF (.not.allocated(wrtobsscale)) THEN
1009 allocate ( wrtobsscale(ngrids) )
1010 dmem(1)=dmem(1)+real(ngrids,r8)
1011 END IF
1012 IF (.not.allocated(wrtrpmod)) THEN
1013 allocate ( wrtrpmod(ngrids) )
1014 dmem(1)=dmem(1)+real(ngrids,r8)
1015 END IF
1016 IF (.not.allocated(wrttlmod)) THEN
1017 allocate ( wrttlmod(ngrids) )
1018 dmem(1)=dmem(1)+real(ngrids,r8)
1019 END IF
1020 IF (.not.allocated(wrtmisfit)) THEN
1021 allocate ( wrtmisfit(ngrids) )
1022 dmem(1)=dmem(1)+real(ngrids,r8)
1023 END IF
1024
1025# if defined I4DVAR_ANA_SENSITIVITY
1026# ifdef OBS_IMPACT
1027 IF (.not.allocated(wrtimpact_tot)) THEN
1028 allocate ( wrtimpact_tot(ngrids) )
1029 dmem(1)=dmem(1)+real(ngrids,r8)
1030 END IF
1031# endif
1032# ifdef OBS_IMPACT_SPLIT
1033 IF (.not.allocated(wrtimpact_ic)) THEN
1034 allocate ( wrtimpact_ic(ngrids) )
1035 dmem(1)=dmem(1)+real(ngrids,r8)
1036 END IF
1037# if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
1038 IF (.not.allocated(wrtimpact_fc)) THEN
1039 allocate ( wrtimpact_fc(ngrids) )
1040 dmem(1)=dmem(1)+real(ngrids,r8)
1041 END IF
1042# endif
1043# if defined ADJUST_BOUNDARY
1044 IF (.not.allocated(wrtimpact_bc)) THEN
1045 allocate ( wrtimpact_bc(ngrids) )
1046 dmem(1)=dmem(1)+real(ngrids,r8)
1047 END IF
1048# endif
1049# endif
1050# endif
1051
1052 IF (.not.allocated(haveadmod)) THEN
1053 allocate ( haveadmod(ngrids) )
1054 dmem(1)=dmem(1)+real(ngrids,r8)
1055 END IF
1056 IF (.not.allocated(havenlmod)) THEN
1057 allocate ( havenlmod(ngrids) )
1058 dmem(1)=dmem(1)+real(ngrids,r8)
1059 END IF
1060 IF (.not.allocated(havetlmod)) THEN
1061 allocate ( havetlmod(ngrids) )
1062 dmem(1)=dmem(1)+real(ngrids,r8)
1063 END IF
1064 IF (.not.allocated(haveobsmeta)) THEN
1065 allocate ( haveobsmeta(ngrids) )
1066 dmem(1)=dmem(1)+real(ngrids,r8)
1067 END IF
1068
1069 IF (.not.allocated(khmin)) THEN
1070 allocate ( khmin(ngrids) )
1071 dmem(1)=dmem(1)+real(ngrids,r8)
1072 END IF
1073 IF (.not.allocated(khmax)) THEN
1074 allocate ( khmax(ngrids) )
1075 dmem(1)=dmem(1)+real(ngrids,r8)
1076 END IF
1077 IF (.not.allocated(kvmin)) THEN
1078 allocate ( kvmin(ngrids) )
1079 dmem(1)=dmem(1)+real(ngrids,r8)
1080 END IF
1081 IF (.not.allocated(kvmax)) THEN
1082 allocate ( kvmax(ngrids) )
1083 dmem(1)=dmem(1)+real(ngrids,r8)
1084 END IF
1085
1086# if defined BALANCE_OPERATOR
1087 IF (.not.allocated(wrtzetaref)) THEN
1088 allocate ( wrtzetaref(ngrids) )
1089 dmem(1)=dmem(1)+real(ngrids,r8)
1090 END IF
1091# endif
1092
1093# ifdef WEAK_CONSTRAINT
1094 IF (.not.allocated(wrtforce)) THEN
1095 allocate ( wrtforce(ngrids) )
1096 dmem(1)=dmem(1)+real(ngrids,r8)
1097 END IF
1098# ifdef RPM_RELAXATION
1099 IF (.not.allocated(lweakrelax)) THEN
1100 allocate ( lweakrelax(ngrids) )
1101 dmem(1)=dmem(1)+real(ngrids,r8)
1102 END IF
1103# endif
1104# ifdef RBL4DVAR_FCT_SENSITIVITY
1105 IF (.not.allocated(lsenfct)) THEN
1106 allocate ( lsenfct(ngrids) )
1107 dmem(1)=dmem(1)+real(ngrids,r8)
1108 END IF
1109# endif
1110# if defined RBL4DVAR_ANA_SENSITIVITY || \
1111 defined rbl4dvar_fct_sensitivity || \
1112 defined r4dvar_ana_sensitivity
1113 IF (.not.allocated(lsen4dvar)) THEN
1114 allocate ( lsen4dvar(ngrids) )
1115 dmem(1)=dmem(1)+real(ngrids,r8)
1116 END IF
1117# ifdef OBS_SPACE
1118 IF (.not.allocated(lobspace)) THEN
1119 allocate ( lobspace(ngrids) )
1120 dmem(1)=dmem(1)+real(ngrids,r8)
1121 END IF
1122# ifndef OBS_IMPACT
1123 IF (.not.allocated(ladjvar)) THEN
1124 allocate ( ladjvar(ngrids) )
1125 dmem(1)=dmem(1)+real(ngrids,r8)
1126 END IF
1127# endif
1128# endif
1129# endif
1130# ifdef RBL4DVAR_FCT_SENSITIVITY
1131 IF (.not.allocated(ntimes_ana)) THEN
1132 allocate ( ntimes_ana(ngrids) )
1133 dmem(1)=dmem(1)+real(ngrids,r8)
1134 END IF
1135 IF (.not.allocated(ntimes_fct)) THEN
1136 allocate ( ntimes_fct(ngrids) )
1137 dmem(1)=dmem(1)+real(ngrids,r8)
1138 END IF
1139# endif
1140 IF (.not.allocated(forcetime)) THEN
1141 allocate ( forcetime(ngrids) )
1142 dmem(1)=dmem(1)+real(ngrids,r8)
1143 END IF
1144# endif
1145!
1146 RETURN

References mod_param::dmem, forcetime, haveadmod, havenlmod, haveobsmeta, havetlmod, khmax, khmin, kvmax, kvmin, load_zobs, lobspace, lsadd, lsen4dvar, lsenfct, lweakrelax, ndatum, nendobs, mod_param::ngrids, nmethod, nobs, nobsvar, nstatevar, nstrobs, nsurvey, ntimes_ana, ntimes_fct, obssurvey, optimality, processobs, rscheme, wrote_zobs, wrtforce, wrtimpact_bc, wrtimpact_fc, wrtimpact_ic, wrtimpact_tot, wrtmisfit, wrtnlmod, wrtobsscale, wrtrpmod, wrttlmod, and wrtzetaref.

◆ deallocate_fourdvar()

subroutine, public mod_fourdvar::deallocate_fourdvar

Definition at line 1149 of file mod_fourdvar.F.

1150!
1151!=======================================================================
1152! !
1153! This routine deallocates only the variables in module associated !
1154! with observation dimension parameters Ndatum, Nsurvey, and Mobs. !
1155! !
1156!=======================================================================
1157
1158# ifdef OBSERVATIONS
1159!
1160!-----------------------------------------------------------------------
1161! Deallocate structure variables for each nested grids.
1162!-----------------------------------------------------------------------
1163!
1164 IF (allocated(fourdvar)) THEN
1165 deallocate (fourdvar)
1166 END IF
1167!
1168!-----------------------------------------------------------------------
1169! Deallocate model and observation variables.
1170!-----------------------------------------------------------------------
1171!
1172# ifdef CURVGRID
1173 IF (allocated(obsangler)) THEN
1174 deallocate (obsangler)
1175 END IF
1176# endif
1177
1178 IF (allocated(obstype)) THEN
1179 deallocate (obstype)
1180 END IF
1181
1182 IF (allocated(obsprov)) THEN
1183 deallocate (obsprov)
1184 END IF
1185
1186 IF (allocated(obserr)) THEN
1187 deallocate (obserr)
1188 END IF
1189
1190 IF (allocated(obsmeta)) THEN
1191 deallocate (obsmeta)
1192 END IF
1193
1194 IF (allocated(obsscale)) THEN
1195 deallocate (obsscale)
1196 END IF
1197
1198# ifndef WEAK_CONSTRAINT
1199 IF (allocated(obsscaleglobal)) THEN
1200 deallocate (obsscaleglobal)
1201 END IF
1202# endif
1203
1204 IF (allocated(obsvetting)) THEN
1205 deallocate (obsvetting)
1206 END IF
1207
1208 IF (allocated(obsval)) THEN
1209 deallocate (obsval)
1210 END IF
1211
1212 IF (allocated(tobs)) THEN
1213 deallocate (tobs)
1214 END IF
1215
1216 IF (allocated(xobs)) THEN
1217 deallocate (xobs)
1218 END IF
1219
1220 IF (allocated(yobs)) THEN
1221 deallocate (yobs)
1222 END IF
1223
1224# ifdef SOLVE3D
1225 IF (allocated(zobs)) THEN
1226 deallocate (zobs)
1227 END IF
1228# endif
1229
1230 IF (allocated(admodval)) THEN
1231 deallocate (admodval)
1232 END IF
1233
1234 IF (allocated(nlmodval)) THEN
1235 deallocate (nlmodval)
1236 END IF
1237
1238 IF (allocated(misfit)) THEN
1239 deallocate (misfit)
1240 END IF
1241
1242 IF (allocated(unvetted)) THEN
1243 deallocate (unvetted)
1244 END IF
1245
1246 IF (allocated(uradial)) THEN
1247 deallocate (uradial)
1248 END IF
1249
1250 IF (allocated(vradial)) THEN
1251 deallocate (vradial)
1252 END IF
1253
1254# ifndef VERIFICATION
1255 IF (allocated(ubgerr)) THEN
1256 deallocate (ubgerr)
1257 END IF
1258
1259 IF (allocated(vbgerr)) THEN
1260 deallocate (vbgerr)
1261 END IF
1262# endif
1263
1264# ifdef FOUR_DVAR
1265 IF (allocated(bgerr)) THEN
1266 deallocate (bgerr)
1267 END IF
1268
1269 IF (allocated(nlincrement)) THEN
1270 deallocate (nlincrement)
1271 END IF
1272
1273 IF (allocated(innovation)) THEN
1274 deallocate (innovation)
1275 END IF
1276
1277 IF (allocated(residual)) THEN
1278 deallocate (residual)
1279 END IF
1280# endif
1281
1282# ifdef BGQC
1283 IF (allocated(bgthresh)) THEN
1284 deallocate (bgthresh)
1285 END IF
1286# endif
1287
1288# ifdef TLM_OBS
1289 IF (allocated(tlmodval)) THEN
1290 deallocate (tlmodval)
1291 END IF
1292
1293# if defined ARRAY_MODES || defined CLIPPING || \
1294 defined rbl4dvar || defined r4dvar || \
1295 defined sensitivity_4dvar || defined tl_rbl4dvar || \
1296 defined tl_r4dvar
1297 IF (allocated(tlmodval_s)) THEN
1298 deallocate (tlmodval_s)
1299 END IF
1300
1301# if defined RPCG && \
1302 (defined sensitivity_4dvar || defined tl_rbl4dvar || \
1303 defined tl_r4dvar)
1304
1305 IF (allocated(tl_tlmodval_s)) THEN
1306 deallocate (tl_tlmodval_s)
1307 END IF
1308
1309 IF (allocated(ad_tlmodval_s)) THEN
1310 deallocate (ad_tlmodval_s)
1311 END IF
1312# endif
1313
1314# endif
1315
1316# if defined SENSITIVITY_4DVAR || defined RBL4DVAR || \
1317 defined tl_rbl4dvar
1318
1319 IF (allocated(bckmodval)) THEN
1320 deallocate (bckmodval)
1321 END IF
1322
1323# ifdef RPCG
1324 IF (allocated(hbk)) THEN
1325 deallocate (hbk)
1326 END IF
1327# if defined SENSITIVITY_4DVAR || defined TL_RBL4DVAR || \
1328 defined tl_r4dvar
1329
1330 IF (allocated(tl_bckmodval)) THEN
1331 deallocate (tl_bckmodval)
1332 END IF
1333 IF (allocated(ad_bckmodval)) THEN
1334 deallocate (ad_bckmodval)
1335 END IF
1336 IF (allocated(tl_hbk)) THEN
1337 deallocate (tl_hbk)
1338 END IF
1339 IF (allocated(ad_hbk)) THEN
1340 deallocate (ad_hbk)
1341 END IF
1342# endif
1343# endif
1344# endif
1345# endif
1346
1347# ifdef WEAK_CONSTRAINT
1348!
1349! Deallocate and initialize weak constraint conjugate gradient vectors.
1350!
1351# if defined SENSITIVITY_4DVAR || defined TL_RBL4DVAR || \
1352 defined tl_r4dvar
1353 IF (allocated(ad_zcglwk)) THEN
1354 deallocate (ad_zcglwk)
1355 END IF
1356
1357 IF (allocated(ad_zgrad0)) THEN
1358 deallocate (ad_zgrad0)
1359 END IF
1360
1361 IF (allocated(ad_cg_innov)) THEN
1362 deallocate (ad_cg_innov)
1363 END IF
1364
1365 IF (allocated(ad_cg_pxsave)) THEN
1366 deallocate (ad_cg_pxsave)
1367 END IF
1368
1369 IF (allocated(tl_zcglwk)) THEN
1370 deallocate (tl_zcglwk)
1371 END IF
1372
1373 IF (allocated(tl_zgrad0)) THEN
1374 deallocate (tl_zgrad0)
1375 END IF
1376
1377 IF (allocated(tl_cg_innov)) THEN
1378 deallocate (tl_cg_innov)
1379 END IF
1380
1381 IF (allocated(tl_cg_pxsave)) THEN
1382 deallocate (tl_cg_pxsave)
1383 END IF
1384
1385 IF (allocated(tl_obsval)) THEN
1386 deallocate (tl_obsval)
1387 END IF
1388
1389# ifdef IMPACT_INNER
1390 IF (allocated(ad_obsval)) THEN
1391 deallocate (ad_obsval)
1392 END IF
1393# else
1394 IF (allocated(ad_obsval)) THEN
1395 deallocate (ad_obsval)
1396 END IF
1397# endif
1398# endif
1399
1400# ifdef RPCG
1401 IF (allocated(zcglwk)) THEN
1402 deallocate (zcglwk)
1403 END IF
1404
1405 IF (allocated(vcglwk)) THEN
1406 deallocate (vcglwk)
1407 END IF
1408
1409 IF (allocated(vcglev)) THEN
1410 deallocate (vcglev)
1411 END IF
1412
1413 IF (allocated(zgrad0)) THEN
1414 deallocate (zgrad0)
1415 END IF
1416
1417 IF (allocated(vgrad0)) THEN
1418 deallocate (vgrad0)
1419 END IF
1420
1421 IF (allocated(vgrad0s)) THEN
1422 deallocate (vgrad0s)
1423 END IF
1424
1425 IF (allocated(cg_innov)) THEN
1426 deallocate (cg_innov)
1427 END IF
1428
1429# if defined RBL4DVAR_ANA_SENSITIVITY || \
1430 defined rbl4dvar_fct_sensitivity || \
1431 defined r4dvar_ana_sensitivity || \
1432 defined tl_rbl4dvar || \
1433 defined tl_r4dvar
1434
1435 IF (allocated(ad_vcglwk)) THEN
1436 deallocate (ad_vcglwk)
1437 END IF
1438
1439 IF (allocated(ad_vgrad0)) THEN
1440 deallocate (ad_vgrad0)
1441 END IF
1442
1443 IF (allocated(ad_vgrad0s)) THEN
1444 deallocate (ad_vgrad0s)
1445 END IF
1446
1447 IF (allocated(ad_jb0)) THEN
1448 deallocate (ad_jb0)
1449 END IF
1450
1451 IF (allocated(ad_obserr)) THEN
1452 deallocate (ad_obserr)
1453 END IF
1454
1455 IF (allocated(ad_obsscale)) THEN
1456 deallocate (ad_obsscale)
1457 END IF
1458
1459 IF (allocated(ad_obsval)) THEN
1460 deallocate (ad_obsval)
1461 END IF
1462
1463 IF (allocated(tl_vcglwk)) THEN
1464 deallocate (tl_vcglwk)
1465 END IF
1466
1467 IF (allocated(tl_vgrad0)) THEN
1468 deallocate (tl_vgrad0)
1469 END IF
1470
1471 IF (allocated(tl_vgrad0s)) THEN
1472 deallocate (tl_vgrad0s)
1473 END IF
1474
1475 IF (allocated(tl_jb0)) THEN
1476 deallocate (tl_jb0)
1477 END IF
1478
1479 IF (allocated(tl_obserr)) THEN
1480 deallocate (tl_obserr)
1481 END IF
1482
1483 IF (allocated(tl_obsscale)) THEN
1484 deallocate (tl_obsscale)
1485 END IF
1486
1487 IF (allocated(tl_obsval)) THEN
1488 deallocate (tl_obsval)
1489 END IF
1490
1491 IF (allocated(tl_cg_gnorm_v)) THEN
1492 deallocate (tl_cg_gnorm_v)
1493 END IF
1494
1495 IF (allocated(ad_cg_gnorm_v)) THEN
1496 deallocate (ad_cg_gnorm_v)
1497 END IF
1498# endif
1499
1500# else
1501
1502 IF (allocated(cg_pxsave)) THEN
1503 deallocate (cg_pxsave)
1504 END IF
1505
1506 IF (allocated(zcglwk)) THEN
1507 deallocate (zcglwk)
1508 END IF
1509
1510 IF (allocated(vcglev)) THEN
1511 deallocate (vcglev)
1512 END IF
1513
1514 IF (allocated(zgrad0)) THEN
1515 deallocate (zgrad0)
1516 END IF
1517
1518 IF (allocated(vgrad0)) THEN
1519 deallocate (vgrad0)
1520 END IF
1521
1522 IF (allocated(vgrad0s)) THEN
1523 deallocate (vgrad0s)
1524 END IF
1525
1526 IF (allocated(gdgw)) THEN
1527 deallocate (gdgw)
1528 END IF
1529
1530 IF (allocated(cg_innov)) THEN
1531 deallocate (cg_innov)
1532 END IF
1533
1534 IF (allocated(cg_pxsave)) THEN
1535 deallocate (cg_pxsave)
1536 END IF
1537
1538# endif
1539# endif
1540# endif
1541!
1542 RETURN

References ad_cg_gnorm_v, ad_cg_innov, ad_cg_pxsave, ad_jb0, ad_obserr, ad_obsscale, ad_obsval, ad_vcglwk, ad_vgrad0, ad_vgrad0s, ad_zcglwk, ad_zgrad0, admodval, bgerr, bgthresh, cg_innov, fourdvar, gdgw, innovation, misfit, nlincrement, nlmodval, obsangler, obserr, obsmeta, obsprov, obsscale, obstype, obsval, obsvetting, residual, tl_cg_gnorm_v, tl_cg_innov, tl_cg_pxsave, tl_jb0, tl_obserr, tl_obsscale, tl_obsval, tl_vcglwk, tl_vgrad0, tl_vgrad0s, tl_zcglwk, tl_zgrad0, tobs, unvetted, uradial, vcglev, vcglwk, vgrad0, vgrad0s, vradial, xobs, yobs, zcglwk, zgrad0, and zobs.

Referenced by mod_arrays::roms_deallocate_arrays(), and roms_kernel_mod::roms_run().

Here is the caller graph for this function:

◆ initialize_fourdvar()

subroutine, public mod_fourdvar::initialize_fourdvar

Definition at line 1545 of file mod_fourdvar.F.

1546!
1547!=======================================================================
1548! !
1549! This routine initializes several variables in module "mod_fourdvar" !
1550! for all nested grids. !
1551! !
1552!=======================================================================
1553!
1554 USE mod_parallel
1555 USE mod_scalars
1556 USE mod_iounits
1557 USE mod_ncparam
1558 USE mod_netcdf
1559# if defined PIO_LIB && defined DISTRIBUTE
1560 USE mod_pio_netcdf
1561# endif
1562!
1563 USE strings_mod, ONLY : founderror
1564 USE strings_mod, ONLY : uppercase
1565!
1566! Local variable declarations.
1567!
1568 integer :: my_NobsVar, i, icount, ng
1569# if defined BALANCE_OPERATOR && defined ZETA_ELLIPTIC
1570 integer :: LBi, UBi, LBj, UBj
1571 integer :: tile
1572
1573 real(r8) :: size2d
1574# endif
1575 real(r8), parameter :: IniVal = 0.0_r8
1576!
1577 character (len=*), parameter :: MyFile = &
1578 & __FILE__//", initialize_fourdvar"
1579
1580 sourcefile=myfile
1581
1582# ifdef OBSERVATIONS
1583!
1584!-----------------------------------------------------------------------
1585! Inquire observations NetCDF and determine the maximum dimension of
1586! several observations arrays.
1587!-----------------------------------------------------------------------
1588!
1589! Inquire about the size of the "datum" unlimitted dimension and
1590! "survey" dimension.
1591!
1592 DO ng=1,ngrids
1593 SELECT CASE (obs(ng)%IOtype)
1594 CASE (io_nf90)
1595 CALL netcdf_get_dim (ng, inlm, trim(obs(ng)%name))
1596# if defined PIO_LIB && defined DISTRIBUTE
1597 CASE (io_pio)
1598 CALL pio_netcdf_get_dim (ng, inlm, trim(obs(ng)%name))
1599# endif
1600 CASE DEFAULT
1601 IF (master) WRITE (stdout,10) obs(ng)%IOtype
1602 exit_flag=3
1603 END SELECT
1604 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1605!
1606 ndatum(ng)=0
1607 nsurvey(ng)=0
1608 DO i=1,n_dim
1609 IF (trim(dim_name(i)).eq.'datum') then
1610 ndatum(ng)=dim_size(i)
1611 ELSE IF (trim(dim_name(i)).eq.'survey') THEN
1612 nsurvey(ng)=dim_size(i)
1613 ELSE IF (trim(dim_name(i)).eq.'state_variable') THEN
1614 my_nobsvar=dim_size(i)
1615 END IF
1616 END DO
1617 IF (ndatum(ng).eq.0) THEN
1618 WRITE (stdout,20) 'datum', trim(obs(ng)%name)
1619 exit_flag=2
1620 RETURN
1621 END IF
1622 IF (nsurvey(ng).eq.0) THEN
1623 WRITE (stdout,20) 'survey', trim(obs(ng)%name)
1624 exit_flag=2
1625 RETURN
1626 END IF
1627 END DO
1628!
1629!-----------------------------------------------------------------------
1630! Allocate module variables.
1631!-----------------------------------------------------------------------
1632!
1633# if defined BALANCE_OPERATOR && defined ZETA_ELLIPTIC
1634# ifdef DISTRIBUTE
1635 tile=myrank
1636# else
1637 tile=0
1638# endif
1639# endif
1640!
1641! Allocate structure (1:Ngrids).
1642!
1643 IF (.not.allocated( fourdvar)) THEN
1644 allocate ( fourdvar(ngrids) )
1645 END IF
1646!
1647! Allocate variables in structure.
1648!
1649 DO ng=1,ngrids
1650
1651# if defined BALANCE_OPERATOR && defined ZETA_ELLIPTIC
1652 lbi=bounds(ng)%LBi(tile)
1653 ubi=bounds(ng)%UBi(tile)
1654 lbj=bounds(ng)%LBj(tile)
1655 ubj=bounds(ng)%UBj(tile)
1656
1657 size2d=real((ubi-lbi+1)*(ubj-lbj+1),r8)
1658# endif
1659!
1660! Number of state variables associated with data assimilation.
1661!
1662! zeta, ubar, vbar, Uvel, Vvel, Tvar(1:MT)
1663! Ustr, Vstr, Tsur(1:MT)
1664!
1665! The rest are ignored.
1666!
1667# ifdef SOLVE3D
1668 nstatevar(ng)=5+nt(ng)
1669# ifdef ADJUST_STFLUX
1670 nstatevar(ng)=nstatevar(ng)+nt(ng)
1671# endif
1672# else
1673 nstatevar(ng)=3
1674# endif
1675# ifdef ADJUST_WSTRESS
1676 nstatevar(ng)=nstatevar(ng)+2
1677# endif
1678
1679 nobsvar(ng)=nstatevar(ng)+nextraobs
1680
1681 IF (.not.allocated(fourdvar(ng) % NobsSurvey)) THEN
1682 allocate ( fourdvar(ng) % NobsSurvey(nsurvey(ng)) )
1683 dmem(ng)=dmem(ng)+real(nsurvey(ng),r8)
1684 fourdvar(ng) % NobsSurvey = 0
1685 END IF
1686
1687 IF (.not.allocated(fourdvar(ng) % SurveyTime)) THEN
1688 allocate ( fourdvar(ng) % SurveyTime(nsurvey(ng)) )
1689 dmem(ng)=dmem(ng)+real(nsurvey(ng),r8)
1690 fourdvar(ng) % SurveyTime = inival
1691 END IF
1692
1693# ifdef FOUR_DVAR
1694 IF (.not.allocated(fourdvar(ng) % BackCost)) THEN
1695 allocate ( fourdvar(ng) % BackCost(0:nstatevar(ng)) )
1696 dmem(ng)=dmem(ng)+real(nstatevar(ng)+1,r8)
1697 fourdvar(ng) % BackCost = inival
1698 END IF
1699
1700 IF (.not.allocated(fourdvar(ng) % Cost0)) THEN
1701 allocate ( fourdvar(ng) % Cost0(nouter) )
1702 dmem(ng)=dmem(ng)+real(nouter,r8)
1703 fourdvar(ng) % Cost0 = inival
1704 END IF
1705
1706 IF (.not.allocated(fourdvar(ng) % CostFun)) THEN
1707 allocate ( fourdvar(ng) % CostFun(0:nobsvar(ng)) )
1708 dmem(ng)=dmem(ng)+real(nobsvar(ng)+1,r8)
1709 fourdvar(ng) % CostFun = inival
1710 END IF
1711
1712 IF (.not.allocated(fourdvar(ng) % CostFunOld)) THEN
1713 allocate ( fourdvar(ng) % CostFunOld(0:nobsvar(ng)) )
1714 dmem(ng)=dmem(ng)+real(nobsvar(ng)+1,r8)
1715 fourdvar(ng) % CostFunOld = inival
1716 END IF
1717
1718 IF (.not.allocated(fourdvar(ng) % CostNorm)) THEN
1719 allocate ( fourdvar(ng) % CostNorm(0:nobsvar(ng)) )
1720 dmem(ng)=dmem(ng)+real(nobsvar(ng)+1,r8)
1721 fourdvar(ng) % CostNorm = inival
1722 END IF
1723
1724 IF (.not.allocated(fourdvar(ng) % ObsCost)) THEN
1725 allocate ( fourdvar(ng) % ObsCost(0:nobsvar(ng)) )
1726 dmem(ng)=dmem(ng)+real(nobsvar(ng)+1,r8)
1727 fourdvar(ng) % ObsCost = inival
1728 END IF
1729
1730 IF (.not.allocated(fourdvar(ng) % DataPenalty)) THEN
1731 allocate ( fourdvar(ng) % DataPenalty(0:nobsvar(ng)) )
1732 dmem(ng)=dmem(ng)+real(nobsvar(ng)+1,r8)
1733 fourdvar(ng) % DataPenalty = inival
1734 END IF
1735
1736# if defined DATALESS_LOOPS || defined RBL4DVAR || \
1737 defined sensitivity_4dvar || defined sp4dvar || \
1738 defined tl_rbl4dvar
1739
1740 IF (.not.allocated(fourdvar(ng) % NLPenalty)) THEN
1741 allocate ( fourdvar(ng) % NLPenalty(0:nobsvar(ng)) )
1742 dmem(ng)=dmem(ng)+real(nobsvar(ng)+1,r8)
1743 fourdvar(ng) % NLPenalty = inival
1744 END IF
1745# endif
1746
1747 IF (.not.allocated(fourdvar(ng) % NLobsCost)) THEN
1748 allocate ( fourdvar(ng) % NLobsCost(0:nobsvar(ng)) )
1749 dmem(ng)=dmem(ng)+real(nobsvar(ng)+1,r8)
1750 fourdvar(ng) % NLobsCost = inival
1751 END IF
1752
1753# if defined BALANCE_OPERATOR && defined ZETA_ELLIPTIC
1754 IF (.not.allocated(fourdvar(ng) % bc_ak)) THEN
1755 allocate ( fourdvar(ng) % bc_ak(nbico(ng)) )
1756 dmem(ng)=dmem(ng)+real(nbico(ng),r8)
1757 fourdvar(ng) % bc_ak = inival
1758 END IF
1759
1760 IF (.not.allocated(fourdvar(ng) % bc_bk)) THEN
1761 allocate ( fourdvar(ng) % bc_bk(nbico(ng)) )
1762 dmem(ng)=dmem(ng)+real(nbico(ng),r8)
1763 fourdvar(ng) % bc_bk = inival
1764 END IF
1765
1766 IF (.not.allocated(fourdvar(ng) % zdf1)) THEN
1767 allocate ( fourdvar(ng) % zdf1(nbico(ng)) )
1768 dmem(ng)=dmem(ng)+real(nbico(ng),r8)
1769 fourdvar(ng) % zdf1 = inival
1770 END IF
1771
1772 IF (.not.allocated(fourdvar(ng) % zdf2)) THEN
1773 allocate ( fourdvar(ng) % zdf2(nbico(ng)) )
1774 dmem(ng)=dmem(ng)+real(nbico(ng),r8)
1775 fourdvar(ng) % zdf2 = inival
1776 END IF
1777
1778 IF (.not.allocated(fourdvar(ng) % zdf3)) THEN
1779 allocate ( fourdvar(ng) % zdf3(nbico(ng)) )
1780 dmem(ng)=dmem(ng)+real(nbico(ng),r8)
1781 fourdvar(ng) % zdf3 = inival
1782 END IF
1783
1784 IF (.not.allocated(fourdvar(ng) % pc_r2d)) THEN
1785 allocate ( fourdvar(ng) % pc_r2d(lbi:ubi,lbj:ubj) )
1786 dmem(ng)=dmem(ng)+size2d
1787 fourdvar(ng) % pc_r2d = inival
1788 END IF
1789
1790 IF (.not.allocated(fourdvar(ng) % rhs_r2d)) THEN
1791 allocate ( fourdvar(ng) % rhs_r2d(lbi:ubi,lbj:ubj) )
1792 dmem(ng)=dmem(ng)+size2d
1793 fourdvar(ng) % rhs_r2d = inival
1794 END IF
1795
1796 IF (.not.allocated(fourdvar(ng) % r2d_ref)) THEN
1797 allocate ( fourdvar(ng) % r2d_ref(lbi:ubi,lbj:ubj) )
1798 dmem(ng)=dmem(ng)+size2d
1799 fourdvar(ng) % r2d_ref = inival
1800 END IF
1801
1802 IF (.not.allocated(fourdvar(ng) % zeta_ref)) THEN
1803 allocate ( fourdvar(ng) % zeta_ref(lbi:ubi,lbj:ubj) )
1804 dmem(ng)=dmem(ng)+size2d
1805 fourdvar(ng) % zeta_ref = inival
1806 END IF
1807
1808 IF (.not.allocated(fourdvar(ng) % p_r2d)) THEN
1809 allocate ( fourdvar(ng) % p_r2d(lbi:ubi,lbj:ubj,nbico(ng)) )
1810 dmem(ng)=dmem(ng)+real(nbico(ng),r8)*size2d
1811 fourdvar(ng) % p_r2d = inival
1812 END IF
1813
1814 IF (.not.allocated(fourdvar(ng) % r_r2d)) THEN
1815 allocate ( fourdvar(ng) % r_r2d(lbi:ubi,lbj:ubj,nbico(ng)) )
1816 dmem(ng)=dmem(ng)+real(nbico(ng),r8)*size2d
1817 fourdvar(ng) % r_r2d = inival
1818 END IF
1819
1820 IF (.not.allocated(fourdvar(ng) % bp_r2d)) THEN
1821 allocate ( fourdvar(ng) % bp_r2d(lbi:ubi,lbj:ubj,nbico(ng)) )
1822 dmem(ng)=dmem(ng)+real(nbico(ng),r8)*size2d
1823 fourdvar(ng) % bp_r2d = inival
1824 END IF
1825
1826 IF (.not.allocated(fourdvar(ng) % br_r2d)) THEN
1827 allocate ( fourdvar(ng) % br_r2d(lbi:ubi,lbj:ubj,nbico(ng)) )
1828 dmem(ng)=dmem(ng)+real(nbico(ng),r8)*size2d
1829 fourdvar(ng) % br_r2d = inival
1830 END IF
1831
1832 IF (.not.allocated(fourdvar(ng) % tl_rhs_r2d)) THEN
1833 allocate ( fourdvar(ng) % tl_rhs_r2d(lbi:ubi,lbj:ubj) )
1834 dmem(ng)=dmem(ng)+size2d
1835 fourdvar(ng) % tl_rhs_r2d = inival
1836 END IF
1837
1838 IF (.not.allocated(fourdvar(ng) % tl_r2d_ref)) THEN
1839 allocate ( fourdvar(ng) % tl_r2d_ref(lbi:ubi,lbj:ubj) )
1840 dmem(ng)=dmem(ng)+size2d
1841 fourdvar(ng) % tl_r2d_ref = inival
1842 END IF
1843
1844 IF (.not.allocated(fourdvar(ng) % tl_zeta_ref)) THEN
1845 allocate ( fourdvar(ng) % tl_zeta_ref(lbi:ubi,lbj:ubj) )
1846 dmem(ng)=dmem(ng)+size2d
1847 fourdvar(ng) % tl_zeta_ref = inival
1848 END IF
1849
1850 IF (.not.allocated(fourdvar(ng) % ad_rhs_r2d)) THEN
1851 allocate ( fourdvar(ng) % ad_rhs_r2d(lbi:ubi,lbj:ubj) )
1852 dmem(ng)=dmem(ng)+size2d
1853 fourdvar(ng) % ad_rhs_r2d = inival
1854 END IF
1855
1856 IF (.not.allocated(fourdvar(ng) % ad_r2d_ref)) THEN
1857 allocate ( fourdvar(ng) % ad_r2d_ref(lbi:ubi,lbj:ubj) )
1858 dmem(ng)=dmem(ng)+size2d
1859 fourdvar(ng) % ad_r2d_ref = inival
1860 END IF
1861
1862 IF (.not.allocated(fourdvar(ng) % ad_zeta_ref)) THEN
1863 allocate ( fourdvar(ng) % ad_zeta_ref(lbi:ubi,lbj:ubj) )
1864 dmem(ng)=dmem(ng)+size2d
1865 fourdvar(ng) % ad_zeta_ref = inival
1866 END IF
1867# endif
1868# endif
1869
1870 IF (.not.allocated(fourdvar(ng) % ObsCount)) THEN
1871 allocate ( fourdvar(ng) % ObsCount(0:nobsvar(ng)) )
1872 dmem(ng)=dmem(ng)+real(nobsvar(ng)+1,r8)
1873 fourdvar(ng) % ObsCount = 0
1874 END IF
1875
1876 IF (.not.allocated(fourdvar(ng) % ObsReject)) THEN
1877 allocate ( fourdvar(ng) % ObsReject(0:nobsvar(ng)) )
1878 dmem(ng)=dmem(ng)+real(nobsvar(ng)+1,r8)
1879 fourdvar(ng) % ObsReject = 0
1880 END IF
1881
1882# ifdef RPCG
1883 IF (.not.allocated(fourdvar(ng) % cg_pxout)) THEN
1884 allocate ( fourdvar(ng) % cg_pxout(nouter) )
1885 dmem(ng)=dmem(ng)+real(nouter)
1886 fourdvar(ng) % cg_pxout = inival
1887 END IF
1888
1889 IF (.not.allocated(fourdvar(ng) % cg_pxsave)) THEN
1890 allocate ( fourdvar(ng) % cg_pxsave(ndatum(ng)+1) )
1891 dmem(ng)=dmem(ng)+real(ndatum(ng)+1,r8)
1892 fourdvar(ng) % cg_pxsave = inival
1893 END IF
1894
1895 IF (.not.allocated(fourdvar(ng) % cg_pxsum)) THEN
1896 allocate ( fourdvar(ng) % cg_pxsum(ndatum(ng)+1) )
1897 dmem(ng)=dmem(ng)+real(ndatum(ng)+1,r8)
1898 fourdvar(ng) % cg_pxsum = inival
1899 END IF
1900
1901 IF (.not.allocated(fourdvar(ng) % cg_Dpxsum)) THEN
1902 allocate ( fourdvar(ng) % cg_Dpxsum(ndatum(ng)+1) )
1903 dmem(ng)=dmem(ng)+real(ndatum(ng)+1,r8)
1904 fourdvar(ng) % cg_Dpxsum = inival
1905 END IF
1906
1907 IF (.not.allocated(fourdvar(ng) % tl_cg_pxsave)) THEN
1908 allocate ( fourdvar(ng) % tl_cg_pxsave(ndatum(ng)+1) )
1909 dmem(ng)=dmem(ng)+real(ndatum(ng)+1,r8)
1910 fourdvar(ng) % tl_cg_pxsave = inival
1911 END IF
1912
1913 IF (.not.allocated(fourdvar(ng) % ad_cg_pxsave)) THEN
1914 allocate ( fourdvar(ng) % ad_cg_pxsave(ndatum(ng)+1) )
1915 dmem(ng)=dmem(ng)+real(ndatum(ng)+1,r8)
1916 fourdvar(ng) % ad_cg_pxsave = inival
1917 END IF
1918# endif
1919
1920 END DO
1921!
1922!-----------------------------------------------------------------------
1923! Read in number of observations available per survey at their times.
1924!-----------------------------------------------------------------------
1925!
1926 DO ng=1,ngrids
1927!
1928! Read in number of observations available per survey and the time of
1929! each observation survey.
1930!
1931 SELECT CASE (obs(ng)%IOtype)
1932 CASE (io_nf90)
1933 CALL netcdf_get_ivar (ng, inlm, trim(obs(ng)%name), &
1934 & trim(vname(1,idnobs)), &
1935 & fourdvar(ng)%NobsSurvey, &
1936 & start = (/1/), &
1937 & total = (/nsurvey(ng)/))
1938 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1939!
1940 CALL netcdf_get_fvar (ng, inlm, trim(obs(ng)%name), &
1941 & trim(vname(1,idoday)), &
1942 & fourdvar(ng)%SurveyTime, &
1943 & start = (/1/), &
1944 & total = (/nsurvey(ng)/))
1945 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1946
1947# if defined PIO_LIB && defined DISTRIBUTE
1948 CASE (io_pio)
1949 CALL pio_netcdf_get_ivar (ng, inlm, trim(obs(ng)%name), &
1950 & trim(vname(1,idnobs)), &
1951 & fourdvar(ng)%NobsSurvey, &
1952 & start = (/1/), &
1953 & total = (/nsurvey(ng)/))
1954 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1955!
1956 CALL pio_netcdf_get_fvar (ng, inlm, trim(obs(ng)%name), &
1957 & trim(vname(1,idoday)), &
1958 & fourdvar(ng)%SurveyTime, &
1959 & start = (/1/), &
1960 & total = (/nsurvey(ng)/))
1961 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1962# endif
1963 END SELECT
1964!
1965! Determine maximum size of observation arrays.
1966!
1967# ifdef WEAK_CONSTRAINT
1968 mobs=maxval(ndatum)
1969# else
1970 mobs=0
1971 DO i=1,nsurvey(ng)
1972 mobs=max(mobs, fourdvar(ng)%NobsSurvey(i))
1973 END DO
1974# endif
1975 END DO
1976!
1977!-----------------------------------------------------------------------
1978! Allocate and initialize model and observation variables.
1979!-----------------------------------------------------------------------
1980!
1981# ifdef CURVGRID
1982 IF (.not.allocated(obsangler)) THEN
1983 allocate ( obsangler(mobs) )
1984 dmem(1)=dmem(1)+real(mobs,r8)
1985 END IF
1986 obsangler = inival
1987# endif
1988
1989 IF (.not.allocated(obstype)) THEN
1990 allocate ( obstype(mobs) )
1991 dmem(1)=dmem(1)+real(mobs,r8)
1992 END IF
1993 obstype = 0
1994
1995 IF (.not.allocated(obsstate2type)) THEN
1996 allocate ( obsstate2type(0:maxval(nobsvar)) )
1997 dmem(1)=dmem(1)+real(maxval(nobsvar)+1,r8)
1998 END IF
1999 obsstate2type = 0
2000
2001 IF (.not.allocated(obsprov)) THEN
2002 allocate ( obsprov(mobs) )
2003 dmem(1)=dmem(1)+real(mobs,r8)
2004 END IF
2005 obsprov = 0
2006
2007 IF (.not.allocated(obserr)) THEN
2008 allocate ( obserr(mobs) )
2009 dmem(1)=dmem(1)+real(mobs,r8)
2010 END IF
2011 obserr = inival
2012
2013 IF (.not.allocated(obsmeta)) THEN
2014 allocate ( obsmeta(mobs) )
2015 dmem(1)=dmem(1)+real(mobs,r8)
2016 END IF
2017 obsmeta = inival
2018
2019 IF (.not.allocated(obsscale)) THEN
2020 allocate ( obsscale(mobs) )
2021 dmem(1)=dmem(1)+real(mobs,r8)
2022 END IF
2023 obsscale = inival
2024
2025# ifndef WEAK_CONSTRAINT
2026 IF (.not.allocated(obsscaleglobal)) THEN
2027 allocate ( obsscaleglobal(maxval(ndatum)) )
2028 dmem(1)=dmem(1)+real(maxval(ndatum),r8)
2029 END IF
2030 obsscaleglobal = inival
2031# endif
2032
2033 IF (.not.allocated(obsvetting)) THEN
2034 allocate ( obsvetting(mobs) )
2035 dmem(1)=dmem(1)+real(mobs,r8)
2036 END IF
2037 obsvetting = inival
2038
2039 IF (.not.allocated(obsval)) THEN
2040 allocate ( obsval(mobs) )
2041 dmem(1)=dmem(1)+real(mobs,r8)
2042 END IF
2043 obsval = inival
2044
2045 IF (.not.allocated(tobs)) THEN
2046 allocate ( tobs(mobs) )
2047 dmem(1)=dmem(1)+real(mobs,r8)
2048 END IF
2049 tobs = 0.0_dp
2050
2051 IF (.not.allocated(xobs)) THEN
2052 allocate ( xobs(mobs) )
2053 dmem(1)=dmem(1)+real(mobs,r8)
2054 END IF
2055 xobs = inival
2056
2057 IF (.not.allocated(yobs)) THEN
2058 allocate ( yobs(mobs) )
2059 dmem(1)=dmem(1)+real(mobs,r8)
2060 END IF
2061 yobs = inival
2062
2063# ifdef SOLVE3D
2064 IF (.not.allocated(zobs)) THEN
2065 allocate ( zobs(mobs) )
2066 dmem(1)=dmem(1)+real(mobs,r8)
2067 END IF
2068 zobs = inival
2069# endif
2070
2071 IF (.not.allocated(admodval)) THEN
2072 allocate ( admodval(mobs) )
2073 dmem(1)=dmem(1)+real(mobs,r8)
2074 END IF
2075 admodval = inival
2076
2077 IF (.not.allocated(nlmodval)) THEN
2078 allocate ( nlmodval(mobs) )
2079 dmem(1)=dmem(1)+real(mobs,r8)
2080 END IF
2081 nlmodval = inival
2082
2083 IF (.not.allocated(misfit)) THEN
2084 allocate ( misfit(mobs) )
2085 dmem(1)=dmem(1)+real(mobs,r8)
2086 END IF
2087 misfit = inival
2088
2089 IF (.not.allocated(unvetted)) THEN
2090 allocate ( unvetted(mobs) )
2091 dmem(1)=dmem(1)+real(mobs,r8)
2092 END IF
2093 unvetted = inival
2094
2095 IF (.not.allocated(uradial)) THEN
2096 allocate ( uradial(mobs) )
2097 dmem(1)=dmem(1)+real(mobs,r8)
2098 END IF
2099 uradial = inival
2100
2101 IF (.not.allocated(vradial)) THEN
2102 allocate ( vradial(mobs) )
2103 dmem(1)=dmem(1)+real(mobs,r8)
2104 END IF
2105 vradial = inival
2106
2107# ifndef VERIFICATION
2108 IF (.not.allocated(ubgerr)) THEN
2109 allocate ( ubgerr(mobs) )
2110 dmem(1)=dmem(1)+real(mobs,r8)
2111 END IF
2112 ubgerr = inival
2113
2114 IF (.not.allocated(vbgerr)) THEN
2115 allocate ( vbgerr(mobs) )
2116 dmem(1)=dmem(1)+real(mobs,r8)
2117 END IF
2118 vbgerr = inival
2119# endif
2120
2121# ifdef SP4DVAR
2122 IF (.not.allocated(admodval_s)) THEN
2123 allocate ( admodval_s(mobs, ninner+4) )
2124 dmem(1)=dmem(1)+real(mobs*(ninner+4),r8)
2125 END IF
2126 admodval_s = inival
2127
2128 IF (.not.allocated(harnoldi)) THEN
2129 allocate ( harnoldi(ninner+1, ninner, nouter) )
2130 dmem(1)=dmem(1)+real((ninner+1)*ninner*nouter,r8)
2131 END IF
2132 harnoldi = inival
2133
2134 IF (.not.allocated(gmze)) THEN
2135 allocate ( gmze(ninner, nouter) )
2136 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2137 END IF
2138 gmze = inival
2139
2140 IF (.not.allocated(cg_beta0)) THEN
2141 allocate ( cg_beta0(nouter) )
2142 dmem(1)=dmem(1)+real(nouter,r8)
2143 END IF
2144 cg_beta0 = inival
2145
2146 IF (.not.allocated(jobs)) THEN
2147 allocate ( jobs(0:ninner,1:nouter) )
2148 dmem(1)=dmem(1)+real((ninner+1)*nouter,r8)
2149 END IF
2150 jobs = inival
2151# endif
2152
2153# ifdef FOUR_DVAR
2154 IF (.not.allocated(bgerr)) THEN
2155 allocate ( bgerr(mobs) )
2156 dmem(1)=dmem(1)+real(mobs,r8)
2157 END IF
2158 bgerr = inival
2159
2160 IF (.not.allocated(nlincrement)) THEN
2161 allocate ( nlincrement(mobs) )
2162 dmem(1)=dmem(1)+real(mobs,r8)
2163 END IF
2164 nlincrement = inival
2165
2166 IF (.not.allocated(innovation)) THEN
2167 allocate ( innovation(mobs) )
2168 dmem(1)=dmem(1)+real(mobs,r8)
2169 END IF
2170 innovation = inival
2171
2172 IF (.not.allocated(residual)) THEN
2173 allocate ( residual(mobs) )
2174 dmem(1)=dmem(1)+real(mobs,r8)
2175 END IF
2176 residual = inival
2177# endif
2178
2179# ifdef BGQC
2180 IF (.not.allocated(bgthresh)) THEN
2181 allocate ( bgthresh(mobs) )
2182 dmem(1)=dmem(1)+real(mobs,r8)
2183 END IF
2184 bgthresh = inival
2185
2186 IF (.not.allocated(jact_s)) THEN
2187 allocate ( jact_s(0:ninner) )
2188 dmem(1)=dmem(1)+real(ninner+1,r8)
2189 END IF
2190 jact_s = inival
2191# endif
2192
2193# ifdef TLM_OBS
2194 IF (.not.allocated(tlmodval)) THEN
2195 allocate ( tlmodval(mobs) )
2196 dmem(1)=dmem(1)+real(mobs,r8)
2197 END IF
2198 tlmodval = inival
2199
2200# if defined ARRAY_MODES || defined CLIPPING || \
2201 defined rbl4dvar || defined r4dvar || \
2202 defined sensitivity_4dvar || defined sp4dvar || \
2203 defined tl_rbl4dvar || defined tl_r4dvar
2204!
2205! The following arrays are only needed and used by the master node.
2206! In order to avoid hogging memory penalties in the other nodes, they
2207! are only allocated by the master node.
2208!
2209 IF (.not.allocated(tlmodval_s)) THEN
2210 allocate ( tlmodval_s(mobs,ninner,nouter) )
2211 dmem(1)=dmem(1)+real(mobs*ninner*nouter,r8)
2212 END IF
2213 tlmodval_s = inival
2214
2215# if defined RPCG && \
2216 (defined sensitivity_4dvar || defined tl_rbl4dvar || \
2217 defined tl_r4dvar)
2218 IF (master) THEN
2219 IF (.not.allocated(tl_tlmodval_s)) THEN
2220 allocate ( tl_tlmodval_s(mobs,ninner,nouter) )
2221 dmem(1)=dmem(1)+real(mobs*ninner*nouter,r8)
2222 END IF
2223 tl_tlmodval_s = inival
2224
2225 IF (.not.allocated(ad_tlmodval_s)) THEN
2226 allocate ( ad_tlmodval_s(mobs,ninner,nouter) )
2227 dmem(1)=dmem(1)+real(mobs*ninner*nouter,r8)
2228 END IF
2229 ad_tlmodval_s = inival
2230 END IF
2231# endif
2232# endif
2233
2234# if defined ARRAY_MODES || defined RBL4DVAR || \
2235 defined sensitivity_4dvar || defined tl_rbl4dvar
2236
2237 IF (.not.allocated(bckmodval)) THEN
2238 allocate ( bckmodval(mobs) )
2239 dmem(1)=dmem(1)+real(mobs,r8)
2240 END IF
2241 bckmodval = inival
2242
2243# ifdef RPCG
2244
2245 IF (.not.allocated(hbk)) THEN
2246 allocate ( hbk(mobs,nouter) )
2247 dmem(1)=dmem(1)+real(mobs*nouter,r8)
2248 END IF
2249 hbk = inival
2250
2251# if defined SENSITIVITY_4DVAR || defined TL_RBL4DVAR || \
2252 defined tl_r4dvar
2253
2254 IF (.not.allocated(tl_bckmodval)) THEN
2255 allocate ( tl_bckmodval(mobs) )
2256 dmem(1)=dmem(1)+real(mobs,r8)
2257 END IF
2258 tl_bckmodval = inival
2259
2260 IF (.not.allocated(tl_hbk)) THEN
2261 allocate ( tl_hbk(mobs) )
2262 dmem(1)=dmem(1)+real(mobs,r8)
2263 END IF
2264 tl_hbk = inival
2265
2266 IF (.not.allocated(ad_bckmodval)) THEN
2267 allocate ( ad_bckmodval(mobs) )
2268 dmem(1)=dmem(1)+real(mobs,r8)
2269 END IF
2270 ad_bckmodval = inival
2271
2272 IF (.not.allocated(ad_hbk)) THEN
2273 allocate ( ad_hbk(mobs) )
2274 dmem(1)=dmem(1)+real(mobs,r8)
2275 END IF
2276 ad_hbk = inival
2277
2278# endif
2279# endif
2280# endif
2281# endif
2282!
2283! Allocate and initialize observation types names and indices.
2284! Notice that a mapping from state-to-type (ObsState2Type) and its
2285! inverse type-to-state (ObsType2State) indices are needed because
2286! the User is allowed to add extra-observation operators with
2287! nonsequential type enumeration.
2288!
2289! Both mapping arrays ObsState2Type and ObsType2State have a zero
2290! array element to allow applications with no extra-observations
2291! to work with their zero associated state index (as initialized
2292! in mod_ncparam.F). For example, if the index "isRadial" is not
2293! redefined below, the following assignment
2294!
2295! ObsState2Type(isRadial)=ObsState2Type(0)=0
2296! ObsType2State(isRadial)=ObsType2State(0)=0
2297!
2298! is still legal with isRadial=0. It avoids a Fortran segmentation
2299! violation (i.e., subscript #1 of the array ObsState2Type has
2300! value 0 which is less than the lower bound of 1). Sorry for
2301! the awkward logic but we need a generic way to specify extra-
2302! observation operators.
2303!
2304 IF (.not.allocated(obsname)) THEN
2305 allocate ( obsname(maxval(nobsvar)) )
2306 dmem(1)=dmem(1)+real(maxval(nobsvar),r8)
2307 END IF
2308 icount=maxval(nstatevar)
2309 obsstate2type(0)=0
2310 DO i=1,icount ! 5+NT
2311 obsstate2type(i)=i
2312 obsname(i)=trim(vname(1,idsvar(i)))
2313 END DO
2314 IF (nextraobs.gt.0) THEN
2315 DO i=1,nextraobs
2316 icount=icount+1
2317 obsname(icount)=trim(extraname(i))
2318 SELECT CASE (trim(uppercase(extraname(i))))
2319 CASE ('RADIAL')
2320 isradial=icount
2321 obsstate2type(icount)=extraindex(i)
2322 END SELECT
2323 END DO
2324 END IF
2325!
2326! Set inverse mapping type-to-state.
2327!
2328 IF (.not.allocated(obstype2state)) THEN
2329 allocate ( obstype2state(0:maxval(obsstate2type)) )
2330 dmem(1)=dmem(1)+real(maxval(obsstate2type)+1,r8)
2331 END IF
2332 obstype2state(0)=0
2333 obstype2state(1:maxval(obsstate2type))=-1
2334 DO i=1,maxval(nstatevar)
2335 obstype2state(i)=i
2336 END DO
2337 IF (nextraobs.gt.0) THEN
2338 DO i=1,nextraobs
2339 icount=extraindex(i)
2340 SELECT CASE (trim(uppercase(extraname(i))))
2341 CASE ('RADIAL')
2342 obstype2state(icount)=isradial
2343 END SELECT
2344 END DO
2345 END IF
2346
2347# ifdef WEAK_CONSTRAINT
2348!
2349! Allocate and initialize weak constraint conjugate gradient vectors.
2350!
2351# if defined RPCG
2352# define MyMobs Mobs+1
2353# else
2354# define MyMobs Mobs
2355# endif
2356
2357 IF (.not.allocated(zcglwk)) THEN
2358 allocate ( zcglwk(mymobs,ninner+1,nouter) )
2359 dmem(1)=dmem(1)+real((mymobs)*(ninner+1)*nouter,r8)
2360 END IF
2361 zcglwk = inival
2362
2363# ifdef RPCG
2364 IF (.not.allocated(vcglwk)) THEN
2365 allocate ( vcglwk(mymobs,ninner+1,nouter) )
2366 dmem(1)=dmem(1)+real((mymobs)*(ninner+1)*nouter,r8)
2367 END IF
2368 vcglwk = inival
2369
2370 IF (.not.allocated(jb0)) THEN
2371 allocate ( jb0(0:nouter) )
2372 dmem(1)=dmem(1)+real(nouter+1,r8)
2373 END IF
2374 jb0 = inival
2375# endif
2376
2377 IF (.not.allocated(vcglev)) THEN
2378 allocate ( vcglev(mymobs,ninner,nouter) )
2379 dmem(1)=dmem(1)+real((mymobs)*ninner*nouter,r8)
2380 END IF
2381 vcglev = inival
2382
2383 IF (.not.allocated(zgrad0)) THEN
2384 allocate ( zgrad0(mymobs,nouter) )
2385 dmem(1)=dmem(1)+real((mymobs)*nouter,r8)
2386 END IF
2387 zgrad0 = inival
2388
2389 IF (.not.allocated(vgrad0)) THEN
2390 allocate ( vgrad0(mymobs) )
2391 dmem(1)=dmem(1)+real(mymobs,r8)
2392 END IF
2393 vgrad0 = inival
2394
2395 IF (.not.allocated(vgrad0s)) THEN
2396 allocate ( vgrad0s(mymobs) )
2397 dmem(1)=dmem(1)+real(mymobs,r8)
2398 END IF
2399 vgrad0s = inival
2400
2401 IF (.not.allocated(gdgw)) THEN
2402 allocate ( gdgw(mobs) )
2403 dmem(1)=dmem(1)+real(mobs,r8)
2404 END IF
2405 gdgw = inival
2406
2407 IF (.not.allocated(vgnorm)) THEN
2408 allocate ( vgnorm(nouter) )
2409 dmem(1)=dmem(1)+real(nouter,r8)
2410 END IF
2411 vgnorm = inival
2412
2413 IF (.not.allocated(cg_innov)) THEN
2414 allocate ( cg_innov(mymobs) )
2415 dmem(1)=dmem(1)+real(mymobs,r8)
2416 END IF
2417 cg_innov = inival
2418
2419 IF (.not.allocated(cg_dla)) THEN
2420 allocate ( cg_dla(ninner,nouter) )
2421 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2422 END IF
2423 cg_dla = inival
2424
2425# ifndef RPCG
2426
2427 IF (.not.allocated(cg_pxsave)) THEN
2428 allocate ( cg_pxsave(mobs) )
2429 dmem(1)=dmem(1)+real(mobs,r8)
2430 END IF
2431 cg_pxsave = inival
2432
2433# endif
2434
2435# if defined SENSITIVITY_4DVAR || defined TL_RBL4DVAR || \
2436 defined tl_r4dvar
2437
2438 IF (.not.allocated(ad_zcglwk)) THEN
2439 allocate ( ad_zcglwk(mymobs,ninner+1) )
2440 dmem(1)=dmem(1)+real((mymobs)*(ninner+1),r8)
2441 END IF
2442 ad_zcglwk = inival
2443
2444 IF (.not.allocated(ad_zgrad0)) THEN
2445 allocate ( ad_zgrad0(mymobs) )
2446 dmem(1)=dmem(1)+real(mymobs,r8)
2447 END IF
2448 ad_zgrad0 = inival
2449
2450# if defined RPCG && \
2451 (defined sensitivity_4dvar || defined tl_rbl4dvar || \
2452 defined tl_r4dvar)
2453
2454 IF (.not.allocated(ad_vcglwk)) THEN
2455 allocate ( ad_vcglwk(mymobs,ninner+1) )
2456 dmem(1)=dmem(1)+real((mymobs)*(ninner+1),r8)
2457 END IF
2458 ad_vcglwk = inival
2459
2460 IF (.not.allocated(ad_vgrad0)) THEN
2461 allocate ( ad_vgrad0(mymobs) )
2462 dmem(1)=dmem(1)+real(mymobs,r8)
2463 END IF
2464 ad_vgrad0 = inival
2465
2466 IF (.not.allocated(ad_vgrad0s)) THEN
2467 allocate ( ad_vgrad0s(mymobs) )
2468 dmem(1)=dmem(1)+real(mymobs,r8)
2469 END IF
2470 ad_vgrad0s = inival
2471
2472 IF (.not.allocated(ad_jb0)) THEN
2473 allocate ( ad_jb0(0:nouter) )
2474 dmem(1)=dmem(1)+real(nouter+1,r8)
2475 END IF
2476 ad_jb0 = inival
2477
2478 IF (.not.allocated(ad_obserr)) THEN
2479 allocate ( ad_obserr(mobs) )
2480 dmem(1)=dmem(1)+real(mobs,r8)
2481 END IF
2482 ad_obserr = inival
2483
2484 IF (.not.allocated(ad_obsscale)) THEN
2485 allocate ( ad_obsscale(mobs) )
2486 dmem(1)=dmem(1)+real(mobs,r8)
2487 END IF
2488 ad_obsscale = inival
2489
2490# endif
2491
2492 IF (.not.allocated(ad_cg_innov)) THEN
2493 allocate ( ad_cg_innov(mobs) )
2494 dmem(1)=dmem(1)+real(mobs,r8)
2495 END IF
2496 ad_cg_innov = inival
2497
2498 IF (.not.allocated(ad_cg_pxsave)) THEN
2499 allocate ( ad_cg_pxsave(mobs) )
2500 dmem(1)=dmem(1)+real(mobs,r8)
2501 END IF
2502 ad_cg_pxsave = inival
2503
2504 IF (.not.allocated(ad_cg_pxout)) THEN
2505 allocate ( ad_cg_pxout(nouter) )
2506 dmem(1)=dmem(1)+real(nouter,r8)
2507 END IF
2508 ad_cg_pxout = inival
2509
2510# ifdef IMPACT_INNER
2511
2512 IF (.not.allocated(ad_obsval)) THEN
2513 allocate ( ad_obsval(mobs,ninner) )
2514 dmem(1)=dmem(1)+real(mobs*ninner,r8)
2515 END IF
2516 ad_obsval = inival
2517
2518# else
2519
2520 IF (.not.allocated(ad_obsval)) THEN
2521 allocate ( ad_obsval(mobs) )
2522 dmem(1)=dmem(1)+real(mobs,r8)
2523 END IF
2524 ad_obsval = inival
2525
2526# endif
2527
2528 IF (.not.allocated(tl_zcglwk)) THEN
2529 allocate ( tl_zcglwk(mymobs,ninner+1) )
2530 dmem(1)=dmem(1)+real((mymobs)*(ninner+1),r8)
2531 END IF
2532 tl_zcglwk = inival
2533
2534 IF (.not.allocated(tl_zgrad0)) THEN
2535 allocate ( tl_zgrad0(mymobs) )
2536 dmem(1)=dmem(1)+real(mymobs,r8)
2537 END IF
2538 tl_zgrad0 = inival
2539
2540# ifdef RPCG
2541
2542 IF (.not.allocated(tl_vcglwk)) THEN
2543 allocate ( tl_vcglwk(mymobs,ninner+1) )
2544 dmem(1)=dmem(1)+real((mymobs)*(ninner+1),r8)
2545 END IF
2546 tl_vcglwk = inival
2547
2548 IF (.not.allocated(tl_vgrad0)) THEN
2549 allocate ( tl_vgrad0(mymobs) )
2550 dmem(1)=dmem(1)+real(mymobs,r8)
2551 END IF
2552 tl_vgrad0 = inival
2553
2554 IF (.not.allocated(tl_vgrad0s)) THEN
2555 allocate ( tl_vgrad0s(mymobs) )
2556 dmem(1)=dmem(1)+real(mymobs,r8)
2557 END IF
2558 tl_vgrad0s = inival
2559
2560 IF (.not.allocated(tl_jb0)) THEN
2561 allocate ( tl_jb0(0:nouter) )
2562 dmem(1)=dmem(1)+real(nouter+1,r8)
2563 END IF
2564 tl_jb0 = inival
2565
2566 IF (.not.allocated(tl_obserr)) THEN
2567 allocate ( tl_obserr(mobs) )
2568 dmem(1)=dmem(1)+real(mobs,r8)
2569 END IF
2570 tl_obserr = inival
2571
2572 IF (.not.allocated(tl_obsscale)) THEN
2573 allocate ( tl_obsscale(mobs) )
2574 dmem(1)=dmem(1)+real(mobs,r8)
2575 END IF
2576 tl_obsscale = inival
2577
2578# endif
2579
2580 IF (.not.allocated(tl_cg_innov)) THEN
2581 allocate ( tl_cg_innov(mobs) )
2582 dmem(1)=dmem(1)+real(mobs,r8)
2583 END IF
2584 tl_cg_innov = inival
2585
2586 IF (.not.allocated(tl_cg_pxsave)) THEN
2587 allocate ( tl_cg_pxsave(mobs) )
2588 dmem(1)=dmem(1)+real(mobs,r8)
2589 END IF
2590 tl_cg_pxsave = inival
2591
2592 IF (.not.allocated(tl_cg_pxout)) THEN
2593 allocate ( tl_cg_pxout(nouter) )
2594 dmem(1)=dmem(1)+real(nouter,r8)
2595 END IF
2596 tl_cg_pxout = inival
2597
2598 IF (.not.allocated(tl_obsval)) THEN
2599 allocate ( tl_obsval(mobs) )
2600 dmem(1)=dmem(1)+real(mobs,r8)
2601 END IF
2602 tl_obsval = inival
2603
2604# endif
2605# endif
2606
2607# else
2608!
2609!-----------------------------------------------------------------------
2610! Allocate other module variables.
2611!-----------------------------------------------------------------------
2612!
2613! Set number of state variables.
2614!
2615 DO ng=1,ngrids
2616# ifdef SOLVE3D
2617 nstatevar(ng)=5+nt(ng)
2618# ifdef ADJUST_STFLUX
2619 nstatevar(ng)=nstatevar(ng)+nt(ng)
2620# endif
2621# else
2622 nstatevar(ng)=3
2623# endif
2624# ifdef ADJUST_WSTRESS
2625 nstatevar(ng)=nstatevar(ng)+2
2626# endif
2627 END DO
2628# endif
2629!
2630! Allocate convolution parameters.
2631!
2632 IF (.not.allocated(nhsteps)) THEN
2633 allocate ( nhsteps(2,mstatevar) )
2634 dmem(1)=dmem(1)+2.0_r8*real(mstatevar,r8)
2635 END IF
2636 nhsteps = 0
2637
2638 IF (.not.allocated(nvsteps)) THEN
2639 allocate ( nvsteps(2,mstatevar) )
2640 dmem(1)=dmem(1)+2.0_r8*real(mstatevar,r8)
2641 END IF
2642 nvsteps = 0
2643
2644 IF (.not.allocated(dtsizeh)) THEN
2645 allocate ( dtsizeh(2,mstatevar) )
2646 dmem(1)=dmem(1)+2.0_r8*real(mstatevar,r8)
2647 END IF
2648 dtsizeh = inival
2649
2650 IF (.not.allocated(dtsizev)) THEN
2651 allocate ( dtsizev(2,mstatevar) )
2652 dmem(1)=dmem(1)+2.0_r8*real(mstatevar,r8)
2653 END IF
2654 dtsizev = inival
2655
2656 IF (.not.allocated(nvstepsb)) THEN
2657 allocate ( nvstepsb(4,mstatevar) )
2658 dmem(1)=dmem(1)+4.0_r8*real(mstatevar,r8)
2659 END IF
2660 nvstepsb = 0
2661
2662 IF (.not.allocated(nhstepsb)) THEN
2663 allocate ( nhstepsb(4,mstatevar) )
2664 dmem(1)=dmem(1)+4.0_r8*real(mstatevar,r8)
2665 END IF
2666 nhstepsb = 0
2667
2668 IF (.not.allocated(dtsizehb)) THEN
2669 allocate ( dtsizehb(4,mstatevar) )
2670 dmem(1)=dmem(1)+4.0_r8*real(mstatevar,r8)
2671 END IF
2672 dtsizehb = inival
2673
2674 IF (.not.allocated(dtsizevb)) THEN
2675 allocate ( dtsizevb(4,mstatevar) )
2676 dmem(1)=dmem(1)+4.0_r8*real(mstatevar,r8)
2677 END IF
2678 dtsizevb = inival
2679
2680
2681# if defined I4DVAR || defined WEAK_CONSTRAINT
2682!
2683! Allocate conjugate gradient variables.
2684!
2685# if defined SENSITIVITY_4DVAR || defined TL_RBL4DVAR || \
2686 defined tl_r4dvar
2687 IF (.not.allocated(ad_cg_beta)) THEN
2688 allocate ( ad_cg_beta(ninner+1) )
2689 dmem(1)=dmem(1)+real(ninner+1,r8)
2690 END IF
2691 ad_cg_beta = inival
2692
2693 IF (.not.allocated(ad_cg_delta)) THEN
2694 allocate ( ad_cg_delta(ninner) )
2695 dmem(1)=dmem(1)+real(ninner,r8)
2696 END IF
2697 ad_cg_delta = inival
2698
2699 IF (.not.allocated(ad_cg_qg)) THEN
2700 allocate ( ad_cg_qg(ninner+1) )
2701 dmem(1)=dmem(1)+real(ninner+1,r8)
2702 END IF
2703 ad_cg_qg = inival
2704
2705 IF (.not.allocated(ad_cg_gnorm)) THEN
2706 allocate ( ad_cg_gnorm(nouter) )
2707 dmem(1)=dmem(1)+real(nouter,r8)
2708 END IF
2709 ad_cg_gnorm = inival
2710
2711 IF (.not.allocated(tl_cg_beta)) THEN
2712 allocate ( tl_cg_beta(ninner+1) )
2713 dmem(1)=dmem(1)+real(ninner+1,r8)
2714 END IF
2715 tl_cg_beta = inival
2716
2717 IF (.not.allocated(tl_cg_delta)) THEN
2718 allocate ( tl_cg_delta(ninner) )
2719 dmem(1)=dmem(1)+real(ninner,r8)
2720 END IF
2721 tl_cg_delta = inival
2722
2723 IF (.not.allocated(tl_cg_qg)) THEN
2724 allocate ( tl_cg_qg(ninner+1) )
2725 dmem(1)=dmem(1)+real(ninner+1,r8)
2726 END IF
2727 tl_cg_qg = inival
2728
2729# ifdef RPCG
2730 IF (.not.allocated(tl_cg_gnorm_v)) THEN
2731 allocate ( tl_cg_gnorm_v(nouter) )
2732 dmem(1)=dmem(1)+real(nouter,r8)
2733 END IF
2734 tl_cg_gnorm_v = inival
2735
2736 IF (.not.allocated(ad_cg_gnorm_v)) THEN
2737 allocate ( ad_cg_gnorm_v(nouter) )
2738 dmem(1)=dmem(1)+real(nouter,r8)
2739 END IF
2740 ad_cg_gnorm_v = inival
2741
2742# else
2743 IF (.not.allocated(tl_cg_gnorm)) THEN
2744 allocate ( tl_cg_gnorm(nouter) )
2745 dmem(1)=dmem(1)+real(nouter,r8)
2746 END IF
2747 tl_cg_gnorm = inival
2748
2749# endif
2750# endif
2751
2752 IF (.not.allocated(cg_alpha)) THEN
2753 allocate ( cg_alpha(0:ninner,nouter) )
2754 dmem(1)=dmem(1)+real((ninner+1)*nouter,r8)
2755 END IF
2756 cg_alpha = inival
2757
2758 IF (.not.allocated(cg_beta)) THEN
2759 allocate ( cg_beta(ninner+1,nouter) )
2760 dmem(1)=dmem(1)+real((ninner+1)*nouter,r8)
2761 END IF
2762 cg_beta = inival
2763
2764 IF (.not.allocated(cg_tau)) THEN
2765 allocate ( cg_tau(0:ninner,nouter) )
2766 dmem(1)=dmem(1)+real((ninner+1)*nouter,r8)
2767 END IF
2768 cg_tau = inival
2769
2770 IF (.not.allocated(ritz)) THEN
2771 allocate ( ritz(ninner+1) )
2772 dmem(1)=dmem(1)+real(ninner+1,r8)
2773 END IF
2774 ritz = inival
2775
2776# if defined POSTERIOR_EOFS || defined POSTERIOR_ERROR_I || \
2777 defined posterior_error_f
2778
2779 IF (.not.allocated(ae_beta)) THEN
2780 allocate ( ae_beta(nposti+1,nouter) )
2781 dmem(1)=dmem(1)+real((nposti+1)*nouter,r8)
2782 END IF
2783 ae_beta = inival
2784
2785 IF (.not.allocated(ae_zv)) THEN
2786 allocate ( ae_zv(nposti,nposti) )
2787 dmem(1)=dmem(1)+real(nposti*nposti,r8)
2788 END IF
2789 ae_zv = inival
2790
2791 IF (.not.allocated(ae_delta)) THEN
2792 allocate ( ae_delta(nposti,nouter) )
2793 dmem(1)=dmem(1)+real(nposti*nouter,r8)
2794 END IF
2795 ae_delta = inival
2796
2797 IF (.not.allocated(ae_trace)) THEN
2798 allocate ( ae_trace(nposti+1) )
2799 dmem(1)=dmem(1)+real(nposti+1,r8)
2800 END IF
2801 ae_trace = inival
2802
2803 IF (.not.allocated(ae_gnorm)) THEN
2804 allocate ( ae_gnorm(nouter) )
2805 dmem(1)=dmem(1)+real(nouter,r8)
2806 END IF
2807 ae_gnorm = inival
2808
2809 IF (.not.allocated(ae_tmatrix)) THEN
2810 allocate ( ae_tmatrix(nposti,3) )
2811 dmem(1)=dmem(1)+3.0_r8*real(nposti,r8)
2812 END IF
2813 ae_tmatrix = inival
2814
2815 IF (.not.allocated(ae_ritz)) THEN
2816 allocate ( ae_ritz(nposti,nouter) )
2817 dmem(1)=dmem(1)+real(nposti*nouter,r8)
2818 END IF
2819 ae_ritz = inival
2820
2821 IF (.not.allocated(ae_ritzerr)) THEN
2822 allocate ( ae_ritzerr(nposti,nouter) )
2823 dmem(1)=dmem(1)+real(nposti*nouter,r8)
2824 END IF
2825 ae_ritzerr = inival
2826
2827# if defined POSTERIOR_ERROR_I || defined POSTERIOR_ERROR_F
2828 IF (.not.allocated(zlanczos_coef)) THEN
2829 allocate ( zlanczos_coef(ninner,ninner) )
2830 dmem(1)=dmem(1)+real(ninner*ninner,r8)
2831 END IF
2832 zlanczos_coef = inival
2833
2834 IF (.not.allocated(zlanczos_inv)) THEN
2835 allocate ( zlanczos_inv(ninner,ninner) )
2836 dmem(1)=dmem(1)+real(ninner*ninner,r8)
2837 END IF
2838 zlanczos_inv = inival
2839
2840 IF (.not.allocated(zlanczos_err)) THEN
2841 allocate ( zlanczos_err(ninner,ninner) )
2842 dmem(1)=dmem(1)+real(ninner*ninner,r8)
2843 END IF
2844 zlanczos_err = inival
2845# endif
2846# endif
2847
2848# ifdef WEAK_CONSTRAINT
2849 IF (.not.allocated(cg_zv)) THEN
2850 allocate ( cg_zv(ninner,ninner,nouter) )
2851 dmem(1)=dmem(1)+real(ninner*ninner*nouter,r8)
2852 END IF
2853# else
2854 IF (.not.allocated(cg_zv)) THEN
2855 allocate ( cg_zv(ninner,ninner) )
2856 dmem(1)=dmem(1)+real(ninner*ninner,r8)
2857 END IF
2858# endif
2859 cg_zv = inival
2860
2861 IF (.not.allocated(cg_delta)) THEN
2862 allocate ( cg_delta(ninner,nouter) )
2863 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2864 END IF
2865 cg_delta = inival
2866
2867 IF (.not.allocated(cg_gamma)) THEN
2868 allocate ( cg_gamma(ninner,nouter) )
2869 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2870 END IF
2871 cg_gamma = inival
2872
2873 IF (.not.allocated(cg_gnorm)) THEN
2874 allocate ( cg_gnorm(nouter) )
2875 dmem(1)=dmem(1)+real(nouter,r8)
2876 END IF
2877 cg_gnorm = inival
2878
2879 IF (.not.allocated(cg_gnorm_v)) THEN
2880 allocate ( cg_gnorm_v(nouter) )
2881 dmem(1)=dmem(1)+real(nouter,r8)
2882 END IF
2883 cg_gnorm_v = inival
2884
2885 IF (.not.allocated(cg_gnorm_y)) THEN
2886 allocate ( cg_gnorm_y(nouter) )
2887 dmem(1)=dmem(1)+real(nouter,r8)
2888 END IF
2889 cg_gnorm_y = inival
2890
2891 IF (.not.allocated(cg_greduc)) THEN
2892 allocate ( cg_greduc(ninner,nouter) )
2893 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2894 END IF
2895 cg_greduc = inival
2896
2897 IF (.not.allocated(cg_qg)) THEN
2898 allocate ( cg_qg(ninner+1,nouter) )
2899 dmem(1)=dmem(1)+real((ninner+1)*nouter,r8)
2900 END IF
2901 cg_qg = inival
2902
2903 IF (.not.allocated(cg_tmatrix)) THEN
2904 allocate ( cg_tmatrix(ninner,3) )
2905 dmem(1)=dmem(1)+3.0_r8*real(ninner,r8)
2906 END IF
2907 cg_tmatrix = inival
2908
2909 IF (.not.allocated(cg_ritz)) THEN
2910 allocate ( cg_ritz(ninner,nouter) )
2911 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2912 END IF
2913 cg_ritz = inival
2914
2915 IF (.not.allocated(cg_ritzerr)) THEN
2916 allocate ( cg_ritzerr(ninner,nouter) )
2917 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2918 END IF
2919 cg_ritzerr = inival
2920
2921 IF (.not.allocated(cg_zu)) THEN
2922 allocate ( cg_zu(ninner,nouter) )
2923 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2924 END IF
2925 cg_zu = inival
2926
2927# endif
2928# if defined I4DVAR_ANA_SENSITIVITY
2929!
2930! Allocate Lanczos algorithm coefficients.
2931!
2932 IF (.not.allocated(cg_beta)) THEN
2933 allocate ( cg_beta(ninner+1,nouter) )
2934 dmem(1)=dmem(1)+real((ninner+1)*nouter,r8)
2935 END IF
2936 cg_beta = inival
2937
2938 IF (.not.allocated(cg_delta)) THEN
2939 allocate ( cg_delta(ninner,nouter) )
2940 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2941 END IF
2942 cg_delta = inival
2943
2944 IF (.not.allocated(cg_zu)) THEN
2945 allocate ( cg_zu(ninner,nouter) )
2946 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2947 END IF
2948 cg_zu = inival
2949# endif
2950# if defined HESSIAN_SV || defined HESSIAN_SO || defined HESSIAN_FSV
2951!
2952! Allocate Lanczos algorithm coefficients.
2953!
2954 IF (.not.allocated(cg_beta)) THEN
2955 allocate ( cg_beta(ninner+1,nouter) )
2956 dmem(1)=dmem(1)+real((ninner+1)*nouter,r8)
2957 END IF
2958 cg_beta = inival
2959
2960 IF (.not.allocated(cg_delta)) THEN
2961 allocate ( cg_delta(ninner,nouter) )
2962 dmem(1)=dmem(1)+real(ninner*nouter,r8)
2963 END IF
2964 cg_delta = inival
2965
2966 IF (.not.allocated(zlanczos_diag)) THEN
2967 allocate ( zlanczos_diag(ninner) )
2968 dmem(1)=dmem(1)+real(ninner,r8)
2969 END IF
2970 zlanczos_diag = inival
2971
2972 IF (.not.allocated(zlanczos_offdiag)) THEN
2973 allocate ( zlanczos_offdiag(ninner) )
2974 dmem(1)=dmem(1)+real(ninner,r8)
2975 END IF
2976 zlanczos_offdiag = inival
2977# ifdef LCZ_FINAL
2978
2979 IF (.not.allocated(gsmatrix)) THEN
2980 allocate ( gsmatrix(ninner,ninner) )
2981 dmem(1)=dmem(1)+real(ninner*ninner,r8)
2982 END IF
2983 gsmatrix = inival
2984
2985 IF (.not.allocated(gsmatinv)) THEN
2986 allocate ( gsmatinv(ninner,ninner) )
2987 dmem(1)=dmem(1)+real(ninner*ninner,r8)
2988 END IF
2989 gsmatinv = inival
2990# endif
2991# endif
2992!
2993!-----------------------------------------------------------------------
2994! Initialize various variables.
2995!-----------------------------------------------------------------------
2996!
2997 DO ng=1,ngrids
2998 load_zobs(ng)=.false.
2999 processobs(ng)=.false.
3000# ifdef SP4DVAR
3001 lsadd(ng)=.false.
3002# endif
3003 haveobsmeta(ng)=.false.
3004 wrote_zobs(ng)=.false.
3005 wrtmisfit(ng)=.false.
3006 wrtnlmod(ng)=.false.
3007 wrtobsscale(ng)=.false.
3008 wrtrpmod(ng)=.false.
3009 wrttlmod(ng)=.false.
3010# ifdef BALANCE_OPERATOR
3011 wrtzetaref(ng)=.false.
3012# endif
3013# ifdef I4DVAR_ANA_SENSITIVITY
3014# ifdef OBS_IMPACT
3015 wrtimpact_tot(ng)=.false.
3016# endif
3017# ifdef OBS_IMPACT_SPLIT
3018 wrtimpact_ic(ng)=.false.
3019# if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
3020 wrtimpact_fc(ng)=.false.
3021# endif
3022# if defined ADJUST_BOUNDARY
3023 wrtimpact_bc(ng)=.false.
3024# endif
3025# endif
3026# endif
3027 khmin(ng)=1.0_r8
3028 khmax(ng)=1.0_r8
3029 kvmin(ng)=1.0_r8
3030 kvmax(ng)=1.0_r8
3031 optimality(ng)=0.0_r8
3032# ifdef WEAK_CONSTRAINT
3033 forcetime(ng)=0.0_r8
3034# endif
3035 END DO
3036
3037# ifdef OBSERVATIONS
3038!
3039 10 FORMAT (/,' INITIALIZE_FOURDVAR - Illegal output type,', &
3040 & ' io_type = ',i0)
3041 20 FORMAT (/,' INITIALIZE_FOURDVAR - error inquiring dimension:', &
3042 & 1x,a,2x,'in input NetCDF file: ',a)
3043# endif
3044
3045 RETURN
type(t_io), dimension(:), allocatable obs
integer stdout
character(len=256) sourcefile
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer, parameter io_pio
Definition mod_ncparam.F:96
character(len=maxlen), dimension(6, 0:nv) vname
integer, dimension(:), allocatable idsvar
integer idoday
integer isradial
integer idnobs
character(len=100), dimension(mdims) dim_name
Definition mod_netcdf.F:168
integer, dimension(mdims) dim_size
Definition mod_netcdf.F:159
integer n_dim
Definition mod_netcdf.F:151
subroutine, public netcdf_get_dim(ng, model, ncname, ncid, dimname, dimsize, dimid)
Definition mod_netcdf.F:330
logical master
subroutine, public pio_netcdf_get_dim(ng, model, ncname, piofile, dimname, dimsize, dimid)
integer ninner
integer nouter
integer mstatevar
integer exit_flag
integer noerror
character(len(sinp)) function, public uppercase(sinp)
Definition strings.F:582
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52

References ad_cg_beta, ad_cg_delta, ad_cg_gnorm, ad_cg_gnorm_v, ad_cg_innov, ad_cg_pxout, ad_cg_pxsave, ad_cg_qg, ad_jb0, ad_obserr, ad_obsscale, ad_obsval, ad_vcglwk, ad_vgrad0, ad_vgrad0s, ad_zcglwk, ad_zgrad0, admodval, admodval_s, ae_beta, ae_delta, ae_gnorm, ae_ritz, ae_ritzerr, ae_tmatrix, ae_trace, ae_zv, bgerr, bgthresh, mod_param::bounds, cg_alpha, cg_beta, cg_beta0, cg_delta, cg_dla, cg_gamma, cg_gnorm, cg_gnorm_v, cg_gnorm_y, cg_greduc, cg_innov, cg_qg, cg_ritz, cg_ritzerr, cg_tau, cg_tmatrix, cg_zu, cg_zv, mod_netcdf::dim_name, mod_netcdf::dim_size, mod_param::dmem, dtsizeh, dtsizehb, dtsizev, dtsizevb, mod_scalars::exit_flag, extraindex, extraname, forcetime, strings_mod::founderror(), fourdvar, gdgw, gmze, gsmatinv, gsmatrix, harnoldi, haveobsmeta, mod_ncparam::idnobs, mod_ncparam::idoday, mod_ncparam::idsvar, mod_param::inlm, innovation, mod_ncparam::io_nf90, mod_ncparam::io_pio, mod_ncparam::isradial, jact_s, jb0, jobs, khmax, khmin, kvmax, kvmin, load_zobs, lsadd, mod_parallel::master, misfit, mobs, mod_scalars::mstatevar, mod_parallel::myrank, mod_netcdf::n_dim, mod_param::nbico, ndatum, mod_netcdf::netcdf_get_dim(), nextraobs, mod_param::ngrids, nhsteps, nhstepsb, mod_scalars::ninner, nlincrement, nlmodval, nobsvar, mod_scalars::noerror, mod_scalars::nouter, nposti, nstatevar, nsurvey, mod_param::nt, nvsteps, nvstepsb, mod_iounits::obs, obsangler, obserr, obsmeta, obsname, obsprov, obsscale, obsstate2type, obstype, obstype2state, obsval, obsvetting, optimality, mod_pio_netcdf::pio_netcdf_get_dim(), processobs, residual, ritz, mod_iounits::sourcefile, mod_iounits::stdout, tl_cg_beta, tl_cg_delta, tl_cg_gnorm_v, tl_cg_innov, tl_cg_pxout, tl_cg_pxsave, tl_cg_qg, tl_jb0, tl_obserr, tl_obsscale, tl_obsval, tl_vcglwk, tl_vgrad0, tl_vgrad0s, tl_zcglwk, tl_zgrad0, tobs, unvetted, strings_mod::uppercase(), uradial, vcglev, vcglwk, vgnorm, vgrad0, vgrad0s, mod_ncparam::vname, vradial, wrote_zobs, wrtimpact_bc, wrtimpact_fc, wrtimpact_ic, wrtimpact_tot, wrtmisfit, wrtnlmod, wrtobsscale, wrtrpmod, wrttlmod, wrtzetaref, xobs, yobs, zcglwk, zgrad0, zlanczos_coef, zlanczos_diag, zlanczos_err, zlanczos_inv, zlanczos_offdiag, and zobs.

Referenced by mod_arrays::roms_initialize_arrays(), and roms_kernel_mod::roms_run().

Here is the call graph for this function:
Here is the caller graph for this function:

Variable Documentation

◆ ad_cg_beta

real(dp), dimension(:), allocatable mod_fourdvar::ad_cg_beta

Definition at line 649 of file mod_fourdvar.F.

649 real(dp), allocatable :: ad_cg_beta(:)

Referenced by initialize_fourdvar().

◆ ad_cg_delta

real(dp), dimension(:), allocatable mod_fourdvar::ad_cg_delta

Definition at line 678 of file mod_fourdvar.F.

678 real(dp), allocatable :: ad_cg_delta(:)

Referenced by initialize_fourdvar().

◆ ad_cg_gnorm

real(dp), dimension(:), allocatable mod_fourdvar::ad_cg_gnorm

Definition at line 680 of file mod_fourdvar.F.

680 real(dp), allocatable :: ad_cg_Gnorm(:)

Referenced by initialize_fourdvar().

◆ ad_cg_gnorm_v

real(dp), dimension(:), allocatable mod_fourdvar::ad_cg_gnorm_v

Definition at line 686 of file mod_fourdvar.F.

686 real(dp), allocatable :: ad_cg_Gnorm_v(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_cg_innov

real(r8), dimension(:), allocatable mod_fourdvar::ad_cg_innov

Definition at line 318 of file mod_fourdvar.F.

318 real(r8), allocatable :: ad_cg_innov(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_cg_pxout

real(r8), dimension(:), allocatable mod_fourdvar::ad_cg_pxout

Definition at line 320 of file mod_fourdvar.F.

320 real(r8), allocatable :: ad_cg_pxout(:)

Referenced by initialize_fourdvar().

◆ ad_cg_pxsave

real(r8), dimension(:), allocatable mod_fourdvar::ad_cg_pxsave

Definition at line 319 of file mod_fourdvar.F.

319 real(r8), allocatable :: ad_cg_pxsave(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_cg_qg

real(dp), dimension(:), allocatable mod_fourdvar::ad_cg_qg

Definition at line 679 of file mod_fourdvar.F.

679 real(dp), allocatable :: ad_cg_QG(:)

Referenced by initialize_fourdvar().

◆ ad_jb0

real(r8), dimension(:), allocatable mod_fourdvar::ad_jb0

Definition at line 314 of file mod_fourdvar.F.

314 real(r8), allocatable :: ad_Jb0(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_obserr

real(r8), dimension(:), allocatable mod_fourdvar::ad_obserr

Definition at line 315 of file mod_fourdvar.F.

315 real(r8), allocatable :: ad_ObsErr(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_obsscale

real(r8), dimension(:), allocatable mod_fourdvar::ad_obsscale

Definition at line 316 of file mod_fourdvar.F.

316 real(r8), allocatable :: ad_ObsScale(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_obsval

real(r8), dimension(:,:), allocatable mod_fourdvar::ad_obsval

Definition at line 322 of file mod_fourdvar.F.

322 real(r8), allocatable :: ad_ObsVal(:,:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), rep_matrix(), and roms_kernel_mod::roms_run().

◆ ad_vcglwk

real(r8), dimension(:,:), allocatable mod_fourdvar::ad_vcglwk

Definition at line 311 of file mod_fourdvar.F.

311 real(r8), allocatable :: ad_vcglwk(:,:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_vgrad0

real(r8), dimension(:), allocatable mod_fourdvar::ad_vgrad0

Definition at line 312 of file mod_fourdvar.F.

312 real(r8), allocatable :: ad_vgrad0(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_vgrad0s

real(r8), dimension(:), allocatable mod_fourdvar::ad_vgrad0s

Definition at line 313 of file mod_fourdvar.F.

313 real(r8), allocatable :: ad_vgrad0s(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_zcglwk

real(r8), dimension(:,:), allocatable mod_fourdvar::ad_zcglwk

Definition at line 306 of file mod_fourdvar.F.

306 real(r8), allocatable :: ad_zcglwk(:,:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ ad_zgrad0

real(r8), dimension(:), allocatable mod_fourdvar::ad_zgrad0

Definition at line 307 of file mod_fourdvar.F.

307 real(r8), allocatable :: ad_zgrad0(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ addotproduct

◆ admodval

◆ admodval_s

real(r8), dimension(:,:), allocatable mod_fourdvar::admodval_s

Definition at line 251 of file mod_fourdvar.F.

251 real(r8), allocatable :: ADmodVal_S(:,:)

Referenced by initialize_fourdvar().

◆ ae_beta

real(r8), dimension(:,:), allocatable mod_fourdvar::ae_beta

Definition at line 739 of file mod_fourdvar.F.

739 real(r8), allocatable :: ae_beta(:,:)

Referenced by initialize_fourdvar(), posterior_mod::lanczos(), and posterior_mod::posterior_tile().

◆ ae_delta

real(r8), dimension(:,:), allocatable mod_fourdvar::ae_delta

Definition at line 738 of file mod_fourdvar.F.

738 real(r8), allocatable :: ae_delta(:,:)

Referenced by posterior_mod::analysis_error(), initialize_fourdvar(), posterior_mod::lanczos(), and posterior_mod::posterior_tile().

◆ ae_gnorm

real(r8), dimension(:), allocatable mod_fourdvar::ae_gnorm

Definition at line 745 of file mod_fourdvar.F.

745 real(r8), allocatable :: ae_Gnorm(:)

Referenced by initialize_fourdvar(), posterior_mod::lanczos(), and posterior_mod::posterior_tile().

◆ ae_ritz

real(r8), dimension(:,:), allocatable mod_fourdvar::ae_ritz

Definition at line 755 of file mod_fourdvar.F.

755 real(r8), allocatable :: ae_Ritz(:,:)

Referenced by initialize_fourdvar(), posterior_mod::posterior_eofs(), and posterior_mod::posterior_tile().

◆ ae_ritzerr

real(r8), dimension(:,:), allocatable mod_fourdvar::ae_ritzerr

Definition at line 756 of file mod_fourdvar.F.

756 real(r8), allocatable :: ae_RitzErr(:,:)

Referenced by initialize_fourdvar(), posterior_mod::posterior_eofs(), and posterior_mod::posterior_tile().

◆ ae_tmatrix

real(r8), dimension(:,:), allocatable mod_fourdvar::ae_tmatrix

Definition at line 749 of file mod_fourdvar.F.

749 real(r8), allocatable :: ae_Tmatrix(:,:)

Referenced by initialize_fourdvar(), and posterior_mod::posterior_tile().

◆ ae_trace

real(r8), dimension(:), allocatable mod_fourdvar::ae_trace

Definition at line 741 of file mod_fourdvar.F.

741 real(r8), allocatable :: ae_trace(:)

Referenced by initialize_fourdvar(), and posterior_mod::posterior_tile().

◆ ae_zu

real(r8), dimension(:,:), allocatable mod_fourdvar::ae_zu

Definition at line 750 of file mod_fourdvar.F.

750 real(r8), allocatable :: ae_zu(:,:)

◆ ae_zv

real(r8), dimension(:,:), allocatable mod_fourdvar::ae_zv

Definition at line 740 of file mod_fourdvar.F.

740 real(r8), allocatable :: ae_zv(:,:)

Referenced by initialize_fourdvar(), posterior_mod::posterior_eofs(), and posterior_mod::posterior_tile().

◆ bgerr

real(r8), dimension(:), allocatable mod_fourdvar::bgerr

◆ bgthresh

real(r8), dimension(:), allocatable mod_fourdvar::bgthresh

Definition at line 594 of file mod_fourdvar.F.

594 real(r8), allocatable :: BgThresh(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_operator(), obs_write_mod::obs_write_nf90(), and obs_write_mod::obs_write_pio().

◆ cg_alpha

real(dp), dimension(:,:), allocatable mod_fourdvar::cg_alpha

Definition at line 653 of file mod_fourdvar.F.

653 real(dp), allocatable :: cg_alpha(:,:)

Referenced by initialize_fourdvar().

◆ cg_beta

◆ cg_beta0

real(r8), dimension(:), allocatable mod_fourdvar::cg_beta0

Definition at line 254 of file mod_fourdvar.F.

254 real(r8), allocatable :: cg_beta0(:)

Referenced by initialize_fourdvar().

◆ cg_delta

◆ cg_dla

◆ cg_gamma

real(dp), dimension(:,:), allocatable mod_fourdvar::cg_gamma

◆ cg_gammam1

real(r8) mod_fourdvar::cg_gammam1 = 0.0_r8

Definition at line 846 of file mod_fourdvar.F.

846 real(r8) :: cg_gammam1 = 0.0_r8

◆ cg_gnorm

◆ cg_gnorm_v

◆ cg_gnorm_y

◆ cg_greduc

◆ cg_innov

real(r8), dimension(:), allocatable mod_fourdvar::cg_innov

Definition at line 299 of file mod_fourdvar.F.

299 real(r8), allocatable :: cg_innov(:)

Referenced by congrad_mod::congrad(), deallocate_fourdvar(), initialize_fourdvar(), and rpcg_lanczos_mod::rpcg_lanczos().

◆ cg_qg

◆ cg_ritz

◆ cg_ritzerr

◆ cg_rnorm

real(r8) mod_fourdvar::cg_rnorm = 0.0_r8

Definition at line 848 of file mod_fourdvar.F.

848 real(r8) :: cg_rnorm = 0.0_r8

◆ cg_sigmam1

real(r8) mod_fourdvar::cg_sigmam1 = 0.0_r8

Definition at line 847 of file mod_fourdvar.F.

847 real(r8) :: cg_sigmam1 = 0.0_r8

◆ cg_tau

real(dp), dimension(:,:), allocatable mod_fourdvar::cg_tau

Definition at line 655 of file mod_fourdvar.F.

655 real(dp), allocatable :: cg_tau(:,:)

Referenced by initialize_fourdvar().

◆ cg_tmatrix

real(dp), dimension(:,:), allocatable mod_fourdvar::cg_tmatrix

◆ cg_zu

◆ cg_zv

◆ dotproduct

real(r8) mod_fourdvar::dotproduct

Definition at line 868 of file mod_fourdvar.F.

868 real(r8) :: DotProduct

Referenced by dotproduct_mod::nl_dotproduct_tile(), and dotproduct_mod::tl_dotproduct_tile().

◆ dtsizeh

real(r8), dimension(:,:), allocatable mod_fourdvar::dtsizeh

◆ dtsizehb

real(r8), dimension(:,:), allocatable mod_fourdvar::dtsizehb

◆ dtsizev

real(r8), dimension(:,:), allocatable mod_fourdvar::dtsizev

◆ dtsizevb

real(r8), dimension(:,:), allocatable mod_fourdvar::dtsizevb

◆ extraindex

integer, dimension(:), allocatable mod_fourdvar::extraindex

Definition at line 383 of file mod_fourdvar.F.

383 integer, allocatable :: ExtraIndex(:)

Referenced by initialize_fourdvar().

◆ extraname

character(len=40), dimension(:), allocatable mod_fourdvar::extraname

Definition at line 385 of file mod_fourdvar.F.

385 character(len=40), allocatable :: ExtraName(:)

Referenced by initialize_fourdvar().

◆ forcetime

◆ fourdvar

◆ g1

real(r8), dimension(1000) mod_fourdvar::g1

Definition at line 876 of file mod_fourdvar.F.

876 real(r8), dimension(1000) :: g1

Referenced by dotproduct_mod::nl_dotproduct_tile().

◆ g2

real(r8), dimension(1000) mod_fourdvar::g2

Definition at line 877 of file mod_fourdvar.F.

877 real(r8), dimension(1000) :: g2

Referenced by dotproduct_mod::tl_dotproduct_tile().

◆ gdgw

real(r8), dimension(:), allocatable mod_fourdvar::gdgw

Definition at line 297 of file mod_fourdvar.F.

297 real(r8), allocatable :: gdgw(:)

Referenced by congrad_mod::congrad(), deallocate_fourdvar(), and initialize_fourdvar().

◆ gmze

real(r8), dimension(:,:), allocatable mod_fourdvar::gmze

Definition at line 253 of file mod_fourdvar.F.

253 real(r8), allocatable :: gmze(:,:)

Referenced by initialize_fourdvar().

◆ graderr

real(r8) mod_fourdvar::graderr

Definition at line 853 of file mod_fourdvar.F.

853 real(r8) :: GradErr

Referenced by wrt_info_mod::wrt_info::wrt_info_nf90(), and wrt_info_mod::wrt_info::wrt_info_pio().

◆ gsmatinv

real(r8), dimension(:,:), allocatable mod_fourdvar::gsmatinv

Definition at line 773 of file mod_fourdvar.F.

773 real(r8), allocatable :: GSmatinv(:,:) ! Inverse Gramm-Schmidt matrix

Referenced by initialize_fourdvar(), and inner2state_mod::tl_inner2state_tile().

◆ gsmatrix

real(r8), dimension(:,:), allocatable mod_fourdvar::gsmatrix

Definition at line 772 of file mod_fourdvar.F.

772 real(r8), allocatable :: GSmatrix(:,:) ! Gramm-Schmidt matrix

Referenced by initialize_fourdvar(), and inner2state_mod::tl_inner2state_tile().

◆ harnoldi

real(r8), dimension(:,:,:), allocatable mod_fourdvar::harnoldi

Definition at line 252 of file mod_fourdvar.F.

252 real(r8), allocatable :: harnoldi(:,:,:)

Referenced by initialize_fourdvar().

◆ haveadmod

logical, dimension(:), allocatable mod_fourdvar::haveadmod

Definition at line 553 of file mod_fourdvar.F.

553 logical, allocatable :: haveADmod(:)

Referenced by allocate_fourdvar().

◆ havenlmod

◆ haveobsmeta

◆ havetlmod

logical, dimension(:), allocatable mod_fourdvar::havetlmod

◆ hevecerr

◆ ig1count

integer mod_fourdvar::ig1count

Definition at line 873 of file mod_fourdvar.F.

873 integer :: ig1count ! counter for g1 dot product

Referenced by dotproduct_mod::nl_dotproduct_tile().

◆ ig2count

integer mod_fourdvar::ig2count

Definition at line 874 of file mod_fourdvar.F.

874 integer :: ig2count ! counter for g2 dot product

Referenced by dotproduct_mod::tl_dotproduct_tile().

◆ innovation

real(r8), dimension(:), allocatable mod_fourdvar::innovation

Definition at line 579 of file mod_fourdvar.F.

579 real(r8), allocatable :: innovation(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_write_nf90(), and obs_write_mod::obs_write_pio().

◆ ioffset

integer, dimension(nweights), parameter mod_fourdvar::ioffset = (/ 0, 1, 0, 1, 0, 1, 0, 1 /)

Definition at line 392 of file mod_fourdvar.F.

392 integer, parameter, dimension(Nweights) :: Ioffset = &
393 & (/ 0, 1, 0, 1, 0, 1, 0, 1 /)

◆ jact_s

real(r8), dimension(:), allocatable mod_fourdvar::jact_s

Definition at line 598 of file mod_fourdvar.F.

598 real(r8), allocatable :: Jact_S(:)

Referenced by initialize_fourdvar(), and rpcg_lanczos_mod::rpcg_lanczos().

◆ jb0

◆ jobs

◆ joffset

integer, dimension(nweights), parameter mod_fourdvar::joffset = (/ 0, 0, 1, 1, 0, 0, 1, 1 /)

Definition at line 394 of file mod_fourdvar.F.

394 integer, parameter, dimension(Nweights) :: Joffset = &
395 & (/ 0, 0, 1, 1, 0, 0, 1, 1 /)

◆ khmax

real(r8), dimension(:), allocatable mod_fourdvar::khmax

Definition at line 635 of file mod_fourdvar.F.

635 real(r8), allocatable :: KhMax(:)

Referenced by allocate_fourdvar(), get_state_mod::get_state_nf90(), get_state_mod::get_state_pio(), initialize_fourdvar(), and metrics_mod::metrics_tile().

◆ khmin

real(r8), dimension(:), allocatable mod_fourdvar::khmin

Definition at line 634 of file mod_fourdvar.F.

634 real(r8), allocatable :: KhMin(:)

Referenced by allocate_fourdvar(), get_state_mod::get_state_nf90(), get_state_mod::get_state_pio(), and initialize_fourdvar().

◆ koffset

integer, dimension(nweights), parameter mod_fourdvar::koffset = (/ 0, 0, 0, 0, 1, 1, 1, 1 /)

Definition at line 396 of file mod_fourdvar.F.

396 integer, parameter, dimension(Nweights) :: Koffset = &
397 & (/ 0, 0, 0, 0, 1, 1, 1, 1 /)

◆ kvmax

real(r8), dimension(:), allocatable mod_fourdvar::kvmax

Definition at line 637 of file mod_fourdvar.F.

637 real(r8), allocatable :: KvMax(:)

Referenced by allocate_fourdvar(), get_state_mod::get_state_nf90(), get_state_mod::get_state_pio(), initialize_fourdvar(), and metrics_mod::metrics_tile().

◆ kvmin

real(r8), dimension(:), allocatable mod_fourdvar::kvmin

Definition at line 636 of file mod_fourdvar.F.

636 real(r8), allocatable :: KvMin(:)

Referenced by allocate_fourdvar(), get_state_mod::get_state_nf90(), get_state_mod::get_state_pio(), and initialize_fourdvar().

◆ laugweak

logical mod_fourdvar::laugweak =.FALSE.

Definition at line 833 of file mod_fourdvar.F.

833 logical :: LaugWeak=.false.

Referenced by convolve_mod::error_covariance(), and rbl4dvar_mod::increment().

◆ lhessianev

◆ lhotstart

logical mod_fourdvar::lhotstart

◆ load_zobs

logical, dimension(:), allocatable mod_fourdvar::load_zobs

Definition at line 353 of file mod_fourdvar.F.

353 logical, allocatable :: Load_Zobs(:)

Referenced by allocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_operator(), obs_read_mod::obs_read_nf90(), and obs_read_mod::obs_read_pio().

◆ lobspace

logical, dimension(:), allocatable mod_fourdvar::lobspace

Definition at line 904 of file mod_fourdvar.F.

904 logical, allocatable :: Lobspace(:)

Referenced by ad_htobs_mod::ad_htobs_tile(), ad_initial(), allocate_fourdvar(), and roms_kernel_mod::roms_run().

◆ lprecond

◆ lritz

◆ lsadd

logical, dimension(:), allocatable mod_fourdvar::lsadd

Definition at line 489 of file mod_fourdvar.F.

489 logical, allocatable :: Lsadd(:)

Referenced by ad_main3d(), allocate_fourdvar(), initialize_fourdvar(), output(), and tl_output().

◆ lsen4dvar

logical, dimension(:), allocatable mod_fourdvar::lsen4dvar

Definition at line 902 of file mod_fourdvar.F.

902 logical, allocatable :: Lsen4DVAR(:)

Referenced by ad_get_data(), ad_initial(), ad_main3d(), ad_set_data_tile(), allocate_fourdvar(), and roms_kernel_mod::roms_run().

◆ lsenfct

logical, dimension(:), allocatable mod_fourdvar::lsenfct

Definition at line 893 of file mod_fourdvar.F.

893 logical, allocatable :: LsenFCT(:)

Referenced by ad_initial(), ad_main3d(), allocate_fourdvar(), and roms_kernel_mod::roms_run().

◆ lweakrelax

logical, dimension(:), allocatable mod_fourdvar::lweakrelax

Definition at line 890 of file mod_fourdvar.F.

890 logical, allocatable :: LweakRelax(:)

Referenced by ad_rhs3d_mod::ad_rhs3d(), allocate_fourdvar(), r4dvar_mod::increment(), and tl_rhs3d_mod::tl_rhs3d().

◆ misfit

real(r8), dimension(:), allocatable mod_fourdvar::misfit

Definition at line 242 of file mod_fourdvar.F.

242 real(r8), allocatable :: misfit(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_write_nf90(), and obs_write_mod::obs_write_pio().

◆ mobs

◆ nconvritz

◆ ndatum

◆ nendobs

◆ nextraobs

integer mod_fourdvar::nextraobs = 0

Definition at line 381 of file mod_fourdvar.F.

381 integer :: NextraObs = 0

Referenced by initialize_fourdvar().

◆ nhsteps

integer, dimension(:,:), allocatable mod_fourdvar::nhsteps

◆ nhstepsb

integer, dimension(:,:), allocatable mod_fourdvar::nhstepsb

◆ nimpact

integer mod_fourdvar::nimpact

◆ nlincrement

real(r8), dimension(:), allocatable mod_fourdvar::nlincrement

◆ nlmodval

◆ nmethod

integer, dimension(:), allocatable mod_fourdvar::nmethod

◆ nobs

integer, dimension(:), allocatable mod_fourdvar::nobs

◆ nobsvar

◆ nposti

◆ nrandom

integer mod_fourdvar::nrandom = 1000

◆ nritzev

◆ nstatevar

integer, dimension(:), allocatable mod_fourdvar::nstatevar

Definition at line 362 of file mod_fourdvar.F.

362 integer, allocatable :: NstateVar(:)

Referenced by ad_def_his_mod::ad_def_his_nf90(), ad_def_his_mod::ad_def_his_pio(), allocate_fourdvar(), r4dvar_mod::analysis(), posterior_mod::analysis_error(), back_cost_mod::back_cost_tile(), cgradient_mod::cgradient_tile(), comp_jb0_mod::comp_jb0_tile(), def_avg_mod::def_avg_nf90(), def_avg_mod::def_avg_pio(), def_dai_mod::def_dai_nf90(), def_dai_mod::def_dai_pio(), def_diags_mod::def_diags_nf90(), def_diags_mod::def_diags_pio(), def_error_mod::def_error_nf90(), def_error_mod::def_error_pio(), def_floats_mod::def_floats_nf90(), def_floats_mod::def_floats_pio(), def_hessian_mod::def_hessian_nf90(), def_hessian_mod::def_hessian_pio(), def_his_mod::def_his_nf90(), def_his_mod::def_his_pio(), def_impulse_mod::def_impulse_nf90(), def_impulse_mod::def_impulse_pio(), def_info_mod::def_info::def_info_nf90(), def_info_mod::def_info::def_info_pio(), def_lanczos_mod::def_lanczos_nf90(), def_lanczos_mod::def_lanczos_pio(), def_norm_mod::def_norm_nf90(), def_norm_mod::def_norm_pio(), def_quick_mod::def_quick_nf90(), def_quick_mod::def_quick_pio(), def_rst_mod::def_rst_nf90(), def_rst_mod::def_rst_pio(), def_state_mod::def_state_nf90(), def_state_mod::def_state_pio(), def_station_mod::def_station_nf90(), def_station_mod::def_station_pio(), cgradient_mod::hessian(), cgradient_mod::hessian_evecs(), r4dvar_mod::increment(), ini_lanczos_mod::ini_lanczos_tile(), initialize_fourdvar(), cgradient_mod::lanczos(), posterior_mod::lanczos(), metrics_mod::metrics_tile(), cgradient_mod::new_cost(), cgradient_mod::new_gradient(), posterior_mod::posterior_eofs(), posterior_mod::posterior_tile(), tl_def_his_mod::tl_def_his_nf90(), tl_def_his_mod::tl_def_his_pio(), tl_def_ini_mod::tl_def_ini_nf90(), tl_def_ini_mod::tl_def_ini_pio(), inner2state_mod::tl_inner2state_tile(), wrt_info_mod::wrt_info::wrt_info_nf90(), and wrt_info_mod::wrt_info::wrt_info_pio().

◆ nstrobs

◆ nsurvey

◆ ntimes_ana

integer, dimension(:), allocatable mod_fourdvar::ntimes_ana

Definition at line 911 of file mod_fourdvar.F.

911 integer, allocatable :: ntimes_ana(:)

Referenced by allocate_fourdvar(), roms_kernel_mod::roms_run(), wrt_info_mod::wrt_info::wrt_info_nf90(), and wrt_info_mod::wrt_info::wrt_info_pio().

◆ ntimes_fct

integer, dimension(:), allocatable mod_fourdvar::ntimes_fct

Definition at line 912 of file mod_fourdvar.F.

912 integer, allocatable :: ntimes_fct(:)

Referenced by allocate_fourdvar(), roms_kernel_mod::roms_run(), wrt_info_mod::wrt_info::wrt_info_nf90(), and wrt_info_mod::wrt_info::wrt_info_pio().

◆ nvct

integer mod_fourdvar::nvct

Definition at line 470 of file mod_fourdvar.F.

470 integer :: Nvct

Referenced by array_modes_mod::rep_check(), array_modes_mod::rep_clip(), and array_modes_mod::rep_eigen().

◆ nvsteps

integer, dimension(:,:), allocatable mod_fourdvar::nvsteps

◆ nvstepsb

integer, dimension(:,:), allocatable mod_fourdvar::nvstepsb

◆ nweights

integer, parameter mod_fourdvar::nweights = 8

Definition at line 390 of file mod_fourdvar.F.

390 integer, parameter :: Nweights = 8

◆ obsangler

real(r8), dimension(:), allocatable mod_fourdvar::obsangler

Definition at line 223 of file mod_fourdvar.F.

223 real(r8), allocatable :: ObsAngler(:)

Referenced by ad_htobs_mod::ad_htobs_tile(), ad_misfit_mod::ad_misfit_tile(), deallocate_fourdvar(), initialize_fourdvar(), and obs_write_mod::obs_operator().

◆ obserr

◆ obsmeta

real(r8), dimension(:), allocatable mod_fourdvar::obsmeta

◆ obsname

◆ obsprov

integer, dimension(:), allocatable mod_fourdvar::obsprov

Definition at line 220 of file mod_fourdvar.F.

220 integer, allocatable :: ObsProv(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_operator(), obs_read_mod::obs_read_nf90(), and obs_read_mod::obs_read_pio().

◆ obsscale

◆ obsstate2type

◆ obssurvey

integer, dimension(:), allocatable mod_fourdvar::obssurvey

◆ obstype

◆ obstype2state

integer, dimension(:), allocatable mod_fourdvar::obstype2state

Definition at line 219 of file mod_fourdvar.F.

219 integer, allocatable :: ObsType2State(:)

Referenced by initialize_fourdvar(), obs_cost(), stats_modobs_mod::stats_modobs_nf90(), and stats_modobs_mod::stats_modobs_pio().

◆ obsval

◆ obsvetting

◆ optimality

real(r8), dimension(:), allocatable mod_fourdvar::optimality

Definition at line 478 of file mod_fourdvar.F.

478 real(r8), allocatable :: Optimality(:)

Referenced by allocate_fourdvar(), i4dvar_mod::increment(), initialize_fourdvar(), tl_wrt_ini_mod::tl_wrt_ini_nf90(), and tl_wrt_ini_mod::tl_wrt_ini_pio().

◆ processobs

◆ residual

real(r8), dimension(:), allocatable mod_fourdvar::residual

Definition at line 583 of file mod_fourdvar.F.

583 real(r8), allocatable :: residual(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_write_nf90(), and obs_write_mod::obs_write_pio().

◆ ritz

◆ ritzmaxerr

◆ rscheme

integer, dimension(:), allocatable mod_fourdvar::rscheme

◆ tl_cg_beta

real(dp), dimension(:), allocatable mod_fourdvar::tl_cg_beta

Definition at line 651 of file mod_fourdvar.F.

651 real(dp), allocatable :: tl_cg_beta(:)

Referenced by initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_cg_delta

real(dp), dimension(:), allocatable mod_fourdvar::tl_cg_delta

Definition at line 682 of file mod_fourdvar.F.

682 real(dp), allocatable :: tl_cg_delta(:)

Referenced by initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_cg_gnorm_v

real(dp), dimension(:), allocatable mod_fourdvar::tl_cg_gnorm_v

Definition at line 685 of file mod_fourdvar.F.

685 real(dp), allocatable :: tl_cg_Gnorm_v(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_cg_innov

real(r8), dimension(:), allocatable mod_fourdvar::tl_cg_innov

Definition at line 337 of file mod_fourdvar.F.

337 real(r8), allocatable :: tl_cg_innov(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ tl_cg_pxout

real(r8), dimension(:), allocatable mod_fourdvar::tl_cg_pxout

Definition at line 339 of file mod_fourdvar.F.

339 real(r8), allocatable :: tl_cg_pxout(:)

Referenced by initialize_fourdvar().

◆ tl_cg_pxsave

real(r8), dimension(:), allocatable mod_fourdvar::tl_cg_pxsave

Definition at line 338 of file mod_fourdvar.F.

338 real(r8), allocatable :: tl_cg_pxsave(:)

Referenced by deallocate_fourdvar(), and initialize_fourdvar().

◆ tl_cg_qg

real(dp), dimension(:), allocatable mod_fourdvar::tl_cg_qg

Definition at line 683 of file mod_fourdvar.F.

683 real(dp), allocatable :: tl_cg_QG(:)

Referenced by initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_jb0

real(r8), dimension(:), allocatable mod_fourdvar::tl_jb0

Definition at line 333 of file mod_fourdvar.F.

333 real(r8), allocatable :: tl_Jb0(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_obserr

real(r8), dimension(:), allocatable mod_fourdvar::tl_obserr

Definition at line 334 of file mod_fourdvar.F.

334 real(r8), allocatable :: tl_ObsErr(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_obsscale

real(r8), dimension(:), allocatable mod_fourdvar::tl_obsscale

Definition at line 335 of file mod_fourdvar.F.

335 real(r8), allocatable :: tl_ObsScale(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_obsval

real(r8), dimension(:), allocatable mod_fourdvar::tl_obsval

Definition at line 340 of file mod_fourdvar.F.

340 real(r8), allocatable :: tl_ObsVal(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_vcglwk

real(r8), dimension(:,:), allocatable mod_fourdvar::tl_vcglwk

Definition at line 330 of file mod_fourdvar.F.

330 real(r8), allocatable :: tl_vcglwk(:,:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_vgrad0

real(r8), dimension(:), allocatable mod_fourdvar::tl_vgrad0

Definition at line 331 of file mod_fourdvar.F.

331 real(r8), allocatable :: tl_vgrad0(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_vgrad0s

real(r8), dimension(:), allocatable mod_fourdvar::tl_vgrad0s

Definition at line 332 of file mod_fourdvar.F.

332 real(r8), allocatable :: tl_vgrad0s(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_zcglwk

real(r8), dimension(:,:), allocatable mod_fourdvar::tl_zcglwk

Definition at line 327 of file mod_fourdvar.F.

327 real(r8), allocatable :: tl_zcglwk(:,:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tl_zgrad0

real(r8), dimension(:), allocatable mod_fourdvar::tl_zgrad0

Definition at line 328 of file mod_fourdvar.F.

328 real(r8), allocatable :: tl_zgrad0(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and tl_rpcg_lanczos().

◆ tobs

◆ unvetted

real(r8), dimension(:), allocatable mod_fourdvar::unvetted

Definition at line 243 of file mod_fourdvar.F.

243 real(r8), allocatable :: unvetted(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_operator(), obs_write_mod::obs_write_nf90(), and obs_write_mod::obs_write_pio().

◆ uradial

real(r8), dimension(:), allocatable mod_fourdvar::uradial

Definition at line 244 of file mod_fourdvar.F.

244 real(r8), allocatable :: uradial(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and obs_write_mod::obs_operator().

◆ vcglev

◆ vcglwk

◆ vgnorm

real(r8), dimension(:), allocatable mod_fourdvar::vgnorm

Definition at line 298 of file mod_fourdvar.F.

298 real(r8), allocatable :: vgnorm(:)

Referenced by initialize_fourdvar().

◆ vgrad0

◆ vgrad0s

real(r8), dimension(:), allocatable mod_fourdvar::vgrad0s

Definition at line 296 of file mod_fourdvar.F.

296 real(r8), allocatable :: vgrad0s(:)

Referenced by congrad_mod::congrad(), deallocate_fourdvar(), initialize_fourdvar(), and rpcg_lanczos_mod::rpcg_lanczos().

◆ vradial

real(r8), dimension(:), allocatable mod_fourdvar::vradial

Definition at line 245 of file mod_fourdvar.F.

245 real(r8), allocatable :: vradial(:)

Referenced by deallocate_fourdvar(), initialize_fourdvar(), and obs_write_mod::obs_operator().

◆ wrote_zobs

logical, dimension(:), allocatable mod_fourdvar::wrote_zobs

◆ wrtforce

logical, dimension(:), allocatable mod_fourdvar::wrtforce

◆ wrtimpact_bc

logical, dimension(:), allocatable mod_fourdvar::wrtimpact_bc

Definition at line 543 of file mod_fourdvar.F.

543 logical, allocatable :: wrtIMPACT_BC(:)

Referenced by allocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_write_nf90(), and obs_write_mod::obs_write_pio().

◆ wrtimpact_fc

logical, dimension(:), allocatable mod_fourdvar::wrtimpact_fc

Definition at line 536 of file mod_fourdvar.F.

536 logical, allocatable :: wrtIMPACT_FC(:)

Referenced by allocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_write_nf90(), and obs_write_mod::obs_write_pio().

◆ wrtimpact_ic

logical, dimension(:), allocatable mod_fourdvar::wrtimpact_ic

Definition at line 530 of file mod_fourdvar.F.

530 logical, allocatable :: wrtIMPACT_IC(:)

Referenced by allocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_write_nf90(), and obs_write_mod::obs_write_pio().

◆ wrtimpact_tot

logical, dimension(:), allocatable mod_fourdvar::wrtimpact_tot

Definition at line 523 of file mod_fourdvar.F.

523 logical, allocatable :: wrtIMPACT_TOT(:)

Referenced by allocate_fourdvar(), initialize_fourdvar(), obs_write_mod::obs_write_nf90(), and obs_write_mod::obs_write_pio().

◆ wrtmisfit

◆ wrtnlmod

◆ wrtobsscale

◆ wrtrpmod

◆ wrttlmod

◆ wrtzetaref

logical, dimension(:), allocatable mod_fourdvar::wrtzetaref

◆ xobs

◆ yobs

◆ zcglwk

◆ zgrad0

◆ zlanczos_coef

real(r8), dimension(:,:), allocatable mod_fourdvar::zlanczos_coef

Definition at line 763 of file mod_fourdvar.F.

763 real(r8), allocatable :: zLanczos_coef(:,:) ! coefficients

Referenced by initialize_fourdvar(), posterior_var_mod::posterior_var_tile(), wrt_error_mod::wrt_error_nf90(), and wrt_error_mod::wrt_error_pio().

◆ zlanczos_diag

real(r8), dimension(:), allocatable mod_fourdvar::zlanczos_diag

Definition at line 769 of file mod_fourdvar.F.

769 real(r8), allocatable :: zLanczos_diag(:) ! diagonal coefs

Referenced by initialize_fourdvar(), and inner2state_mod::tl_inner2state_tile().

◆ zlanczos_err

real(r8), dimension(:,:), allocatable mod_fourdvar::zlanczos_err

Definition at line 765 of file mod_fourdvar.F.

765 real(r8), allocatable :: zLanczos_err(:,:) ! error, identity

Referenced by initialize_fourdvar(), posterior_var_mod::posterior_var_tile(), wrt_error_mod::wrt_error_nf90(), and wrt_error_mod::wrt_error_pio().

◆ zlanczos_inv

real(r8), dimension(:,:), allocatable mod_fourdvar::zlanczos_inv

Definition at line 764 of file mod_fourdvar.F.

764 real(r8), allocatable :: zLanczos_inv(:,:) ! inverse

Referenced by initialize_fourdvar(), posterior_var_mod::posterior_var_tile(), wrt_error_mod::wrt_error_nf90(), and wrt_error_mod::wrt_error_pio().

◆ zlanczos_offdiag

real(r8), dimension(:), allocatable mod_fourdvar::zlanczos_offdiag

Definition at line 770 of file mod_fourdvar.F.

770 real(r8), allocatable :: zLanczos_offdiag(:) ! off-diagonal coefs

Referenced by initialize_fourdvar(), and inner2state_mod::tl_inner2state_tile().

◆ zobs