Add the Himmelblau function, as defined in Applied Nonlinear Programming.

diff --git a/optimization/test_functions/himmelblau.m b/optimization/test_functions/himmelblau.m
new file mode 100644 (file)
index 0000000..c973a86
--- /dev/null
+++ b/optimization/test_functions/himmelblau.m
@@ -0,0 +1,17 @@
+function f = himmelblau(x1, x2)
+  ##
+  ## The eponymous function defined by Himmelblau in his Applied
+  ## Nonlinear Programming, 1972.
+  ##
+  ## It has four identical local minima,
+  ##
+  ##   * f(3,2)                                    == 0
+  ##
+  ##   * f(-2.805118086952745,  3.131312518250573) == 0
+  ##
+  ##   * f(-3.779310253377747, -3.283185991286170) == 0
+  ##
+  ##   * f(3.584428340330492, -1.848126526964404)  == 0
+  ##
+  f = (x1^2 + x2 - 11)^2 + (x1 + x2^2 - 7)^2;
+end

diff --git a/optimization/test_functions/himmelblau1.m b/optimization/test_functions/himmelblau1.m
new file mode 100644 (file)
index 0000000..c975af2
--- /dev/null
+++ b/optimization/test_functions/himmelblau1.m
@@ -0,0 +1,12 @@
+function f = himmelblau1(x)
+  ##
+  ## A version of the Himmelblau function which takes a column
+  ## 2-vector instead of two distinct arguments. See himmelblau.m for
+  ## more information.
+  ##
+  if (length(x) == 2)
+    f = himmelblau(x(1), x(2));
+  else
+    f = NA;
+  end
+end

diff --git a/optimization/test_functions/himmelblau_gradient.m b/optimization/test_functions/himmelblau_gradient.m
new file mode 100644 (file)
index 0000000..b50e7b6
--- /dev/null
+++ b/optimization/test_functions/himmelblau_gradient.m
@@ -0,0 +1,9 @@
+function g = himmelblau_gradient(x1,x2)
+  ##
+  ## The gradient of the Himmelblau function. See himmelblau.m for
+  ## more information.
+  ##
+  f_x1 = 4*(x1^2 + x2 - 11)*x1 + 2*x2^2 + 2*x1 - 14;
+  f_x2 = 4*(x2^2 + x1 - 7)*x2 + 2*x1^2 + 2*x2 - 22;
+  g = [f_x1; f_x2];
+end

diff --git a/optimization/test_functions/himmelblau_gradient1.m b/optimization/test_functions/himmelblau_gradient1.m
new file mode 100644 (file)
index 0000000..91db28d
--- /dev/null
+++ b/optimization/test_functions/himmelblau_gradient1.m
@@ -0,0 +1,12 @@
+function g = himmelblau_gradient1(x)
+  ##
+  ## A version of the himmelblau_gradient() function which takes a
+  ## column 2-vector instead of two distinct arguments. See
+  ## himmelblau_gradient.m for more information.
+  ##
+  if (length(x) == 2)
+    g = himmelblau_gradient(x(1), x(2));
+  else
+    g = NA;
+  end
+end

diff --git a/optimization/test_functions/himmelblau_hessian.m b/optimization/test_functions/himmelblau_hessian.m
new file mode 100644 (file)
index 0000000..5a8752f
--- /dev/null
+++ b/optimization/test_functions/himmelblau_hessian.m
@@ -0,0 +1,11 @@
+function H = himmelblau_hessian(x1, x2)
+  ##
+  ## The Hessian of the Himmelblau function. See himmelblau.m for more
+  ## information.
+  ##
+  H = zeros(2,2);
+  H(1,1) = 12*x1^2 + 4*x2 - 42;
+  H(1,2) = 4*x1 + 4*x2;
+  H(2,1) = 4*x1 + 4*x2;
+  H(2,2) = 12*x2^2 + 4*x1 - 26;
+end

diff --git a/optimization/test_functions/himmelblau_hessian1.m b/optimization/test_functions/himmelblau_hessian1.m
new file mode 100644 (file)
index 0000000..84efb31
--- /dev/null
+++ b/optimization/test_functions/himmelblau_hessian1.m
@@ -0,0 +1,12 @@
+function H = himmelblau_hessian1(x)
+  ##
+  ## A version of the himmelblau_hessian() function which takes a column
+  ## 2-vector instead of two distinct arguments. See himmelblau_hessian.m
+  ## for more information.
+  ##
+  if (length(x) == 2)
+    H = himmelblau_hessian(x(1), x(2));
+  else
+    H = NA;
+  end
+end

diff --git a/tests/himmelblau_gradient_tests.m b/tests/himmelblau_gradient_tests.m
new file mode 100644 (file)
index 0000000..e76562d
--- /dev/null
+++ b/tests/himmelblau_gradient_tests.m
@@ -0,0 +1,26 @@
+## The gradient should be zero at the optimal points.
+
+unit_test_equals("himmelblau_gradient(3,2) == 0", ...
+                0, ...
+                himmelblau_gradient(3,2));
+
+x1 = -2.805118086952745;
+x2 = 3.131312518250573;
+msg = sprintf("himmelblau_gradient(%f, %f) == 0", x1, x2);
+unit_test_equals(msg, ...
+                true, ...
+                norm(himmelblau_gradient(x1, x2)) < 1e-13);
+
+x1 = -3.779310253377747;
+x2 = -3.283185991286170;
+msg = sprintf("himmelblau_gradient(%f, %f) == 0", x1, x2);
+unit_test_equals(msg, ...
+                true, ...
+                norm(himmelblau_gradient(x1, x2)) < 1e-13);
+
+x1 = 3.584428340330492;
+x2 = -1.848126526964404;
+msg = sprintf("himmelblau_gradient(%f, %f) == 0", x1, x2);
+unit_test_equals(msg, ...
+                true, ...
+                norm(himmelblau_gradient(x1, x2)) < 1e-13);

diff --git a/tests/himmelblau_tests.m b/tests/himmelblau_tests.m
new file mode 100644 (file)
index 0000000..75a145f
--- /dev/null
+++ b/tests/himmelblau_tests.m
@@ -0,0 +1,26 @@
+## Test the optimal points.
+
+unit_test_equals("himmelblau(3,2) == 0", ...
+                0, ...
+                himmelblau(3,2));
+
+x1 = -2.805118086952745;
+x2 = 3.131312518250573;
+msg = sprintf("himmelblau(%f, %f) == 0", x1, x2);
+unit_test_equals(msg, ...
+                true, ...
+                norm(himmelblau(x1, x2)) < 1e-12);
+
+x1 = -3.779310253377747;
+x2 = -3.283185991286170;
+msg = sprintf("himmelblau(%f, %f) == 0", x1, x2);
+unit_test_equals(msg, ...
+                true, ...
+                norm(himmelblau(x1, x2)) < 1e-12);
+
+x1 = 3.584428340330492;
+x2 = -1.848126526964404;
+msg = sprintf("himmelblau(%f, %f) == 0", x1, x2);
+unit_test_equals(msg, ...
+                true, ...
+                norm(himmelblau(x1, x2)) < 1e-12);