← Back to branch summary

~fluidity-core/fluidity/sediment_restructure

« back to all changes in this revision

Viewing changes to tests/shallow_water_optimisation_bounds/adjoint_C.swml

merge with branch including remap_surface_to_field but then removed this change and approached it a different way which seems to work

Show diffs side-by-side

added added

removed removed

Lines of Context:
tests/shallow_water_optimisation_bounds/adjoint_C.swml
1
 
<?xml version='1.0' encoding='utf8'?>
2
 
<shallow_water_options>
3
 
    <simulation_name>
4
 
        <string_value lines="1">wave_C</string_value>
5
 
    </simulation_name>
6
 
    <geometry>
7
 
        <dimension>
8
 
            <integer_value rank="0">3</integer_value>
9
 
        </dimension>
10
 
        <mesh name="CoordinateMesh">
11
 
            <from_file file_name="src/mesh_C">
12
 
                <format name="triangle"/>
13
 
                <dimension>
14
 
                    <integer_value rank="0">1</integer_value>
15
 
                </dimension>
16
 
                <stat>
17
 
                    <include_in_stat/>
18
 
                </stat>
19
 
            </from_file>
20
 
        </mesh>
21
 
        <mesh name="VelocityMesh">
22
 
            <from_mesh>
23
 
                <mesh name="PeriodicMesh"/>
24
 
                <mesh_continuity>
25
 
                    <string_value>discontinuous</string_value>
26
 
                </mesh_continuity>
27
 
                <stat>
28
 
                    <exclude_from_stat/>
29
 
                </stat>
30
 
            </from_mesh>
31
 
        </mesh>
32
 
        <mesh name="PressureMesh">
33
 
            <from_mesh>
34
 
                <mesh name="PeriodicMesh"/>
35
 
                <mesh_shape>
36
 
                    <polynomial_degree>
37
 
                        <integer_value rank="0">2</integer_value>
38
 
                    </polynomial_degree>
39
 
                </mesh_shape>
40
 
                <stat>
41
 
                    <exclude_from_stat/>
42
 
                </stat>
43
 
            </from_mesh>
44
 
        </mesh>
45
 
        <mesh name="PeriodicMesh">
46
 
            <from_mesh>
47
 
                <mesh name="CoordinateMesh"/>
48
 
                <periodic_boundary_conditions name="periodicity">
49
 
                    <physical_boundary_ids>
50
 
                        <integer_value rank="1" shape="1">1</integer_value>
51
 
                    </physical_boundary_ids>
52
 
                    <aliased_boundary_ids>
53
 
                        <integer_value rank="1" shape="1">2</integer_value>
54
 
                    </aliased_boundary_ids>
55
 
                    <coordinate_map>
56
 
                        <string_value lines="20" type="python">def val(X,t):
57
 
  result = list(X)
58
 
  result[0]=result[0]-1.0
59
 
  return result</string_value>
60
 
                    </coordinate_map>
61
 
                </periodic_boundary_conditions>
62
 
                <stat>
63
 
                    <exclude_from_stat/>
64
 
                </stat>
65
 
            </from_mesh>
66
 
        </mesh>
67
 
        <quadrature>
68
 
            <degree>
69
 
                <integer_value rank="0">4</integer_value>
70
 
            </degree>
71
 
        </quadrature>
72
 
    </geometry>
73
 
    <io>
74
 
        <dump_format>
75
 
            <string_value>vtk</string_value>
76
 
        </dump_format>
77
 
        <dump_period_in_timesteps>
78
 
            <constant>
79
 
                <integer_value rank="0">1</integer_value>
80
 
            </constant>
81
 
        </dump_period_in_timesteps>
82
 
        <output_mesh name="CoordinateMesh"/>
83
 
    </io>
84
 
    <timestepping>
85
 
        <current_time>
86
 
            <real_value rank="0">0</real_value>
87
 
        </current_time>
88
 
        <timestep>
89
 
            <real_value rank="0">0.0625</real_value>
90
 
        </timestep>
91
 
        <nonlinear_iterations>
92
 
            <integer_value rank="0">1</integer_value>
93
 
        </nonlinear_iterations>
94
 
        <finish_time>
95
 
            <real_value rank="0">1</real_value>
96
 
        </finish_time>
97
 
    </timestepping>
98
 
    <physical_parameters>
99
 
        <gravity>
100
 
            <magnitude>
101
 
                <real_value rank="0">1</real_value>
102
 
            </magnitude>
103
 
            <vector_field name="GravityDirection" rank="1">
104
 
                <prescribed>
