Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I created a simple numeric-symbolic integration. Here you can use symbolical and numerical techniques at the same time. You can also interpolate numerical integrals.

The problem with my function is the use of TimeConstrained for testing if the integral can be done fully symbolically. This wastes time; without it the result is much much faster.

How can it be optimized? Maybe there is some other implementation like this I haven't found? Note: It can easily be extended to Sum-s, just exchange Integrate for Sum :)

Options@NSIntegrate = {
    InterpolationVariable -> {}
};

(* Linearity *)
NSIntegrate[f_ + g__, s_, options___] := NSIntegrate[f, s, options] + NSIntegrate[+g, s, options]
NSIntegrate[α_*f_, {x_, limits__}, options___] /; FreeQ[α, x] := α*NSIntegrate[f, {x, limits}, options]

(* Constant Function *)
NSIntegrate[α_, {x_, limits__}, options___] /; FreeQ[α, x] \[And] α =!= 1 := α*NSIntegrate[1, {x, limits}, options]

(* Exponentials *)
NSIntegrate[α_^(f_+g__), {x_, limits__}, options___] /; FreeQ[f, x] := α^f*NSIntegrate[α^g, {x, limits}, options]
NSIntegrate[α_^(f_+g__)*h_, {x_, limits__}, options___] /; FreeQ[f, x] := α^f*NSIntegrate[α^g*h, {x, limits}, options]

(* Some Basic Integrals - Polynomials *)
NSIntegrate[1, {x_, a_, b_}, ___] := b - a
NSIntegrate[x_, {x_, a_, b_}, ___] := (b^2 - a^2)/2
NSIntegrate[x_^-1, {x_, a_, b_}, ___] /; If[NumericQ@a \[And] NumericQ@b, Sign@a === Sign@b, True] \[And] a =!= b := Log[b/a]
NSIntegrate[x_^n_, {x_, a_, b_}, ___] /; n =!= -1 := (b^(n+1) - a^(n+1))/(n + 1)

(* Some Basic Integrals - Exponentials *)
NSIntegrate[α_^(A_*x_), {x_, a_, b_}, ___] /; α =!= 0 \[And] α =!= 1 := (α^(A*b) - α^(A*a))/(A*Log@α)
NSIntegrate[α_^x_, {x_, a_, b_}, ___] /; α =!= 0 \[And] α =!= 1 := (α^b - α^a)/Log@α
NSIntegrate[0^x_, s__, ___] := 0

(* Scaling and offset *)
NSIntegrate[f_[α_*x_ + β_, args___], {x_, a_, b_}, options___] /; (And @@ (FreeQ[#, x] & /@ Flatten@{args, α, β})) := 
    NSIntegrate[f[x, args]/α, {x, a*α + β, b*α + β}, options]
NSIntegrate[f_[α_*x_, args___], {x_, a_, b_}, options___] /; (And @@ (FreeQ[#, x] & /@ Flatten@{args, α})) := 
    NSIntegrate[f[x, args]/α, {x, a*α, b*α}, options]
NSIntegrate[f_[x_ + β_, args___], {x_, a_, b_}, options___] /; (And @@ (FreeQ[#, x] & /@ Flatten@{args, β})) := 
NSIntegrate[f[x, args], {x, a + β, b + β}, options]

(* Interpolation *)
NSIntegrate[f_, {x_, a_, b_}, options1___, InterpolationVariable -> varRange_List, options2___] /; Length@varRange != 0 := Module[{
    vars, range = varRange, table
    },
    If[Depth@varRange == 2, vars = varRange[[{1}]]; range = List@range, vars = varRange[[All, 1]]];

    table = Table[{vars, NSIntegrate[f, {x, a, b}, options1, options2]}, ##1] & @@ range;
    InterpolatingPolynomial[table, vars] /. s_Real :> Chop[s, 10^-6] // Expand
]

(* If not in database, use NIntegrate *)
NSIntegrate[f_, s:{x_, a_, b_}, opts___] /; ExpandAll@f =!= f := NSIntegrate[ExpandAll@f, s, opts]
NSIntegrate[f_, s:{x_, a_?NumericQ, b_?NumericQ}, opts___] := NIntegrate[f, s, opts]

*Usage:

NSIntegrate[Cos[a*x], {x, 0, 1}, InterpolationVariable -> {a, 1, 5, 1/2}]

InterpolationVariable work just like in a Table. In case or more than 1 interpolation use in the form {{a, 1, 5}, {b, 1, 10, 0.3}}

** Added (* Scaling and offset *) section.

share|improve this question
    
"The problem with my function is the use of TimeConstrained for testing if the integral can be done fully symbolically." - unless you want to go through all the trouble of reimplementing the Risch algorithm, I'm not sure you can do better... –  J. M. Oct 24 '12 at 2:59
    
Yeay, I thought so too, to re-implement all Mathematica is waste of time, it's better just to use interpolation and add some simple rules. I wanted also to make NSIntegrate multidimensional, but this seams a litte bit hard –  Thales Fernandes Oct 24 '12 at 14:47

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.