CODEWITHSUNDEEP
functionmatlabsymbolic-integration
SJa
SJa

Reputation: 515

solve (limit) integral expression into numerical values in Matlab

My expression (after final evaluation) turned out to be a very complicated integral expression. It involves inverse of a function, within function.

I am using symbolic function to solve this. Unfortunately i could not get the answer in terms of expression that I need. Can someone guide me with that ?

Following is my code;

    dn=2;
    pfa = .1;
    dref = 1;

    syms x z
    lo_x_lim = erfinv(pfa-1)*sqrt(2);
    hi_x_lim = inf;

    lo_z_lim = 0
    hi_z_lim = dn

    %% Expression

    Q(x,z) = .5 +.5*erf(x/sqrt(2))
    Qinv_ = erfinv(Q - pfa/2) %Note that this need to be changed
    f_step_2 = (Qinv_ - x)/(1-z/dref)
    f_step_3 = Q(f_step_2,0)

    g(x,z) = 2*z/dn^2
    h(x,z) = exp(-x^2/2)

    %combining three equations to form one
    fin = f_step_3*g*h

    % Doing double integral wrt x and z in two steps
    %integrate wrt x
    int_x = int(fin,x,-0.2,hi_x_lim)

    %integrate wrt z

    int_z = int(int_x,z,lo_z_lim,hi_z_lim)

Following is the final answer that I am getting. (Note that I need it in terms floating-point format)

int_z =

int(int((z*exp(-x^2/2)*(erf((2^(1/2)*(x - erfinv(erf((2^(1/2)*x)/2)/2 + 9/20)))/(2*(z - 1)))/2 + 1/2))/2, x, -.2, Inf), z, 0, 2)

Upvotes: 0

Views: 67

Answers (1)

TroyHaskin
TroyHaskin

Reputation: 8401

Symbolic computations typically strive to maintain as exact an answer as possible. For a definite integral, this means the engine needs to find the anti-derivative of the integrand and then evaluate it at the integration bounds. If it cannot find the anti-derivative, the engine may throw an error/warning or, in the case of the Matlab Symbolic Toolbox and some other engines, "it just returns int(f)". However, you can force int to return a numerical approximation to the definite integral using double. Changing the final line to 

int_z = double(int(int_x,z,lo_z_lim,hi_z_lim));

will produce a numerical approximation. However, due to the complexities of the integrand, it takes a fairly long time. As I repeatedly like to soapbox about, if you're in search of a number, you probably want to use numbers at the end. Consider the direct numerical approximation using integral2:

fin_num   = matlabFunction(fin,'Vars',[x,z]);
int_z_num = integral2(fin_num,-0.2,hi_x_lim,lo_z_lim,hi_z_lim);

Running the modified code outputs

int_z =
    0.6957
Elapsed time is 841.910756 seconds.
int_z_num =
    0.6957
Elapsed time is 0.219571 seconds.

Upvotes: 2

Related Questions