105
 
                    <mesh name="CoordinateMesh"/>
106
 
                    <value name="WholeMesh">
107
 
                        <constant>
108
 
                            <real_value rank="1" shape="3" dim1="dim">0 0 -1</real_value>
109
 
                        </constant>
110
 
                    </value>
111
 
                    <output/>
112
 
                    <stat>
113
 
                        <include_in_stat/>
114
 
                    </stat>
115
 
                    <detectors>
116
 
                        <exclude_from_detectors/>
117
 
                    </detectors>
118
 
                    <adjoint_storage>
119
 
                        <exists_in_forward/>
120
 
                    </adjoint_storage>
121
 
                </prescribed>
122
 
            </vector_field>
123
 
        </gravity>
124
 
    </physical_parameters>
125
 
    <material_phase name="Fluid">
126
 
        <vector_field name="Velocity" rank="1">
127
 
            <prognostic>
128
 
                <mesh name="VelocityMesh"/>
129
 
                <equation name="ShallowWater"/>
130
 
                <spatial_discretisation>
131
 
                    <discontinuous_galerkin>
132
 
                        <advection_scheme>
133
 
                            <none/>
134
 
                        </advection_scheme>
135
 
                    </discontinuous_galerkin>
136
 
                    <conservative_advection>
137
 
                        <real_value rank="0">0</real_value>
138
 
                    </conservative_advection>
139
 
                </spatial_discretisation>
140
 
                <solver>
141
 
                    <iterative_method name="preonly"/>
142
 
                    <preconditioner name="lu"/>
143
 
                    <relative_error>
144
 
                        <real_value rank="0">1e-07</real_value>
145
 
                    </relative_error>
146
 
                    <max_iterations>
147
 
                        <integer_value rank="0">500</integer_value>
148
 
                    </max_iterations>
149
 
                    <never_ignore_solver_failures/>
150
 
                    <diagnostics>
151
 
                        <monitors/>
152
 
                    </diagnostics>
153
 
                </solver>
154
 
                <initial_condition name="WholeMesh">
155
 
                    <python>
156
 
                        <string_value lines="20" type="python">def val(X, t):
157
 
  import shallow_water_optimisation_bounds as constants
158
 
  return constants.u_exact(X[0], t)</string_value>
159
 
                    </python>
160
 
                </initial_condition>
161
 
                <output/>
162
 
                <stat>
163
 
                    <include_in_stat/>
164
 
                    <previous_time_step>
165
 
                        <exclude_from_stat/>
166
 
                    </previous_time_step>
167
 
                    <nonlinear_field>
168
 
                        <exclude_from_stat/>
169
 
                    </nonlinear_field>
170
 
                </stat>
171
 
                <convergence>
172
 
                    <include_in_convergence/>
173
 
                </convergence>
174
 
                <detectors>
175
 
                    <include_in_detectors/>
176
 
                </detectors>
177
 
                <steady_state>
178
 
                    <include_in_steady_state/>
179
 
                </steady_state>
180
 
                <consistent_interpolation/>
181
 
            </prognostic>
182
 
        </vector_field>
183
 
        <scalar_field name="LayerThickness" rank="0">
184
 
            <prognostic>
185
 
                <mesh name="PressureMesh"/>
186
 
                <spatial_discretisation>
187
 
                    <continuous_galerkin>
188
 
                        <advection_terms>
189
 
                            <exclude_advection_terms/>
190
 
                        </advection_terms>
191
 
                    </continuous_galerkin>
192
 
                    <conservative_advection>
193
 
                        <real_value rank="0">0</real_value>
194
 
                    </conservative_advection>
195
 
                </spatial_discretisation>
196
 
                <temporal_discretisation>
197
 
                    <theta>
198
 
                        <real_value rank="0">0.5</real_value>
199
 
                    </theta>
200
 
                    <relaxation>
201
 
                        <real_value rank="0">1</real_value>
202
 
                    </relaxation>
203
 
                </temporal_discretisation>
204
 
                <solver>
205
 
                    <iterative_method name="preonly"/>
206
 
                    <preconditioner name="lu"/>
207
 
                    <relative_error>
208
 
                        <real_value rank="0">1e-07</real_value>
209
 
                    </relative_error>
210
 
                    <max_iterations>
211
 
                        <integer_value rank="0">500</integer_value>
212
 
                    </max_iterations>
213
 
                    <never_ignore_solver_failures/>
214
 
                    <cache_solver_context/>
215
 
                    <diagnostics>
216
 
                        <monitors/>
217
 
                    </diagnostics>
