main.f90
Go to the documentation of this file.
1 !*************************************************************************
2 ! This file is part of HERCULES.
3 !
4 ! HERCULES is free software; you can redistribute it and/or modify
5 ! it under the terms of the GNU General Public License as published by
6 ! the Free Software Foundation; either version 3 of the License, or
7 ! (at your option) any later version.
8 !
9 ! HERCULES is distributed in the hope that it will be useful,
10 ! but WITHOUT ANY WARRANTY; without even the implied warranty of
11 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 ! GNU General Public License for more details.
13 !
14 ! You should have received a copy of the GNU General Public License
15 ! along with this program. If not, see <http://www.gnu.org/licenses/>.
16 !
17 ! Copyright 2015, Ping He, MEAS, NCSU
18 !*************************************************************************
19 
20 
21 
22 program main
23 
26  use mpi
27  use decomp_2d
28 
29  implicit none
30 
31  ! read and initiate all the parameters
33 
34  ! initiate mpi
35  call initiate_mpi
36 
37  ! create a folder to store the output results
38  if (myid .eq. 0) then
39 
40  call system("mkdir -p ./results")
41 
42  if (isxy2d .eq. 1 .or. isxz2d .eq. 1 .or. isyz2d .eq. 1) then
43  call system("mkdir -p ./results/slices")
44  endif
45 
46  endif
47 
48  ! initiate the 2decomp lib
49  call decomp_2d_init(nx,ny,nz,p_row,p_col)
50 
51  ! initiate some indices, e.g., istr3, iend3, i_offset, etc.
52  ! we need to do this after decomp_2d_init
53  call initiate_indices
54 
55  ! welcome info
56  if (myid .eq. 0) then
57  call welcome_info
58  endif
59 
60  ! main subroutine for HERCULES
61  call drive
62 
63  ! goodbye info
64  if (myid .eq. 0) then
65  call goodby_info
66  endif
67 
68  ! finalize mpi
69  call finalize_mpi
70  call decomp_2d_finalize
71 
72  stop
73 
74 
75 
76 
77 contains
78 
79 
80  ! this function prints some information on the screen
81  subroutine welcome_info
82 
83  use mpi
85 
86  implicit none
87 
88  ! these are start date, time, etc
89  character(len=8) :: date1
90  character(len=10) :: time1
91  character(len=5) :: zone1
92  integer,dimension(8) :: d1
93 
94  print *, '!**********************************************************!'
95  print *, '! HERCULES V1.1 !'
96  print *, '!**********************************************************!'
97  print *, ''
98 
99  ! channel configuration
100  if (ochannel .eq. 0) then
101  print *, 'Configurations: Closed Channel'
102  elseif (ochannel .eq. 1) then
103  print *, 'Configurations: Open Channel'
104  else
105  print *,'ochannel error'
106  stop
107  endif
108 
109  ! domain sizes
110  write(*,'(A8,F8.3,A14,F8.3,A14,F8.3)'), &
111  ' Lx =',lx,' Ly =',ly,' Lz =',lz
112 
113  ! grid numbers
114  write(*,'(A8,I8,A14,I8,A14,I8)'), &
115  ' Nx =',nx,' Ny =',ny,' Nz =',nz
116 
117  ! Reynolds number etc
118  if (dp_opt .eq. 1) then ! the constant friction case
119 
120  write(*,'(A8,F8.1,A14,F8.1,A14,F8.2)'), &
121  ' Re_t =',1.d0/nu,' Ri_t =',got,' Pr =',nu/alpha
122 
123  elseif (dp_opt .eq. 2) then ! the constant mass flow rate case
124 
125  write(*,'(A8,F8.1,A14,F8.3,A14,F8.2)'), &
126  ' Re_b =',1.d0/nu,' Ri_b =',got,' Pr =',nu/alpha
127 
128  elseif (dp_opt .eq. 3) then ! the Ekman case
129 
130  write(*,'(A8,F8.5,A14,F8.2,A14,F8.2)'), &
131  ' nu =',nu,' g/t =',got,' Pr =',nu/alpha
132 
133  write(*,'(A8,F8.5,A14,F8.3,A14,F8.3)'), &
134  ' f =',fc,' Ug =',ug,' Vg =',vg
135 
136  endif
137 
138  ! temperature boundary conditions
139  write(*,'(A8,F8.2,A14,F8.2,A14,F8.2)'), &
140  ' Tb =',tbot,' Tt =',ttop,' T_ref =',t_ref
141 
142  ! simulation setup
143  write(*,'(A8,I8,A14,F8.4)'), ' Its =',imax,' dt =',dt
144  write(*,'(A11,I5)'), ' IsScalar: ',isscalar
145  write(*,'(A11,I5)'), ' Isdamp : ',isdamp
146  write(*,'(A11,F5.2)'),' zstretch: ',zstretch
147  write(*,'(A11,F5.2)'), ' u_mrf : ',u_mrf
148  write(*,'(A11,I5)'), ' Nprc : ',nprc
149 
150  ! spatial schemes
151  select case (cds)
152  case (1)
153  write(*,'(A36)'), ' Spatial Derivatives : 2nd Order CDS'
154  case (2)
155  write(*,'(A36)'), ' Spatial Derivatives : 4th Order CDS'
156  case (3)
157  write(*,'(A36)'), ' Spatial Derivatives : Fourier '
158  case default
159  print *, 'cds cheme error'
160  stop
161  end select
162 
163  ! print if the skew-symmetric form is used for nonlinear terms
164  if (cds .eq. 3 .and. issk .eq. 1) then
165  write(*,'(A25)'), ' Skew-Symmetric : On'
166  elseif (cds .eq. 3 .and. issk .ne. 1) then
167  write(*,'(A26)'), ' Skew-Symmetric : Off'
168  endif
169 
170  ! temporal scheme
171  select case (dts)
172  case (1)
173  write(*,'(A33)'), ' Time Scheme : AB2 and CN'
174  case (2)
175  write(*,'(A33)'), ' Time Scheme : RK3 and CN'
176  case default
177  print *, 'dts scheme error'
178  stop
179  end select
180 
181  ! how to drive the flow
182  select case (dp_opt)
183  case (1)
184  print *, 'dP Option : Constant Friction'
185  case (2)
186  print *, 'dP Option : Constant Mass Flow Rate'
187  case (3)
188  print *, 'dP Option : Constant Geostrophic Wind'
189  case default
190  print *, 'dp_opt error'
191  stop
192  end select
193 
194  if (is_ri_var .eq. 1) then
195  print *, 'Cooling Surface : On'
196  write(*,'(A22,F8.4)'), &
197  ' Cooling Rate : ',(ri_end-ri_str)*1.d0/imax/dt
198  elseif (is_ri_var .eq. 0) then
199  print *, 'Cooling Surface : Off'
200  else
201  print *, 'is_ri_var error'
202  stop
203  endif
204 
205  print *, ''
206 
207  ! print the start time
208  call date_and_time(date1,time1,zone1,d1)
209  write(*,'(A21,A8,A2,I2,A1,I2,A1,I2)'), &
210  ' Simulation Started: ',date1, ', ',d1(5),':',d1(6),':',d1(7)
211  print *, ''
212 
213 
214  return
215 
216  end subroutine
217 
218 
219  ! this function prints some information on the screen
220  subroutine goodby_info
222  use mpi
224 
225  implicit none
226 
227  ! these are start date, time, etc
228  character(len=8) :: date1
229  character(len=10) :: time1
230  character(len=5) :: zone1
231  integer,dimension(8) :: d1
232 
233  ! print the end time
234  call date_and_time(date1,time1,zone1,d1)
235  write(*,'(A22,A8,A2,I2,A1,I2,A1,I2)'), ' Simulation finished: ',date1, &
236  ', ',d1(5),':',d1(6),':',d1(7)
237 
238  print *, ''
239  print *, '!***********************************************************!'
240  print *, '! Simulation Finished!! !'
241  print *, '!***********************************************************!'
242 
243  return
244 
245  end subroutine
246 
247 
248 end program
real(mytype), save, protected ri_str
integer, save, protected is_ri_var
integer, save, protected myid
integer, save, protected isxy2d
integer, save, protected ochannel
integer, save, protected nprc
integer, save, protected ny
integer, save, protected isyz2d
real(mytype), save, protected fc
integer, save, protected issk
real(mytype), save, protected dt
integer, save, protected dts
real(mytype), save, protected nu
real(mytype), save, protected ttop
integer, save, protected nz
subroutine goodby_info
Definition: main.f90:221
real(mytype), save, protected ly
real(mytype), save, protected ug
integer, save, protected p_row
integer, save, protected isscalar
real(mytype), save, protected alpha
real(mytype), save, protected vg
program main
Definition: main.f90:22
real(mytype), save, protected t_ref
integer, save, protected isxz2d
subroutine welcome_info
Definition: main.f90:82
real(mytype), save, protected u_mrf
integer, save, protected p_col
integer, save, protected cds
integer, save, protected isdamp
real(mytype), save, protected lz
integer, save, protected imax
real(mytype), save, protected lx
real(mytype), save, protected zstretch
integer, save, protected nx
real(mytype), save, protected ri_end
subroutine initiate_parameters
real(mytype), save, protected got
real(mytype), save, protected tbot
integer, save, protected dp_opt