64 class(*),
target,
intent(inout) :: data
65 type(
mesh_t),
pointer :: msh
66 real(kind=
dp),
pointer :: params(:)
67 character(len=3),
pointer :: cbc(:,:)
68 integer,
allocatable :: curve_type(:,:)
69 logical,
allocatable :: curve_element(:)
70 character(len=1) :: chtemp
71 integer :: ndim, nparam, nskip, nlogic, ncurve
72 integer :: nelgs, nelgv, i, j, ierr, l
73 integer :: el_idx, pt_idx
74 logical :: read_param, read_bcs, read_map
75 real(kind=
dp) :: xc(8), yc(8), zc(8),
curve(5)
76 real(kind=
dp),
allocatable :: bc_data(:,:,:)
77 real(kind=
dp),
allocatable :: curve_data(:,:,:)
81 character(len=1024) :: re2_fname, map_fname
82 integer :: start_el, end_el, nel, edge
86 integer :: sym_facet, pids(4), p_el_idx, p_facet
88 integer,
parameter,
dimension(6) :: facet_map = (/3, 2, 4, 1, 5, 6/)
89 logical :: curve_skip = .false.
90 character(len=LOG_SIZE) :: log_buf
92 call this%check_exists()
110 if (read_param .and. read_bcs .and.
pe_size .gt. 1)
then
111 call neko_error(
'Reading NEKTON session data only implemented in serial')
115 open(unit=9,
file=trim(this%fname), status=
'old', iostat=ierr)
116 call neko_log%message(
'Reading NEKTON file ' // this%fname)
123 if (.not. read_param)
then
129 allocate(params(nparam))
150 read(9, *) nelgs,ndim, nelgv
151 if (nelgs .lt. 0)
then
152 re2_fname = trim(this%fname(1:scan(trim(this%fname), &
153 '.', back=.true.)))//
're2'
157 write(log_buf,1) ndim, nelgv
158 1
format(
'gdim = ', i1,
', nelements =', i7)
162 inquire(
file=map_fname, exist=read_map)
168 call neko_log%warning(
'No NEKTON map file found')
173 nel = dist%num_local()
174 start_el = dist%start_idx() + 1
175 end_el = dist%end_idx() + 1
177 call msh%init(ndim, dist)
179 call htp%init(2**ndim * nel, ndim)
185 if (ndim .eq. 2)
then
186 read(9, *) (xc(j),j=1,4)
187 read(9, *) (yc(j),j=1,4)
188 if (i .ge. start_el .and. i .le. end_el)
then
194 call msh%add_element(el_idx, p(1), p(2), p(4), p(3))
196 else if (ndim .eq. 3)
then
197 read(9, *) (xc(j),j=1,4)
198 read(9, *) (yc(j),j=1,4)
199 read(9, *) (zc(j),j=1,4)
200 read(9, *) (xc(j),j=5,8)
201 read(9, *) (yc(j),j=5,8)
202 read(9, *) (zc(j),j=5,8)
203 if (i .ge. start_el .and. i .le. end_el)
then
209 call msh%add_element(el_idx, &
210 p(1), p(2), p(4), p(3), p(5), p(6), p(8), p(7))
213 if (i .ge. start_el .and. i .le. end_el)
then
222 allocate(curve_data(5,8,nelgv))
223 allocate(curve_element(nelgv))
224 allocate(curve_type(8,nelgv))
226 curve_element(i) = .false.
230 curve_data(l,j,i) = 0d0
235 read(9, *) edge, el_idx, (
curve(j),j=1,5), chtemp
237 curve_data(j,edge,el_idx) =
curve(j)
239 curve_element(el_idx) = .true.
240 select case(trim(chtemp))
242 curve_type(edge,el_idx) = 1
245 curve_type(edge,el_idx) = 2
248 curve_type(edge,el_idx) = 3
250 curve_type(edge,el_idx) = 4
254 call neko_log%warning(
'Curve type: s, e are not supported, treating mesh as non-curved.')
257 if (curve_element(el_idx))
then
258 call msh%mark_curve_element(el_idx, &
259 curve_data(1,1,el_idx), curve_type(1,el_idx))
263 deallocate(curve_data)
264 deallocate(curve_element)
265 deallocate(curve_type)
270 if (.not. read_bcs)
then
271 allocate(cbc(6,nelgv))
272 allocate(bc_data(6,2*ndim,nelgv))
275 if (nelgv .lt. 1000) off = 1
277 if (i .ge. start_el .and. i .le. end_el)
then
278 el_idx = i - start_el + 1
280 read(9, *) cbc(j, i), (bc_data(l,j,i),l=1,6)
281 sym_facet = facet_map(j)
282 select case(trim(cbc(j,i)))
284 call msh%mark_wall_facet(sym_facet, el_idx)
286 call msh%mark_inlet_facet(sym_facet, el_idx)
288 call msh%mark_outlet_facet(sym_facet, el_idx)
290 call msh%mark_outlet_normal_facet(sym_facet, el_idx)
292 call msh%mark_sympln_facet(sym_facet, el_idx)
294 p_el_idx = int(bc_data(2+off,j,i))
295 p_facet = facet_map(int(bc_data(3+off,j,i)))
296 call msh%get_facet_ids(sym_facet, el_idx, pids)
297 call msh%mark_periodic_facet(sym_facet, el_idx, &
298 p_facet, p_el_idx, pids)
304 if (i .ge. start_el .and. i .le. end_el)
then
305 el_idx = i - start_el + 1
307 sym_facet = facet_map(j)
308 select case(trim(cbc(j,i)))
310 p_el_idx = int(bc_data(2+off,j,i))
311 p_facet = facet_map(int(bc_data(3+off,j,i)))
312 call msh%create_periodic_ids(sym_facet, el_idx, &
319 if (i .ge. start_el .and. i .le. end_el)
then
320 el_idx = i - start_el + 1
322 sym_facet = facet_map(j)
323 select case(trim(cbc(j,i)))
325 p_el_idx = int(bc_data(2+off,j,i))
326 p_facet = facet_map(int(bc_data(3+off,j,i)))
327 call msh%create_periodic_ids(sym_facet, el_idx, &
334 if (i .ge. start_el .and. i .le. end_el)
then
335 el_idx = i - start_el + 1
337 sym_facet = facet_map(j)
338 select case(trim(cbc(j,i)))
340 p_el_idx = int(bc_data(2+off,j,i))
341 p_facet = facet_map(int(bc_data(3+off,j,i)))
342 call msh%create_periodic_ids(sym_facet, el_idx, &
351 allocate(cbc(6,nelgv))
354 read(9,
'(a1, a3)') chtemp, cbc(j, i)
369 class(*),
target,
intent(in) :: data
370 real(kind=
rp),
intent(in),
optional :: t
375 type(
point_t),
intent(inout) :: p
376 integer,
intent(inout) :: idx
379 if (htp%get(p, tmp) .gt. 0)
then
type(mpi_comm) neko_comm
MPI communicator.
integer pe_size
MPI size of communicator.
Defines a domain as a subset of facets in a mesh.
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 dp
integer, parameter, public rp
Global precision used in computations.
NEKTON mesh data in re2 format.
NEKTON session data reader.
subroutine rea_file_read(this, data)
Load NEKTON session data from an ascii file.
subroutine rea_file_add_point(htp, p, idx)
subroutine rea_file_write(this, data, t)
subroutine, public rea_free(r)
Free a NEKTON session data.
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 .
Interface for NEKTON re2 files.
NEKTON session data struct.
Interface for NEKTON ascii files.