--- /dev/null
+function f = wood(x1, x2, x3, x4)
+ ##
+ ## The Wood function. See Dennis & Schnabel, Appendix B, problem #5.
+ ## This function has a global minimum at x=(1,1,1,1) with f(x) == 0.
+ ## The suggested starting point is x0=(-3, -1, -3, -1).
+ ##
+ f = 100*(x1^2 - x2)^2 + (x1 - 1)^2 + (x3 - 1)^2;
+ f = f + 90*(x3^2 - x4)^2;
+ f = f + 10.1*((x2 - 1)^2 + (x4 - 1)^2);
+ f = f + 19.8*(x2 - 1)*(x4 - 1);
+end
--- /dev/null
+function f = wood1(x)
+ ##
+ ## A version of the Wood function which takes a column 4-vector
+ ## instead of four distinct arguments. See wood.m for more
+ ## information.
+ ##
+ f = wood(x(1), x(2), x(3), x(4));
+end
--- /dev/null
+function g = wood_gradient(x1, x2, x3, x4)
+ ##
+ ## The gradient of the Wood function. See wood.m for more information.
+ ##
+ f_x1 = 400*(x1^2 - x2)*x1 + 2*x1 - 2;
+ f_x2 = -200*x1^2 + 220.2*x2 + 19.8*x4 - 40;
+ f_x3 = 360*(x3^2 - x4)*x3 + 2*x3 - 2;
+ f_x4 = -180*x3^2 + 19.8*x2 + 200.2*x4 - 40;
+
+ g = [f_x1; f_x2; f_x3; f_x4];
+end
--- /dev/null
+function g = wood_gradient1(x)
+ ##
+ ## A version of the wood_gradient() function which takes a column
+ ## 4-vector instead of four distinct arguments. See wood_gradient.m
+ ## for more information.
+ ##
+ g = wood_gradient(x(1), x(2), x(3), x(4));
+end
--- /dev/null
+function H = wood_hessian(x1, x2, x3, x4)
+ ##
+ ## The Hessian of the Wood function. See wood.m for more
+ ## information.
+ ##
+ H = zeros(4,4);
+ H(1,1) = 1200*x(1)^2 - 400*x(2) + 2;
+ H(1,2) = -400*x(1);
+ H(1,3) = 0;
+ H(1,4) = 0;
+ H(2,1) = H(1,2);
+ H(2,2) = 220.2;
+ H(2,3) = 0;
+ H(2,4) = 19.8;
+ H(3,1) = H(1,3);
+ H(3,2) = H(2,3);
+ H(3,3) = 1080*x(3)^2 - 360*x(4) + 2;
+ H(3,4) = -360*x(3);
+ H(4,1) = H(1,4);
+ H(4,2) = H(2,4);
+ H(4,3) = H(3,4);
+ H(4,4) = 200.2;
+end
--- /dev/null
+function H = wood_hessian1(x)
+ ##
+ ## A version of the wood_hessian() function which takes a column
+ ## 4-vector instead of four distinct arguments. See wood_hessian.m
+ ## for more information.
+ ##
+ H = zeros(4,4);
+ H(1,1) = 1200*x(1)^2 - 400*x(2) + 2;
+ H(1,2) = -400*x(1);
+ H(1,3) = 0;
+ H(1,4) = 0;
+ H(2,1) = H(1,2);
+ H(2,2) = 220.2;
+ H(2,3) = 0;
+ H(2,4) = 19.8;
+ H(3,1) = H(1,3);
+ H(3,2) = H(2,3);
+ H(3,3) = 1080*x(3)^2 - 360*x(4) + 2;
+ H(3,4) = -360*x(3);
+ H(4,1) = H(1,4);
+ H(4,2) = H(2,4);
+ H(4,3) = H(3,4);
+ H(4,4) = 200.2;
+end
--- /dev/null
+## The gradient should be zero at the optimal point.
+
+g = wood_gradient(1,1,1,1);
+unit_test_equals("wood_gradient(1,1,1,1) == 0", ...
+ 0, ...
+ wood_gradient(1,1,1,1));
--- /dev/null
+## Test the optimal point.
+
+f = wood(1,1,1,1);
+unit_test_equals("wood(1,1,1,1) == 0", ...
+ 0, ...
+ wood(1,1,1,1));