47 use json_module,
only : json_file, json_core, json_value
79 type(
field_t),
pointer :: indicator_unfiltered
88 procedure,
public, pass(this) :: init => &
115 type(json_file),
intent(inout) :: json
117 type(
coef_t),
intent(in),
target :: coef
118 character(len=*),
intent(in) :: variable_name
119 real(kind=
rp) :: start_time, end_time
121 character(len=:),
allocatable :: filter_type
122 real(kind=
rp) :: filter_radius
123 real(kind=
rp),
dimension(:),
allocatable :: brinkman_limits
124 real(kind=
rp) :: brinkman_penalty
126 type(json_value),
pointer :: json_object_list
127 type(json_core) :: core
129 character(len=:),
allocatable :: object_type
130 type(json_file) :: object_settings
141 call json_get(json,
'brinkman.limits', brinkman_limits)
142 call json_get(json,
'brinkman.penalty', brinkman_penalty)
144 if (
size(brinkman_limits) .ne. 2)
then
145 call neko_error(
'brinkman_limits must be a 2 element array of reals')
149 call this%init_base(fields, coef, start_time, end_time)
161 this%indicator_unfiltered => &
168 call json%info(
'objects', n_children = n_regions)
174 select case (object_type)
175 case (
'boundary_mesh')
176 call this%init_boundary_mesh(object_settings)
178 call this%init_point_zone(object_settings)
181 call object_settings%print()
182 call neko_error(
'Brinkman source term objects require a region type')
184 call neko_error(
'Brinkman source term unknown region type')
193 select case (filter_type)
196 call this%indicator_unfiltered%init(coef%dof)
202 call this%filter%init(json, coef)
205 call field_copy(this%indicator_unfiltered, this%indicator)
208 call this%filter%apply(this%indicator, this%indicator_unfiltered)
212 call output%fields%assign_to_field(1, this%indicator_unfiltered)
213 call output%fields%assign_to_field(2, this%indicator)
214 call output%fields%assign_to_field(3, this%brinkman)
219 call output%fields%assign_to_field(1, this%indicator)
220 call output%fields%assign_to_field(2, this%brinkman)
223 call neko_error(
'Brinkman source term unknown filter type')
229 this%brinkman = this%indicator
231 brinkman_limits(1), brinkman_limits(2), brinkman_penalty)
234 call output%sample(0.0_rp)
242 nullify(this%indicator)
243 nullify(this%indicator_unfiltered)
244 nullify(this%brinkman)
246 if (
allocated(this%filter))
then
247 call this%filter%free()
248 deallocate(this%filter)
250 call this%free_base()
259 type(
field_t),
pointer :: u, v, w, fu, fv, fw
262 n = this%fields%item_size(1)
268 fu => this%fields%get(1)
269 fv => this%fields%get(2)
270 fw => this%fields%get(3)
284 type(json_file),
intent(inout) :: json
287 character(len=:),
allocatable :: mesh_file_name
288 character(len=:),
allocatable :: distance_transform
289 character(len=:),
allocatable :: filter_type
290 character(len=:),
allocatable :: mesh_transform
295 real(kind=
rp) :: scalar_r
296 real(kind=
dp) :: scalar_d
297 logical :: cache, cache_exist
298 character(len=:),
allocatable :: cache_filename
299 type(
file_t) :: cache_file
307 real(kind=
dp),
dimension(:),
allocatable :: box_min, box_max
308 logical :: keep_aspect_ratio
309 real(kind=
dp),
dimension(3) :: scaling
310 real(kind=
dp),
dimension(3) :: translation
312 type(
aabb_t) :: mesh_box, target_box
314 character(len=LOG_SIZE) :: log_msg
319 call json_get(json,
'name', mesh_file_name)
323 call json_get(json,
'distance_transform.type', distance_transform)
328 call json_get(json,
'cache_file', cache_filename)
330 inquire(
file=trim(cache_filename) //
"0.nek5000", exist=cache_exist)
331 write(log_msg,
'(A)')
"Checking for Brinkman source term cache."
334 if (cache_exist)
then
335 write(log_msg,
'(A)')
"Loading Brinkman source term from cache."
338 call cache_data%init()
339 call temp_field%init(this%coef%dof)
341 call cache_file%init(cache_filename //
"0.fld")
342 call cache_file%set_counter(0)
343 call cache_file%read(cache_data)
351 if (cache_data%glb_nelv .ne. temp_field%msh%glb_nelv .or. &
352 cache_data%gdim .ne. temp_field%msh%gdim)
then
353 call neko_error(
"The fld file must match the current mesh! " // &
354 "Use 'interpolate': 'true' to enable interpolation.")
358 call prev_xh%init(
gll, cache_data%lx, cache_data%ly, cache_data%lz)
359 call space_interp%init(temp_field%Xh, prev_xh)
360 call space_interp%map_host(temp_field%x, cache_data%p%x, &
361 cache_data%nelv, temp_field%Xh)
362 call space_interp%free()
375 call cache_data%free()
376 call temp_field%free()
377 call cache_file%free()
385 call mesh_file%init(mesh_file_name)
386 call mesh_file%read(boundary_mesh)
388 if (boundary_mesh%nelv .eq. 0)
then
389 call neko_error(
'No elements in the boundary mesh')
396 mesh_transform,
'none')
398 select case (mesh_transform)
401 case (
'bounding_box')
402 call json_get(json,
'mesh_transform.box_min', box_min)
403 call json_get(json,
'mesh_transform.box_max', box_max)
405 keep_aspect_ratio, .true.)
407 if (
size(box_min) .ne. 3 .or.
size(box_max) .ne. 3)
then
409 &box_min and box_max must be 3 element arrays of reals')
412 call target_box%init(box_min, box_max)
416 scaling = target_box%get_diagonal() / mesh_box%get_diagonal()
417 if (keep_aspect_ratio)
then
418 scaling = minval(scaling)
421 translation = - scaling * mesh_box%get_min() + target_box%get_min()
423 do idx_p = 1, boundary_mesh%mpts
424 boundary_mesh%points(idx_p)%x = &
425 scaling * boundary_mesh%points(idx_p)%x + translation
429 write(log_msg,
'(A)')
"The following transformation was applied:"
431 write(log_msg,
'(A, 3F12.6)')
"Scaling: ", scaling
433 write(log_msg,
'(A, 3F12.6)')
"Translation: ", translation
448 call temp_field%init(this%coef%dof)
451 select case (distance_transform)
453 call json_get(json,
'distance_transform.value', scalar_d)
454 scalar_r =
real(scalar_d, kind=
rp)
461 call json_get(json,
'distance_transform.value', scalar_d)
462 scalar_r =
real(scalar_d, kind=
rp)
473 write(log_msg,
'(A)')
"Writing Brinkman source term to cache."
475 call cache_output%init(
dp, cache_filename, 1)
476 call cache_output%fields%assign_to_field(1, temp_field)
477 call cache_output%sample(0.0_rp)
483 call temp_field%free()
490 type(json_file),
intent(inout) :: json
493 character(len=:),
allocatable :: zone_name
502 call json_get(json,
'name', zone_name)
505 call temp_field%init(this%coef%dof)
511 zone%mask%get_d(), zone%size)
513 call cfill_mask(temp_field%x, 1.0_rp, temp_field%size(), &
514 zone%mask%get(), zone%size)
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.
Axis Aligned Bounding Box (aabb) implementation in Fortran.
type(aabb_t) function, public get_aabb(object, padding)
Construct the aabb of a predefined object.
Implements the brinkman_source_term_t type.
subroutine init_point_zone(this, json)
Initializes the source term from a point zone.
subroutine init_boundary_mesh(this, json)
Initializes the source term from a boundary mesh.
subroutine brinkman_source_term_init_from_json(this, json, fields, coef, variable_name)
The common constructor using a JSON object.
subroutine brinkman_source_term_free(this)
Destructor.
subroutine brinkman_source_term_compute(this, time)
Computes the source term and adds the result to fields.
subroutine, public device_cfill_mask(a_d, c, n, mask_d, n_mask, strm)
Fill a constant to a masked vector. .
subroutine, public device_pwmax2(a_d, b_d, n, strm)
Compute the point-wise maximum of two vectors .
Device abstraction, common interface for various accelerators.
integer, parameter, public host_to_device
subroutine, public field_pwmax2(a, b, n)
Point-wise max operation .
subroutine, public field_copy(a, b, n)
Copy a vector .
subroutine, public field_subcol3(a, b, c, n)
Returns .
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.
Filter to be applied to a scalar field.
Simple module to handle fld file series. Provides an interface to the different fields sotred in a fl...
Implements fld_file_output_t.
Implements global_interpolation given a dofmap.
Routines to interpolate between different spaces.
Utilities for retrieving parameters from the case files.
integer, parameter, public neko_log_debug
Debug log level.
type(log_t), public neko_log
Global log stream.
integer, parameter, public log_size
A module containing mapping functions and subroutines. These functions are used to modify fields in a...
subroutine, public permeability_field(f, k_0, k_1, q)
Apply a permeability function to a field.
subroutine, public smooth_step_field(f, edge0, edge1)
Apply a smooth step function to a field.
subroutine, public step_function_field(f, x0, value0, value1)
Apply a step function to a field.
subroutine, public pwmax2(a, b, n)
Point-wise maximum of two vectors .
subroutine, public cfill_mask(a, c, n, mask, n_mask)
Fill a constant to a masked vector. .
integer, parameter neko_bcknd_device
integer, parameter, public dp
integer, parameter, public sp
integer, parameter, public rp
Global precision used in computations.
type(point_zone_registry_t), target, public neko_point_zone_registry
Global point_zone registry.
subroutine, public profiler_start_region(name, region_id)
Started a named (name) profiler region.
subroutine, public profiler_end_region(name, region_id)
End the most recently started profiler region.
Module containing Signed Distance Functions.
subroutine, public signed_distance_field(field_data, object, max_distance)
Signed distance field.
Implements the source_term_t type and a wrapper source_term_wrapper_t.
Defines a function space.
integer, parameter, public gll
Module with things related to the simulation time.
Defines a triangular surface mesh.
Axis Aligned Bounding Box (aabb) data structure.
A Brinkman source term. The region and strength are controlled by assigning regions types and brinkma...
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
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...
Base abstract class for filter.
A simple output saving a list of fields to a .fld file.
Implements global interpolation for arbitrary points in the domain.
Interpolation between two space::space_t.
A PDE based filter mapping $\rho \mapsto \tilde{\rho}$, see Lazarov & O. Sigmund 2010,...
Base abstract type for point zones.
Base abstract type for source terms.
The function space for the SEM solution fields.
A struct that contains all info about the time, expand as needed.