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 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
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
60 simulation_adjoint_step, simulation_adjoint_finalize
66 type(case_t),
public :: neko_case
70 class(fluid_scheme_incompressible_t),
public,
pointer :: fluid => null()
72 type(scalars_t),
public,
pointer :: scalars => null()
79 type(fld_file_output_t),
public :: output_forward
82 type(fld_file_output_t),
public :: output_adjoint
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
103 subroutine simulation_init(this, parameters)
105 type(json_file),
intent(inout) :: parameters
106 integer :: i, n_scalars
109 call neko_init(this%neko_case)
113 select type (fluid => this%neko_case%fluid)
114 type is (fluid_pnpn_t)
119 select type (adjoint_fluid => this%adjoint_case%fluid_adj)
121 this%adjoint_fluid => adjoint_fluid
125 if (
allocated(this%neko_case%scalars))
then
126 this%scalars => this%neko_case%scalars
129 if (
allocated(this%adjoint_case%adjoint_scalars))
then
130 this%adjoint_scalars => this%adjoint_case%adjoint_scalars
137 if (
allocated(this%neko_case%scalars))
then
138 n_scalars =
size(this%neko_case%scalars%scalar_fields)
140 call this%output_forward%init(sp,
'forward_fields', 4 + n_scalars)
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)
148 if (
allocated(this%neko_case%scalars))
then
150 call this%output_forward%fields%assign(4 + i, &
151 this%scalars%scalar_fields(i)%s)
156 if (
allocated(this%adjoint_case%adjoint_scalars))
then
157 n_scalars =
size(this%adjoint_case%adjoint_scalars%adjoint_scalar_fields)
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)
166 if (
allocated(this%adjoint_case%adjoint_scalars))
then
168 call this%output_adjoint%fields%assign(4 + i, &
169 this%adjoint_scalars%adjoint_scalar_fields(i)%s_adj)
173 end subroutine simulation_init
176 subroutine simulation_free(this)
179 call adjoint_free(this%adjoint_case)
180 call neko_finalize(this%neko_case)
182 end subroutine simulation_free
185 subroutine simulation_run_forward(this)
189 this%neko_case%time%t = 0.0_rp
190 this%neko_case%time%tstep = 0
193 call neko_solve(this%neko_case)
195 end subroutine simulation_run_forward
198 subroutine simulation_run_backward(this)
200 type(time_step_controller_t) :: dt_controller
201 real(kind=dp) :: tstep_loop_start_time
204 call dt_controller%init(this%neko_case%params)
206 call simulation_adjoint_init(this%adjoint_case, dt_controller)
209 cfl = this%adjoint_case%fluid_adj%compute_cfl(this%adjoint_case%time%dt)
210 tstep_loop_start_time = mpi_wtime()
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)
219 call simulation_adjoint_finalize(this%adjoint_case)
221 end subroutine simulation_run_backward
224 subroutine simulation_reset(this)
226 integer :: i, n_scalars
228 call reset(this%neko_case)
233 this%adjoint_case%time%t = 0.0_dp
234 this%adjoint_case%time%tstep = 0
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)
243 this%adjoint_case%adjoint_scalars%adjoint_scalar_fields(i)%s_adj)
247 end subroutine simulation_reset
250 subroutine simulation_write(this, idx)
252 integer,
intent(in) :: idx
254 call this%output_forward%sample(real(idx, kind=rp))
255 call this%output_adjoint%sample(real(idx, kind=rp))
257 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.
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.