ExeClim / Isca

Idealized GCM from the University of Exeter
https://execlim.github.io/IscaWebsite
GNU General Public License v3.0
139 stars 125 forks source link

Held-Suarez setup with input file #269

Open jmartinezclaros opened 1 month ago

jmartinezclaros commented 1 month ago

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

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