[Dart-dev] DART/branches Revision: 12964
dart at ucar.edu
dart at ucar.edu
Thu Jan 31 16:43:08 MST 2019
nancy at ucar.edu
2019-01-31 16:43:08 -0700 (Thu, 31 Jan 2019)
904
remove the allocated space for obs_seq.final from all
tasks other than the one which is going to collect the
data and write that file. for large ensemble sizes where
all the obs values are being output, this is a big
memory savings.
also add the option to avoid computing the posterior
forward operator values, and don't allocate space in
the obs_seq.final for posteriors. this should be
both a time and space savings.
add a check so users can't ask for both posterior
inflation and to not compute posteriors.
add some code in filter as a preparation for streamlining
the trace and timing messages to try to make the code
easier to read.
update the dopplerfold version of filter to stay in sync
with the changes to filter, even though this is only
used (as far as i know) with the wrf model. still,
this is a proposed version to be migrated to the
rma_trunk and on to Manhattan, so update it now.
Modified: DART/branches/recam/assimilation_code/modules/assimilation/adaptive_inflate_mod.f90
===================================================================
--- DART/branches/recam/assimilation_code/modules/assimilation/adaptive_inflate_mod.f90 2019-01-29 17:54:38 UTC (rev 12963)
+++ DART/branches/recam/assimilation_code/modules/assimilation/adaptive_inflate_mod.f90 2019-01-31 23:43:08 UTC (rev 12964)
@@ -364,7 +364,8 @@
subroutine validate_inflate_options(inf_flavor, inf_damping, inf_initial_from_restart, &
inf_sd_initial_from_restart, inf_deterministic, inf_sd_max_change, &
- do_prior_inflate, do_posterior_inflate, output_inflation)
+ do_prior_inflate, do_posterior_inflate, output_inflation, &
+ compute_posterior)
integer, intent(in) :: inf_flavor(2)
real(r8), intent(inout) :: inf_damping(2)
@@ -375,6 +376,7 @@
logical, intent(out) :: do_prior_inflate
logical, intent(out) :: do_posterior_inflate
logical, intent(out) :: output_inflation
+logical, intent(in) :: compute_posterior
integer :: i
character(len=32) :: string(2)
@@ -411,6 +413,13 @@
if(inf_flavor(1) == 4) call error_handler(E_ERR, 'validate_inflate_options', &
'RTPS inflation (type 4) only supported for Posterior inflation', source, revision, revdate)
+! Cannot select posterior options if not computing posterior
+if(.not. compute_posterior .and. inf_flavor(2) > 0) then
+ write(msgstring, *) 'cannot enable posterior inflation if not computing posterior values'
+ call error_handler(E_ERR,'validate_inflate_options', msgstring, source, revision, revdate, &
+ text2='"compute_posterior" is false; posterior inflation flavor must be 0')
+endif
+
! RTPS needs a single parameter from namelist: inf_initial(2).
! Do not read in any files. Also, no damping. but warn the user if they try to set different
! values in the namelist.
@@ -669,7 +678,7 @@
call linear_bayes(dist_2, sigma_p_2, sigma_o_2, lambda_mean, lambda_sd_2, gamma_corr, &
new_cov_inflate)
-! Bail out to save cost when lower bound is reached on lambda standard deviation
+ ! Bail out to save cost when lower bound is reached on lambda standard deviation
! The original test to see if lambda_sd was less than the lower bound
! would sometimes return false because of roundoff error and the computation
! would go through the expensive part of the code when the minimum was
@@ -678,44 +687,47 @@
! and double precision.) In any case, the test was changed to return if
! the value is within TINY of the limit.
if(abs(lambda_sd - sd_lower_bound_in) <= TINY(0.0_r8)) then
- new_cov_inflate_sd = lambda_sd
+ new_cov_inflate_sd = lambda_sd
return
-else
- ! Compute by forcing a Gaussian fit at one positive SD
-! First compute the new_max value for normalization purposes
- new_max = compute_new_density(dist_2, sigma_p_2, sigma_o_2, lambda_mean, lambda_sd, &
- gamma_corr, new_cov_inflate)
-
-! Find value at a point one OLD sd above new mean value
+ else
+ ! Compute by forcing a Gaussian fit at one positive SD
+ ! First compute the new_max value for normalization purposes
+ new_max = compute_new_density(dist_2, sigma_p_2, sigma_o_2, lambda_mean, lambda_sd, &
+ gamma_corr, new_cov_inflate)
+
+ ! Find value at a point one OLD sd above new mean value
new_1_sd = compute_new_density(dist_2, sigma_p_2, sigma_o_2, lambda_mean, lambda_sd, gamma_corr, &
- new_cov_inflate + lambda_sd)
-
- ! If either the numerator or denominator of the following computation
- ! of 'ratio' is going to be zero (or almost so), return the original incoming
- ! inflation value. The computation would have resulted in either Inf or NaN.
- if (abs(new_max) <= TINY(0.0_r8) .or. abs(new_1_sd) <= TINY(0.0_r8)) then
- new_cov_inflate_sd = lambda_sd
- return
+ new_cov_inflate + lambda_sd)
+
+ ! If either the numerator or denominator of the following computation
+ ! of 'ratio' is going to be zero (or almost so), return the original incoming
+ ! inflation value. The computation would have resulted in either Inf or NaN.
+ if (abs(new_max) <= TINY(0.0_r8) .or. abs(new_1_sd) <= TINY(0.0_r8)) then
+ new_cov_inflate_sd = lambda_sd
+ return
+ endif
+
+ ratio = new_1_sd / new_max
+
+ ! Another error for numerical issues; if ratio is larger than 0.99, bail out
+ if(ratio > 0.99) then
+ new_cov_inflate_sd = lambda_sd
+ return
+ endif
+
+ ! Can now compute the standard deviation consistent with this as
+ ! sigma = sqrt(-x^2 / (2 ln(r)) where r is ratio and x is lambda_sd (distance from mean)
+ new_cov_inflate_sd = sqrt( -1.0_r8 * lambda_sd_2 / (2.0_r8 * log(ratio)))
+
+ ! Prevent an increase in the sd of lambda???
+ ! For now, this is mostly countering numerical errors in this computation
+ if(new_cov_inflate_sd > lambda_sd) then
+ new_cov_inflate_sd = lambda_sd
More information about the Dart-dev
mailing list