Commit 05f92d89 authored by Francois Tessier's avatar Francois Tessier

Add S3D-IO benchmark

parent a481d569
Copyright 2003-2013 Northwestern University
Portions of this software were developed by the Sandia National Laboratory.
Access and use of this software shall impose the following obligations
and understandings on the user. The user is granted the right, without
any fee or cost, to use, copy, modify, alter, enhance and distribute
this software, and any derivative works thereof, and its supporting
documentation for any purpose whatsoever, provided that this entire
notice appears in all copies of the software, derivative works and
supporting documentation. Further, Northwestern University requests
that the user credit Northwestern University in any publications that
result from the use of this software or in any product that includes
this software. The name Northwestern University, however, may not be
used in any advertising or publicity to endorse or promote any
products or commercial entity unless specific written permission is
obtained from Northwestern University. The user also understands that
Northwestern University is not obligated to provide the user with
any support, consulting, training or assistance of any kind with regard
to the use, operation and performance of this software nor to provide
the user with any updates, revisions, new versions or "bug fixes."
THIS SOFTWARE IS PROVIDED BY NORTHWESTERN UNIVERSITY "AS IS" AND ANY
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL NORTHWESTERN UNIVERSITY BE
LIABLE FOR ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY
DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS,
WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION,
ARISING OUT OF OR IN CONNECTION WITH THE ACCESS, USE OR PERFORMANCE
OF THIS SOFTWARE.
#
# Copyright (C) 2013, Northwestern University
# See COPYRIGHT notice in top-level directory.
#
# $Id: Makefile 3485 2015-12-27 00:06:31Z wkliao $
#
#
# Please change the following variables:
# MPIF90 -- MPI Fortran compiler
# FCFLAGS -- Compile flag
# PNETCDF_DIR -- PnetCDF library installation directory
#
MPIF90 = mpif90
FCFLAGS = -Wall -g
PNETCDF_DIR = $(HOME)/PnetCDF
PNETCDF_DIR = $(HOME)
COMPILE_F90 = $(MPIF90) $(FCFLAGS) $(INC) -c
LINK = $(MPIF90) $(FCFLAGS)
INC = -I$(PNETCDF_DIR)/include
LIBS = -L$(PNETCDF_DIR)/lib -lpnetcdf
SRCS = runtime_m.f90 \
param_m.f90 \
topology_m.f90 \
variables_m.f90 \
io_profiling_m.f90 \
pnetcdf_m.f90 \
init_field.f90 \
io.f90 \
random_number.f90 \
solve_driver.f90 \
main.f90
OBJS = $(SRCS:.f90=.o)
MODS = $(SRCS:.f90=.mod)
TARGET = s3d_io.x
all: $(TARGET)
%.o:%.f90
$(COMPILE_F90) $<
$(TARGET): $(OBJS)
$(LINK) $(OBJS) -o $(TARGET) $(LIBS)
PACKAGE_NAME = s3d-io-pnetcdf-1.1
PACKING_LIST = $(SRCS) Makefile README COPYRIGHT RELEASE_NOTE
dist:
/bin/rm -rf $(PACKAGE_NAME) $(PACKAGE_NAME).tar.gz
mkdir $(PACKAGE_NAME)
cp $(PACKING_LIST) $(PACKAGE_NAME)
tar -cf $(PACKAGE_NAME).tar $(PACKAGE_NAME)
gzip $(PACKAGE_NAME).tar
/bin/rm -rf $(PACKAGE_NAME)
clean:
/bin/rm -f $(OBJS) $(MODS) $(TARGET)
distclean: clean
/bin/rm -rf $(PACKAGE_NAME).tar.gz $(PACKAGE_NAME)
#
# Copyright (C) 2013, Northwestern University
# See COPYRIGHT notice in top-level directory.
#
# $Id: README 3457 2015-11-21 23:07:56Z wkliao $
This benchmark programs is the I/O kernel of S3D combustion simulation code.
http://exactcodesign.org/ There are several I/O methods implemented in S3D.
This software only contains the method of Parallel NetCDF.
S3D is a continuum scale first principles direct numerical simulation code
which solves the compressible governing equations of mass continuity, momenta,
energy and mass fractions of chemical species including chemical reactions.
Readers are referred to the published paper below. J. Chen, A. Choudhary, B.
de Supinski, M. DeVries, E. Hawkes, S. Klasky, W. Liao, K. Ma, J. Crummey, N.
Podhorszki, R. Sankaran, S. Shende, and C. Yoo. Teras-cale Direct Numerical
Simulations of Turbulent Combustion Using S3D. In Computational Science and
Discovery Volume 2, January 2009.
I/O pattern:
A checkpoint is performed at regular intervals, and its data consist of 8-byte
three-dimensional arrays. At each checkpoint, four global arrays, representing
mass, velocity, pressure, and temperature, respectively, are written to a newly
created file in the canonical order. Mass and velocity are four-dimensional
arrays while pressure and temperature are three-dimensional arrays. All four
arrays share the same size for the lowest three spatial dimensions X, Y, and Z,
which are partitioned among MPI processes in a block-block-block fashion. For
the mass and velocity arrays, the length of the fourth dimension is 11 and 3,
respectively. The fourth dimension, the most significant one, is not
partitioned. As the number of MPI processes increases, the aggregate I/O
amount proportionally increases as well.
For more detailed description of the data partitioning and I/O patterns,
please refer to the following paper.
W. Liao and A. Choudhary. Dynamically Adapting File Domain Partitioning
Methods for Collective I/O Based on Underlying Parallel File System
Locking Protocols. In the Proceedings of International Conference for
High Performance Computing, Networking, Storage and Analysis, Austin,
Texas, November 2008.
To compile:
Edit Makefile and set/change variables:
MPIF90 - MPI Fortran 90 compiler
FCFLAGS - compile flags
PNETCDF_DIR - the path of PnetCDF library
(1.4.0 and higher is required)
For example:
MPIF90 = mpif90
FCFLAGS = -O2
PNETCDF_DIR = ${HOME}/PnetCDF
To run:
Usage: s3d_io.x nx_g ny_g nz_g npx npy npz dir_path
There are 9 command-line arguments:
nx_g - GLOBAL grid size along X dimension
ny_g - GLOBAL grid size along Y dimension
nz_g - GLOBAL grid size along Z dimension
npx - number of MPI processes along X dimension
npy - number of MPI processes along Y dimension
npz - number of MPI processes along Z dimension
method - 0: using PnetCDF blocking APIs, 1: nonblocking APIs
restart - restart from reading a previous written file (True/False)
dir_path - the directory name to store the output files
To change the number of checkpoint dumps (default is set to 5), edit
file param_m.f90 and set a different value for i_time_end:
i_time_end = 5 ! number of checkpoints (also number of output files)
The contents of all variables written to files are set to random numbers.
This setting can be disabled by comment out the line below in file
solve_driver.f90
call random_set
Example run command:
For a test run with small data size and a short return time, here is an
example command for running on 4 MPI processes.
mpiexec -n 4 ./s3d_io.x 10 10 10 2 2 1 1 F .
The command below runs on 4096 MPI processes with the global array
of size 800x800x800 and local array of size 50x50x50, output directory
/scratch1/scratchdirs/wkliao/FS_1M_96 using nonblocking APIs, and without
restart.
mpiexec -l -n 512 ./s3d_io.x 800 800 800 16 16 16 1 F /scratch1/scratchdirs/wkliao/FS_1M_96
Example output from stdout:
++++ I/O is done through PnetCDF ++++
I/O method : nonblocking APIs
Run with restart : False
No. MPI processes : 4096
Global array size : 800 x 800 x 800
output file path : /scratch1/scratchdirs/wkliao/FS_1M_96
file striping count : 96
file striping size : 1048576 bytes
-----------------------------------------------
Time for open : 0.11 sec
Time for read : 0.00 sec
Time for write : 18.04 sec
Time for close : 0.02 sec
no. read calls : 0 per process
no. write calls : 20 per process
total read amount : 0.00 GiB
total write amount : 305.18 GiB
read bandwidth : 0.00 MiB/s
write bandwidth : 17318.78 MiB/s
-----------------------------------------------
total I/O amount : 305.18 GiB
total I/O time : 18.17 sec
I/O bandwidth : 17201.53 MiB/s
Questions/Comments:
email: wkliao@eecs.northwestern.edu
Release note
version 1.1
* Add command-line option to choose PnetCDF blocking APIs or nonblocking APIs.
* Add command-line option for whether to run from a restart file.
* Reorganize performance result output.
!
! Copyright (C) 2013, Northwestern University
! See COPYRIGHT notice in top-level directory.
!
! $Id: init_field.f90 2192 2013-11-14 19:48:08Z wkliao $
!----< initialize_field >----------------------------------------
subroutine initialize_field
! initializes primitive and solution variables amoung other things
use variables_m, only: allocate_variables_arrays
use runtime_m, only: restart, time_save, time_save_inc
use topology_m, only: gcomm
implicit none
integer err
! allocate checkpointing variable arrays
call allocate_variables_arrays(1)
time_save_inc = 1.0e+5
time_save = time_save_inc
! Restart code from prevous data files.
if (restart) then
call MPI_Barrier(gcomm,err)
call read_savefile
endif
end subroutine initialize_field
!
! Copyright (C) 2013, Northwestern University
! See COPYRIGHT notice in top-level directory.
!
! $Id: io.f90 2192 2013-11-14 19:48:08Z wkliao $
!----< write_savefile >------------------------------------------
subroutine write_savefile
! writes data file for post-processing and restarting
use mpi
use runtime_m, only: time, time_ref, run_title
use io_profiling_m, only: dir_path
use pnetcdf_m, only: pnetcdf_write
implicit none
! local variables
character*10 time_ext
character*256 filename
! set file extention string
write(time_ext,'(1pe10.3)') time*time_ref
! set file name
filename = trim(dir_path)//'/'//trim(run_title)//'.'// &
trim(adjustl(time_ext))//'.field.nc'
! call PnetCDF APIs to write data to file
call pnetcdf_write(filename)
end subroutine write_savefile
!----< read_savefile >-------------------------------------------
subroutine read_savefile
! reads restart file
use mpi
use runtime_m, only: time, time_ref, run_title
use io_profiling_m, only: dir_path
use pnetcdf_m, only: pnetcdf_read
use topology_m, only : myid
implicit none
! local variables
character*10 time_ext
character*256 filename
logical exist
! set file extention strings
write(time_ext,'(1pe10.3)') time*time_ref
! set file name
filename = trim(dir_path)//'/'//trim(run_title)//'.'// &
trim(adjustl(time_ext))//'.field.nc'
! inquire about file existence
inquire(file=trim(filename),exist=exist)
if (.NOT. exist) then
if (myid .EQ. 0) then
print*, 'restart file does not exist ', trim(filename)
print*, 'skip reading restart files'
endif
return
endif
! call PnetCDF APIs to read data from file
call pnetcdf_read(filename)
end subroutine read_savefile
!
! Copyright (C) 2013, Northwestern University
! See COPYRIGHT notice in top-level directory.
!
! $Id: io_profiling_m.f90 2713 2014-08-18 22:51:07Z wkliao $
module io_profiling_m
! module for Read-Write Restart files
use mpi
implicit none
character*256 :: dir_path ! directory name for output files
integer file_info, info_used
integer cb_nodes ! collective buffering nodes
integer read_num ! number of read calls
integer write_num ! number of write calls
double precision read_amount ! total read amount
double precision write_amount ! total write amount
double precision io_amount ! total I/O amount
double precision openT, writeT, readT, closeT
! time for open, write, read, and close
contains
! ----------------------------------------------------------------
! print the MPI info objects to stdout
! ----------------------------------------------------------------
subroutine print_io_hints(info)
implicit none
integer, intent(in) :: info
! local variables
character*(MPI_MAX_INFO_VAL) key, value
integer i, nkeys, valuelen, err
logical flag
1001 FORMAT(' ',A32,' = ',A)
call MPI_Info_get_nkeys(info, nkeys, err)
print *, '---- MPI file info used ----'
do i=0, nkeys-1
key(:) = ' '
call MPI_Info_get_nthkey(info, i, key, err)
call MPI_Info_get(info, key, MPI_MAX_INFO_VAL, value, flag, err)
call MPI_Info_get_valuelen(info, key, valuelen, flag, err)
value(valuelen+1:) = ' '
if (key(len_trim(key):len_trim(key)) .EQ. char(0)) &
key(len_trim(key):) = ' '
print 1001, trim(key), trim(value)
enddo
print *
end subroutine print_io_hints
! ----------------------------------------------------------------
! get the file striping information from the MPI info objects
! ----------------------------------------------------------------
subroutine get_file_striping(info, striping_factor, striping_unit)
implicit none
integer, intent(in) :: info
integer, intent(out) :: striping_factor
integer, intent(out) :: striping_unit
! local variables
character*(MPI_MAX_INFO_VAL) key, value
integer i, nkeys, valuelen, err
logical flag
call MPI_Info_get_nkeys(info, nkeys, err)
do i=0, nkeys-1
key(:) = ' '
call MPI_Info_get_nthkey(info, i, key, err)
call MPI_Info_get(info, key, MPI_MAX_INFO_VAL, value, flag, err)
call MPI_Info_get_valuelen(info, key, valuelen, flag, err)
value(valuelen+1:) = ' '
if (key(len_trim(key):len_trim(key)) .EQ. char(0)) &
key(len_trim(key):) = ' '
if (trim(key) .EQ. 'striping_factor') &
read(value, '(i10)') striping_factor
if (trim(key) .EQ. 'striping_unit') &
read(value, '(i10)') striping_unit
enddo
end subroutine get_file_striping
!----< set_io_hints() >-------------------------------------------
subroutine set_io_hints(flag)
implicit none
integer, intent(in) :: flag
! local variables
integer err
character*16 int_str
! free up info and file type
if (flag .EQ. -1) then
if (file_info .NE. MPI_INFO_NULL) then
call MPI_Info_free(file_info, err)
file_info = MPI_INFO_NULL
endif
return
endif
read_amount = 0.0
write_amount = 0.0
io_amount = 0.0
read_num = 0
write_num = 0
openT = 0.0
writeT = 0.0
readT = 0.0
closeT = 0.0
info_used = MPI_INFO_NULL
! set MPI I/O hints for performance enhancement
call MPI_Info_create(file_info, err)
! disable ROMIO data sieving
call MPI_Info_set(file_info, 'romio_ds_write', 'disable', err)
call MPI_Info_set(file_info, 'romio_ds_read', 'disable', err)
call MPI_Info_set(file_info, 'romio_no_indep_rw', 'true', err)
! set the number of aggregate I/O nodes (for advanced users)
write(int_str,'(I16)') cb_nodes
! call MPI_Info_set(file_info, 'cb_nodes', int_str, err)
! set PnetCDF hints
! call MPI_Info_set (file_info, 'nc_header_align_size', '512', err)
! call MPI_Info_set (file_info, 'nc_var_align_size', '512', err)
! call MPI_Info_set (file_info, 'nc_header_read_chunk_size', '262144', err)
end subroutine set_io_hints
!----< print_io_performance() >-----------------------------------
subroutine print_io_performance(io_time)
use topology_m, only : gcomm, myid, npes
use param_m, only : nx_g, ny_g, nz_g
use runtime_m, only : method, restart
implicit none
double precision, intent(inout) :: io_time
! local variables
double precision d_tmp, read_bandwidth, write_bandwidth, io_bandwidth
integer striping_factor, striping_unit, err
call MPI_Reduce(openT, d_tmp, 1, MPI_DOUBLE_PRECISION, MPI_MAX, 0, gcomm, err)
openT = d_tmp
call MPI_Reduce(writeT, d_tmp, 1, MPI_DOUBLE_PRECISION, MPI_MAX, 0, gcomm, err)
writeT = d_tmp
call MPI_Reduce(readT, d_tmp, 1, MPI_DOUBLE_PRECISION, MPI_MAX, 0, gcomm, err)
readT = d_tmp
call MPI_Reduce(closeT, d_tmp, 1, MPI_DOUBLE_PRECISION, MPI_MAX, 0, gcomm, err)
closeT = d_tmp
call MPI_Reduce(io_time, d_tmp, 1, MPI_DOUBLE_PRECISION, MPI_MAX, 0, gcomm, err)
io_time = d_tmp
call MPI_Reduce(read_amount, d_tmp, 1, MPI_DOUBLE_PRECISION, MPI_SUM, 0, gcomm, err)
read_amount = d_tmp
call MPI_Reduce(write_amount, d_tmp, 1, MPI_DOUBLE_PRECISION, MPI_SUM, 0, gcomm, err)
write_amount = d_tmp
io_amount = read_amount + write_amount
if (myid == 0) then
io_bandwidth = 0.0
read_bandwidth = 0.0
write_bandwidth = 0.0
if (io_time > 0) io_bandwidth = io_amount / io_time
if (readT > 0) read_bandwidth = read_amount / readT
if (writeT > 0) write_bandwidth = write_amount / writeT
io_bandwidth = io_bandwidth / 1048576.0
io_amount = io_amount / 1073741824.0
read_bandwidth = read_bandwidth / 1048576.0
read_amount = read_amount / 1073741824.0
write_bandwidth = write_bandwidth / 1048576.0
write_amount = write_amount / 1073741824.0
call print_io_hints(info_used)
striping_factor = 0
striping_unit = 0
call get_file_striping(info_used, striping_factor, striping_unit)
2000 FORMAT(A)
2001 FORMAT(A, A)
2002 FORMAT(A ,I7)
2003 FORMAT(A ,I7, A)
2004 FORMAT(A ,F10.2, A)
2005 FORMAT(A, i7,' x',i7,' x',i7)
write(6,*) '++++ I/O is done through PnetCDF ++++'
if (method .EQ. 0) then
write(6,*) 'I/O method : blocking APIs'
else
write(6,*) 'I/O method : nonblocking APIs'
endif
if (restart) then
write(6,*) 'Run with restart : True'
else
write(6,*) 'Run with restart : False'
endif
write(6, 2002) ' No. MPI processes : ', npes
write(6, 2005) ' Global array size : ', nx_g, ny_g, nz_g
write(6, 2001) ' output file path : ',trim(dir_path)
write(6, 2002) ' file striping count : ',striping_factor
write(6, 2003) ' file striping size : ',striping_unit, ' bytes'
write(6, 2000) ' -----------------------------------------------'
write(6, 2004) ' Time for open : ',openT, ' sec'
write(6, 2004) ' Time for read : ',readT, ' sec'
write(6, 2004) ' Time for write : ',writeT, ' sec'
write(6, 2004) ' Time for close : ',closeT, ' sec'
write(6, 2003) ' no. read calls : ',read_num, ' per process'
write(6, 2003) ' no. write calls : ',write_num,' per process'
write(6, 2004) ' total read amount : ',read_amount , ' GiB'
write(6, 2004) ' total write amount : ',write_amount , ' GiB'
write(6, 2004) ' read bandwidth : ',read_bandwidth , ' MiB/s'
write(6, 2004) ' write bandwidth : ',write_bandwidth , ' MiB/s'
write(6, 2000) ' -----------------------------------------------'
write(6, 2004) ' total I/O amount : ',io_amount , ' GiB'
write(6, 2004) ' total I/O time : ',io_time , ' sec'
write(6, 2004) ' I/O bandwidth : ',io_bandwidth, ' MiB/s'
endif
if (info_used .NE. MPI_INFO_NULL) &
call MPI_Info_free(info_used, err)
end subroutine print_io_performance
end module io_profiling_m
!
! Copyright (C) 2013, Northwestern University
! See COPYRIGHT notice in top-level directory.
!
! $Id: main.f90 3457 2015-11-21 23:07:56Z wkliao $
!----< main >-----------------------------------------------------
program main
use mpi
use param_m, only: npx, npy, npz
use param_m, only: initialize_param
use topology_m, only: gcomm, npes, myid, initialize_topology
implicit none
integer err
logical isArgvRight
call MPI_Init(err)
call MPI_Comm_rank(MPI_COMM_WORLD, myid, err)
call MPI_Comm_size(MPI_COMM_WORLD, npes, err)
call MPI_Comm_dup (MPI_COMM_WORLD, gcomm,err)
call read_command_line_arg(isArgvRight)
if (.NOT. isArgvRight) goto 999
! initialize parameters: nx, ny, nz, nsc, n_spec
call initialize_param(myid, gcomm)
! intialize MPI process topology
call initialize_topology(npx,npy,npz)
! main computation task is here
call MPI_Barrier(gcomm,err)
call solve_driver
999 call MPI_Comm_free(gcomm,err)
call MPI_Finalize(err)
end program main
!----< read_command_line_arg >------------------------------------
subroutine read_command_line_arg(isArgvRight)
use mpi
use param_m, only: nx_g, ny_g, nz_g, npx, npy, npz
use param_m, only: initialize_param
use runtime_m, only: method, restart
use topology_m, only: gcomm, npes, myid, initialize_topology
use io_profiling_m, only: dir_path
implicit none
character(len=128) executable
logical isArgvRight
! declare external functions
integer IARGC
! local variables for reading command-line arguments
character(len = 256) :: argv(9)
integer i, argc, int_argv(7), err
! Only root process reads command-line arguments
if (myid .EQ. 0) then
isArgvRight = .TRUE.
call getarg(0, executable)
argc = IARGC()
if (argc .NE. 9) then
print *, 'Usage: ',trim(executable), &
' nx_g ny_g nz_g npx npy npz method restart dir_path'
isArgvRight = .FALSE.
else
do i=1, argc-2
call getarg(i, argv(i))
read(argv(i), FMT='(I16)') int_argv(i)
enddo
call getarg(argc-1, argv(argc-1))
read(argv(argc-1), FMT='(L)') restart
call getarg(argc, argv(argc))
dir_path = argv(argc)
nx_g = int_argv(1)
ny_g = int_argv(2)
nz_g = int_argv(3)
npx = int_argv(4)
npy = int_argv(5)
npz = int_argv(6)
method = int_argv(7)
endif
endif
! broadcast if arguments are valid
call MPI_Bcast(isArgvRight, 1, MPI_LOGICAL, 0, gcomm, err)
if (.NOT. isArgvRight) return
call MPI_Bcast(nx_g, 1, MPI_INTEGER, 0, gcomm, err)
call MPI_Bcast(ny_g, 1, MPI_INTEGER, 0, gcomm, err)
call MPI_Bcast(nz_g, 1, MPI_INTEGER, 0, gcomm, err)
call MPI_Bcast(npx, 1, MPI_INTEGER, 0, gcomm, err)
call MPI_Bcast(npy, 1, MPI_INTEGER, 0, gcomm, err)
call MPI_Bcast(npz, 1, MPI_INTEGER, 0, gcomm, err)
call MPI_Bcast(method, 1, MPI_INTEGER, 0, gcomm, err)
call MPI_Bcast(restart, 1, MPI_LOGICAL, 0, gcomm, err)
call MPI_Bcast(dir_path, 256, MPI_CHARACTER, 0, gcomm, err)
end subroutine read_command_line_arg
!
! Copyright (C) 2013, Northwestern University
! See COPYRIGHT notice in top-level directory.
!
! $Id: param_m.f90 3485 2015-12-27 00:06:31Z wkliao $
module param_m
implicit none
integer nx_g ! global number of grid points in x-direction
integer ny_g ! global number of grid points in y-direction
integer nz_g ! global number of grid points in z-direction
integer nx ! local number of grid points in x-direction
integer ny ! local number of grid points in y-direction
integer nz ! local number of grid points in z-direction
integer npx ! number of processors in x-direction
integer npy ! number of processors in y-direction
integer npz ! number of processors in z-direction
integer nsc ! number of chemical species (excluding N2)
integer n_spec ! number of chemical species including N2
contains
!----< initialize_param() ---------------------------------------
subroutine initialize_param(myid,gcomm)
! sets various parameters
use mpi
use runtime_m, only: i_time_end, run_title
implicit none
! declarations passed in
integer myid, gcomm
! local declarations
integer err, iflag
i_time_end = 5 ! number of checkpoints (also number of output files)
run_title = 'pressure_wave_test' ! prefix of output file names
! set iflag
iflag=0
if (myid .EQ. 0) then
! error trapping for number of grid points in x-direction
if (mod(nx_g,npx) .NE. 0 ) then
print*, ' Grid Pts in X dimension ', nx_g, &
' are not exactly divisible by ', npx , &
' number of PEs in the X dimension '
iflag=1
endif
! error trapping for nu