ARCHER2 fixes from MOSRS (includes all changes)

ScFull_128.mcf

@@ -0,0 +1,161 @@
+# Global configuration
+global_configuration=global_config
+
+logging=4
+
+# Override global component defaults
+cfltest_enabled=.true.
+checkpointer_enabled=.true.
+fftsolver_enabled=.false.
+iterativesolver_enabled=.true.
+iterativesolver_single_prec_enabled=.false.
+petsc_solver_enabled=.false.
+pw_advection_enabled=.true.
+simplesetup_enabled=.true.
+smagorinsky_enabled=.true.
+lower_bc_enabled=.true.
+setfluxlook_enabled=.true. #This must be set to true if running with lower_bc
+viscosity_enabled=.true.
+diffusion_enabled=.true.
+simplecloud_enabled=.true.
+coriolis_enabled=.true.
+damping_enabled=.true.
+forcing_enabled=.true. 
+galilean_transformation=.true. # Needs debugging
+randomnoise_enabled=.true.
+mean_profiles_enabled=.true. #This must be set to true if running with damping
+th_advection_enabled=.false.
+lwrad_exponential_enabled=.true.
+iobridge_enabled=.true.
+scalar_diagnostics_enabled=.true.
+profile_diagnostics_enabled=.true.
+subgrid_profile_diagnostics_enabled=.true.
+flux_budget_enabled=.true.
+
+# Control configuration
+
+display_synopsis_frequency=1
+nn_consolidation=1
+termination_time=17500.
+dtm=0.4
+
+# IO server configuration
+ioserver_configuration_file="io/io_cfg_files/data_write_1file.xml"
+diagnostic_file="diagnostic_files/MONC_ScFull_128_diagnostics.nc"
+moncs_per_io_server=15
+sampling_frequency=20
+3d_sampling_frequency=4500
+mm=1800.0
+mm1=20.0
+diag_write_freq=7260.0
+
+# Checkpoint configuration
+checkpoint_frequency=0
+checkpoint_file="checkpoint_files/MONC_ScFull_128_dump.nc"
+
+check_walltime_frequency=20
+walltime_limit=00:05:00
+
+# Advection choices
+advection_flow_fields=pw
+advection_theta_field=tvd
+advection_q_fields=tvd
+
+# CFL configuration
+cfl_frequency=1000
+cfl_cvismax=0.4
+cfl_cvelmax=0.4
+cfl_dtmmax=0.4
+cfl_dtmmin=0.001
+
+# Simple setup configuration
+# We'll want to change this reference profile later
+thref0=288.0
+surface_pressure=100000.
+surface_reference_pressure=100000.
+x_size=128
+y_size=128
+z_size=50
+dxx=50
+dyy=50
+zztop=1600.0
+kgd=10, 30, 40, 50
+hgd=300., 700., 1000., 1600.
+nsmth=40
+rmlmax=6.9
+enable_theta=.true.
+number_q_fields=2
+use_anelastic_equations=.false.
+origional_vertical_grid_setup=.true.
+passive_th=.false.
+passive_q=.false.
+backscatter=.false.
+use_viscosity_and_diffusion=.true.
+
+# Initialization of fields
+l_init_pl_theta=.true.
+z_init_pl_theta=0.0, 600.0, 601.0, 1600.0
+f_init_pl_theta=287.2, 287.2, 297.2, 300.2
+l_init_pl_u=.true.
+z_init_pl_u=0.0, 1601.0
+f_init_pl_u=5.0, 5.0
+l_init_pl_v=.false.
+l_init_pl_q=.true.
+names_init_pl_q=vapour
+z_init_pl_q=0.0, 600.0, 601.0, 1600.0
+f_init_pl_q=9.2e-3, 9.2e-3, 6.2e-3, 6.2e-3
+
+l_matchthref=.false.
+
+# Random noise
+l_rand_pl_theta=.true.
+z_rand_pl_theta=0.0, 600.0, 601.0, 1600.
+f_rand_pl_theta=0.1, 0.1, 0.0001, 0.0001
+
+# Simple cloud
+max_height_cloud=1600.
+
+# physical constants
+z0=0.002
+z0th=0.0002
+
+# Coriolis
+fcoriol=0.00007
+geostrophic_wind_rate_of_change_in_x=0.0
+geostrophic_wind_rate_of_change_in_y=0.0
+surface_geostrophic_wind_x=5.0
+surface_geostrophic_wind_y=0.0
+
+# Damping configuration
+dmptim=0.001
+zdmp=1000.0
+hdmp=650.0
+
+# Subsidence profile
+l_subs_pl_theta=.true.
+z_subs_pl=0.0, 1600.0
+f_subs_pl=0.000005,0.000005 
+l_subs_pl_q=.true.
+
+#SUBSIDENCE=1, DIVERGENCE=0
+subsidence_input_type=0
+subsidence_local_theta=.true.
+subsidence_local_q=.true.
+
+# surface flux config
+# type_of_surface_boundary_conditions=PRESCRIBED_FLUX=0
+# type_of_surface_boundary_conditions=PRESCRIBED_SURFACE_VALUE=1
+use_surface_boundary_conditions=.true.
+use_time_varying_surface_values= .false.
+type_of_surface_boundary_conditions = 1 
+surface_boundary_input_times = 0
+surface_temperatures  = 288.2
+#surface_latent_heat_flux   = 130.052
+#surface_sensible_heat_flux = 8.04
+
+# Radiation settings - variables and switches for exponential radiation (LWEXP) 
+longwave_exp_decay = 130.0 # longwave rad exponential decay factor
+cltop_longwave_flux = 70.0 # cloud top longwave flux (Wm-2)
+clbase_longwave_flux = 0.0 # cloud base longwave flux (Wm-2)
+
+

