X-Git-Url: http://gitweb.michael.orlitzky.com/?p=octave.git;a=blobdiff_plain;f=poisson_matrix.m;h=3ed1c0c673948c5ae44217cb78d0d45559335956;hp=696185ff0e6418649aecb60d47e05d274e9b33f3;hb=c40bbb2f7d9073193f107bbe252d4e9c8509a203;hpb=6de25e4c4273789fd2f98e309db4ef3f0902b185

diff --git a/poisson_matrix.m b/poisson_matrix.m
index 696185f..3ed1c0c 100644
--- a/poisson_matrix.m
+++ b/poisson_matrix.m
@@ -1,37 +1,37 @@
 function A = poisson_matrix(integerN, x0, xN)
-  ##
-  ## In the numerical solution of the poisson equation,
-  ##
-  ##   -u''(x) = f(x)
-  ##
-  ## in one dimension, subject to the boundary conditions,
-  ##
-  ##   u(x0) = 0
-  ##   u(xN) = 1
-  ##
-  ## over the interval [x0,xN], we need to compute a matrix A which
-  ## is then multiplied by the vector of u(x0), ..., u(xN). The entries
-  ## of A are determined from the second order finite difference formula,
-  ## i.e. the second order forward Euler method.
-  ##
-  ## INPUTS:
-  ##
-  ##   * ``integerN`` - An integer representing the number of
-  ##     subintervals we should use to approximate `u`. Must be
-  ##     greater than or equal to 2, since we have at least two
-  ##     values for u(x0) and u(xN).
-  ##
-  ##   * ``f`` - The function on the right hand side of the poisson
-  ##     equation.
-  ##
-  ##   * ``x0`` - The initial point.
-  ##
-  ##   * ``xN`` - The terminal point.
-  ##
-  ## OUTPUTS:
-  ##
-  ##   * ``A`` - The (N+1)x(N+1) matrix of coefficients for u(x0),
-  ##             ..., u(xN).
+  %
+  % In the numerical solution of the poisson equation,
+  %
+  %   -u''(x) = f(x)
+  %
+  % in one dimension, subject to the boundary conditions,
+  %
+  %   u(x0) = 0
+  %   u(xN) = 1
+  %
+  % over the interval [x0,xN], we need to compute a matrix A which
+  % is then multiplied by the vector of u(x0), ..., u(xN). The entries
+  % of A are determined from the second order finite difference formula,
+  % i.e. the second order forward Euler method.
+  %
+  % INPUTS:
+  %
+  %   * ``integerN`` - An integer representing the number of
+  %     subintervals we should use to approximate `u`. Must be
+  %     greater than or equal to 2, since we have at least two
+  %     values for u(x0) and u(xN).
+  %
+  %   * ``f`` - The function on the right hand side of the poisson
+  %     equation.
+  %
+  %   * ``x0`` - The initial point.
+  %
+  %   * ``xN`` - The terminal point.
+  %
+  % OUTPUTS:
+  %
+  %   * ``A`` - The (N+1)x(N+1) matrix of coefficients for u(x0),
+  %             ..., u(xN).
 
   if (integerN < 2)
     A = NA
@@ -40,26 +40,27 @@ function A = poisson_matrix(integerN, x0, xN)
 
   [xs,h] = partition(integerN, x0, xN);
 
-  ## We cannot evaluate u_xx at the endpoints because our
-  ## differentiation algorithm relies on the points directly to the left
-  ## and right of `x`.
+  % We cannot evaluate u_xx at the endpoints because our
+  % differentiation algorithm relies on the points directly to the left
+  % and right of `x`.
   differentiable_points = xs(2:end-1);
 
-  ## These are the coefficient vectors for the u(x0) and u(xn)
-  ## constraints. There should be N zeros and a single 1.
+  % These are the coefficient vectors for the u(x0) and u(xn)
+  % constraints. There should be N zeros and a single 1.
   the_rest_zeros = zeros(1, integerN);
   u_x0_coeffs = cat(2, 1, the_rest_zeros);
   u_xN_coeffs = cat(2, the_rest_zeros, 1);
 
-  ## Start with the u(x0) row.
+  % Start with the u(x0) row.
   A = u_x0_coeffs;
 
   for x = differentiable_points
-    ## Append each row obtained from the forward Euler method to A.
+    % Append each row obtained from the forward Euler method to A.
     u_row = forward_euler(2, xs, x);
     A = cat(1, A, u_row);
   end
 
-  ## Finally, append the last row for xN.
-  A = cat(1, A, u_xN_coeffs);
+  % Finally, append the last row for xN and negate the whole thing (see
+  % the definition of the poisson problem).
+  A = - cat(1, A, u_xN_coeffs);
 end