Back to home page

MITgcm

 
 

    


File indexing completed on 2024-06-06 05:11:45 UTC

view on githubraw file Latest commit af61e5eb on 2024-06-06 03:30:35 UTC
af61e5eb16 Mart*0001 #!/usr/bin/env python
                0002 # -*- coding: iso-8859-15 -*-
                0003 ######################## -*- coding: utf-8 -*-
                0004 """Usage: gendata.py
                0005 generates extra input files of input_ad.obcs test
                0006 """
                0007 
                0008 import numpy as np
                0009 import sys
                0010 
                0011 nx=80
                0012 ny=42
                0013 nr=1
                0014 nt=1
                0015 
                0016 def writefield(fname,data):
                0017     """Call signatures::
                0018 
                0019     writefield(filename, numpy.ndarray)
                0020 
                0021     Write unblocked binary data.
                0022     """
                0023 
                0024     if sys.byteorder == 'little': data.byteswap(True)
                0025 
                0026     fid = open(fname,"wb")
                0027     data.tofile(fid)
                0028     fid.close()
                0029 
                0030     # go back to original byte ordering
                0031     if sys.byteorder == 'little': data.byteswap(True)
                0032 
                0033 def readfield(fname,dims,datatype):
                0034     """Call signatures::
                0035 
                0036     readfield(filename, dims, numpy.datatype)
                0037 
                0038     Read unblocked binary data with dimentions "dims".
                0039     """
                0040     fid = open(fname,"rb")
                0041     v   = np.fromfile(fid, datatype)
                0042     fid.close()
                0043 
                0044     if sys.byteorder == 'little': v.byteswap(True)
                0045 
                0046     if   len(v) == np.prod(dims):     v = v.reshape(dims)
                0047     elif len(v) == np.prod(dims[1:]): v = v.reshape(dims[1:])
                0048     else:
                0049         errstr = (  "dimensions do not match: \n len(data) = " + str(len(v))
                0050                   + ", but prod(dims) = " + str(np.prod(dims)) )
                0051         raise RuntimeError(errstr)
                0052 
                0053     return v
                0054 
                0055 # use default topography and modify it
                0056 fld = readfield('../input/bathy_3c.bin',[ny,nx],'float64')
                0057 
                0058 fmin = fld.min()
                0059 
                0060 bnew = np.copy(fld)
                0061 bnew[0,:17] = fmin
                0062 bnew[0,63:] = fmin
                0063 
                0064 writefield('./bathy_3c.obcs',bnew)
                0065 
                0066 # use initial conditions (of constant ocean) to derive obcs conditions
                0067 # for ocean
                0068 ufld = readfield('../input/uVel_3c0.bin',[ny,nx],'float64')
                0069 vfld = readfield('../input/vVel_3c0.bin',[ny,nx],'float64')
                0070 
                0071 uS = ufld[1,:]
                0072 vS = vfld[2,:]
                0073 uN = ufld[-1,:]
                0074 vN = vfld[-1,:]
                0075 
                0076 uW = ufld[:,1]
                0077 vW = vfld[:,2]
                0078 uE = ufld[:,-1]
                0079 vE = vfld[:,-1]
                0080 
                0081 writefield('OBSu.bin',uS)
                0082 writefield('OBSv.bin',vS)
                0083 writefield('OBNu.bin',uN)
                0084 writefield('OBNv.bin',vN)
                0085 writefield('OBEu.bin',uE)
                0086 writefield('OBEv.bin',vE)
                0087 writefield('OBWu.bin',uW)
                0088 writefield('OBWv.bin',vW)
                0089 
                0090 # initial conditions for sea ice variables
                0091 iceC0 = 1.
                0092 iceConc = np.ones((ny,nx))*iceC0
                0093 iceConc[bnew==0] = 0.
                0094 writefield('ice0_area.obcs',iceConc)
                0095 
                0096 iceH0=0.2
                0097 iceVol=iceConc*iceH0
                0098 writefield('ice0_heff.obcs',iceVol)
                0099 
                0100 wind0 = 10.
                0101 windy = np.ones((ny,nx))*wind0
                0102 windy[:,:nx//2] = -wind0
                0103 writefield('windy.bin',windy)