How to solve the following four integration in MATLAB

Question

chaaru datta 2024-12-23

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2172419-how-to-solve-the-following-four-integration-in-matlab

评论： chaaru datta 2025-1-8

Hello all, I am trying to solve the following expression involving four integrations in MATLAB but note getting it correctly.

This is how I had tried to code:

integrand = @(y, h, z1, z2) (1 / gamma(m_0)) * ...
    gammainc( ((u_1 ./ ...
    (a_1 * zeta_1 * y .* g_1_abs_sq .* h .* P_h1_2 - ...
    u_1 * a_2 * zeta_1 * y .*  g_1_abs_sq .* h .*P_h1_2 - ...
    u_1 * a_1 * zeta_1 *  g_1_abs_sq .* h .* z1 - ...
    u_1 * a_1 * zeta_1 *P_h1_2 * z2 .* g_1_abs_sq .* h - ...
    u_1 * a_2 * zeta_1 * g_1_abs_sq .* h .* z1 - ...
    u_1 * a_2 * zeta_1 * P_h1_2 * z2 .* g_1_abs_sq .* h)) / ...
    (Omega_0 / m_0)),m_0) .* ...
    (1/gamma(m_0))*((m_0/Omega_0)^m_0)*y^(m_0-1)*exp(-m_0*y/Omega_0).* 1/(2*pi*A_1).*lambda_1*exp(-lambda_1*z1).*lambda_2*exp(-lambda_2*z2);
outer_integral = @(z2) arrayfun(@(z2_val) integral(@(y,h,z1) integrand(y,h,z1,z2_val ), 0, y_max(z2_val)), 0:1000);
% Perform the integration
op = integral(outer_integral, 0, 1000);

I am not getting why I am getting such errors:

Not enough input arguments.

Error in Analytical_2>@(z1,z2)(P_h1_2*(a_1-u_1*a_2))./(u_1*a_1*z1+u_1*a_1*P_h1_2*z2+u_1*a_2*z1+u_1*a_2*P_h1_2*z2) (line 96)

