KTH framework for Nek5000 toolboxes; testing version  0.0.1
pwit.f
Go to the documentation of this file.
1 
6 !=======================================================================
10  subroutine stepper_register()
11  implicit none
12 
13  include 'SIZE'
14  include 'INPUT'
15  include 'FRAMELP'
16  include 'TSTPRD'
17  include 'PWITD'
18 
19  ! local variables
20  integer lpmid, il
21  real ltim
22  character*2 str
23 
24  ! functions
25  real dnekclock
26 !-----------------------------------------------------------------------
27  ! timing
28  ltim = dnekclock()
29 
30  ! check if the current module was already registered
31  call mntr_mod_is_name_reg(lpmid,pwit_name)
32  if (lpmid.gt.0) then
33  call mntr_warn(lpmid,
34  $ 'module ['//trim(pwit_name)//'] already registered')
35  return
36  endif
37 
38  ! find parent module
39  call mntr_mod_is_name_reg(lpmid,tstpr_name)
40  if (lpmid.le.0) then
41  lpmid = 1
42  call mntr_abort(lpmid,
43  $ 'parent module ['//trim(tstpr_name)//'] not registered')
44  endif
45 
46  ! register module
47  call mntr_mod_reg(pwit_id,lpmid,pwit_name,
48  $ 'Power iterations for time stepper')
49 
50  ! register timers
51  ! initialisation
52  call mntr_tmr_reg(pwit_tmr_ini_id,tstpr_tmr_ini_id,pwit_id,
53  $ 'PWIT_INI','Power iteration initialisation time',.true.)
54  ! submodule operation
55  call mntr_tmr_reg(pwit_tmr_evl_id,tstpr_tmr_evl_id,pwit_id,
56  $ 'PWIT_EVL','Power iteration evolution time',.true.)
57 
58  ! register and set active section
59  call rprm_sec_reg(pwit_sec_id,pwit_id,'_'//adjustl(pwit_name),
60  $ 'Runtime paramere section for power iteration module')
61  call rprm_sec_set_act(.true.,pwit_sec_id)
62 
63  ! register parameters
64  call rprm_rp_reg(pwit_l2n_id,pwit_sec_id,'L2N',
65  $ 'Vector initial norm',rpar_real,0,1.0,.false.,' ')
66 
67  ! set initialisation flag
68  pwit_ifinit=.false.
69 
70  ! timing
71  ltim = dnekclock() - ltim
72  call mntr_tmr_add(pwit_tmr_ini_id,1,ltim)
73 
74  return
75  end subroutine
76 !=======================================================================
80  subroutine stepper_init()
81  implicit none
82 
83  include 'SIZE'
84  include 'SOLN' ! V[XYZ]P, TP
85  include 'MASS' ! BM1
86  include 'FRAMELP'
87  include 'TSTPRD'
88  include 'PWITD'
89 
90  ! local variables
91  integer itmp, il, set_in
92  real rtmp, ltim, lnorm
93  logical ltmp
94  character*20 ctmp
95 
96  ! to get checkpoint runtime parameters
97  integer ierr, lmid, lsid, lrpid
98 
99  ! functions
100  real dnekclock, cnht_glsc2_wt
101  logical chkpts_is_initialised
102 !-----------------------------------------------------------------------
103  ! check if the module was already initialised
104  if (pwit_ifinit) then
105  call mntr_warn(pwit_id,
106  $ 'module ['//trim(pwit_name)//'] already initialised.')
107  return
108  endif
109 
110  ! timing
111  ltim = dnekclock()
112 
113  ! get runtime parameters
114  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,pwit_l2n_id,rpar_real)
115  pwit_l2n = rtmp
116 
117  ! check the restart flag
118  ! check if checkpointing module was registered and take parameters
119  ierr = 0
120  call mntr_mod_is_name_reg(lmid,'CHKPT')
121  if (lmid.gt.0) then
122  call rprm_sec_is_name_reg(lsid,lmid,'_CHKPT')
123  if (lsid.gt.0) then
124  ! restart flag
125  call rprm_rp_is_name_reg(lrpid,lsid,'READCHKPT',rpar_log)
126  if (lrpid.gt.0) then
127  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,lrpid,rpar_log)
128  pwit_ifrst = ltmp
129  else
130  ierr = 1
131  goto 30
132  endif
133  if (pwit_ifrst) then
134  ! checkpoint set number
135  call rprm_rp_is_name_reg(lrpid,lsid,'CHKPFNUMBER',
136  $ rpar_int)
137  if (lrpid.gt.0) then
138  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,lrpid,rpar_int)
139  pwit_fnum = itmp
140  else
141  ierr = 1
142  goto 30
143  endif
144  endif
145  else
146  ierr = 1
147  endif
148  else
149  ierr = 1
150  endif
151 
152  30 continue
153 
154  ! check for errors
155  call mntr_check_abort(pwit_id,ierr,
156  $ 'Error reading checkpoint parameters')
157 
158  ! read checkpoint file
159  if (pwit_ifrst) then
160  if(.not.chkpts_is_initialised()) call mntr_abort(pwit_id,
161  $ 'Checkpointing module not initialised')
162  set_in = pwit_fnum -1
163  call stepper_read(set_in)
164  endif
165 
166  ! initial growth rate
167  pwit_grw = 0.0
168 
169  ! normalise vector
170  lnorm = cnht_glsc2_wt(vxp,vyp,vzp,tp,vxp,vyp,vzp,tp,bm1)
171  lnorm = sqrt(pwit_l2n/lnorm)
172  call cnht_opcmult (vxp,vyp,vzp,tp,lnorm)
173 
174  if (tstpr_pr.ne.0) then
175  ! normalise pressure ???????????????????????????
176  itmp = nx2*ny2*nz2*nelv
177  call cmult(prp,lnorm,itmp)
178  endif
179 
180  ! make sure the velocity and temperature fields are continuous at
181  ! element faces and edges
182  call tstpr_dssum
183 
184  ! save intial vector
185  call cnht_opcopy (pwit_vx,pwit_vy,pwit_vz,pwit_t,vxp,vyp,vzp,tp)
186 
187  ! stamp log file
188  call mntr_log(pwit_id,lp_prd,'POWER ITERATIONS initialised')
189  call mntr_logr(pwit_id,lp_prd,'L2NORM = ',pwit_l2n)
190 
191  ! everything is initialised
192  pwit_ifinit=.true.
193 
194  ! timing
195  ltim = dnekclock() - ltim
196  call mntr_tmr_add(pwit_tmr_ini_id,1,ltim)
197 
198  return
199  end subroutine
200 !=======================================================================
204  logical function stepper_is_initialised()
205  implicit none
206 
207  include 'SIZE'
208  include 'PWITD'
209 !-----------------------------------------------------------------------
210  stepper_is_initialised = pwit_ifinit
211 
212  return
213  end function
214 !=======================================================================
219  subroutine stepper_vsolve
220  implicit none
221 
222  include 'SIZE' ! NIO
223  include 'TSTEP' ! TIME, LASTEP, NSTEPS
224  include 'INPUT' ! IFHEAT
225  include 'MASS' ! BM1
226  include 'SOLN' ! V[XYZ]P, TP
227  include 'FRAMELP'
228  include 'TSTPRD'
229  include 'PWITD'
230 
231  ! scratch space
232  real TA1 (LPX1*LPY1*LPZ1*LELV), TA2 (LPX1*LPY1*LPZ1*LELV),
233  $ TA3 (LPX1*LPY1*LPZ1*LELV), TAT (LPX1*LPY1*LPZ1*LELT)
234  COMMON /scruz/ ta1, ta2, ta3, tat
235 
236  ! local variables
237  integer itmp
238  real lnorm, grth_old, ltim
239 
240  ! functions
241  real dnekclock, cnht_glsc2_wt
242 !-----------------------------------------------------------------------
243  ! timing
244  ltim=dnekclock()
245 
246  ! normalise vector
247  lnorm = cnht_glsc2_wt(vxp,vyp,vzp,tp,vxp,vyp,vzp,tp,bm1)
248  lnorm = sqrt(pwit_l2n/lnorm)
249  call cnht_opcmult (vxp,vyp,vzp,tp,lnorm)
250 
251  if (tstpr_pr.ne.0) then
252  ! normalise pressure ???????????????????????????
253  itmp = nx2*ny2*nz2*nelv
254  call cmult(prp,lnorm,itmp)
255  endif
256 
257  ! make sure the velocity and temperature fields are continuous at
258  ! element faces and edges
259  call tstpr_dssum
260 
261  ! compare current and prevoius growth rate
262  grth_old = pwit_grw
263  pwit_grw = 1.0/lnorm
264  grth_old = pwit_grw - grth_old
265 
266  ! get L2 norm of the update
267  call cnht_opsub3 (ta1,ta2,ta3,tat,pwit_vx,pwit_vy,pwit_vz,pwit_t,
268  $ vxp,vyp,vzp,tp)
269  lnorm = cnht_glsc2_wt(ta1,ta2,ta3,tat,ta1,ta2,ta3,tat,bm1)
270  lnorm = sqrt(lnorm)
271 
272  ! log stamp
273  call mntr_log(pwit_id,lp_prd,'POWER ITERATIONS: convergence')
274  call mntr_logr(pwit_id,lp_prd,'||V-V_old|| = ',lnorm)
275  call mntr_logr(pwit_id,lp_prd,'Growth ',pwit_grw)
276 
277  itmp = 0
278  if (ifheat) itmp = 1
279 
280  !write down current field
281  call outpost2(vxp,vyp,vzp,prp,tp,itmp,'PWI')
282 
283  ! write down field difference
284  call outpost2(ta1,ta2,ta3,prp,tat,itmp,'VDF')
285 
286  ! check convergence
287  if(lnorm.lt.tstpr_tol.and.grth_old.lt.tstpr_tol) then
288  call mntr_log(pwit_id,lp_prd,'Reached stopping criteria')
289  ! mark the last step
290  lastep = 1
291  else
292  ! save current vector and restart stepper
293  call cnht_opcopy(pwit_vx,pwit_vy,pwit_vz,pwit_t,vxp,vyp,vzp,tp)
294  endif
295 
296  ! save checkpoint
297  if (lastep.eq.1.or.tstpr_cmax.eq.tstpr_vstep) then
298  call stepper_write()
299 
300  ! mark the last step
301  lastep = 1
302  endif
303 
304  ! timing
305  ltim = dnekclock() - ltim
306  call mntr_tmr_add(pwit_tmr_evl_id,1,ltim)
307 
308  if (lastep.eq.1) then
309  ! final log stamp
310  call mntr_log(pwit_id,lp_prd,'POWER ITERATIONS finalised')
311  call mntr_logr(pwit_id,lp_prd,'||V-V_old|| = ',lnorm)
312  call mntr_logr(pwit_id,lp_prd,'Growth ',pwit_grw)
313  endif
314 
315  return
316  end subroutine
317 !=======================================================================
321  subroutine stepper_read(set_in)
322  implicit none
323 
324  include 'SIZE' ! NIO
325  include 'TSTEP' ! TIME, LASTEP, NSTEPS
326  include 'INPUT' ! IFMVBD, IFREGUO
327  include 'FRAMELP'
328  include 'CHKPTD'
329  include 'CHKPTMSD'
330  include 'PWITD'
331 
332  ! argument list
333  integer set_in
334 
335  ! local variables
336  integer ifile, step_cnt, fnum
337  character*132 fname(chkptms_fmax)
338  logical ifreguol
339 !-----------------------------------------------------------------------
340  ! no regular mesh
341  ifreguol= ifreguo
342  ifreguo = .false.
343 
344  call mntr_log(pwit_id,lp_inf,'Reading checkpoint snapshot')
345 
346  ! initialise I/O data
347  call io_init
348 
349  ! get set of file names in the snapshot
350  ifile = 1
351  call chkptms_set_name(fname, fnum, set_in, ifile)
352 
353  ! read files
354  call chkptms_restart_read(fname, fnum)
355 
356  ! put parameters back
357  ifreguo = ifreguol
358 
359  return
360  end subroutine
361 !=======================================================================
364  subroutine stepper_write
365  implicit none
366 
367  include 'SIZE' ! NIO
368  include 'TSTEP' ! TIME, LASTEP, NSTEPS
369  include 'INPUT' ! IFMVBD, IFREGUO
370  include 'FRAMELP'
371  include 'CHKPTD'
372  include 'CHKPTMSD'
373  include 'PWITD'
374 
375  ! local variables
376  integer ifile, step_cnt, set_out, fnum
377  character*132 fname(chkptms_fmax)
378  logical ifcoord
379  logical ifreguol
380 !-----------------------------------------------------------------------
381  ! no regular mesh
382  ifreguol= ifreguo
383  ifreguo = .false.
384 
385  call mntr_log(pwit_id,lp_inf,'Writing checkpoint snapshot')
386 
387  ! initialise I/O data
388  call io_init
389 
390  ! get set of file names in the snapshot
391  ifile = 1
392  call chkpt_get_fset(step_cnt, set_out)
393  call chkptms_set_name(fname, fnum, set_out, ifile)
394 
395  ifcoord = .true.
396  ! write down files
397  call chkptms_restart_write(fname, fnum, ifcoord)
398 
399  ! put parameters back
400  ifreguo = ifreguol
401 
402  return
403  end subroutine
404 !=======================================================================
405 
stepper_register
subroutine stepper_register()
Register Arnoldi ARPACK module.
Definition: arna.f:24
stepper_vsolve
subroutine stepper_vsolve
Create Krylov space, get Ritz values and restart stepper phase.
Definition: arna.f:354
stepper_is_initialised
logical function stepper_is_initialised()
Check if module was initialised.
Definition: arna.f:338
stepper_init
subroutine stepper_init()
Initilise Arnoldi ARPACK module.
Definition: arna.f:97
chkpt_get_fset
subroutine chkpt_get_fset(step_cnt, set_out)
Get step count to the checkpoint and a set number.
Definition: chkpt.f:256
chkptms_restart_read
subroutine chkptms_restart_read(fname, fnum)
Read checkpoint snapshot.
Definition: chkptms.f:604
chkptms_restart_write
subroutine chkptms_restart_write(fname, fnum, ifcoord)
Write checkpoint snapshot.
Definition: chkptms.f:512
chkptms_set_name
subroutine chkptms_set_name(fname, fnum, nset, ifile)
Generate set of restart file names in snapshot.
Definition: chkptms.f:384
cnht_opcopy
subroutine cnht_opcopy(a1, a2, a3, a4, b1, b2, b3, b4)
Copy vectors A=B (velocity and temperature)
Definition: cnht_tools.f:279
cnht_opsub3
subroutine cnht_opsub3(a1, a2, a3, a4, b1, b2, b3, b4, c1, c2, c3, c4)
Subtract vectors A = B-C (velocity and temperature)
Definition: cnht_tools.f:374
cnht_opcmult
subroutine cnht_opcmult(a1, a2, a3, a4, const)
Multiply vector by constant A = c*A (single coeff. for velocity and temperature)
Definition: cnht_tools.f:405
mntr_logr
subroutine mntr_logr(mid, priority, logs, rvar)
Write log message adding single real.
Definition: mntrlog.f:731
mntr_warn
subroutine mntr_warn(mid, logs)
Write warning message.
Definition: mntrlog.f:803
mntr_tmr_add
subroutine mntr_tmr_add(mid, icount, time)
Check if timer id is registered. This operation is performed locally.
Definition: mntrtmr.f:237
mntr_mod_is_name_reg
subroutine mntr_mod_is_name_reg(mid, mname)
Check if module name is registered and return its id.
Definition: mntrlog.f:459
mntr_abort
subroutine mntr_abort(mid, logs)
Abort simulation.
Definition: mntrlog.f:837
mntr_log
subroutine mntr_log(mid, priority, logs)
Write log message.
Definition: mntrlog.f:600
mntr_mod_reg
subroutine mntr_mod_reg(mid, pmid, mname, mdscr)
Register new module.
Definition: mntrlog.f:346
mntr_tmr_reg
subroutine mntr_tmr_reg(mid, pmid, modid, mname, mdscr, ifsum)
Register new timer.
Definition: mntrtmr.f:16
mntr_check_abort
subroutine mntr_check_abort(mid, ierr, logs)
Abort simulation.
Definition: mntrlog.f:856
stepper_read
subroutine stepper_read(set_in)
Read restart files.
Definition: pwit.f:322
stepper_write
subroutine stepper_write
Write restart files.
Definition: pwit.f:365
rprm_rp_is_name_reg
subroutine rprm_rp_is_name_reg(rpid, mid, pname, ptype)
Check if parameter name is registered and return its id. Check flags as well.
Definition: rprm.f:658
rprm_rp_get
subroutine rprm_rp_get(ipval, rpval, lpval, cpval, rpid, ptype)
Get runtime parameter form active section. This operation is performed locally.
Definition: rprm.f:883
rprm_sec_is_name_reg
subroutine rprm_sec_is_name_reg(rpid, mid, pname)
Check if section name is registered and return its id. Check mid as well.
Definition: rprm.f:284
rprm_rp_reg
subroutine rprm_rp_reg(rpid, mid, pname, pdscr, ptype, ipval, rpval, lpval, cpval)
Register new runtime parameter.
Definition: rprm.f:484
rprm_sec_set_act
subroutine rprm_sec_set_act(ifact, rpid)
Set section's activation flag. Master value is broadcasted.
Definition: rprm.f:422
rprm_sec_reg
subroutine rprm_sec_reg(rpid, mid, pname, pdscr)
Register new parameter section.
Definition: rprm.f:165
tstpr_dssum
subroutine tstpr_dssum
Average velocity and temperature at element faces.
Definition: tstpr.f:387
cmult
subroutine cmult(a, const, n)
Definition: math.f:315
outpost2
subroutine outpost2(v1, v2, v3, vp, vt, nfldt, name3)
Definition: prepost.f:1394
io_init
subroutine io_init
Definition: prepost.f:1024