asubmonc.sb

@@ -0,0 +1,56 @@
+#!/bin/bash
+#SBATCH --job-name=MONC_ScFull_128
+#SBATCH --output=%x.o%j
+  # %x gives job-name (SLURM_JOB_NAME)
+  # %j gives jobid (individual SLURM_JOB_ID)
+  # %A gives jobid (master     SLURM_ARRAY_JOB_ID)
+  # %a gives array task id number
+  #  https://slurm.schedmd.com/sbatch.html
+#SBATCH --open-mode=append
+#SBATCH --nodes=1
+#SBATCH --ntasks-per-node=128
+#SBATCH --ntasks=128
+#SBATCH --cpus-per-task=1
+#SBATCH --time=00:20:00
+#SBATCH --account=n02
+#SBATCH --partition=standard
+#SBATCH --qos=standard
+
+### Push to development queue
+## #SBATCH --reservation=shortqos
+## #SBATCH --qos=short
+
+# MODULES
+module load atp
+export ATP_ENABLED=1
+. /work/n02/n02/grenvill/monc/gnu-mod
+module switch gcc/10.1.0 gcc/9.3.0
+module load petsc/3.13.3
+module list
+pwd
+
+ulimit -c unlimited
+
+# set variables for submission command----------------------
+export SUBMISSION_SCRIPT_NAME=asubmonc.sb
+export MONC_EXEC=./build/bin/monc_driver.exe
+
+export TESTCASE=ScFull_128.mcf
+export STDOUT_DIR=monc_stdout
+export CP_DIR=checkpoint_files
+export RUN_NAME=${SLURM_JOB_NAME}_dump_
+export NPES=${SLURM_NTASKS}
+export MAX_CONTINUATION_RUNS=200
+# ----------------------------------------------------------
+
+echo -e "\nSubmission time: $(date)\n"
+
+. misc/continuation.sbatch.sh
+
+run_monc
+
+# output job statisitcs to .o (%x.o%j)
+echo -e "\nCompletion time: $(date)\n"
+scontrol show job $SLURM_JOB_ID
+sstat $SLURM_JOB_ID --format="AveRSS,MaxRSS"
+#  Run after batch job concludes: seff $SLURM_JOB_ID

components/buoyancy/src/buoyancy.F90

