yacas 1.3.6-2+b1 / usr / share / yacas / scripts / integrate.rep / code.ys

RuleBase("IntegrateMultiplicative",{var,from,a,b});
UnFence("IntegrateMultiplicative",4);

10# (Integrate(_var)(expr_IsList)) 
    <-- Map("Integrate",{FillList(var,Length(expr)),expr});
20 # (Integrate(_var)(_expr)) <-- IntSub(var,expr,AntiDeriv(var,IntClean(var,expr)));


10 # IntSub(_var,_expr,Integrate(_var)(_expr2)) <--
     `Hold(Integrate(@var)(@expr));
20 # IntSub(_var,_expr,_result) <-- result;	// + UniqueConstant();

////////////////////////////////////////////////
//
// Integrate over a range
//
////////////////////////////////////////////////
10# (Integrate(_var,_from,_to)(expr_IsList)) 
    <-- Map("Integrate",{FillList(var,Length(expr)),
                         FillList(from,Length(expr)),
			 FillList(to,Length(expr)),
                         expr});

20 # (Integrate(_var,_from,_to)(_expr))
    <-- indefIntegrate(var,from,to,expr,a,b);

////////////////////////////////////////////////
//
// separate rules can be added here for specific integrals
// to indefIntegrate
//
////////////////////////////////////////////////

10 # indefIntegrate(_var,_from,_to,_expr,_a,_b)_(from = -to And IsOddFunction(expr,var)) <-- 0;

// We need to define this case (integrating from 0 to 0 over an even function) explicitly, otherwise
// the integration ends up going in to infinite recursion. Extended it a little bit more, since if
// you are integrating from A to A, then the result is obviously zero. There are perhaps situations
// where this does not work, where we need to simplify (to-from) first. A naive implementation caused
// a test to fail.
10 # indefIntegrate(_var,_from,_from,_expr,_a,_b) <-- 0;

12 # indefIntegrate(_var,_from,_to,_expr,_a,_b)_(from = -to And IsEvenFunction(expr,var)) <--
     2*indefIntegrate(var,0,to,expr,a,b);

100 # indefIntegrate(_var,_from,_to,_expr,_a,_b)_(Type(AntiDeriv(var,IntClean(var,expr))) != "AntiDeriv") 
    <--
    IntegrateRange(var,expr,from,to,AntiDeriv(var,IntClean(var,expr)));
101 # indefIntegrate(_var,_from,_to,_expr,_a,_b)
    <-- `Hold(Integrate(@var,@from,@to)(@expr));
//    <-- IntegrateRange(var,expr,from,to,AntiDeriv(var,expr));



////////////////////////////////////////////////
//
// No anti-derivative found, return unavaluated.
//
////////////////////////////////////////////////
10 # IntegrateRange(_var,_expr,_from,_to,Integrate(_var)_expr2)
    <-- `Hold(Integrate(@var,@from,@to)@expr);
    
////////////////////////////////////////////////
//
// Anti-derivative found, return result.
//
////////////////////////////////////////////////
20 # IntegrateRange(_var,_expr,_from,_to,_antideriv)
    <-- `(@antideriv Where @var == @to) - `(@antideriv Where @var == @from);

////////////////////////////////////////////////
//
// IntClean cleans up an expression before passing
// it on to integration. This function normalizes
// an expression in a way desirable for integration.
// TrigSimpCombine, for instance, expands expressions
// containing trigonometric functions so that they are
// additive as opposed to multiplicative.
//
// If the expression doesn't contain the variable,
// just return it as-is. This fixes: 
// In> Integrate(x) z^100
////////////////////////////////////////////////
10 # IntClean(_var,_expr) <-- 
[
	if( IsFreeOf(var,expr))[
		expr;
	] else if ( HasFunc(expr,Sin) Or HasFunc(expr,Cos) )[
		Simplify(TrigSimpCombine(expr));
	] else [
		Simplify(expr);
	];
];

////////////////////////////////////////////////
//
// Anti-derivative of a univariate polynomial
//
////////////////////////////////////////////////
5  # AntiDeriv(_var, poly_CanBeUni(var) )
     <-- NormalForm(AntiDeriv(var,`MakeUni(@poly,@var)));
