problem with put_att (maximum number of attributes)
Wei-keng Liao
wkliao at ece.northwestern.edu
Fri Feb 20 15:27:33 CST 2009
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 at 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
More information about the parallel-netcdf
mailing list