在 Matlab 中将积分表达式求解(限制)为数值
solve (limit) integral expression into numerical values in Matlab
我的表达式(经过最终评估)原来是一个非常复杂的积分表达式。它涉及函数的反函数,在函数内。
我正在使用符号函数来解决这个问题。不幸的是,我无法得到我需要的表达方式的答案。有人可以指导我吗?
以下是我的代码;
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)
以下是我得到的最终答案。 (请注意,我需要它的浮点格式)
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)
符号计算通常会努力保持尽可能准确的答案。对于定积分,这意味着引擎需要找到被积函数的反导数,然后在积分边界对其求值。如果找不到反导数,引擎可能会抛出 error/warning,或者在 Matlab Symbolic Toolbox 和其他一些引擎的情况下,抛出 "it just returns int(f)
"。但是,您可以使用 double
强制 int
到 return 定积分的数值近似。将最后一行更改为
int_z = double(int(int_x,z,lo_z_lim,hi_z_lim));
将产生一个数值近似值。但是,由于被积函数的复杂性,需要相当长的时间。正如我反复喜欢讲的那样,如果您正在寻找一个数字,您可能希望在最后使用数字。考虑使用 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);
运行修改后的代码输出
int_z =
0.6957
Elapsed time is 841.910756 seconds.
int_z_num =
0.6957
Elapsed time is 0.219571 seconds.
我的表达式(经过最终评估)原来是一个非常复杂的积分表达式。它涉及函数的反函数,在函数内。
我正在使用符号函数来解决这个问题。不幸的是,我无法得到我需要的表达方式的答案。有人可以指导我吗?
以下是我的代码;
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)
以下是我得到的最终答案。 (请注意,我需要它的浮点格式)
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)
符号计算通常会努力保持尽可能准确的答案。对于定积分,这意味着引擎需要找到被积函数的反导数,然后在积分边界对其求值。如果找不到反导数,引擎可能会抛出 error/warning,或者在 Matlab Symbolic Toolbox 和其他一些引擎的情况下,抛出 "it just returns int(f)
"。但是,您可以使用 double
强制 int
到 return 定积分的数值近似。将最后一行更改为
int_z = double(int(int_x,z,lo_z_lim,hi_z_lim));
将产生一个数值近似值。但是,由于被积函数的复杂性,需要相当长的时间。正如我反复喜欢讲的那样,如果您正在寻找一个数字,您可能希望在最后使用数字。考虑使用 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);
运行修改后的代码输出
int_z =
0.6957
Elapsed time is 841.910756 seconds.
int_z_num =
0.6957
Elapsed time is 0.219571 seconds.