Neko  0.9.0
A portable framework for high-order spectral element flow simulations
advection_fctry.f90
Go to the documentation of this file.
1 ! Copyright (c) 2024, 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 !
34 submodule(advection) advection_fctry
36 
37  ! Advection and derivatives
38  use adv_dealias, only : adv_dealias_t
40  use adv_oifs, only : adv_oifs_t
41 
42 
43 contains
44 
56  module subroutine advection_factory(object, json, coef, ulag, vlag, wlag, &
57  dtlag, tlag, time_scheme, slag)
58  class(advection_t), allocatable, intent(inout) :: object
59  type(json_file), intent(inout) :: json
60  type(coef_t), intent(inout), target :: coef
61  type(field_series_t), intent(in), target :: ulag, vlag, wlag
62  real(kind=rp), intent(in), target :: dtlag(10)
63  real(kind=rp), intent(in), target :: tlag(10)
64  type(time_scheme_controller_t), intent(in), target :: time_scheme
65  type(field_series_t), target, optional :: slag
66 
67  logical :: dealias, oifs
68  real(kind=rp) :: ctarget
69  integer :: lxd, order
70 
71  ! Read the parameters from the json file
72  call json_get(json, 'case.numerics.dealias', dealias)
73  call json_get(json, 'case.numerics.polynomial_order', order)
74  call json_get_or_default(json, 'case.numerics.oifs', oifs, .false.)
75 
76  call json_get_or_default(json, 'case.numerics.dealiased_polynomial_order', &
77  lxd, ( 3 * (order + 1) ) / 2)
78 
79  call json_get_or_default(json, 'case.numerics.target_cfl', ctarget, 1.9_rp)
80 
81  ! Free allocatables if necessary
82  if (allocated(object)) then
83  call object%free
84  deallocate(object)
85  end if
86 
87  if (oifs) then
88  allocate(adv_oifs_t::object)
89  else
90  if (dealias) then
91  allocate(adv_dealias_t::object)
92  else
93  allocate(adv_no_dealias_t::object)
94  end if
95  end if
96 
97  select type (adv => object)
98  type is (adv_dealias_t)
99  call adv%init(lxd, coef)
100  type is (adv_no_dealias_t)
101  call adv%init(coef)
102  type is (adv_oifs_t)
103  if (present(slag)) then
104  call adv%init(lxd, coef, ctarget, ulag, vlag, wlag, &
105  dtlag, tlag, time_scheme, slag)
106  else
107  call adv%init(lxd, coef, ctarget, ulag, vlag, wlag, &
108  dtlag, tlag, time_scheme)
109  end if
110  end select
111 
112  end subroutine advection_factory
113 
114 
115 end submodule advection_fctry
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
Definition: json_utils.f90:54
Retrieves a parameter by name or throws an error.
Definition: json_utils.f90:45
Subroutines to add advection terms to the RHS of a transport equation.
Definition: adv_dealias.f90:34
Subroutines to add advection terms to the RHS of a transport equation.
Subroutines to add advection terms to the RHS of a transport equation.
Definition: adv_oifs.f90:34
Subroutines to add advection terms to the RHS of a transport equation.
Definition: advection.f90:34
Utilities for retrieving parameters from the case files.
Definition: json_utils.f90:34
Base class for time integration schemes.
Definition: time_scheme.f90:61
Type encapsulating advection routines with dealiasing.
Definition: adv_dealias.f90:52
Type encapsulating advection routines with no dealiasing applied.