[MITgcm-support] Problem creating Binary Forcing Files for EXF

[Stephen Hunter [EAR]]

Hi Stanislav and Martin,

Thank you both for your quick replies.

Stanislav - Regarding the saving of the data to binary. When I first tried running my model setup I wrote the array in one go similar to the way that you suggest but hit the same problem. As a sanity check I've just run the model with new forcing files generated with the array written in one go and I can confirm that I get the same model behaviour.

I forgot to mention that I am starting the ocean with zero moment (I only define initial temperature and salinity), and so the 1500 years is so that my simulation approaches a state of equilibrium. I look at the final 50 yr of climatology just to get a snapshot of the ocean state. I have begun testing with just running the model from rest for 20 year, and my problem is present by then, so I can confirm that it persists through the entire length of the integration.

Martin- thanks for the heads up regarding these diagnostics. If I implement some of these into my model I should be able to get some insight into my problem.  I will give this a go.
I’ll keep the forum updated

Cheers

Stephen











On 18 May 2020, at 13:31, Martin Losch <Martin.Losch at awi.de<mailto:Martin.Losch at awi.de>> wrote:

Hi Stephen,

Stanislav’s description of how to generate time dependent 2D forcing fields to be used and interpolated by pkg/exf is correct. I just want to add, that sometimes the forcing products (re-analysis data) have a different ordering with (i,j)=(1,1) at the NW corner, but MITgcm/pkg/exf expects them at the SW corner of the domain.

If you don’t want to use exf-interpolation, then you do the same as for bathymetry (I guess a (ni=6*96,nj=96) field) for each of nt time slices and also save them like that, i.e. before saving your field would have these dimension: interpolated_field(1:ni,1:nj,1:nt), or if that’s too large for your memory you can also write 2D slices one after another to the same file. For the vector data (taux, tauy) it’s far more tricky, because the directions on the cube faces change from face to face, so using exf-interpolation is recommended.

Your second question:
2) Currently I have MITgcm output monthly means of ocean variables but to help me understand my problem,  is it possible to have MITgcm write out the monthly averages of the Input Forcing Fields? (or at least the surface wind vectors)
Easy with the diagnostics package: in availabe_diagnostics.log, there’s a long list of EXF related fields that you can save at any frequency you like, e.g. I have this for specific configuration that uses EXF:

------------------------------------------------------------------------------------
  Num  |<-Name->|Levs|  mate |<- code ->|<--  Units   -->|<- Tile (max=80c)
------------------------------------------------------------------------------------
   199 |EXFhs   |  1 |       |SM      U1|W/m^2           |Sensible heat flux into ocean, >0 increases theta
   200 |EXFhl   |  1 |       |SM      U1|W/m^2           |Latent heat flux into ocean, >0 increases theta
   201 |EXFlwnet|  1 |       |SM      U1|W/m^2           |Net upward longwave radiation, >0 decreases theta
   202 |EXFswnet|  1 |       |SM      U1|W/m^2           |Net upward shortwave radiation, >0 decreases theta
   203 |EXFlwdn |  1 |       |SM      U1|W/m^2           |Downward longwave radiation, >0 increases theta
   204 |EXFswdn |  1 |       |SM      U1|W/m^2           |Downward shortwave radiation, >0 increases theta
   205 |EXFqnet |  1 |       |SM      U1|W/m^2           |Net upward heat flux (turb+rad), >0 decreases theta
   206 |EXFtaux |  1 |       |UU      U1|N/m^2           |zonal surface wind stress, >0 increases uVel
   207 |EXFtauy |  1 |       |VV      U1|N/m^2           |meridional surface wind stress, >0 increases vVel
   208 |EXFuwind|  1 |       |UM      U1|m/s             |zonal 10-m wind speed, >0 eastward
   209 |EXFvwind|  1 |       |VM      U1|m/s             |meridional 10-m wind speed, >0 northward
   210 |EXFwspee|  1 |       |SM      U1|m/s             |10-m wind speed modulus ( >= 0 )
   211 |EXFatemp|  1 |       |SM      U1|degK            |surface (2-m) air temperature
   212 |EXFaqh  |  1 |       |SM      U1|kg/kg           |surface (2-m) specific humidity
   213 |EXFevap |  1 |       |SM      U1|m/s             |evaporation, > 0 increases salinity
   214 |EXFpreci|  1 |       |SM      U1|m/s             |precipitation, > 0 decreases salinity
   215 |EXFsnow |  1 |       |SM      U1|m/s             |snow precipitation, > 0 decreases salinity
   216 |EXFempmr|  1 |       |SM      U1|m/s             |net upward freshwater flux, > 0 increases salinity
   217 |EXFpress|  1 |       |SM      U1|N/m^2           |atmospheric pressure field
   218 |EXFroff |  1 |       |SM      U1|m/s             |river runoff, > 0 decreases salinity
   219 |EXFroft |  1 |       |SM      U1|deg C           |river runoff temperature
   220 |EXFsalfx|  1 |       |SM      U1|psu.kg/m^2/s    |upward salt flux, > 0 decreases salinity


I find this very useful in debugging forcing fields.

I am wondering why your problem only appears in the last 50 years of 1500 year run, if you forcing is cycle and repeats every 20 years. It’s not clear what should go wrong after 1500/20 = 75 successful cycles where you don’t see any problems.

