Neko-TOP
A portable framework for high-order spectral element flow toplogy optimization.
Loading...
Searching...
No Matches
neko_ext.f90
Go to the documentation of this file.
1
6
9module neko_ext
10 use case, only: case_t
11 use json_utils, only: json_get, json_get_or_default, json_extract_object, &
12 json_extract_item
13 use num_types, only: rp
14 use simcomp_executor, only: neko_simcomps
15 use flow_ic, only: set_flow_ic
16 use scalar_ic, only: set_scalar_ic
17 use field, only: field_t
18 use chkp_output, only: chkp_output_t
19 use output_controller, only: output_controller_t
20 ! for vector/field math
21 use math, only: copy
22 use device_math, only: device_copy
23 use neko_config, only : neko_bcknd_device
24 use vector, only: vector_t
25 use field, only: field_t
26 use utils, only: neko_error
27 use json_module, only : json_file
28 use scalars, only: scalars_t
30
31 implicit none
32
33 ! ========================================================================= !
34 ! Module interface
35 ! ========================================================================= !
36 private
39
40contains
41
42 ! ========================================================================= !
43 ! Public routines
44 ! ========================================================================= !
45
54 subroutine reset(neko_case)
55 type(case_t), intent(inout) :: neko_case
56 real(kind=rp) :: t
57 integer :: i, n_scalars
58 character(len=:), allocatable :: string_val
59 logical :: has_scalar, freezeflow
60 type(field_t), pointer :: u, v, w, p, s
61 type(json_file) :: json_subdict, scalar_params
62
63 ! ------------------------------------------------------------------------ !
64 ! Setup shorthand notation
65 ! ------------------------------------------------------------------------ !
66
67 u => neko_case%fluid%u
68 v => neko_case%fluid%v
69 w => neko_case%fluid%w
70 p => neko_case%fluid%p
71 if (allocated(neko_case%scalars)) then
72 s => neko_case%scalars%scalar_fields(1)%s
73 else
74 nullify(s)
75 end if
76
77 ! ------------------------------------------------------------------------ !
78 ! Reset the timing parameters
79 ! ------------------------------------------------------------------------ !
80
81 t = 0.0_rp
82 neko_case%time%t = t
83 neko_case%time%tstep = 0
84
85 ! Setup lagged time step parameters
86 neko_case%time%tlag = t
87 neko_case%time%dtlag = neko_case%time%dt
88 do i = 1, size(neko_case%time%tlag)
89 neko_case%time%tlag(i) = t - i*neko_case%time%dtlag(i)
90 end do
91
92 ! Reset the time step counter
93 call neko_case%output_controller%set_counter(neko_case%time)
94
95 ! Restart the fields
96 call neko_case%fluid%restart(neko_case%chkp)
97 if (allocated(neko_case%scalars)) then
98 call neko_case%scalars%restart(neko_case%chkp)
99 end if
100
101 ! Reset the external BDF coefficients
102 do i = 1, size(neko_case%time%dtlag)
103 call neko_case%fluid%ext_bdf%set_coeffs(neko_case%time%dtlag)
104 end do
105
106 ! Restart the simulation components
107 call neko_simcomps%restart(neko_case%time)
108
109 ! ------------------------------------------------------------------------ !
110 ! Reset the fluid field to the initial condition
111 ! ------------------------------------------------------------------------ !
112
113 call json_get(neko_case%params, &
114 'case.fluid.initial_condition.type', string_val)
115 call json_extract_object(neko_case%params, 'case.fluid.initial_condition', &
116 json_subdict)
117
118 if (trim(string_val) .ne. 'user') then
119 call set_flow_ic(u, v, w, p, &
120 neko_case%fluid%c_Xh, neko_case%fluid%gs_Xh, &
121 string_val, json_subdict)
122 else
123 call set_flow_ic(u, v, w, p, &
124 neko_case%fluid%c_Xh, neko_case%fluid%gs_Xh, &
125 neko_case%user%initial_conditions, neko_case%fluid%name)
126 end if
127
128 ! ------------------------------------------------------------------------ !
129 ! Reset the scalar field to the initial condition
130 ! ------------------------------------------------------------------------ !
131
132 call json_get_or_default(neko_case%params, &
133 'case.scalar.enabled', has_scalar, .false.)
134
135 if (has_scalar) then
136 if (neko_case%params%valid_path('case.adjoint_scalar')) then
137 ! we shouldn't fallback to the primal here.
138 call json_get(neko_case%params, &
139 'case.adjoint_scalar.initial_condition.type', string_val)
140 call json_extract_object(neko_case%params, &
141 'case.adjoint_scalar.initial_condition', json_subdict)
142
143 !call neko_log%section("Adjoint scalar initial condition ")
144
145 if (trim(string_val) .ne. 'user') then
146 call set_scalar_ic(neko_case%scalars%scalar_fields(1)%s, &
147 neko_case%fluid%c_Xh, neko_case%fluid%gs_Xh, string_val, &
148 json_subdict)
149 else
150 call neko_error("user defined ICs not implemented for " // &
151 "adjoint scalar")
152 ! call set_scalar_ic(this%adjoint_scalars%s_adj, &
153 ! this%adjoint_scalars%c_Xh, this%adjoint_scalars%gs_Xh, &
154 ! this%usr%scalar_user_ic, neko_case%params)
155 end if
156
157 ! call neko_log%end_section()
158 else
159
160 ! Handle multiple scalars
161 call neko_case%params%info('case.scalars', n_children = n_scalars)
162
163 do i = 1, n_scalars
164 call json_extract_item(neko_case%params, 'case.adjoint_scalars', &
165 i, scalar_params)
166 call json_get(scalar_params, 'initial_condition.type', string_val)
167 call json_extract_object(scalar_params, 'initial_condition', &
168 json_subdict)
169
170 if (trim(string_val) .ne. 'user') then
171 call set_scalar_ic(neko_case%scalars%scalar_fields(i)%s, &
172 neko_case%scalars%scalar_fields(i)%c_Xh, &
173 neko_case%scalars%scalar_fields(i)%gs_Xh, string_val, &
174 json_subdict)
175 else
176 call neko_error("user defined ICs not implemented for " // &
177 "adjoint scalar")
178 end if
179 end do
180 end if
181 end if
182
183 ! ------------------------------------------------------------------------ !
184 ! Reset the "freeze" parameter of the flow
185 ! ------------------------------------------------------------------------ !
186
187 call json_get_or_default(neko_case%params, &
188 'case.fluid.freeze_flow', freezeflow, .false.)
189
190 neko_case%fluid%freeze = freezeflow
191
192 end subroutine reset
193
202 subroutine setup_iteration(neko_case, iter)
203 type(case_t), intent(inout) :: neko_case
204 integer, intent(in) :: iter
205
206 character(len=:), allocatable :: dirname
207 character(len=80) :: file_name
208
209 if (iter .ne. 1) then
210 call reset(neko_case)
211 end if
212
213 call json_get_or_default(neko_case%params, &
214 'case.output_directory', dirname, './')
215
216 write (file_name, '(a,a,i5.5,a)') &
217 trim(adjustl(dirname)), '/topopt_', iter, '_.fld'
218
219 call neko_case%f_out%output_t%file_%init(trim(file_name))
220 call neko_case%output_controller%execute(neko_case%time, .true.)
221
222 end subroutine setup_iteration
223
231 subroutine vector_to_field(field, vector)
232 type(field_t), intent(inout) :: field
233 type(vector_t), intent(in) :: vector
234
235 ! first check they're the same size
236 if (field%size() .ne. vector%size()) then
237 call neko_error("vector and field are not the same size")
238 end if
239
240 if (neko_bcknd_device .eq. 1) then
241 call device_copy(field%x_d, vector%x_d, field%size())
242 else
243 call copy(field%x, vector%x, field%size())
244 end if
245
246 end subroutine vector_to_field
247
255 subroutine field_to_vector(vector, field)
256 type(vector_t), intent(inout) :: vector
257 type(field_t), intent(in) :: field
258
259 ! first check they're the same size
260 if (field%size() .ne. vector%size()) then
261 call neko_error("vector and field are not the same size")
262 end if
263
264 if (neko_bcknd_device .eq. 1) then
265 call device_copy(vector%x_d, field%x_d, field%size())
266 else
267 call copy(vector%x, field%x, field%size())
268 end if
269
270 end subroutine field_to_vector
271
280 subroutine get_scalar_indicies(i_primal, i_adjoint, scalars, &
281 adjoint_scalars, primal_name)
282 integer, intent(out) :: i_primal
283 integer, intent(out) :: i_adjoint
284 type(scalars_t), intent(inout) :: scalars
285 type(adjoint_scalars_t), intent(inout) :: adjoint_scalars
286 character(len=*), intent(in) :: primal_name
287 integer :: i, n_primal_scalars, n_adjoint_scalars
288
289 i_primal = -1
290 i_adjoint = -1
291 n_adjoint_scalars = size(adjoint_scalars%adjoint_scalar_fields)
292 n_primal_scalars = size(scalars%scalar_fields)
293
294 if ((n_adjoint_scalars .eq. 1) .and. (n_primal_scalars .eq. 1)) then
295 i_primal = 1
296 i_adjoint = 1
297 return
298 end if
299
300 do i = 1, n_adjoint_scalars
301 if (adjoint_scalars%adjoint_scalar_fields(i)%primal_name &
302 .eq. primal_name) then
303 i_adjoint = i
304 exit
305 end if
306 end do
307
308 do i = 1, n_primal_scalars
309 if (scalars%scalar_fields(i)%name .eq. primal_name) then
310 i_primal = i
311 exit
312 end if
313 end do
314
315 if (i_primal .le. 0 .or. i_adjoint .le. 0) then
316 call neko_error('could not find matching primal and adjoint' // &
317 ' scalar fields')
318 end if
319
320 end subroutine get_scalar_indicies
321
322end module neko_ext
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 field_to_vector(vector, field)
Field to vector.
Definition neko_ext.f90:256
subroutine, public reset(neko_case)
Reset the case data structure.
Definition neko_ext.f90:55
subroutine, public vector_to_field(field, vector)
Vector to field.
Definition neko_ext.f90:232
subroutine, public setup_iteration(neko_case, iter)
Setup the iteration.
Definition neko_ext.f90:203
subroutine, public get_scalar_indicies(i_primal, i_adjoint, scalars, adjoint_scalars, primal_name)
get scalar indices
Definition neko_ext.f90:282
Type to manage multiple adjoint scalar transport equations.