Scheduled Downtime
On Tuesday 24 October 2023 @ 5pm MT the forums will be in read only mode in preparation for the downtime. On Wednesday 25 October 2023 @ 5am MT, this website will be down for maintenance and expected to return online later in the morning.
Normal Operations
The forums are back online with normal operations. If you notice any issues or errors related to the forums, please reach out to help@ucar.edu

Time varying surface forcing from netcdf file for an idealized Southern Ocean channel model

npoumaere

Nelson Poumaere
New Member
Hello everyone,

Hope you're all doing well.

SETUP: main branch of MOM6 (downloaded and built 2 months ago from main repo), no modification to the code

I'm rather new to MOM6, and I am trying to force a Southern Ocean zonally reentrant channel model with time varying surface wind and buoyancy profiles, based on this configuration from H. Khatri and G. MacGilchrist, in which forcing profiles are time-independent.

This time-independent forcing simulation runs without problem on Archer2 (UK) computing facility.

My forcing profiles are varying monthly, so that for example, my zonal wind stress in the wind.nc (see end of post for wind.nc and input.nml) forcing file has dimensions taux(Time, lath, lonh) = taux(12,272,320).

I use a "NOLEAP" calendar, with 360 days long years. Looking into ocean_only examples "single_column" and "global", I guessed that the time evolution/periodicity of forcing is either specified through:

# 1) average_T1, average_T2 and average_DT data variables,
and/or
# 2) modulo attribute of Time coordinate.

I wrote the forcing files with both of these specified.

In MOM_input I'm specifying "VARIABLE_WINDS = True"

PROBLEM: At runtime, MOM6 directly raises the following error:

####
320 272 1 320 272 12

FATAL from PE 6: fms_io(read_data_3d_new), field taux in file INPUT/wind.nc: field size mismatch 1
####

which comes from this bit of code from fms_io.F90:

####
5421 subroutine read_data_3d_new(filename,fieldname,data,domain,timelevel, &
1 no_domain, scalar_or_1d, position, tile_count)

................

5517 if (siz_in(1) /= gxsize .or. siz_in(2) /= gysize .or. siz_in(3) /= size(data,3)) then
1 PRINT *, gxsize, gysize, size(data, 3), siz_in(1), siz_in(2), siz_in(3)
2 call mpp_error(FATAL,'fms_io(read_data_3d_new), field '//trim(fieldname)// &
3 ' in file '//trim(filename)//': field size mismatch 1')
####

QUESTION: Do you have any idea how to solve this issue? I don't see how this happens. and more generally, I'm really not sure how time periodicity/evolution of file-specified forcing fields are understood by MOM6.

Thanks in advance for any info on the matter!

All the best,
Nelson

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

ncdump -h wind.nc :

dimensions:
Time = 12 ;
lath = 272 ;
lonh = 320 ;
variables:
double Time(Time) ;
Time:_FillValue = 0. ;
Time:units = "days since 0001-01-01 00:00:00" ;
Time:calendar = "NOLEAP" ;
Time:cartesian_axis = "T" ;
Time:modulo = 360. ;
double lath(lath) ;
lath:_FillValue = 0. ;
lath:units = "degrees north" ;
lath:standard_name = "Latitude" ;
lath:cartesian_axis = "Y" ;
double lonh(lonh) ;
lonh:_FillValue = 0. ;
lonh:units = "degrees east" ;
lonh:standard_name = "Longitude" ;
lonh:cartesian_axis = "X" ;
double taux(Time, lath, lonh) ;
taux:_FillValue = -1.e+34 ;
taux:units = "N/m2" ;
taux:standard_name = "Zonal wind stress" ;
taux:missing_value = -1.e+34 ;
double tauy(Time, lath, lonh) ;
tauy:_FillValue = -1.e+34 ;
tauy:units = "N/m2" ;
tauy:standard_name = "Meridional wind stress" ;
tauy:missing_value = -1.e+34 ;
double average_T1(Time) ;
average_T1:_FillValue = 1.e+20 ;
average_T1:long_name = "Start time for average period" ;
average_T1:units = "days since 0001-01-01 00:00:00" ;
average_T1:missing_value = 1.e+20 ;
double average_T2(Time) ;
average_T2:_FillValue = 1.e+20 ;
average_T2:long_name = "End time for average period" ;
average_T2:units = "days since 0001-01-01 00:00:00" ;
average_T2:missing_value = 1.e+20 ;
double average_DT(Time) ;
average_DT:_FillValue = 1.e+20 ;
average_DT:long_name = "Length of average period" ;
average_DT:units = "days" ;
average_DT:missing_value = 1.e+20 ;
}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

