36 use case,
only: case_t
37 use neko,
only: neko_init, neko_finalize, neko_solve
39 use fluid_scheme_incompressible,
only: fluid_scheme_incompressible_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
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
70 simulation_adjoint_step, simulation_adjoint_finalize
71 use simulation,
only: simulation_init, simulation_step, simulation_finalize, &
80 type(case_t),
public :: neko_case
84 class(fluid_scheme_incompressible_t),
public,
pointer :: fluid => null()
86 type(scalars_t),
public,
pointer :: scalars => null()
93 type(fld_file_output_t),
public :: output_forward
96 type(fld_file_output_t),
public :: output_adjoint
98 logical :: have_scalar = .false.
99 integer :: n_timesteps = 0
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
126 subroutine simulation_initialize(this, parameters)
128 type(json_file),
intent(inout) :: parameters
129 type(json_file) :: checkpoint_params
130 integer :: i, n_scalars
133 call neko_init(this%neko_case)
140 select type (fluid => this%neko_case%fluid)
141 type is (fluid_pnpn_t)
145 select type (adjoint_fluid => this%adjoint_case%fluid_adj)
147 this%adjoint_fluid => adjoint_fluid
150 if (
allocated(this%neko_case%scalars))
then
151 this%scalars => this%neko_case%scalars
154 if (
allocated(this%adjoint_case%adjoint_scalars))
then
155 this%adjoint_scalars => this%adjoint_case%adjoint_scalars
162 if (
allocated(this%neko_case%scalars))
then
163 n_scalars =
size(this%neko_case%scalars%scalar_fields)
165 call this%output_forward%init(sp,
'forward_fields', 4 + n_scalars)
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)
173 if (
allocated(this%neko_case%scalars))
then
175 call this%output_forward%fields%assign(4 + i, &
176 this%scalars%scalar_fields(i)%s)
181 if (
allocated(this%adjoint_case%adjoint_scalars))
then
182 n_scalars =
size(this%adjoint_case%adjoint_scalars%adjoint_scalar_fields)
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)
191 if (
allocated(this%adjoint_case%adjoint_scalars))
then
193 call this%output_adjoint%fields%assign(4 + i, &
194 this%adjoint_scalars%adjoint_scalar_fields(i)%s_adj)
198 if (
"checkpoints" .in. parameters)
then
199 call json_get(parameters,
'checkpoints', checkpoint_params)
200 call this%checkpoint%init(this%neko_case, checkpoint_params)
203 end subroutine simulation_initialize
206 subroutine simulation_free(this)
212 call this%checkpoint%free()
213 call adjoint_free(this%adjoint_case)
214 call neko_finalize(this%neko_case)
216 end subroutine simulation_free
219 subroutine simulation_run_forward(this)
221 type(time_step_controller_t) :: dt_controller
222 real(kind=dp) :: loop_start
224 call dt_controller%init(this%neko_case%params)
226 call simulation_init(this%neko_case, dt_controller)
228 call profiler_start_region(
"Forward simulation")
229 loop_start = mpi_wtime()
231 do while (this%neko_case%time%t .lt. this%neko_case%time%end_time)
232 this%n_timesteps = this%n_timesteps + 1
234 call simulation_step(this%neko_case, dt_controller, loop_start)
236 call this%checkpoint%save(this%neko_case)
238 call profiler_end_region(
"Forward simulation")
240 call simulation_finalize(this%neko_case)
242 end subroutine simulation_run_forward
245 subroutine simulation_run_backward(this)
247 type(time_step_controller_t) :: dt_controller
248 real(kind=dp) :: loop_start
252 call dt_controller%init(this%neko_case%params)
254 call simulation_adjoint_init(this%adjoint_case, dt_controller)
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)
262 call simulation_adjoint_step(this%adjoint_case, dt_controller, cfl, &
265 call profiler_end_region(
"Adjoint simulation")
267 call simulation_adjoint_finalize(this%adjoint_case)
269 end subroutine simulation_run_backward
272 subroutine simulation_reset(this)
275 call reset(this%neko_case)
277 call this%checkpoint%reset()
279 end subroutine simulation_reset
282 subroutine simulation_write(this, idx)
284 integer,
intent(in) :: idx
286 call this%output_forward%sample(real(idx, kind=rp))
287 call this%output_adjoint%sample(real(idx, kind=rp))
289 end subroutine simulation_write
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.
subroutine, public reset(neko_case)
Reset the case data structure.
subroutine, public reset_adjoint(adjoint_case, neko_case)
Reset the adjoint case data structure.
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...
Base type of all fluid formulations.
Type to manage multiple adjoint scalar transport equations.