[Dart-dev] DART/branches Revision: 12921

[dart at ucar.edu]

nancy at ucar.edu
2018-10-29 13:49:41 -0600 (Mon, 29 Oct 2018)
198
replace < with abs(diff) < TINY for testing to see if
the inflation standard deviation is at the lower threshold.
see comments in the code for reasons.  this is already
committed to the rma_trunk.




Modified: DART/branches/recam/assimilation_code/modules/assimilation/adaptive_inflate_mod.f90
===================================================================
--- DART/branches/recam/assimilation_code/modules/assimilation/adaptive_inflate_mod.f90	2018-10-28 23:37:16 UTC (rev 12920)
+++ DART/branches/recam/assimilation_code/modules/assimilation/adaptive_inflate_mod.f90	2018-10-29 19:49:41 UTC (rev 12921)
@@ -670,8 +670,16 @@
       new_cov_inflate)
 
 ! Bail out to save cost when lower bound is reached on lambda standard deviation
-if(lambda_sd <= sd_lower_bound_in) then
+   ! 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
+   ! really reached.  (sd_lower_bound comes from a namelist and precision
+   ! errors may be because of the conversion between ascii, single precision
+   ! 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
+      return
 else
    ! Compute by forcing a Gaussian fit at one positive SD
 ! First compute the new_max value for normalization purposes
@@ -718,9 +726,10 @@
                     gamma_corr, ens_size, rate, new_cov_inflate)
 
    ! Bail out to save cost when lower bound is reached on lambda standard deviation
-   if(lambda_sd <= sd_lower_bound_in) then
+   ! See comment in Anderson case for why we use abs and TINY for this comparison.
+   if(abs(lambda_sd - sd_lower_bound_in) <= TINY(0.0_r8)) then
       new_cov_inflate_sd = lambda_sd
-   
+      return 
    else
       ! Compute the shape parameter of the prior IG
       ! This comes from the assumption that the mode of the IG is the mean/mode of the input Gaussian