MathGroup Archive 1996

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Re: Solving a numerical integration

  • To: mathgroup at smc.vnet.net
  • Subject: [mg5502] Re: [mg5476] Re: [mg5440] Solving a numerical integration
  • From: Allan Hayes <hay at haystack.demon.co.uk>
  • Date: Sat, 14 Dec 1996 19:26:06 -0500
  • Sender: owner-wri-mathgroup at wolfram.com

fransm at win.tue.nl (Frans Martens)
in [mg5476] Re: [mg5440] Solving a numerical integration

Gave the attached solution.
This can be speeded up (and maybe simplified ?) by using NDSolve  
instead of NIntegrate - see below.
We could to automate the search for an approximate zero  
economically by gradually extending the solution and looking for a  
change of sign.

I compare the times of analous parts of the two solutions.

New solution using NDSolve:

f[x_]:= Sqrt[1+x +x^2] -3

soln = y/.First[NDSolve[{y'[x]==f[x], y[0]==0},y,{x,0,6}]];//Timing

	{0.133333 Second, Null}

Plot[soln[x],{x,0,6}];//Timing

	{0.2 Second, Null}

FindRoot[soln[x]==0,{x, 4.8}, Jacobian -> f[x]]//Timing
	{0.0166667 Second, {x -> 4.66231}}

Frans' Solution using NIntegrate:

Clear[int]

int[x_]:= NIntegrate[f[t],{t,0,x},
    AccuracyGoal -> 6];//Timing

	{0.0166667 Second, Null}

Plot[int[x],{x,0,6}];//Timing

	{3.8 Second, Null}

FindRoot[int[x]==0,{x,4.8}, Jacobian -> f[x]]//Timing

	{0.783333 Second, {x -> 4.66231}}


Allan Hayes
hay at haystack.demon.cc.uk
http://www.haystack.demon.co.uk

******************
Frans' posting

Pere Llosas wrote:


> I would like to solve an equation of this kind
>
> NIntegrate[f[n],{n,0,x}]==0    (f[x_]=Sqrt[1+x+x^2...)
> where x is the searched value, and f cannot be integrated

> analytically.
>
> NSolve[NIntegrate[f[n],{n,0,x}]==0,x], tries to evaluate
> NIntegrate[f[n],{n,0,x}] before assignin a numerical value to x and
> returns an error.
>
> How could this calculation be done without having to write a

> program that searches the root?

The equation NIntegrate[f[n],{n,0,x}]==0 has root x = 0 and the
NSolve function tries to compute the inverse of the function x |->  
NIntegrate[f[n],{n,0,x}] .

The function FindRoot is more suitable. Here is an example with a
second root in the neigbourhood of x = 4.8 .

In[25]:=
  Clear[int,f]
  int[x_]:=NIntegrate[f[t],{t,0,x},
    AccuracyGoal -> 6];
  f[x_]:=Sqrt[1+x+x^2]-3;
In[28]:=
  FindRoot[int[x]==0,{x,4.8},
    Jacobian -> f[x]]
Out[28]=
  {x -> 4.66231}

>>>>>     The whole above without messages   <<<<<<<<<

There are two precautions:

1)The option AccuracyGoal in NIntegrate is set to 6 because the

  integral int[x] equals zero for the root x.
2)FindRoot uses the method of Newton-Raphson and setting the option 

  Jacobian prevents the symbolic computation of int'[x] . Note that 

  int'[x] equals f[x].


You must have a global idea of the roots of the original equation.

Frans Martens
Eindhoven
The Netherlands





  • Prev by Date: Re: InterpolatingFunction::dmwarn: Warning
  • Next by Date: Re: RSA and PowerMod-Function
  • Previous by thread: Re: Solving a numerical integration
  • Next by thread: Re: How can I Flatten from the inside out, not the outside