KTH framework for Nek5000 toolboxes; testing version  0.0.1
spnb.f
Go to the documentation of this file.
1 
6 !=======================================================================
10  subroutine spnb_register()
11  implicit none
12 
13  include 'SIZE'
14  include 'INPUT'
15  include 'FRAMELP'
16  include 'SPNBD'
17 
18  ! local variables
19  integer lpmid
20  real ltim
21 
22  ! functions
23  real dnekclock
24 !-----------------------------------------------------------------------
25  ! timing
26  ltim = dnekclock()
27 
28  ! check if the current module was already registered
29  call mntr_mod_is_name_reg(lpmid,spnb_name)
30  if (lpmid.gt.0) then
31  call mntr_warn(lpmid,
32  $ 'module ['//trim(spnb_name)//'] already registered')
33  return
34  endif
35 
36  ! find parent module
37  call mntr_mod_is_name_reg(lpmid,'FRAME')
38  if (lpmid.le.0) then
39  lpmid = 1
40  call mntr_abort(lpmid,
41  $ 'Parent module ['//'FRAME'//'] not registered')
42  endif
43 
44  ! register module
45  call mntr_mod_reg(spnb_id,lpmid,spnb_name,
46  $ 'Sponge/fringe for rectangular domain')
47 
48  ! register timer
49  call mntr_tmr_is_name_reg(lpmid,'FRM_TOT')
50  call mntr_tmr_reg(spnb_tmr_id,lpmid,spnb_id,
51  $ 'spnb_INI','Sponge calculation initialisation time',.false.)
52 
53  ! register and set active section
54  call rprm_sec_reg(spnb_sec_id,spnb_id,'_'//adjustl(spnb_name),
55  $ 'Runtime paramere section for sponge box module')
56  call rprm_sec_set_act(.true.,spnb_sec_id)
57 
58  ! register parameters
59  call rprm_rp_reg(spnb_str_id,spnb_sec_id,'STRENGTH',
60  $ 'Sponge strength',rpar_real,0,0.0,.false.,' ')
61 
62  call rprm_rp_reg(spnb_wl_id(1),spnb_sec_id,'WIDTHLX',
63  $ 'Sponge left section width; dimension X ',
64  $ rpar_real,0,0.0,.false.,' ')
65 
66  call rprm_rp_reg(spnb_wl_id(2),spnb_sec_id,'WIDTHLY',
67  $ 'Sponge left section width; dimension Y ',
68  $ rpar_real,0,0.0,.false.,' ')
69 
70  if (if3d) call rprm_rp_reg(spnb_wl_id(ndim),spnb_sec_id,
71  $ 'WIDTHLZ','Sponge left section width; dimension Z ',
72  $ rpar_real,0,0.0,.false.,' ')
73 
74  call rprm_rp_reg(spnb_wr_id(1),spnb_sec_id,'WIDTHRX',
75  $ 'Sponge right section width; dimension X ',
76  $ rpar_real,0,0.0,.false.,' ')
77 
78  call rprm_rp_reg(spnb_wr_id(2),spnb_sec_id,'WIDTHRY',
79  $ 'Sponge right section width; dimension Y ',
80  $ rpar_real,0,0.0,.false.,' ')
81 
82  if (if3d) call rprm_rp_reg(spnb_wr_id(ndim),spnb_sec_id,
83  $ 'WIDTHRZ','Sponge right section width; dimension Z ',
84  $ rpar_real,0,0.0,.false.,' ')
85 
86  call rprm_rp_reg(spnb_dl_id(1),spnb_sec_id,'DROPLX',
87  $ 'Sponge left drop/rise section width; dimension X ',
88  $ rpar_real,0,0.0,.false.,' ')
89 
90  call rprm_rp_reg(spnb_dl_id(2),spnb_sec_id,'DROPLY',
91  $ 'Sponge left drop/rise section width; dimension Y ',
92  $ rpar_real,0,0.0,.false.,' ')
93 
94  if (if3d) call rprm_rp_reg(spnb_dl_id(ndim),spnb_sec_id,
95  $ 'DROPLZ','Sponge left drop/rise section width; dimension Z ',
96  $ rpar_real,0,0.0,.false.,' ')
97 
98  call rprm_rp_reg(spnb_dr_id(1),spnb_sec_id,'DROPRX',
99  $ 'Sponge right drop/rise section width; dimension X ',
100  $ rpar_real,0,0.0,.false.,' ')
101 
102  call rprm_rp_reg(spnb_dr_id(2),spnb_sec_id,'DROPRY',
103  $ 'Sponge right drop/rise section width; dimension Y ',
104  $ rpar_real,0,0.0,.false.,' ')
105 
106  if (if3d) call rprm_rp_reg(spnb_dr_id(ndim),spnb_sec_id,
107  $ 'DROPRZ','Sponge right drop/rise section width; dimension Z ',
108  $ rpar_real,0,0.0,.false.,' ')
109 
110  ! timing
111  ltim = dnekclock() - ltim
112  call mntr_tmr_add(spnb_tmr_id,1,ltim)
113 
114  return
115  end subroutine
116 !=======================================================================
122  subroutine spnb_init(lvx,lvy,lvz)
123  implicit none
124 
125  include 'SIZE'
126  include 'INPUT'
127  include 'GEOM'
128  include 'FRAMELP'
129  include 'SPNBD'
130 
131  ! argument list
132  real lvx(LX1*LY1*LZ1*LELV),lvy(LX1*LY1*LZ1*LELV),
133  $ lvz(LX1*LY1*LZ1*LELV)
134 
135  ! local variables
136  integer ierr, nhour, nmin
137  integer itmp
138  real rtmp, ltim
139  logical ltmp
140  character*20 ctmp
141 
142  integer ntot, il, jl
143  real bmin(LDIM), bmax(LDIM)
144 
145  real xxmax, xxmax_c, xxmin, xxmin_c, arg
146  real lcoord(LX1*LY1*LZ1*LELV)
147  common /scruz/ lcoord
148 
149  ! functions
150  real dnekclock, glmin, glmax, math_stepf
151 !-----------------------------------------------------------------------
152  ! check if the module was already initialised
153  if (spnb_ifinit) then
154  call mntr_warn(spnb_id,
155  $ 'module ['//trim(spnb_name)//'] already initiaised.')
156  return
157  endif
158 
159  ! timing
160  ltim = dnekclock()
161 
162  ! get runtime parameters
163  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_str_id,rpar_real)
164  spnb_str = rtmp
165 
166  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_wl_id(1),rpar_real)
167  spnb_wl(1) = rtmp
168 
169  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_wl_id(2),rpar_real)
170  spnb_wl(2) = rtmp
171 
172  if (if3d) then
173  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_wl_id(ndim),
174  $ rpar_real)
175  spnb_wl(ndim) = rtmp
176  endif
177 
178  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_wr_id(1),rpar_real)
179  spnb_wr(1) = rtmp
180 
181  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_wr_id(2),rpar_real)
182  spnb_wr(2) = rtmp
183 
184  if (if3d) then
185  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_wr_id(ndim),
186  $ rpar_real)
187  spnb_wr(ndim) = rtmp
188  endif
189 
190  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_dl_id(1),rpar_real)
191  spnb_dl(1) = rtmp
192 
193  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_dl_id(2),rpar_real)
194  spnb_dl(2) = rtmp
195 
196  if (if3d) then
197  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_dl_id(ndim),
198  $ rpar_real)
199  spnb_dl(ndim) = rtmp
200  endif
201 
202  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_dr_id(1),rpar_real)
203  spnb_dr(1) = rtmp
204 
205  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_dr_id(2),rpar_real)
206  spnb_dr(2) = rtmp
207 
208  if (if3d) then
209  call rprm_rp_get(itmp,rtmp,ltmp,ctmp,spnb_dr_id(ndim),
210  $ rpar_real)
211  spnb_dr(ndim) = rtmp
212  endif
213 
214  ! initialise sponge variables
215 
216  ! get box size
217  ntot = nx1*ny1*nz1*nelv
218  bmin(1) = glmin(xm1,ntot)
219  bmax(1) = glmax(xm1,ntot)
220  bmin(2) = glmin(ym1,ntot)
221  bmax(2) = glmax(ym1,ntot)
222  if(if3d) then
223  bmin(ndim) = glmin(zm1,ntot)
224  bmax(ndim) = glmax(zm1,ntot)
225  endif
226 
227  ! zero spnb_fun
228  call rzero(spnb_fun,ntot)
229 
230 !
231 ! spnb_fun
232 ! *********
233 !
234 ! bmin xxmin_c xxmax_c bmax
235 ! | | xxmin xxmax | |
236 ! |_____| | | |____|
237 ! ^ \ /
238 ! | \ x-> /
239 ! spnb_str| \ /
240 ! v \___________/
241 ! <--------> <-------->
242 ! wl wr
243 ! <--> <-->
244 ! dl dr
245 
246  if(spnb_str.gt.0.0) then
247  call mntr_log(spnb_id,lp_inf,"Sponge turned on")
248 
249  ! save reference field
250  call copy(spnb_vr(1,1),lvx, ntot)
251  call copy(spnb_vr(1,2),lvy, ntot)
252  if (if3d) call copy(spnb_vr(1,ndim),lvz, ntot)
253 
254  ! for every dimension
255  do il=1,ndim
256 
257  if (spnb_wl(il).gt.0.0.or.spnb_wr(il).gt.0.0) then
258  if (spnb_wl(il).lt.spnb_dl(il).or.
259  $ spnb_wr(il).lt.spnb_dr(il)) then
260  call mntr_abort(spnb_id,"Wrong sponge parameters")
261  endif
262 
263  ! sponge beginning (rise at xmax; right)
264  xxmax = bmax(il) - spnb_wr(il)
265  ! end (drop at xmin; left)
266  xxmin = bmin(il) + spnb_wl(il)
267  ! beginnign of constant part (right)
268  xxmax_c = xxmax + spnb_dr(il)
269  ! beginnign of constant part (left)
270  xxmin_c = xxmin - spnb_dl(il)
271 
272  ! get spnb_FUN
273  if (xxmax.le.xxmin) then
274  call mntr_abort(spnb_id,"Sponge too wide")
275  else
276  ! this should be done by pointers, but for now I avoid it
277  if (il.eq.1) then
278  call copy(lcoord,xm1, ntot)
279  elseif (il.eq.2) then
280  call copy(lcoord,ym1, ntot)
281  elseif (il.eq.3) then
282  call copy(lcoord,zm1, ntot)
283  endif
284 
285  do jl=1,ntot
286  rtmp = lcoord(jl)
287  if(rtmp.lt.xxmin_c) then ! constant; xmin
288  rtmp=spnb_str
289  elseif(rtmp.lt.xxmin) then ! fall; xmin
290  arg = (xxmin-rtmp)/spnb_dl(il)
291  rtmp = spnb_str*math_stepf(arg)
292  elseif (rtmp.le.xxmax) then ! zero
293  rtmp = 0.0
294  elseif (rtmp.lt.xxmax_c) then ! rise
295  arg = (rtmp-xxmax)/spnb_dr(il)
296  rtmp = spnb_str*math_stepf(arg)
297  else ! constant
298  rtmp = spnb_str
299  endif
300  spnb_fun(jl)=max(spnb_fun(jl),rtmp)
301  enddo
302 
303  endif ! xxmax.le.xxmin
304 
305  endif ! spnb_w(il).gt.0.0
306  enddo
307 
308 
309  endif
310 
311 #ifdef DEBUG
312  ! for debugging
313  ltmp = ifto
314  ifto = .true.
315  call outpost2(spnb_vr,spnb_vr(1,2),spnb_vr(1,ndim),spnb_fun,
316  $ spnb_fun,1,'spg')
317  ifto = ltmp
318 #endif
319 
320  ! is everything initialised
321  spnb_ifinit=.true.
322 
323  ! timing
324  ltim = dnekclock() - ltim
325  call mntr_tmr_add(spnb_tmr_id,1,ltim)
326 
327  return
328  end subroutine
329 !=======================================================================
333  logical function spnb_is_initialised()
334  implicit none
335 
336  include 'SIZE'
337  include 'SPNBD'
338 !-----------------------------------------------------------------------
339  spnb_is_initialised = spnb_ifinit
340 
341  return
342  end function
343 !=======================================================================
349  subroutine spnb_forcing(ffx,ffy,ffz,ix,iy,iz,ieg)
350  implicit none
351 
352  include 'SIZE' !
353  include 'INPUT' ! IF3D
354  include 'PARALLEL' ! GLLEL
355  include 'SOLN' ! JP
356  include 'SPNBD'
357 
358  ! argument list
359  real ffx, ffy, ffz
360  integer ix,iy,iz,ieg
361 
362  ! local variables
363  integer iel, ip
364 !-----------------------------------------------------------------------
365  iel=gllel(ieg)
366  if (spnb_str.gt.0.0) then
367  ip=ix+nx1*(iy-1+ny1*(iz-1+nz1*(iel-1)))
368 
369  if (jp.eq.0) then
370  ! dns
371  ffx = ffx + spnb_fun(ip)*(spnb_vr(ip,1) - vx(ix,iy,iz,iel))
372  ffy = ffy + spnb_fun(ip)*(spnb_vr(ip,2) - vy(ix,iy,iz,iel))
373  if (if3d) ffz = ffz + spnb_fun(ip)*
374  $ (spnb_vr(ip,ndim) - vz(ix,iy,iz,iel))
375  else
376  ! perturbation
377  ffx = ffx - spnb_fun(ip)*vxp(ip,jp)
378  ffy = ffy - spnb_fun(ip)*vyp(ip,jp)
379  if(if3d) ffz = ffz - spnb_fun(ip)*vzp(ip,jp)
380  endif
381 
382  endif
383 
384  return
385  end subroutine
386 !=======================================================================
gllel
integer function gllel(ieg)
Definition: dprocmap.f:183
mntr_tmr_is_name_reg
subroutine mntr_tmr_is_name_reg(mid, mname)
Check if timer name is registered and return its id.
Definition: mntrtmr.f:146
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
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_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
spnb_init
subroutine spnb_init(lvx, lvy, lvz)
Initilise sponge box module.
Definition: spnb.f:123
spnb_forcing
subroutine spnb_forcing(ffx, ffy, ffz, ix, iy, iz, ieg)
Get sponge forcing.
Definition: spnb.f:350
spnb_register
subroutine spnb_register()
Register sponge box module.
Definition: spnb.f:11
spnb_is_initialised
logical function spnb_is_initialised()
Check if module was initialised.
Definition: spnb.f:334
copy
subroutine copy(a, b, n)
Definition: math.f:260
rzero
subroutine rzero(a, n)
Definition: math.f:208
outpost2
subroutine outpost2(v1, v2, v3, vp, vt, nfldt, name3)
Definition: prepost.f:1394