From 66221076e7a636590c3c5299aa468e6b0857fcc7 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Thu, 1 Feb 2024 15:13:11 -0700 Subject: [PATCH] addressed code review --- doc/ChangeLog | 4 +- musica/micm/micm.F90 | 375 +++++++++++++++++++++--------------------- musica/micm/micm.meta | 2 +- 3 files changed, 190 insertions(+), 191 deletions(-) diff --git a/doc/ChangeLog b/doc/ChangeLog index 9305d071..47f37983 100644 --- a/doc/ChangeLog +++ b/doc/ChangeLog @@ -4,12 +4,12 @@ Tag name:atmos_phys0_02_001 Originator(s): boulderdaze Date: 29 Jan 2024 One-line Summary: Update MICM wrapper for species unit conversions -Github PR URL: https://github.com/ESCOMP/atmospheric_physics/issues/73 +Github PR URL: https://github.com/ESCOMP/atmospheric_physics/pull/78 Purpose of changes (include the issue number and title text for each relevant GitHub issue): The main purpose of this PR is to add unit conversion functions for species concentrations -from CAM-SIMA units (kg kg-1) to MICM units (mol m-3) and vise versa. +to convert CAM-SIMA units (kg kg-1) to MICM units (mol m-3) and vise versa. This PR fixes the following NCAR/atmospheric_physics Github issues: diff --git a/musica/micm/micm.F90 b/musica/micm/micm.F90 index 4e2afa78..fc6e0523 100644 --- a/musica/micm/micm.F90 +++ b/musica/micm/micm.F90 @@ -1,201 +1,200 @@ module micm - ! Wrapper for MICM functionality - use iso_c_binding + use iso_c_binding - ! Note: "micm_core" is included in an external pre-built MICM library that the host - ! model is responsible for linking to during compilation - use micm_core, only: micm_t - use ccpp_kinds, only: kind_phys + ! Note: "micm_core" is included in an external pre-built MICM library that the host + ! model is responsible for linking to during compilation + use micm_core, only: micm_t + use ccpp_kinds, only: kind_phys - implicit none + implicit none - public :: micm_init, micm_run, micm_final - private :: convert_to_mol_per_cubic_meter, convert_to_mass_mixing_ratio + public :: micm_init, micm_run, micm_final + private :: convert_to_mol_per_cubic_meter, convert_to_mass_mixing_ratio - type(micm_t), allocatable :: micm_obj + type(micm_t), allocatable :: micm_obj contains - !> \section arg_table_micm_init Argument Table - !! \htmlinclude micm_init.html - subroutine micm_init(config_path, iulog, errcode, errmsg) - ! Initializes MICM objects by creating solvers from configure files - - character(len=*), intent(in) :: config_path - integer, intent(in) :: iulog - integer, intent(out) :: errcode - character(len=512), intent(out) :: errmsg - - errcode = 0 - errmsg = '' - - ! Constructs MICM object - allocate(micm_obj) - micm_obj = micm_t(config_path) - - ! Creates solver - errcode = micm_obj%create_solver() - - if (errcode /= 0) then - !TODO(jiwon) - which file? - errmsg = "INIT MICM: FATAL: Failed in creating MICM solver because parsing & - configuration files failed. Please look over at the other file." - else - write(iulog,*) " [INIT MICM]: INFO: Created MICM solver" - endif + !> \section arg_table_micm_init Argument Table + !! \htmlinclude micm_init.html + subroutine micm_init(config_path, iulog, errcode, errmsg) + character(len=*), intent(in) :: config_path + integer, intent(in) :: iulog + integer, intent(out) :: errcode + character(len=512), intent(out) :: errmsg + + errcode = 0 + errmsg = '' + + allocate(micm_obj) + micm_obj = micm_t(config_path) + + errcode = micm_obj%create_solver() + + if (errcode /= 0) then + errmsg = "[fatal] [micm] Failed to creat MICM solver. Parsing configuration failed. & + Please look over at the other file for further information." + write(iulog,*) trim(errmsg) + return + endif + + write(iulog,*) "[info] [micm] Created MICM solver." end subroutine micm_init - !> \section arg_table_micm_run Argument Table - !! \htmlinclude micm_run.html - subroutine micm_run(time_step, temperature, pressure, dry_air_density, constituent_props, & - constituents, iulog, errcode, errmsg) - use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t - - real(kind_phys), intent(in) :: time_step ! s - real(kind_phys), intent(in) :: temperature(:,:) ! K - real(kind_phys), intent(in) :: pressure(:,:) ! Pa - real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3 - type(ccpp_constituent_prop_ptr_t), intent(in) :: constituent_props(:) - real(kind_phys), intent(inout) :: constituents(:,:,:) ! kg kg-1 - integer, intent(in) :: iulog - integer, intent(out) :: errcode - character(len=512), intent(out) :: errmsg - - ! local variables - real(c_double) :: c_time_step - real(c_double), dimension(size(temperature, dim=1), & - size(temperature, dim=2)) :: c_temperature - real(c_double), dimension(size(pressure, dim=1), & - size(pressure, dim=2)) :: c_pressure - real(c_double), dimension(size(constituents, dim=1), & - size(constituents, dim=2), & - size(constituents, dim=3)) :: c_constituents - - real(kind_phys), dimension(size(constituents, dim=3)) :: molar_mass_arr ! kg mol-1 - - integer :: num_columns, num_layers, num_constituents - integer :: i_column, i_layer, i_elem - - num_columns = size(constituents, dim=1) - num_layers = size(constituents, dim=2) - num_constituents = size(constituents, dim=3) - - errcode = 0 - errmsg = '' - - ! Get the molar_mass that is set in the call to instantiate() - do i_elem = 1, num_constituents - call constituent_props(i_elem)%molar_mass(molar_mass_arr(i_elem), errcode, errmsg) - - if (errcode /= 0) then - write(iulog,*) " [RUN MICM]: FATAL: Unable to get molar mass at n-th element: ", & - i_elem, ". Error msg: ", errmsg - stop 3 - end if - end do - - ! TODO(jiwon) Checks the denominator is non zero - ! this code needs to go when ccpp framework does the check - do i_elem = 1, num_constituents - if (molar_mass_arr(i_elem) == 0) then - write(iulog,*) " [RUN MICM]: FATAL: molar mass cannot be zero" - errcode = 1 - stop 3 - end if - end do - - ! Convert CAM-SIMA unit to MICM unit (kg kg-1 -> mol m-3) - call convert_to_mol_per_cubic_meter(dry_air_density, molar_mass_arr, constituents) - - c_temperature = real(temperature, c_double) - c_pressure = real(pressure, c_double) - c_time_step = real(time_step, c_double) - c_constituents = real(constituents, c_double) - - write(iulog,*) " [RUN MICM]: INFO: Running MICM solver..." - do i_column = 1, num_columns - do i_layer = 1, num_layers - call micm_obj%solve(c_temperature(i_column, i_layer), c_pressure(i_column, i_layer), & - c_time_step, num_constituents, c_constituents(i_column, i_layer, :)) - end do - end do - write(iulog,*) " [RUN MICM]: INFO: MICM solver has finished." - - constituents = real(c_constituents, kind_phys) - - ! Convert MICM unit back to CAM-SIMA unit (mol m-3 -> kg kg-1) - call convert_to_mass_mixing_ratio(dry_air_density, molar_mass_arr, constituents) - - end subroutine micm_run - - !> \section arg_table_micm_final Argument Table - !! \htmlinclude micm_final.html - subroutine micm_final(iulog, errcode, errmsg) - integer, intent(in) :: iulog - integer, intent(out) :: errcode - character(len=512), intent(out) :: errmsg - - errcode = 0 - errmsg = '' - - write(iulog,*) " [FINAL MICM]: INFO: Deallocating MICM object..." - if (allocated(micm_obj)) deallocate(micm_obj) - - end subroutine micm_final - - ! Convert CAM-SIMA unit to MICM unit (kg kg-1 -> mol m-3) - subroutine convert_to_mol_per_cubic_meter(dry_air_density, molar_mass_arr, constituents) - real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3 - real(kind_phys), intent(in) :: molar_mass_arr(:) ! kg mol-1 - real(kind_phys), intent(inout) :: constituents(:,:,:) ! in: kg kg-1 | out: mol m-3 - - integer :: num_columns, num_layers, num_constituents - integer :: i_column, i_layer, i_elem - - real(c_double) :: val - - num_columns = size(constituents, dim=1) - num_layers = size(constituents, dim=2) - num_constituents = size(constituents, dim=3) - - do i_column = 1, num_columns - do i_layer = 1, num_layers - do i_elem = 1, num_constituents - val = constituents(i_column, i_layer, i_elem) * dry_air_density(i_column, i_layer) & - / molar_mass_arr(i_elem) - constituents(i_column, i_layer, i_elem) = val - end do - end do - end do - - end subroutine convert_to_mol_per_cubic_meter - - ! Convert MICM unit to CAM-SIMA unit (mol m-3 -> kg kg-1) - subroutine convert_to_mass_mixing_ratio(dry_air_density, molar_mass_arr, constituents) - real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3 - real(kind_phys), intent(in) :: molar_mass_arr(:) ! kg mol-1 - real(kind_phys), intent(inout) :: constituents(:,:,:) ! in: mol m-3 | out: kg kg-1 - - integer :: num_columns, num_layers, num_constituents - integer :: i_column, i_layer, i_elem - - real(c_double) :: val - - num_columns = size(constituents, dim=1) - num_layers = size(constituents, dim=2) - num_constituents = size(constituents, dim=3) - - do i_column = 1, num_columns - do i_layer = 1, num_layers - do i_elem = 1, num_constituents - val = constituents(i_column, i_layer, i_elem) / dry_air_density(i_column, i_layer) & - * molar_mass_arr(i_elem) - constituents(i_column, i_layer, i_elem) = val - end do - end do - end do - - end subroutine convert_to_mass_mixing_ratio + !> \section arg_table_micm_run Argument Table + !! \htmlinclude micm_run.html + subroutine micm_run(time_step, temperature, pressure, dry_air_density, constituent_props, & + constituents, iulog, errcode, errmsg) + use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t + + real(kind_phys), intent(in) :: time_step ! s + real(kind_phys), intent(in) :: temperature(:,:) ! K + real(kind_phys), intent(in) :: pressure(:,:) ! Pa + real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3 + type(ccpp_constituent_prop_ptr_t), intent(in) :: constituent_props(:) + real(kind_phys), intent(inout) :: constituents(:,:,:) ! kg kg-1 + integer, intent(in) :: iulog + integer, intent(out) :: errcode + character(len=512), intent(out) :: errmsg + + ! local variables + real(c_double) :: c_time_step + real(c_double), dimension(size(temperature, dim=1), & + size(temperature, dim=2)) :: c_temperature + real(c_double), dimension(size(pressure, dim=1), & + size(pressure, dim=2)) :: c_pressure + real(c_double), dimension(size(constituents, dim=1), & + size(constituents, dim=2), & + size(constituents, dim=3)) :: c_constituents + + real(kind_phys), dimension(size(constituents, dim=3)) :: molar_mass_arr ! kg mol-1 + + integer :: num_columns, num_layers, num_constituents + integer :: i_column, i_layer, i_elem + + num_columns = size(constituents, dim=1) + num_layers = size(constituents, dim=2) + num_constituents = size(constituents, dim=3) + + errcode = 0 + errmsg = '' + + ! Get the molar_mass that is set in the call to instantiate() + do i_elem = 1, num_constituents + call constituent_props(i_elem)%molar_mass(molar_mass_arr(i_elem), errcode, errmsg) + + if (errcode /= 0) then + errmsg = "[error] [micm] Unable to get molar mass." + write(iulog,*) trim(errmsg) + return + end if + end do + + ! TODO(jiwon) Check molar mass is non zero as it becomes a denominator for unit converison + ! this code needs to go when ccpp framework does the check + do i_elem = 1, num_constituents + if (molar_mass_arr(i_elem) == 0) then + errcode = 1 + errmsg = "[error] [micm] Molar mass must be a non zero value." + write(iulog,*) trim(errmsg) + return + end if + end do + + ! Convert CAM-SIMA unit to MICM unit (kg kg-1 -> mol m-3) + call convert_to_mol_per_cubic_meter(dry_air_density, molar_mass_arr, constituents) + + c_time_step = real(time_step, c_double) + c_temperature = real(temperature, c_double) + c_pressure = real(pressure, c_double) + c_constituents = real(constituents, c_double) + + write(iulog,*) "[info] [micm] Running MICM solver..." + + do i_column = 1, num_columns + do i_layer = 1, num_layers + call micm_obj%solve(c_temperature(i_column, i_layer), c_pressure(i_column, i_layer), & + c_time_step, num_constituents, c_constituents(i_column, i_layer, :)) + end do + end do + + write(iulog,*) "[info] [micm] MICM solver has finished." + + constituents = real(c_constituents, kind_phys) + + ! Convert MICM unit back to CAM-SIMA unit (mol m-3 -> kg kg-1) + call convert_to_mass_mixing_ratio(dry_air_density, molar_mass_arr, constituents) + + end subroutine micm_run + + !> \section arg_table_micm_final Argument Table + !! \htmlinclude micm_final.html + subroutine micm_final(iulog, errcode, errmsg) + integer, intent(in) :: iulog + integer, intent(out) :: errcode + character(len=512), intent(out) :: errmsg + + errcode = 0 + errmsg = '' + + write(iulog,*) "[debug] [micm] Deallocating MICM object..." + if (allocated(micm_obj)) deallocate(micm_obj) + + end subroutine micm_final + + ! Convert CAM-SIMA unit to MICM unit (kg kg-1 -> mol m-3) + subroutine convert_to_mol_per_cubic_meter(dry_air_density, molar_mass_arr, constituents) + real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3 + real(kind_phys), intent(in) :: molar_mass_arr(:) ! kg mol-1 + real(kind_phys), intent(inout) :: constituents(:,:,:) ! in: kg kg-1 | out: mol m-3 + + integer :: num_columns, num_layers, num_constituents + integer :: i_column, i_layer, i_elem + + real(kind_phys) :: val + + num_columns = size(constituents, dim=1) + num_layers = size(constituents, dim=2) + num_constituents = size(constituents, dim=3) + + do i_column = 1, num_columns + do i_layer = 1, num_layers + do i_elem = 1, num_constituents + val = constituents(i_column, i_layer, i_elem) * dry_air_density(i_column, i_layer) & + / molar_mass_arr(i_elem) + constituents(i_column, i_layer, i_elem) = val + end do + end do + end do + + end subroutine convert_to_mol_per_cubic_meter + + ! Convert MICM unit to CAM-SIMA unit (mol m-3 -> kg kg-1) + subroutine convert_to_mass_mixing_ratio(dry_air_density, molar_mass_arr, constituents) + real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3 + real(kind_phys), intent(in) :: molar_mass_arr(:) ! kg mol-1 + real(kind_phys), intent(inout) :: constituents(:,:,:) ! in: mol m-3 | out: kg kg-1 + + integer :: num_columns, num_layers, num_constituents + integer :: i_column, i_layer, i_elem + + real(kind_phys) :: val + + num_columns = size(constituents, dim=1) + num_layers = size(constituents, dim=2) + num_constituents = size(constituents, dim=3) + + do i_column = 1, num_columns + do i_layer = 1, num_layers + do i_elem = 1, num_constituents + val = constituents(i_column, i_layer, i_elem) / dry_air_density(i_column, i_layer) & + * molar_mass_arr(i_elem) + constituents(i_column, i_layer, i_elem) = val + end do + end do + end do + + end subroutine convert_to_mass_mixing_ratio end module micm \ No newline at end of file diff --git a/musica/micm/micm.meta b/musica/micm/micm.meta index ff066caf..a54f4071 100644 --- a/musica/micm/micm.meta +++ b/musica/micm/micm.meta @@ -52,7 +52,7 @@ dimensions = (horizontal_loop_extent,vertical_layer_dimension) intent = in [ dry_air_density ] - standard_name = density_of_dry_air + standard_name = dry_air_density units = kg m-3 dimensions = (horizontal_loop_extent, vertical_layer_dimension) type = real | kind = kind_phys