# Speeding up integral3 command is taking too long to execute

2 visualizzazioni (ultimi 30 giorni)
DM il 26 Feb 2015
Modificato: Mike Hosea il 28 Feb 2015
I am using MATLAB to execute a triple integral using integral3 and it is running very slow. I was wondering if there ways to speed it. Below is my program. Note that some of the variables in my code are nonnegative constants: A1, a1, B1, b2, lambda1, lambda2, delta1, delta2, theta..
clear all
syms gamma1
syms gamma2
syms z
syms v
Nt=16; sigmanoise=10^(-7.9); c3=0.129; c1=(1-c3)/2;a2=0;b2=0;
a1=0.0030; b1= 0.0030; A1= 1.5625e-04,A2=0; B1= 7.8125e-05;B2=0;
theta= 3.1623;lambda1= 4.9736e-05;lambda2=0;p1=1;p2=0; alpha1=2; alpha2=4;delta1=2/alpha1; delta2=2/alpha2;beta1=0.025; beta2=0.025;
a= gamma1^-1+gamma2^-1+2*gamma1^(-0.5)*gamma2^(-0.5);
laplacesgi=(exp(+2*pi*j.*z*a)-1)./(2*pi*j*z);
laplacesgi=matlabFunction(laplacesgi);
laplacenoi=exp(-2*pi*j.*z*theta*sigmanoise/Nt);
laplacenoi=matlabFunction(laplacenoi);
interfere= @(gamma1,gamma2,v,z)( (1 -2*c1-c3./(1+2*pi*j*z*theta*v.^(-1))).*(A1.*(v).^(delta1-1).*exp(-a1.*(v).^ (delta1./2))+B1.*(v).^(delta2-1) .*(1-exp(-b1.*(v).^ (delta2./2)))));
gscalar =@(gamma1,gamma2,z)integral(@(v)(interfere(gamma1,gamma2,v,z)),gamma2,inf);
g = @(gamma1,gamma2,z)arrayfun(gscalar,gamma1,gamma2,z);
lp= A1*(gamma1)^(delta1-1)*exp(-a1*(gamma1)^ (delta1/2))+B1*(gamma1)^(delta2-1)*(1-exp(-b1*(gamma1)^ (delta2/2)))+A2*gamma1^(delta1-1)*exp(-a2*gamma1^(delta1/2))+ B2*gamma1^(delta2-1)*(1-exp(-b2*gamma1^ (delta2/2)));%;
dk1=((2*pi*lambda1))/(beta1^2)*(1-exp(-a1*(gamma2)^(delta1/2))*(1+(gamma2)^(delta1/2)*a1))+ pi*lambda1*gamma2^(delta2)*p1^delta2-((2*pi*lambda1)/(beta1^2))*(1-exp(-b1*(gamma2)^(delta2/2))*(1+(gamma2)^(delta2/2)*b1));
dk2=((2*pi*lambda2))/(beta2^2)*(1-exp(-a2*(gamma2)^(delta1/2))*(1+(gamma2)^(delta1/2)*a2))+ pi*lambda2*gamma2^(delta2)*p2^delta2-((2*pi*lambda2)/(beta2^2))*(1-exp(-b2*(gamma2)^(delta2/2))*(1+(gamma2)^(delta2/2)*b2));
dk=dk1+dk2;
lcp= A1*(gamma2)^(delta1-1)*exp(-a1*(gamma2)^ (delta1/2))+B1*(gamma2)^(delta2-1)*(1-exp(-b1*(gamma2)^ (delta2/2)))+A2*gamma2^(delta1-1)*exp(-a2*gamma2^ (delta1/2))+ B2*gamma2^(delta2-1)*(1-exp(-b2*gamma2^(delta2/2)));%;
pdflast=lp*lcp*exp(-dk);
pdflast=matlabFunction(pdflast);
pdflast= @(gamma1,gamma2)arrayfun(pdflast,gamma1,gamma2);
gamma2min=@(gamma1)gamma1;
warning('off','MATLAB:integral:MinStepSize');
T = integral3(@(gamma1,gamma2,z)(laplacenoi(z).*laplacesgi(gamma1,gamma2,z).*pdflast(gamma1,gamma2).*exp(-g(gamma1,gamma2,z))),0,inf,@(gamma2)gamma2,inf,0.05,1000,'abstol',1e-3)
I appreciate any ideas or suggestions.
##### 0 CommentiMostra -2 commenti meno recentiNascondi -2 commenti meno recenti

Accedi per commentare.

### Risposte (1)

Mike Hosea il 26 Feb 2015
These two articles might help.
Try different variable orderings. Make sure your tolerances are sane (do this first). Also, consider whether it is possible to speed up the integrand function.
##### 8 CommentiMostra 6 commenti meno recentiNascondi 6 commenti meno recenti
Mike Hosea il 28 Feb 2015
Yeah, it should work, but trying things out it looked like compiling it with MATLAB Coder speeded it up about 10x. I'm running a slightly modified version now that uses PARFOR for the loop computing quadgk.
Mike Hosea il 28 Feb 2015
Modificato: Mike Hosea il 28 Feb 2015
BTW, just FYI (since you don't have MATLAB Coder), the loop doing quadgk can be a "parfor" in MATLAB Coder. You don't need the Parallel Computing Toolbox for that since it's only multithreading (i.e., not distributed). I do not think it would be a good idea to make it a parfor in MATLAB, but I don't have that much experience with the Parallel Computing Toolbox, so I don't have a good feeling for how much overhead there will be.
But with Coder's multithreading, I get (I've got several spare cores, so YMMV)
>> tic; T = integral3(@integrand_mex,0,inf,@(gamma1)gamma1,inf,0.05,1000,'AbStol',1e-5,'RelTol',1e-3), toc
T =
5.329614020381315e+00 - 2.196183908954074e+02i
Elapsed time is 1434.025262 seconds.

Accedi per commentare.

### Categorie

Scopri di più su Linear Algebra in Help Center e File Exchange

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by