← Back to branch summary

~reducedmodelling/fluidity/ROM_Non-intrusive-ann

« back to all changes in this revision

Viewing changes to tests/burgers_mms_steady_adjoint_gradient_init/mms_E.bml

  • Committer: Patrick Farrell
  • Date: 2011-07-08 09:23:50 UTC
  • mto: This revision was merged to the branch mainline in revision 3521.
  • Revision ID: patrick.farrell06@imperial.ac.uk-20110708092350-jnpl7o3hwc2rxrx0
Add another test. The functional is the L2norm at t=0. The gradient is
computed with the adjoint, which only kicks in for the first forward equation.
As expected, the derivative test converges at second order.

However, if we change the functional to the L2norm at t=1, the
gradient is wrong. I'm still investigating.

Show diffs side-by-side

added added

removed removed

Lines of Context:
tests/burgers_mms_steady_adjoint_gradient_init/mms_E.bml
 
1
<?xml version='1.0' encoding='utf-8'?>
 
2
<burgers_equation>
 
3
  <simulation_name>
 
4
    <string_value lines="1">mms_adjoint_E</string_value>
 
5
  </simulation_name>
 
6
  <geometry>
 
7
    <dimension>
 
8
      <integer_value rank="0">1</integer_value>
 
9
    </dimension>
 
10
    <mesh name="CoordinateMesh">
 
11
      <from_file file_name="mms_E">
 
12
        <format name="triangle"/>
 
13
        <stat>
 
14
          <include_in_stat/>
 
15
        </stat>
 
16
      </from_file>
 
17
    </mesh>
 
18
    <mesh name="VelocityMesh">
 
19
      <from_mesh>
 
20
        <mesh name="CoordinateMesh"/>
 
21
        <stat>
 
22
          <exclude_from_stat/>
 
23
        </stat>
 
24
      </from_mesh>
 
25
    </mesh>
 
26
    <quadrature>
 
27
      <degree>
 
28
        <integer_value rank="0">4</integer_value>
 
29
      </degree>
 
30
    </quadrature>
 
31
  </geometry>
 
32
  <io>
 
33
    <dump_format>vtk</dump_format>
 
34
    <dump_period_in_timesteps>
 
35
      <constant>
 
36
        <integer_value rank="0">1</integer_value>
 
37
      </constant>
 
38
    </dump_period_in_timesteps>
 
39
    <output_mesh name="VelocityMesh"/>
 
40
  </io>
 
41
  <timestepping>
 
42
    <current_time>
 
43
      <real_value rank="0">0</real_value>
 
44
    </current_time>
 
45
    <timestep>
 
46
      <real_value rank="0">1</real_value>
 
47
    </timestep>
 
48
    <finish_time>
 
49
      <real_value rank="0">1</real_value>
 
50
    </finish_time>
 
51
    <nonlinear_iterations>
 
52
      <integer_value rank="0">1</integer_value>
 
53
    </nonlinear_iterations>
 
54
  </timestepping>
 
55
  <material_phase name="Fluid">
 
56
    <scalar_field name="Velocity" rank="0">
 
57
      <prognostic>
 
58
        <mesh name="VelocityMesh"/>
 
59
        <temporal_discretisation>
 
60
          <theta>
 
61
            <real_value rank="0">0.5</real_value>
 
62
          </theta>
 
63
          <relaxation>
 
64
            <real_value rank="0">0.5</real_value>
 
65
          </relaxation>
 
66
        </temporal_discretisation>
 
67
        <solver>
 
68
          <iterative_method name="preonly"/>
 
69
          <preconditioner name="lu"/>
 
70
          <relative_error>
 
71
            <real_value rank="0">1e-16</real_value>
 
72
          </relative_error>
 
73
          <absolute_error>
 
74
            <real_value rank="0">1e-12</real_value>
 
75
          </absolute_error>
 
76
          <max_iterations>
 
77
            <integer_value rank="0">10000</integer_value>
 
78
          </max_iterations>
 
79
          <never_ignore_solver_failures/>
 
80
          <diagnostics>
 
81
            <monitors/>
 
82
          </diagnostics>
 
83
        </solver>
 
84
        <remove_advection_term/>
 
85
        <initial_condition name="WholeMesh">
 
86
          <python>
 
87
            <string_value lines="20" type="python">def val(X, t):
 
88
  from math import sin, cos
 
89
  x = X[0]
 
90
  return sin(x) + cos(x)</string_value>
 
91
          </python>
 
92
        </initial_condition>
 
93
        <boundary_conditions name="bc">
 
94
          <surface_ids>
 
95
            <integer_value shape="2" rank="1">1 2</integer_value>
 
96
          </surface_ids>
 
97
          <type name="dirichlet">
 
98
            <python>
 
99
              <string_value lines="20" type="python">def val(X, t):
 
100
  from math import sin, cos
 
101
  x = X[0]
 
102
  return sin(x) + cos(x)</string_value>
 
103
            </python>
 
104
          </type>
 
105
        </boundary_conditions>
 
106
        <viscosity>
 
107
          <real_value rank="0">1</real_value>
 
108
        </viscosity>
 
109
        <stat/>
 
110
        <adjoint_storage>
 
111
          <exists_in_both/>
 
112
        </adjoint_storage>
 
113
        <scalar_field name="Source" rank="0">
 
114
          <prescribed>
 
115
            <value name="WholeMesh">
 
116
              <python>
 
