Held-Suarez setup with input file

[jmartinezclaros]

Description

(This could be a bug.) I have been setting up a Held-Suarez experiment on Isca, using an input file that contains a local heating source (K/s units). I ran the model with two setups: one with input file, and a control setup without any input file. Then I calculated the difference between the zonal-mean time-averaged zonal winds of both setups. The difference is exactly ZERO. I am not sure if this an accurate output, but it likely isn't. The internal Fortran script (interpolator.F90) does see the input file (any small change in the file creates an array of familiar fatal error in the model). However, I suspect that my script is just interpolating the file and then doing nothing afterwards. I am including the python script. No change is done on either hs_forcing.F90 or interpolator.F90, which are the main Fortran programs it runs under.

Screenshots

REPRODUCIBLE EXAMPLE

Very little is changed from the Held-Suarez test case. Notable changes are the call to

input file option in the namelist parameters, and telling Isca where the file is

import numpy as np
import time
import os
import sys
import glob

from isca import DryCodeBase, DiagTable, Experiment, Namelist, GFDL_BASE

NCORES = 16
RESOLUTION = 'T85', 25

cb = DryCodeBase.from_directory(GFDL_BASE)

cb.compile()

create an Experiment object to handle the configuration of model parameters

and output diagnostics

exp_name = 'jets_ctrl_bam'
exp = Experiment(exp_name, codebase=cb)

#Tell model how to write diagnostics
diag = DiagTable()
diag.add_file('ctrl_bam_local_heating_added', 30, 'days', time_units='days')

#Tell model which diagnostics to write
diag.add_field('dynamics', 'ps', time_avg=True)
diag.add_field('dynamics', 'bk')
diag.add_field('dynamics', 'pk')
diag.add_field('dynamics', 'ucomp', time_avg=True)
diag.add_field('dynamics', 'vcomp', time_avg=True)
diag.add_field('dynamics', 'temp', time_avg=True)
diag.add_field('dynamics', 'vor', time_avg=True)
diag.add_field('dynamics', 'div', time_avg=True)

exp.diag_table = diag

#local_heating='from_file'
#local_heating_file='a_file_I_just_invented'

define namelist values as python dictionary

wrapped as a namelist object.

namelist = Namelist({
'main_nml': {
'dt_atmos': 600,
'days': 3650,
'calendar': 'thirty_day',
'current_date': [2000,1,1,0,0,0]
},

'atmosphere_nml': {
    'idealized_moist_model': False  # False for Newtonian Cooling.  True for Isca/Frierson
},

'spectral_dynamics_nml': {
    'damping_order'           : 4,                      # default: 2
    'water_correction_limit'  : 200.e2,                 # default: 0
    'reference_sea_level_press': 1.0e5,                  # default: 101325
    'valid_range_t'           : [100., 800.],           # default: (100, 500)
    'initial_sphum'           : 0.0,                  # default: 0
    'vert_coord_option'       : 'uneven_sigma',         # default: 'even_sigma'
    'scale_heights': 6.0,
    'exponent': 7.5,
    'surf_res': 0.5
},

'hs_forcing_nml': {
    'local_heating_option':'from_file',
    'local_heating_file':'local_heating',
    't_zero': 315.,    # temperature at reference pressure at equator (default 315K)
    't_strat': 200.,   # stratosphere temperature (default 200K)
    'delh': 60.,       # equator-pole temp gradient (default 60K)
    'delv': 10.,       # lapse rate (default 10K)
    'eps': 0.,         # stratospheric latitudinal variation (default 0K)
    'sigma_b': 0.7,    # boundary layer friction height (default p/ps = sigma = 0.7)

    # negative sign is a flag indicating that the units are days
    'ka':   -40.,      # Constant Newtonian cooling timescale (default 40 days)
    'ks':    -4.,      # Boundary layer dependent cooling timescale (default 4 days)
    'kf':   -1.,       # BL momentum frictional timescale (default 1 days)

    'do_conserve_energy':   True,  # convert dissipated momentum into heat (default True)
},

'diag_manager_nml': {
    'mix_snapshot_average_fields': False
},

'fms_nml': {
    'domains_stack_size': 600000                        # default: 0
},

'fms_io_nml': {
    'threading_write': 'single',                         # default: multi
    'fileset_write': 'single',                           # default: multi
}

})

exp.namelist = namelist
exp.set_resolution(*RESOLUTION)
exp.inputfiles = [os.path.join('/data/jmartinezclaros/local_heating.nc')]

#Run
if name == 'main':
st = time.time()
exp.run(1, num_cores=NCORES, use_restart=False)
for i in range(2, 2):
exp.run(i, num_cores=NCORES) # use the restart i-1 by default
elapsed_time = time.time() - st
print('Execution time:', time.strftime("%H:%M:%S", time.gmtime(elapsed_time)))

END OF EXAMPLE

########################################################################################