Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
scalar_ic.f90
Go to the documentation of this file.
1! Copyright (c) 2021-2025, The Neko Authors
2! All rights reserved.
3!
4! Redistribution and use in source and binary forms, with or without
5! modification, are permitted provided that the following conditions
6! are met:
7!
8! * Redistributions of source code must retain the above copyright
9! notice, this list of conditions and the following disclaimer.
10!
11! * Redistributions in binary form must reproduce the above
12! copyright notice, this list of conditions and the following
13! disclaimer in the documentation and/or other materials provided
14! with the distribution.
15!
16! * Neither the name of the authors nor the names of its
17! contributors may be used to endorse or promote products derived
18! from this software without specific prior written permission.
19!
20! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31! POSSIBILITY OF SUCH DAMAGE.
32!
35 use gather_scatter, only : gs_t, gs_op_add
37 use num_types, only : rp
38 use device_math, only : device_col2
40 use field, only : field_t
43 use coefs, only : coef_t
44 use math, only : col2, cfill, cfill_mask
46 use json_module, only : json_file
48 use point_zone, only: point_zone_t
51 use logger, only: neko_log, log_size
53 use fld_file, only: fld_file_t
54 use checkpoint, only: chkp_t
55 use file, only: file_t
58 use space, only: space_t, gll
59 use field_list, only : field_list_t
60 implicit none
61 private
62
63 interface set_scalar_ic
64 module procedure set_scalar_ic_int, set_scalar_ic_usr
65 end interface set_scalar_ic
66
67 public :: set_scalar_ic
68
69contains
70
83 subroutine set_scalar_ic_int(s, coef, gs, type, params, i)
84 type(field_t), intent(inout) :: s
85 type(coef_t), intent(in) :: coef
86 type(gs_t), intent(inout) :: gs
87 character(len=*) :: type
88 type(json_file), intent(inout) :: params
89 integer, intent(in) :: i
90
91 ! Variables for retrieving JSON parameters
92 real(kind=rp) :: ic_value
93 character(len=:), allocatable :: read_str
94 character(len=NEKO_FNAME_LEN) :: fname, mesh_fname
95 real(kind=rp) :: zone_value, tol
96 logical :: interpolate
97
98 if (trim(type) .eq. 'uniform') then
99
100 call json_get(params, 'value', ic_value)
101 call set_scalar_ic_uniform(s, ic_value)
102
103 else if (trim(type) .eq. 'point_zone') then
104
105 call json_get(params, 'base_value', ic_value)
106 call json_get(params, 'zone_name', read_str)
107 call json_get(params, 'zone_value', zone_value)
108
109 call set_scalar_ic_point_zone(s, ic_value, read_str, zone_value)
110
111 else if (trim(type) .eq. 'field') then
112
113 call json_get(params, 'file_name', read_str)
114 fname = trim(read_str)
115 call json_get_or_default(params, 'interpolate', interpolate, &
116 .false.)
117 call json_get_or_default(params, 'tolerance', tol, 0.000001_rp)
118 call json_get_or_default(params, 'mesh_file_name', read_str, &
119 "none")
120 mesh_fname = trim(read_str)
121
122 call set_scalar_ic_fld(s, fname, interpolate, tol, mesh_fname, i)
123
124 else
125 call neko_error('Invalid initial condition')
126 end if
127
128 call set_scalar_ic_common(s, coef, gs)
129
130 end subroutine set_scalar_ic_int
131
139 subroutine set_scalar_ic_usr(scheme_name, s, coef, gs, user_proc)
140 character(len=*), intent(in) :: scheme_name
141 type(field_t), target, intent(inout) :: s
142 type(coef_t), intent(in) :: coef
143 type(gs_t), intent(inout) :: gs
144 procedure(user_initial_conditions_intf) :: user_proc
145 type(field_list_t) :: fields
146
147 call neko_log%message("Type: user")
148
149 call fields%init(1)
150 call fields%assign_to_field(1, s)
151
152 call user_proc(scheme_name, fields)
153 call set_scalar_ic_common(s, coef, gs)
154
155 end subroutine set_scalar_ic_usr
156
163 subroutine set_scalar_ic_common(s, coef, gs)
164 type(field_t), intent(inout) :: s
165 type(coef_t), intent(in) :: coef
166 type(gs_t), intent(inout) :: gs
167 integer :: n
168
169 n = s%dof%size()
170 if (neko_bcknd_device .eq. 1) then
171 call device_memcpy(s%x, s%x_d, n, host_to_device, sync = .false.)
172 end if
173
174 ! Ensure continuity across elements for initial conditions
175 call gs%op(s%x, n, gs_op_add)
176
177 if (neko_bcknd_device .eq. 1) then
178 call device_col2(s%x_d, coef%mult_d, n)
179 else
180 call col2(s%x, coef%mult, n)
181 end if
182
183 end subroutine set_scalar_ic_common
184
189 subroutine set_scalar_ic_uniform(s, ic_value)
190 type(field_t), intent(inout) :: s
191 real(kind=rp), intent(in) :: ic_value
192 integer :: n
193 character(len=LOG_SIZE) :: log_buf
194
195 call neko_log%message("Type : uniform")
196 write (log_buf, '(A,ES12.6)') "Value: ", ic_value
197 call neko_log%message(log_buf)
198
199 s = ic_value
200 n = s%dof%size()
201 if (neko_bcknd_device .eq. 1) then
202 call cfill(s%x, ic_value, n)
203 end if
204
205 end subroutine set_scalar_ic_uniform
206
214 subroutine set_scalar_ic_point_zone(s, base_value, zone_name, zone_value)
215 type(field_t), intent(inout) :: s
216 real(kind=rp), intent(in) :: base_value
217 character(len=*), intent(in) :: zone_name
218 real(kind=rp), intent(in) :: zone_value
219
220 ! Internal variables
221 character(len=LOG_SIZE) :: log_buf
222 class(point_zone_t), pointer :: zone
223 integer :: size
224
225 call neko_log%message("Type : point_zone")
226 write (log_buf, '(A,ES12.6)') "Base value: ", base_value
227 call neko_log%message(log_buf)
228 call neko_log%message("Zone name : " // trim(zone_name))
229 write (log_buf, '(A,ES12.6)') "Zone value: ", zone_value
230 call neko_log%message(log_buf)
231
232 size = s%dof%size()
233 zone => neko_point_zone_registry%get_point_zone(trim(zone_name))
234
235 call set_scalar_ic_uniform(s, base_value)
236 call cfill_mask(s%x, zone_value, size, zone%mask%get(), zone%size)
237
238 end subroutine set_scalar_ic_point_zone
239
255 subroutine set_scalar_ic_fld(s, file_name, &
256 interpolate, tolerance, mesh_file_name, i)
257 type(field_t), intent(inout) :: s
258 character(len=*), intent(in) :: file_name
259 logical, intent(in) :: interpolate
260 real(kind=rp), intent(in) :: tolerance
261 character(len=*), intent(inout) :: mesh_file_name
262 integer, intent(in) :: i
263
264 character(len=LOG_SIZE) :: log_buf
265 integer :: sample_idx, sample_mesh_idx
266 integer :: last_index
267 type(fld_file_data_t) :: fld_data
268 type(file_t) :: f
269 logical :: mesh_mismatch
270
271 ! ---- For the mesh to mesh interpolation
272 type(global_interpolation_t) :: global_interp
273 ! -----
274
275 ! ---- For space to space interpolation
276 type(space_t) :: prev_Xh
277 type(interpolator_t) :: space_interp
278 ! ----
279
280 call neko_log%message("Type : field")
281 call neko_log%message("File name : " // trim(file_name))
282 write (log_buf, '(A,L1)') "Interpolation : ", interpolate
283 call neko_log%message(log_buf)
284
285 ! Extract sample index from the file name
286 sample_idx = extract_fld_file_index(file_name, -1)
287
288 if (sample_idx .eq. -1) then
289 call neko_error("Invalid file name for the initial condition. The " // &
290 "file format must be e.g. 'mean0.f00001'")
291 end if
292
293 ! Change from "field0.f000*" to "field0.fld" for the fld reader
294 call filename_chsuffix(file_name, file_name, 'fld')
295
296 call fld_data%init
297 call f%init(trim(file_name))
298
299 if (interpolate) then
300
301 ! If no mesh file is specified, use the default file name
302 if (mesh_file_name .eq. "none") then
303 mesh_file_name = trim(file_name)
304 sample_mesh_idx = sample_idx
305 else
306
307 ! Extract sample index from the mesh file name
308 sample_mesh_idx = extract_fld_file_index(mesh_file_name, -1)
309
310 if (sample_mesh_idx .eq. -1) then
311 call neko_error("Invalid file name for the initial condition." // &
312 " The file format must be e.g. 'mean0.f00001'")
313 end if
314
315 write (log_buf, '(A,ES12.6)') "Tolerance : ", tolerance
316 call neko_log%message(log_buf)
317 write (log_buf, '(A,A)') "Mesh file : ", &
318 trim(mesh_file_name)
319 call neko_log%message(log_buf)
320
321 end if ! if mesh_file_name .eq. none
322
323 ! Read the mesh coordinates if they are not in our fld file
324 if (sample_mesh_idx .ne. sample_idx) then
325 call f%set_counter(sample_mesh_idx)
326 call f%read(fld_data)
327 end if
328
329 end if
330
331 ! Read the field file containing (u,v,w,p)
332 call f%set_counter(sample_idx)
333 call f%read(fld_data)
334
335 !
336 ! Check that the data in the fld file matches the current case.
337 ! Note that this is a safeguard and there are corner cases where
338 ! two different meshes have the same dimension and same # of elements
339 ! but this should be enough to cover obvious cases.
340 !
341 mesh_mismatch = (fld_data%glb_nelv .ne. s%msh%glb_nelv .or. &
342 fld_data%gdim .ne. s%msh%gdim)
343
344 if (mesh_mismatch .and. .not. interpolate) then
345 call neko_error("The fld file must match the current mesh! " // &
346 "Use 'interpolate': 'true' to enable interpolation.")
347 else if (.not. mesh_mismatch .and. interpolate) then
348 call neko_log%warning("You have activated interpolation but you " // &
349 "might still be using the same mesh.")
350 end if
351
352
353 ! Mesh interpolation if specified
354 if (interpolate) then
355 ! Issue a warning if the mesh is in single precision
356 select type (ft => f%file_type)
357 type is (fld_file_t)
358 if (.not. ft%dp_precision) then
359 call neko_warning("The coordinates read from the field file " // &
360 "are in single precision.")
361 call neko_log%message("It is recommended to use a mesh in " // &
362 "double precision for better interpolation results.")
363 call neko_log%message("If the interpolation does not work, " // &
364 "you can try to increase the tolerance.")
365 end if
366 class default
367 end select
368
369 ! Copy all fld data to device since the reader loads everything on the host
370 call fld_data%x%copy_from(host_to_device, .false.)
371 call fld_data%y%copy_from(host_to_device, .false.)
372 call fld_data%z%copy_from(host_to_device, .false.)
373 call fld_data%t%copy_from(host_to_device, .true.)
374
375 ! Generates an interpolator object and performs the point search
376 call fld_data%generate_interpolator(global_interp, s%dof, s%msh, &
377 tolerance)
378
379 ! Evaluate scalar
380
381 ! i == 0 means it's the temperature field
382 if (i .ne. 0) then
383 call global_interp%evaluate(s%x, fld_data%s(i)%x, .false.)
384 else
385 call global_interp%evaluate(s%x, fld_data%t%x, .false.)
386 end if
387
388 call global_interp%free
389
390 ! Copy back to the host for set_scalar_ic_common
391 call fld_data%t%copy_from(device_to_host, .true.)
392
393 else ! No interpolation, just potentially from different spaces
394
395 ! Build a space_t object from the data in the fld file
396 call prev_xh%init(gll, fld_data%lx, fld_data%ly, fld_data%lz)
397 call space_interp%init(s%Xh, prev_xh)
398
399 ! Do the space-to-space interpolation
400 ! i == 0 means it's the temperature field
401 if (i .ne. 0) then
402 call space_interp%map_host(s%x, fld_data%s(i)%x, fld_data%nelv, s%Xh)
403 else
404 call space_interp%map_host(s%x, fld_data%t%x, fld_data%nelv, s%Xh)
405 end if
406
407 call space_interp%free
408
409 end if
410
411 call fld_data%free
412
413 end subroutine set_scalar_ic_fld
414
415end module scalar_ic
Copy data between host and device (or device and device)
Definition device.F90:71
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
Retrieves a parameter by name or throws an error.
Abstract interface for user defined initial conditions.
Definition user_intf.f90:65
Defines a checkpoint.
Coefficients.
Definition coef.f90:34
subroutine, public device_col2(a_d, b_d, n, strm)
Vector multiplication .
Device abstraction, common interface for various accelerators.
Definition device.F90:34
integer, parameter, public host_to_device
Definition device.F90:47
integer, parameter, public device_to_host
Definition device.F90:47
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Defines a field.
Definition field.f90:34
Module for file I/O operations.
Definition file.f90:34
Simple module to handle fld file series. Provides an interface to the different fields sotred in a fl...
NEKTON fld file format.
Definition fld_file.f90:35
Gather-scatter.
Implements global_interpolation given a dofmap.
Routines to interpolate between different spaces.
Utilities for retrieving parameters from the case files.
Logging routines.
Definition log.f90:34
type(log_t), public neko_log
Global log stream.
Definition log.f90:70
integer, parameter, public log_size
Definition log.f90:46
Definition math.f90:60
subroutine, public cfill(a, c, n)
Set all elements to a constant c .
Definition math.f90:487
subroutine, public col2(a, b, n)
Vector multiplication .
Definition math.f90:854
subroutine, public cfill_mask(a, c, n, mask, n_mask)
Fill a constant to a masked vector. .
Definition math.f90:397
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
type(point_zone_registry_t), target, public neko_point_zone_registry
Global point_zone registry.
Scalar initial condition.
Definition scalar_ic.f90:34
subroutine set_scalar_ic_uniform(s, ic_value)
Uniform initial condition.
subroutine set_scalar_ic_point_zone(s, base_value, zone_name, zone_value)
Point zone initial condition.
subroutine set_scalar_ic_common(s, coef, gs)
Set scalar initial condition (common)
subroutine set_scalar_ic_int(s, coef, gs, type, params, i)
Set scalar initial condition (builtin)
Definition scalar_ic.f90:84
subroutine set_scalar_ic_fld(s, file_name, interpolate, tolerance, mesh_file_name, i)
Set the initial condition of the scalar based on a field. @detail The field is read from an fld file....
subroutine set_scalar_ic_usr(scheme_name, s, coef, gs, user_proc)
Set scalar intial condition (user defined)
Defines a function space.
Definition space.f90:34
integer, parameter, public gll
Definition space.f90:49
Interfaces for user interaction with NEKO.
Definition user_intf.f90:34
Utilities.
Definition utils.f90:35
integer function, public extract_fld_file_index(fld_filename, default_index)
Extracts the index of a field file. For example, "myfield.f00045" will return 45. If the suffix of th...
Definition utils.f90:157
integer, parameter, public neko_fname_len
Definition utils.f90:40
subroutine, public neko_warning(warning_msg)
Reports a warning to standard output.
Definition utils.f90:346
subroutine, public filename_chsuffix(fname, new_fname, new_suffix)
Change a filename's suffix.
Definition utils.f90:140
subroutine, public filename_suffix(fname, suffix)
Extract a filename's suffix.
Definition utils.f90:108
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:55
field_list_t, To be able to group fields together
A wrapper around a polymorphic generic_file_t that handles its init. This is essentially a factory fo...
Definition file.f90:55
Interface for NEKTON fld files.
Definition fld_file.f90:64
Implements global interpolation for arbitrary points in the domain.
Interpolation between two space::space_t.
Base abstract type for point zones.
The function space for the SEM solution fields.
Definition space.f90:63