5 # AntiDeriv(_var,UniVariate(_var,_first,_coefs)) <--
[
  Local(result,i);
  result:=FlatCopy(coefs);
  For(i:=1,i<=Length(result),i++)
  [
    result[i]:= result[i]/(first+i);
  ];
  UniVariate(var,first+1,result);
];



////////////////////////////////////////////////
//
// Standard additive properties of integration.
//
////////////////////////////////////////////////
10 # AntiDeriv(_var,_x + _y) <-- AntiDeriv(var,x) + AntiDeriv(var,y);
10 # AntiDeriv(_var,_x - _y) <-- AntiDeriv(var,x) - AntiDeriv(var,y);
10 # AntiDeriv(_var,   - _y) <--            - AntiDeriv(var,y);

10 # AntiDeriv(_var,_x/c_IsFreeOf(var) )_(HasExpr(x,var)) <-- AntiDeriv(var,x)/c;
10 # AntiDeriv(_var,c_IsFreeOf(var)/_x )_(HasExpr(x,var) And c!= 1)
	<-- c*AntiDeriv(var,1/x);


////////////////////////////////////////////////
//
// Multiplying a polynomial with another (integrable)
// function, Integrate by parts.
//
////////////////////////////////////////////////
1570 # IntegrateMultiplicative(_var,(exy_CanBeUni(var)) * _exx,_dummy1,_dummy2)
     <-- IntByParts(var,exy*exx,AntiDeriv(var,exx));
1570 # IntegrateMultiplicative(_var,_exx * (exy_CanBeUni(var)),_dummy1,_dummy2)
     <-- IntByParts(var,exy*exx,AntiDeriv(var,exx));
10 # IntByParts(_var,_exy * _exx,Integrate(_var)(_something)) <--
     `Hold(AntiDeriv(@var,((@exy)*(@exx))));
20 # IntByParts(_var,_exy * _exx,_anti)_(Not IsFreeOf(anti,exx)) <--
     `Hold(AntiDeriv(@var,((@exy)*(@exx))));
30 # IntByParts(_var,_exy * _exx,_anti) <--
     [
       Local(cf);
       cf:=anti*Deriv(var)exy;
//  Echo({exy*anti,exy*exx,cf});
       exy*anti - `(AntiDeriv(@var,@cf));
     ];

////////////////////////////////////////////////
//
// Rational functions: f(x)/g(x) where f and g are
// polynomials.
//
////////////////////////////////////////////////
1570 # IntegrateMultiplicative(_var,(exy_CanBeUni(var)) / (exx_CanBeUni(var)),_dummy1,_dummy2) <--
     IntRat(var,exy/exx,MakeUni(exy,var),MakeUni(exx,var));

10 # IntRat(_var,_exy / _exx,_exyu,_exxu)_
     (Degree(exyu) > Degree(exxu) Or Degree(Gcd(exyu,exxu)) > 0) <--
     [
     Local(gcd);
     gcd:=Gcd(exxu,exyu);
     exyu:=Div(exyu,gcd);
     exxu:=Div(exxu,gcd);
     AntiDeriv(var,NormalForm(Div(exyu,exxu))) +
       AntiDeriv(var,NormalForm(Mod(exyu,exxu))/NormalForm(exxu));
     ];

11 # IntRat(_var,_exy / _exx,_exyu,_exxu)_
        (Degree(exxu,var) > 1 And LeadingCoef(exxu)=1 And
	 IsNumericList(Coef(exxu,var,0 .. Degree(exxu)))) <--
[
  Local(ee);
  ee:=Apart(exy/exx,var);
  `AntiDeriv(@var,@ee);
];


20 # IntRat(_var,_exy / _exx,_exyu,_exxu) <--
     `Hold(AntiDeriv(@var,((@exy)/(@exx))));


30 # AntiDeriv(_var,Deriv(_var)(_expr)) <-- expr;

////////////////////////////////////////////////
//
// No simple form, try something else
//
////////////////////////////////////////////////
100 # AntiDeriv(_var,_exp) <--
[
  IntegrateMultiplicative(var,exp,a,b);
];


////////////////////////////////////////////////
//
// Special anti-derivatives can be added here.
//
////////////////////////////////////////////////

// integrating expressions containing if:
10 # IntegrateMultiplicative(_var,if(_cond)(_body),_a,_b)
     <--
     [
       body := AntiDeriv(var,body);
       `Hold(if(@cond)(@body));
     ];
