[622] | 1 | import dynamico.unstructured as unst |
---|
| 2 | import dynamico.xios as xios |
---|
[694] | 3 | from dynamico.meshes import radian, MPAS_Format, Unstructured_PMesh as PMesh, Local_Mesh as Mesh |
---|
[622] | 4 | |
---|
[663] | 5 | from mpi4py import MPI |
---|
| 6 | comm = MPI.COMM_WORLD |
---|
| 7 | mpi_rank, mpi_size = comm.Get_rank(), comm.Get_size() |
---|
| 8 | print '%d/%d starting'%(mpi_rank,mpi_size) |
---|
| 9 | |
---|
| 10 | import numpy as np |
---|
[682] | 11 | import cProfile |
---|
| 12 | import re |
---|
[663] | 13 | |
---|
[622] | 14 | #----------------------------- read MPAS mesh -------------------------------- |
---|
| 15 | |
---|
[682] | 16 | grid, llm, nqdyn, nqtot = 2562, 4,1,1 # 2562, 10242, 40962 |
---|
[622] | 17 | Omega, radius, g, gh0 = 2.*np.pi/86400., 6.4e6, 1., 2.94e4 |
---|
| 18 | N, T, courant = 40, 10800., 1.2 # simulation length = N*T |
---|
| 19 | print 'Omega, planetary PV', Omega, 2*Omega/gh0 |
---|
| 20 | |
---|
| 21 | def f(lon,lat): return 2*Omega*np.sin(lat) # Coriolis parameter |
---|
| 22 | print 'Reading MPAS mesh ...' |
---|
[680] | 23 | meshfile = MPAS_Format('grids/x1.%d.grid.nc'%grid) |
---|
[682] | 24 | pmesh = PMesh(comm,meshfile) |
---|
| 25 | mesh = Mesh(pmesh, llm, nqdyn, radius, f) |
---|
[622] | 26 | print '...Done' |
---|
| 27 | |
---|
[694] | 28 | #--------------------------------- write some data ---------------------------------------- |
---|
[622] | 29 | |
---|
[696] | 30 | context=xios.XIOS_Context(pmesh,mesh,nqtot, 3600) |
---|
[622] | 31 | |
---|
[694] | 32 | lat_i = radian*mesh.lat_i |
---|
[696] | 33 | lat_ik, junk = np.meshgrid(lat_i, np.arange(llm), indexing='ij') |
---|
[682] | 34 | |
---|
[622] | 35 | no_error=True |
---|
| 36 | try: |
---|
| 37 | for i in range(100): |
---|
[694] | 38 | context.update_calendar(i) |
---|
[696] | 39 | print 'send_field', i, lat_ik.shape |
---|
[694] | 40 | context.send_field_primal('ps', lat_i) |
---|
| 41 | context.send_field_primal('theta', lat_ik) |
---|
[622] | 42 | except: |
---|
| 43 | print 'An error has occured, trying to close XIOS properly' |
---|
| 44 | no_error=False |
---|
| 45 | |
---|
| 46 | print 'xios.context_finalize()' |
---|
[694] | 47 | context.finalize() |
---|
[622] | 48 | print 'xios.finalize()' |
---|
| 49 | xios.finalize() |
---|
| 50 | print 'Done' |
---|
| 51 | |
---|
| 52 | if not no_error : raise |
---|