Martin


On 18. May 2020, at 13:16, Stanislav Martyanov <martyanov.sd at gmail.com<mailto:martyanov.sd at gmail.com>> wrote:

Hello, Stephen!

When I create my EXF input data, I just write 2D XY fields of meteo forcing for each required variable sequentially for each time-step. Say, I have an array of 2m air temperature 't2m' prepared and stored in the matlab's current workspace and I want to write it into a binary file for use in MITgcm exf package. If I use daily-mean fields and 1deg horizontal resolution of the forcing data, I have the array of 360x180x366 elements for an ordinary year and 360x180x367 for a leap year. The additional record (366th or 367th) will be used for the internal time-interpolation implemented in the model when the model's time step in less than 24 hours. When such array is correctly prepared in matlab, I write it:

fout = fopen('exf_hourly_t2m.bin','wb','ieee-be');
fwrite(fout,t2m,'float32');
fclose(fout);

And so on for each variable..

Also, the corresponding parameters should be correctly set in data.exf and EXF_OPTIONS.h if you use spatial interpolation.

Regards,
Stanislav Martyanov


пн, 18 мая 2020 г. в 13:44, Stephen Hunter [EAR] <S.Hunter at leeds.ac.uk<mailto:S.Hunter at leeds.ac.uk>>:

Firstly my apologies for the long query (my 2 questions are summarised at the bottom)

I have set up a CS96 model for a deep-time palaeo experiment. I am using 20 years of monthly surface climatology in a cyclical manner to run the model for 1500 years. The surface climatology is from a global high-resolution atmosphere-only model. I am using Matlab to generate the MITgcm boundary conditions (bathymetry) and forcing files.  I believe that I have MITgcm correctly set-up as I can run the model from these boundary conditions and forcing files for 1500 years without problems. Globally integrated temperature on ocean levels comes into equilibrium and I see a drift in globally integrated salinity (as I expect as I haven’t implemented a river scheme yet).

I generate my bathymetry binary on a CS96 grid using the lat-lon specification within the grid definitions. I am satisfied that this works as I expect because I can use the reverse of my regridding process to unwrap the MITgcm climatological average back to a regular lat-lon grid and the imprinted land sea mask comes out as expected.

Whilst the model is happy to run for 1500 years, the processed climatology from the final 50 years appears incorrect. Whilst the surface ocean temperature appears visually OK (latitudinal distribution of temperatures appears fine). It is when I look closer, things appear to be wrong. For example when I look at the surface ocean currents I have a lot of diverging currents where I would expect zonal equatorial currents. I also don’t have flow through ocean gateways and channels. I do get some western boundary currents and some upwelling, but these are inconsistent in some regions.

I suspect that this mishmash is because I am creating the Forcing binaries incorrectly.  I would therefore like to check if I am creating them correctly. Here is a simplified account of what I am doing….

For each climatological field I generate a binary file for each year of monthly data. I do this in Matlab by building a vector of the data and then saving this as a float64 binary.  Consider a year of climatological data is an array of size (360,180,12), in which longitude runs from -179.5 to 179.5 in 1-degree cells, and latitude runs from 89.5 to -89.5 in 1-degree cells. I then create a vector from this array with the following recipe (simplified pseudo code):

vec=zeros(360*180*12,1);
counter=1;
for month=1:1:12
  for lat=180:-1:1
   for lon=1:1:360
    vec(counter,1)=data(lon,lat,month)
    counter=counter+1
  end
end
end

I then write vec to a 64 bit binary. With regards to CPP options for EXF I am using configuration A1,B3 but currently without runoff. In data.exf I define Ustress, Vstress, atemp, aqh, precip, wspeed, swdown and lwdown to be read in from yearly binary files with lon0 and lat0 as -179.5 and -89.5. My non scalar fields are on slightly different grid but these are specified accordingly.

I cannot find a description of the ordering of cell data within binary forcing files for the EXF package. I have tried writing array data directly to binary file, transposing my input data, changing the ordering that I fill vec, but these appear wrong also. I have tried to write my forcing data to the CS96 grid (like I do for bathymetry) but I didn’t know how to order the data wrt to the time dimension. I would prefer to have EXF read and interpolate my lat-lon climatology.


My questions are therefore:

1) What is the ordering of the data within binary forcing files that EXF reads?
2) Currently I have MITgcm output monthly means of ocean variables but to help me understand my problem,  is it possible to have MITgcm write out the monthly averages of the Input Forcing Fields? (or at least the surface wind vectors)


Kind regards
Stephen




_______________________________________________
MITgcm-support mailing list
MITgcm-support at mitgcm.org<mailto:MITgcm-support at mitgcm.org>
http://mailman.mitgcm.org/mailman/listinfo/mitgcm-support
_______________________________________________
MITgcm-support mailing list
MITgcm-support at mitgcm.org
http://mailman.mitgcm.org/mailman/listinfo/mitgcm-support

_______________________________________________
MITgcm-support mailing list
MITgcm-support at mitgcm.org<mailto:MITgcm-support at mitgcm.org>
http://mailman.mitgcm.org/mailman/listinfo/mitgcm-support

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.mitgcm.org/pipermail/mitgcm-support/attachments/20200518/360e2c12/attachment-0001.html>