param_m.f90 9.45 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62
!=========================================================================================
  module param_m
!=========================================================================================
! module for param variables

  implicit none
!-----------------------------------------------------------------------------------------
! integers

  character run_mode*20        !mode in which to run

  integer numdim      !number of dimensions of problem

  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 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 nxyz        !miscelaneous for compatibility checks
  integer nxyz2       !miscelaneous for compatibility checks
  integer nxm         !miscelaneous for compatibility checks
  integer nym         !miscelaneous for compatibility checks
  integer nzm         !miscelaneous for compatibility checks

  integer nxyz_g      !miscelaneous for compatibility checks
  integer nxyz2_g     !miscelaneous for compatibility checks
  integer nxm_g       !miscelaneous for compatibility checks
  integer nym_g       !miscelaneous for compatibility checks
  integer nzm_g       !miscelaneous for compatibility checks

  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 nvar_tot    !total number of variables (equations) = 5 + nsc

  integer n_elem      !number of elements in chemical mechanism
  integer n_spec      !number of chemical species including N2
  integer n_reac      !total number of reactions

  integer ntr         !number of reactions with third-body efficiencies  

  integer n_reg       !number of registers used in runge-kutta integration

  integer periodic_x  !switch for periodicity in x-direction
  integer periodic_y  !switch for periodicity in y-direction
  integer periodic_z  !switch for periodicity in z-direction

  integer vary_in_x   !switch to turn on x-direction
  integer vary_in_y   !switch to turn on y-direction
  integer vary_in_z   !switch to turn on z-direction

  integer iorder      !order of derivatives
  integer iforder     !order of filter

  integer io_method

63
  character*100 output_dir
64

65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289
  character*8  dat_1, dat_2     !for start and end date (wall clock)
  character*10 tim_1, tim_2     !for start and end time (wall clock)

!-----------------------------------------------------------------------------------------
  contains
!========================================================================================
  subroutine initialize_param(io,myid,ierr,gcomm)
!========================================================================================
! routine sets various parameters
!----------------------------------------------------------------------------------------
  implicit none
!-----------------------------------------------------------------------------------------
  include 'mpif.h'
!----------------------------------------------------------------------------------------
! declarations passed in

  integer io
  integer myid, ierr, gcomm

! local declarations

!!$  integer nxnynz(3)
  integer i, iflag
!----------------------------------------------------------------------------------------
! write header

  if(myid.eq.0) then
    write(io,*) 'initializing param module...'
    write(io,*)
  endif
!----------------------------------------------------------------------------------------
! 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

      if( myid.eq.0 ) then
        write(io,*) ' Grid Pts in X dimension ',  nx_g,  &
                    ' are not exactly divisible by ', npx ,  &
                    ' number of PEs in the X dimension '
      endif

      iflag=1

    endif

!   error trapping for number of grid points in x-direction

    if(mod(ny_g,npy).ne.0 ) then

      if(myid.eq.0) then
        write(io,*) ' Grid Pts in Y dimension ',  ny_g,  &
                    ' are not exactly divisible by ', npy ,  &
                    ' number of PEs in the Y dimension '
      endif

      iflag=1

    endif

!   error trapping for number of grid points in x-direction

    if(mod(nz_g,npz).ne.0 ) then

      if( myid.eq.0 ) then
        write(io,*) ' Grid Pts in Z dimension ',  nz_g,  &
                    ' are not exactly divisible by ', npz ,  &
                    ' number of PEs in the Z dimension '
      endif

      iflag=1

    endif

!   check for domain size conflict with stencils size in x-direction

    if((vary_in_x.eq.1).and.(nx_g.lt.9)) then

      if(myid.eq.0) then
        write(io,*) 'input error: nx_g < 9 and vary_in_x = 1'
      endif

      iflag=1

    endif

!   check for domain size conflict with stencils size in y-direction

    if((vary_in_y.eq.1).and.(ny_g.lt.9)) then

      if(myid.eq.0) then
        write(io,*) 'input error: ny_g < 9 and vary_in_y = 1'
      endif

      iflag=1

    endif

!   check for domain size conflict with stencils size in z-direction

    if((vary_in_z.eq.1).and.(nz_g.lt.9)) then

      if(myid.eq.0) then
        write(io,*) 'input error: nz_g < 9 and vary_in_z = 1'
      endif

      iflag=1

    endif