117
                <string_value lines="20" type="python">def val(X, t):
 
118
  from math import sin, cos
 
119
  x = X[0]
 
120
  y = sin(x) + cos(x)
 
121
  return y</string_value>
 
122
              </python>
 
123
            </value>
 
124
          </prescribed>
 
125
        </scalar_field>
 
126
      </prognostic>
 
127
    </scalar_field>
 
128
    <scalar_field name="AnalyticalSolution" rank="0">
 
129
      <prescribed>
 
130
        <mesh name="VelocityMesh"/>
 
131
        <value name="WholeMesh">
 
132
          <python>
 
133
            <string_value lines="20" type="python">def val(X, t):
 
134
  from math import sin, cos
 
135
  x = X[0]
 
136
  return sin(x) + cos(x)</string_value>
 
137
          </python>
 
138
        </value>
 
139
        <adjoint_storage>
 
140
          <exists_in_forward/>
 
141
        </adjoint_storage>
 
142
      </prescribed>
 
143
    </scalar_field>
 
144
    <scalar_field name="Error" rank="0">
 
145
      <diagnostic>
 
146
        <mesh name="VelocityMesh"/>
 
147
        <algorithm name="scalar_python_diagnostic" material_phase_support="single">
 
148
          <string_value lines="20" type="python">soln = state.scalar_fields["Velocity"]
 
149
exact = state.scalar_fields["AnalyticalSolution"]
 
150
 
 
151
for i in range(field.node_count):
 
152
  field.set(i, abs(soln.node_val(i) - exact.node_val(i)))</string_value>
 
153
        </algorithm>
 
154
        <stat/>
 
155
        <adjoint_storage>
 
156
          <exists_in_forward/>
 
157
        </adjoint_storage>
 
158
      </diagnostic>
 
159
    </scalar_field>
 
160
    <scalar_field name="funct_derivative" rank="0">
 
161
      <diagnostic>
 
162
        <mesh name="VelocityMesh"/>
 
163
        <algorithm name="scalar_python_diagnostic" material_phase_support="single">
 
164
          <string_value lines="20" type="python">field.val[:] = 0.0
 
165
if time == 1.0:
 
166
  coord = state.vector_fields['Coordinate']
 
167
  adj_u = state.scalar_fields['AdjointVelocity']
 
168
  #print "Adjoint velocity(t=1.0): ", adj_u.val
 
169
  for ele in range(coord.element_count):
 
170
    t = Transform(ele, coord)
 
171
    shape = adj_u.ele_shape(ele)
 
172
    mass = t.shape_shape(shape, shape)
 
173
    field.addto(adj_u.ele_nodes(ele), numpy.dot(mass, adj_u.ele_val(ele)))
 
174
    
 
175
#print "Funct_derivative (diagnostic): ", field.val</string_value>
 
176
        </algorithm>
 
177
        <stat/>
 
178
        <adjoint_storage>
 
179
          <exists_in_adjoint/>
 
180
        </adjoint_storage>
 
181
      </diagnostic>
 
182
    </scalar_field>
 
183
  </material_phase>
 
184
  <adjoint>
 
185
    <functional name="time_integral_ad">
 
186
      <functional_value>
 
187
        <algorithm name="functional_value">
 
188
          <string_value lines="20" type="python">import numpy
 
189
 
 
190
J = 0.0
 
191
 
 
192
if time == 0.0:
 
193
  coord = states[0]["Fluid"].vector_fields["Coordinate"]
 
194
  u = states[0]["Fluid"].scalar_fields["Velocity"]
 
195
  
 
196
  #print "u(t=0): ", states[0]["Fluid"].scalar_fields["Velocity"].val
 
197
  #print "u(t=1): ", states[1]["Fluid"].scalar_fields["Velocity"].val
 
198
 
 
199
  for ele in range(coord.element_count):
 
200
    t = Transform(ele, coord)
 
201
    shape = u.ele_shape(ele)
 
202
    mass = t.shape_shape(shape, shape)
 
203
    J = J + 0.5 * numpy.dot(u.ele_val(ele), numpy.dot(mass, u.ele_val(ele)))
 
204
#  print "Got functional value: ", J
 
205
    
 
206
#if hasattr(J, "nominal_value"):
 
207
#  print "u.val: ", u.val
 
208
#  print "delJ/delu: ", [J.derivatives[x] for x in u.val]
 
209
#  print "J.nominal_value: ", J.nominal_value</string_value>
 
210
        </algorithm>
 
211
        <reduction>
 
212
          <sum/>
 
213
        </reduction>
 
214
      </functional_value>
 
215
      <functional_dependencies>
 
216
        <algorithm name="functional_dependencies">
 
217
          <string_value lines="20" type="python">def dependencies(times, timestep):
 
218
  if times[0] &lt;= 0.0 &lt; times[1]:
 
219
    return {"Fluid::Coordinate": [0],
 
220
            "Fluid::Velocity": [0]}
 
221
  else:
 
222
    return {}</string_value>
 
223
        </algorithm>
 
224
      </functional_dependencies>
 
225
    </functional>
 
226
    <controls>
 
227
      <control name="InitialCondition">
 
228
        <type field_name="Fluid::Velocity" name="initial_condition"/>
 
229
      </control>
 
230
    </controls>
 
231
    <debug>
 
232
      <replay_forward_run/>
 
233
      <check_action_transposes/>
 
234
    </debug>
 
235
  </adjoint>
 
236
</burgers_equation>