Changeset 794 for codes/icosagcm/devel/Python
- Timestamp:
- 12/13/18 18:01:37 (6 years ago)
- Location:
- codes/icosagcm/devel/Python
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/devel/Python/dynamico/LAM.py
r791 r794 10 10 N1, N2 = self.N1, self.N2 11 11 m = zeros(c.size) 12 c1 ,c2 = c0-N1*delta, c0-(N1+N2)*delta12 c1 = c0+N2*delta 13 13 for i in range(c.size): 14 14 ci=c[i] 15 m[i] = (1.+cos((ci-c 2)*pi/(N2*delta)))/2.016 if ci < c 2: m[i]=1.15 m[i] = (1.+cos((ci-c0)*pi/(N2*delta)))/2.0 16 if ci < c0 : m[i]=1. 17 17 if ci > c1 : m[i]=0. 18 18 return m -
codes/icosagcm/devel/Python/test/py/Baroclinic_3D_ullrich.py
r792 r794 36 36 # Baroclinic instability test based on Ullrich et al. 2015, QJRMS 37 37 38 39 def baroclinic_3D(Lx,nx,Ly,ny,llm,ztop=25000.): 38 def create_pmesh(nx,ny): 39 filename = 'cart_%03d_%03d.nc'%(nx,ny) 40 logging.info('Reading Cartesian mesh ...') 41 meshfile = meshes.DYNAMICO_Format(filename) 42 pmesh = meshes.Unstructured_PMesh(comm,meshfile) 43 pmesh.partition_curvilinear(args.mpi_ni,args.mpi_nj) 44 return pmesh 45 46 def baroclinic_3D(pmesh, Lx,Ly,llm,ztop=25000.): 40 47 Rd = 287.0 # Gas constant for dryy air (j kg^-1 K^-1) 41 48 T0 = 288.0 # Reference temperature (K) … … 75 82 def coriolis(x,y): return f0+beta0*y # here y is 0 at domain center 76 83 77 filename = 'cart_%03d_%03d.nc'%(nx,ny)78 logging.info('Reading Cartesian mesh ...')79 meshfile = meshes.DYNAMICO_Format(filename)80 pmesh = meshes.Unstructured_PMesh(comm,meshfile)81 pmesh.partition_curvilinear(args.mpi_ni,args.mpi_nj)82 84 nqdyn, radius = 1, None 83 mesh 85 mesh = meshes.Local_Mesh(pmesh, llm, nqdyn, radius, coriolis) 84 86 85 87 alpha_k = (np.arange(llm) +.5)/llm … … 158 160 g, Lx, Ly = 9.81, 4e7, 6e6 159 161 nx, ny, llm = args.nx, args.ny, args.llm 160 dx,dy = Lx/nx,Ly/ny 162 dx = Lx/nx 163 dy = dx 164 165 if True: # physical domain excludes relaxation zone 166 Dy = Ly + dy*2*(args.Davies_N1+args.Davies_N2) 167 print('Dy=%f, Ly = %f'%(Dy,Ly)) 168 else: # physical domain includes relaxation zone 169 Dy=Ly 170 # print('WARNING: Davies_N1 not equal to Davies_N2; average value will be used') 171 161 172 162 173 unst.setvar('g',g) 163 164 thermo, mesh, params, flow0, gas0 = baroclinic_3D(Lx,nx,Ly,ny,llm) 174 175 pmesh = create_pmesh(nx,ny) 176 thermo, mesh, params, flow0, gas0 = baroclinic_3D(pmesh, Lx,Ly,llm) 165 177 166 178 mass_bl,mass_dak,mass_dbk = meshes.compute_hybrid_coefs(flow0[0]) … … 173 185 dt = T/nt 174 186 logging.info( 'Time step : %d x %g = %g s' % (nt,dt,nt*dt)) 187 175 188 176 189 caldyn_thermo, caldyn_eta = unst.thermo_entropy, unst.eta_mass … … 197 210 caldyn_step.mass[:,:], caldyn_step.theta_rhodz[:,:], caldyn_step.u[:,:] = m,S,u 198 211 caldyn_step.geopot[:,:], caldyn_step.W[:,:] = Phi,W 199 kappa = dx**4/432000. 212 #Relaxation zone inside the domain 213 kappa = 3.26e+15 #nonhyd=>#Y3.3#Y3.27#Y3.26#Y3.255#dblchkd;N3.25#N3.2#N3.0#Y3.4#dx**4/432000. 214 #Relaxation zone outside the domain 215 #kappa = 3.6e+15 #nonhyd=>#Y3.3#Y3.27#Y3.26#Y3.255#dblchkd;N3.25#N3.2#N3.0#Y3.4#dx**4/432000. 216 #print('Kappa=',kappa) 200 217 for i in range(nt): 201 218 caldyn_step.next()
Note: See TracChangeset
for help on using the changeset viewer.