Skip to content

Commit

Permalink
Merge pull request #310 from jcapriot/assert_tests
Browse files Browse the repository at this point in the history 
Assert tests
  • Loading branch information
@jcapriot
jcapriot authored Apr 13, 2023
2 parents d5a27bf + c3c2a54 commit 53002e1
Show file tree
Hide file tree
Showing 3 changed files with 243 additions and 104 deletions.
252 changes: 173 additions & 79 deletions discretize/tests.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -359,54 +359,142 @@ def orderTest(self):
"expectedOrders must have the same length as the meshTypes"
)

for ii_meshType, mesh_type in enumerate(self.meshTypes):
def test_func(n_cells):
max_h = self.setupMesh(n_cells)
err = self.getError()
return err, max_h

for mesh_type, order, tolerance in zip(
self.meshTypes, self.tolerance, self.expectedOrders
):
self._meshType = mesh_type
self._tolerance = self.tolerance[ii_meshType]
self._expectedOrder = self.expectedOrders[ii_meshType]

order = []
err_old = 0.0
max_h_old = 0.0
for ii, nc in enumerate(self.meshSizes):
# Leave these as setupMesh and getError for deprecated classes that might have extended these two methods
max_h = self.setupMesh(nc)
err = self.getError()
if ii == 0:
print("")
print(self._meshType + ": " + self.name)
print("_____________________________________________")
print(" h | error | e(i-1)/e(i) | order")
print("~~~~~~|~~~~~~~~~~~~~|~~~~~~~~~~~~~|~~~~~~~~~~")
print("{0:4d} | {1:8.2e} |".format(nc, err))
else:
order.append(np.log(err / err_old) / np.log(max_h / max_h_old))
print(
"{0:4d} | {1:8.2e} | {2:6.4f} | {3:6.4f}".format(
nc, err, err_old / err, order[-1]
)
)
err_old = err
max_h_old = max_h
print("---------------------------------------------")
passTest = np.mean(np.array(order)) > self._tolerance * self._expectedOrder
if passTest:
print(happiness[np.random.randint(len(happiness))])
else:
print("Failed to pass test on " + self._meshType + ".")
print(sadness[np.random.randint(len(sadness))])
print("")
self.assertGreater(
np.mean(np.array(order)), self._tolerance * self._expectedOrder
assert_expected_order(
test_func,
self.meshSizes,
expected_order=order,
rtol=np.abs(1 - tolerance),
test_type="mean_at_least",
)

# expectedOrders = deprecate_property("expectedOrders", "expectedOrders", removal_version="1.0.0")
# meshSizes = deprecate_property("meshSizes", "meshSizes", removal_version="1.0.0")
# meshTypes = deprecate_property("meshTypes", "meshTypes", removal_version="1.0.0")
# meshDimension = deprecate_property("meshDimension", "meshDimension", removal_version="1.0.0")
# setupMesh = deprecate_method("setupMesh", "setupMesh", removal_version="1.0.0")
# getError = deprecate_method("getError", "getError", removal_version="1.0.0")
# orderTest = deprecate_method("orderTest", "orderTest", removal_version="1.0.0")
# _meshType = deprecate_property("_meshType", "_meshType", removal_version="1.0.0")

def assert_expected_order(
func, n_cells, expected_order=2.0, rtol=0.15, test_type="mean"
):
"""Perform an order test.
For number of cells specified in `mesh_sizes` call `func`
and prints mesh size, error, ratio between current and previous error,
and estimated order of convergence.
Parameters
----------
func : callable
Function which should accept an integer representing the number of
discretizations on the domain and return a tuple of the error and
the discretization widths.
n_cells : array_like of int
List of number of discretizations to pass to func.
expected_order : float, optional
The expected order of accuracy for you test
rtol : float, optional
The relative tolerance of the order test.
test_type : {'mean', 'min', 'last', 'all', 'mean_at_least'}
Which property of the list of calculated orders to test.
Returns
-------
numpy.ndarray
The calculated order values on success
Raises
------
AssertionError
Notes
-----
For the different ``test_type`` arguments, different properties of the
order is tested:
- `mean`: the mean value of all calculated orders is tested for
approximate equality with the expected order.
- `min`: The minimimum value of calculated orders is tested for
approximate equality with the expected order.
- `last`: The last calculated order is tested for approximate equality
with the expected order.
- `all`: All calculated orders are tested for approximate equality with
the expected order.
- `mean_at_least`: The mean is tested to be at least approximately the
expected order. This is the default test for the previous ``OrderTest``
class in older versions of `discretize`.
Examples
--------
Testing the convergence order of an central difference operator
>>> from discretize.tests import assert_expected_order
>>> func = lambda y: np.cos(y)
>>> func_deriv = lambda y: -np.sin(y)
Define the function that returns the error and cell width for
a given number of discretizations.
>>> def deriv_error(n):
... # grid points
... nodes = np.linspace(0, 1, n+1)
... cc = 0.5 * (nodes[1:] + nodes[:-1])
... dh = nodes[1]-nodes[0]
... # evaluate the function on nodes
... node_eval = func(nodes)
... # calculate the numerical derivative
... num_deriv = (node_eval[1:] - node_eval[:-1]) / dh
... # calculate the true derivative
... true_deriv = func_deriv(cc)
... # compare the L-inf norm of the error vector
... err = np.linalg.norm(num_deriv - true_deriv, ord=np.inf)
... return err, dh
Then run the expected order test.
>>> assert_expected_order(deriv_error, [10, 20, 30, 40, 50])
"""
n_cells = np.asarray(n_cells, dtype=int)
if test_type not in ["mean", "min", "last", "all", "mean_at_least"]:
raise ValueError
orders = []
# Do first values:
nc = n_cells[0]
err_last, h_last = func(nc)

print("_______________________________________________________")
print(" nc | h | error | e(i-1)/e(i) | order ")
print("~~~~~~|~~~~~~~~~|~~~~~~~~~~~~~|~~~~~~~~~~~~~|~~~~~~~~~~")
print(f"{nc:^6d}|{h_last:^9.2e}|{err_last:^13.3e}| |")

for nc in n_cells[1:]:
err, h = func(nc)
order = np.log(err / err_last) / np.log(h / h_last)
print(f"{nc:^6d}|{h:^9.2e}|{err:^13.3e}|{err_last / err:^13.4f}|{order:^10.4f}")
err_last = err
h_last = h
orders.append(order)

print("-------------------------------------------------------")

try:
if test_type == "mean":
np.testing.assert_allclose(np.mean(orders), expected_order, rtol=rtol)
elif test_type == "mean_at_least":
assert np.mean(orders) > expected_order * (1 - rtol)
elif test_type == "min":
np.testing.assert_allclose(np.min(orders), expected_order, rtol=rtol)
elif test_type == "last":
np.testing.assert_allclose(orders[-1], expected_order, rtol=rtol)
elif test_type == "all":
np.testing.assert_allclose(orders, expected_order, rtol=rtol)
print(np.random.choice(happiness))
except AssertionError as err:
print(np.random.choice(sadness))
raise err

return orders


def rosenbrock(x, return_g=True, return_H=True):
Expand Down Expand Up @@ -561,52 +649,58 @@ def l2norm(x):
)
)

