EDA365电子论坛网

标题: matlab实现曲面积分（surf_integral函数） [打印本页]

作者: ulppknot 时间: 2021-2-1 10:03
标题: matlab实现曲面积分（surf_integral函数）

MATLAB语言并未直接提供曲面积分的现成函数，因此，此处给出计算曲面积分的函数。

目录
函数说明
应用举例
第一类曲面积分
第二类曲面积分
函数实现

函数说明

function I = surf_integral(f,vars,t,a,b)
%surf_integral
%第一类曲面积分
% I = surf_integral(f, z, [x,y], [x_m,x_M], [y_m,y_M])
% I = surf_integral(f, [x,y,z], [u,v], [u_m,u_M], [v_m,v_M])
% Examples:
% 计算int_int(x^2*y+z*y^2)dS, 积分曲面如下：
% x=ucosv, y=usinv, z=v, 0<=u<=a, 0<=v<=2*pi
% MATLAB求解语句
% syms u v; syms a positive;
% x=u*cos(v); y=u*sin(v); z=v; f=x^2*y+z*y^2;
% I = surf_integral(f,[x,y,z],[u,v],[0,a],[0,2*pi])
%
%第二类曲面积分
% I = surf_integral([P,Q,R], z, [x,y], [x_m,x_M], [y_m,y_M])
% I = surf_integral([P,Q,R], [x,y,z], [u,v], [u_m,u_M], [v_m,v_M])
% 注意：I = int_int_S(P*dydz+Q*dxdz+R*dxdy)
% Examples:
% 计算曲面积分int_int(x*y+z)dxdy, 积分曲面如下：
% (x/a)^2+(y/b)^2+(z/c)^2=1的上半部，且积分沿椭球面的上面。
% 引入参数方程：x=a*sin(u)*cos(v)，y=b*sin(u)*sin(v)，z=c*cos(u);
% 且0<=u<=pi/2, 0<=v<=2*pi 。
% MATLAB求解语句
% syms u v; syms a b c positive;
% x=a*sin(u)*cos(v); y=b*sin(u)*sin(v); z=c*cos(u);
% I = surf_integral([0, 0, x*y+z],[x,y,z],[u,v],[0,pi/2],[0,2*pi])
- Y7 t0 I6 m8 g5 |, E

应用举例

第一类曲面积分
1. 计算

, 积分曲面如下：

求解方法

syms u v; syms a positive;
x=u*cos(v); y=u*sin(v); z=v; f=x^2*y+z*y^2;
I = surf_integral(f,[x,y,z],[u,v],[0,a],[0,2*pi])3 |. q9 D& w& d Q5 g$ K6 \

第二类曲面积分

2. 计算曲面积分

, 积分曲面如下：

的上半部，且积分沿椭球面的上面。
注：引入参数方程：

; 且

求解方法

syms u v; syms a b c positive;
x=a*sin(u)*cos(v); y=b*sin(u)*sin(v); z=c*cos(u);
I = surf_integral([0, 0, x*y+z],[x,y,z],[u,v],[0,pi/2],[0,2*pi]). t z8 Y$ W/ v; s% p$ N2 [7 x5 }

函数实现

function I = surf_integral(f,vars,t,a,b)
if length(f)==1
if length(vars)~=1
E = simplify(sum(diff(vars,t(1)).^2));
F = sum(diff(vars,t(1)).*diff(vars,t(2)));
G = simplify(sum(diff(vars,t(2)).^2));
else
E = simplify(1+diff(vars,t(1))^2);
F = diff(vars,t(1))*diff(vars,t(2));
G = simplify(1+diff(vars,t(2))^2);
end
I = int(int(simplify(f*sqrt(E*G-F^2)),t(1),a(1),a(2)),t(2),b(1),b(2));
else
if length(vars)~=1
A = diff(vars(2),t(1))*diff(vars(3),t(2)) - diff(vars(3),t(1))*diff(vars(2),t(2));
B = diff(vars(3),t(1))*diff(vars(1),t(2)) - diff(vars(1),t(1))*diff(vars(3),t(2));
C = diff(vars(1),t(1))*diff(vars(2),t(2)) - diff(vars(2),t(1))*diff(vars(1),t(2));
else
A = - diff(vars,t(1));
B = - diff(vars,t(2));
C = 1;
end
f = f(：); abc = [A, B, C];
I = int(int(simplify(abc*f),t(1),a(1),a(2)),t(2),b(1),b(2));
end
; H: ^: \$ H' E1 U) l

作者: yin123 时间: 2021-2-1 11:07
matlab实现曲面积分（surf_integral函数）

欢迎光临 EDA365电子论坛网 (https://bbs.eda365.com/)