计算量不小 区间不宜过大
zlabel('f5')
title('f5')
......................
...........................
matlab求二重积分
g = @(x,y)
1.*((x-2.5).^2+(y-1).^2<=1)+0.*((x-2.5).^2+(y-1).^2>1);
f = @(x,y)
1.*((x-2).^2+(y-1).^2<=1)+0.*((x-2).^2+(y-1).^2>1);
z = @(x,y) f(x,y)+g(x,y)
Q = dblquad(z,0,4,-1,3)
...........
% 以下代码在7.1版以上均可运行。
r=1;
a=2; % 输入a的值
Phi=3; % 输入Phi的值
f1 = @(Theta_2,Phi_2) sin(Theta_2).*sin(Theta_2).*cos(Phi_2);
f2 = @(Theta_2,Phi_2) sqrt(r^2+a^2-2*r*a.*sin(Theta_2).*cos(Phi-Phi_2));
f3 = @(Theta_2,Phi_2) f1(Theta_2,Phi_2)./f2(Theta_2,Phi_2);
f5 = dblquad(f3,0,pi,0,pi)
......................
...........................
..................
dblquad-
Numerically evaluate double
integral over rectangle
Syntax
q = dblquad(fun,xmin,xmax,ymin,ymax)
q = dblquad(fun,xmin,xmax,ymin,ymax,tol)
q = dblquad(fun,xmin,xmax,ymin,ymax,tol,method)
Description
q = dblquad(fun,xmin,xmax,ymin,ymax)
callsthe quad
function to evaluatethe double integral fun(x,y) over the
rectangle xmin<= x <= xmax,
ymin <=y <= ymax.
fun isa function handle. See Function
Handles in the MATLAB Programming documentationfor more
information. fun(x,y) must accept a vector x anda
scalar y and return a vector of values of
theintegrand.
Parameterizing
Functions, in the MATLAB Mathematicsdocumentation, explains how
to provide additional parameters to thefunction fun, if
necessary.
q = dblquad(fun,xmin,xmax,ymin,ymax,tol)
usesa tolerance tol instead of the default, which is
1.0e-6.
q =
dblquad(fun,xmin,xmax,ymin,ymax,tol,method) usesthe
quadrature function specified as method, insteadof the
default quad. Valid values for method are
@quadl orthe function handle of a user-defined quadrature
method that has thesame calling sequence as quad and
quadl.
Examples
Pass function handle @integrnd to dblquad:
Q = dblquad(@integrnd,pi,2*pi,0,pi);
where the function integrnd.m is:
function z = integrnd(x, y)
z = y*sin(x)+x*cos(y);
Pass anonymous function handle F to
dblquad:
F = @(x,y)y*sin(x)+x*cos(y);
Q = dblquad(F,pi,2*pi,0,pi);
The integrnd function integrates
y*sin(x)+x*cos(y) overthe square pi <=
x <= 2*pi, 0 <= y
<= pi. Note that the integrand can be evaluated
with a vector x anda scalar y.
Nonsquare regions can be handled by setting the integrand tozero
outside of the region. For example, the volume of a
hemisphereis:
dblquad(@(x,y)sqrt(max(1-(x.^2+y.^2),0)), -1, 1, -1, 1)
or
dblquad(@(x,y)sqrt(1-(x.^2+y.^2)).*(x.^2+y.^2<=1), -1, 1, -1, 1)
See Also
...............
quad2d-
Numerically evaluate double
integral over planar region
Syntax
q = quad2d(fun,a,b,c,d)
[q,errbnd] = quad2d(...)
q = quad2d(fun,a,b,c,d,param1,val1,param2,val2,...)
Description
q = quad2d(fun,a,b,c,d) approximates
theintegral of fun(x,y) over the planar region
and
. fun is afunction handle,
c and d mayeach be a scalar or a function
handle.
All input functions must be vectorized. The function
Z=fun(X,Y) mustaccept 2-D matrices X and
Y ofthe same size and return a matrix Z of
correspondingvalues. The functions ymin=c(X) and
ymax=d(X) mustaccept matrices and return matrices of the
same size with correspondingvalues.
[q,errbnd] = quad2d(...). errbnd
isan approximate upper bound on the absolute error, |Q -
I|,where I denotes the exact value of the
integral.
q =
quad2d(fun,a,b,c,d,param1,val1,param2,val2,...)
performsthe integration as above with specified values of optional
parameters:
AbsTol
absolute error tolerance
RelTol
relative error tolerance
quad2d
attempts to satisfy ERRBND<=
max(AbsTol,RelTol*|Q|). This is absolute error controlwhen
|Q| is sufficiently small and relative errorcontrol when
|Q| is larger. A default tolerancevalue is used when a
tolerance is not specified. The default valueof AbsTol is
1e-5. The default value of RelTol is
100*eps(class(Q)).This is also the minimum value of
RelTol. Smaller RelTol valuesare automatically
increased to the default value.
MaxFunEvals
Maximum allowed number of evaluations of
fun reached.
The MaxFunEvals parameter limits the numberof
vectorized calls to fun. The default is 2000.
FailurePlot
Generate a plot if MaxFunEvals is
reached.
Setting FailurePlot to true generatesa
graphical representation of the regions needing further
refinementwhen MaxFunEvals is reached. No plot is
generatedif the integration succeeds before reaching
MaxFunEvals.These (generally) 4-sided regions are mapped
to rectangles internally.Clusters of small regions indicate the
areas of difficulty. The defaultis false.
Singular
Problem may have boundary singularities
With Singular set to true, quad2d
will employ transformations to weaken boundary singularities for
better performance. The defaultis true. Setting
Singular to false willturn these transformations
off, which may provide a performance benefiton some smooth
problems.
Examples
Example 1
Integrate over
,
. The true value of the integralis
.
Q = quad2d(@(x,y) y.*sin(x)+x.*cos(y),pi,2*pi,0,pi)
Example 2
Integrate over the triangle
and
. The integrand is infinite
at(0,0). The true value of the integral is .
fun = @(x,y) 1./(sqrt(x + y) .* (1 + x + y).^2 )
In Cartesian coordinates:
ymax = @(x) 1 - x;
Q = quad2d(fun,0,1,0,ymax)
In polar coordinates:
polarfun = @(theta,r) fun(r.*cos(theta),r.*sin(theta)).*r;
rmax = @(theta) 1./(sin(theta) + cos(theta));
Q = quad2d(polarfun,0,pi/2,0,rmax)
Limitations
quad2d
begins by mappingthe region of integration to a rectangle.
Consequently, it may havetrouble integrating over a region that
does not have four sides orhas a side that cannot be mapped
smoothly to a straight line. Ifthe integration is unsuccessful,
some helpful tactics are leaving Singular set to its
default value of true, changing betweenCartesian and polar
coordinates, or breaking the region of integrationinto pieces and
adding the results of integration over the pieces.
For example:
fun = @(x,y)abs(x.^2 + y.^2 - 0.25);
c = @(x)-sqrt(1 - x.^2);
d = @(x)sqrt(1 - x.^2);
quad2d(fun,-1,1,c,d,'AbsTol',1e-8,...
'FailurePlot',true,'Singular',false)
Warning: Reached the maximum number of function ...
evaluations (2000). The result fails the ...
global error test.
The failure plot shows twoareas of difficulty, near the points
(-1,0) and (1,0) andnear the circle
:
Changing the value of Singular to true
willcope with the geometric singularities at (-1,0) and
(1,0).The larger shaded areas may need refinement but are
probably not areasof difficulty.
Q = quad2d(fun,-1,1,c,d,'AbsTol',1e-8, ...
'FailurePlot',true,'Singular',true)
Warning: Reached the maximum number of function ...
evaluations (2000). The result passes the ...
global error test.
From here youcan take advantage of symmetry:
Q = 4*quad2d(fun,0,1,0,d,'Abstol',1e-8,...
'Singular',true, 'FailurePlot',true)
However, the code is still working very hard near the
singularity.It may not be able to provide higher accuracy:
Q = 4*quad2d(fun,0,1,0,d,'Abstol',1e-10,...
'Singular',true,'FailurePlot',true)
Warning: Reached the maximum number of function ...
evaluations (2000). The result passes the ...
global error test.
At higher accuracy, a change in coordinates may work better.
polarfun = @(theta,r) fun(r.*cos(theta),r.*sin(theta)).*r;
Q = 4*quad2d(polarfun,0,pi/2,0,1,'AbsTol',1e-10)
It is best to put the singularity on the boundary by
splittingthe region of integration into two parts:
Q1 = 4*quad2d(polarfun,0,pi/2,0,0.5,'AbsTol',5e-11);
Q2 = 4*quad2d(polarfun,0,pi/2,0.5,1,'AbsTol',5e-11);
Q = Q1 + Q2
References
[1] L.F. Shampine, "Matlab Program for Quadraturein
2D."Applied Mathematics and Computation. Vol.
202,Issue 1, 2008, pp. 266–274.
See Also