@requires({"matplotlib": matplotlib})
def _plot_it(axes, passed):
if axes is None:
axes = plt.subplot(111)
axes.loglog(h, E0, "b")
axes.loglog(h, E1, "g--")
axes.set_title(
"Check Derivative - {0!s}".format(("PASSED :)" if passed else "FAILED :("))
)
axes.set_xlabel("h")
axes.set_ylabel("Error")
leg = axes.legend(
[r"$\mathcal{O}(h)$", r"$\mathcal{O}(h^2)$"],
loc="best",
title=r"$f(x + h\Delta x) - f(x) - h g(x) \Delta x - \mathcal{O}(h^2) = 0$",
frameon=False,
)
plt.setp(leg.get_title(), fontsize=15)
plt.show()

# Ensure we are about precision
order0 = order0[E0[1:] > eps]
order1 = order1[E1[1:] > eps]
belowTol = order1.size == 0 and order0.size >= 0
# Make sure we get the correct order
correctOrder = order1.size > 0 and np.mean(order1) > tolerance * expectedOrder

passTest = belowTol or correctOrder

if passTest:
# belowTol = order1.size == 0 and order0.size >= 0
# # Make sure we get the correct order
# correctOrder = order1.size > 0 and np.mean(order1) > tolerance * expectedOrder
#
# passTest = belowTol or correctOrder
try:
if order1.size == 0:
# This should happen if the original function was a linear function
# Thus it has no higher order derivatives.
assert order0.size >= 0
else:
assert np.mean(order1) > tolerance * expectedOrder
print("{0!s} PASS! {1!s}".format("=" * 25, "=" * 25))
print(happiness[np.random.randint(len(happiness))] + "\n")
else:
print(np.random.choice(happiness) + "\n")
if plotIt:
_plot_it(ax, True)
except AssertionError as err:
print(
"{0!s}\n{1!s} FAIL! {2!s}\n{3!s}".format(
"*" * 57, "<" * 25, ">" * 25, "*" * 57
)
)
print(sadness[np.random.randint(len(sadness))] + "\n")

