9 integer,
parameter :: max_boundary_conds = 10
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(:, :, :, :)
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)
28 integer,
public,
protected,
allocatable ::
bc_data_ix(:, :)
43 integer :: i, iw, ib, n_files, n_bc
44 double precision :: xmax(3)
61 xmax = bc%origin + (bc%n_points-1) * bc%dx
66 call mpistop(
"bc_data_init: only some files are time varying")
70 call mpistop(
"bc_data_init: 1D case not supported")
75 xmax(1:
ndim), bc%values(:, :, 1:1, 1))
79 xmax(1), bc%values(:, 1, 1:1, 1))
85 xmax(1:
ndim), bc%values(:, :, :, 1:1))
89 xmax(1:
ndim-1), bc%values(:, :, 1:1, 1))
101 character(len=*),
intent(in) :: files(:)
106 do n = 1,
size(files)
107 open(
unitpar, file=trim(files(n)), status=
"old")
108 read(
unitpar, bd_list,
end=111)
115 integer,
intent(in) :: n_bc
116 double precision,
intent(in) :: x1, x2, qt
117 double precision :: val
120 val = lt3_get_col(lt_3d(n_bc), 1, x1, x2, qt)
122 val = lt2_get_col(lt_2d(n_bc), 1, x1, x2)
127 integer,
intent(in) :: n_bc
128 double precision,
intent(in) :: x1, qt
129 double precision :: val
132 val = lt2_get_col(lt_2d(n_bc), 1, x1, qt)
134 val = lt_get_col(lt_1d(n_bc), 1, x1)
141 integer,
intent(in) :: ixi^
l, ixo^
l, ib
142 double precision,
intent(in) :: qt, x(ixi^s,1:
ndim)
143 double precision,
intent(inout) :: w(ixi^s,1:nw)
144 double precision :: tmp(ixo^s)
145 integer :: i, ix, iw, n_bc
152 if (n_bc == -1) cycle
155 x(ixomin1, ixomin2:ixomax2, 2), qt)
165 do i = 0, ixomax1-ixomin1
166 w(ixomin1+i, ixomin2:ixomax2, iw) = &
167 2 * tmp(ixomin1, ixomin2:ixomax2) - &
168 w(ix, ixomin2:ixomax2, iw)
174 if (n_bc == -1) cycle
177 x(ixomin1:ixomax1, ixomin2, 1), qt)
187 do i = 0, ixomax2-ixomin2
188 w(ixomin1:ixomax1, ixomin2+i, iw) = &
189 2 * tmp(ixomin1:ixomax1, ixomin2) - &
190 w(ixomin1:ixomax1, ix, iw)
194 call mpistop(
"bc_data_set: unknown iB")
203 if (n_bc == -1) cycle
205 tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3) =
bc_data_get_3d(n_bc, &
206 x(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3, 2), &
207 x(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3, 3), qt)
217 do i = 0, ixomax1-ixomin1
218 w(ixomin1+i, ixomin2:ixomax2, ixomin3:ixomax3, iw) = &
219 2 * tmp(ixomin1, ixomin2:ixomax2, ixomin3:ixomax3) - &
220 w(ix, ixomin2:ixomax2, ixomin3:ixomax3, iw)
226 if (n_bc == -1) cycle
228 tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3) =
bc_data_get_3d(n_bc, &
229 x(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3, 1), &
230 x(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3, 3), qt)
240 do i = 0, ixomax2-ixomin2
241 w(ixomin1:ixomax1, ixomin2+i, ixomin3:ixomax3, iw) = &
242 2 * tmp(ixomin1:ixomax1, ixomin2, ixomin3:ixomax3) - &
243 w(ixomin1:ixomax1, ix, ixomin3:ixomax3, iw)
249 if (n_bc == -1) cycle
251 tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3) =
bc_data_get_3d(n_bc, &
252 x(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3, 1), &
253 x(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3, 2), qt)
263 do i = 0, ixomax3-ixomin3
264 w(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3+i, iw) = &
265 2 * tmp(ixomin1:ixomax1, ixomin2:ixomax2, ixomin3) - &
266 w(ixomin1:ixomax1, ixomin2:ixomax2, ix, iw)
270 call mpistop(
"bc_data_set: unknown iB")
277 character(len=*),
intent(in) :: fname
278 type(bc_data_t),
intent(out) :: bc
279 double precision,
allocatable :: tmp_data(:, :, :, :)
280 integer,
parameter :: max_variables = 10
281 character(len=40) :: tmp_names(max_variables)
282 character(len=256) :: line
283 character(len=40) :: word, typename
284 integer,
parameter :: my_unit = 123
285 integer :: n, n_points_total
286 integer :: n_components
288 open(my_unit, file=trim(fname), status=
"old", action=
"read")
291 read(my_unit,
"(A)") line
294 read(my_unit,
"(A)") line
297 read(my_unit,
"(A)") line
299 if (line /=
"ASCII")
then
300 print *,
"line: ", trim(line)
301 error stop
"read_vtk: not an ASCII file"
305 read(my_unit,
"(A)") line
307 if (line /=
"DATASET STRUCTURED_POINTS")
then
308 print *,
"line: ", trim(line)
309 error stop
"read_vtk must have: DATASET STRUCTURED_POINTS"
313 read(my_unit,
"(A)") line
314 read(line, *) word, bc%n_points
316 if (word /=
"DIMENSIONS")
then
317 print *,
"line: ", trim(line)
318 error stop
"read_vtk expects: DIMENSIONS"
322 read(my_unit, *) word, bc%dx
324 if (word /=
"SPACING")
then
325 print *,
"line: ", trim(line)
326 error stop
"read_vtk expects: SPACING"
330 read(my_unit, *) word, bc%origin
331 if (word /=
"ORIGIN")
then
332 print *,
"line: ", trim(line)
333 error stop
"read_vtk expects: ORIGIN"
337 read(my_unit, *) word, n_points_total
339 if (word /=
"POINT_DATA")
then
340 print *,
"line: ", trim(line)
341 error stop
"read_vtk expects: POINT_DATA n_points"
344 if (n_points_total /= product(bc%n_points)) &
345 error stop
"read_vtk: n_points not equal to NX*NY*NZ"
347 allocate(tmp_data(bc%n_points(1), bc%n_points(2), bc%n_points(3), &
351 do n = 1, max_variables
354 read(my_unit, *,
end=900) word, tmp_names(n), typename, n_components
356 if (word /=
"SCALARS")
then
357 print *,
"line: ", trim(line)
358 error stop
"read_vtk expects: SCALARS name type ncomponents"
361 if (n_components /= 1) error stop
"read_vtk: ncomponents should be 1"
364 read(my_unit, *) word, typename
366 if (word /=
"LOOKUP_TABLE")
then
367 print *,
"line: ", trim(line)
368 error stop
"read_vtk expects: LOOKUP_TABLE name"
372 read(my_unit, *) tmp_data(:, :, :, n)
380 if (n == max_variables + 1) &
381 error stop
"read_vtk: increase max_variables"
385 bc%values = tmp_data(:, :, :, 1:n-1)
386 bc%names = tmp_names(1:n-1)
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
Module to set boundary conditions from user data.
subroutine read_vtk_structured_points(fname, bc)
character(len=std_len), public, protected boundary_data_file_name
data file name
elemental double precision function, public bc_data_get_2d(n_bc, x1, qt)
subroutine bc_read_params(files)
Read this module"s parameters from a file.
elemental double precision function, public bc_data_get_3d(n_bc, x1, x2, qt)
subroutine, public bc_data_set(qt, ixIL, ixOL, iB, w, x)
Set boundary conditions according to user data.
subroutine, public bc_data_init()
integer, dimension(:, :), allocatable, public, protected bc_data_ix
Integer array for indexing lookup tables per variable per direction.
logical, public, protected bc_data_time_varying
Whether boundary condition data is time varying.
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.
integer, parameter bc_data
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.