diff --git a/YOG_convection/SAM_consts.F90 b/YOG_convection/SAM_consts.F90 index 99a1700..dc654b3 100644 --- a/YOG_convection/SAM_consts.F90 +++ b/YOG_convection/SAM_consts.F90 @@ -85,7 +85,10 @@ module SAM_consts_mod ! --------------------------- integer, parameter :: input_ver_dim = 30 !! The number of cells in a SAM atmospheric column on which the neural net was trained - + integer, parameter :: num_sam_cells = 30 + !! number of SAM cells + integer, parameter :: sam_sounding = 48 + !! SAM sounding data ! Outputs from NN are supplied at lowest 30 half-model levels for sedimentation fluxes, ! and at 29 levels for fluxes (as flux at bottom boundary is zero). integer, parameter :: nrf = 30 @@ -97,6 +100,18 @@ module SAM_consts_mod !! SAM timestep in seconds !--------------------------------------------------------------------- + + ! --------------------------- + ! CAM column variables, number of cells + ! --------------------------- + integer, parameter :: num_cols = 4 + !! Number of columns + integer, parameter :: num_cells = 3 + !! Number of cells + integer, parameter :: num_cam_cells_fine = 90 + !! Number of fine CAM cells + integer, parameter :: num_cam_cells_coarse = 10 + !! Number of coarse CAM cells ! Functions and Subroutines contains diff --git a/tests/test_CAM_interface/test_cam_interface.F90 b/tests/test_CAM_interface/test_cam_interface.F90 index 37986d6..f73e238 100644 --- a/tests/test_CAM_interface/test_cam_interface.F90 +++ b/tests/test_CAM_interface/test_cam_interface.F90 @@ -6,6 +6,8 @@ module cam_tests use netcdf use precision, only: dp + use SAM_consts_mod, only: nrf, num_cols, sam_sounding, num_cells, num_cam_cells_fine, num_sam_cells, & + num_cam_cells_coarse use nn_interface_CAM, only: nn_convection_flux_CAM, nn_convection_flux_CAM_init, nn_convection_flux_CAM_finalize, & interp_to_sam, interp_to_cam, fetch_sam_data, SAM_var_conversion, CAM_var_conversion use test_utils, only: assert_array_equal @@ -14,68 +16,75 @@ module cam_tests character(len=15) :: pass = char(27)//'[32m'//'PASSED'//char(27)//'[0m' character(len=15) :: fail = char(27)//'[31m'//'FAILED'//char(27)//'[0m' - integer, parameter :: nrf = 30 integer, parameter :: n_nn_out = 148 real(dp), dimension(n_nn_out) :: nn_out_ones contains subroutine test_interp_to_sam_match(test_name) - !! Check interpolation to SAM grid by interpolating variable of value 1.0 - !! Define a CAM grid from 1111.0 to 10.0 + !! Check interpolation to SAM grid by defining an idential CAM grid and + !! interpolating a variable equal to the pressure + !! Define a CAM grid consiting of 4 atmospheric columns + !! from 1111.0 to 10.0 character(len=*), intent(in) :: test_name integer :: i - real(dp), dimension(4, 48) :: p_cam, p_int_cam - real(dp), dimension(4, 48) :: var_cam - real(dp), dimension(4) :: var_cam_surface - real(dp), dimension(4) :: ps_cam - real(dp), dimension(4, 48) :: var_sam, var_sam_exp + real(dp), dimension(num_cols, nrf) :: p_cam, p_int_cam + real(dp), dimension(num_cols, nrf) :: var_cam + real(dp), dimension(num_cols) :: var_cam_surface + real(dp), dimension(num_cols) :: ps_cam + real(dp), dimension(num_cols, nrf) :: var_sam, var_sam_exp - real(dp), dimension(48) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam + real(dp), dimension(sam_sounding) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam !! Data from the SAM soundings used in tests ! Fetch SAM grid data call fetch_sam_data(pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam) - do i=1,4 - p_cam(i, 1:48) = pres_sam(1:48) - p_int_cam(i, 1:47) = presi_sam(2:48) + do i = 1, num_cols + p_cam(i, :) = pres_sam(1:nrf) + p_int_cam(i, :) = presi_sam(1:nrf) ! Set SAM variable equal to cell size (density 1.0) - var_cam(i, :) = pres_sam(1:48) + var_cam(i, :) = pres_sam(1:nrf) enddo - ! Set interface of top of CAM grid - p_int_cam(:,48) = p_cam(:, 48) + (p_int_cam(:, 47)-p_cam(:, 48)) ps_cam(:) = presi_sam(1) var_cam_surface(:) = presi_sam(1) call interp_to_sam(p_cam, ps_cam, var_cam, var_sam, var_cam_surface) - ! Set anything above 30 elems to zero as per interpolation routine + ! Compare the results of the interpolation scheme to expected output + ! Set anything above 30 elems to zero as the parameterization and interpolation + ! code only uses the bottom 30 cells on the SAM grid var_sam_exp = var_cam - var_sam_exp(:, 31:48) = 0.0 call assert_array_equal(var_sam, var_sam_exp, test_name) end subroutine test_interp_to_sam_match subroutine test_interp_to_sam_one(test_name) - !! Check interpolation to SAM grid by interpolating variable of value 1.0 - !! Define a CAM grid from 1111.0 to 10.0 + !! Check interpolation to SAM grid by setting up a coarse CAM grid with every variable being 1.0 + !! interpolating to a new grid should also have of value 1.0 everywhere. character(len=*), intent(in) :: test_name - real(dp), dimension(4, 3) :: p_cam - real(dp), dimension(4, 3) :: var_cam - real(dp), dimension(4) :: var_cam_surface - real(dp), dimension(4) :: ps_cam - real(dp), dimension(4, 30) :: var_sam, var_sam_exp + integer :: i - p_cam = reshape((/ 1000.0, 1000.0, 1000.0, 1000.0, 500.0, 500.0, 500.0, 500.0, 10.0, 10.0, 10.0, 10.0 /), (/ 4, 3 /)) + real(dp), dimension(num_cols, num_cells) :: p_cam + real(dp), dimension(num_cols, num_cells) :: var_cam + real(dp), dimension(num_cols) :: var_cam_surface + real(dp), dimension(num_cols) :: ps_cam + real(dp), dimension(num_cols, num_sam_cells) :: var_sam, var_sam_exp + + ! Set up a coarse CAM grid of 4 columns of pressures [1000, 500, 10] hPa with surface pressure 1111 hPa + do i = 1, num_cols + p_cam(i, 1) = 1000.0 + p_cam(i, 2) = 500.0 + p_cam(i, 3) = 10.0 + end do ps_cam = (/ 1111.0, 1111.0, 1111.0, 1111.0/) - + var_cam = 1.0 var_cam_surface = 1.0 var_sam_exp = 1.0 @@ -88,41 +97,44 @@ end subroutine test_interp_to_sam_one subroutine test_interp_to_sam_pres(test_name) !! Check interpolation to SAM grid by interpolating pressure to pressure - !! Set top of CAM to 1.0d-4 - !! => should match pres from SAM + !! Use a coarse CAM grid of 3 cells and 4 columns + !! => expected variable on SAM grid should be pressure at that point character(len=*), intent(in) :: test_name integer :: i - real(dp), dimension(4, 3) :: p_cam - real(dp), dimension(4, 3) :: var_cam - real(dp), dimension(4) :: var_cam_surface - real(dp), dimension(4) :: ps_cam - real(dp), dimension(4, 30) :: var_sam, var_sam_exp + real(dp), dimension(num_cols, num_cells) :: p_cam + real(dp), dimension(num_cols, num_cells) :: var_cam + real(dp), dimension(num_cols) :: var_cam_surface + real(dp), dimension(num_cols) :: ps_cam + real(dp), dimension(num_cols, num_sam_cells) :: var_sam, var_sam_exp - real(dp), dimension(48) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam + real(dp), dimension(sam_sounding) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam !! Data from the SAM soundings used in tests ! Fetch SAM grid data call fetch_sam_data(pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam) - p_cam(:, 1) = pres_sam(5) - p_cam(:, 2) = pres_sam(10) - p_cam(:, 3) = 1.0d-4 + ! Set up coarse CAM grid with surface pressure equal to SAM surface, + ! a top pressure of 1.0d-4 above the SAM grid, and other at 10 and 20 from the SAM grid ps_cam = presi_sam(1) - var_cam(:, 1) = pres_sam(5) - var_cam(:, 2) = pres_sam(10) - var_cam(:, 3) = 1.0d-4 - var_cam_surface = presi_sam(1) + p_cam(:, 1) = pres_sam(10) + p_cam(:, 2) = pres_sam(20) + p_cam(:, 3) = 1.0d-4 - do i = 1,4 - var_sam_exp(i, :) = pres_sam(1:30) + ! Set the variable on the CAM grid equal to the pressure + var_cam_surface = ps_cam + var_cam = p_cam + + ! Expected variable value on the SAM grid will be equal to pressure at that point + do i = 1, num_cols + var_sam_exp(i, :) = pres_sam(1:num_sam_cells) end do call interp_to_sam(p_cam, ps_cam, var_cam, var_sam, var_cam_surface) - call assert_array_equal(var_sam_exp, var_sam, test_name) + call assert_array_equal(var_sam_exp, var_sam, test_name, rtol_opt = 1.0d-14) end subroutine test_interp_to_sam_pres @@ -135,43 +147,104 @@ subroutine test_interp_to_cam_match(test_name) integer :: i - real(dp), dimension(4, 30) :: p_cam, var_cam, rho_cam, rho_cam_exp - real(dp), dimension(4, 31) :: p_int_cam - real(dp), dimension(4) :: ps_cam, var_cam_surface - real(dp), dimension(4, 30) :: var_sam + real(dp), dimension(num_cols, nrf) :: p_cam, var_cam, rho_cam, rho_cam_exp + real(dp), dimension(num_cols, nrf + 1) :: p_int_cam + real(dp), dimension(num_cols) :: ps_cam, var_cam_surface + real(dp), dimension(num_cols, nrf) :: var_sam - real(dp), dimension(48) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam + real(dp), dimension(sam_sounding) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam !! Data from the SAM soundings used in tests ! Fetch SAM grid data call fetch_sam_data(pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam) - do i=1,4 - p_cam(i, 1:30) = pres_sam(1:30) - p_int_cam(i, 1:31) = presi_sam(1:31) + do i = 1, num_cols + ! Set up a CAM grid that matches the lower nrf cells of the SAM grid + p_cam(i, :) = pres_sam(1 : num_sam_cells) + p_int_cam(i, :) = presi_sam(1 : nrf + 1) ! Set SAM variable equal to cell size ("density" 1.0) - var_sam(i, :) = (presi_sam(1:30) - presi_sam(2:31)) - ! var_sam(i, :) = 1.0 + var_sam(i, :) = (presi_sam(1 : nrf) - presi_sam(2 : nrf + 1)) enddo call interp_to_cam(p_cam, p_int_cam, p_int_cam(:, 1), var_sam, var_cam) - do i=1,30 - rho_cam(:, i) = var_cam(:, i) / (p_int_cam(:, i)-p_int_cam(:, i+1)) + ! Calculate resulting density by dividing by cell size + do i = 1, nrf + rho_cam(:, i) = var_cam(:, i) / (p_int_cam(:, i) - p_int_cam(:, i + 1)) end do rho_cam_exp = 1.0 call assert_array_equal(rho_cam_exp, rho_cam, test_name) - ! Could do an additional check on the sums of the variables to check - ! variable has been conserved - ! do i=1,4 - ! write(*,*) sum(var_sam(i,:)), sum(var_cam(i,:)) - ! enddo - + ! Here in this test, we have conservation + !do i = 1, num_cols + ! write(*,*) sum(var_sam(i, :)), sum(var_cam(i, :)) + !enddo + end subroutine test_interp_to_cam_match + subroutine test_interp_to_cam_coarse(test_name) + !! Check conservative regridding to CAM coarse grid + !! => should conserve density + + character(len=*), intent(in) :: test_name + integer :: i, j + + real(dp), dimension(num_cols, num_cam_cells_coarse) :: p_cam, var_cam, rho_cam, rho_cam_exp + real(dp), dimension(num_cols, num_cam_cells_coarse - 1) :: p_int_cam + real(dp), dimension(num_cols) :: ps_cam, sum_sam, sum_cam, var_cam_surface + real(dp), dimension(num_cols, num_sam_cells) :: var_sam + + real(dp), dimension(sam_sounding) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam + !! Data from the SAM soundings used in tests + + ! Fetch SAM grid data + call fetch_sam_data(pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam) + + ! Define CAM grid coarser than SAM grid + do j = 1, num_cam_cells_coarse - 1 + p_int_cam(:,j) = presi_sam(3 * (j - 1) + 1) + enddo + + !write(*,*) p_int_cam(:,1) + + ! Get the CAM pressures from average interface pressures + do j = 1, num_cam_cells_coarse - 1 + p_cam(:, j) = (p_int_cam(:, j + 1) + p_int_cam(:, j)) / 2.0 + end do + + !write(*,*) p_cam(:,1) + !write(*,*) p_cam(:,2) + + do j = 1, num_cols + ! Set SAM variable equal to cell size + var_sam(j, :) = (presi_sam(1 : num_sam_cells) - presi_sam(2 : num_sam_cells + 1)) + end do + + call interp_to_cam(p_cam, p_int_cam, p_int_cam(:, 1), var_sam, var_cam) + + ! Loop over cam coarse cells + do i = 1, num_cam_cells_coarse - 1 + rho_cam(:, i) = var_cam(:, i) / (p_int_cam(:, i) - p_int_cam(:, i + 1)) + enddo + + write(*,*) rho_cam(1,:) + write(*,*) var_cam(1,:) + write(*,*) p_int_cam(1,:) + + ! Expected density + rho_cam_exp = 1.0 + + call assert_array_equal(rho_cam, rho_cam_exp, test_name, rtol_opt = 2.0D-5) + + ! checking conservation here (not yet conserved) + do i = 1, num_cols + write(*,*) sum(var_sam(i, :)), sum(var_cam(i, :)) + enddo + + end subroutine test_interp_to_cam_coarse + subroutine test_interp_to_cam_fine(test_name) !! Check conservative regridding to CAM grid by interpolating constant density !! => should conserve density @@ -181,65 +254,55 @@ subroutine test_interp_to_cam_fine(test_name) integer :: i, j - real(dp), dimension(4, 90) :: p_cam, var_cam, rho_cam, rho_cam_exp - real(dp), dimension(4, 91) :: p_int_cam - real(dp), dimension(4) :: var_cam_surface - real(dp), dimension(4) :: ps_cam, sum_sam, sum_cam - real(dp), dimension(4, 30) :: var_sam + real(dp), dimension(num_cols, num_cam_cells_fine) :: p_cam, var_cam, rho_cam, rho_cam_exp + real(dp), dimension(num_cols, num_cam_cells_fine + 1) :: p_int_cam + real(dp), dimension(num_cols) :: var_cam_surface + real(dp), dimension(num_cols) :: ps_cam, sum_sam, sum_cam + real(dp), dimension(num_cols, num_sam_cells) :: var_sam - real(dp), dimension(48) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam + real(dp), dimension(sam_sounding) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam !! Data from the SAM soundings used in tests ! Fetch SAM grid data call fetch_sam_data(pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam) ! CAM grid finer than SAM grid - do j=0, 29 - do i=1,4 + do j = 0, num_sam_cells - 1 + do i = 1, num_cols p_int_cam(i, 3*j+1) = presi_sam(j+1) - p_int_cam(i, 3*j+2) = presi_sam(j+1) + (presi_sam(j+2)-presi_sam(j+1))/3 - p_int_cam(i, 3*j+3) = presi_sam(j+1) + 2*(presi_sam(j+2)-presi_sam(j+1))/3 + p_int_cam(i, 3*j+2) = presi_sam(j+1) + (presi_sam(j+2) - presi_sam(j+1)) / 3 + p_int_cam(i, 3*j+3) = presi_sam(j+1) + 2*(presi_sam(j+2) - presi_sam(j+1)) / 3 + enddo enddo - enddo - p_int_cam(:, 91) = presi_sam(31) - - do j=1, 90 - p_cam(:, j) = (p_int_cam(:, j+1)+p_int_cam(:, j)) / 2.0 + p_int_cam(:, num_cam_cells_fine + 1) = presi_sam(num_sam_cells + 1) + + !write(*,*) p_int_cam(1, 1:9) + !write(*,*) presi_sam(1:3) + do j = 1, num_cam_cells_fine + p_cam(:, j) = (p_int_cam(:, j+1) + p_int_cam(:, j)) / 2.0 end do - do j=1, 4 + do j = 1, num_cols ! Set SAM variable equal to cell size - var_sam(j, :) = (presi_sam(1:30) - presi_sam(2:31)) + var_sam(j, :) = (presi_sam(1 : num_sam_cells) - presi_sam(2 : num_sam_cells + 1)) end do + !write(*,*) var_sam + call interp_to_cam(p_cam, p_int_cam, p_int_cam(:, 1), var_sam, var_cam) - do i=1,92 - rho_cam(:, i) = var_cam(:, i) / (p_int_cam(:, i)-p_int_cam(:, i+1)) + do i = 1, num_cam_cells_fine + 2 + rho_cam(:, i) = var_cam(:, i) / (p_int_cam(:, i) - p_int_cam(:, i+1)) end do rho_cam_exp = 1.0 -! ! Compare individual cells of regridded data -! do i = 1,30 -! write(*,*) i, var_sam(1,i), presi_sam(i) -! ! write(*,*) i, p_int_cam(1,i) -! ! write(*,*) i, p_cam(1,i) -! enddo -! write(*,*) presi_sam(31) -! do i = 1,90 -! write(*,*) i, var_cam(1,i), rho_cam(1,i), rho_cam_exp(1,i), p_int_cam(1, i) -! ! write(*,*) i, p_int_cam(1,i) -! ! write(*,*) i, p_cam(1,i) -! enddo -! write(*,*) p_int_cam(1, 91) -! do j=0, 30 -! write(*,*) var_sam(1, j+1), sum(var_cam(1, 3*j+1:3*j+3)), var_cam(1, 3*j+1:3*j+3) -! enddo - call assert_array_equal(rho_cam, rho_cam_exp, test_name, rtol_opt=2.0D-5) + call assert_array_equal(rho_cam, rho_cam_exp, test_name, rtol_opt = 2.0D-5) + + !write(*,*) rho_cam ! Compare sums of the variables to check conservation between grids - do i=1,4 + do i = 1, num_cols sum_sam(i) = sum(var_sam(i,:)) sum_cam(i) = sum(var_cam(i,:)) enddo @@ -252,10 +315,10 @@ subroutine test_var_conv_sam_zero(test_name) character(len=*), intent(in) :: test_name - real(dp), dimension(4, 48) :: t, q, tabs, qv, qc, qi - real(dp), dimension(4, 48) :: tabs_exp, qv_exp, qc_exp, qi_exp + real(dp), dimension(num_cols, sam_sounding) :: t, q, tabs, qv, qc, qi + real(dp), dimension(num_cols, sam_sounding) :: tabs_exp, qv_exp, qc_exp, qi_exp integer :: i - real(dp), dimension(48) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam + real(dp), dimension(sam_sounding) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam !! Data from the SAM soundings used in tests ! Fetch SAM grid data @@ -277,7 +340,7 @@ subroutine test_var_conv_sam_zero(test_name) qi = qi + 1.0 ! Check that, without moisture, tabs is (t - gamaz) - do i=1,48 + do i = 1, sam_sounding tabs_exp(:,i) = 0.0 - gamaz_sam(i) end do qv_exp = 1.0 @@ -296,10 +359,10 @@ subroutine test_rev_var_conv_sat(test_name) character(len=*), intent(in) :: test_name - real(dp), dimension(4, 48) :: t, q, tabs, qv, qc, qi - real(dp), dimension(4, 48) :: tabs_exp, qv_exp, qc_exp, qi_exp + real(dp), dimension(num_cols, sam_sounding) :: t, q, tabs, qv, qc, qi + real(dp), dimension(num_cols, sam_sounding) :: tabs_exp, qv_exp, qc_exp, qi_exp integer :: i - real(dp), dimension(48) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam + real(dp), dimension(sam_sounding) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam real(dp), parameter :: rgas = 287.0 !! Gas constant for dry air used in SAM [j / kg / K] @@ -308,7 +371,7 @@ subroutine test_rev_var_conv_sat(test_name) call fetch_sam_data(pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam) ! Set up dry saturated air - do i = 1,48 + do i = 1, sam_sounding tabs(:, i) = 100.0 * pres_sam(i) / (rho_sam(i) * rgas) end do qv = 0.0 @@ -344,11 +407,11 @@ subroutine test_rev_var_conv_moist(test_name) character(len=*), intent(in) :: test_name - real(dp), dimension(4, 48) :: t, q, tabs, qv, qc, qi - real(dp), dimension(4, 48) :: r, p_sat, q_sat - real(dp), dimension(4, 48) :: tabs_exp, qv_exp, qc_exp, qi_exp + real(dp), dimension(num_cols, sam_sounding) :: t, q, tabs, qv, qc, qi + real(dp), dimension(num_cols, sam_sounding) :: r, p_sat, q_sat + real(dp), dimension(num_cols, sam_sounding) :: tabs_exp, qv_exp, qc_exp, qi_exp integer :: i - real(dp), dimension(48) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam + real(dp), dimension(sam_sounding) :: pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam real(dp), parameter :: rgas = 287.0 !! Gas constant for dry air used in SAM [j / kg / K] @@ -359,7 +422,7 @@ subroutine test_rev_var_conv_moist(test_name) call fetch_sam_data(pres_sam, presi_sam, gamaz_sam, rho_sam, z_sam) ! Set up moist air - do i = 1,48 + do i = 1, sam_sounding qv(:, i) = 0.01 * exp(-5.0d-4 * z_sam(i)) r(:, i) = qv(:, i) / (1.0 - qv(:, i)) @@ -410,10 +473,10 @@ program run_cam_tests implicit none - real(dp), dimension(4) :: a = [1.0, 2.0, 3.0, 4.0] + real(dp), dimension(num_cols) :: a = [1.0, 2.0, 3.0, 4.0] real(dp), dimension(5) :: b = [1.0, 2.0, 3.0, 4.0, 6.0] - integer, dimension(4) :: aa = [1, 2, 3, 4] - integer, dimension(4) :: bb = [1, 2, 3, 4] + integer, dimension(num_cols) :: aa = [1, 2, 3, 4] + integer, dimension(num_cols) :: bb = [1, 2, 3, 4] character(len=1024) :: nn_file = "../../YOG_convection/NN_weights_YOG_convection.nc" character(len=1024) :: sounding_file = "../../YOG_convection/resources/SAM_sounding.nc" @@ -430,6 +493,7 @@ program run_cam_tests ! SAM to CAM grid interpolation call test_interp_to_cam_match("Test interpolation to CAM mapping density 1:1 match grid") call test_interp_to_cam_fine("Test interpolation to CAM mapping density 1:1 fine grid") + call test_interp_to_cam_coarse("Test interpolation to CAM mapping density 1:1 coarse grid") ! Variable conversion call test_var_conv_sam_zero("Test variable conversion SAM->CAM for 0.0")