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 fluid_pnpn, only: fluid_pnpn_t
47 use time_step_controller, only: time_step_controller_t
48 use fld_file_output, only: fld_file_output_t
49 use simcomp_executor, only: neko_simcomps
50 use neko_ext, only: reset
51 use field_math, only: field_rzero
52 use json_file_module, only: json_file
53 use json_utils, only: json_extract_item
54 use num_types, only: rp, sp, dp
55 use logger, only: log_size
56 use mpi_f08, only: mpi_wtime
57 use jobctrl, only: jobctrl_time_limit
58 use profiler, only: profiler_start, profiler_stop
59 use simulation_adjoint, only: simulation_adjoint_init, &
60 simulation_adjoint_step, simulation_adjoint_finalize
61 implicit none
62 private
63
66 type(case_t), public :: neko_case
68 type(adjoint_case_t), public :: adjoint_case
70 class(fluid_scheme_incompressible_t), public, pointer :: fluid => null()
72 type(scalars_t), public, pointer :: scalars => null()
74 class(adjoint_fluid_scheme_t), public, pointer :: adjoint_fluid => null()
76 type(adjoint_scalars_t), public, pointer :: adjoint_scalars => null()
79 type(fld_file_output_t), public :: output_forward
82 type(fld_file_output_t), public :: output_adjoint
83
84 contains
86 procedure, pass(this) :: init => simulation_init
88 procedure, pass(this) :: free => simulation_free
90 procedure, pass(this) :: run_forward => simulation_run_forward
92 procedure, pass(this) :: run_backward => simulation_run_backward
94 procedure, pass(this) :: reset => simulation_reset
96 procedure, pass(this) :: write => simulation_write
97 end type simulation_t
98
99 public :: simulation_t
100contains
101
103 subroutine simulation_init(this, parameters)
104 class(simulation_t), intent(inout), target :: this
105 type(json_file), intent(inout) :: parameters
106 integer :: i, n_scalars
107
108 ! initialize the primal
109 call neko_init(this%neko_case)
110 ! initialize the adjoint
111 call adjoint_init(this%adjoint_case, this%neko_case)
112
113 select type (fluid => this%neko_case%fluid)
114 type is (fluid_pnpn_t)
115 this%fluid => fluid
116
117 end select
118
119 select type (adjoint_fluid => this%adjoint_case%fluid_adj)
120 type is (adjoint_fluid_pnpn_t)
121 this%adjoint_fluid => adjoint_fluid
122
123 end select
124
125 if (allocated(this%neko_case%scalars)) then
126 this%scalars => this%neko_case%scalars
127 end if
128
129 if (allocated(this%adjoint_case%adjoint_scalars)) then
130 this%adjoint_scalars => this%adjoint_case%adjoint_scalars
131 end if
132
133 ! init the sampler
134 !---------------------------------------------------------
135 ! Allocate the output type
136 n_scalars = 0
137 if (allocated(this%neko_case%scalars)) then
138 n_scalars = size(this%neko_case%scalars%scalar_fields)
139 end if
140 call this%output_forward%init(sp, 'forward_fields', 4 + n_scalars)
141
142 call this%output_forward%fields%assign(1, this%fluid%p)
143 call this%output_forward%fields%assign(2, this%fluid%u)
144 call this%output_forward%fields%assign(3, this%fluid%v)
145 call this%output_forward%fields%assign(4, this%fluid%w)
146
147 ! Assign all scalar fields
148 if (allocated(this%neko_case%scalars)) then
149 do i = 1, n_scalars
150 call this%output_forward%fields%assign(4 + i, &
151 this%scalars%scalar_fields(i)%s)
152 end do
153 end if
154
155 n_scalars = 0
156 if (allocated(this%adjoint_case%adjoint_scalars)) then
157 n_scalars = size(this%adjoint_case%adjoint_scalars%adjoint_scalar_fields)
158 end if
159 call this%output_adjoint%init(sp, 'adjoint_fields', 4 + n_scalars)
160 call this%output_adjoint%fields%assign(1, this%adjoint_fluid%p_adj)
161 call this%output_adjoint%fields%assign(2, this%adjoint_fluid%u_adj)
162 call this%output_adjoint%fields%assign(3, this%adjoint_fluid%v_adj)
163 call this%output_adjoint%fields%assign(4, this%adjoint_fluid%w_adj)
164
165 ! Assign all scalar fields
166 if (allocated(this%adjoint_case%adjoint_scalars)) then
167 do i = 1, n_scalars
168 call this%output_adjoint%fields%assign(4 + i, &
169 this%adjoint_scalars%adjoint_scalar_fields(i)%s_adj)
170 end do
171 end if
172
173 end subroutine simulation_init
174
176 subroutine simulation_free(this)
177 class(simulation_t), intent(inout) :: this
178
179 call adjoint_free(this%adjoint_case)
180 call neko_finalize(this%neko_case)
181
182 end subroutine simulation_free
183
185 subroutine simulation_run_forward(this)
186 class(simulation_t), intent(inout) :: this
187
188 ! set forward time to zero
189 this%neko_case%time%t = 0.0_rp
190 this%neko_case%time%tstep = 0
191
192 ! run the primal
193 call neko_solve(this%neko_case)
194
195 end subroutine simulation_run_forward
196
198 subroutine simulation_run_backward(this)
199 class(simulation_t), intent(inout) :: this
200 type(time_step_controller_t) :: dt_controller
201 real(kind=dp) :: tstep_loop_start_time
202 real(kind=rp) :: cfl
203
204 call dt_controller%init(this%neko_case%params)
205
206 call simulation_adjoint_init(this%adjoint_case, dt_controller)
207
208 call profiler_start
209 cfl = this%adjoint_case%fluid_adj%compute_cfl(this%adjoint_case%time%dt)
210 tstep_loop_start_time = mpi_wtime()
211
212 do while (this%adjoint_case%time%t .lt. this%adjoint_case%time%end_time &
213 .and. (.not. jobctrl_time_limit()))
214 call simulation_adjoint_step(this%adjoint_case, dt_controller, cfl, &
215 tstep_loop_start_time)
216 end do
217 call profiler_stop
218
219 call simulation_adjoint_finalize(this%adjoint_case)
220
221 end subroutine simulation_run_backward
222
224 subroutine simulation_reset(this)
225 class(simulation_t), intent(inout) :: this
226 integer :: i, n_scalars
227
228 call reset(this%neko_case)
229
230 ! TODO
231 ! reset for the adjoint
232 ! call reset(this%adjoint_case)
233 this%adjoint_case%time%t = 0.0_dp
234 this%adjoint_case%time%tstep = 0
235
236 call field_rzero(this%adjoint_case%fluid_adj%u_adj)
237 call field_rzero(this%adjoint_case%fluid_adj%v_adj)
238 call field_rzero(this%adjoint_case%fluid_adj%w_adj)
239 if (allocated(this%adjoint_case%adjoint_scalars)) then
240 n_scalars = size(this%adjoint_case%adjoint_scalars%adjoint_scalar_fields)
241 do i = 1, n_scalars
242 call field_rzero(&
243 this%adjoint_case%adjoint_scalars%adjoint_scalar_fields(i)%s_adj)
244 end do
245 end if
246
247 end subroutine simulation_reset
248
250 subroutine simulation_write(this, idx)
251 class(simulation_t), intent(inout) :: this
252 integer, intent(in) :: idx
253
254 call this%output_forward%sample(real(idx, kind=rp))
255 call this%output_adjoint%sample(real(idx, kind=rp))
256
257 end subroutine simulation_write
258
259end 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:55
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.