Neko-TOP
A portable framework for high-order spectral element flow toplogy optimization.
Loading...
Searching...
No Matches
adjoint_mixing_scalar_source_term.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! this is a such a dumb name
36 use num_types, only : rp
37 use field_list, only : field_list_t
38 use field, only: field_t
39 use field_registry, only: neko_field_registry
40 use scratch_registry, only: neko_scratch_registry
41 use json_module, only : json_file
42 use time_state, only: time_state_t
43 use json_utils, only: json_get, json_get_or_default
44 use source_term, only : source_term_t
45 use coefs, only : coef_t
46 use field_math, only: field_add2s2, field_copy, field_cadd
47 use mask_ops, only: mask_exterior_const, compute_masked_volume
48 use point_zone, only: point_zone_t
49 implicit none
50 private
51
52 ! this will be the adjoint forcing from Casper's objective function
53 ! TODO
54 ! it actually isn't this. Infact, his is a surface integral on the outlet
55 ! Which means we get strange BC's in our adjoint problem.
56 ! This source term would be if we had a certain volume that we wanted more
57 ! mixed
58 ! the forcing is of the form:
59 ! \f$ \phi - \phi_ref \f$
60 ! ie, difference between it and the average.
61 type, public, extends(source_term_t) :: adjoint_mixing_scalar_source_term_t
63 type(field_t), pointer :: s => null()
65 real(kind=rp) :: obj_scale
67 real(kind=rp) :: phi_ref
69 class(point_zone_t), pointer :: mask => null()
71 logical :: if_mask
73 real(kind=rp) :: mask_volume
74 contains
76 procedure, pass(this) :: init => &
77 adjoint_mixing_scalar_source_term_init_from_json
79 procedure, pass(this) :: init_from_components => &
80 adjoint_mixing_scalar_source_term_init_from_components
82 procedure, pass(this) :: free => adjoint_mixing_scalar_source_term_free
84 procedure, pass(this) :: compute_ => &
85 adjoint_mixing_scalar_source_term_compute
87
88contains
95 subroutine adjoint_mixing_scalar_source_term_init_from_json(this, &
96 json, fields, coef, variable_name)
97 class(adjoint_mixing_scalar_source_term_t), intent(inout) :: this
98 type(json_file), intent(inout) :: json
99 type(field_list_t), intent(in), target :: fields
100 type(coef_t), intent(in), target :: coef
101 character(len=*), intent(in) :: variable_name
102
103
104 end subroutine adjoint_mixing_scalar_source_term_init_from_json
105
106
107 subroutine adjoint_mixing_scalar_source_term_init_from_components(this, &
108 f_s, s, obj_scale, phi_ref, mask, if_mask, coef)
109 class(adjoint_mixing_scalar_source_term_t), intent(inout) :: this
110 type(field_t), pointer, intent(in) :: f_s
111 real(kind=rp), intent(in) :: phi_ref
112 type(field_t), intent(in), target :: s
113 class(point_zone_t), intent(in), target :: mask
114 type(field_list_t) :: fields
115 type(coef_t) :: coef
116 real(kind=rp) :: start_time
117 real(kind=rp) :: end_time
118 real(kind=rp) :: obj_scale
119
120 logical :: if_mask
121
122 start_time = 0.0_rp
123 end_time = 100000000.0_rp
124
125
126 call this%free()
127
128 ! this is copying the fluid source term init
129 ! We package the fields for the source term to operate on in a field list.
130 call fields%init(1)
131 call fields%assign(1, f_s)
132
133 call this%init_base(fields, coef, start_time, end_time)
134
135 ! point everything in the correct places
136 this%s => s
137 this%obj_scale = obj_scale
138 this%phi_ref = phi_ref
139 this%if_mask = if_mask
140
141 if (this%if_mask) then
142 this%mask => mask
143 this%mask_volume = compute_masked_volume(this%mask, coef)
144 else
145 this%mask_volume = coef%volume
146 end if
147
148 end subroutine adjoint_mixing_scalar_source_term_init_from_components
149
151 subroutine adjoint_mixing_scalar_source_term_free(this)
152 class(adjoint_mixing_scalar_source_term_t), intent(inout) :: this
153
154 call this%free_base()
155 end subroutine adjoint_mixing_scalar_source_term_free
156
160 subroutine adjoint_mixing_scalar_source_term_compute(this, time)
161 class(adjoint_mixing_scalar_source_term_t), intent(inout) :: this
162 type(time_state_t), intent(in) :: time
163 type(field_t), pointer :: fs
164 type(field_t), pointer :: work
165 integer :: temp_indices(1)
166
167
168 fs => this%fields%get(1)
169
170 call neko_scratch_registry%request_field(work, temp_indices(1))
171 ! \phi
172 call field_copy(work, this%s)
173 ! \phi - \phi_ref
174 call field_cadd(work, -this%phi_ref)
175 ! mask
176 if (this%if_mask) then
177 call mask_exterior_const(work, this%mask, 0.0_rp)
178 end if
179
180 ! append to RHS with scaling and mask volume
181 call field_add2s2(fs, work, this%obj_scale / this%mask_volume)
182 call neko_scratch_registry%relinquish_field(temp_indices)
183
184
185 end subroutine adjoint_mixing_scalar_source_term_compute
186
Implements the adjoint_mixing_scalar_source_term type.
Some common Masking operations we may need.
Definition mask_ops.f90:34