Merge pull request #149 from dolfin-adjoint/dolci/fix_tests_for_L2_riesz

Fix tests for L2 riesz maps

Author: dham

pyadjoint/control.py

@@ -1,5 +1,5 @@
 from .overloaded_type import OverloadedType, create_overloaded_object
-import warnings
+import logging
 
 
 class Control(object):
@@ -47,13 +47,13 @@ def tape_value(self):
 
     def get_derivative(self, options={}):
         if self.block_variable.adj_value is None:
-            warnings.warn("Adjoint value is None, is the functional independent of the control variable?")
+            logging.warning("Adjoint value is None, is the functional independent of the control variable?")
             return self.control._ad_convert_type(0., options=options)
         return self.control._ad_convert_type(self.block_variable.adj_value, options=options)
 
     def get_hessian(self, options={}):
         if self.block_variable.adj_value is None:
-            warnings.warn("Hessian value is None, is the functional independent of the control variable?")
+            logging.warning("Hessian value is None, is the functional independent of the control variable?")
             return self.control._ad_convert_type(0., options=options)
         return self.control._ad_convert_type(self.block_variable.hessian_value, options=options)
 

tests/firedrake_adjoint/test_assemble.py

@@ -3,7 +3,7 @@
 pytest.importorskip("firedrake")
 
 from numpy.random import rand
-from numpy.testing import assert_approx_equal
+from numpy.testing import assert_allclose
 
 from firedrake import *
 from firedrake.adjoint import *
@@ -22,8 +22,8 @@ def test_assemble_0_forms():
     # where stored as a Function instead of Vector()
     s = a1 + a2 + 2.0 * a3
     rf = ReducedFunctional(s, Control(u))
-    # derivative is: (1+2*u+6*u**2)*dx - summing is equivalent to testing with 1
-    assert_approx_equal(rf.derivative().vector().sum(), 1. + 2. * 4 + 6 * 16.)
+    dJdm = rf.derivative()
+    assert_allclose(dJdm.dat.data_ro, 1. + 2. * 4. + 6. * 16.)
 
 
 def test_assemble_0_forms_mixed():
@@ -41,7 +41,9 @@ def test_assemble_0_forms_mixed():
     s -= a3  # this is done deliberately to end up with an adj_input of 0.0 for the a3 AssembleBlock
     rf = ReducedFunctional(s, Control(u))
     # derivative is: (1+4*u)*dx - summing is equivalent to testing with 1
-    assert_approx_equal(rf.derivative().vector().sum(), 1. + 4. * 7)
+    dJdm = rf.derivative()
+    assert_allclose(dJdm.dat.data_ro[0], 1. + 4. * 7)
+    assert_allclose(dJdm.dat.data_ro[1], 0.0)
 
 
 def test_assemble_1_forms_adjoint():

tests/firedrake_adjoint/test_hessian.py

@@ -43,8 +43,8 @@ def test_simple_solve():
     tape.evaluate_adj()
 
     m = f.copy(deepcopy=True)
-    dJdm = Jhat.derivative().vector().inner(h.vector())
-    Hm = Jhat.hessian(h).vector().inner(h.vector())
+    dJdm = assemble(inner(Jhat.derivative(), h)*dx)
+    Hm = assemble(inner(Jhat.hessian(h), h)*dx)
     assert taylor_test(Jhat, m, h, dJdm=dJdm, Hm=Hm) > 2.9
 
 
@@ -133,15 +133,13 @@ def test_function():
 
     # Total derivative
     dJdc, dJdf = compute_gradient(J, [control_c, control_f])
-    dJdm = dJdc.vector().inner(h_c) + dJdf.vector().inner(h_f)
+    dJdm = assemble(dJdc * h_c * dx + dJdf * h_f * dx)
 
     # Hessian
     Hcc, Hff = compute_hessian(J, [control_c, control_f], [h_c, h_f])
-    Hm = Hff.vector().inner(h_f.vector()) + Hcc.vector().inner(h_c.vector())
-
+    Hm = assemble(Hcc * h_c * dx + Hff * h_f * dx)
     assert taylor_test(Jhat, [c, f], [h_c, h_f], dJdm=dJdm, Hm=Hm) > 2.9
 