218
 
                </solver>
219
 
                <initial_condition name="WholeMesh">
220
 
                    <python>
221
 
                        <string_value lines="20" type="python">def val(X, t):
222
 
  import shallow_water_optimisation_bounds as constants
223
 
  #return constants.eta_exact(X[0], t)
224
 
  # Instead of prescribing the exact eta we will optimise to it
225
 
  return 0.0</string_value>
226
 
                    </python>
227
 
                </initial_condition>
228
 
                <mean_layer_thickness>
229
 
                    <real_value rank="0">1</real_value>
230
 
                </mean_layer_thickness>
231
 
                <output/>
232
 
                <stat/>
233
 
                <consistent_interpolation/>
234
 
            </prognostic>
235
 
        </scalar_field>
236
 
        <scalar_field name="AnalyticalLayerThickness" rank="0">
237
 
            <prescribed>
238
 
                <mesh name="PressureMesh"/>
239
 
                <value name="WholeMesh">
240
 
                    <python>
241
 
                        <string_value lines="20" type="python">def val(X,t):
242
 
 import shallow_water_optimisation_bounds as constants
243
 
 x=X[0]
244
 
 return constants.eta_exact(x, t)</string_value>
245
 
                    </python>
246
 
                </value>
247
 
                <output/>
248
 
                <stat/>
249
 
                <detectors>
250
 
                    <exclude_from_detectors/>
251
 
                </detectors>
252
 
                <adjoint_storage>
253
 
                    <exists_in_forward/>
254
 
                </adjoint_storage>
255
 
            </prescribed>
256
 
        </scalar_field>
257
 
        <scalar_field name="LayerThicknessError" rank="0">
258
 
            <diagnostic>
259
 
                <algorithm name="scalar_difference" source_field_2_type="scalar" source_field_1_name="LayerThickness" source_field_2_name="AnalyticalLayerThickness" material_phase_support="single" source_field_1_type="scalar">
260
 
                    <absolute_difference/>
261
 
                </algorithm>
262
 
                <mesh name="PressureMesh"/>
263
 
                <output/>
264
 
                <stat/>
265
 
                <convergence>
266
 
                    <include_in_convergence/>
267
 
                </convergence>
268
 
                <detectors>
269
 
                    <include_in_detectors/>
270
 
                </detectors>
271
 
                <steady_state>
272
 
                    <include_in_steady_state/>
273
 
                </steady_state>
274
 
                <adjoint_storage>
275
 
                    <exists_in_forward/>
276
 
                </adjoint_storage>
277
 
            </diagnostic>
278
 
        </scalar_field>
279
 
        <scalar_field name="InitEta_UpperBound" rank="0">
280
 
            <prescribed>
281
 
                <mesh name="PressureMesh"/>
282
 
                <value name="WholeMesh">
283
 
                    <constant>
284
 
                        <real_value rank="0">1.1</real_value>
285
 
                    </constant>
286
 
                </value>
287
 
                <output/>
288
 
                <stat/>
289
 
                <detectors>
290
 
                    <exclude_from_detectors/>
291
 
                </detectors>
292
 
                <adjoint_storage>
293
 
                    <exists_in_forward/>
294
 
                </adjoint_storage>
295
 
            </prescribed>
296
 
        </scalar_field>
297
 
        <scalar_field name="InitEta_LowerBound" rank="0">
298
 
            <prescribed>
299
 
                <mesh name="PressureMesh"/>
300
 
                <value name="WholeMesh">
301
 
                    <constant>
302
 
                        <real_value rank="0">-1.1</real_value>
303
 
                    </constant>
304
 
                </value>
305
 
                <output/>
306
 
                <stat/>
307
 
                <detectors>
308
 
                    <exclude_from_detectors/>
309
 
                </detectors>
310
 
                <adjoint_storage>
311
 
                    <exists_in_forward/>
312
 
                </adjoint_storage>
313
 
            </prescribed>
314
 
        </scalar_field>
315
 
        <vector_field name="AnalyticalVelocity" rank="1">
316
 
            <prescribed>
317
 
                <mesh name="VelocityMesh"/>
318
 
                <value name="WholeMesh">
319
 
                    <python>
320
 
                        <string_value lines="20" type="python">def val(X,t):
321
 
 import shallow_water_optimisation_bounds as constants
322
 
 x = X[0]
323
 
 return constants.u_exact(x, t)</string_value>
324
 
                    </python>
325
 
                </value>
326
 
                <output/>
327
 
                <stat>
328
 
                    <include_in_stat/>
329
 
                </stat>