cat input.nml

&MOM_input_nml
output_directory = 'OUTPUT_0/',
input_filename='n'
restart_input_dir = 'RESTART_IN/',
restart_output_dir = 'RESTART_OUT/',
parameter_filename = 'MOM_input',
'MOM_override' /

&diag_manager_nml
/

&fms_nml
domains_stack_size = 5000000,
stack_size = 0 /

&ocean_solo_nml
date_init = 1,1,1,0,0,0
calendar = 'thirty_day' /
 

npoumaere

Nelson Poumaere
New Member
Hi everyone,

Just a quick update on the situation,

1) I made all the calendar attributes coherent with the "thirty_day_month" given in input.nml, should have noticed that earlier but this didn't improve the situation.

2) Then I came to the conclusion that for temporally varying forcing files I needed to go through data_override (as in the single_column example).
I thus turned on

BUOY_CONFIG = "data_override"
and
WIND_CONFIG = "data_override"

in MOM_input and modified data_table from single_column to the one given at the end of the post (I believe this is fine).

I had to take the 'grid_spec.nc' and 'ocean_mosaic.nc' from the global example to orient MOM6 to my own ocean_hgrid.nc.
This is a bit hacky but seems to work fine.

3) There was still a problem with the dimension of taux so I modified this by adding a DEPTH variable (of length 1) so that taux now has dimensions (lonh,lath,DEPTH,Time) as seems to be the default layout required by FMS routines.

PROBLEM: Running MOM6 now gives me

%%%%
NOTE: callTree: loop Main loop, MOM_driver.F90 1
NOTE: callTree: ---> set_forcing, MOM_surface_forcing.F90
NOTE: callTree: ---> wind_forcing_by_data_override, MOM_surface_forcing.F90

FATAL: external field taux does not have an associated record dimension (REQUIRED)
[...]
%%%%

which is raised by function init_external_field in file "FMS2/time_interp/time_interp_external.F90" line 378, this being called by data_override_3d in file "./FMS2/data_override/data_override.F90", itself called by surface_forcing module.

So there's still a problem with the time layout of my data.
Note that a (Time,lath,lonh) layout leads to the same error.
My issue is that in single_column, the forcing data is indeed in a (Time,lath,lonh) layout and works perfectly fine.

Any advice on that? I'm still digging but it might be that I'm missing something completely obvious.

All the best,
Nelson

P.S.: Discovered that MOM6/FMS do not like trailing empty lines in data_table, and will crash if there are.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5

cat data_table

# gridname, fieldname_code, fieldname_file, file_name, ongrid, factor
"OCN" , "taux", "taux", "./INPUT/wind.nc", "none", 1.0
"OCN" , "tauy", "tauy", "./INPUT/wind.nc", "none", 1.0
"OCN" , "liq_precip", "liq_precip", "./INPUT/buoyancy.nc", "none", 1.0
"OCN" , "froz_precip", "froz_precip", "./INPUT/buoyancy.nc", "none", 1.0
"OCN" , "sens", "sensible", "./INPUT/buoyancy.nc", "none", 1.0
"OCN" , "lw", "", "" , "none", 0.0
"OCN" , "sw", "", "" , "none", 0.0
"OCN" , "rain", "", "" , "none", 0.0
"OCN" , "evap", "", "" , "none", 0.0
"OCN" , "snow", "", "" , "none", 0.0
"OCN" , "runoff", "", "" , "none", 0.0
"OCN" , "calving", "", "" , "none", 0.0

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

ncdump -h wind.nc

