Message Boards

GROUPS:

Staff Picks Mathematics Calculus High-Performance Computing Wolfram Language Numerical Computation

Eric Brown

A parallel adaptive sampling strategy (integration/plotting)

Eric Brown, RCG Global Services

Posted 2 months ago

510 Views

1 Reply

6 Total Likes

Follow this post

A parallel adaptive sampling strategy (integration/plotting)
by Eric Brown

Hello Wolfram Community,

I am interested in the parallel evaluation integrals where the integrand is very expensive to compute. Such a feature can be found in some implementations of “cubature” where the whole vector of points needed for the current recursive subdivision are made available so that they can be handled in whatever batch seems appropriate.

I have not found such a capability in Wolfram Language/Mathematica, but I have attached a notebook where I demonstrate it is possible “trick” NIntegrate into yielding its points by supplying a fake/dummy function that is rapid to evaluate, such as f[x]:=0. Those are Reap/Sowed, and then sent to ParallelMap. (Note that previous recursion must be successively provided to this trick so that it knows how to subdivide the current level!)

After memoizing those values, they can be supplied to a second invocation of NIntegrate, but this time using the memoized variable as the integrand.

(Here, I use “LocalAdaptive” because I believe that it optimally-reuses lower levels of recursion.)

Disclaimer: this is a rough cut of an idea, and not well tested. But I am finding it useful in higher dimensions, and the technique is applicable to other functions which can be expensive an adaptive, like ContourPlot3D as well.

In[]:=

LaunchKernels[]

Out[]=

KernelObject

Name: Local kernel

KernelID: 1

,KernelObject

Name: Local kernel

KernelID: 2

,KernelObject

Name: Local kernel

KernelID: 3

,KernelObject

Name: Local kernel

KernelID: 4



In[]:=

f[x_?NumericQ]:=

Cos[Log[x]/x];

In[]:=

Plot[f[x],{x,0.1,0.2}]

Out[]=

NIntegrate Reference Values

Establish reference values for integration using standard way:

MaxRecursion  0

In[]:=

NIntegrate[f[x],{x,0.1,0.2},Method->{"LocalAdaptive","SymbolicProcessing"->0},MinRecursion->0,MaxRecursion0]

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.2}. NIntegrate obtained -0.0955463 and 0.437466 for the integral and error estimates.

Out[]=

-0.0955463

MaxRecursion  1

In[]:=

NIntegrate[f[x],{x,0.1,0.2},Method->{"LocalAdaptive","SymbolicProcessing"->0},MinRecursion->1,MaxRecursion1]

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.15}. NIntegrate obtained -0.0260387 and 0.0690225 for the integral and error estimates.

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.2}. NIntegrate obtained -0.0744665 and 0.000155576 for the integral and error estimates.

Out[]=

-0.100505

MaxRecursion  2

In[]:=

NIntegrate[f[x],{x,0.1,0.2},Method->{"LocalAdaptive","SymbolicProcessing"->0},MinRecursion->2,MaxRecursion2]

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.125}. NIntegrate obtained 0.00566124 and 0.00275659 for the integral and error estimates.

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.15}. NIntegrate obtained -0.0316706 and 0.000201692 for the integral and error estimates.

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.175}. NIntegrate obtained 0.0281529 and 0.0000213722 for the integral and error estimates.

General

:Further output of NIntegrate::ncvb will be suppressed during this calculation.

Out[]=

-0.100476

MaxRecursion  3

In[]:=

NIntegrate[f[x],{x,0.1,0.2},Method->{"LocalAdaptive","SymbolicProcessing"->0},MinRecursion->3,MaxRecursion3]

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.1125}. NIntegrate obtained -0.0429909 and 0.0000300396 for the integral and error estimates.

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.125}. NIntegrate obtained 0.0486525 and 8.68573×

-6

for the integral and error estimates.

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.1375}. NIntegrate obtained -0.0755929 and 7.38119×

-7

for the integral and error estimates.

General

:Further output of NIntegrate::ncvb will be suppressed during this calculation.

Out[]=

-0.100476

ParallelNIntegrate Prototype

In[]:=