@requires({"matplotlib": matplotlib})
def plot_it(ax):
print(np.random.choice(sadness) + "\n")
if plotIt:
if ax is None:
ax = plt.subplot(111)
ax.loglog(h, E0, "b")
ax.loglog(h, E1, "g--")
ax.set_title(
"Check Derivative - {0!s}".format(
("PASSED :)" if passTest else "FAILED :(")
)
)
ax.set_xlabel("h")
ax.set_ylabel("Error")
leg = ax.legend(
[r"$\mathcal{O}(h)$", r"$\mathcal{O}(h^2)$"],
loc="best",
title=r"$f(x + h\Delta x) - f(x) - h g(x) \Delta x - \mathcal{O}(h^2) = 0$",
frameon=False,
)
plt.setp(leg.get_title(), fontsize=15)
plt.show()

plot_it(ax)
_plot_it(ax, False)
raise err

return passTest
return True


def get_quadratic(A, b, c=0):
Expand Down
71 changes: 70 additions & 1 deletion tests/base/test_tests.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,8 @@
import discretize
import subprocess
import numpy as np
from discretize.tests import assert_isadjoint
import scipy.sparse as sp
from discretize.tests import assert_isadjoint, checkDerivative, assert_expected_order


class TestAssertIsAdjoint:
Expand Down Expand Up @@ -67,6 +68,74 @@ def test_complex_clinear(self):
)


class TestCheckDerivative:
def test_simplePass(self):
def simplePass(x):
return np.sin(x), sp.diags(np.cos(x))

checkDerivative(simplePass, np.random.randn(5), plotIt=False)

def test_simpleFunction(self):
def simpleFunction(x):
return np.sin(x), lambda xi: np.cos(x) * xi

checkDerivative(simpleFunction, np.random.randn(5), plotIt=False)

def test_simpleFail(self):
def simpleFail(x):
return np.sin(x), -sp.diags(np.cos(x))

with pytest.raises(AssertionError):
checkDerivative(simpleFail, np.random.randn(5), plotIt=False)


@pytest.mark.parametrize("test_type", ["mean", "min", "last", "all", "mean_at_least"])
def test_expected_order_pass(test_type):
func = lambda y: np.cos(y)
func_deriv = lambda y: -np.sin(y)

def deriv_error(n):
# The l-inf norm of the average error vector does
# follow second order convergence.
nodes = np.linspace(0, 1, n + 1)
cc = 0.5 * (nodes[1:] + nodes[:-1])
dh = nodes[1] - nodes[0]
node_eval = func(nodes)
num_deriv = (node_eval[1:] - node_eval[:-1]) / dh
true_deriv = func_deriv(cc)
err = np.linalg.norm(num_deriv - true_deriv, ord=np.inf)
return err, dh

assert_expected_order(deriv_error, [10, 20, 30, 40, 50], test_type=test_type)


@pytest.mark.parametrize("test_type", ["mean", "min", "last", "all", "mean_at_least"])
def test_expected_order_failed(test_type):
func = lambda y: np.cos(y)
func_deriv = lambda y: -np.sin(y)

def deriv_error(n):
# The l2 norm of the average error vector does not
# follow second order convergence.
nodes = np.linspace(0, 1, n + 1)
cc = 0.5 * (nodes[1:] + nodes[:-1])
dh = nodes[1] - nodes[0]
node_eval = func(nodes)
num_deriv = (node_eval[1:] - node_eval[:-1]) / dh
true_deriv = func_deriv(cc)
err = np.linalg.norm(num_deriv - true_deriv)
return err, dh

with pytest.raises(AssertionError):
assert_expected_order(deriv_error, [10, 20, 30, 40, 50], test_type=test_type)


def test_expected_order_bad_test_type():
# should fail fast if given a bad test_type
with pytest.raises(ValueError):
assert_expected_order(None, [10, 20, 30, 40, 50], test_type="not_a_test")


@pytest.mark.skipif(not sys.platform.startswith("linux"), reason="Not Linux.")
def test_import_time():
# Relevant for the CLI: How long does it take to import?
Expand Down
Loading

0 comments on commit 53002e1

Please sign in to comment.