46 use json_module,
only : json_file
82 type(
field_t),
intent(inout) :: s
83 type(
coef_t),
intent(in) :: coef
84 type(
gs_t),
intent(inout) :: gs
85 character(len=*) :: type
86 type(json_file),
intent(inout) :: params
89 real(kind=
rp) :: ic_value
90 character(len=:),
allocatable :: read_str
91 character(len=NEKO_FNAME_LEN) :: fname, mesh_fname
92 real(kind=
rp) :: zone_value, tol
93 logical :: interpolate
95 if (trim(type) .eq.
'uniform')
then
97 call json_get(params,
'case.scalar.initial_condition.value', ic_value)
100 else if (trim(type) .eq.
'point_zone')
then
102 call json_get(params,
'case.scalar.initial_condition.base_value', &
104 call json_get(params,
'case.scalar.initial_condition.zone_name', &
106 call json_get(params,
'case.scalar.initial_condition.zone_value', &
111 else if (trim(type) .eq.
'field')
then
113 call json_get(params,
'case.scalar.initial_condition.file_name', &
115 fname = trim(read_str)
117 'case.scalar.initial_condition.interpolate', interpolate, &
120 'case.scalar.initial_condition.tolerance', tol, 0.000001_rp)
122 'case.scalar.initial_condition.mesh_file_name', read_str, &
124 mesh_fname = trim(read_str)
144 type(
field_t),
intent(inout) :: s
145 type(
coef_t),
intent(in) :: coef
146 type(
gs_t),
intent(inout) :: gs
148 type(json_file),
intent(inout) :: params
151 call usr_ic(s, params)
164 type(
field_t),
intent(inout) :: s
165 type(
coef_t),
intent(in) :: coef
166 type(
gs_t),
intent(inout) :: gs
176 call gs%op(s%x, n, gs_op_add)
181 call col2(s%x, coef%mult, n)
191 type(
field_t),
intent(inout) :: s
192 real(kind=
rp),
intent(in) :: ic_value
194 character(len=LOG_SIZE) :: log_buf
196 call neko_log%message(
"Type : uniform")
197 write (log_buf,
'(A,ES12.6)')
"Value: ", ic_value
203 call cfill(s%x, ic_value, n)
216 type(
field_t),
intent(inout) :: s
217 real(kind=
rp),
intent(in) :: base_value
218 character(len=*),
intent(in) :: zone_name
219 real(kind=
rp),
intent(in) :: zone_value
222 character(len=LOG_SIZE) :: log_buf
226 call neko_log%message(
"Type : point_zone")
227 write (log_buf,
'(A,ES12.6)')
"Base value: ", base_value
229 call neko_log%message(
"Zone name : " // trim(zone_name))
230 write (log_buf,
'(A,ES12.6)')
"Zone value: ", zone_value
237 call cfill_mask(s%x, zone_value,
size, zone%mask, zone%size)
256 interpolate, tolerance, mesh_file_name)
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
263 character(len=LOG_SIZE) :: log_buf
264 integer :: sample_idx, sample_mesh_idx
265 integer :: last_index
268 logical :: mesh_mismatch
279 call neko_log%message(
"Type : field")
280 call neko_log%message(
"File name : " // trim(file_name))
281 write (log_buf,
'(A,L1)')
"Interpolation : ", interpolate
283 if (interpolate)
then
289 if (sample_idx .eq. -1) &
290 call neko_error(
"Invalid file name for the initial condition. The&
291 & file format must be e.g. 'mean0.f00001'")
297 f =
file_t(trim(file_name))
299 if (interpolate)
then
302 if (mesh_file_name .eq.
"none")
then
303 mesh_file_name = trim(file_name)
304 sample_mesh_idx = sample_idx
310 if (sample_mesh_idx .eq. -1) &
311 call neko_error(
"Invalid file name for the initial condition. &
312 &The file format must be e.g. 'mean0.f00001'")
314 write (log_buf,
'(A,ES12.6)')
"Tolerance : ", tolerance
316 write (log_buf,
'(A,A)')
"Mesh file : ", &
323 if (sample_mesh_idx .ne. sample_idx)
then
324 call f%set_counter(sample_mesh_idx)
325 call f%read(fld_data)
331 call f%set_counter(sample_idx)
332 call f%read(fld_data)
340 mesh_mismatch = (fld_data%glb_nelv .ne. s%msh%glb_nelv .or. &
341 fld_data%gdim .ne. s%msh%gdim)
343 if (mesh_mismatch .and. .not. interpolate)
then
344 call neko_error(
"The fld file must match the current mesh! &
345 &Use 'interpolate': 'true' to enable interpolation.")
346 else if (.not. mesh_mismatch .and. interpolate)
then
347 call neko_log%warning(
"You have activated interpolation but you might &
348 &still be using the same mesh.")
353 if (interpolate)
then
355 select type (ft => f%file_type)
357 if (.not. ft%dp_precision)
then
358 call neko_warning(
"The coordinates read from the field file are &
359 &in single precision.")
360 call neko_log%message(
"It is recommended to use a mesh in double &
361 &precision for better interpolation results.")
362 call neko_log%message(
"If the interpolation does not work, you&
363 &can try to increase the tolerance.")
369 global_interp = fld_data%generate_interpolator(s%dof, s%msh, tolerance)
372 call global_interp%evaluate(s%x, fld_data%t%x)
373 call global_interp%free
378 call prev_xh%init(
gll, fld_data%lx, fld_data%ly, fld_data%lz)
379 call space_interp%init(s%Xh, prev_xh)
382 call space_interp%map_host(s%x, fld_data%t%x, fld_data%nelv, s%Xh)
384 call space_interp%free
Copy data between host and device (or device and device)
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 scalar initial conditions.
subroutine, public device_col2(a_d, b_d, n)
Vector multiplication .
Device abstraction, common interface for various accelerators.
integer, parameter, public host_to_device
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Module for file I/O operations.
Simple module to handle fld file series. Provides an interface to the different fields sotred in a fl...
Implements global_interpolation given a dofmap.
Routines to interpolate between different spaces.
Utilities for retrieving parameters from the case files.
type(log_t), public neko_log
Global log stream.
integer, parameter, public log_size
subroutine, public cfill(a, c, n)
Set all elements to a constant c .
subroutine, public col2(a, b, n)
Vector multiplication .
subroutine, public cfill_mask(a, c, size, mask, mask_size)
Fill a constant to a masked vector. .
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
type(point_zone_registry_t), target, public neko_point_zone_registry
Global point_zone registry.
Scalar initial condition.
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_usr(s, coef, gs, usr_ic, params)
Set scalar intial condition (user defined)
subroutine set_scalar_ic_fld(s, file_name, interpolate, tolerance, mesh_file_name)
Set the initial condition of the scalar based on a field. @detail The field is read from an fld file....
subroutine set_scalar_ic_int(s, coef, gs, type, params)
Set scalar initial condition (builtin)
Defines a function space.
integer, parameter, public gll
Interfaces for user interaction with NEKO.
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...
integer, parameter, public neko_fname_len
subroutine, public neko_warning(warning_msg)
Reports a warning to standard output.
subroutine, public filename_chsuffix(fname, new_fname, new_suffix)
Change a filename's suffix.
subroutine, public filename_suffix(fname, suffix)
Extract a filename's suffix.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
A wrapper around a polymorphic generic_file_t that handles its init. This is essentially a factory fo...
Interface for NEKTON fld files.
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.