(u_1 * a_1 * z1 + u_1 * a_1 * P_h1_2 * z2 + ...

Error in Analytical_2>@(z2_val)integral(@(y,h,z1)integrand(y,h,z1,z2_val),0,y_max(z2_val)) (line 132)

outer_integral = @(z2) arrayfun(@(z2_val) integral(@(y,h,z1) integrand(y,h,z1,z2_val ), 0, y_max(z2_val)), 0:1000);

Error in Analytical_2>@(z2)arrayfun(@(z2_val)integral(@(y,h,z1)integrand(y,h,z1,z2_val),0,y_max(z2_val)),0:1000) (line 132)

outer_integral = @(z2) arrayfun(@(z2_val) integral(@(y,h,z1) integrand(y,h,z1,z2_val ), 0, y_max(z2_val)), 0:1000);

Error in integralCalc/iterateScalarValued (line 314)

fx = FUN(t);

Error in integralCalc/vadapt (line 132)

[q,errbnd] = iterateScalarValued(u,tinterval,pathlen);

Error in integralCalc (line 75)

[q,errbnd] = vadapt(@AtoBInvTransform,interval);

Error in integral (line 87)

Q = integralCalc(fun,a,b,opstruct);

Error in Analytical_2 (line 136)

op = integral(outer_integral, 0, 1000);

Please note that I had intentionally not given other part of code which is used to produce the values of elements that are used in this integrand. Any help in this regard will be highly appreciated.

4 个评论
显示 2更早的评论隐藏 2更早的评论

chaaru datta 2025-1-4

Thank u sir for ur response, but I am not getting why MATLAB cant solve such 4 integrals.

Torsten 2025-1-4

编辑：Torsten 2025-1-4

In MATLAB, there exist functions that can directly be applied to 1d, 2d and 3d integrals (namely integral, integral2 and integral3). If you have a 4d integral, the existing MATLAB routines can be combined (e.g. 3d with 1d, 2d with 2d, 2d with two times 1d and four times 1d) to solve the higher-dimensional integral. And that's exactly what "integralN" does.

If you want to use your personal method to solve your integral (e.g. by taking "integral" four times), you can do so. But as you see from your code, you don't succeed. I recommended "integralN" because it simplifies the set-up of such calculations.

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Walter Roberson 2025-1-4

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2172419-how-to-solve-the-following-four-integration-in-matlab#answer_1556917

在 MATLAB Online 中打开

integral expects a function handle to a function that accepts a single input parameter. It normally passes a vector of values to the function handle, and expects a vector of results of the same size. The size of the vector that it passes in is variable.

You are first passing integral() a function handle to a function that accepts a single parameter, which is good. Your code does

outer_integral = @(z2) arrayfun(@(z2_val) integral(@(y,h,z1) integrand(y,h,z1,z2_val ), 0, y_max(z2_val)), 0:1000);

The arrayfun() part of it is good. But inside the arrayfun() you are passing integral() a function handle of a function that expects up to three input parameters. integral() cannot handle functions that expect three input parameters. integral3 on the other hand is able to accepts function handles that expect three input parameters. Thus, proper code would more likely be

outer_integral = @(z2) arrayfun(@(z2_val) integral3(@(y,h,z1) integrand(y,h,z1,z2_val ), 0, y_max(z2_val)), 0:1000);

17 个评论
显示 15更早的评论隐藏 15更早的评论

chaaru datta 2025-1-6

在 MATLAB Online 中打开

Thank u for ur response sir...but still its not clearing my doubt...In my previous research work , the analytical expression was

where f_z(z) and f_Y(y) are PDF of exponential random variable and Z_lim =

and I had written the MATLAB code for this double integration as follows:

% Define Z_lim a s a function of y
Z_lim = @(y) (gamma_th * ((A + beta_0^2 * y.^2) * (P_u1_lin / jj1) + 1 / jj1)) / (beta_0^2 * y * A_11);
% Define the integrand
integrand = @(y, z) (1-exp(-lambda*(gamma_th*((A + beta_0^2 * y.^2) * (P_u1_lin / jj1) + 1 / jj1)-beta_0^2 * y * z *A_11)/A_11)).* lambda .* exp(-lambda * y) .* lambda .* exp(-lambda * z);
% Define the outer integral with respect to y from 0 to infinity
outer_integral = @(y) arrayfun(@(y_val) integral(@(z) integrand(y_val, z), 0, Z_lim(y_val)), y);
% Perform the double integration
op = integral(outer_integral, 0, inf);

My query is if the above MATLAB code of double integration is working correctly and giving me proper results then why the MATLAB code for four integrals is not working.

Torsten 2025-1-6

编辑：Torsten 2025-1-6

在 MATLAB Online 中打开

@chaaru datta

For example:

f = @(x,y,z,u) x + y + z + u;
value_integral_num = integral(@(x)integral3(@(y,z,u)f(x,y,z,u),0,1,0,1,0,1),0,1,'ArrayValued',1)
value_integral_num = 2.0000
syms x y z u
f =  x + y + z + u;
value_integral_sym = int(int(int(int(f,x,0,1),y,0,1),z,0,1),u,0,1)
value_integral_sym = 2

both compute the four-fold integral of f over [0,1]^4.

But you should try the numerical variant. The symbolic approach is slow and should only be used if there are analytical antiderivatives of your function with respect to the integration variables.

Are you sure that the expression in the denominator has no zeros in the region of integration (which would cause problems for the integration process) ?

    (a_1 * zeta_1 * y .* g_1_abs_sq .* h .* P_h1_2 - ...
    u_1 * a_2 * zeta_1 * y .*  g_1_abs_sq .* h .*P_h1_2 - ...
    u_1 * a_1 * zeta_1 *  g_1_abs_sq .* h .* z1 - ...
    u_1 * a_1 * zeta_1 *P_h1_2 * z2 .* g_1_abs_sq .* h - ...
    u_1 * a_2 * zeta_1 * g_1_abs_sq .* h .* z1 - ...
    u_1 * a_2 * zeta_1 * P_h1_2 * z2 .* g_1_abs_sq .* h)

chaaru datta 2025-1-7

在 MATLAB Online 中打开

I am getting confused now, I will once again share my query here:

I am working on two seperate analytical expression let us call them expression A and expression B.

Expression A is:

where

and

are PDFs of exponential random variable and

.

The MATLAB code (which is working correctly) for evaluating expression (A) is as:

% Define Z_lim a s a function of y
Z_lim = @(y) (gamma_th * ((A + beta_0^2 * y.^2) * (P_u1_lin / jj1) + 1 / jj1)) / (beta_0^2 * y * A_11);
% Define the integrand
integrand = @(y, z) (1-exp(-lambda*(gamma_th*((A + beta_0^2 * y.^2) * (P_u1_lin / jj1) + 1 / jj1)-beta_0^2 * y * z *A_11)/A_11)).* lambda .* exp(-lambda * y) .* lambda .* exp(-lambda * z);  
% Define the outer integral with respect to y from 0 to infinity
outer_integral = @(y) arrayfun(@(y_val) integral(@(z) integrand(y_val, z), 0, Z_lim(y_val)), y);
% Perform the double integration
op = integral(outer_integral, 0, inf);

Next, the expression B is:

---- (B)

where

and

My query is based on the style of expression (A), how to write MATLAB code for expression (B) ?

Torsten 2025-1-7

编辑：Torsten 2025-1-7

在 MATLAB Online 中打开

My query is based on the style of expression (A), how to write MATLAB code for expression (B) ?

Why are you getting confused ?

I gave you the answer previously (some * are replaced by .* and / are replaced by ./ ):

integrand = @(y, h, z1, z2) (1 ./ gamma(m_0)) .* ...
    gammainc( ((u_1 ./ ...
    (a_1 .* zeta_1 .* y .* g_1_abs_sq .* h .* P_h1_2 - ...
    u_1 .* a_2 .* zeta_1 .* y .*  g_1_abs_sq .* h .* P_h1_2 - ...
    u_1 .* a_1 .* zeta_1 .*  g_1_abs_sq .* h .* z1 - ...
    u_1 .* a_1 .* zeta_1 .* P_h1_2 .* z2 .* g_1_abs_sq .* h - ...
    u_1 .* a_2 .* zeta_1 .* g_1_abs_sq .* h .* z1 - ...
    u_1 .* a_2 .* zeta_1 .* P_h1_2 .* z2 .* g_1_abs_sq .* h)) ./ ...
    (Omega_0 ./ m_0)),m_0) .* ...
    (1./gamma(m_0)).*((m_0./Omega_0).^m_0).*y.^(m_0-1).*exp(-m_0.*y./Omega_0).* 1./(2.*pi.*A_1).*lambda_1.*exp(-lambda_1.*z1).*lambda_2.*exp(-lambda_2.*z2);
y_max = @(z1,z2) P_h1_2 .* (a_1 - u_1 .* a_2) ./ ...
                 (u_1 .* a_1 .* z1 + u_1 .* a_1 .* P_h1_2 .* z2 + ...
                  u_1 .* a_2 .* z1 + u_1 .* a_2 .* P_h1_2 .* z2);
op = integral(@(h)integral3(@(z1,z2,y)integrand(y,h,z1,z2),0, inf, 0, inf, 0, y_max),0, Inf, 'ArrayValued',1)
% or
%outer_integral = @(h) arrayfun(@(h_val) integral3(@(z1,z2,y) integrand(y,h_val,z1,z2 ), 0, inf, 0, inf, 0, y_max), h);
%op = integral(outer_integral, 0, inf);

I assumed z1 = Z1 and z2 = Z2 in your mathematical description.

But as written, it will not work because h = 0 gives a division by zero because the denominator of the argument to the gammainc function becomes 0.

chaaru datta 2025-1-8

Ok sir...

请先登录，再进行评论。

How to solve the following four integration in MATLAB

4 个评论
显示 2更早的评论隐藏 2更早的评论

回答（1 个）

17 个评论
显示 15更早的评论隐藏 15更早的评论

另请参阅

类别

标签

Community Treasure Hunt

How to solve the following four integration in MATLAB

4 个评论 显示 2更早的评论隐藏 2更早的评论

回答（1 个）

17 个评论 显示 15更早的评论隐藏 15更早的评论

另请参阅

类别

标签

Community Treasure Hunt

4 个评论
显示 2更早的评论隐藏 2更早的评论

17 个评论
显示 15更早的评论隐藏 15更早的评论