run-tests.m

   1 #!/usr/bin/octave --silent
   2
   3 unit_init(1, {});
   4
   5 unit_test_equals("sin[0] == 0", ...
   6                  0, ...
   7                  divided_difference(@sin, 0));
   8
   9 unit_test_equals("sin[0, pi] == 0", ...
  10                  0, ...
  11                  divided_difference(@sin, [0,pi]));
  12
  13 unit_test_equals("sin[0, pi, 2*pi] == 0", ...
  14                  0, ...
  15                  divided_difference(@sin, [0,pi,2*pi]));
  16
  17 unit_test_equals("zero order divided_difference_coefficients", ...
  18                  [1], ...
  19                  divided_difference_coefficients([0]));
  20
  21 unit_test_equals("first order divided_difference_coefficients", ...
  22                  [-1, 1] / pi, ...
  23                  divided_difference_coefficients([0, pi]));
  24
  25 unit_test_equals("second order divided_difference_coefficients", ...
  26                  [1, -2, 1] / (2*pi^2), ...
  27                  divided_difference_coefficients([0, pi, 2*pi]));
  28
  29
  30 unit_test_equals("1 is odd", ...
  31                  true, ...
  32                  odd(1));
  33
  34 unit_test_equals("1 is not even", ...
  35                  false, ...
  36                  even(1));
  37
  38 unit_test_equals("2 is not odd", ...
  39                  false, ...
  40                  odd(2));
  41
  42 unit_test_equals("2 is even", ...
  43                  true, ...
  44                  even(2));
  45
  46 expected_A = [1, 0, 0, 0, 0; ...
  47               16, -32,  16,  0,  0; ...
  48               0,   16, -32,  16, 0; ...
  49               0,   0,   16, -32, 16; ...
  50               0,   0,   0,   0,  1];
  51 unit_test_equals("Homework #1 problem #1 Poisson matrix is correct", ...
  52                  true, ...
  53                  expected_A == poisson_matrix(4, 0, 1));
  54
  55
  56 g = @(x) 1 + atan(x);
  57 expected_fp = 2.1323;
  58 tol = 1 / 10^10;
  59 x0 = 2.4;
  60 unit_test_equals("Homework #2 problem #5 fixed point is correct", ...
  61                  expected_fp, ...
  62                  fixed_point_method(g, tol, x0));
  63
  64
  65 h = 0.5;
  66 g1 = @(u) 1 + h*exp(-u(1)^2)/(1+u(2)^2);
  67 g2 = @(u) 0.5 + h*atan(u(1)^2 + u(2)^2);
  68 my_g = @(u) [g1(u), g2(u)];
  69 tol = 1 / 10^9;
  70 u0 = [1,1];
  71 expected_fp = [1.0729, 1.0821];
  72 unit_test_equals("Homework #3 problem #3i fixed point is correct", ...
  73                  expected_fp, ...
  74                  fixed_point_method(my_g, tol, u0));
  75
  76
  77 f = @(x) x^6 - x - 1;
  78 f_prime = @(x) 6*x^5 - 1;
  79 tol = 1/1000000;
  80 x0 = 2;
  81 expected_root = 1.1347;
  82 unit_test_equals("Newton's method agrees with Haskell", ...
  83                  expected_root, ...
  84                  newtons_method(f, f_prime, tol, x0));
  85
  86
  87
  88 f1 = @(u) u(1)^2 + u(1)*u(2)^3 - 9;
  89 f2 = @(u) 3*u(1)^2*u(2) - u(2)^3 - 4;
  90 f = @(u) [f1(u); f2(u)];
  91 ## The partials for the Jacobian.
  92 f1x = @(u) 2*u(1) + u(2)^3;
  93 f1y = @(u) 3*u(1)*u(2)^2;
  94 f2x = @(u) 6*u(1)*u(2);
  95 f2y = @(u) 3*u(1)^2 - 3*u(2)^2;
  96 ## f_prime == Jacobian.
  97 f_prime = @(u) [ f1x(u), f1y(u); f2x(u), f2y(u) ];
  98 tol = 1 / 10^12;
  99 u0 = [1.2; 2.5];
 100 expected_root = [1.33635; 1.75424];
 101 [actual_root, iterations] = newtons_method(f, f_prime, tol, u0);
 102 unit_test_equals("Homework #3 problem #4 root is correct", ...
 103                  expected_root, ...
 104                  actual_root);