Neko-TOP
A portable framework for high-order spectral element flow toplogy optimization.
Loading...
Searching...
No Matches
simulation.f90
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!
34! Here, we simply march forward to steady state solutions
36 use case, only: case_t
37 use neko, only: neko_init, neko_finalize, neko_solve
38 use adjoint_case, only: adjoint_case_t, adjoint_init, adjoint_free
39 use fluid_scheme_incompressible, only: fluid_scheme_incompressible_t
40 use adjoint_fluid_scheme, only: adjoint_fluid_scheme_t
42 use scalar_pnpn, only: scalar_pnpn_t
45 use scalars, only: scalars_t
46 use scalar_scheme, only: scalar_scheme_t
47 use fluid_pnpn, only: fluid_pnpn_t
48 use time_step_controller, only: time_step_controller_t
49 use time_state, only: time_state_t
50 use fld_file_output, only: fld_file_output_t
51 use chkp_output, only: chkp_output_t
52 use simcomp_executor, only: neko_simcomps
53 use neko_ext, only: reset, reset_adjoint
54 use field, only: field_t
55 use field_registry, only: neko_field_registry
56 use field_math, only: field_rzero, field_copy
57 use checkpoint, only: chkp_t
58 use file, only: file_t
59 use utils, only: neko_warning, neko_error
60 use comm, only: pe_rank
61 use json_file_module, only: json_file
62 use json_utils, only: json_get, json_get_or_default
63 use num_types, only: rp, sp, dp
64 use logger, only: log_size, neko_log
65 use mpi_f08, only: mpi_wtime
66 use jobctrl, only: jobctrl_time_limit
67 use profiler, only: profiler_start, profiler_stop, &
68 profiler_start_region, profiler_end_region
69 use simulation_adjoint, only: simulation_adjoint_init, &
70 simulation_adjoint_step, simulation_adjoint_finalize
71 use simulation, only: simulation_init, simulation_step, simulation_finalize, &
72 simulation_restart
73 use simulation_checkpoint, only: simulation_checkpoint_t
74 implicit none
75 private
76
78
80 type(case_t), public :: neko_case
82 type(adjoint_case_t), public :: adjoint_case
84 class(fluid_scheme_incompressible_t), public, pointer :: fluid => null()
86 type(scalars_t), public, pointer :: scalars => null()
88 class(adjoint_fluid_scheme_t), public, pointer :: adjoint_fluid => null()
90 type(adjoint_scalars_t), public, pointer :: adjoint_scalars => null()
93 type(fld_file_output_t), public :: output_forward
96 type(fld_file_output_t), public :: output_adjoint
97
98 logical :: have_scalar = .false.
99 integer :: n_timesteps = 0
100
101 ! ----------------------------------------------------------------------- !
102 ! Checkpoint system
103
105 type(simulation_checkpoint_t) :: checkpoint
106
107 contains
109 procedure, pass(this) :: init => simulation_initialize
111 procedure, pass(this) :: free => simulation_free
113 procedure, pass(this) :: run_forward => simulation_run_forward
115 procedure, pass(this) :: run_backward => simulation_run_backward
117 procedure, pass(this) :: reset => simulation_reset
119 procedure, pass(this) :: write => simulation_write
120
121 end type simulation_t
122 public :: simulation_t
123contains
124
126 subroutine simulation_initialize(this, parameters)
127 class(simulation_t), intent(inout), target :: this
128 type(json_file), intent(inout) :: parameters
129 type(json_file) :: checkpoint_params
130 integer :: i, n_scalars
131
132 ! initialize the primal
133 call neko_init(this%neko_case)
134 ! initialize the adjoint
135 call adjoint_init(this%adjoint_case, this%neko_case)
136
137 ! Start the profiler
138 call profiler_start
139
140 select type (fluid => this%neko_case%fluid)
141 type is (fluid_pnpn_t)
142 this%fluid => fluid
143 end select
144
145 select type (adjoint_fluid => this%adjoint_case%fluid_adj)
146 type is (adjoint_fluid_pnpn_t)
147 this%adjoint_fluid => adjoint_fluid
148 end select
149
150 if (allocated(this%neko_case%scalars)) then
151 this%scalars => this%neko_case%scalars
152 end if
153
154 if (allocated(this%adjoint_case%adjoint_scalars)) then
155 this%adjoint_scalars => this%adjoint_case%adjoint_scalars
156 end if
157
158 ! init the sampler
159 !---------------------------------------------------------
160 ! Allocate the output type
161 n_scalars = 0
162 if (allocated(this%neko_case%scalars)) then
163 n_scalars = size(this%neko_case%scalars%scalar_fields)
164 end if
165 call this%output_forward%init(sp, 'forward_fields', 4 + n_scalars)
166
167 call this%output_forward%fields%assign(1, this%fluid%p)
168 call this%output_forward%fields%assign(2, this%fluid%u)
169 call this%output_forward%fields%assign(3, this%fluid%v)
170 call this%output_forward%fields%assign(4, this%fluid%w)
171
172 ! Assign all scalar fields
173 if (allocated(this%neko_case%scalars)) then
174 do i = 1, n_scalars
175 call this%output_forward%fields%assign(4 + i, &
176 this%scalars%scalar_fields(i)%s)
177 end do
178 end if
179
180 n_scalars = 0
181 if (allocated(this%adjoint_case%adjoint_scalars)) then
182 n_scalars = size(this%adjoint_case%adjoint_scalars%adjoint_scalar_fields)
183 end if
184 call this%output_adjoint%init(sp, 'adjoint_fields', 4 + n_scalars)
185 call this%output_adjoint%fields%assign(1, this%adjoint_fluid%p_adj)
186 call this%output_adjoint%fields%assign(2, this%adjoint_fluid%u_adj)
187 call this%output_adjoint%fields%assign(3, this%adjoint_fluid%v_adj)
188 call this%output_adjoint%fields%assign(4, this%adjoint_fluid%w_adj)
189
190 ! Assign all scalar fields
191 if (allocated(this%adjoint_case%adjoint_scalars)) then
192 do i = 1, n_scalars
193 call this%output_adjoint%fields%assign(4 + i, &
194 this%adjoint_scalars%adjoint_scalar_fields(i)%s_adj)
195 end do
196 end if
197
198 if ("checkpoints" .in. parameters) then
199 call json_get(parameters, 'checkpoints', checkpoint_params)
200 call this%checkpoint%init(this%neko_case, checkpoint_params)
201 end if
202
203 end subroutine simulation_initialize
204
206 subroutine simulation_free(this)
207 class(simulation_t), intent(inout) :: this
208
209 ! Stop the profiler
210 call profiler_stop
211
212 call this%checkpoint%free()
213 call adjoint_free(this%adjoint_case)
214 call neko_finalize(this%neko_case)
215
216 end subroutine simulation_free
217
219 subroutine simulation_run_forward(this)
220 class(simulation_t), intent(inout) :: this
221 type(time_step_controller_t) :: dt_controller
222 real(kind=dp) :: loop_start
223
224 call dt_controller%init(this%neko_case%params)
225
226 call simulation_init(this%neko_case, dt_controller)
227
228 call profiler_start_region("Forward simulation")
229 loop_start = mpi_wtime()
230 this%n_timesteps = 0
231 do while (this%neko_case%time%t .lt. this%neko_case%time%end_time)
232 this%n_timesteps = this%n_timesteps + 1
233
234 call simulation_step(this%neko_case, dt_controller, loop_start)
235
236 call this%checkpoint%save(this%neko_case)
237 end do
238 call profiler_end_region("Forward simulation")
239
240 call simulation_finalize(this%neko_case)
241
242 end subroutine simulation_run_forward
243
245 subroutine simulation_run_backward(this)
246 class(simulation_t), intent(inout) :: this
247 type(time_step_controller_t) :: dt_controller
248 real(kind=dp) :: loop_start
249 real(kind=rp) :: cfl
250 integer :: i
251
252 call dt_controller%init(this%neko_case%params)
253
254 call simulation_adjoint_init(this%adjoint_case, dt_controller)
255
256 call profiler_start_region("Adjoint simulation")
257 cfl = this%adjoint_case%fluid_adj%compute_cfl(this%adjoint_case%time%dt)
258 loop_start = mpi_wtime()
259 do i = this%n_timesteps, 1, -1
260 call this%checkpoint%restore(this%neko_case, i)
261
262 call simulation_adjoint_step(this%adjoint_case, dt_controller, cfl, &
263 loop_start)
264 end do
265 call profiler_end_region("Adjoint simulation")
266
267 call simulation_adjoint_finalize(this%adjoint_case)
268
269 end subroutine simulation_run_backward
270
272 subroutine simulation_reset(this)
273 class(simulation_t), intent(inout) :: this
274
275 call reset(this%neko_case)
276 call reset_adjoint(this%adjoint_case, this%neko_case)
277 call this%checkpoint%reset()
278
279 end subroutine simulation_reset
280
282 subroutine simulation_write(this, idx)
283 class(simulation_t), intent(inout) :: this
284 integer, intent(in) :: idx
285
286 call this%output_forward%sample(real(idx, kind=rp))
287 call this%output_adjoint%sample(real(idx, kind=rp))
288
289 end subroutine simulation_write
290
291end module simulation_m
Adjoint Pn/Pn formulation.
Contains the adjoint_scalar_pnpn_t type.
Contains the adjoint_scalars_t type that manages multiple scalar fields.
Contains extensions to the neko library required to run the topology optimization code.
Definition neko_ext.f90:9
subroutine, public reset(neko_case)
Reset the case data structure.
Definition neko_ext.f90:56
subroutine, public reset_adjoint(adjoint_case, neko_case)
Reset the adjoint case data structure.
Definition neko_ext.f90:194
Adjoint simulation driver.
Implements the steady_problem_t type.
Adjoint case type. Todo: This should Ideally be a subclass of case_t, however, that is not yet suppor...
Type to manage multiple adjoint scalar transport equations.