function coefficients = forward_euler(integer_order, xs, x)
- ##
- ## Return the coefficients of u(x0), u(x1), ..., u(xn) as a vector.
- ## Take for example a first order approximation, with,
- ##
- ## xs = [x0,x1,x2,x3,x4]
- ##
- ## f'(x1) ~= [f(x2)-f(x1)]/(x2-x1)
- ##
- ## This would return [0, -1/(x2-x1), 2/(x2-x1), 0, 0]. This aids the
- ## solution of linear systems.
- ##
- ##
- ## INPUTS:
- ##
- ## * ``integer_order`` - The order of the derivative which we're
- ## approximating.
- ##
- ## * ``xs`` - The vector of x-coordinates.
- ##
- ## * ``x`` - The point `x` at which you'd like to evaluate the
- ## derivative of the specified `integer_order`. This should be an
- ## element of `xs`.
- ##
- ##
- ## OUTPUTS:
- ##
- ## * ``coefficients`` - The vector of coefficients, in order, of
- ## f(x0), f(x1), ..., f(xn).
- ##
+ %
+ % Return the coefficients of u(x0), u(x1), ..., u(xn) as a vector.
+ % Take for example a first order approximation, with,
+ %
+ % xs = [x0,x1,x2,x3,x4]
+ %
+ % f'(x1) ~= [f(x2)-f(x1)]/(x2-x1)
+ %
+ % This would return [0, -1/(x2-x1), 2/(x2-x1), 0, 0]. This aids the
+ % solution of linear systems.
+ %
+ %
+ % INPUTS:
+ %
+ % * ``integer_order`` - The order of the derivative which we're
+ % approximating.
+ %
+ % * ``xs`` - The vector of x-coordinates.
+ %
+ % * ``x`` - The point `x` at which you'd like to evaluate the
+ % derivative of the specified `integer_order`. This should be an
+ % element of `xs`.
+ %
+ %
+ % OUTPUTS:
+ %
+ % * ``coefficients`` - The vector of coefficients, in order, of
+ % f(x0), f(x1), ..., f(xn).
+ %
if (integer_order < 0)
- ## You have made a grave mistake.
+ % You have made a grave mistake.
coefficients = NA;
return;
end
end
if (length(xs) < 2)
- ## You can't approximate a derivative of order greater than zero
- ## with zero or one points!
+ % You can't approximate a derivative of order greater than zero
+ % with zero or one points!
coefficients = NA
return;
end
offset_b = integer_order / 2;
offset_f = offset_b;
else
- ## When the order is odd, we need one more "forward" point than we
- ## do "backward" points.
+ % When the order is odd, we need one more "forward" point than we
+ % do "backward" points.
offset_b = (integer_order - 1) / 2;
offset_f = offset_b + 1;
end
- ## Zero out the coefficients for terms that won't appear. We compute
- ## where `x` is, and we just computed how far back/forward we need to
- ## look from `x`, so we just need to make the rest zeros.
+ % Zero out the coefficients for terms that won't appear. We compute
+ % where `x` is, and we just computed how far back/forward we need to
+ % look from `x`, so we just need to make the rest zeros.
x_idx = find(xs == x);
first_nonzero_idx = x_idx - offset_b;
last_nonzero_idx = x_idx + offset_f;
targets = xs(first_nonzero_idx : last_nonzero_idx);
- # The multiplier comes from the Taylor expansion.
+ % The multiplier comes from the Taylor expansion.
multiplier = factorial(integer_order);
cs = divided_difference_coefficients(targets) * multiplier;