38 use json_module,
only : json_file
62 real(kind=
rp) :: ref_value = 0
71 procedure, pass(this) :: init_from_compenents => &
86 type(json_file),
intent(inout) :: json
88 type(
coef_t),
intent(inout) :: coef
89 real(kind=
rp),
allocatable :: values(:)
90 real(kind=
rp) :: start_time, end_time, ref_value
91 character(len=:),
allocatable :: scalar_name
92 real(kind=
rp),
allocatable :: g(:)
95 if (.not.
size(fields%fields) == 3)
then
96 call neko_error(
"Boussinesq term expects 3 fields to work on.")
105 if (.not.
size(g) == 3)
then
106 call neko_error(
"The gravity vector should have 3 components")
109 call json_get(json,
"reference_value", ref_value)
113 ref_value, g, beta, coef, start_time, end_time)
127 scalar_name, ref_value, g, beta, coef, start_time, end_time)
130 character(len=*),
intent(in) :: scalar_name
131 real(kind=
rp),
intent(in) :: ref_value
132 real(kind=
rp),
intent(in) :: g(3)
133 real(kind=
rp),
intent(in) :: beta
135 real(kind=
rp),
intent(in) :: start_time
136 real(kind=
rp),
intent(in) :: end_time
139 call this%init_base(fields, coef, start_time, end_time)
146 this%ref_value = ref_value
155 call this%free_base()
164 real(kind=
rp),
intent(in) :: t
165 integer,
intent(in) :: tstep
166 integer :: n_fields, i, n
168 n_fields =
size(this%fields%fields)
169 n = this%fields%fields(1)%f%dof%size()
173 this%ref_value, this%g, this%beta)
176 this%ref_value, this%g, this%beta)
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.
Implements the cpu kernel for the boussinesq_source_term_t type.
subroutine, public boussinesq_source_term_compute_cpu(fields, s, ref_value, g, beta)
Computes the Boussinesq source term on the cpu.
Implements the device kernel for the boussinesq_source_term_t type.
subroutine, public boussinesq_source_term_compute_device(fields, s, ref_value, g, beta)
Computes the Boussinesq source term on the device.
Implements the boussinesq_source_term_t type.
subroutine boussinesq_source_term_compute(this, t, tstep)
Computes the source term and adds the result to fields.
subroutine boussinesq_source_term_init_from_json(this, json, fields, coef)
The common constructor using a JSON object.
subroutine boussinesq_source_term_init_from_components(this, fields, scalar_name, ref_value, g, beta, coef, start_time, end_time)
The constructor from type components.
subroutine boussinesq_source_term_free(this)
Destructor.
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Utilities for retrieving parameters from the case files.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Implements the source_term_t type and a wrapper source_term_wrapper_t.
Bouyancy source term accroding to the Boussinesq approximation.
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
Base abstract type for source terms.