!   set local number of grid points

    nx = nx_g / npx
    ny = ny_g / npy
    nz = nz_g / npz

!   set some other stuff

    nxyz = max(nx, ny, nz);      nxyz_g = max(nx_g, ny_g, nz_g)
    nxyz2= 2*nxyz         ;      nxyz2_g= 2*nxyz_g      
    nxm  = max(1,nx - 1)  ;      nxm_g  = max(1,nx_g - 1)   
    nym  = max(1,ny - 1)  ;      nym_g  = max(1,ny_g - 1)   
    nzm  = max(1,nz - 1)  ;      nzm_g  = max(1,nz_g - 1)   

!   set chemistry parameters based on chemkin initialization

    call set_number_elem_spec_reac(n_elem,n_spec,n_reac,myid,io)
    call set_number_third_body_reactions(ntr,io)

    n_reac=n_reac*2   !convert n_reac for DNS/getrates purposes (2x)
    nsc=n_spec-1      !set number of chemical species for DNS purposes
    nvar_tot=nsc+5    !number of species + 5 (density, energy, momentum)

  endif
!----------------------------------------------------------------------------------------
! check status of error

  call MPI_Bcast(iflag,1,MPI_INTEGER,0,gcomm,ierr)

  if(iflag.eq.1) call terminate_run(io,0)  !must be called by all processors
!----------------------------------------------------------------------------------------
! broadcast parameters
!----------------------------------------------------------------------------------------
! broadcast local grid parameters

  call MPI_Bcast(nx       ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(ny       ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nz       ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nxyz     ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nxyz2    ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nxm      ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nym      ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nzm      ,1, MPI_INTEGER, 0, gcomm, ierr)

! broadcast global grid parameters

  call MPI_Bcast(nxyz_g  ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nxyz2_g ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nxm_g   ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nym_g   ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nzm_g   ,1, MPI_INTEGER, 0, gcomm, ierr)

! broadcast chemistry parameters

  call MPI_Bcast(nsc     ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(nvar_tot,1, MPI_INTEGER, 0, gcomm, ierr)

  call MPI_Bcast(n_elem  ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(n_spec  ,1, MPI_INTEGER, 0, gcomm, ierr)
  call MPI_Bcast(n_reac  ,1, MPI_INTEGER, 0, gcomm, ierr)

  call MPI_Bcast(ntr     ,1, MPI_INTEGER, 0, gcomm, ierr)

! broadcast runge-kutta parameters

  call MPI_Bcast(n_reg   ,1, MPI_INTEGER, 0, gcomm, ierr)

! sync processors

  call MPI_Barrier( gcomm,ierr )
!----------------------------------------------------------------------------------------
! miscellaneous initializations to be set on all processors
!----------------------------------------------------------------------------------------
! set number of dimensions

  numdim = 0
  if(vary_in_x.eq.1) numdim = numdim + 1
  if(vary_in_y.eq.1) numdim = numdim + 1
  if(vary_in_z.eq.1) numdim = numdim + 1

!!$  nxnynz(1)=nx
!!$  nxnynz(2)=ny
!!$  nxnynz(3)=nz
!!$
!!$! set number of waves (nw) and nh
!!$
!!$  if (numdim == 3) then
!!$    nw=int(min(nx_g,ny_g,nz_g)/2)-1
!!$  elseif(numdim == 2) then
!!$    if( vary_in_x==1 .and. vary_in_y==1 ) then
!!$      nw = int(min(nx_g,ny_g)/2)-1
!!$    elseif (vary_in_x==1 .and. vary_in_z==1) then
!!$      nw = int(min(nx_g,nz_g)/2)-1
!!$    elseif (vary_in_y==1 .and. vary_in_z==1) then
!!$      nw = int(min(ny_g,nz_g)/2)-1
!!$    endif
!!$  endif
!!$  nh = 1 + nw / 2

!----------------------------------------------------------------------------------------
! write header

  if(myid.eq.0) then
    call write_header(io,'-')
  endif
!----------------------------------------------------------------------------------------
  return
  end subroutine initialize_param
!----------------------------------------------------------------------------------------
  end module param_m