35 use gather_scatter,
only : gs_t, gs_op_add
36 use neko_config,
only : neko_bcknd_device
37 use num_types,
only : rp
38 use device_math,
only : device_col2
39 use device,
only : device_memcpy, host_to_device
40 use field,
only : field_t
41 use utils,
only : neko_error, filename_chsuffix, filename_suffix, &
42 neko_warning, neko_fname_len, extract_fld_file_index
43 use coefs,
only : coef_t
44 use math,
only : col2, cfill, cfill_mask
45 use json_module,
only : json_file
46 use json_utils,
only: json_get, json_get_or_default
47 use point_zone,
only: point_zone_t
48 use point_zone_registry,
only: neko_point_zone_registry
49 use field_registry,
only: neko_field_registry
50 use logger,
only: neko_log, log_size
51 use fld_file_data,
only: fld_file_data_t
52 use fld_file,
only: fld_file_t
53 use checkpoint,
only: chkp_t
54 use file,
only: file_t
55 use global_interpolation,
only: global_interpolation_t
56 use interpolation,
only: interpolator_t
57 use space,
only: space_t, gll
62 module procedure set_optimization_ic_int
83 subroutine set_optimization_ic_int(fld, coef, gs, params)
84 type(field_t),
intent(inout) :: fld
85 type(coef_t),
intent(in) :: coef
86 type(gs_t),
intent(inout) :: gs
87 type(json_file),
intent(inout) :: params
90 character(len=:),
allocatable :: type
91 real(kind=rp) :: ic_value
92 character(len=:),
allocatable :: read_str
93 character(len=:),
allocatable :: fname, mesh_fname
94 real(kind=rp) :: zone_value, tol
95 logical :: interpolate
97 call json_get(params,
"type", type)
99 select case (trim(type))
101 call json_get(params,
'value', ic_value)
102 call set_optimization_ic_uniform(fld, ic_value)
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)
109 call set_optimization_ic_point_zone(fld, ic_value, read_str, zone_value)
112 call json_get(params,
'file_name', fname)
113 call json_get_or_default(params,
'interpolate', interpolate, .false.)
114 call json_get_or_default(params,
'tolerance', tol, 0.000001_rp)
115 call json_get_or_default(params,
'mesh_file_name', mesh_fname,
"none")
117 call set_optimization_ic_fld(fld, fname, interpolate, tol, mesh_fname)
120 call neko_error(
'Invalid initial condition')
124 call set_optimization_ic_common(fld, coef, gs)
126 end subroutine set_optimization_ic_int
134 subroutine set_optimization_ic_common(fld, coef, gs)
135 type(field_t),
intent(inout) :: fld
136 type(coef_t),
intent(in) :: coef
137 type(gs_t),
intent(inout) :: gs
141 if (neko_bcknd_device .eq. 1)
then
142 call device_memcpy(fld%x, fld%x_d, n, &
143 host_to_device, sync = .false.)
147 call gs%op(fld%x, n, gs_op_add)
149 if (neko_bcknd_device .eq. 1)
then
150 call device_col2(fld%x_d, coef%mult_d, n)
152 call col2(fld%x, coef%mult, n)
155 end subroutine set_optimization_ic_common
161 subroutine set_optimization_ic_uniform(fld, ic_value)
162 type(field_t),
intent(inout) :: fld
163 real(kind=rp),
intent(in) :: ic_value
165 character(len=LOG_SIZE) :: log_buf
167 call neko_log%message(
"Type : uniform")
168 write (log_buf,
'(A,ES12.6)')
"Value: ", ic_value
169 call neko_log%message(log_buf)
173 if (neko_bcknd_device .eq. 1)
then
174 call cfill(fld%x, ic_value, n)
177 end subroutine set_optimization_ic_uniform
186 subroutine set_optimization_ic_point_zone(fld, base_value, zone_name, &
188 type(field_t),
intent(inout) :: fld
189 real(kind=rp),
intent(in) :: base_value
190 character(len=*),
intent(in) :: zone_name
191 real(kind=rp),
intent(in) :: zone_value
194 character(len=LOG_SIZE) :: log_buf
195 class(point_zone_t),
pointer :: zone
198 call neko_log%message(
"Type : point_zone")
199 write (log_buf,
'(A,ES12.6)')
"Base value: ", base_value
200 call neko_log%message(log_buf)
201 call neko_log%message(
"Zone name : " // trim(zone_name))
202 write (log_buf,
'(A,ES12.6)')
"Zone value: ", zone_value
203 call neko_log%message(log_buf)
205 size = fld%dof%size()
206 zone => neko_point_zone_registry%get_point_zone(trim(zone_name))
208 call set_optimization_ic_uniform(fld, base_value)
209 call cfill_mask(fld%x, zone_value,
size, zone%mask%get(), zone%size)
211 end subroutine set_optimization_ic_point_zone
227 subroutine set_optimization_ic_fld(fld, file_name, &
228 interpolate, tolerance, mesh_file_name)
229 type(field_t),
intent(inout) :: fld
230 character(len=*),
intent(in) :: file_name
231 logical,
intent(in) :: interpolate
232 real(kind=rp),
intent(in) :: tolerance
233 character(len=*),
intent(inout) :: mesh_file_name
235 character(len=LOG_SIZE) :: log_buf
236 integer :: sample_idx, sample_mesh_idx
237 type(fld_file_data_t) :: fld_data
239 logical :: mesh_mismatch
242 type(global_interpolation_t) :: global_interp
246 type(space_t) :: prev_Xh
247 type(interpolator_t) :: space_interp
250 call neko_log%message(
"Type : field")
251 call neko_log%message(
"File name : " // trim(file_name))
252 write (log_buf,
'(A,L1)')
"Interpolation: ", interpolate
253 call neko_log%message(log_buf)
256 sample_idx = extract_fld_file_index(file_name, -1)
258 if (sample_idx .eq. -1)
then
259 call neko_error(
"Invalid file name for the initial condition. The " // &
260 "file format must be e.g. 'mean0.f00001'")
264 call filename_chsuffix(file_name, trim(file_name),
'fld')
267 call f%init(trim(file_name))
269 if (interpolate)
then
272 if (mesh_file_name .eq.
"none")
then
273 mesh_file_name = trim(file_name)
274 sample_mesh_idx = sample_idx
278 sample_mesh_idx = extract_fld_file_index(mesh_file_name, -1)
280 if (sample_mesh_idx .eq. -1)
then
281 call neko_error(
"Invalid file name for the initial condition." // &
282 " The file format must be e.g. 'mean0.f00001'")
285 write (log_buf,
'(A,ES12.6)')
"Tolerance : ", tolerance
286 call neko_log%message(log_buf)
287 write (log_buf,
'(A,A)')
"Mesh file : ", &
289 call neko_log%message(log_buf)
294 if (sample_mesh_idx .ne. sample_idx)
then
295 call f%set_counter(sample_mesh_idx)
296 call f%read(fld_data)
302 call f%set_counter(sample_idx)
303 call f%read(fld_data)
311 mesh_mismatch = (fld_data%glb_nelv .ne. fld%msh%glb_nelv .or. &
312 fld_data%gdim .ne. fld%msh%gdim)
314 if (mesh_mismatch .and. .not. interpolate)
then
315 call neko_error(
"The fld file must match the current mesh! " // &
316 "Use 'interpolate': 'true' to enable interpolation.")
317 else if (.not. mesh_mismatch .and. interpolate)
then
318 call neko_log%warning(
"You have activated interpolation but you " // &
319 "might still be using the same mesh.")
324 if (interpolate)
then
326 select type (ft => f%file_type)
328 if (.not. ft%dp_precision)
then
329 call neko_warning(
"The coordinates read from the field file " // &
330 "are in single precision.")
331 call neko_log%message(
"It is recommended to use a mesh in " // &
332 "double precision for better interpolation results.")
333 call neko_log%message(
"If the interpolation does not work, " // &
334 "you can try to increase the tolerance.")
340 call fld_data%generate_interpolator(global_interp, fld%dof, fld%msh, &
346 call global_interp%evaluate(fld%x, fld_data%u%x)
347 call global_interp%free()
352 call prev_xh%init(gll, fld_data%lx, fld_data%ly, fld_data%lz)
353 call space_interp%init(fld%Xh, prev_xh)
356 call space_interp%map_host(fld%x, fld_data%t%x, fld_data%nelv, fld%Xh)
358 call space_interp%free()
364 end subroutine set_optimization_ic_fld
Optimization initial condition.