330
 
                <detectors>
331
 
                    <exclude_from_detectors/>
332
 
                </detectors>
333
 
                <adjoint_storage>
334
 
                    <exists_in_forward/>
335
 
                </adjoint_storage>
336
 
            </prescribed>
337
 
        </vector_field>
338
 
        <vector_field name="VelocityError" rank="1">
339
 
            <diagnostic>
340
 
                <algorithm name="vector_difference" source_field_2_type="vector" source_field_1_name="Velocity" source_field_2_name="AnalyticalVelocity" material_phase_support="single" source_field_1_type="vector">
341
 
                    <absolute_difference/>
342
 
                </algorithm>
343
 
                <mesh name="VelocityMesh"/>
344
 
                <output/>
345
 
                <stat>
346
 
                    <include_in_stat/>
347
 
                </stat>
348
 
                <convergence>
349
 
                    <include_in_convergence/>
350
 
                </convergence>
351
 
                <detectors>
352
 
                    <include_in_detectors/>
353
 
                </detectors>
354
 
                <steady_state>
355
 
                    <include_in_steady_state/>
356
 
                </steady_state>
357
 
                <adjoint_storage>
358
 
                    <exists_in_forward/>
359
 
                </adjoint_storage>
360
 
            </diagnostic>
361
 
        </vector_field>
362
 
    </material_phase>
363
 
    <adjoint>
364
 
        <functional name="integral_eta_t1">
365
 
            <functional_value>
366
 
                <algorithm name="functional_value">
367
 
                    <string_value lines="20" type="python">J = 0.0
368
 
T = 1.0 # the time at which to evaluate
369
 
if time &lt; T &lt;= time+dt:
370
 
  import numpy
371
 
  eta_prev = states[n-1]["Fluid"].scalar_fields["LayerThickness"]
372
 
  eta      = states[n]["Fluid"].scalar_fields["LayerThickness"]
373
 
  eta_exact_prev = states[n-1]["Fluid"].scalar_fields["AnalyticalLayerThickness"]
374
 
  eta_exact      = states[n]["Fluid"].scalar_fields["AnalyticalLayerThickness"]
375
 
    
376
 
  # We want to temporally interpolate to evaluate eta at t=1.0
377
 
  alpha = (time + dt - T) / dt
378
 
  assert 0 &lt;= alpha &lt; 1
379
 
  tmp_eta = alpha * (eta_prev.val - eta_exact_prev.val) + (1-alpha) * (eta.val - eta_exact.val)
380
 
  
381
 
  # Now we want to integrate that over space
382
 
  coord = states[0]["Fluid"].vector_fields["Coordinate"]
383
 
  assert eta.element_count == eta_prev.element_count == coord.element_count
384
 
  for ele in range(coord.element_count):
385
 
    t = Transform(ele, coord)
386
 
    shape = eta_prev.ele_shape(ele)
387
 
    mass = t.shape_shape(shape, shape)
388
 
    nodes = eta_prev.ele_nodes(ele)
389
 
    J = J + numpy.dot(tmp_eta[nodes], numpy.dot(mass, tmp_eta[nodes]))</string_value>
390
 
                </algorithm>
391
 
                <reduction>
392
 
                    <sum/>
393
 
                </reduction>
394
 
            </functional_value>
395
 
            <functional_dependencies>
396
 
                <algorithm name="functional_dependencies">
397
 
                    <string_value lines="20" type="python">def dependencies(times, timestep):
398
 
  if times[0] &lt; 1.0 &lt;= times[1]:
399
 
    return {"Fluid::Coordinate": [0],
400
 
            "Fluid::LayerThickness": [timestep-1, timestep],
401
 
            "Fluid::AnalyticalLayerThickness": [timestep-1, timestep]
402
 
            }
403
 
  else:
404
 
    return {}</string_value>
405
 
                </algorithm>
406
 
            </functional_dependencies>
407
 
        </functional>
408
 
        <controls>
409
 
            <control name="InitEta">
410
 
                <type name="initial_condition" field_name="Fluid::LayerThickness"/>
411
 
                <bounds>
412
 
                    <upper_bound field_name="Fluid::InitEta_UpperBound"/>
413
 
                    <lower_bound field_name="Fluid::InitEta_LowerBound"/>
414
 
                </bounds>
415
 
            </control>
416
 
            <load_controls/>
417
 
        </controls>
418
 
        <debug>
419
 
            <html_output/>
420
 
        </debug>
421
 
    </adjoint>
422
 
</shallow_water_options>
 
 
b'\\ No newline at end of file'