netcdf wind {
dimensions:
lonh = 320 ;
lath = 272 ;
DEPTH = 1 ;
Time = 12 ;
variables:
double lonh(lonh) ;
lonh:_FillValue = 0. ;
lonh:units = "degrees east" ;
lonh:standard_name = "Longitude" ;
lonh:cartesian_axis = "X" ;
double lath(lath) ;
lath:_FillValue = 0. ;
lath:units = "degrees north" ;
lath:standard_name = "Latitude" ;
lath:cartesian_axis = "Y" ;
double DEPTH(DEPTH) ;
DEPTH:_FillValue = 0. ;
DEPTH:units = "m" ;
DEPTH:standard_name = "Z coordinate" ;
DEPTH:cartesian_axis = "Z" ;
double Time(Time) ;
Time:_FillValue = 0. ;
Time:long_name = "Time" ;
Time:units = "days since 0001-01-01 00:00:00" ;
Time:cartesian_axis = "T" ;
Time:calendar_type = "THIRTY_DAY_MONTHS" ;
Time:calendar = "THIRTY_DAY_MONTHS" ;
Time:bounds = "Time_bnds" ;
Time:modulo = 360. ;
double taux(lonh, lath, DEPTH, Time) ;
taux:_FillValue = -1.e+34 ;
taux:units = "N/m2" ;
taux:standard_name = "Zonal wind stress" ;
taux:missing_value = -1.e+34 ;
double tauy(lonh, lath, DEPTH, Time) ;
tauy:_FillValue = -1.e+34 ;
tauy:units = "N/m2" ;
tauy:standard_name = "Meridional wind stress" ;
tauy:missing_value = -1.e+34 ;
}
 

npoumaere

Nelson Poumaere
New Member
Hi again,

Following recommendations from Theresa Morrison, I enabled the "UNLIMITED" property for the Time variable and set the Time:modulo attribute to "''. I also reverted to file-specified forcing.

The model is now running! I still need to check that everything is functioning properly, but this is very encouraging.

In the end, my files have basically the same layout as the "forcing_monthly.nc" files in the single_column example, I guess I just missed the importance of the "UNLIMITED" property, though I don't know at the moment how that solved the issue. I attach the layout of the forcing file for reference.

All the best,
Nelson

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5

ncdump -h wind.nc

netcdf wind {
dimensions:
Time = UNLIMITED ; // (12 currently)
lath = 272 ;
lonh = 320 ;
nv = 2 ;
variables:
double Time(Time) ;
Time:_FillValue = 0. ;
Time:long_name = "Time" ;
Time:units = "days since 0001-01-01 00:00:00" ;
Time:cartesian_axis = "T" ;
Time:calendar_type = "THIRTY_DAY_MONTHS" ;
Time:calendar = "THIRTY_DAY_MONTHS" ;
Time:bounds = "Time_bnds" ;
Time:modulo = "" ;
double lath(lath) ;
lath:_FillValue = 0. ;
lath:units = "degrees north" ;
lath:standard_name = "Latitude" ;
lath:cartesian_axis = "Y" ;
double lonh(lonh) ;
lonh:_FillValue = 0. ;
lonh:units = "degrees east" ;
lonh:standard_name = "Longitude" ;
lonh:cartesian_axis = "X" ;
double taux(Time, lath, lonh) ;
taux:_FillValue = -1.e+34 ;
taux:units = "N/m2" ;
taux:standard_name = "Zonal wind stress" ;
taux:missing_value = -1.e+34 ;
double tauy(Time, lath, lonh) ;
tauy:_FillValue = -1.e+34 ;
tauy:units = "N/m2" ;
tauy:standard_name = "Meridional wind stress" ;
tauy:missing_value = -1.e+34 ;
double average_T1(Time) ;
average_T1:_FillValue = 1.e+20 ;
average_T1:long_name = "Start time for average period" ;
average_T1:units = "days since 0001-01-01 00:00:00" ;
average_T1:missing_value = 1.e+20 ;
double average_T2(Time) ;
average_T2:_FillValue = 1.e+20 ;
average_T2:long_name = "End time for average period" ;
average_T2:units = "days since 0001-01-01 00:00:00" ;
average_T2:missing_value = 1.e+20 ;
double average_DT(Time) ;
average_DT:_FillValue = 1.e+20 ;
average_DT:long_name = "Length of average period" ;
average_DT:units = "days" ;
average_DT:missing_value = 1.e+20 ;
double nv(nv) ;
nv:_FillValue = NaN ;
nv:long_name = "vertex number" ;
nv:units = "none" ;
nv:cartesian_axis = "N" ;
double Time_bnds(Time, nv) ;
Time_bnds:_FillValue = 1.e+20 ;
Time_bnds:long_name = "Time axis boundaries" ;
Time_bnds:units = "days" ;
Time_bnds:missing_value = 1.e+20 ;
}
 
Top