MPI-AMRVAC 3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
Loading...
Searching...
No Matches
mod_bc_data.t
Go to the documentation of this file.
1!> Module to set boundary conditions from user data
4 use mod_global_parameters, only: std_len
5
6 implicit none
7 private
8
9 integer, parameter :: max_boundary_conds = 10
10
11 type bc_data_t
12 integer :: n_variables
13 double precision :: origin(3)
14 double precision :: dx(3)
15 integer :: n_points(3)
16 character(len=40), allocatable :: names(:)
17 double precision, allocatable :: values(:, :, :, :)
18 end type bc_data_t
19
20 type(lt_t) :: lt_1d(max_boundary_conds)
21 type(lt2_t) :: lt_2d(max_boundary_conds)
22 type(lt3_t) :: lt_3d(max_boundary_conds)
23
24 !> Whether boundary condition data is time varying
25 logical, public, protected :: bc_data_time_varying = .false.
26
27
28 !> Integer array for indexing lookup tables per variable per direction
29 integer, public, protected, allocatable :: bc_data_ix(:, :)
30
31 !> data file name
32 character(len=std_len), public, protected :: boundary_data_file_name
33 logical, public, protected :: interp_phy_first_row=.true.
34 logical, public, protected :: interp_phy_first_row_same=.false.
35 logical, public, protected :: bc_phy_first_row=.false.
36 logical, public, protected :: boundary_data_primitive=.false.
37
38 public :: bc_data_init
39 public :: bc_data_set
40 public :: bc_data_get_2d
41 public :: bc_data_get_3d
42 public :: lt_3d
43contains
44
45 subroutine bc_data_init()
47 use mod_comm_lib, only: mpistop
48
49 integer :: i, iw, ib, n_files, n_bc
50 double precision :: xmax(3)
51 type(bc_data_t) :: bc
52
53 call bc_read_params(par_files)
54
55 allocate(bc_data_ix(nwfluxbc, 2*ndim))
56
57 bc_data_ix(:, :) = -1
58 n_bc = 0
59
60 if(any(typeboundary(1:nwfluxbc,1:2*ndim)==bc_data ) .or. any(typeboundary(1:nwfluxbc,1:2*ndim)==bc_icarus)) then
61 call read_vtk_structured_points(trim(boundary_data_file_name), bc)
62 xmax = bc%origin + (bc%n_points-1) * bc%dx
63 bc_data_time_varying = (bc%n_points(ndim) > 1)
64
65 endif
66
67 do ib = 1, 2 * ndim
68 do iw = 1, nwfluxbc
69 if (typeboundary(iw, ib)==bc_data .or. typeboundary(iw, ib)==bc_icarus) then
70 n_bc = n_bc + 1
71 bc_data_ix(iw, ib) = n_bc
72
73
74 {^ifoned
75 call mpistop("bc_data_init: 1D case not supported")
76 }
77 {^iftwod
78 if (bc_data_time_varying) then
79 lt_2d(n_bc) = lt2_create_from_data(bc%origin(1:ndim), &
80 xmax(1:ndim), bc%values(:, :, 1:1, n_bc))
81 else
82 ! Use first point in time
83 lt_1d(n_bc) = lt_create_from_data(bc%origin(1), &
84 xmax(1), bc%values(:, 1, 1:1, n_bc))
85 end if
86 }
87 {^ifthreed
88 if (bc_data_time_varying) then
89 lt_3d(n_bc) = lt3_create_from_data(bc%origin(1:ndim), &
90 xmax(1:ndim), bc%values(:, :, :, n_bc:n_bc))
91 else
92 ! Use first point in time
93 lt_2d(n_bc) = lt2_create_from_data(bc%origin(1:ndim-1), &
94 xmax(1:ndim-1), bc%values(:, :, 1:1, n_bc))
95 end if
96 }
97 end if
98 end do
99 end do
100
101 end subroutine bc_data_init
102
103 !> Read this module"s parameters from a file
104 subroutine bc_read_params(files)
106 character(len=*), intent(in) :: files(:)
107 integer :: n
108
111
112 do n = 1, size(files)
113 open(unitpar, file=trim(files(n)), status="old")
114 read(unitpar, bd_list, end=111)
115111 close(unitpar)
116 end do
117
118 end subroutine bc_read_params
119
120 elemental function bc_data_get_3d(n_bc, x1, x2, qt) result(val)
121 integer, intent(in) :: n_bc
122 double precision, intent(in) :: x1, x2, qt
123 double precision :: val
124
125 if (bc_data_time_varying) then
126 val = lt3_get_col(lt_3d(n_bc), 1, x1, x2, qt)
127 else
128 val = lt2_get_col(lt_2d(n_bc), 1, x1, x2)
129 end if
130 end function bc_data_get_3d
131
132 elemental function bc_data_get_2d(n_bc, x1, qt) result(val)
133 integer, intent(in) :: n_bc
134 double precision, intent(in) :: x1, qt
135 double precision :: val
136
137 if (bc_data_time_varying) then
138 val = lt2_get_col(lt_2d(n_bc), 1, x1, qt)
139 else
140 val = lt_get_col(lt_1d(n_bc), 1, x1)
141 end if
142 end function bc_data_get_2d
143
144 !> Set boundary conditions according to user data
145 subroutine bc_data_set(qt,ixI^L,ixO^L,iB,w,x)
148 use mod_comm_lib, only: mpistop
149
150 integer, intent(in) :: ixi^l, ixo^l, ib
151 double precision, intent(in) :: qt, x(ixi^s,1:ndim)
152 double precision, intent(inout) :: w(ixi^s,1:nw)
153 double precision :: tmp(ixo^s)
154 integer :: i, ix, iw, n_bc
155
156 integer :: ixoo^l
157
158 ixoo^l=ixo^l;
159
160 {^iftwod
161 select case (ib)
162 case (1, 2)
163 if(interp_phy_first_row) then
164 if (ib == 1) then
165 ix = ixomax1+1
166 else
167 ix = ixomin1-1
168 end if
170 call phys_to_primitive(ixi^l,ix,ixomin2,ix,ixomax2,w,x)
171 endif
172 endif
173 ! the reason for this is that not all bc for all the vars
174 ! might be set with bc_data
175 if(bc_phy_first_row)then
176 if (ib == 1) then
177 ixoomax1=ixoomax1+1
178 else
179 ixoomin1=ixoomin1-1
180 endif
182 if (ib == 1) then
183 ixoomin1=ixoomax1
184 else
185 ixoomax1=ixoomin1
186 endif
187 call phys_to_primitive(ixi^l,ixoo^l,w,x)
188 if (ib == 1) then
189 ixoomin1=ixomin1 ! to be used later in to_conserved
190 else
191 ixoomax1=ixomax1
192 endif
193 endif
194 endif
195 do iw = 1, nwfluxbc
196 n_bc = bc_data_ix(iw, ib)
197 if (n_bc == -1) cycle
198
199 tmp(ixomin1, ixomin2:ixomax2) = bc_data_get_2d(n_bc, &
200 x(ixomin1, ixomin2:ixomax2, 2), qt)
201 if(interp_phy_first_row) then
202 ! Approximate boundary value by linear interpolation to first ghost
203 ! cell, rest of ghost cells contains the same value
205 do i = 0, ixoomax1-ixoomin1
206 w(ixoomin1+i, ixomin2:ixomax2, iw) = &
207 2 * tmp(ixomin1, ixomin2:ixomax2) - &
208 w(ix, ixomin2:ixomax2, iw)
209 end do
210 else
211 if (ib==1) then
212 do i = 0, ixoomax1-ixoomin1
213 if (i==0) then
214 w(ixoomin1+i, ixomin2:ixomax2, iw) = &
215 4 * tmp(ixomin1, ixomin2:ixomax2) - &
216 3 * w(ix, ixomin2:ixomax2, iw)
217 elseif (i==1)then
218 w(ixoomin1+i, ixomin2:ixomax2, iw) = &
219 2 * tmp(ixomin1, ixomin2:ixomax2) - &
220 w(ix, ixomin2:ixomax2, iw)
221 else
222 w(ixoomin1+i, ixomin2:ixomax2, iw) = &
223 tmp(ixomin1, ixomin2:ixomax2)
224 endif
225
226 end do
227 else
228 do i = 0, ixoomax1-ixoomin1
229 if (i==ixoomax1-ixoomin1) then
230 w(ixoomin1+i, ixomin2:ixomax2, iw) = &
231 4 * tmp(ixomin1, ixomin2:ixomax2) - &
232 3 * w(ix, ixomin2:ixomax2, iw)
233 elseif (i==ixoomax1-ixoomin1-1)then
234 w(ixoomin1+i, ixomin2:ixomax2, iw) = &
235 2 * tmp(ixomin1, ixomin2:ixomax2) - &
236 w(ix, ixomin2:ixomax2, iw)
237 else
238 w(ixoomin1+i, ixomin2:ixomax2, iw) = &
239 tmp(ixomin1, ixomin2:ixomax2)
240 endif
241
242 end do
243 endif
244
245 endif
246 else
247 do i = ixoomin1, ixoomax1
248 w(i, ixomin2:ixomax2, iw) = &
249 tmp(ixomin1, ixomin2:ixomax2)
250 end do
251
252 endif
253 end do
254 if(interp_phy_first_row) then
256 call phys_to_conserved(ixi^l,ix,ixomin2,ix,ixomax2,w,x)
257 endif
258 endif
259 case (3, 4)
260 if(interp_phy_first_row) then
261 if (ib == 3) then
262 ix = ixomax2+1
263 else
264 ix = ixomin2-1
265 end if
267 call phys_to_primitive(ixi^l,ixomin1,ix,ixomax1,ix,w,x)
268 endif
269 endif
270 if(bc_phy_first_row)then
271 if (ib == 3) then
272 ixoomax2=ixoomax2+1
273 else
274 ixoomin2=ixoomin2-1
275 endif
277 if (ib == 3) then
278 ixoomin2=ixoomax2
279 else
280 ixoomax2=ixoomin2
281 endif
282 call phys_to_primitive(ixi^l,ixoo^l,w,x)
283 if (ib == 3) then
284 ixoomin2=ixomin2 ! to be used later in to_conserved
285 else
286 ixoomax2=ixomax2
287 endif
288 endif
289 endif
290 do iw = 1, nwfluxbc
291 n_bc = bc_data_ix(iw, ib)
292 if (n_bc == -1) cycle
293
294 tmp(ixomin1:ixomax1, ixomin2) = bc_data_get_2d(n_bc, &
295 x(ixomin1:ixomax1, ixomin2, 1), qt)
296
297 if(interp_phy_first_row) then
298
299 ! Approximate boundary value by linear interpolation to first ghost
300 ! cell, rest of ghost cells contains the same value
302 do i = 0, ixoomax2-ixoomin2
303 w(ixomin1:ixomax1, ixoomin2+i, iw) = &
304 2 * tmp(ixomin1:ixomax1, ixomin2) - &
305 w(ixomin1:ixomax1, ix, iw)
306 end do
307 else
308 if (ib==3) then
309 do i = 0, ixoomax2-ixoomin2
310 if (i==0) then
311 w(ixomin1:ixomax1, ixoomin2+i, iw) = &
312 4 * tmp(ixomin1:ixomax1, ixomin2) - &
313 3 * w(ixomin1:ixomax1, ix, iw)
314 elseif (i==1)then
315 w(ixomin1:ixomax1, ixoomin2+i, iw) = &
316 2 * tmp(ixomin1:ixomax1, ixomin2) - &
317 w(ixomin1:ixomax1, ix, iw)
318 else
319 w(ixomin1:ixomax1, ixoomin2+i, iw) = &
320 tmp(ixomin1:ixomax1, ixomin2)
321 endif
322
323 end do
324 else
325 do i = 0, ixoomax2-ixoomin2
326 if (i==ixoomax2-ixoomin2) then
327 w(ixomin1:ixomax1, ixoomin2+i, iw) = &
328 4 * tmp(ixomin1:ixomax1, ixomin2) - &
329 3 * w(ixomin1:ixomax1, ix, iw)
330 elseif (i==ixoomax2-ixoomin2-1)then
331 w(ixomin1:ixomax1, ixoomin2+i, iw) = &
332 2 * tmp(ixomin1:ixomax1, ixomin2) - &
333 w(ixomin1:ixomax1, ix, iw)
334 else
335 w(ixomin1:ixomax1, ixoomin2+i, iw) = &
336 tmp(ixomin1:ixomax1, ixomin2)
337 endif
338
339 end do
340 endif
341
342 endif
343 else
344 do i = ixoomin2, ixoomax2
345 w(ixomin1:ixomax1, i, iw) = &
346 tmp(ixomin1:ixomax1, ixomin2)
347 end do
348
349 endif
350 end do
351 if(interp_phy_first_row) then
353 call phys_to_conserved(ixi^l,ixomin1,ix,ixomax1,ix,w,x)
354 endif
355 endif
356 case default
357 call mpistop("bc_data_set: unknown iB")
358 end select
359 }
360
361 {^ifthreed
362 select case (ib)
363 case (1, 2)
364 if(interp_phy_first_row) then
365 if (ib == 1) then
366 ix = ixomax1+1
367 else
368 ix = ixomin1-1
369 end if
371 call phys_to_primitive(ixi^l,ix,ixomin2,ixomin3,ix,ixomax2,ixomax3,w,x)
372 endif
373 endif
374 if(bc_phy_first_row)then
375 if (ib == 1) then
376 ixoomax1=ixoomax1+1
377 else
378 ixoomin1=ixoomin1-1
379 endif
381 if (ib == 1) then
382 ixoomin1=ixoomax1
383 else
384 ixoomax1=ixoomin1
385 endif
386 call phys_to_primitive(ixi^l,ixoo^l,w,x)
387 if (ib == 1) then
388 ixoomin1=ixomin1 ! to be used later in to_conserved
389 else
390 ixoomax1=ixomax1
391 endif
392 endif
393 endif
394 do iw = 1, nwfluxbc
395 n_bc = bc_data_ix(iw, ib)
396 if (n_bc == -1) cycle
397
398 tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3) = bc_data_get_3d(n_bc, &
399 x(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3, 2), &
400 x(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3, 3), qt)
401
402 if(interp_phy_first_row) then
403
404 ! Approximate boundary value by linear interpolation to first ghost
405 ! cell, rest of ghost cells contains the same value
407 do i = 0, ixoomax1-ixoomin1
408 w(ixoomin1+i, ixomin2:ixomax2, ixomin3:ixomax3, iw) = &
409 2 * tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3) - &
410 w(ix, ixomin2:ixomax2, ixomin3:ixomax3, iw)
411 end do
412 else
413 if (ib==1) then
414 do i = 0, ixoomax1-ixoomin1
415 if (i==0) then
416 w(ixoomin1+i, ixomin2:ixomax2, ixomin3:ixomax3, iw) = &
417 4 * tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3) - &
418 3*w(ix, ixomin2:ixomax2, ixomin3:ixomax3, iw)
419 elseif (i==1)then
420 w(ixoomin1+i, ixomin2:ixomax2, ixomin3:ixomax3, iw) = &
421 2 * tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3) - &
422 w(ix, ixomin2:ixomax2, ixomin3:ixomax3, iw)
423 else
424 w(ixoomin1+i, ixomin2:ixomax2, ixomin3:ixomax3, iw) = &
425 tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3)
426 endif
427
428 end do
429 else
430 do i = 0, ixoomax1-ixoomin1
431 if (i==ixoomax1-ixoomin1) then
432 w(ixoomin1+i, ixomin2:ixomax2, ixomin3:ixomax3, iw) = &
433 4 * tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3) - &
434 3*w(ix, ixomin2:ixomax2, ixomin3:ixomax3, iw)
435 elseif (i==ixoomax1-ixoomin1-1)then
436 w(ixoomin1+i, ixomin2:ixomax2, ixomin3:ixomax3, iw) = &
437 2 * tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3) - &
438 w(ix, ixomin2:ixomax2, ixomin3:ixomax3, iw)
439 else
440 w(ixoomin1+i, ixomin2:ixomax2, ixomin3:ixomax3, iw) = &
441 tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3)
442 endif
443
444 end do
445 endif
446
447 endif
448 else
449 do i = ixoomin1,ixoomax1
450 w(i, ixomin2:ixomax2, ixomin3:ixomax3, iw) = &
451 tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3)
452 end do
453 end if
454
455 end do
456 if(interp_phy_first_row) then
458 call phys_to_conserved(ixi^l,ix,ixomin2,ixomin3,ix,ixomax2,ixomax3,w,x)
459 endif
460 endif
461 case (3, 4)
462 if(interp_phy_first_row) then
463 if (ib == 3) then
464 ix = ixomax2+1
465 else
466 ix = ixomin2-1
467 end if
469 call phys_to_primitive(ixi^l,ixomin1,ix,ixomin3,ixomax1,ix,ixomax3,w,x)
470 endif
471 endif
472 if(bc_phy_first_row)then
473 if (ib == 3) then
474 ixoomax2=ixoomax2+1
475 else
476 ixoomin2=ixoomin2-1
477 endif
479 if (ib == 3) then
480 ixoomin2=ixoomax2
481 else
482 ixoomax2=ixoomin2
483 endif
484 call phys_to_primitive(ixi^l,ixoo^l,w,x)
485 if (ib == 3) then
486 ixoomin2=ixomin2 ! to be used later in to_conserved
487 else
488 ixoomax2=ixomax2
489 endif
490 endif
491 endif
492 do iw = 1, nwfluxbc
493 n_bc = bc_data_ix(iw, ib)
494 if (n_bc == -1) cycle
495
496 tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3) = bc_data_get_3d(n_bc, &
497 x(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3, 1), &
498 x(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3, 3), qt)
499
500 if(interp_phy_first_row) then
501
503 ! Approximate boundary value by linear interpolation to first ghost
504 ! cell, rest of ghost cells contains the same value
505 do i = 0, ixoomax2-ixoomin2
506 w(ixomin1:ixomax1, ixoomin2+i, ixomin3:ixomax3, iw) = &
507 2 * tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3) - &
508 w(ixomin1:ixomax1, ix, ixomin3:ixomax3, iw)
509 end do
510 else
511 if (ib==3) then
512 do i = 0, ixoomax2-ixoomin2
513 if (i==0) then
514 w(ixomin1:ixomax1, ixoomin2+i, ixomin3:ixomax3, iw) = &
515 4 * tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3) - &
516 3* w(ixomin1:ixomax1, ix, ixomin3:ixomax3, iw)
517 elseif (i==1)then
518 w(ixomin1:ixomax1, ixoomin2+i, ixomin3:ixomax3, iw) = &
519 2 * tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3) - &
520 w(ixomin1:ixomax1, ix, ixomin3:ixomax3, iw)
521 else
522 w(ixomin1:ixomax1, ixoomin2+i, ixomin3:ixomax3, iw) = &
523 tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3)
524 endif
525
526 end do
527 else
528 do i = 0, ixoomax2-ixoomin2
529 if (i==ixoomax2-ixoomin2) then
530 w(ixomin1:ixomax1, ixoomin2+i, ixomin3:ixomax3, iw) = &
531 4 * tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3) - &
532 3* w(ixomin1:ixomax1, ix, ixomin3:ixomax3, iw)
533 elseif (i==ixoomax2-ixoomin2-1)then
534 w(ixomin1:ixomax1, ixoomin2+i, ixomin3:ixomax3, iw) = &
535 2 * tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3) - &
536 w(ixomin1:ixomax1, ix, ixomin3:ixomax3, iw)
537 else
538 w(ixomin1:ixomax1, ixoomin2+i, ixomin3:ixomax3, iw) = &
539 tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3)
540 endif
541
542 end do
543 endif
544
545 endif
546 else
547 do i = ixoomin2,ixoomax2
548 w(ixomin1:ixomax1, i, ixomin3:ixomax3, iw) = &
549 tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3)
550 end do
551 end if
552 end do
553 if(interp_phy_first_row) then
555 call phys_to_conserved(ixi^l,ixomin1,ix,ixomin3,ixomax1,ix,ixomax3,w,x)
556 endif
557 endif
558 case (5, 6)
559 if(interp_phy_first_row) then
560 if (ib == 5) then
561 ix = ixomax3+1
562 else
563 ix = ixomin3-1
564 end if
566 call phys_to_primitive(ixi^l,ixomin1,ixomax1,ixomin2,ixomax2,ix,ix,w,x)
567 endif
568 endif
569 if(bc_phy_first_row)then
570 if (ib == 5) then
571 ixoomax3=ixoomax3+1
572 else
573 ixoomin3=ixoomin3-1
574 endif
576 if (ib == 5) then
577 ixoomin3=ixoomax3
578 else
579 ixoomax3=ixoomin3
580 endif
581 call phys_to_primitive(ixi^l,ixoo^l,w,x)
582 if (ib == 5) then
583 ixoomin3=ixomin3 ! to be used later in to_conserved
584 else
585 ixoomax3=ixomax3
586 endif
587 endif
588 endif
589 do iw = 1, nwfluxbc
590 n_bc = bc_data_ix(iw, ib)
591 if (n_bc == -1) cycle
592
593 tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3) = bc_data_get_3d(n_bc, &
594 x(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3, 1), &
595 x(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3, 2), qt)
596
597 if(interp_phy_first_row) then
598
600 ! Approximate boundary value by linear interpolation to first ghost
601 ! cell, rest of ghost cells contains the same value
602 do i = 0, ixoomax3-ixoomin3
603 w(ixomin1:ixomax1, ixomin2:ixomax2, ixoomin3+i, iw) = &
604 2 * tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3) - &
605 w(ixomin1:ixomax1, ixomin2:ixomax2, ix, iw)
606 end do
607 else
608 if (ib==5) then
609 do i = 0, ixoomax3-ixoomin3
610 if (i==0) then
611 w(ixomin1:ixomax1, ixomin2:ixomax2, ixoomin3+i, iw) = &
612 4 * tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3) - &
613 3 * w(ixomin1:ixomax1, ixomin2:ixomax2, ix, iw)
614 elseif (i==1)then
615 w(ixomin1:ixomax1, ixomin2:ixomax2, ixoomin3+i, iw) = &
616 2 * tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3) - &
617 w(ixomin1:ixomax1, ixomin2:ixomax2, ix, iw)
618 else
619 w(ixomin1:ixomax1, ixomin2:ixomax2, ixoomin3+i, iw) = &
620 tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3)
621 endif
622
623 end do
624 else
625 do i = 0, ixoomax3-ixoomin3
626 if (i==ixoomax3-ixoomin3) then
627 w(ixomin1:ixomax1, ixomin2:ixomax2, ixoomin3+i, iw) = &
628 4 * tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3) - &
629 3 * w(ixomin1:ixomax1, ixomin2:ixomax2, ix, iw)
630 elseif (i==ixoomax3-ixoomin3-1)then
631 w(ixomin1:ixomax1, ixomin2:ixomax2, ixoomin3+i, iw) = &
632 2 * tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3) - &
633 w(ixomin1:ixomax1, ixomin2:ixomax2, ix, iw)
634 else
635 w(ixomin1:ixomax1, ixomin2:ixomax2, ixoomin3+i, iw) = &
636 tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3)
637 endif
638
639 end do
640 endif
641
642 endif
643 else
644 do i = ixoomin3,ixoomax3
645 w(ixomin1:ixomax1, ixomin2:ixomax2, i, iw) = &
646 tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3)
647 end do
648 end if
649 end do
650 if(interp_phy_first_row) then
652 call phys_to_conserved(ixi^l,ixomin1,ixomax1,ixomin2,ixomax2,ix,ix,w,x)
653 endif
654 endif
655 case default
656 call mpistop("bc_data_set: unknown iB")
657 end select
658 }
659
661 call phys_to_conserved(ixi^l,ixoo^l,w,x)
662 endif
663
664 end subroutine bc_data_set
665
666 subroutine read_vtk_structured_points(fname, bc)
667 character(len=*), intent(in) :: fname
668 type(bc_data_t), intent(out) :: bc
669 double precision, allocatable :: tmp_data(:, :, :, :)
670 integer, parameter :: max_variables = 10
671 character(len=40) :: tmp_names(max_variables)
672 character(len=256) :: line
673 character(len=40) :: word, typename
674 integer, parameter :: my_unit = 123
675 integer :: n, n_points_total
676 integer :: n_components
677
678 open(my_unit, file=trim(fname), status="old", action="read")
679
680 ! Header, e.g. # vtk DataFile Version 2.0
681 read(my_unit, "(A)") line
682
683 ! Dataset name
684 read(my_unit, "(A)") line
685
686 ! ASCII / BINARY
687 read(my_unit, "(A)") line
688
689 if (line /= "ASCII") then
690 print *, "line: ", trim(line)
691 error stop "read_vtk: not an ASCII file"
692 end if
693
694 ! DATASET STRUCTURED_POINTS
695 read(my_unit, "(A)") line
696
697 if (line /= "DATASET STRUCTURED_POINTS") then
698 print *, "line: ", trim(line)
699 error stop "read_vtk must have: DATASET STRUCTURED_POINTS"
700 end if
701
702 ! DIMENSIONS NX NY NZ
703 read(my_unit, "(A)") line
704 read(line, *) word, bc%n_points
705
706 if (word /= "DIMENSIONS") then
707 print *, "line: ", trim(line)
708 error stop "read_vtk expects: DIMENSIONS"
709 end if
710
711 ! SPACING DX DY DZ
712 read(my_unit, *) word, bc%dx
713
714 if (word /= "SPACING") then
715 print *, "line: ", trim(line)
716 error stop "read_vtk expects: SPACING"
717 end if
718
719 ! ORIGIN 0 0 0
720 read(my_unit, *) word, bc%origin
721 if (word /= "ORIGIN") then
722 print *, "line: ", trim(line)
723 error stop "read_vtk expects: ORIGIN"
724 end if
725
726 ! POINT_DATA NPOINTS
727 read(my_unit, *) word, n_points_total
728
729 if (word /= "POINT_DATA") then
730 print *, "line: ", trim(line)
731 error stop "read_vtk expects: POINT_DATA n_points"
732 end if
733
734 if (n_points_total /= product(bc%n_points)) &
735 error stop "read_vtk: n_points not equal to NX*NY*NZ"
736
737 allocate(tmp_data(bc%n_points(1), bc%n_points(2), bc%n_points(3), &
738 max_variables))
739
740 ! Read all scalar variables
741 do n = 1, max_variables
742
743 ! SCALARS name type ncomponents
744 read(my_unit, *, end=900) word, tmp_names(n), typename, n_components
745
746 if (word /= "SCALARS") then
747 print *, "line: ", trim(line)
748 error stop "read_vtk expects: SCALARS name type ncomponents"
749 end if
750
751 if (n_components /= 1) error stop "read_vtk: ncomponents should be 1"
752
753 ! LOOKUP_TABLE default
754 read(my_unit, *) word, typename
755
756 if (word /= "LOOKUP_TABLE") then
757 print *, "line: ", trim(line)
758 error stop "read_vtk expects: LOOKUP_TABLE name"
759 end if
760
761 ! Read list of data values
762 read(my_unit, *) tmp_data(:, :, :, n)
763 end do
764
765 ! Done reading variables
766900 continue
767
768 close(my_unit)
769
770 if (n == max_variables + 1) &
771 error stop "read_vtk: increase max_variables"
772
773 ! Loop index is one higher than number of variables
774 bc%n_variables = n-1
775 bc%values = tmp_data(:, :, :, 1:n-1)
776 bc%names = tmp_names(1:n-1)
777 end subroutine read_vtk_structured_points
778
779end module mod_bc_data
Module to set boundary conditions from user data.
Definition mod_bc_data.t:2
logical, public, protected interp_phy_first_row_same
Definition mod_bc_data.t:34
subroutine, public bc_data_set(qt, ixil, ixol, ib, w, x)
Set boundary conditions according to user data.
character(len=std_len), public, protected boundary_data_file_name
data file name
Definition mod_bc_data.t:32
elemental double precision function, public bc_data_get_2d(n_bc, x1, qt)
logical, public, protected bc_phy_first_row
Definition mod_bc_data.t:35
type(lt_t), dimension(max_boundary_conds) lt_1d
Definition mod_bc_data.t:20
type(lt3_t), dimension(max_boundary_conds), public lt_3d
Definition mod_bc_data.t:22
logical, public, protected boundary_data_primitive
Definition mod_bc_data.t:36
elemental double precision function, public bc_data_get_3d(n_bc, x1, x2, qt)
subroutine, public bc_data_init()
Definition mod_bc_data.t:46
integer, dimension(:, :), allocatable, public, protected bc_data_ix
Integer array for indexing lookup tables per variable per direction.
Definition mod_bc_data.t:29
logical, public, protected bc_data_time_varying
Whether boundary condition data is time varying.
Definition mod_bc_data.t:25
logical, public, protected interp_phy_first_row
Definition mod_bc_data.t:33
subroutine, public mpistop(message)
Exit MPI-AMRVAC with an error message.
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, parameter unitpar
file handle for IO
integer, dimension(:, :), allocatable typeboundary
Array indicating the type of boundary condition per variable and per physical boundary.
integer, parameter ndim
Number of spatial dimensions for grid variables.
character(len=std_len), dimension(:), allocatable par_files
Which par files are used as input.
A Fortran 90 module for creating 1D and 2D lookup tables. These tables can be used to efficiently int...
type(lt2_t) function, public lt2_create_from_data(x_min, x_max, spaced_data)
This function returns a new lookup table from regularly spaced data.
type(lt_t) function, public lt_create_from_data(x_min, x_max, spaced_data)
This function returns a new lookup table from regularly spaced data.
type(lt3_t) function, public lt3_create_from_data(x_min, x_max, spaced_data)
This function returns a new lookup table from regularly spaced data.
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition mod_physics.t:4
procedure(sub_convert), pointer phys_to_primitive
Definition mod_physics.t:55
procedure(sub_convert), pointer phys_to_conserved
Definition mod_physics.t:54
The 2D lookup table type.
The 3D lookup table type.
The lookup table type. There can be one or more columns, for which values can be looked up for a give...