-
 def test_nonlinear():
     tape = Tape()
     set_working_tape(tape)

tests/firedrake_adjoint/test_shape_derivatives.py

@@ -21,8 +21,10 @@ def test_sin_weak_spatial():
     computed = Jhat.derivative().vector().get_local()
 
     V = TestFunction(S)
-    dJV = div(V)*sin(x[0])*dx + V[0]*cos(x[0])*dx
-    actual = assemble(dJV).vector().get_local()
+    # Derivative (Cofunction)
+    dJV = assemble(div(V)*sin(x[0])*dx + V[0]*cos(x[0])*dx)
+    # Apply L2 riesz representation to obtain the gradient.
+    actual = dJV.riesz_representation().vector().get_local()
     assert np.allclose(computed, actual, rtol=1e-14)
 
 
@@ -71,7 +73,7 @@ def test_shape_hessian():
     s = Function(S,name="deform")
 
     mesh.coordinates.assign(mesh.coordinates + s)
-    J = assemble(sin(x[1])* dx(domain=mesh))
+    J = assemble(cos(x[1])* dx(domain=mesh))
     c = Control(s)
     Jhat = ReducedFunctional(J, c)
 
@@ -80,12 +82,13 @@ def test_shape_hessian():
     h.interpolate(as_vector((cos(x[2]), A*cos(x[1]), A*x[1])))
 
     # Second order taylor
-    dJdm = Jhat.derivative().vector().inner(h.vector())
-    Hm = compute_hessian(J, c, h).vector().inner(h.vector())
+    dJdm = assemble(inner(Jhat.derivative(), h)*dx)
+    Hm = assemble(inner(compute_hessian(J, c, h), h)*dx)
     r2 = taylor_test(Jhat, s, h, dJdm=dJdm, Hm=Hm)
+    print(r2)
     assert(r2 > 2.9)
     Jhat(s)
-    dJdmm_exact = derivative(derivative(sin(x[1])* dx(domain=mesh),x,h), x, h)
+    dJdmm_exact = derivative(derivative(cos(x[1]) * dx(domain=mesh), x, h), x, h)
     assert(np.isclose(assemble(dJdmm_exact), Hm))
 
 
@@ -108,8 +111,8 @@ def test_PDE_hessian_neumann():
     l = f*v*dx
     u = Function(V)
     solve(a==l, u, solver_parameters={'ksp_type':'preonly', 'pc_type':'lu',
-                                          "mat_type": "aij",
-                                          "pc_factor_mat_solver_type": "mumps"})
+                                        "mat_type": "aij",
+                                        "pc_factor_mat_solver_type": "mumps"})
     J = assemble(u*dx(domain=mesh))
     c = Control(s)
     Jhat = ReducedFunctional(J, c)
@@ -136,8 +139,8 @@ def test_PDE_hessian_neumann():
     Jhat(s)
 
     # # Second order taylor
-    dJdm = Jhat.derivative().vector().inner(h.vector())
-    Hm = compute_hessian(J, c, h).vector().inner(h.vector())
+    dJdm = assemble(inner(Jhat.derivative(), h)*dx)
+    Hm = assemble(inner(compute_hessian(J, c, h), h)*dx)
     r2 = taylor_test(Jhat, s, h, dJdm=dJdm, Hm=Hm)
     assert(r2>2.95)
 
@@ -190,8 +193,8 @@ def test_PDE_hessian_dirichlet():
     Jhat(s)
 
     # # Second order taylor
-    dJdm = Jhat.derivative().vector().inner(h.vector())
-    Hm = compute_hessian(J, c, h).vector().inner(h.vector())
+    dJdm = assemble(inner(Jhat.derivative(), h)*dx)
+    Hm = assemble(inner(compute_hessian(J, c, h), h)*dx)
     r2 = taylor_test(Jhat, s, h, dJdm=dJdm, Hm=Hm)
     assert(r2>2.95)