// integrating expressions containing else
10 # IntegrateMultiplicative(_var,(_left) else (_right),_a,_b)
     <--
     [
       left  := AntiDeriv(var,left);
       right := AntiDeriv(var,right);
       `Hold( (@left) else (@right) );
     ];


////////////////////////////////////////////////
//
// Could not find anti-derivative, return unsimplified
//
////////////////////////////////////////////////
1600 # IntegrateMultiplicative(_var,_exp,_a,_b) <-- `Hold(Integrate(@var)(@exp));

////////////////////////////////////////////////
//
// IntFunc declares the anti-derivative of a function
// that has one argument.
// Calling sequence: IntFunc(variable,from,to);
// Example: IntFunc(x,Cos(_x),Sin(x));
//
////////////////////////////////////////////////
LocalSymbols(intpred)
[
  intpred := 50;
  IntFunc(_vr,_from,_to) <--
  [
    `((@intpred) # IntegrateMultiplicative(_var,@from,_dummy1,_dummy2)_MatchLinear(var,@vr) <-- (@to)/Matched'a());
    intpred++;
  ];
];


IntPureSquare(_vr,_from,_sign2,_sign0,_to) <--
[
  `(50 # IntegrateMultiplicative(_var,@from,_dummy1,_dummy2)_MatchPureSquared(var,@sign2,@sign0,@vr) <-- (@to));
];




////////////////////////////////////////////////
//
// Declaration of the anti-derivatives of a few analytic functions
//
////////////////////////////////////////////////


IntFunc(x,Sqrt(_x),(2*Sqrt(x)^(3))/3);
IntFunc(x,1/Sqrt(_x),2*Sqrt(x));
IntFunc(x,1/_x^(_n),x^(1-n)/(1-n) );
IntFunc(x,Sin(_x),-Cos(x));
IntFunc(x,1/Sin(_x), Ln( 1/Sin(x) - Cos(x)/Sin(x) ) );
IntFunc(x,Cos(_x),Sin(x));
IntFunc(x,1/Cos(_x),Ln(1/Cos(x)+Tan(x)));
IntFunc(x,Tan(_x),-Ln(Cos(x)));
IntFunc(x,1/Tan(_x),Ln(Sin(x)) );
IntFunc(x,Cos(_x)/Sin(_x),Ln(Sin(x)));
IntFunc(x,Exp(_x),Exp(x));
IntFunc(x,(C_IsFreeOf(var))^(_x),C^x/Ln(C));
// we don't need Ln(Abs(x))
IntFunc(x,num_IsFreeOf(var) / (_x),num*Ln(x));
IntFunc(x,Ln(_x),x*Ln(x)-x);
// where did these 1+1's come from?
IntFunc(x,(_x)*Ln(_x),(1/(1+1))*x^(1+1)*Ln(x) - (1/(1+1)^2)*x^(1+1) );
IntFunc(x,Ln(_x)*(_x),(1/(1+1))*x^(1+1)*Ln(x) - (1/(1+1)^2)*x^(1+1) );

IntFunc(x,1/Sin(_x)^2,-Cos(x)/Sin(x) );
IntFunc(x,1/Cos(_x)^2,Tan(x) );
IntFunc(x,1/(Sin(_x)*Tan(_x)),-1/Sin(x));
IntFunc(x,Tan(_x)/Cos(_x),1/Cos(x));
IntFunc(x,1/Sinh(_x)^2,-1/Tanh(x));
IntFunc(x,1/Cosh(_x)^2,Tanh(x));
IntFunc(x,1/(Sinh(_x)*Tan(_x)),-1/Sinh(x));
IntFunc(x,Tanh(_x)/Cosh(_x),-1/Cosh(x));

IntFunc(x,1/Sqrt(m_IsFreeOf(x)-_x^2),ArcSin(x/Sqrt(m)) );

IntFunc(x,Exp(n_IsNumber*_x)*Sin(m_IsNumber*_x),Exp(n*x)*(n*Sin(m*x)- m*Cos(m*x))/(m^2+n^2) );

// n>0
IntFunc(x,Ln(_x)*(_x)^n_IsNumber,(1/(n+1))*x^(n+1)*Ln(x) - (1/(n+1)^2)*x^(n+1) );

// n>0
IntFunc(x,Ln(A_IsNumber*_x)*(_x)^n_IsNumber,(1/(n+1))*x^(n+1)*Ln(A*x) - (1/(n+1)^2)*x^(n+1) );

IntFunc(x,Sin(Ln(_x)),x*Sin(Ln(x))/2 - x*Cos(Ln(x))/2 );
IntFunc(x,Cos(Ln(_x)),x*Sin(Ln(x))/2 - x*Cos(Ln(x))/2 );

IntFunc(x,1/((_x)*Ln(_x)),Ln(Ln(x)));

IntFunc(x,(_x)^(-1),Ln(x));

IntFunc(x,(_x)^(n_IsFreeOf(x)),x^(n+1)/(n+1));
IntFunc(x,Sinh(_x),Cosh(x));
IntFunc(x,Sinh(_x)^2,Sinh(2*x)/4 - x/2);
IntFunc(x,1/Sinh(_x),Ln(Tanh(x/2)));
IntFunc(x,Cosh(_x),Sinh(x));
IntFunc(x,Cosh(_x)^2,Sinh(2*x)/4 + x/2);
IntFunc(x,1/Cosh(_x),ArcTan(Sinh(x)));
IntFunc(x,Tanh(_x),Ln(Cosh(x)));
IntFunc(x,Tanh(_x)/Cosh(_x),-1/Cosh(x));
IntFunc(x,1/Cosh(_x)^2,Tanh(x));
//IntFunc(x,1/Sech(_x)*Coth(_x),-1/Sinh(x));
IntFunc(x,1/Tanh(_x),Ln(Sinh(x)));

IntFunc(x,Abs(_x),Abs(x)*x/2);	// not 2*a

IntFunc(x,ArcTan(_x),x*ArcTan(x) - Ln(x^2 + 1)/2);
//IntFunc(x,ArcSin(_x),(x*ArcSin(x)) + Sqrt(1-x^2) );
IntFunc(x,ArcCos(_x),x*ArcCos(x) - Sqrt(1-x^2) );

IntFunc(x,ArcTanh(_x),x*ArcTanh(x) + Ln(1-x^2)/2 );
IntFunc(x,ArcSinh(_x),x*ArcSinh(x) - Sqrt(x^2 + 1) );
IntFunc(x,ArcCosh(_x),x*ArcCosh(x) - Sqrt(x-1)*Sqrt(x+1) );


// n^2 > x^2 
//IntFunc(x,num_IsFreeOf(var)/(-(_x)^2 + n_IsNumber),num*ArcTanh(x/Sqrt(n))/n);

// x^2 > n^2
//IntFunc(x,num_IsFreeOf(var)/((_x)^2 - n_IsNumber),num * -ArcCoth(x/Sqrt(n))/Sqrt(n));

// n^2 > x^2
//IntFunc(x,num_IsFreeOf(var)/Sqrt(n_IsNumber - (_x)^2),num*ArcSin(x/Sqrt(n)));

// previous code is killing this....
IntFunc(x,num_IsFreeOf(var)/(A_IsNumber + B_IsNumber*(_x))^2,-num/(A*b + B^2*x));

// Code works now?
IntFunc(x,num_IsFreeOf(var)/(n_IsNumber + m_IsNumber*Exp(p_IsNumber*(_x))),num*x/n - num*Ln(n + m*Exp(p*x))/(n*p));
IntFunc(x,num_IsFreeOf(var)/(m_IsNumber*Exp(p_IsNumber*(_x)) + n_IsNumber),num*x/n - num*Ln(n + m*Exp(p*x))/(n*p));

IntPureSquare(x,num_IsFreeOf(var)/(_x),1,1,(num/(Sqrt(Matched'b()*Matched'a())))*ArcTan(var/Sqrt(Matched'b()/Matched'a())));



///// Integrating Special Functions
IntFunc(x,Erf(_x), x*Erf(x)+ 1/(Exp(x^2)*Sqrt(Pi)) );