@@ -29,8 +29,6 @@ module buoyancy_mod
   real(kind=DEFAULT_PRECISION), dimension(:), allocatable :: tend_pr_tot_w
   logical :: l_tend_pr_tot_w
 
-  integer :: diagnostic_generation_frequency
-
   public buoyancy_get_descriptor
 
 contains
@@ -170,9 +168,6 @@ subroutine initialisation_callback(current_state)
       allocate( tend_pr_tot_w(current_state%local_grid%size(Z_INDEX)) )
     endif
 
-    ! Save the sampling_frequency to force diagnostic calculation on select time steps
-    diagnostic_generation_frequency=options_get_integer(current_state%options_database, "sampling_frequency")
-
   end subroutine initialisation_callback  
 
 
@@ -197,6 +192,10 @@ subroutine timestep_callback(current_state)
 
     integer :: k, n
     integer :: current_x_index, current_y_index, target_x_index, target_y_index
+    logical :: calculate_diagnostics
+
+    calculate_diagnostics = current_state%diagnostic_sample_timestep &
+                            .and. .not. current_state%halo_column
 
     current_x_index=current_state%column_local_x
     current_y_index=current_state%column_local_y
@@ -210,10 +209,8 @@ subroutine timestep_callback(current_state)
       endif
     endif  ! zero totals
 
-    if (mod(current_state%timestep, diagnostic_generation_frequency) == 0 .and. .not. current_state%halo_column) then
-      call save_precomponent_tendencies(current_state, current_x_index, current_y_index, target_x_index, target_y_index)
-    end if
-
+    if (calculate_diagnostics) &
+        call save_precomponent_tendencies(current_state, current_x_index, current_y_index, target_x_index, target_y_index)
 
 #ifdef W_ACTIVE
     if (.not. current_state%passive_th .and. current_state%th%active) then
@@ -252,9 +249,8 @@ subroutine timestep_callback(current_state)
     end if
 #endif
 
-    if (mod(current_state%timestep, diagnostic_generation_frequency) == 0 .and. .not. current_state%halo_column) then
-      call compute_component_tendencies(current_state, current_x_index, current_y_index, target_x_index, target_y_index)
-    end if
+    if (calculate_diagnostics) &
+        call compute_component_tendencies(current_state, current_x_index, current_y_index, target_x_index, target_y_index)
 
   end subroutine timestep_callback
 

components/casim/src/casim.F90

@@ -74,7 +74,15 @@ module casim_mod
      , l_pssub & ! sublimation of snow
      , l_pgsub & ! sublimation of graupel
      , l_pisub & ! sublimation of ice
-     , l_pimlt   ! ice melting
+     , l_pimlt & ! ice melting
+     ! New switches for sedimentation (these are sort-of temporary)
+     , l_gamma_online & ! when true use standard vn0.3.3 sed, when false use precalced gamma
+     , l_subseds_maxv & ! Use a CFL criteria based on max terminal velocity
+                    ! and sed_1M_2M 
+     , l_sed_eulexp & ! switch for eulexp sed based on UM. Default is false
+                      ! so standard casim sed used
+     , cfl_vt_max & ! cfl limit for sedimentation (default = 1.0)
+     , l_kfsm
 
   use micro_main, only: shipway_microphysics
   use generic_diagnostic_variables, ONLY: casdiags, allocate_diagnostic_space, &
@@ -93,7 +101,7 @@ module casim_mod
      , nc(:,:,:), qr(:,:,:), nr(:,:,:), m3r(:,:,:),rho(:,:,:) &
      , exner(:,:,:), w(:,:,:), tke(:,:,:)                               &
      , qi(:,:,:), ni(:,:,:), qs(:,:,:), ns(:,:,:), m3s(:,:,:) &
-     , qg(:,:,:), ng(:,:,:), m3g(:,:,:) 
+     , qg(:,:,:), ng(:,:,:), m3g(:,:,:), cfliq(:,:,:), cfice(:,:,:)
 
   REAL(wp), allocatable :: AccumSolMass(:,:,:), AccumSolNumber(:,:,:) ! Accumulation mode aerosol
   REAL(wp), allocatable :: ActiveSolLiquid(:,:,:)                      ! Activated aerosol
