Skip to content

Commit

Permalink
Adds extra bounds on w for the cases
Browse files Browse the repository at this point in the history 
when the uniformly distributed point is near 0 or 1.

Also adds a comment about outputting SILHS subcolumns
to the README.

For #1096.
  • Loading branch information
@vlarson
vlarson committed Jul 12, 2023
1 parent a073701 commit 128a9f8
Show file tree
Hide file tree
Showing 2 changed files with 95 additions and 4 deletions.
15 changes: 15 additions & 0 deletions README
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1163,6 +1163,21 @@ the zt file, the rest (e.g. variance, flux) all occur in the zm file.
Note that the scalars can be used in a host model as well. See SAM-CLUBB for
an example of how to do this.

----------------------------------------------------------------------
- (4.2) SILHS
---------------------------------------------------------------------

SILHS is code that draws samples from CLUBB's PDF and
feeds the samples into microphysics.

An example of how to configure a SILHS run can be
seen in input/case_setups/arm_97_model.in.

In order to output individual SILHS subcolumns to disk,
set l_output_2D_lognormal_dist = .true. in
src/SILHS/latin_hypercube_driver_module.F90.


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% CHAPTER 5: CODE CONTRIBUTED BY EXTERNAL RESEARCH GROUPS
Expand Down
84 changes: 80 additions & 4 deletions src/SILHS/transform_to_pdf_module.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -212,7 +212,10 @@ subroutine sample_w_using_invrs_cdf &
iiPDF_w ! Variable(s)

use constants_clubb, only: &
one
one, &
two, &
fstderr


use clubb_precision, only: &
core_rknd
Expand Down Expand Up @@ -252,9 +255,15 @@ subroutine sample_w_using_invrs_cdf &
integer :: i, k, sample, p ! Loop counters

real( kind = core_rknd ) :: &
U_pt, & ! values from uniformly distributed sample array, X_u_all_levs
mixt_frac, & ! value of mixture fraction from pdf_params%mixt_frac
interp_coef ! interpolation coefficient between low and high values of w sample
U_pt, & ! values from uniformly distributed sample array, X_u_all_levs
mixt_frac, & ! value of mixture fraction from pdf_params%mixt_frac
interp_coef, & ! interpolation coefficient between low and high values of w sample
U_approx ! back-solved value of U, used to check accuracy of sampling

real( kind = core_rknd ), dimension(1,1,1,1) :: &
U_pt_array, & ! Modified sample
w_n_max, & ! normalized version of w_max
w_n_min ! normalized version of w_min

! Flag to clip sample point values of chi in extreme situations.
logical, parameter :: &
Expand Down Expand Up @@ -285,6 +294,38 @@ subroutine sample_w_using_invrs_cdf &
+ std_normal(iiPDF_w,1:ngrdcol,1:nz,sample) * sqrt( pdf_params%varnce_w_1(1:ngrdcol,1:nz) )
end do

! Overwrite w_max if we can find a better upper bound
do k = 1, nz
do sample = 1, num_samples
do i = 1, ngrdcol
if ( ( X_u_all_levs(i,sample,k,iiPDF_w) < 0.9_core_rknd * ( one - pdf_params%mixt_frac(i,k) ) ) &
.and. ( pdf_params%mixt_frac(i,k) < 0.5_core_rknd ) ) then
U_pt_array = X_u_all_levs(i,sample,k,iiPDF_w)/(one-pdf_params%mixt_frac(i,k))
call cdfnorminv( 1, 1, 1, 1, U_pt_array, & ! In
w_n_max )
w_max(i,k,sample) = pdf_params%w_2(i,k) + w_n_max(1,1,1,1) * sqrt( pdf_params%varnce_w_2(i,k) )
end if
end do
end do
end do

! Overwrite w_min if we can find a better lower bound
do k = 1, nz
do sample = 1, num_samples
do i = 1, ngrdcol
if ( ( X_u_all_levs(i,sample,k,iiPDF_w) > ( one - 0.9_core_rknd * pdf_params%mixt_frac(i,k) ) ) &
.and. ( pdf_params%mixt_frac(i,k) > 0.5_core_rknd ) ) then
U_pt_array = ( X_u_all_levs(i,sample,k,iiPDF_w) - (one-pdf_params%mixt_frac(i,k)) ) / pdf_params%mixt_frac(i,k)
call cdfnorminv( 1, 1, 1, 1, &
U_pt_array, & ! In
w_n_min )
w_min(i,k,sample) = pdf_params%w_1(i,k) + w_n_min(1,1,1,1) * sqrt( pdf_params%varnce_w_1(i,k) )
end if
end do
end do
end do


do k = 1, nz
do sample = 1, num_samples
do i = 1, ngrdcol
Expand All @@ -298,6 +339,41 @@ subroutine sample_w_using_invrs_cdf &
end do
end do

! Assertion check: Are w samples accurate?
do k = 1, nz
do sample = 1, num_samples
do i = 1, ngrdcol
U_pt = X_u_all_levs(i,sample,k,iiPDF_w)
mixt_frac = pdf_params%mixt_frac(i,k)
interp_coef = ( U_pt - ( one - mixt_frac ) * U_pt ) &
/ ( ( one - mixt_frac ) * ( one - U_pt ) + mixt_frac * U_pt )
U_approx = mixt_frac * 0.5_core_rknd * ( one &
+ erf( ( X_nl_all_levs(i,sample,k,iiPDF_w) - pdf_params%w_1(i,k) ) &
/ sqrt( pdf_params%varnce_w_1(i,k)*two + 1e-9_core_rknd ) ) ) &
+ ( one - mixt_frac ) * 0.5_core_rknd * ( one &
+ erf( ( X_nl_all_levs(i,sample,k,iiPDF_w) - pdf_params%w_2(i,k) ) &
/ sqrt( pdf_params%varnce_w_2(i,k)*two + 1e-9_core_rknd ) ) )
if ( abs( U_approx - U_pt ) > 0.2_core_rknd &
.and. ( abs( pdf_params%w_1(i,k) ) > 0.02_core_rknd &
.or. abs( pdf_params%w_2(i,k) ) > 0.02_core_rknd ) &
) then
!write(fstderr, *) "abs( U_approx - U ) > 0.1 in subroutine sample_w_using_invrs_cdf!"
write(fstderr, *) "------------------"
write(fstderr, *) "k = ", k
write(fstderr, *) "U_pt = ", U_pt
write(fstderr, *) "U_approx = ", U_approx
write(fstderr, *) "1 - mixt_frac = ", one - mixt_frac
write(fstderr, *) "mixt_frac = ", mixt_frac
write(fstderr, *) "interp_coef = ", interp_coef
write(fstderr, *) "w_sample = ", X_nl_all_levs(i,sample,k,iiPDF_w)
write(fstderr, *) "w_min, w_2 = ", w_min(i,k,sample), pdf_params%w_2(i,k)
write(fstderr, *) "w_1, w_max = ", pdf_params%w_1(i,k), w_max(i,k,sample)
!write(fstderr, *) "w_2(i,k) = ", pdf_params%w_2(i,k)
!write(fstderr, *) "w_1(i,k) = ", pdf_params%w_1(i,k)
end if
end do
end do
end do

end subroutine sample_w_using_invrs_cdf

Expand Down

0 comments on commit 128a9f8

Please sign in to comment.