Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
les_simcomp.f90
Go to the documentation of this file.
1! Copyright (c) 2023, 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!
33!
35
37 use num_types, only : rp
38 use json_module, only : json_file
40 use case, only : case_t
41 use time_state, only : time_state_t
42 use les_model, only : les_model_t, les_model_factory
45 use utils, only : neko_error
46 implicit none
47 private
48
51 type, public, extends(simulation_component_t) :: les_simcomp_t
53 class(les_model_t), allocatable :: les_model
55 type(field_writer_t) :: writer
56 contains
58 procedure, pass(this) :: init => les_simcomp_init_from_json
60 procedure, pass(this) :: free => les_simcomp_free
62 procedure, pass(this) :: preprocess_ => les_simcomp_compute
64 procedure, pass(this) :: restart_ => les_simcomp_restart
65 end type les_simcomp_t
66
67contains
68
70 subroutine les_simcomp_init_from_json(this, json, case)
71 class(les_simcomp_t), intent(inout), target :: this
72 type(json_file), intent(inout) :: json
73 class(case_t), intent(inout), target :: case
74 character(len=:), allocatable :: name
75 character(len=:), allocatable :: nut_field
76 character(len=20) :: fields(2)
77
78 call this%free()
79
80 ! Add fields keyword to the json so that the field_writer picks it up.
81 ! Will also add fields to the registry if missing.
82 call json_get_or_default(json, "nut_field", nut_field, "nut")
83 fields(1) = "les_delta"
84 fields(2) = nut_field
85
86 call json%add("fields", fields)
87
88 call this%init_base(json, case)
89 call this%writer%init(json, case)
90
91 call json_get(json, "model", name)
92
93 call les_model_factory(this%les_model, name, case%fluid, json)
94 end subroutine les_simcomp_init_from_json
95
97 subroutine les_simcomp_free(this)
98 class(les_simcomp_t), intent(inout) :: this
99 call this%free_base()
100 call this%writer%free()
101
102 if (allocated(this%les_model)) then
103 call this%les_model%free()
104 deallocate(this%les_model)
105 end if
106 end subroutine les_simcomp_free
107
110 subroutine les_simcomp_compute(this, time)
111 class(les_simcomp_t), intent(inout) :: this
112 type(time_state_t), intent(in) :: time
113
114 call this%les_model%compute(time%t, time%tstep)
115 end subroutine les_simcomp_compute
116
119 subroutine les_simcomp_restart(this, time)
120 class(les_simcomp_t), intent(inout) :: this
121 type(time_state_t), intent(in) :: time
122
123 call this%les_model%compute(time%t, 0)
124 end subroutine les_simcomp_restart
125
126end module les_simcomp
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.
Defines a simulation case.
Definition case.f90:34
Implements the field_writer_t type.
Utilities for retrieving parameters from the case files.
Implements les_model_t.
Definition les_model.f90:35
Implements the les_simcomp_t type.
subroutine les_simcomp_restart(this, time)
Compute the les_simcomp field when restart.
subroutine les_simcomp_init_from_json(this, json, case)
Constructor from json.
subroutine les_simcomp_compute(this, time)
Compute the les_simcomp field.
subroutine les_simcomp_free(this)
Destructor.
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Simulation components are objects that encapsulate functionality that can be fit to a particular comp...
subroutine restart_(this, time)
Dummy restart function.
subroutine preprocess_(this, time)
Dummy preprocessing function.
Module with things related to the simulation time.
Utilities.
Definition utils.f90:35
A simulation component that writes a 3d field to a file.
Base abstract type for LES models based on the Boussinesq approximation.
Definition les_model.f90:64
A simulation component that drives the computation of the SGS viscosity.
Base abstract class for simulation components.
A struct that contains all info about the time, expand as needed.