66 class(*),
target,
intent(inout) :: data
67 type(
mesh_t),
pointer :: msh
68 character(len=5) :: hdr_ver
69 character(len=54) :: hdr_str
70 character(len=80) :: hdr_full
71 integer :: nel, ndim, nelv, ierr, nBCre2
72 type(mpi_status) :: status
74 integer (kind=MPI_OFFSET_KIND) :: mpi_offset
77 integer :: ncurv, nbcs
81 character(len=1024) :: map_fname
83 integer :: re2_data_xy_size
84 integer :: re2_data_xyz_size
85 integer :: re2_data_cv_size
86 integer :: re2_data_bc_size
88 character(len=LOG_SIZE) :: log_buf
90 call this%check_exists()
100 open(unit=9,
file=trim(this%fname), status=
'old', iostat=ierr)
101 call neko_log%message(
'Reading binary NEKTON file ' // this%fname)
103 read(9,
'(a80)') hdr_full
104 read(hdr_full,
'(a5)') hdr_ver
106 if (hdr_ver .eq.
'#v004')
then
107 read(hdr_full,
'(a5,i16,i3,i16,i4,a36)') hdr_ver, nel, ndim, nelv, nbcre2, hdr_str
109 else if (hdr_ver .eq.
'#v002' .or. hdr_ver .eq.
'#v003')
then
110 read(hdr_full,
'(a5,i9,i3,i9,a54)') hdr_ver, nel, ndim, nelv, hdr_str
112 else if (hdr_ver .eq.
'#v001')
then
113 read(hdr_full,
'(a5,i9,i3,i9,a54)') hdr_ver, nel, ndim, nelv, hdr_str
128 write(log_buf,1) ndim, nelv
129 1
format(
'gdim = ', i1,
', nelements =', i7)
135 inquire(
file=map_fname, exist=read_map)
141 call neko_log%warning(
'No NEKTON map file found')
144 call mpi_file_open(
neko_comm, trim(this%fname), &
145 mpi_mode_rdonly, mpi_info_null, fh, ierr)
147 if (ierr .ne. 0)
then
148 call neko_log%error(
"Can't open binary NEKTON file ")
153 call msh%init(ndim, dist)
158 call mpi_file_read_at_all(fh, mpi_offset, test, 1, mpi_real, status, ierr)
162 call neko_error(
'Invalid endian of re2 file, byte swap not implemented yet')
166 mpi_offset, re2_data_xy_size, re2_data_xyz_size, v2_format)
171 if (ndim .eq. 2)
then
172 mpi_offset = mpi_offset + int(dist%num_global(),
i8) * int(re2_data_xy_size,
i8)
174 mpi_offset = mpi_offset + int(dist%num_global(),
i8) * int(re2_data_xyz_size,
i8)
180 call mpi_file_read_at_all(fh, mpi_offset, t2, 1, mpi_double_precision, status, ierr)
184 mpi_offset = mpi_offset + int(ncurv,
i8) * int(re2_data_cv_size,
i8)
185 call mpi_file_read_at_all(fh, mpi_offset, t2, 1, mpi_double_precision, status, ierr)
191 call mpi_file_read_at_all(fh, mpi_offset, ncurv, 1, mpi_integer, status, ierr)
194 mpi_offset = mpi_offset + int(ncurv,
i8) * int(re2_data_cv_size,
i8)
195 call mpi_file_read_at_all(fh, mpi_offset, nbcs, 1, mpi_integer, status, ierr)
202 call mpi_file_close(fh, ierr)
213 class(*),
target,
intent(in) :: data
214 real(kind=
rp),
intent(in),
optional :: t
215 type(
re2v1_xy_t),
allocatable :: re2_data_xy(:)
217 type(
mesh_t),
pointer :: msh
218 character(len=5),
parameter :: RE2_HDR_VER =
'#v001'
219 character(len=54),
parameter :: RE2_HDR_STR =
'RE2 exported by NEKO'
220 integer :: i, j, ierr, nelgv
221 type(mpi_status) :: status
223 integer (kind=MPI_OFFSET_KIND) :: mpi_offset
224 integer :: element_offset
225 integer :: re2_data_xy_size
226 integer :: re2_data_xyz_size
227 character(len=LOG_SIZE) :: log_buf
238 call mpi_reduce(msh%nelv, nelgv, 1, &
239 mpi_integer, mpi_sum, 0,
neko_comm, ierr)
241 call mpi_exscan(msh%nelv, element_offset, 1, &
244 call neko_log%message(
'Writing data as a binary NEKTON file ' // this%fname)
247 open(unit=9,
file=trim(this%fname), status=
'new', iostat=ierr)
248 write(9,
'(a5,i9,i3,i9,a54)') re2_hdr_ver, nelgv, msh%gdim,&
254 call mpi_file_open(
neko_comm, trim(this%fname), &
255 mpi_mode_wronly + mpi_mode_create, mpi_info_null, fh, ierr)
259 mpi_real, status, ierr)
262 if (msh%gdim .eq. 2)
then
263 allocate(re2_data_xy(msh%nelv))
265 re2_data_xy(i)%rgroup = 1.0
267 re2_data_xy(i)%x(j) =
real(msh%elements(i)%e%pts(j)%p%x(1))
268 re2_data_xy(i)%y(j) =
real(msh%elements(i)%e%pts(j)%p%x(2))
271 mpi_offset = mpi_offset + element_offset * re2_data_xy_size
272 call mpi_file_write_at(fh, mpi_offset, &
275 deallocate(re2_data_xy)
277 else if (msh%gdim .eq. 3)
then
278 allocate(re2_data_xyz(msh%nelv))
280 re2_data_xyz(i)%rgroup = 1.0
282 re2_data_xyz(i)%x(j) =
real(msh%elements(i)%e%pts(j)%p%x(1))
283 re2_data_xyz(i)%y(j) =
real(msh%elements(i)%e%pts(j)%p%x(2))
284 re2_data_xyz(i)%z(j) =
real(msh%elements(i)%e%pts(j)%p%x(3))
287 mpi_offset = mpi_offset + element_offset * re2_data_xyz_size
288 call mpi_file_write_at(fh, mpi_offset, &
291 deallocate(re2_data_xyz)
296 call mpi_file_close(fh, ierr)
304 mpi_offset, re2_data_xy_size, re2_data_xyz_size, v2_format)
305 type(
mesh_t),
intent(inout) :: msh
306 integer (kind=MPI_OFFSET_KIND) :: mpi_offset
307 integer,
intent(inout) :: ndim
308 integer,
intent(inout) :: nel
309 type(mpi_file),
intent(inout) :: fh
310 integer,
intent(in) :: re2_data_xy_size
311 integer,
intent(in) :: re2_data_xyz_size
312 logical,
intent(in) :: v2_format
314 integer :: element_offset
315 type(
re2v1_xy_t),
allocatable :: re2v1_data_xy(:)
316 type(
re2v1_xyz_t),
allocatable :: re2v1_data_xyz(:)
317 type(
re2v2_xy_t),
allocatable :: re2v2_data_xy(:)
318 type(
re2v2_xyz_t),
allocatable :: re2v2_data_xyz(:)
319 type(mpi_status) :: status
322 integer :: pt_idx, nelv
323 integer :: i, j, ierr
326 nelv = dist%num_local()
327 element_offset = dist%start_idx()
329 call htp%init(2*nel, ndim)
331 if (ndim .eq. 2)
then
332 mpi_offset = mpi_offset + element_offset * re2_data_xy_size
333 if (.not. v2_format)
then
334 allocate(re2v1_data_xy(nelv))
335 call mpi_file_read_at_all(fh, mpi_offset, &
340 real(re2v1_data_xy(i)%y(j),dp), 0.0d0)
344 if(mod(i,nelv/10) .eq. 0)
write(*,*) i,
'elements read'
347 call msh%add_element(i, p(1), p(2), p(4), p(3))
349 deallocate(re2v1_data_xy)
351 allocate(re2v2_data_xy(nelv))
352 call mpi_file_read_at_all(fh, mpi_offset, &
356 p(j) =
point_t(re2v2_data_xy(i)%x(j), &
357 re2v2_data_xy(i)%y(j), 0.0d0)
361 if(mod(i,nelv/10) .eq. 0)
write(*,*) i,
'elements read'
364 call msh%add_element(i, p(1), p(2), p(4), p(3))
366 deallocate(re2v2_data_xy)
368 else if (ndim .eq. 3)
then
369 mpi_offset = mpi_offset + element_offset * re2_data_xyz_size
370 if (.not. v2_format)
then
371 allocate(re2v1_data_xyz(nelv))
372 call mpi_file_read_at_all(fh, mpi_offset, &
377 real(re2v1_data_xyz(i)%y(j),dp),&
378 real(re2v1_data_xyz(i)%z(j),dp))
382 if(mod(i,nelv/100) .eq. 0)
write(*,*) i,
'elements read'
385 call msh%add_element(i, &
386 p(1), p(2), p(4), p(3), p(5), p(6), p(8), p(7))
388 deallocate(re2v1_data_xyz)
390 allocate(re2v2_data_xyz(nelv))
391 call mpi_file_read_at_all(fh, mpi_offset, &
395 p(j) =
point_t(re2v2_data_xyz(i)%x(j), &
396 re2v2_data_xyz(i)%y(j),&
397 re2v2_data_xyz(i)%z(j))
401 if(mod(i,nelv/100) .eq. 0)
write(*,*) i,
'elements read'
404 call msh%add_element(i, &
405 p(1), p(2), p(4), p(3), p(5), p(6), p(8), p(7))
407 deallocate(re2v2_data_xyz)
415 type(
mesh_t),
intent(inout) :: msh
416 integer (kind=MPI_OFFSET_KIND) :: mpi_offset
417 integer,
intent(inout) :: ncurve
419 type(mpi_file),
intent(inout) :: fh
420 logical,
intent(in) :: v2_format
421 type(mpi_status) :: status
422 integer :: i, j, l, ierr, el_idx, id
425 real(kind=
dp),
allocatable :: curve_data(:,:,:)
426 integer,
allocatable :: curve_type(:,:)
427 logical,
allocatable :: curve_element(:)
428 character(len=1) :: chtemp
429 logical :: curve_skip = .false.
431 allocate(curve_data(5,12,msh%nelv))
432 allocate(curve_element(msh%nelv))
433 allocate(curve_type(12,msh%nelv))
435 curve_element(i) = .false.
439 curve_data(l,j,i) = 0d0
443 write(*,*)
'reading curved data'
444 if (.not. v2_format)
then
445 allocate(re2v1_data_curve(ncurve))
446 call mpi_file_read_at_all(fh, mpi_offset, re2v1_data_curve, ncurve, &
449 allocate(re2v2_data_curve(ncurve))
450 call mpi_file_read_at_all(fh, mpi_offset, re2v2_data_curve, ncurve, &
456 el_idx = re2v2_data_curve(i)%elem - dist%start_idx()
457 id = re2v2_data_curve(i)%zone
458 chtemp = re2v2_data_curve(i)%type
460 curve_data(j,id, el_idx) = re2v2_data_curve(i)%point(j)
463 el_idx = re2v1_data_curve(i)%elem - dist%start_idx()
464 id = re2v1_data_curve(i)%zone
465 chtemp = re2v1_data_curve(i)%type
467 curve_data(j,id, el_idx) =
real(re2v1_data_curve(i)%point(j),
dp)
471 curve_element(el_idx) = .true.
473 select case(trim(chtemp))
475 curve_type(id,el_idx) = 1
476 call neko_log%warning(
'curve type s not supported, treating mesh as non-curved')
480 curve_type(id,el_idx) = 2
481 call neko_log%warning(
'curve type e not supported, treating mesh as non-curved')
485 curve_type(id,el_idx) = 3
487 curve_type(id,el_idx) = 4
489 write(*,*) chtemp,
'curve type not supported yet, treating mesh as non-curved',id, el_idx
496 deallocate(re2v2_data_curve)
498 deallocate(re2v1_data_curve)
500 if (.not. curve_skip)
then
501 do el_idx = 1, msh%nelv
502 if (curve_element(el_idx))
then
503 call msh%mark_curve_element(el_idx, &
504 curve_data(1,1,el_idx), curve_type(1,el_idx))
509 deallocate(curve_data)
510 deallocate(curve_element)
511 deallocate(curve_type)
516 type(
mesh_t),
intent(inout) :: msh
517 integer (kind=MPI_OFFSET_KIND) :: mpi_offset
518 integer,
intent(inout) :: nbcs
520 type(mpi_file),
intent(inout) :: fh
521 logical,
intent(in) :: v2_format
522 type(mpi_status) :: status
524 integer :: sym_facet, label
525 integer :: p_el_idx, p_facet
526 integer :: i, j, ierr, el_idx
527 integer,
parameter,
dimension(6) :: facet_map = (/3, 2, 4, 1, 5, 6/)
529 type(
re2v1_bc_t),
allocatable :: re2v1_data_bc(:)
530 type(
re2v2_bc_t),
allocatable :: re2v2_data_bc(:)
532 if (.not. v2_format)
then
533 allocate(re2v1_data_bc(nbcs))
534 call mpi_file_read_at_all(fh, mpi_offset, re2v1_data_bc, nbcs, &
537 allocate(re2v2_data_bc(nbcs))
538 call mpi_file_read_at_all(fh, mpi_offset, re2v2_data_bc, nbcs, &
547 el_idx = int(re2v2_data_bc(i)%elem) - dist%start_idx()
548 sym_facet = facet_map(int(re2v2_data_bc(i)%face))
550 select case(trim(re2v2_data_bc(i)%type))
552 call msh%mark_wall_facet(sym_facet, el_idx)
554 call msh%mark_inlet_facet(sym_facet, el_idx)
556 call msh%mark_outlet_facet(sym_facet, el_idx)
558 call msh%mark_outlet_normal_facet(sym_facet, el_idx)
560 call msh%mark_sympln_facet(sym_facet, el_idx)
563 p_el_idx = int(re2v2_data_bc(i)%bc_data(1))
564 p_facet = facet_map(int(re2v2_data_bc(i)%bc_data(2)))
565 call msh%get_facet_ids(sym_facet, el_idx, pids)
566 call msh%mark_periodic_facet(sym_facet, el_idx, &
567 p_facet, p_el_idx, pids)
568 case (
'MSH',
'msh',
'EXO',
'exo')
569 label = int(re2v2_data_bc(i)%bc_data(5))
571 call neko_error(
'Invalid label id (valid range [1,...,20])')
573 call msh%mark_labeled_facet(sym_facet, el_idx, label)
575 write (*,*) re2v2_data_bc(i)%type,
'bc type not supported yet'
576 write (*,*) re2v2_data_bc(i)%bc_data
586 el_idx = re2v2_data_bc(i)%elem - dist%start_idx()
587 sym_facet = facet_map(int(re2v2_data_bc(i)%face))
588 select case(trim(re2v2_data_bc(i)%type))
590 p_el_idx = int(re2v2_data_bc(i)%bc_data(1))
591 p_facet = facet_map(int(re2v2_data_bc(i)%bc_data(2)))
592 call msh%create_periodic_ids(sym_facet, el_idx, &
598 deallocate(re2v2_data_bc)
602 el_idx = re2v1_data_bc(i)%elem - dist%start_idx()
603 sym_facet = facet_map(re2v1_data_bc(i)%face)
604 select case(trim(re2v1_data_bc(i)%type))
606 call msh%mark_wall_facet(sym_facet, el_idx)
608 call msh%mark_inlet_facet(sym_facet, el_idx)
610 call msh%mark_outlet_facet(sym_facet, el_idx)
612 call msh%mark_sympln_facet(sym_facet, el_idx)
615 p_el_idx = int(re2v1_data_bc(i)%bc_data(1))
616 p_facet = facet_map(int(re2v1_data_bc(i)%bc_data(2)))
617 call msh%get_facet_ids(sym_facet, el_idx, pids)
618 call msh%mark_periodic_facet(sym_facet, el_idx, &
619 p_facet, p_el_idx, pids)
621 label = int(re2v1_data_bc(i)%bc_data(5))
622 call msh%mark_labeled_facet(sym_facet, el_idx, label)
624 write (*,*) re2v1_data_bc(i)%type,
'bc type not supported yet'
625 write (*,*) re2v1_data_bc(i)%bc_data
636 el_idx = re2v1_data_bc(i)%elem - dist%start_idx()
637 sym_facet = facet_map(re2v1_data_bc(i)%face)
638 select case(trim(re2v1_data_bc(i)%type))
640 p_el_idx = int(re2v1_data_bc(i)%bc_data(1))
641 p_facet = facet_map(int(re2v1_data_bc(i)%bc_data(2)))
642 call msh%create_periodic_ids(sym_facet, el_idx, &
649 deallocate(re2v1_data_bc)
657 type(
point_t),
intent(inout) :: p
658 integer,
intent(inout) :: idx
661 if (htp%get(p, tmp) .gt. 0)
then
type(mpi_comm) neko_comm
MPI communicator.
integer pe_size
MPI size of communicator.
Defines practical data distributions.
Module for file I/O operations.
Implements a hash table ADT.
type(log_t), public neko_log
Global log stream.
integer, parameter, public neko_msh_max_zlbls
Max num. zone labels.
integer, public mpi_double_precision_size
Size of MPI type double precision.
type(mpi_datatype), public mpi_re2v2_data_xyz
MPI derived type for 3D NEKTON re2 data.
type(mpi_datatype), public mpi_re2v2_data_xy
MPI derived type for 2D NEKTON re2 data.
integer, public mpi_character_size
Size of MPI type character.
integer, public mpi_real_size
Size of MPI type real.
integer, public mpi_integer_size
Size of MPI type integer.
type(mpi_datatype), public mpi_re2v1_data_cv
MPI derived type for NEKTON re2 cv data.
type(mpi_datatype), public mpi_re2v2_data_cv
MPI derived type for NEKTON re2 cv data.
type(mpi_datatype), public mpi_re2v1_data_xy
MPI derived type for 2D NEKTON re2 data.
type(mpi_datatype), public mpi_re2v1_data_bc
MPI derived type for NEKTON re2 bc data.
type(mpi_datatype), public mpi_re2v2_data_bc
MPI derived type for NEKTON re2 bc data.
type(mpi_datatype), public mpi_re2v1_data_xyz
MPI derived type for 3D NEKTON re2 data.
integer, parameter, public i8
integer, parameter, public dp
integer, parameter, public sp
integer, parameter, public rp
Global precision used in computations.
NEKTON mesh data in re2 format.
subroutine re2_file_read(this, data)
Load a binary NEKTON mesh from a re2 file.
subroutine re2_file_write(this, data, t)
subroutine re2_file_read_bcs(msh, nbcs, dist, fh, mpi_offset, v2_format)
subroutine re2_file_add_point(htp, p, idx)
subroutine re2_file_read_curve(msh, ncurve, dist, fh, mpi_offset, v2_format)
subroutine re2_file_read_points(msh, ndim, nel, dist, fh, mpi_offset, re2_data_xy_size, re2_data_xyz_size, v2_format)
real(kind=sp), parameter re2_endian_test
NEKTION re2 endian test.
integer, parameter re2_hdr_size
NEKTON re2 header size.
subroutine filename_chsuffix(fname, new_fname, new_suffix)
Change a filename's suffix.
Load-balanced linear distribution .
Interface for NEKTON map files.
A point in with coordinates .
NEKTON re2 bc data (version 1)
NEKTON re2 curve data (version 1)
NEKTON re2 element data (2d) (version 1)
NEKTON re2 element data (3d) (version 1)
NEKTON re2 bc data (version 2)
NEKTON re2 curve data (version 2)
NEKTON re2 element data (2d) (version 2)
NEKTON re2 element data (3d) (version 2)
Interface for NEKTON re2 files.