/*
Functions for the calculation of lattice Green functions in Maxima.

Copyright (C) 2005 Exstrom Laboratories LLC
This work is freely redistributable under the terms of the GNU General Public License as
published by the Free Software Foundation. This work comes with ABSOLUTELY NO WARRANTY.
*/

s(n,m):=
block([f,t1,t2,t3],
if n < 0 then return(0),
if evenp(n) then f:1 else f:(n+1)*(n+1+m),
t3:f,
for j:fix(n/2) thru 1 step -1 do
(
  t1:n - j + 1,
  t2:t1 - j,
  f:f * 4 * j * t1 * (2 * t1 + m) * (m + 2 * j) / (t2 * (t2 + 1))^2,
  t3:t3 + f
),
return(t3)
)$

hodd(n):=
block([j],
if n=0 then return(0),
if n=1 then return(1),
if n>1 then return(sum(1/(2*j-1),j,1,n)),
if n<0 then return(1-sum(1/(2*j-1),j,1,1-n))
)$

/* g(n,m) calculates the square lattice Green function value. Created July 18, 2005 */
g(n,m):=
block([t1,t2],
if n=m then return(hodd(n)/%PI) else (
t1: s(n-m-1,2*m) / 4,
t2: sum(((-1)^j*s(n-m-1,2*m-2*j+1)+(-1)^m)/(2*j-1),j,1,m)/%PI,
t3: sum((-1)^j*s(n+m-2*j,2*j-2*m-1)/(2*j-1),j,m+1,fix((n+m)/2))/%PI,
t4: if m=0 then t3 else sum((-1)^j*s(n-m-2*j,2*j+2*m-1)/(2*j-1),j,1,fix((n-m)/2))/%PI,
return((-1)^m*(t1+t2+t3+t4))
)
)$