ParallelNIntegrate[f_,{xmin_,xmax_},maxr_]:=Module[{memoize,mtemp,res,fake,mem,xnew,known,minr,integrand},(*keepalistofknownxvalues*)known={};(*ensurenumericevaluationinsideNIntegrate*)mem[x_?NumericQ]:=memoize[x];fake[x_?NumericQ]:=0;(*incrementMaxRecursionfrom0tomaxr*)Do[(*FirstInvocation:"trick"NIntegratesothatiteitherlooksupvaluesfrompreviouslowerrecursivelevels,orquicklyyieldsafakevalue*)(*Print["Recursion r = "<>ToString[r]];Print["---------"];*)xnew=Quiet[Reap[NIntegrate[If[MemberQ[known,x],mem[x],fake[x]],{x,xmin,xmax},Method->{"LocalAdaptive","SymbolicProcessing"->0},MinRecursion->r,MaxRecursionr,EvaluationMonitor:>((*Print[{"1st: ",x,If[MemberQ[known,x],mem[x],fake[x]],MemberQ[known,x]}];*)Sow[x])]][[2,1]]];(*returnthepreviousSecondInvocationresultifnonewdata*)If[Length[xnew]==0,Return[res]];(*forvaluesthatwehavenotseen,evaluateinparallel*)xnew=Complement[xnew,known];mtemp=ParallelMap[f[#]&,xnew,Method->"FinestGrained"];(*memoizethosevaluesthatwerecomputedinparallel*)Do[memoize[xnew[[i]]]=mtemp[[i]],{i,Length[mtemp]}];(*appendthosex-valuestotheknownlist*)known=Join[known,xnew];(*SecondInvocation:NIntegratewithmemoizedvalues*)res=Quiet[NIntegrate[mem[x],{x,xmin,xmax},Method->{"LocalAdaptive","SymbolicProcessing"->0},MinRecursion->r,MaxRecursionr(*,EvaluationMonitor:>Print[{"2nd: ",x,mem[x]}]*)]];(*incrementuptomaxr*),{r,0,maxr}];res];

In[]:=

ParallelNIntegrate[f,{0.1,0.2},0]

Out[]=

-0.0955463

In[]:=

ParallelNIntegrate[f,{0.1,0.2},1]

Out[]=

-0.100505

In[]:=

ParallelNIntegrate[f,{0.1,0.2},2]

Out[]=

-0.100476

In[]:=

ParallelNIntegrate[f,{0.1,0.2},3]

Out[]=

-0.100476

Sleep Test

Insert a Pause statement to simulate a long running integration

In[]:=

f[x_?NumericQ]:=Pause[5];

Cos[Log[x]/x];

In[]:=

AbsoluteTiming[NIntegrate[f[x],{x,0.1,0.2},Method->{"LocalAdaptive","SymbolicProcessing"->0},MinRecursion->1,MaxRecursion1]]

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.15}. NIntegrate obtained -0.0260387 and 0.0690225 for the integral and error estimates.

NIntegrate

:NIntegrate failed to converge to prescribed accuracy after 0 recursive bisections in x near {x} = {0.2}. NIntegrate obtained -0.0744665 and 0.000155576 for the integral and error estimates.

Out[]=

{95.0225,-0.100505}

In[]:=

AbsoluteTiming[ParallelNIntegrate[f,{0.1,0.2},1]]

Out[]=

{35.0655,-0.100505}

Best regards,
Eric

POSTED BY: Eric Brown

Answer

1 Reply

Sort By:

Moderation Team

Moderation Team, WOLFRAM

Posted 2 months ago

-- you have earned *Featured Contributor Badge* Your exceptional post has been selected for our editorial column *Staff Picks* http://wolfr.am/StaffPicks and Your Profile is now distinguished by a *Featured Contributor Badge* and is displayed on the Featured Contributor Board. Thank you!

POSTED BY: Moderation Team

Answer

Reply to this discussion

Reply Preview

Attachments

Remove Add a file to this post

Follow this discussion

or Discard

Group Abstract

Feedback

Enable JavaScript to interact with content and submit forms on Wolfram websites. Learn how »