@@ -277,6 +285,8 @@ subroutine initialisation_callback(current_state)
     allocate(qg(kte,1,1))
     allocate(ng(kte,1,1))
     allocate(m3g(kte,1,1))
+    allocate(cfliq(kte,1,1))
+    allocate(cfice(kte,1,1))
 
     allocate(AccumSolMass(kte,1,1))
     allocate(AccumSolNumber(kte,1,1))
@@ -581,6 +591,7 @@ subroutine timestep_callback(current_state)
       iqx = iql
       qc(:,1,1) = current_state%zq(iqx)%data(:,jcol,icol)
       dqc(:,1,1) = current_state%sq(iqx)%data(:,jcol,icol)
+      cfliq(:,1,1) = 1.0
     end IF
     IF (nq_r > 0)then
       iqx = iqr
@@ -608,6 +619,7 @@ subroutine timestep_callback(current_state)
       iqx = iqi
       qi(:,1,1) = current_state%zq(iqx)%data(:,jcol,icol)
       dqi(:,1,1) = current_state%sq(iqx)%data(:,jcol,icol)
+      cfice(:,1,1) = 1.0
     end IF
     IF (nq_s > 0)then
       iqx = iqs
@@ -712,7 +724,7 @@ subroutine timestep_callback(current_state)
        pressure, rho,                              &
        w, tke,                                     &
        z_half, z_centre,                           &
-       dz,                                         &
+       dz, cfliq, cfice,                           &
                                 ! in/out
        dqv, dqc, dqr, dnc, dnr, dm3r,              &
        dqi, dqs, dqg, dni, dns, dng, dm3s, dm3g,   &
@@ -838,7 +850,7 @@ subroutine timestep_callback(current_state)
     ! and surface
     ! snow rate (precip_s), which is the sum of ice, snow and graupel (See micromain.F90 in casim for
     ! calculation).
-    if (l_warm) then
+    if (l_warm .or. .not. casdiags % l_surface_snow ) then
        surface_precip(target_y_index,target_x_index) = &
             casdiags % SurfaceRainR(1,1)
     else
@@ -947,6 +959,11 @@ subroutine read_configuration(current_state)
     l_pgsub         = options_get_logical(current_state%options_database, 'l_pgsub')
     l_pisub         = options_get_logical(current_state%options_database, 'l_pisub')
     l_pimlt         = options_get_logical(current_state%options_database, 'l_pimlt')
+    l_gamma_online  = options_get_logical(current_state%options_database, 'l_gamma_online')
+    l_subseds_maxv  = options_get_logical(current_state%options_database, 'l_subseds_maxv')
+    l_sed_eulexp  = options_get_logical(current_state%options_database, 'l_sed_eulexp')
+    cfl_vt_max      = options_get_real(current_state%options_database, 'cfl_vt_max')
+    l_kfsm          = options_get_logical(current_state%options_database, 'l_kfsm')
 
   end subroutine read_configuration
 

components/casim/src/casim_monc_diagnostics/casim_monc_dgs_space.F90

@@ -354,9 +354,12 @@ subroutine populate_casim_monc_dg(current_state, casdiags )
          casdiags % dqr(1,1,:)
 
     if (.not. l_warm) then
-       if ( casdiags % l_precip ) & 
-            casim_monc_dgs % precip(target_y_index,target_x_index) = &
+       if ( casdiags % l_precip .and. casdiags % l_surface_snow ) & 
+         casim_monc_dgs % precip(target_y_index,target_x_index) = &
             casdiags % SurfaceRainR(1,1) + casdiags % SurfaceSnowR(1,1)
+      if ( casdiags % l_precip .and. .not. casdiags % l_surface_snow ) & 
+         casim_monc_dgs % precip(target_y_index,target_x_index) = &
+            casdiags % SurfaceRainR(1,1)
        if ( casdiags % l_surface_snow ) & 
             casim_monc_dgs % SurfaceSnowR(target_y_index,target_x_index) = &
             casdiags % SurfaceSnowR(1,1)