[Maxima] Help with compile_file()

[Joshua Stults]

Thanks; using your fix of replacing ev() with subst() and adding f to
either the argument list (as in your example), or to the local
variable list of the block fixes things, no warnings or errors.  I
guess it's a matter of preference where you declare f, I tend not to
know my f's analytically, so no need to pass it as an argument.

explicit_fde(steps, deriv) := block([A, a_f, b, coef_list, fde,
func_vals, f_coefs, f_taylor, M, n, points, i, f],
  n : length(steps),
  f_taylor : taylor("f"(x), x, x0, n),
  points : [x=x0+steps[1]*h],
  for i : 2 thru n do
  (
    points : append(points, [x=x0+steps[i]*h])
    ),
  /* right-hand side vector (based on the derivative requested): */
  b : zeromatrix(n,1),
  b[deriv+1,1] : 1,
  /* function values vector: */
  func_vals : matrix([at(f(x), points[1])]),
  for i : 2 thru n do
  (
    func_vals : addcol(func_vals, [at(f(x), points[i])])
    ),
  /* coefficients list */
  coef_list : [at(f(x),x=x0)],
  for i : 2 thru n do
  (
    coef_list : append(coef_list, [ 'at('diff(f(x),x,i-1), x=x0) ] )
    ),
  f_coefs : submatrix(factor(augcoefmatrix([f_taylor], coef_list)), n+1),
  /* system matrix: */
  A : addcol(transpose(ev(f_coefs, points[1]))),
  for i : 2 thru n do
  (
    A : addcol(A, transpose(subst(points[i],f_coefs)))
    ),
  M : lu_factor(A),
  /* solve for the coefficients: */
  a_f : fullratsimp( lu_backsub(M, b) ),
  /* finite difference expression: */
  fde : fullratsimp(func_vals.a_f),
  return(fde)
  );

2009/5/22 Martin Schönecker <ms_usenet at gmx.de>:
> Joshua Stults schrieb:
>>
>> Hello,
>>
>> I'm trying to learn how to use compile_file(), but am having a bit of
>> trouble.
>
> I have no experience with compliling either, but let's see: Firstly, it
> seems that the compiler cannot handle the undefined function f, however
> doing as maxima proposes:
>
>> Warning: f has a function or macro call which has not been translated
>> properly.
>> The function was totaly undefined. Maybe you want to quote it.
>
> ...doesn't seem to work.  So we try to give the function name as a third
> argument (There may be more elegant approaches?).  Then, i does not seem to
> be evaluated to the number it is ought to have in the for loop.  This seems
> to be due to the evaluation partialities of ev() (of which I'm unaware).  It
> seems that the correct position of points should be substituted into
> f_coefs, so we'll try to use subst().
>
> Now the following compiles and seems to work properly (om Maxima 5.18.1):
>
> explicit_fde(steps, deriv, f) := block([A, a_f, b, coef_list, fde,
> func_vals, f_coefs, f_taylor, M, n, points, i],
>  n : length(steps),
>  f_taylor : taylor(f(x), x, x0, n),
>  points : [x=x0+steps[1]*h],
>  for i : 2 thru n do
>  (
>    points : append(points, [x=x0+steps[i]*h])
>    ),
>  /* right-hand side vector (based on the derivative requested): */
>  b : zeromatrix(n,1),
>  b[deriv+1,1] : 1,
>  /* function values vector: */
>  func_vals : matrix([at(f(x), points[1])]),
>  for i : 2 thru n do
>  (
>    func_vals : addcol(func_vals, [at(f(x), points[i])])
>    ),
>  /* coefficients list */
>  coef_list : [at(f(x),x=x0)],
>  for i : 2 thru n do
>  (
>    coef_list : append(coef_list, [ 'at('diff(f(x),x,i-1), x=x0) ] )
>    ),
>  f_coefs : submatrix(factor(augcoefmatrix([f_taylor], coef_list)), n+1),
>  print(f_coefs),
>  /* system matrix: */
>  A : addcol(transpose(ev(f_coefs, points[1]))),
>  for i : 2 thru n do
>  (
>    A : addcol(A, transpose(subst(points[i], f_coefs)))
>    ),
>  M : lu_factor(A),
>  /* solve for the coefficients: */
>  a_f : fullratsimp( lu_backsub(M, b) ),
>  /* finite difference expression: */
>  fde : fullratsimp(func_vals.a_f),
>  return(fde)
>  );
>
>
> (%i1) kill(all);
> (%o0) done
> (%i1) compile_file("m:/maxima/explicit_fde.mac");
> Translation begun on #pm:/maxima/explicit_fde.mac.
> See the `unlisp' file for possible optimizations.
> Compiling m:/maxima/explicit_fde.LISP.
> End of Pass 1.
> End of Pass 2.
> OPTIMIZE levels: Safety=2, Space=3, Speed=3
> Finished compiling NIL.
> (%o1)
> [m:/maxima/explicit_fde.mac,m:/maxima/explicit_fde.LISP,m:/maxima/explicit_fde.UNLISP,#pm:/maxima/explicit_fde.o]
> (%i2) explicit_fde([-2, -1,0,1,2], 1, f);
> matrix([1,-(x0-x),(x0-x)^2/2,-(x0-x)^3/6,(x0-x)^4/24])
> Is  h  positive, negative, or zero?p;
> (%o2) -(f(x0+2*h)-8*f(x0+h)+8*f(x0-h)-f(x0-2*h))/(12*h)
>



-- 
Joshua Stults
Website: http://j-stults.blogspot.com