/*  -*-  Mode: Maxima -*- */

/*
**
** Copyright (C) 2015 Marek Rychlik <rychlik@email.arizona.edu>
**
** This program is free software; you can redistribute it and/or modify
** it under the terms of the GNU General Public License as published by
** the Free Software Foundation; either version 2 of the License, or
** (at your option) any later version.
**
** This program is distributed in the hope that it will be useful,
** but WITHOUT ANY WARRANTY; without even the implied warranty of
** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
** GNU General Public License for more details.
**
** You should have received a copy of the GNU General Public License
** along with this program; if not, write to the Free Software
** Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
**
*/

kill(all);
load("linode.mac");

/* This function verifies that linode solves the equation correctly */
linode_check(de,y,x,[bconds]):=block([sub:apply(linode,append([de,y,x],bconds)),e],
  local(y),
  e:subst(sub,de),
  e:ev(e,[diff]),
  e:ratsimp(exponentialize(expand(e))),
  if not is(e) then error("Equation",de,"not solved correctly.") else print("[OK]")
  );


/* Examples of handled equations */
define_variable(linode_examples,[],any_check);

linode_examples:[
  ['diff(y,x,2)-3*'diff(y,x)+'y=0,y,x],
  ['diff(y,x,2)+'y=0,y,x],
  ['diff(y,x,2)+'y=cos('x),y,x],
  ['diff(y,x,3)+3*'diff(y,x,2)+3*'diff(y,x)+'y=0,y,x],
  ['diff(y,x,3)-4*'diff(y,x,2)+5*'diff(y,x)-2*'y=0,y,x],
  ['diff(y,x,2)-a*'diff(y,x)+'y=0,y,x],
  [m*'diff(y,x,2)+b*'diff(y,x)+'k*'y=0,y,x],
  [x^2*'diff(y,x,2)-3*'x*'diff(y,x)+y=0,y,x],
  [x^3*'diff(y,x,3)-3*x^2*'diff(y,x,2)+6*x*'diff(y,x)-6*'y=x^4*log(x),y,x],
  [x*'diff(y,x,4)+'diff(y,x,3)=63*x^4,y,x],
  [16*x^2*'diff(y,x,3)+48*x*'diff(y,x,2)+15*'diff(y,x)=0,y,x],
  [x^2*'diff(y,x,2)-x*'diff(y,x)+y=x^3*log(x),y,x],
  ['diff(y,x,3)+'y='x,y,x],
  /* Simple IVP */
  ['diff(y,x,2)+'y=0,y,x,x=0,[y=1,'diff(y,x)=2]],
  /* A more complex IVP */
  ['diff(y,x,3)+y='x,y,x,x=0,[y=1,'diff(y,x)=2,'diff(y,x,2)=3]],
  /* Check use of alternative variables */
  ['diff(z,t,2)+2*'diff(z,t)+z=cos(t),z,t],
  ['diff(y,x,2)-2*'diff(y,x)+y=x^100*exp(x),y,x],
  ['diff(y,x,2)-2*'diff(y,x)+y=x^100*exp(2*x),y,x],
  ['diff(y,x,3)-3*'diff(y,x,2)+3*'diff(y,x)-y=x^100*exp(2*x),y,x],
  /* A simple Boundary value problem */
  ['diff(y,x,2)+'y=0,y,x,x=0,[y=1],x=1,[y=2]],
  /* Boundary value problem for a deflected beam */
  ['diff(y,x,4)=-1,y,x,x=0,[y=0,'diff(y,x)=0],x=1,[y=0,'diff(y,x)=0]],
  ['diff(y,x,4)=-1,y,x,x=0,[y=0,'diff(y,x,2)=0],x=1,[y=0,'diff(y,x,2)=0]],
  /* Test non-zero coefficient on high derivative */
  [EI*'diff(y,x,4)=f0,y,x,x=0,[y=0,'diff(y,x,2)=0],x=1,[y=0,'diff(y,x,2)=0]],
  /* Test symbolic right-hand side */
  [EI*'diff(y,x,4)=f(x),y,x,x=0,[y=0,'diff(y,x,2)=0],x=1,[y=0,'diff(y,x,2)=0]],
  /* Using y(x) in place of y - compatible with 'desolve' */
  ['diff(y(x),x,3)+3*'diff(y(x),x,2)+3*'diff(y(x),x)+y(x)=-30*exp(-x),y(x),x, x=0,[y(x)=3,'diff(y(x),x)=-3,'diff(y(x),x,2)=-47]]
  ];

test_linode():=block([z,t],
  local(z,t),
  /* Do the list of examples */
  for example in linode_examples do (
    apply(linode_check,example)
    ),
  return('done));

test_linode();