New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
make_INITICE.py in NEMO/trunk/tests/ICE_ADV2D/EXPREF – NEMO

source: NEMO/trunk/tests/ICE_ADV2D/EXPREF/make_INITICE.py @ 10515

Last change on this file since 10515 was 10515, checked in by clem, 5 years ago

add a demo configuration for a diagonal advection of a patch of ice in a square basin

  • Property svn:executable set to *
File size: 3.2 KB
RevLine 
[10515]1#!/usr/bin/python
2
3import os,sys
4from netCDF4 import Dataset as netcdf
5import numpy as np
6import matplotlib.pyplot as plt
7from math import exp
8from math import ceil
9
10resname=''
11
12# input file
13fcoord='mesh_mask.nc'
14
15# output file
16fflx='initice.nc'
17
18print '   creating init ice file  ' +fflx
19
20# Reading coordinates file
21nccoord=netcdf(fcoord,'r')
22nav_lon=nccoord.variables['nav_lon']
23nav_lat=nccoord.variables['nav_lat']
24time_counter=1
25LON1= nav_lon.shape[1]
26LAT1= nav_lon.shape[0]
27print 'nav_lon.shape[1]' ,nav_lon.shape[1]
28print 'LON1 ', LON1
29print 'LAT1 ', LAT1
30
31# Creating INITICE netcdf file
32nc=netcdf(fflx,'w')
33nc.createDimension('y',LAT1)
34nc.createDimension('x',LON1)
35nc.createDimension('time_counter',None)    # Setting dimension size to 0 or None makes it unlimited.
36
37cdflon=nc.createVariable('nav_lon','f',('y','x'))
38cdflat=nc.createVariable('nav_lat','f',('y','x'))
39cdftimecounter=nc.createVariable('time_counter','f',('time_counter'))
40
41# ati : Fraction of open waters in sea ice - units %
42# hti : Sea ice thickness - units m
43# hts : Snow thickness - units m
44# smi : Sea ice salinity:
45# tmi : Sea ice internal temperature - units K
46# tsu : Sea ice surface temperature - units K
47#
48# Take constant values from namelist &namiceini of NEMO
49rn_hti_ini=2.0
50rn_hts_ini=0.2            #  initial real snow thickness (m)
51rn_ati_ini=0.9            #  initial ice concentration   (-)
52rn_smi_ini=6.3            #  initial ice salinity     (g/kg)
53rn_tmi_ini=270.           #  initial ice/snw temperature (K)
54rn_tsu_ini=270.           #  initial sea ice temperature (K)
55#
56cdfati=nc.createVariable('ati','f',('time_counter','y','x'))
57cdfati.units='Percentage'
58cdfati.long_name='Sea ice concentration'
59cdfhti=nc.createVariable('hti','f',('time_counter','y','x'))
60cdfhti.long_name='Sea ice thickness'
61cdfhti.units='m'
62cdfhts=nc.createVariable('hts','f',('time_counter','y','x'))
63cdfhts.long_name='Snow thickness'
64cdfhts.units='m'
65cdfsmi=nc.createVariable('smi','f',('time_counter','y','x'))
66cdfsmi.long_name='Sea ice salinity'
67cdfsmi.units='pss'
68cdftmi=nc.createVariable('tmi','f',('time_counter','y','x'))
69cdftmi.long_name='Sea ice internal temperature'
70cdftmi.units='Kelvin'
71cdftsu=nc.createVariable('tsu','f',('time_counter','y','x'))
72cdftsu.long_name='Sea ice surface temperature'
73cdftsu.units='Kelvin'
74
75cdflon[:,:]=nav_lon[:,:]
76cdflat[:,:]=nav_lat[:,:]
77cdftimecounter[0]=1
78
79# Fill fields
80#print 'cdfati[:,1]', cdfati[:,1]  -> 32 values
81
82# Add a gaussian for sea ice thickness here
83cdfhti[:,:,:]=0.
84cdfhts[:,:,:]=0.
85cdfati[:,:,:]=0.
86cdfsmi[:,:,:]=0.
87cdftmi[:,:,:]=rn_tmi_ini
88cdftsu[:,:,:]=rn_tsu_ini
89
90# --------------------------------------
91# for basin=99x99km with dx=1km ; dy=3km
92#sigx=-0.04
93#sigy=-0.04*9.
94#xshift=50.-1.
95#yshift=17.-1.
96#dlat=7
97#dlon=21
98
99# --------------------------------------
100# for basin=99x99km with dx=3km ; dy=3km
101sigx=-0.04
102sigy=-0.04
103xshift=50.-1.
104yshift=50.-1.
105dlat=21
106dlon=21
107
108# --- gaussian and square experiment ---
109for y in np.arange(dlat,LAT1-dlat,1) :
110    for x in np.arange(dlon,LON1-dlon,1) :
111        cdfhti[:,y,x] = rn_hti_ini*exp(sigx*(x-xshift)**2)*exp(sigy*(y-yshift)**2)
112        cdfhts[:,y,x] = rn_hts_ini*exp(sigx*(x-xshift)**2)*exp(sigy*(y-yshift)**2)
113        cdfati[:,y,x] = rn_ati_ini
114        cdfsmi[:,y,x] = rn_smi_ini
115
116nc.close()
117nccoord.close()
118
119#sys.exit()
Note: See TracBrowser for help on using the repository browser.