problem with put_att (maximum number of attributes)
Rob Ross
rross at mcs.anl.gov
Tue Feb 24 14:19:15 CST 2009
I'm not familiar enough with Fortran compilation to understand what
went amiss. -- Rob
On Feb 24, 2009, at 11:30 AM, Wei-keng Liao wrote:
> Lex,
> I think it may be the Fortran compiler you used.
> I am using gfortran, gcc version 4.1.2 and have no problem even if
> I used (/istart,jstart/),(/nx_mpi,ny_mpi/) in the arguments.
>
> Rob, can you explain why this is happening? I wonder if it is because
> our Fortran interfaces could not handle such arguments.
>
> Lex got the MPI error (see below) when he used
>
> nfmpi_put_vara_real(NCIDP,array_varIDP,(/istart,jstart/),(/
> nx_mpi,ny_mpi/), buf)
>
> but no error if he used array arguments, starts(2) and counts(2)
>
> integer(MPI_OFFSET_KIND) :: starts(2), counts(2)
> ....
> starts = (/istart,jstart/)
> counts = (/nx_mpi,ny_mpi/)
> nfmpi_put_vara_real_all(NCIDP,array_varIDP,starts,counts,buf)
>
> Wei-keng
>
>
> On Feb 24, 2009, at 11:15 AM, Alex Hoffmann wrote:
>
>> Yes, this works, also if I replace in the old version the vectors for
>> start and count by corresponding vector variables. Thank you.
>> Is this a function interface issue?
>> Lex
>>
>> Wei-keng Liao wrote:
>>> Hi, Lex,
>>>
>>> Can you try the attached program? I removed the serial netCDF part.
>>> Since I don't have any problem of running it, I would suggest you to
>>> run mpiexec in the debug mode and find the exact error location in
>>> pnetCDF.
>>> Please let us know.
>>>
>>> Wei-keng
>>>
>>>
>>>
>>>
>>>
>>>
>>> On Feb 24, 2009, at 3:33 AM, Alex Hoffmann wrote:
>>>
>>>> Dear Wei-keng,
>>>>
>>>> The following are the MPI errors I get, also if I subdivide the
>>>> line
>>>> (surely it can't be necessary to divide lines, though?). I never
>>>> had
>>>> any compilation problems...
>>>>
>>>> The job aborts at the nfmpi_put_vara_real_all() call, so my serial
>>>> netCDF subfiles are still assigned and the parallel netCDF file
>>>> still
>>>> created.
>>>>
>>>> Cheers,
>>>>
>>>> Lex
>>>>
>>>> aborting job:
>>>> Fatal error in MPI_Waitall: Other MPI error, error stack:
>>>> MPI_Waitall(273): MPI_Waitall(count=5, req_array=0xa11220,
>>>> status_array=0xa1dc70) failed
>>>> MPIDI_CH3_Progress_wait(207): an error occurred while handling an
>>>> event
>>>> returned by MPIDU_Sock_Wait()
>>>> MPIDI_CH3I_Progress_handle_sock_event(540): ch3|sock|immedread
>>>> 0xa07e50
>>>> 0xa116e0 0xa0cd70
>>>> MPIDU_Sock_readv(331): the supplied buffer contains invalid memory
>>>> (set=0,sock=3,errno=14:Bad address)
>>>> aborting job:
>>>> Fatal error in MPI_Alltoall: Other MPI error, error stack:
>>>> MPI_Alltoall(814): MPI_Alltoall(sbuf=0xa1b180, scount=1, MPI_INT,
>>>> rbuf=0xa1b160, rcount=1, MPI_INT, comm=0x84000001) failed
>>>> MPI_Waitall(273): MPI_Waitall(count=8, req_array=0xa1b240,
>>>> status_array=0xa1b270) failed
>>>> MPIDI_CH3_Progress_wait(207): an error occurred while handling an
>>>> event
>>>> returned by MPIDU_Sock_Wait()
>>>> MPIDI_CH3I_Progress_handle_sock_event(492):
>>>> connection_recv_fail(1728):
>>>> MPIDU_Socki_handle_read(590): connection closed by peer
>>>> (set=0,sock=1)
>>>> aborting job:
>>>> Fatal error in MPI_Allreduce: Other MPI error, error stack:
>>>> MPI_Allreduce(1058): MPI_Allreduce(sbuf=0x7fff2f18140c,
>>>> rbuf=0x7fff2f181408, count=1, MPI_INT, MPI_MAX, comm=0x84000000)
>>>> failed
>>>> MPIR_Allreduce(317):
>>>> MPIC_Sendrecv(161):
>>>> MPIC_Wait(308):
>>>> MPIDI_CH3_Progress_wait(207): an error occurred while handling an
>>>> event
>>>> returned by MPIDU_Sock_Wait()
>>>> MPIDI_CH3I_Progress_handle_sock_event(492):
>>>> connection_recv_fail(1728):
>>>> MPIDU_Socki_handle_read(614): connection failure
>>>> (set=0,sock=2,errno=104:Connection reset by peer)
>>>> rank 0 in job 37 Herzog131_50855 caused collective abort of
>>>> all ranks
>>>> exit status of rank 0: return code 13
>>>>
>>>>
>>>> Wei-keng Liao wrote:
>>>>> Hi, Lex,
>>>>>
>>>>> I got a compile error of "Syntax error in array constructor"
>>>>> from the
>>>>> line calling
>>>>> nfmpi_put_vara_real_all(). But after divide this line into two
>>>>> lines:
>>>>> call
>>>>> handle_err_nc_p(nfmpi_put_vara_real_all(NCIDP,array_varIDP,(/
>>>>> istart,jstart/),
>>>>>
>>>>> &
>>>>>
>>>>> (/nx_mpi,ny_mpi/),reshape(array,(/nx_mpi,ny_mpi/))),10)
>>>>>
>>>>> Everything compiles, runs OK.
>>>>> Please tell us what MPI errors you have got.
>>>>>
>>>>> The followings is from the command-line output and ncmpidump:
>>>>>
>>>>> ----------------------------------------------------------------------------------------
>>>>>
>>>>>
>>>>> % mpiexec -n 4 mpi_setup_netCDFout
>>>>> Myrank 0: 0
>>>>> MPI initialized!
>>>>> Number of processes [nprocs]: 4
>>>>> Coords of process 0 : 0 0
>>>>> whole 4*4 array on each process:
>>>>> 11.00000 12.00000 13.00000 14.00000
>>>>> 21.00000 22.00000 23.00000 24.00000
>>>>> 31.00000 32.00000 33.00000 34.00000
>>>>> 41.00000 42.00000 43.00000 44.00000
>>>>>
>>>>> subgrid array on process with coords 0 0 :
>>>>> shape: 2 2
>>>>> 11.00000 12.00000
>>>>> 21.00000 22.00000
>>>>>
>>>>>
>>>>> array in nc_out_p subroutine
>>>>> 11.00000 21.00000 12.00000 22.00000
>>>>>
>>>>> (/nx,ny/) 2 2
>>>>> (/istart,jstart/) 1 1
>>>>>
>>>>> pnetCDF var dimensions in xy
>>>>> 4 4
>>>>> ntx_mpi 4
>>>>> ntx 4
>>>>> nx_mpi 2
>>>>> nx 2
>>>>>
>>>>>
>>>>> -------------------------------------------------------------------------------------------------------
>>>>>
>>>>>
>>>>> % ncmpidump out_glb.nc
>>>>> netcdf out_glb {
>>>>> dimensions:
>>>>> x = 4 ;
>>>>> y = 4 ;
>>>>> variables:
>>>>> float x(x) ;
>>>>> float y(y) ;
>>>>> float array(y, x) ;
>>>>> data:
>>>>>
>>>>> x = 0, 0, 0, 0 ;
>>>>>
>>>>> y = 0, 0, 0, 0 ;
>>>>>
>>>>> array =
>>>>> 11, 21, 31, 41,
>>>>> 12, 22, 32, 42,
>>>>> 13, 23, 33, 43,
>>>>> 14, 24, 34, 44 ;
>>>>> }
>>>>>
>>>>> --------------------------------------------------------------------------------------------------------
>>>>>
>>>>>
>>>>>
>>>>> On Feb 23, 2009, at 5:57 AM, Alex Hoffmann wrote:
>>>>>
>>>>>> Dear Wei-keng,
>>>>>>
>>>>>> below is the script section in collective data mode:
>>>>>>
>>>>>> -------
>>>>>> ! this should have been possible in collective data mode, but did
>>>>>> not seem to work !!!
>>>>>> ! call handle_err_nc_p(nfmpi_begin_indep_data(NCIDP),9)
>>>>>>
>>>>>> ! call
>>>>>> handle_err_nc_p(nfmpi_put_vara_real(NCIDP,array_varIDP,(/
>>>>>> istart,jstart/),(/nx_mpi,ny_mpi/),reshape(array,(/
>>>>>> nx_mpi,ny_mpi/))),10)
>>>>>>
>>>>>>
>>>>>>
>>>>>> call
>>>>>> handle_err_nc_p(nfmpi_put_vara_real_all(NCIDP,array_varIDP,(/
>>>>>> istart,jstart/),(/nx_mpi,ny_mpi/),reshape(array,(/
>>>>>> nx_mpi,ny_mpi/))),10)
>>>>>>
>>>>>>
>>>>>>
>>>>>> ! call handle_err_nc_p(nfmpi_end_indep_data(NCIDP),11)
>>>>>>
>>>>>> call handle_err_nc_p(nfmpi_close(NCIDP),12)
>>>>>> -------
>>>>>>
>>>>>> Basically, I replaced the independent call by the collective
>>>>>> call for
>>>>>> all processes, and commented out the call into/out of
>>>>>> independent data
>>>>>> mode. This produces MPI error messages...
>>>>>>
>>>>>> Cheers,
>>>>>>
>>>>>> Lex
>>>>>>
>>>>>> Wei-keng Liao wrote:
>>>>>>> Hi, Lex,
>>>>>>>
>>>>>>> Have you tried nfmpi_put_vara_real_all() in collective mode?
>>>>>>> Collective mode is the default mode.
>>>>>>>
>>>>>>> If yes, can you show me how you call it?
>>>>>>>
>>>>>>> Wei-keng
>>>>>>>
>>>>>>>
>>>>>>> On Feb 20, 2009, at 1:24 PM, Alex Hoffmann wrote:
>>>>>>>
>>>>>>>> Thank you for the clarifications.
>>>>>>>>
>>>>>>>> I included the mpi setup and netCDF output test code and the
>>>>>>>> corresponding Makefile as attachments. Maybe it can be of
>>>>>>>> use as a
>>>>>>>> simple example.
>>>>>>>>
>>>>>>>> I am currently running my code with the bugs in the output file
>>>>>>>> fixed,
>>>>>>>> and it seems to be working just fine.
>>>>>>>>
>>>>>>>> Thanks again,
>>>>>>>>
>>>>>>>> Lex
>>>>>>>>
>>>>>>>> Wei-keng Liao wrote:
>>>>>>>>> Hi, Lex, Please read below.
>>>>>>>>>
>>>>>>>>> On Feb 20, 2009, at 5:13 AM, Alex Hoffmann wrote:
>>>>>>>>>
>>>>>>>>>> Many thanks, Wei-keng, the issue seems to have been tracked
>>>>>>>>>> down...
>>>>>>>>>>
>>>>>>>>>> For reasons that are still not completely clear to me, using
>>>>>>>>>> trim() for
>>>>>>>>>> trimming the coordinate variable names and units (lines
>>>>>>>>>> 480/81,
>>>>>>>>>> 485/86)
>>>>>>>>>> seemed to leave a strange character (such as y diacritic or Y
>>>>>>>>>> accent)
>>>>>>>>>> trailing on e.g. the units (such as degrees_east). Though
>>>>>>>>>> this
>>>>>>>>>> did not
>>>>>>>>>> affect the run on its own, any subsequent attribute
>>>>>>>>>> definition
>>>>>>>>>> made the
>>>>>>>>>> code crash. I fixed this by removing the trim().
>>>>>>>>>
>>>>>>>>> If you called trim() on a string argument in a function, the
>>>>>>>>> tailing
>>>>>>>>> character '\0' will not be passed to the function. So if the
>>>>>>>>> string
>>>>>>>>> length
>>>>>>>>> argument has a larger value, the function will include
>>>>>>>>> undefined
>>>>>>>>> tailing
>>>>>>>>> characters as part of the string. Glad to know you solve the
>>>>>>>>> problem.
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>> If I may, I'd like to clarify one further issue with PnetCDF.
>>>>>>>>>> I initially implemented the code (also a small test
>>>>>>>>>> program) in
>>>>>>>>>> collective data mode. This did not work (the different
>>>>>>>>>> processes/Cartesian sub domains) did not write their data
>>>>>>>>>> into the
>>>>>>>>>> netCDF file. I then changed the data write statements into
>>>>>>>>>> independent
>>>>>>>>>> data mode, with which it worked immediately. Am I right to
>>>>>>>>>> use
>>>>>>>>>> independent data mode for this purpose?
>>>>>>>>>
>>>>>>>>> I am not sure what you meant here.
>>>>>>>>> Can you post your small test program to the list?
>>>>>>>>>
>>>>>>>>> Wei-keng
>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> Thank you very much for your support.
>>>>>>>>>>
>>>>>>>>>> Cheers
>>>>>>>>>>
>>>>>>>>>> Lex
>>>>>>>>>>
>>>>>>>>>> Wei-keng Liao wrote:
>>>>>>>>>>> Hi, Lex,
>>>>>>>>>>>
>>>>>>>>>>> If this error message "NC definations on multiprocesses
>>>>>>>>>>> conflict"
>>>>>>>>>>> is the only one you got before every crash, it is most
>>>>>>>>>>> likely the
>>>>>>>>>>> attribute values are not the same across all processes.
>>>>>>>>>>> Please check that first, so you can rule out that
>>>>>>>>>>> possibility.
>>>>>>>>>>>
>>>>>>>>>>> Also, can you check the attributes in the output file for
>>>>>>>>>>> the case
>>>>>>>>>>> that ran successfully? (using ncmpidump -c filename) This
>>>>>>>>>>> is just
>>>>>>>>>>> to see if the attributes are the ones you expected. Also,
>>>>>>>>>>> please
>>>>>>>>>>> make sure the number of attributes created is less than 4K
>>>>>>>>>>> global
>>>>>>>>>>> or per variable, as Rob indicated.
>>>>>>>>>>>
>>>>>>>>>>> Wei-keng
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> On Feb 19, 2009, at 5:53 PM, Alex Hoffmann wrote:
>>>>>>>>>>>
>>>>>>>>>>>> Dear Wei-keng,
>>>>>>>>>>>>
>>>>>>>>>>>> I understand that the attribute define arguments must be
>>>>>>>>>>>> the same
>>>>>>>>>>>> across
>>>>>>>>>>>> all processes. I am also almost certain that the
>>>>>>>>>>>> arguments you
>>>>>>>>>>>> cited
>>>>>>>>>>>> are identical across processes but I will perform a check.
>>>>>>>>>>>>
>>>>>>>>>>>> Be this as it may, I get the same arrow message if I
>>>>>>>>>>>> uncomment
>>>>>>>>>>>> one of
>>>>>>>>>>>> the text attributes specifying the unit or one of the real
>>>>>>>>>>>> attributes
>>>>>>>>>>>> specifying the fill value, for instance. There is no
>>>>>>>>>>>> possibility
>>>>>>>>>>>> that
>>>>>>>>>>>> these are different across processes as they are introduced
>>>>>>>>>>>> straight
>>>>>>>>>>>> into the function call. Any check I performed seemed to
>>>>>>>>>>>> indicate
>>>>>>>>>>>> that
>>>>>>>>>>>> the problem arises from the number of arguments I
>>>>>>>>>>>> defined, and
>>>>>>>>>>>> not
>>>>>>>>>>>> from
>>>>>>>>>>>> the arguments themselves, since they work all fine
>>>>>>>>>>>> individually.
>>>>>>>>>>>>
>>>>>>>>>>>> I was actually thinking to write the attributes ad-hoc
>>>>>>>>>>>> once the
>>>>>>>>>>>> dataset
>>>>>>>>>>>> is defined using regular netCDF functions. This however
>>>>>>>>>>>> seemed to
>>>>>>>>>>>> be a
>>>>>>>>>>>> little elegant solution defeating the purpose of using
>>>>>>>>>>>> PnetCDF.
>>>>>>>>>>>>
>>>>>>>>>>>> Thanks,
>>>>>>>>>>>>
>>>>>>>>>>>> Lex
>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>> Wei-keng Liao a écrit :
>>>>>>>>>>>>> Hi, Lex,
>>>>>>>>>>>>>
>>>>>>>>>>>>> PnetCDF requires the attribute define arguments be the
>>>>>>>>>>>>> same
>>>>>>>>>>>>> across
>>>>>>>>>>>>> all processes. This is why you got the error message for
>>>>>>>>>>>>> defining
>>>>>>>>>>>>> attribute "history" in line 400, if its value
>>>>>>>>>>>>> "conflicts" among
>>>>>>>>>>>>> processes.
>>>>>>>>>>>>> I think you have got this.
>>>>>>>>>>>>>
>>>>>>>>>>>>> I checked your codes. Starting from line 562, you have a
>>>>>>>>>>>>> few
>>>>>>>>>>>>> text
>>>>>>>>>>>>> attributes
>>>>>>>>>>>>> from des_tgas_picture(i), des_trac_picture(i),
>>>>>>>>>>>>> des_tpas_picture(i),
>>>>>>>>>>>>> and des.
>>>>>>>>>>>>> Similar for out_type == 2. Are you sure their values are
>>>>>>>>>>>>> all the
>>>>>>>>>>>>> same
>>>>>>>>>>>>> across
>>>>>>>>>>>>> all the processes? May be you can use MPI_Bcast to check
>>>>>>>>>>>>> if
>>>>>>>>>>>>> anyone has
>>>>>>>>>>>>> different value.
>>>>>>>>>>>>>
>>>>>>>>>>>>> Wei-keng
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> On Feb 19, 2009, at 11:14 AM, Alex Hoffmann wrote:
>>>>>>>>>>>>>
>>>>>>>>>>>>>> Dear Rob
>>>>>>>>>>>>>> Thanks a lot for your answer.
>>>>>>>>>>>>>> Isolating the PnetCDF output module from the code as a
>>>>>>>>>>>>>> stand-alone
>>>>>>>>>>>>>> would
>>>>>>>>>>>>>> require a substantial effort, and I do not encounter
>>>>>>>>>>>>>> this error
>>>>>>>>>>>>>> in my
>>>>>>>>>>>>>> small test programs. I attached the output module for
>>>>>>>>>>>>>> reference.
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> The attribute section that causes trouble is between
>>>>>>>>>>>>>> lines
>>>>>>>>>>>>>> 369-405
>>>>>>>>>>>>>> (global) and 468-613 (variables). I am aware that the
>>>>>>>>>>>>>> attribute
>>>>>>>>>>>>>> history
>>>>>>>>>>>>>> (line 400) may be problematic as every process may have a
>>>>>>>>>>>>>> slightly
>>>>>>>>>>>>>> different value (for the time element taken from the
>>>>>>>>>>>>>> fortran
>>>>>>>>>>>>>> function
>>>>>>>>>>>>>> date_and_time). Yet this is not the only issue.
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> If I prematurely exit at line 492, uncommenting any
>>>>>>>>>>>>>> further
>>>>>>>>>>>>>> attribute
>>>>>>>>>>>>>> definition between 468 and 489 causes a crash with the
>>>>>>>>>>>>>> multiprocesses
>>>>>>>>>>>>>> conflict.
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Note also that the attribute definitions between 558
>>>>>>>>>>>>>> and 613
>>>>>>>>>>>>>> only
>>>>>>>>>>>>>> produce no crash currently because of the print
>>>>>>>>>>>>>> statements on
>>>>>>>>>>>>>> lines 560
>>>>>>>>>>>>>> and 591, which made me think that maybe there is an
>>>>>>>>>>>>>> optimization
>>>>>>>>>>>>>> problem. What is strange is that if these are not
>>>>>>>>>>>>>> included,
>>>>>>>>>>>>>> the
>>>>>>>>>>>>>> crash
>>>>>>>>>>>>>> occurs only when I start writing the first file (movie
>>>>>>>>>>>>>> file,
>>>>>>>>>>>>>> out_type=2), AND if I leave the first block (if
>>>>>>>>>>>>>> out_type==1)
>>>>>>>>>>>>>> uncommented, although this block is not even entered in
>>>>>>>>>>>>>> that
>>>>>>>>>>>>>> particular
>>>>>>>>>>>>>> case, and vice versa. Hope this is more or less clear...
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Cheers,
>>>>>>>>>>>>>> Lex
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Robert Latham wrote:
>>>>>>>>>>>>>>> On Thu, Feb 19, 2009 at 10:09:37AM +0000, Alex
>>>>>>>>>>>>>>> Hoffmann wrote:
>>>>>>>>>>>>>>>> Hello,
>>>>>>>>>>>>>>>> I am getting the same error message as had already been
>>>>>>>>>>>>>>>> inquired
>>>>>>>>>>>>>>>> on in
>>>>>>>>>>>>>>>> this list, i.e.:
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> NC definations on multiprocesses conflict.
>>>>>>>>>>>>>>>> NOTE: Definitions across all process
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> and this seems to be related somehow randomly to the
>>>>>>>>>>>>>>>> total
>>>>>>>>>>>>>>>> number of
>>>>>>>>>>>>>>>> attributes (global or for each variable) that I
>>>>>>>>>>>>>>>> intend to
>>>>>>>>>>>>>>>> write
>>>>>>>>>>>>>>>> during
>>>>>>>>>>>>>>>> the the definition phase (and, though I am not sure
>>>>>>>>>>>>>>>> about
>>>>>>>>>>>>>>>> this, to
>>>>>>>>>>>>>>>> the
>>>>>>>>>>>>>>>> number of processes used to run my F90 parallel code).
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> More precisely, if I uncomment one put_att to replace
>>>>>>>>>>>>>>>> another,
>>>>>>>>>>>>>>>> the
>>>>>>>>>>>>>>>> program still runs, if I do not comment out another
>>>>>>>>>>>>>>>> one, the
>>>>>>>>>>>>>>>> program
>>>>>>>>>>>>>>>> crashes with the above error message.
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Is there any way of identifying where this problem
>>>>>>>>>>>>>>>> comes from
>>>>>>>>>>>>>>>> or is
>>>>>>>>>>>>>>>> there a maximum amount of attributes (though I don't
>>>>>>>>>>>>>>>> write an
>>>>>>>>>>>>>>>> excessive
>>>>>>>>>>>>>>>> amount of them)? There is no problem running it on one
>>>>>>>>>>>>>>>> process(or).
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Maybe you can show us with a small test program how
>>>>>>>>>>>>>>> you are
>>>>>>>>>>>>>>> encountering this error? The maximum number of
>>>>>>>>>>>>>>> attributes
>>>>>>>>>>>>>>> is 4k
>>>>>>>>>>>>>>> (the
>>>>>>>>>>>>>>> constant NC_MAX_ATTRS).
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> ==rob
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> --
>>>>>>>>>>>>>> ______________________________________________________________________
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Alex HOFFMANN PhD Candidate
>>>>>>>>>>>>>> Centre for Atmospheric Science
>>>>>>>>>>>>>> Department of Geography University of Cambridge
>>>>>>>>>>>>>> Downing Place, Cambridge CB2 3EN, UK
>>>>>>>>>>>>>> e-mail: ah519 at cam.ac.uk tel: +44 (0)1223 766581
>>>>>>>>>>>>>> www.geog.cam.ac.uk/people/hoffmann/
>>>>>>>>>>>>>> ______________________________________________________________________
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !-*- F90 -*-
>>>>>>>>>>>>>> module atham_pnetcdf
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !--------------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! author: Alex
>>>>>>>>>>>>>> HOFFMANN !
>>>>>>>>>>>>>> ! email:
>>>>>>>>>>>>>> ah519
>>>>>>>>>>>>>> @cam.ac.uk !
>>>>>>>>>>>>>> ! date: February
>>>>>>>>>>>>>> 2009 !
>>>>>>>>>>>>>> ! version:
>>>>>>>>>>>>>> v0.1
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> ! collection of routines for ATHAM PnetCDF
>>>>>>>>>>>>>> output !
>>>>>>>>>>>>>> !--------------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! module is to replace existing atham_netcdf.mod (GL)
>>>>>>>>>>>>>> to write
>>>>>>>>>>>>>> ATHAM
>>>>>>>>>>>>>> output as netCDF files
>>>>>>>>>>>>>> ! same as atham_netcdf.mod (AH) but adapted to run with
>>>>>>>>>>>>>> the
>>>>>>>>>>>>>> Parallel
>>>>>>>>>>>>>> NetCDF API,
>>>>>>>>>>>>>> ! developed by Northwestern University and Argonne
>>>>>>>>>>>>>> National
>>>>>>>>>>>>>> Laboratory,
>>>>>>>>>>>>>> ! which is built on top of MPI-IO (tested with MPICH
>>>>>>>>>>>>>> 2.1)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use precision, only : kint, kreal, real4, my_real
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : nx, ny, nz, ntx, nty, npx, npy
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : ntgas_picture, ntrac_picture,
>>>>>>>>>>>>>> ntpas_picture, &
>>>>>>>>>>>>>> ntgas_movie, ntrac_movie, ntpas_movie
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> implicit none
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> private
>>>>>>>>>>>>>> public :: netcdf_data
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! 'pnetcdf.inc' file includes cpp directives that were
>>>>>>>>>>>>>> not
>>>>>>>>>>>>>> understood
>>>>>>>>>>>>>> at compile-time
>>>>>>>>>>>>>> ! this is to be changed and hard-wired link to be
>>>>>>>>>>>>>> replaced
>>>>>>>>>>>>>> by the
>>>>>>>>>>>>>> include 'pnetcdf.inc' command
>>>>>>>>>>>>>> !include 'pnetcdf.inc'
>>>>>>>>>>>>>> #include "/usr/local/parallel-netcdf-1.0.1/include/
>>>>>>>>>>>>>> pnetcdf.inc"
>>>>>>>>>>>>>> include 'mpif.h'
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! ######### tried saving netCDF variable IDs here so
>>>>>>>>>>>>>> that they
>>>>>>>>>>>>>> must
>>>>>>>>>>>>>> not be retrieved by inquiry
>>>>>>>>>>>>>> ! ######### if this creates conflict, especially between
>>>>>>>>>>>>>> picture
>>>>>>>>>>>>>> and
>>>>>>>>>>>>>> movie file, this must be changed
>>>>>>>>>>>>>> integer(kint), save :: x_dimID,
>>>>>>>>>>>>>> y_dimID,
>>>>>>>>>>>>>> z_dimID, t_dimID
>>>>>>>>>>>>>> integer(kint), save :: x_varID,
>>>>>>>>>>>>>> y_varID,
>>>>>>>>>>>>>> z_varID, t_varID
>>>>>>>>>>>>>> integer(kint), save :: u_varID,
>>>>>>>>>>>>>> v_varID,
>>>>>>>>>>>>>> w_varID
>>>>>>>>>>>>>> integer(kint), save :: pnew_varID,
>>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>>> tempnew_varID, density_varID
>>>>>>>>>>>>>> integer(kint), save :: turbhor_varID,
>>>>>>>>>>>>>> turbver_varID, turblen_varID
>>>>>>>>>>>>>> integer(kint), dimension(99), save :: tgas_varID,
>>>>>>>>>>>>>> trac_varID,
>>>>>>>>>>>>>> tpas_varID, d_trac_varID
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> integer(MPI_OFFSET_KIND), save :: tcounter
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> real(real4), parameter :: undef =
>>>>>>>>>>>>>> -0.99e33_real4 !
>>>>>>>>>>>>>> make sure precision of _FillValue and output
>>>>>>>>>>>>>> (NF90_FLOAT) is
>>>>>>>>>>>>>> the
>>>>>>>>>>>>>> same
>>>>>>>>>>>>>> !!!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> contains
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =========================================================
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> subroutine netcdf_data(outfile,out_type,timetot)
>>>>>>>>>>>>>> ! IF: outfile = picture_file or movie_file (path/
>>>>>>>>>>>>>> filename for
>>>>>>>>>>>>>> output created in atham_module)
>>>>>>>>>>>>>> ! out_type = 1 (picture_file) or 2 (movie_file)
>>>>>>>>>>>>>> ! timetot = total simulation time since simulation
>>>>>>>>>>>>>> start
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : unew, vnew, wnew, pnew, tetnew,
>>>>>>>>>>>>>> tempnew, &
>>>>>>>>>>>>>> density, turbhor, turbver, turblen
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : tgasnew, tracnew, tpasnew,
>>>>>>>>>>>>>> tracdep
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : itgas_picture, itrac_picture,
>>>>>>>>>>>>>> itpas_picture, &
>>>>>>>>>>>>>> itgas_movie, itrac_movie, itpas_movie
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : mycoor, my_cart, myrank
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> character(*), intent(in ) :: outfile
>>>>>>>>>>>>>> integer(kint), intent(in ) :: out_type
>>>>>>>>>>>>>> real(kreal), intent(in ) :: timetot
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !------------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! local
>>>>>>>>>>>>>> variables
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> !------------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> integer(MPI_OFFSET_KIND) :: nx_mpi,ny_mpi,
>>>>>>>>>>>>>> nz_mpi,
>>>>>>>>>>>>>> y_cl_mpi
>>>>>>>>>>>>>> integer(MPI_OFFSET_KIND) :: istart,jstart
>>>>>>>>>>>>>> integer(kint) :: STATUS
>>>>>>>>>>>>>> character(len=50) :: OP ! debugging
>>>>>>>>>>>>>> variable
>>>>>>>>>>>>>> corresponding to the called netCDF OPeration
>>>>>>>>>>>>>> integer(kint) :: NCID
>>>>>>>>>>>>>> integer(kint) :: i
>>>>>>>>>>>>>> logical :: DynMovie_flag,
>>>>>>>>>>>>>> bounds_control
>>>>>>>>>>>>>> integer(kint) :: y_st, y_sp, y_cl
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! replace DynMovie precompiler flag by program flag for
>>>>>>>>>>>>>> more
>>>>>>>>>>>>>> convenient code structure
>>>>>>>>>>>>>> DynMovie_flag = .false.
>>>>>>>>>>>>>> #ifdef DynMovie
>>>>>>>>>>>>>> DynMovie_flag = .true.
>>>>>>>>>>>>>> #endif
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! re-assign dimension lengths as an MPI_Offset type
>>>>>>>>>>>>>> (necessary
>>>>>>>>>>>>>> for
>>>>>>>>>>>>>> the PnetCDF interface)
>>>>>>>>>>>>>> ! cut-off subdomain borders
>>>>>>>>>>>>>> nx_mpi = nx-2
>>>>>>>>>>>>>> ny_mpi = ny-2
>>>>>>>>>>>>>> nz_mpi = nz
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! assign coordinates of first element of each process
>>>>>>>>>>>>>> in the
>>>>>>>>>>>>>> netCDF file
>>>>>>>>>>>>>> istart = mycoor(1)*nx_mpi + 1_MPI_OFFSET_KIND
>>>>>>>>>>>>>> jstart = mycoor(2)*ny_mpi + 1_MPI_OFFSET_KIND
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> bounds_control = .true. ! debugging flag
>>>>>>>>>>>>>> tcounter = tcounter+1_MPI_OFFSET_KIND
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! === open netCDF dataset: enter data mode
>>>>>>>>>>>>>> ------------------------!
>>>>>>>>>>>>>> ! if dataset does not yet exist (STATUS = ), create
>>>>>>>>>>>>>> file and
>>>>>>>>>>>>>> enter data mode
>>>>>>>>>>>>>> ! else continue
>>>>>>>>>>>>>> OP = 'open existing netCDF file'
>>>>>>>>>>>>>> STATUS = nfmpi_open(my_cart, trim(outfile), NF_WRITE,
>>>>>>>>>>>>>> MPI_INFO_NULL, NCID)
>>>>>>>>>>>>>> if ( STATUS /= nf_noerr ) then
>>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>>> print *, ""
>>>>>>>>>>>>>> print *, "P-NetCDF operation: ",trim(OP)
>>>>>>>>>>>>>> print *, trim(nfmpi_strerror(STATUS))
>>>>>>>>>>>>>> print *, "atham_pnetcdf netcdf_data(): new netCDF file
>>>>>>>>>>>>>> created: ", outfile
>>>>>>>>>>>>>> print *, ""
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>> call netcdf_define(outfile,out_type,timetot,NCID)
>>>>>>>>>>>>>> tcounter = 1_MPI_OFFSET_KIND ! reset time counter for
>>>>>>>>>>>>>> output if
>>>>>>>>>>>>>> new file is created
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! --- write data values
>>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>>> ! unlimited record dimension time (only one process
>>>>>>>>>>>>>> writes the
>>>>>>>>>>>>>> whole dimension, i.e. accross all processes)
>>>>>>>>>>>>>> OP = 'write data values - unlimited dimension t'
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call handle_err(nfmpi_begin_indep_data(NCID))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID, t_varID,
>>>>>>>>>>>>>> (/tcounter/), (/1_MPI_OFFSET_KIND/),
>>>>>>>>>>>>>> real(timetot,real4)),OP)
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call handle_err(nfmpi_end_indep_data(NCID))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! --- write data values
>>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>>> ! discriminate between 2D and 3D case for output
>>>>>>>>>>>>>> generation
>>>>>>>>>>>>>> if (nty>4) then
>>>>>>>>>>>>>> y_st = 2 ! y index start
>>>>>>>>>>>>>> y_sp = ny-1 ! y index stop
>>>>>>>>>>>>>> y_cl = ny-2 ! y number of columns
>>>>>>>>>>>>>> else
>>>>>>>>>>>>>> y_st = ny-1 ! 3 in 2D case
>>>>>>>>>>>>>> y_sp = ny-1 ! 3 in 2D case
>>>>>>>>>>>>>> y_cl = 1
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>> y_cl_mpi = y_cl
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! dynamics data
>>>>>>>>>>>>>> OP = 'write data values - data'
>>>>>>>>>>>>>> if (out_type == 1 .OR. DynMovie_flag == .true.) then
>>>>>>>>>>>>>> call write_data(real(unew,real4),u_varID)
>>>>>>>>>>>>>> bounds_control = .false.
>>>>>>>>>>>>>> call write_data(real(vnew, real4),v_varID)
>>>>>>>>>>>>>> call write_data(real(wnew, real4),w_varID)
>>>>>>>>>>>>>> call write_data(real(pnew, real4),pnew_varID)
>>>>>>>>>>>>>> call write_data(real(tetnew, real4),tetnew_varID)
>>>>>>>>>>>>>> call write_data(real(tempnew,real4),tempnew_varID)
>>>>>>>>>>>>>> call write_data(real(density,real4),density_varID)
>>>>>>>>>>>>>> call write_data(real(turbhor,real4),turbhor_varID)
>>>>>>>>>>>>>> call write_data(real(turbver,real4),turbver_varID)
>>>>>>>>>>>>>> call write_data(real(turblen,real4),turblen_varID)
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! tracer data
>>>>>>>>>>>>>> if (out_type == 1) then
>>>>>>>>>>>>>> do i=1,ntgas_picture
>>>>>>>>>>>>>> call
>>>>>>>>>>>>>> write_data
>>>>>>>>>>>>>> (real
>>>>>>>>>>>>>> (tgasnew(:,:,:,itgas_picture(i)),real4),tgas_varID(i))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>> call
>>>>>>>>>>>>>> write_data
>>>>>>>>>>>>>> (real
>>>>>>>>>>>>>> (tracnew(:,:,:,itrac_picture(i)),real4),trac_varID(i))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntpas_picture
>>>>>>>>>>>>>> call
>>>>>>>>>>>>>> write_data
>>>>>>>>>>>>>> (real
>>>>>>>>>>>>>> (tpasnew(:,:,:,itpas_picture(i)),real4),tpas_varID(i))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>> call handle_err(nfmpi_begin_indep_data(NCID))
>>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID,
>>>>>>>>>>>>>> d_trac_varID(i),(/istart,jstart,tcounter/), &
>>>>>>>>>>>>>> (/nx_mpi,y_cl_mpi,1_MPI_OFFSET_KIND/),
>>>>>>>>>>>>>> real
>>>>>>>>>>>>>> (reshape
>>>>>>>>>>>>>> (tracdep(2:nx-1,y_st:y_sp,itrac_picture(i)),shape=(/
>>>>>>>>>>>>>> nx-2,y_cl,1/)),real4)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_end_indep_data(NCID))
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> else if (out_type == 2) then
>>>>>>>>>>>>>> do i=1,ntgas_movie
>>>>>>>>>>>>>> call
>>>>>>>>>>>>>> write_data
>>>>>>>>>>>>>> (real(tgasnew(:,:,:,itgas_movie(i)),real4),tgas_varID(i))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntrac_movie
>>>>>>>>>>>>>> call
>>>>>>>>>>>>>> write_data
>>>>>>>>>>>>>> (real(tracnew(:,:,:,itrac_movie(i)),real4),trac_varID(i))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntpas_movie
>>>>>>>>>>>>>> call
>>>>>>>>>>>>>> write_data
>>>>>>>>>>>>>> (real(tpasnew(:,:,:,itpas_movie(i)),real4),tpas_varID(i))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> endif
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! === close netCDF dataset: exit data mode & save
>>>>>>>>>>>>>> -----------------!
>>>>>>>>>>>>>> OP = 'close netCDF file'
>>>>>>>>>>>>>> call handle_err(nfmpi_close(NCID))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>>> print *, ''
>>>>>>>>>>>>>> print *, 'wrote netCDF output to file ', outfile
>>>>>>>>>>>>>> print *, ''
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> contains
>>>>>>>>>>>>>> !----------------------------------------------------------------
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> subroutine write_data(array_in,varID)
>>>>>>>>>>>>>> ! IF: array_in = data array (3D) to be written to netCDF
>>>>>>>>>>>>>> ! varID = data variable netCDF ID
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : iflgs ! topography mask:
>>>>>>>>>>>>>> topography
>>>>>>>>>>>>>> (iflgs=0) remains undef
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : myrank
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> real(real4),dimension(nx,ny,nz), intent (in ) ::
>>>>>>>>>>>>>> array_in
>>>>>>>>>>>>>> integer(kint), intent(in ) ::
>>>>>>>>>>>>>> varID
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! local
>>>>>>>>>>>>>> variables !
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> real(real4), dimension(nx,ny,nz) ::
>>>>>>>>>>>>>> array_masked
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! dimensions & bounds control
>>>>>>>>>>>>>> integer(MPI_OFFSET_KIND) :: dimlenx,
>>>>>>>>>>>>>> dimleny, dimlenz, dimlent
>>>>>>>>>>>>>> integer(kint), dimension(NF_max_var_dims) :: dimIDs
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! dimensions & bounds control
>>>>>>>>>>>>>> if (myrank == 0) then
>>>>>>>>>>>>>> if (bounds_control) then
>>>>>>>>>>>>>> print *, ''
>>>>>>>>>>>>>> print *, 'atham var in xyz
>>>>>>>>>>>>>> ',shape(array_in)
>>>>>>>>>>>>>> print *, 'atham var reshaped to xyzt
>>>>>>>>>>>>>> ',shape(reshape(array_in(2:nx-1,y_st:y_sp,1:nz),shape=(/
>>>>>>>>>>>>>> nx-2,y_cl,nz,1/)))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call handle_err(nfmpi_inq_vardimid(NCID, u_varID,
>>>>>>>>>>>>>> dimIDs)) ! taken unew as ref
>>>>>>>>>>>>>> call handle_err(nfmpi_inq_dimlen(NCID, dimIDs(1),
>>>>>>>>>>>>>> dimlenx))
>>>>>>>>>>>>>> call handle_err(nfmpi_inq_dimlen(NCID, dimIDs(2),
>>>>>>>>>>>>>> dimleny))
>>>>>>>>>>>>>> call handle_err(nfmpi_inq_dimlen(NCID, dimIDs(3),
>>>>>>>>>>>>>> dimlenz))
>>>>>>>>>>>>>> call handle_err(nfmpi_inq_dimlen(NCID, dimIDs(4),
>>>>>>>>>>>>>> dimlent))
>>>>>>>>>>>>>> print *, 'netCDF var dimensions in xyzt
>>>>>>>>>>>>>> ',dimlenx,dimleny,dimlenz,dimlent
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! mask topography
>>>>>>>>>>>>>> where (iflgs == 1)
>>>>>>>>>>>>>> array_masked = array_in
>>>>>>>>>>>>>> elsewhere
>>>>>>>>>>>>>> array_masked = undef
>>>>>>>>>>>>>> endwhere
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! mask control
>>>>>>>>>>>>>> !print *, ''
>>>>>>>>>>>>>> !print *, 'topography mask'
>>>>>>>>>>>>>> !print *, iflgs(1:5,3,3)
>>>>>>>>>>>>>> !print *, iflgs(1:5,3,2)
>>>>>>>>>>>>>> !print *, iflgs(1:5,3,1)
>>>>>>>>>>>>>> !print *, ''
>>>>>>>>>>>>>> !print *, 'masked variable array'
>>>>>>>>>>>>>> !print *, array_masked(:,3,3)
>>>>>>>>>>>>>> !print *, array_masked(:,3,2)
>>>>>>>>>>>>>> !print *, array_masked(:,3,1)
>>>>>>>>>>>>>> !print *, ''
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call handle_err(nfmpi_begin_indep_data(NCID))
>>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID,
>>>>>>>>>>>>>> varID,(/istart,jstart,1_MPI_OFFSET_KIND,tcounter/), &
>>>>>>>>>>>>>> (/nx_mpi,y_cl_mpi,nz_mpi,1_MPI_OFFSET_KIND/),
>>>>>>>>>>>>>> reshape(array_masked(2:nx-1,y_st:y_sp,1:nz),shape=(/
>>>>>>>>>>>>>> nx-2,y_cl,nz,1/))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_end_indep_data(NCID))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> end subroutine write_data
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> ----------------------------------------------------------------
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> end subroutine netcdf_data
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =========================================================
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> subroutine
>>>>>>>>>>>>>> netcdf_define(outfile,out_type,timetot,ncid_def)
>>>>>>>>>>>>>> ! IF: outfile = picture_file or movie_file (path/
>>>>>>>>>>>>>> filename for
>>>>>>>>>>>>>> output created in atham_module)
>>>>>>>>>>>>>> ! out_type = 1 (picture_file) or 2 (movie_file)
>>>>>>>>>>>>>> ! timetot = simulation time since simulation start
>>>>>>>>>>>>>> (from
>>>>>>>>>>>>>> atham/atham_out)
>>>>>>>>>>>>>> ! ncid_def = netCDF ID assigned to file during
>>>>>>>>>>>>>> creation
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : var_tgas_picture,
>>>>>>>>>>>>>> var_tgas_movie, &
>>>>>>>>>>>>>> des_tgas_picture,
>>>>>>>>>>>>>> des_tgas_movie, &
>>>>>>>>>>>>>> var_trac_picture,
>>>>>>>>>>>>>> var_trac_movie, &
>>>>>>>>>>>>>> des_trac_picture,
>>>>>>>>>>>>>> des_trac_movie, &
>>>>>>>>>>>>>> var_tpas_picture,
>>>>>>>>>>>>>> var_tpas_movie, &
>>>>>>>>>>>>>> des_tpas_picture, des_tpas_movie
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : iindex, spinup, nrep, periodt
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : volcano_setup,
>>>>>>>>>>>>>> coignimbrite_setup, &
>>>>>>>>>>>>>> nuclear_setup, convection_setup,
>>>>>>>>>>>>>> procsconfig_setup
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : x, y, z, xg, yg
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : grid2geo, sim_d_YYYY,
>>>>>>>>>>>>>> sim_d_MM,
>>>>>>>>>>>>>> sim_d_DD, &
>>>>>>>>>>>>>> sim_t_hh, sim_t_mm, sim_t_ss
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : my_cart, myrank
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> character(*), intent(in ) :: outfile
>>>>>>>>>>>>>> integer(kint), intent(in ) :: out_type
>>>>>>>>>>>>>> real(kreal), intent(in ) :: timetot
>>>>>>>>>>>>>> integer(kint), intent( out) :: ncid_def
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !------------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! local
>>>>>>>>>>>>>> variables
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> !------------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> integer(MPI_OFFSET_KIND) ::
>>>>>>>>>>>>>> ntx_mpi,nty_mpi,nz_mpi
>>>>>>>>>>>>>> character(len=50) :: OP ! debugging
>>>>>>>>>>>>>> variable
>>>>>>>>>>>>>> corresponding to the called netCDF OPeration
>>>>>>>>>>>>>> integer(kint) :: NCID
>>>>>>>>>>>>>> integer(kint) :: i
>>>>>>>>>>>>>> character(len=50) :: var, des
>>>>>>>>>>>>>> real(kreal) :: toffset
>>>>>>>>>>>>>> real(real4) :: coord_scale,
>>>>>>>>>>>>>> coord_x_offset, coord_y_offset, coord_t_offset
>>>>>>>>>>>>>> logical :: DynMovie_flag
>>>>>>>>>>>>>> character(len=3) :: coord_x_name,
>>>>>>>>>>>>>> coord_y_name, coord_z_name
>>>>>>>>>>>>>> character(len=4) :: coord_t_name
>>>>>>>>>>>>>> character(len=9) :: coord_x_lname,
>>>>>>>>>>>>>> coord_y_lname
>>>>>>>>>>>>>> character(len=13) :: coord_x_units,
>>>>>>>>>>>>>> coord_y_units
>>>>>>>>>>>>>> character(len=33) :: coord_t_units
>>>>>>>>>>>>>> character(len=8) :: date
>>>>>>>>>>>>>> character(len=10) :: time
>>>>>>>>>>>>>> character(len=19) :: simdatetime
>>>>>>>>>>>>>> character(len=10) :: simdate
>>>>>>>>>>>>>> character(len=20) :: login
>>>>>>>>>>>>>> character(len=63) :: hist_string
>>>>>>>>>>>>>> integer(kint) :: isetup
>>>>>>>>>>>>>> character(len=12) :: model_setup(6)
>>>>>>>>>>>>>> data model_setup /'volcano ', 'coignimbrite',
>>>>>>>>>>>>>> 'nuclear ',
>>>>>>>>>>>>>> 'convection ', 'procsconfig ','unknown '/
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! replace DynMovie precompiler flag by program flag for
>>>>>>>>>>>>>> more
>>>>>>>>>>>>>> convenient code structure
>>>>>>>>>>>>>> DynMovie_flag = .false.
>>>>>>>>>>>>>> #ifdef DynMovie
>>>>>>>>>>>>>> DynMovie_flag = .true.
>>>>>>>>>>>>>> #endif
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! re-assign dimension lengths as an MPI_Offset type
>>>>>>>>>>>>>> (necessary
>>>>>>>>>>>>>> for
>>>>>>>>>>>>>> the PnetCDF interface)
>>>>>>>>>>>>>> ! cut-off domain borders
>>>>>>>>>>>>>> ntx_mpi = ntx-2
>>>>>>>>>>>>>> nty_mpi = nty-2
>>>>>>>>>>>>>> nz_mpi = nz
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! convert date and time integers into string 'YYYYMMDD
>>>>>>>>>>>>>> hh:mm:ss'
>>>>>>>>>>>>>> write (simdatetime,10)
>>>>>>>>>>>>>> sim_d_YYYY,sim_d_MM,sim_d_DD,sim_t_hh,sim_t_mm,sim_t_ss
>>>>>>>>>>>>>> 10 format (I4,'-',I2,'-',I2,1X,I2,':',I2,':',I2)
>>>>>>>>>>>>>> simdate = simdatetime(1:10)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! switch between cartesian axes/time (default) and
>>>>>>>>>>>>>> geographic
>>>>>>>>>>>>>> axes/date (linear approx)
>>>>>>>>>>>>>> coord_scale = 1._real4
>>>>>>>>>>>>>> coord_x_offset = 0._real4
>>>>>>>>>>>>>> coord_x_name = 'x '
>>>>>>>>>>>>>> coord_x_lname = 'x_coor '
>>>>>>>>>>>>>> coord_x_units = 'm '
>>>>>>>>>>>>>> coord_y_offset = 0._real4
>>>>>>>>>>>>>> coord_y_name = 'y '
>>>>>>>>>>>>>> coord_y_lname = 'y_coor '
>>>>>>>>>>>>>> coord_y_units = 'm '
>>>>>>>>>>>>>> coord_z_name = 'z '
>>>>>>>>>>>>>> coord_t_offset = 0._real4
>>>>>>>>>>>>>> coord_t_name = 'time'
>>>>>>>>>>>>>> coord_t_units = 'sec '
>>>>>>>>>>>>>> if (grid2geo) then
>>>>>>>>>>>>>> call
>>>>>>>>>>>>>> gridconvert
>>>>>>>>>>>>>> (coord_scale
>>>>>>>>>>>>>> ,coord_x_offset,coord_x_name,coord_x_lname,coord_x_units,
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> &
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> coord_y_offset
>>>>>>>>>>>>>> ,coord_y_name,coord_y_lname,coord_y_units, &
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> coord_z_name,coord_t_offset,coord_t_name,coord_t_units)
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! === create netCDF dataset: enter define mode
>>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> OP = 'enter define mode'
>>>>>>>>>>>>>> call handle_err(nfmpi_create(my_cart, trim(outfile), 0,
>>>>>>>>>>>>>> MPI_INFO_NULL,NCID),OP)
>>>>>>>>>>>>>> ncid_def = NCID
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! === create netCDF dataset: assign metadata
>>>>>>>>>>>>>> ----------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! --- set global attributes
>>>>>>>>>>>>>> ----------------!
>>>>>>>>>>>>>> OP = 'set global attributes'
>>>>>>>>>>>>>> if (out_type == 1) then
>>>>>>>>>>>>>> toffset = timetot - spinup
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>>> 'title',
>>>>>>>>>>>>>> 13_MPI_OFFSET_KIND, 'atham_picture'),OP)
>>>>>>>>>>>>>> else
>>>>>>>>>>>>>> toffset = timetot - spinup - nrep*periodt
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>>> 'title',
>>>>>>>>>>>>>> 11_MPI_OFFSET_KIND, 'atham_movie'),OP)
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>>> 'content',
>>>>>>>>>>>>>> 223_MPI_OFFSET_KIND, &
>>>>>>>>>>>>>> 'This data corresponds to a simulation run of the
>>>>>>>>>>>>>> mesoscale-beta/gamma type non-hydrostatic Active Tracer
>>>>>>>>>>>>>> High-resolution Atmospheric Model ATHAM ' &
>>>>>>>>>>>>>> // '(Oberhuber et al., 1998, Herzog, 1998), for a setup
>>>>>>>>>>>>>> as
>>>>>>>>>>>>>> specified under history'),OP) ! ######### add further
>>>>>>>>>>>>>> description as
>>>>>>>>>>>>>> necessary
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>>> 'version',
>>>>>>>>>>>>>> 21_MPI_OFFSET_KIND, 'no current versioning'),OP) !
>>>>>>>>>>>>>> #########
>>>>>>>>>>>>>> change
>>>>>>>>>>>>>> once available
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call date_and_time(date, time)
>>>>>>>>>>>>>> call getlog(login)
>>>>>>>>>>>>>> if (volcano_setup) then
>>>>>>>>>>>>>> isetup=1
>>>>>>>>>>>>>> else if (coignimbrite_setup) then
>>>>>>>>>>>>>> isetup=2
>>>>>>>>>>>>>> else if (nuclear_setup) then
>>>>>>>>>>>>>> isetup=3
>>>>>>>>>>>>>> else if (convection_setup) then
>>>>>>>>>>>>>> isetup=4
>>>>>>>>>>>>>> else if (procsconfig_setup) then
>>>>>>>>>>>>>> isetup=5
>>>>>>>>>>>>>> else
>>>>>>>>>>>>>> isetup=6
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>> hist_string = trim(date//' '//time//' '//
>>>>>>>>>>>>>> trim(login)//' atham
>>>>>>>>>>>>>> '//model_setup(isetup))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>>> 'history',
>>>>>>>>>>>>>> 63_MPI_OFFSET_KIND, hist_string),OP) ! 1 line/
>>>>>>>>>>>>>> modification: -
>>>>>>>>>>>>>> cmd
>>>>>>>>>>>>>> arg
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>>> 'Conventions', 28_MPI_OFFSET_KIND, 'COARDS,
>>>>>>>>>>>>>> //cf-pcmdi.llnl.gov/'),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>>> 'sim_onset_date_time', 19_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> simdatetime),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, NF_GLOBAL,
>>>>>>>>>>>>>> 'sim_time_offset', NF_FLOAT, 1_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> real(toffset,real4)),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, NF_GLOBAL,
>>>>>>>>>>>>>> 'x_zoom_location', NF_FLOAT, 1_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> real(x(max(1,iindex)),real4)),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! --- define dimensions
>>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>>> OP = 'define dimensions'
>>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, trim(coord_x_name),
>>>>>>>>>>>>>> ntx_mpi, x_dimID),OP)
>>>>>>>>>>>>>> if (nty>4) then
>>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, trim(coord_y_name),
>>>>>>>>>>>>>> nty_mpi, y_dimID),OP)
>>>>>>>>>>>>>> else
>>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, trim(coord_y_name),
>>>>>>>>>>>>>> 1_MPI_OFFSET_KIND, y_dimID),OP)
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, trim(coord_z_name),
>>>>>>>>>>>>>> nz_mpi, z_dimID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, coord_t_name,
>>>>>>>>>>>>>> nfmpi_unlimited, t_dimID),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! --- define variables
>>>>>>>>>>>>>> ---------------------!
>>>>>>>>>>>>>> ! dimensions
>>>>>>>>>>>>>> OP = 'define variables - dimensions'
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, trim(coord_x_name),
>>>>>>>>>>>>>> NF_FLOAT, 1, x_dimID, x_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, trim(coord_y_name),
>>>>>>>>>>>>>> NF_FLOAT, 1, y_dimID, y_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, trim(coord_z_name),
>>>>>>>>>>>>>> NF_FLOAT, 1, z_dimID, z_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, coord_t_name,
>>>>>>>>>>>>>> NF_FLOAT, 1, t_dimID, t_varID),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! dynamics data
>>>>>>>>>>>>>> OP = 'define variables - data'
>>>>>>>>>>>>>> if (out_type == 1 .OR. DynMovie_flag == .true.) then
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'u',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), u_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'v',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), v_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'w',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), w_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'pnew',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), pnew_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'tetnew',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), tetnew_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'tempnew',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), tempnew_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'density',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), density_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'turbhor',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), turbhor_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'turbver',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), turbver_varID),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'turblen',
>>>>>>>>>>>>>> NF_FLOAT, 4,
>>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), turblen_varID),OP)
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! tracer data
>>>>>>>>>>>>>> if (out_type == 1) then
>>>>>>>>>>>>>> do i=1,ntgas_picture
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, var_tgas_picture(i),
>>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>>> tgas_varID(i)),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, var_trac_picture(i),
>>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>>> trac_varID(i)),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntpas_picture
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, var_tpas_picture(i),
>>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>>> tpas_varID(i)),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>> var='d_'//trim(var_trac_picture(i))
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, var,
>>>>>>>>>>>>>> NF_FLOAT, 3, (/x_dimID,y_dimID,t_dimID/),
>>>>>>>>>>>>>> d_trac_varID(i)),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> else if (out_type == 2) then
>>>>>>>>>>>>>> do i=1,ntgas_movie
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, var_tgas_movie(i),
>>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>>> tgas_varID(i)),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntrac_movie
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, var_trac_movie(i),
>>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>>> trac_varID(i)),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntpas_movie
>>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, var_tpas_movie(i),
>>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>>> tpas_varID(i)),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> endif
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! --- assign attributes
>>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>>> ! common attributes: long_name, units, valid_range,
>>>>>>>>>>>>>> scale_factor,
>>>>>>>>>>>>>> add_offset, _FillValue, history, conventions
>>>>>>>>>>>>>> ! dimensions
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> OP = 'assign attributes - dimensions'
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, x_varID,
>>>>>>>>>>>>>> 'long_name', 9_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> trim(coord_x_lname)),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, x_varID,
>>>>>>>>>>>>>> 'units', 13_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> trim(coord_x_units)),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, x_varID,
>>>>>>>>>>>>>> 'scale_factor', NF_FLOAT, 1_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> coord_scale),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, x_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> coord_x_offset),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, y_varID,
>>>>>>>>>>>>>> 'long_name', 9_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> trim(coord_y_lname)),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, y_varID,
>>>>>>>>>>>>>> 'units', 13_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> trim(coord_y_units)),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, y_varID,
>>>>>>>>>>>>>> 'scale_factor', NF_FLOAT, 1_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> coord_scale),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, y_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> coord_y_offset),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, z_varID,
>>>>>>>>>>>>>> 'long_name', 6_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'height'),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, z_varID,
>>>>>>>>>>>>>> 'units', 1_MPI_OFFSET_KIND, 'm'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, z_varID,
>>>>>>>>>>>>>> 'positive', 2_MPI_OFFSET_KIND, 'up'),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, t_varID,
>>>>>>>>>>>>>> 'long_name', 4_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> coord_t_name),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, t_varID,
>>>>>>>>>>>>>> 'units', 33_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> trim(coord_t_units)),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, t_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> coord_t_offset),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! OP = 'exit define mode NOW for testing
>>>>>>>>>>>>>> ****************************************************' !
>>>>>>>>>>>>>> ah_testing
>>>>>>>>>>>>>> ********************************
>>>>>>>>>>>>>> ! call
>>>>>>>>>>>>>> handle_err(nfmpi_enddef(NCID),OP)
>>>>>>>>>>>>>> ! ah_testing ********************************
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! dynamics data
>>>>>>>>>>>>>> OP = 'assign attributes - data'
>>>>>>>>>>>>>> if (out_type == 1 .OR. DynMovie_flag == .true.) then
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, u_varID,
>>>>>>>>>>>>>> 'long_name', 6_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'x_wind'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, u_varID,
>>>>>>>>>>>>>> 'units', 7_MPI_OFFSET_KIND, 'm
>>>>>>>>>>>>>> sec-1'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, u_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, u_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, u_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, v_varID,
>>>>>>>>>>>>>> 'long_name', 6_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'y_wind'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, v_varID,
>>>>>>>>>>>>>> 'units', 7_MPI_OFFSET_KIND, 'm
>>>>>>>>>>>>>> sec-1'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, v_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, v_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, v_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, w_varID,
>>>>>>>>>>>>>> 'long_name', 19_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'upward_air_velocity'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID, w_varID,
>>>>>>>>>>>>>> 'units', 7_MPI_OFFSET_KIND, 'm
>>>>>>>>>>>>>> sec-1'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, w_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, w_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID, w_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>>> 'long_name', 20_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'air_pressure_anomaly'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>>> 'units', 2_MPI_OFFSET_KIND, 'Pa'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>>> 'long_name', 25_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'air_potential_temperature'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>>> 'units', 1_MPI_OFFSET_KIND, 'K'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>>> 'long_name', 15_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'air_temperature'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>>> 'units', 7_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'Celsius'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, -273.15),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>>> 'long_name', 11_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'air_density'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>>> 'units', 6_MPI_OFFSET_KIND, 'kg
>>>>>>>>>>>>>> m-3'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>>> 'long_name', 15_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'hor_turb_energy'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>>> 'units', 8_MPI_OFFSET_KIND, 'm2
>>>>>>>>>>>>>> sec-2'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>>> 'long_name', 15_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'ver_turb_energy'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>>> 'units', 8_MPI_OFFSET_KIND, 'm2
>>>>>>>>>>>>>> sec-2'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>>> 'long_name', 17_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> 'turb_length_scale'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_text(NCID,
>>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>>> 'units', 1_MPI_OFFSET_KIND, 'm'),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> ! call handle_err(nfmpi_put_att_real(NCID,
>>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! tracer data
>>>>>>>>>>>>>> if (out_type == 1) then
>>>>>>>>>>>>>> print *, 'NetCDF output type (1-picture,2-movie)
>>>>>>>>>>>>>> ',out_type !
>>>>>>>>>>>>>> ##### DO NOT REMOVE!!! optimization problem, trace
>>>>>>>>>>>>>> source
>>>>>>>>>>>>>> do i=1,ntgas_picture
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, tgas_varID(i),
>>>>>>>>>>>>>> 'long_name', 25_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> des_tgas_picture(i)),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, tgas_varID(i),
>>>>>>>>>>>>>> 'units', 9_MPI_OFFSET_KIND, 'gram
>>>>>>>>>>>>>> kg-1'),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, trac_varID(i),
>>>>>>>>>>>>>> 'long_name', 25_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> des_trac_picture(i)),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, trac_varID(i),
>>>>>>>>>>>>>> 'units', 9_MPI_OFFSET_KIND, 'gram
>>>>>>>>>>>>>> kg-1'),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntpas_picture
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, tpas_varID(i),
>>>>>>>>>>>>>> 'long_name', 25_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> des_tpas_picture(i)),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, tpas_varID(i),
>>>>>>>>>>>>>> 'units', 1_MPI_OFFSET_KIND, '%'),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 100.),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>> des='deposited_'//trim(des_trac_picture(i))
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, d_trac_varID(i),
>>>>>>>>>>>>>> 'long_name', 35_MPI_OFFSET_KIND, des),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, d_trac_varID(i),
>>>>>>>>>>>>>> 'units', 9_MPI_OFFSET_KIND, 'gram
>>>>>>>>>>>>>> m-2'),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, d_trac_varID(i),
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, d_trac_varID(i),
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, d_trac_varID(i),
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> else if (out_type == 2) then
>>>>>>>>>>>>>> print *, 'NetCDF output type (1-picture,2-movie)
>>>>>>>>>>>>>> ',out_type !
>>>>>>>>>>>>>> ##### DO NOT REMOVE!!! optimization problem, trace
>>>>>>>>>>>>>> source
>>>>>>>>>>>>>> do i=1,ntgas_movie
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, tgas_varID(i),
>>>>>>>>>>>>>> 'long_name', 25_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> des_tgas_movie(i)),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, tgas_varID(i),
>>>>>>>>>>>>>> 'units', 9_MPI_OFFSET_KIND, 'gram
>>>>>>>>>>>>>> kg-1'),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntrac_movie
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, trac_varID(i),
>>>>>>>>>>>>>> 'long_name', 25_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> des_trac_movie(i)),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, trac_varID(i),
>>>>>>>>>>>>>> 'units', 9_MPI_OFFSET_KIND, 'gram
>>>>>>>>>>>>>> kg-1'),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> do i=1,ntpas_movie
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, tpas_varID(i),
>>>>>>>>>>>>>> 'long_name', 25_MPI_OFFSET_KIND,
>>>>>>>>>>>>>> des_tpas_movie(i)),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, tpas_varID(i),
>>>>>>>>>>>>>> 'units', 1_MPI_OFFSET_KIND, '%'),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>>> 'scale_factor',NF_FLOAT, 1_MPI_OFFSET_KIND, 100.),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>>> 'add_offset', NF_FLOAT, 1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>>> '_FillValue', NF_FLOAT, 1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>> enddo
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! === create netCDF dataset: exit define mode, enter
>>>>>>>>>>>>>> data mode
>>>>>>>>>>>>>> ----!
>>>>>>>>>>>>>> OP = 'exit define mode'
>>>>>>>>>>>>>> call handle_err(nfmpi_enddef(NCID),OP)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! --- write data values
>>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>>> ! dimensions (only one process writes the whole
>>>>>>>>>>>>>> dimension, i.e.
>>>>>>>>>>>>>> accross all processes)
>>>>>>>>>>>>>> OP = 'write data values - dimensions'
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call handle_err(nfmpi_begin_indep_data(NCID))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID, x_varID,
>>>>>>>>>>>>>> (/1_MPI_OFFSET_KIND/), (/
>>>>>>>>>>>>>> ntx_mpi/),real(xg(2:ntx-1),real4)),OP)
>>>>>>>>>>>>>> if (nty>4) then
>>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID, y_varID,
>>>>>>>>>>>>>> (/1_MPI_OFFSET_KIND/), (/nty_mpi/),
>>>>>>>>>>>>>> real(yg(2:nty-1),real4)),OP)
>>>>>>>>>>>>>> else
>>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID, y_varID,
>>>>>>>>>>>>>> (/1_MPI_OFFSET_KIND/), (/1_MPI_OFFSET_KIND/),
>>>>>>>>>>>>>> real(yg(nty-1),real4)),OP)
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID, z_varID,
>>>>>>>>>>>>>> (/1_MPI_OFFSET_KIND/), (/nz_mpi/), real(z,real4)),OP)
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> call handle_err(nfmpi_end_indep_data(NCID))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> contains
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> ----------------------------------------------------------------
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> subroutine
>>>>>>>>>>>>>> gridconvert
>>>>>>>>>>>>>> (cd_scale,cd_x_offset,cd_x_name,cd_x_lname,cd_x_units,
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> &
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> cd_y_offset,cd_y_name,cd_y_lname,cd_y_units, &
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> cd_z_name,cd_t_offset,cd_t_name,cd_t_units)
>>>>>>>>>>>>>> ! IF: cd_* = coordinate parameters
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! linear conversion from cartesian grid into lon/lat
>>>>>>>>>>>>>> system !
>>>>>>>>>>>>>> ! only valid as first approxiamtion over small
>>>>>>>>>>>>>> domains !
>>>>>>>>>>>>>> ! to be replaced once a dyn core in a lat/lon system
>>>>>>>>>>>>>> available !
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : deglat_start,
>>>>>>>>>>>>>> deglon_start
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use phys_constants, only : pi
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> real(real4), intent(inout) :: cd_scale,
>>>>>>>>>>>>>> cd_x_offset,
>>>>>>>>>>>>>> cd_y_offset, cd_t_offset
>>>>>>>>>>>>>> character(len=3), intent(inout) :: cd_x_name,
>>>>>>>>>>>>>> cd_y_name,
>>>>>>>>>>>>>> cd_z_name
>>>>>>>>>>>>>> character(len=4), intent(inout) :: cd_t_name
>>>>>>>>>>>>>> character(len=9), intent(inout) :: cd_x_lname,
>>>>>>>>>>>>>> cd_y_lname
>>>>>>>>>>>>>> character(len=13), intent(inout) :: cd_x_units,
>>>>>>>>>>>>>> cd_y_units
>>>>>>>>>>>>>> character(len=33), intent(inout) :: cd_t_units
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! local
>>>>>>>>>>>>>> variables !
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> real(kreal) :: r
>>>>>>>>>>>>>> real(kreal), parameter :: a =
>>>>>>>>>>>>>> 6378137._kreal !
>>>>>>>>>>>>>> (WGS 84)
>>>>>>>>>>>>>> real(kreal), parameter :: b =
>>>>>>>>>>>>>> 6356752.3142_kreal !
>>>>>>>>>>>>>> (WGS 84)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ! ellipsoid radius at ATHAM domain start
>>>>>>>>>>>>>> r = a*b / sqrt((b*b*(cos(deglat_start*pi/180._kreal))**2
>>>>>>>>>>>>>> )+(a*a*(sin(deglat_start*pi/180._kreal))**2 ))
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> cd_scale = real((1/r) * (180._kreal/pi),real4)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> cd_x_offset = deglon_start
>>>>>>>>>>>>>> cd_x_name = 'lon'
>>>>>>>>>>>>>> cd_x_lname = 'longitude'
>>>>>>>>>>>>>> cd_x_units = 'degrees_east '
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> cd_y_offset = deglat_start
>>>>>>>>>>>>>> cd_y_name = 'lat'
>>>>>>>>>>>>>> cd_y_lname = 'latitude '
>>>>>>>>>>>>>> cd_y_units = 'degrees_north'
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> cd_z_name = 'lev'
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> coord_t_offset = real((sim_t_hh*3600 + sim_t_mm*60
>>>>>>>>>>>>>> +sim_t_ss),real4)
>>>>>>>>>>>>>> coord_t_name = 'date'
>>>>>>>>>>>>>> coord_t_units = 'seconds since '//simdate//' 00:00:00
>>>>>>>>>>>>>> [local]'
>>>>>>>>>>>>>> ! COARDS conventions require time zone from UTC to be
>>>>>>>>>>>>>> specified
>>>>>>>>>>>>>> ! not necessary for local area version of ATHAM
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> end subroutine gridconvert
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> ----------------------------------------------------------------
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> end subroutine netcdf_define
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =========================================================
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =========================================================
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> subroutine handle_err(status,op)
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> use atham_module, only : myrank
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> integer, intent (in ) :: status
>>>>>>>>>>>>>> character(len=*), optional, intent (in ) :: op
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> if ( status /= nf_noerr ) then
>>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>>> print *, ""
>>>>>>>>>>>>>> print *, "P-NetCDF operation: ",trim(op)
>>>>>>>>>>>>>> print *, trim(nfmpi_strerror(status))
>>>>>>>>>>>>>> print *, "atham_pnetcdf error: program execution stopped"
>>>>>>>>>>>>>> print *, ""
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> stop
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>> end if
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> return
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> end subroutine handle_err
>>>>>>>>>>>>>> !
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =
>>>>>>>>>>>>>> =========================================================
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> end module atham_pnetcdf
>>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> --
>>>>>>>>>> ______________________________________________________________________
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> Alex HOFFMANN PhD Candidate
>>>>>>>>>> Centre for Atmospheric Science
>>>>>>>>>> Department of Geography University of Cambridge
>>>>>>>>>> Downing Place, Cambridge CB2 3EN, UK
>>>>>>>>>> e-mail: ah519 at cam.ac.uk tel: +44 (0)1223 766581
>>>>>>>>>> www.geog.cam.ac.uk/people/hoffmann/
>>>>>>>>>> ______________________________________________________________________
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>> --
>>>>>>>> ______________________________________________________________________
>>>>>>>>
>>>>>>>>
>>>>>>>> Alex HOFFMANN PhD Candidate
>>>>>>>> Centre for Atmospheric Science
>>>>>>>> Department of Geography University of Cambridge
>>>>>>>> Downing Place, Cambridge CB2 3EN, UK
>>>>>>>> e-mail: ah519 at cam.ac.uk tel: +44 (0)1223 766581
>>>>>>>> www.geog.cam.ac.uk/people/hoffmann/
>>>>>>>> ______________________________________________________________________
>>>>>>>>
>>>>>>>> !-*- F90 -*-
>>>>>>>> program mpi_setup_netCDFout
>>>>>>>>
>>>>>>>> !--------------------------------------------------------------------!
>>>>>>>>
>>>>>>>> ! author: Alex
>>>>>>>> HOFFMANN !
>>>>>>>> ! email:
>>>>>>>> ah519 at cam.ac.uk !
>>>>>>>> ! date: February
>>>>>>>> 2009 !
>>>>>>>> ! version:
>>>>>>>> development !
>>>>>>>> !
>>>>>>>> !
>>>>>>>> ! simple program to test mpi runs & netCDF output (serial &
>>>>>>>> parallel)!
>>>>>>>> !--------------------------------------------------------------------!
>>>>>>>>
>>>>>>>>
>>>>>>>> ! test code designed to be used with the variables from the
>>>>>>>> mesoscale
>>>>>>>> ! atmospheric model ATHAM
>>>>>>>> ! re-structured as stand-alone
>>>>>>>> ! reads and writes as netCDF an index matrix for a 4x4
>>>>>>>> matrix in
>>>>>>>> real4
>>>>>>>>
>>>>>>>> ! use precision
>>>>>>>> ! use atham_module, only : nx,ny,ntx,nty,npx,npy
>>>>>>>>
>>>>>>>> implicit none
>>>>>>>>
>>>>>>>> !private
>>>>>>>> !public
>>>>>>>>
>>>>>>>> ! --- mpi interface for Fortran 77 only, use header
>>>>>>>> files-----------
>>>>>>>> include 'mpif.h'
>>>>>>>>
>>>>>>>>
>>>>>>>> ! --- precision
>>>>>>>> ----------------------------------------------------
>>>>>>>> integer, parameter :: real4 =
>>>>>>>> selected_real_kind(6), &
>>>>>>>> real8 =
>>>>>>>> selected_real_kind(12), &
>>>>>>>> kreal = real8
>>>>>>>> integer, parameter :: int4 =
>>>>>>>> selected_int_kind(7), &
>>>>>>>> int8 =
>>>>>>>> selected_int_kind(13), &
>>>>>>>> kint = int4
>>>>>>>> ! --- atham_module
>>>>>>>> -------------------------------------------------
>>>>>>>> integer(kint) ::
>>>>>>>> nx,ny,ntx,nty,npx,npy
>>>>>>>> integer(kint) :: my_ranks,
>>>>>>>> nprocs, &
>>>>>>>> my_cart, myrank,
>>>>>>>> mycoor(2)
>>>>>>>> logical :: mpi_initialized,
>>>>>>>> mpi
>>>>>>>> ! --- other
>>>>>>>> --------------------------------------------------------
>>>>>>>> integer(kint) :: ierror
>>>>>>>> real(real4), allocatable, dimension(:,:) :: array_g, array
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> mpi_initialized=.false.
>>>>>>>> ! ---for initial
>>>>>>>> testing--------------------------------------------
>>>>>>>> nx = 4_kint ! user-specified number of grid points along x
>>>>>>>> ny = 4_kint
>>>>>>>> npx = 2_kint ! number of processes along x
>>>>>>>> npy = 2_kint
>>>>>>>>
>>>>>>>> ntx = nx ! re-adjusted total number of grid points along x
>>>>>>>> nty = ny
>>>>>>>> nx = nx/npx ! number of grid points along x per process
>>>>>>>> ny = ny/npy
>>>>>>>>
>>>>>>>>
>>>>>>>> ! --- MAIN
>>>>>>>> ! --- mpi setup calls
>>>>>>>> ----------------------------------------------
>>>>>>>> call initialize_mpi
>>>>>>>> call setup_mpi
>>>>>>>>
>>>>>>>> ! --- netCDF out calls
>>>>>>>> ---------------------------------------------
>>>>>>>> call nc_out_s
>>>>>>>> call nc_out_p
>>>>>>>>
>>>>>>>> ! --- mpi finalize call
>>>>>>>> --------------------------------------------
>>>>>>>> call mpi_Finalize(ierror)
>>>>>>>> ! --- END MAIN
>>>>>>>>
>>>>>>>> contains
>>>>>>>> !
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> ===============================================================
>>>>>>>>
>>>>>>>> subroutine initialize_mpi()
>>>>>>>>
>>>>>>>> ! use atham_module, only : mpi, mpi_initialized
>>>>>>>> ! use atham_module, only : myrank, nprocs
>>>>>>>>
>>>>>>>> integer(kint) :: ierror, errorclass
>>>>>>>> character(len=100) :: message
>>>>>>>>
>>>>>>>> ! assign variables for serial version
>>>>>>>> ierror=0
>>>>>>>> myrank=0
>>>>>>>> nprocs=1
>>>>>>>> mpi=.false.
>>>>>>>> ! initialize MPI
>>>>>>>> mpi=.true.
>>>>>>>> if (.not. mpi_initialized) call mpi_init(ierror)
>>>>>>>> call MPI_Comm_rank(mpi_comm_world,myrank,ierror)
>>>>>>>> call MPI_Comm_size(mpi_comm_world,nprocs,ierror)
>>>>>>>> mpi=.true.
>>>>>>>> mpi_initialized=.true.
>>>>>>>>
>>>>>>>> if (myrank==0) then
>>>>>>>> print *, 'Myrank 0: ',myrank
>>>>>>>> print *, 'MPI initialized!'
>>>>>>>> print *, 'Number of processes [nprocs]: ',nprocs
>>>>>>>> endif
>>>>>>>>
>>>>>>>> if (ierror>0) then
>>>>>>>> call MPI_Error_string(ierror,message,100)
>>>>>>>> call MPI_Error_class(ierror,errorclass)
>>>>>>>> call exec_stop(message,errorclass)
>>>>>>>> endif
>>>>>>>>
>>>>>>>> end subroutine initialize_mpi
>>>>>>>> !
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> ===============================================================
>>>>>>>>
>>>>>>>> subroutine setup_mpi()
>>>>>>>>
>>>>>>>> ! use atham_module, only : mpi, mpi_initialized
>>>>>>>> ! use atham_module, only : myrank,nprocs,my_cart,mycoor
>>>>>>>>
>>>>>>>> integer(kint) :: ierror, errorclass
>>>>>>>> character(len=100) :: message
>>>>>>>> integer(kint) :: idims(2),shape_array(2)
>>>>>>>> logical :: reorder
>>>>>>>> integer(kint) :: i,j,istart,jstart
>>>>>>>>
>>>>>>>> if (nprocs/=npx*npy) then
>>>>>>>> print *, 'nprocs /= npx*npy - execution aborted'
>>>>>>>> call mpi_abort(mpi_comm_world,ierror)
>>>>>>>> stop
>>>>>>>> endif
>>>>>>>> ! create virtual topology
>>>>>>>> idims(1:2) = (/npx,npy/)
>>>>>>>> reorder = .true.
>>>>>>>>
>>>>>>>> call
>>>>>>>> MPI_Cart_create(mpi_comm_world,2,idims,
>>>>>>>> (/.false.,.false./),reorder,my_cart,ierror)
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> call MPI_Comm_rank(my_cart,myrank,ierror) ! ranks on new
>>>>>>>> communicator
>>>>>>>> call MPI_Cart_coords(my_cart,myrank,2,mycoor,ierror)
>>>>>>>>
>>>>>>>> if(myrank==0) then
>>>>>>>> print *, 'Coords of process ',myrank,' :',mycoor
>>>>>>>> endif
>>>>>>>>
>>>>>>>> if (ierror>0) then
>>>>>>>> call MPI_Error_string(ierror,message,100)
>>>>>>>> call MPI_Error_class(ierror,errorclass)
>>>>>>>> call exec_stop(message,errorclass)
>>>>>>>> endif
>>>>>>>>
>>>>>>>> ! allocate and assign global grid on each process
>>>>>>>> allocate(array_g(ntx,nty))
>>>>>>>>
>>>>>>>> array_g(1,:) = (/ 11.0,12.0,13.0,14.0 /)
>>>>>>>> array_g(2,:) = (/ 21.0,22.0,23.0,24.0 /)
>>>>>>>> array_g(3,:) = (/ 31.0,32.0,33.0,34.0 /)
>>>>>>>> array_g(4,:) = (/ 41.0,42.0,43.0,44.0 /)
>>>>>>>> if (myrank==0) then
>>>>>>>> print *, 'whole 4*4 array on each process:'
>>>>>>>> print *, array_g(1,:)
>>>>>>>> print *, array_g(2,:)
>>>>>>>> print *, array_g(3,:)
>>>>>>>> print *, array_g(4,:)
>>>>>>>> endif
>>>>>>>>
>>>>>>>> ! allocate and assign subgrid portions on various processes
>>>>>>>> allocate(array(nx,ny))
>>>>>>>> istart = mycoor(1)*nx
>>>>>>>> jstart = mycoor(2)*ny
>>>>>>>>
>>>>>>>> do i=1,nx
>>>>>>>> do j=1,ny
>>>>>>>> array(i,j) = array_g(istart+i,jstart+j)
>>>>>>>> enddo
>>>>>>>> enddo
>>>>>>>>
>>>>>>>> ! output results from one process
>>>>>>>> if (myrank==0) then
>>>>>>>> print *, ''
>>>>>>>> print *, 'subgrid array on process with coords ',mycoor,':'
>>>>>>>> shape_array = shape(array)
>>>>>>>> print *, 'shape: ',shape_array
>>>>>>>> do i=1,shape_array(1)
>>>>>>>> print *, array(i,:)
>>>>>>>> enddo
>>>>>>>> print *, ''
>>>>>>>> endif
>>>>>>>>
>>>>>>>> end subroutine setup_mpi
>>>>>>>> !
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> ===============================================================
>>>>>>>>
>>>>>>>> subroutine nc_out_s()
>>>>>>>>
>>>>>>>> use netcdf
>>>>>>>>
>>>>>>>> ! use atham_module, only : makename
>>>>>>>> ! use atham_module, only : mycoor
>>>>>>>> ! use atham_module, only : nx,ny
>>>>>>>>
>>>>>>>> integer(kint) :: NCID
>>>>>>>> character(len=10) :: name
>>>>>>>> integer(kint) ::
>>>>>>>> x_dimID,y_dimID,x_varID,y_varID,array_varID
>>>>>>>>
>>>>>>>>
>>>>>>>> call makename('out',mycoor,name)
>>>>>>>> name = trim(name)//'.nc'
>>>>>>>>
>>>>>>>> call handle_err_nc_s(nf90_create(path = './'//trim(name),
>>>>>>>> cmode =
>>>>>>>> 0, ncid = NCID))
>>>>>>>> call handle_err_nc_s(nf90_def_dim(ncid = NCID, name =
>>>>>>>> 'x', len
>>>>>>>> = nx, dimid = x_dimID))
>>>>>>>> call handle_err_nc_s(nf90_def_dim(ncid = NCID, name =
>>>>>>>> 'y', len
>>>>>>>> = ny, dimid = y_dimID))
>>>>>>>> call handle_err_nc_s(nf90_def_var(ncid = NCID, name = 'x',
>>>>>>>> xtype = NF90_FLOAT, dimids = x_dimID, varid = x_varID))
>>>>>>>> call handle_err_nc_s(nf90_def_var(ncid = NCID, name = 'y',
>>>>>>>> xtype = NF90_FLOAT, dimids = y_dimID, varid = y_varID))
>>>>>>>> call handle_err_nc_s(nf90_def_var(ncid = NCID, name = 'array',
>>>>>>>> xtype = NF90_FLOAT, dimids = (/x_dimID,y_dimID/), varid =
>>>>>>>> array_varID))
>>>>>>>> call handle_err_nc_s(nf90_enddef(ncid = NCID))
>>>>>>>> call handle_err_nc_s(nf90_put_var(ncid = NCID, varid =
>>>>>>>> array_varID,
>>>>>>>> values = array))
>>>>>>>> call handle_err_nc_s(nf90_close(ncid = NCID))
>>>>>>>>
>>>>>>>> end subroutine nc_out_s
>>>>>>>> !
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> ===============================================================
>>>>>>>>
>>>>>>>> subroutine nc_out_p
>>>>>>>>
>>>>>>>> ! use atham_module, only : mycoor,my_cart,myrank
>>>>>>>> ! use atham_module, only : nx,ny
>>>>>>>>
>>>>>>>> ! precompiler problem: include file had to be hardwired into
>>>>>>>> code
>>>>>>>> ! include 'pnetcdf.inc'
>>>>>>>> #include "/usr/local/parallel-netcdf-1.0.1/include/pnetcdf.inc"
>>>>>>>>
>>>>>>>>
>>>>>>>> integer(kint) :: NCIDP
>>>>>>>> character(len=10) :: name
>>>>>>>> integer(kint) ::
>>>>>>>> x_dimIDP,y_dimIDP,x_varIDP,y_varIDP,array_varIDP
>>>>>>>> integer(kint) :: dimIDPs(2)
>>>>>>>> integer(MPI_OFFSET_KIND) :: dimlenx,dimleny
>>>>>>>> integer(MPI_OFFSET_KIND) :: istart,jstart
>>>>>>>> integer(MPI_OFFSET_KIND) :: ntx_mpi,nty_mpi,nx_mpi,ny_mpi
>>>>>>>>
>>>>>>>> ntx_mpi = ntx
>>>>>>>> nty_mpi = nty
>>>>>>>> nx_mpi = nx
>>>>>>>> ny_mpi = ny
>>>>>>>>
>>>>>>>> name = 'out_glb.nc'
>>>>>>>> istart = mycoor(1)*nx_mpi + 1
>>>>>>>> jstart = mycoor(2)*ny_mpi + 1
>>>>>>>>
>>>>>>>>
>>>>>>>> call
>>>>>>>> handle_err_nc_p(nfmpi_create(my_cart,'./'//trim(name),
>>>>>>>> 0,MPI_INFO_NULL,NCIDP),0)
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> call
>>>>>>>> handle_err_nc_p(nfmpi_def_dim(NCIDP,'x',ntx_mpi,x_dimIDP),1)
>>>>>>>>
>>>>>>>> call
>>>>>>>> handle_err_nc_p(nfmpi_def_dim(NCIDP,'y',nty_mpi,y_dimIDP),2)
>>>>>>>>
>>>>>>>> call
>>>>>>>> handle_err_nc_p(nfmpi_def_var(NCIDP,'x',NF_FLOAT,
>>>>>>>> 1,x_dimIDP,x_varIDP),3)
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> call
>>>>>>>> handle_err_nc_p(nfmpi_def_var(NCIDP,'y',NF_FLOAT,
>>>>>>>> 1,y_dimIDP,y_varIDP),4)
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> call
>>>>>>>> handle_err_nc_p(nfmpi_def_var(NCIDP,'array',NF_FLOAT,2,(/
>>>>>>>> x_dimIDP,y_dimIDP/),array_varIDP),5)
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> call handle_err_nc_p(nfmpi_enddef(NCIDP),6)
>>>>>>>>
>>>>>>>> if (myrank==0) then
>>>>>>>> print *, ''
>>>>>>>> print *, 'array in nc_out_p subroutine'
>>>>>>>> print *, array
>>>>>>>> print *, ''
>>>>>>>> print *, '(/nx,ny/) ',(/nx,ny/)
>>>>>>>> print *, '(/istart,jstart/) ',(/istart,jstart/)
>>>>>>>> endif
>>>>>>>>
>>>>>>>> if (myrank==0) then
>>>>>>>> print *, ''
>>>>>>>> call handle_err_nc_p(nfmpi_inq_dimlen(NCIDP,x_dimIDP,dimlenx),
>>>>>>>> 7)
>>>>>>>> call handle_err_nc_p(nfmpi_inq_dimlen(NCIDP,y_dimIDP,dimleny),
>>>>>>>> 8)
>>>>>>>> print *, 'pnetCDF var dimensions in xy ',dimlenx,dimleny
>>>>>>>> print *, 'ntx_mpi ',ntx_mpi
>>>>>>>> print *, 'ntx ',ntx
>>>>>>>> print *, 'nx_mpi ',nx_mpi
>>>>>>>> print *, 'nx ',nx
>>>>>>>> endif
>>>>>>>>
>>>>>>>> ! this should have been possible in collective data mode, but
>>>>>>>> did
>>>>>>>> not seem to work !!!
>>>>>>>> call handle_err_nc_p(nfmpi_begin_indep_data(NCIDP),9)
>>>>>>>>
>>>>>>>> call
>>>>>>>> handle_err_nc_p(nfmpi_put_vara_real(NCIDP,array_varIDP,(/
>>>>>>>> istart,jstart/),(/nx_mpi,ny_mpi/),reshape(array,(/
>>>>>>>> nx_mpi,ny_mpi/))),10)
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> call handle_err_nc_p(nfmpi_end_indep_data(NCIDP),11)
>>>>>>>>
>>>>>>>> call handle_err_nc_p(nfmpi_close(NCIDP),12)
>>>>>>>>
>>>>>>>> end subroutine nc_out_p
>>>>>>>> !
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> ===============================================================
>>>>>>>>
>>>>>>>> subroutine handle_err_nc_s(status)
>>>>>>>>
>>>>>>>> use netcdf
>>>>>>>>
>>>>>>>> integer, intent (in ) :: status
>>>>>>>>
>>>>>>>> if ( status /= nf90_noerr ) then
>>>>>>>> print *, ""
>>>>>>>> print *, trim(nf90_strerror(status))
>>>>>>>> print *, "atham_netcdf error: program execution stopped"
>>>>>>>> print *, ""
>>>>>>>> stop
>>>>>>>> end if
>>>>>>>> return
>>>>>>>> end subroutine handle_err_nc_s
>>>>>>>> !
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> ===============================================================
>>>>>>>>
>>>>>>>> subroutine handle_err_nc_p(status,callID)
>>>>>>>>
>>>>>>>> ! use atham_module, only : myrank
>>>>>>>>
>>>>>>>> #include "/usr/local/parallel-netcdf-1.0.1/include/pnetcdf.inc"
>>>>>>>>
>>>>>>>> integer, intent (in ) :: status
>>>>>>>> integer, intent (in ) :: callID
>>>>>>>>
>>>>>>>> if ( status /= nf_noerr ) then
>>>>>>>> if (myrank==0) then
>>>>>>>> print *, ""
>>>>>>>> print *, trim(nfmpi_strerror(status))
>>>>>>>> print *, "atham_pnetcdf error: program execution stopped
>>>>>>>> ",callID
>>>>>>>> print *, ""
>>>>>>>> ! stop
>>>>>>>> end if
>>>>>>>> endif
>>>>>>>> return
>>>>>>>> end subroutine handle_err_nc_p
>>>>>>>> !
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> ===============================================================
>>>>>>>>
>>>>>>>> subroutine exec_stop(MESSAGE,ERRORCLASS)
>>>>>>>>
>>>>>>>> character(len=*), intent(in ) :: MESSAGE
>>>>>>>> integer(kint), intent(in ) :: ERRORCLASS
>>>>>>>>
>>>>>>>> integer(kint) :: ierror
>>>>>>>>
>>>>>>>> print *, 'ERROR: MPI execution stopped'
>>>>>>>> print *, 'error class: ', ERRORCLASS
>>>>>>>> print *, 'message: ', trim(MESSAGE)
>>>>>>>> call mpi_abort(mpi_comm_world,ierror)
>>>>>>>> stop
>>>>>>>>
>>>>>>>> end subroutine exec_stop
>>>>>>>> !
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> ===============================================================
>>>>>>>>
>>>>>>>> subroutine makename(prefix,icoor,name)
>>>>>>>> !------------------------------------------------------------------!
>>>>>>>> ! create character
>>>>>>>> name=prefix_(icoor(1)+1)_(icoor(2)+1) !
>>>>>>>> !------------------------------------------------------------------!
>>>>>>>> character(len=*), intent(in) :: prefix
>>>>>>>> integer(kint), dimension(2), intent(in) :: icoor
>>>>>>>> character(len=*), intent(out) :: name
>>>>>>>> !------------------------------------------------------------------!
>>>>>>>> ! local
>>>>>>>> variables !
>>>>>>>> !------------------------------------------------------------------!
>>>>>>>> character(len=25) :: numx,numy
>>>>>>>> integer(kint) :: lenx,leny,lenf
>>>>>>>>
>>>>>>>> write(numx,'(i5)') icoor(1)+1
>>>>>>>> numx=adjustl(numx)
>>>>>>>> lenx=len_trim(numx)
>>>>>>>>
>>>>>>>> write(numy,'(i5)') icoor(2)+1
>>>>>>>> numy=adjustl(numy)
>>>>>>>> leny=len_trim(numy)
>>>>>>>>
>>>>>>>> name=adjustl(prefix)
>>>>>>>> lenf=len_trim(name)
>>>>>>>>
>>>>>>>> name=name(1:lenf)//'_'//numx(1:lenx)//'_'//numy(1:leny)
>>>>>>>>
>>>>>>>> end subroutine makename
>>>>>>>> !
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> =
>>>>>>>> ===============================================================
>>>>>>>>
>>>>>>>>
>>>>>>>> end program mpi_setup_netCDFout
>>>>>>>> # MAKEFILE used for development module testing
>>>>>>>> # Alex Hoffmann
>>>>>>>>
>>>>>>>> # set netCDF path
>>>>>>>> NETCDF = /usr/local/netcdf-3.6.2-ifort
>>>>>>>> PNETCDF= /usr/local/parallel-netcdf-1.0.1
>>>>>>>>
>>>>>>>> # set compiler
>>>>>>>> F90 = mpif90 # mpif90 refers to a shell script that
>>>>>>>> includes
>>>>>>>> links
>>>>>>>> to the libraries
>>>>>>>>
>>>>>>>>
>>>>>>>> # set compiler flags
>>>>>>>> OFLAGS = -g
>>>>>>>> # set include flags
>>>>>>>> IFLAGS = -I$(NETCDF)/include/ -I$(PNETCDF)/include/
>>>>>>>> # set pre-processing flags
>>>>>>>> #PFLAGS = -fpp
>>>>>>>> # concatenate flags
>>>>>>>> FLAGS = $(OFLAGS) $(IFLAGS) $(PFLAGS)
>>>>>>>> # set library flags
>>>>>>>> LFLAGS = -L$(NETCDF)/lib -lnetcdf -L$(PNETCDF)/lib -lpnetcdf
>>>>>>>>
>>>>>>>> RM = rm -f
>>>>>>>>
>>>>>>>> #MAIN =
>>>>>>>> MAIN = mpi_setup_netCDFout
>>>>>>>> OBJECTS =
>>>>>>>>
>>>>>>>>
>>>>>>>> all: $(MAIN)
>>>>>>>>
>>>>>>>>
>>>>>>>> $(MAIN): $(MAIN).F90 $(OBJECTS)
>>>>>>>> $(F90) $(FLAGS) -o $@ $^ $(LFLAGS)
>>>>>>>>
>>>>>>>> clean:
>>>>>>>> $(RM) *.o *.mod
>>>>>>>
>>>>>>
>>>>>> --
>>>>>> ______________________________________________________________________
>>>>>>
>>>>>> Alex HOFFMANN PhD Candidate
>>>>>> Centre for Atmospheric Science
>>>>>> Department of Geography University of Cambridge
>>>>>> Downing Place, Cambridge CB2 3EN, UK
>>>>>> e-mail: ah519 at cam.ac.uk tel: +44 (0)1223 766581
>>>>>> www.geog.cam.ac.uk/people/hoffmann/
>>>>>> ______________________________________________________________________
>>>>>>
>>>>>
>>>>
>>>> --
>>>> ______________________________________________________________________
>>>>
>>>> Alex HOFFMANN PhD Candidate
>>>> Centre for Atmospheric Science
>>>> Department of Geography University of Cambridge
>>>> Downing Place, Cambridge CB2 3EN, UK
>>>> e-mail: ah519 at cam.ac.uk tel: +44 (0)1223 766581
>>>> www.geog.cam.ac.uk/people/hoffmann/
>>>> ______________________________________________________________________
>>>>
>>>
>>
>> --
>> ______________________________________________________________________
>>
>> Alex HOFFMANN PhD Candidate
>> Centre for Atmospheric Science
>> Department of Geography University of Cambridge
>> Downing Place, Cambridge CB2 3EN, UK
>> e-mail: ah519 at cam.ac.uk tel: +44 (0)1223 766581
>> www.geog.cam.ac.uk/people/hoffmann/
>> ______________________________________________________________________
>>
>
More information about the parallel-netcdf
mailing list