diff --git a/curvilinear_and_surface_integaration.py b/curvilinear_and_surface_integaration.py new file mode 100644 index 0000000..5157d2e --- /dev/null +++ b/curvilinear_and_surface_integaration.py @@ -0,0 +1,290 @@ +#曲线积分与曲面积分:包括两类曲线积分计算、两类曲面积分计算、环量流量计算、散度旋度计算等 + +import unittest +import sympy as sp +from convert_formula import * +import math + +def first_kind_curve_integration(r, t_range, function): + """ + 计算第一类曲线积分 + :param r: 曲线的参数方程,以元组形式给出 (x(t), y(t)) 或 (x(t), y(t), z(t)) + :param t_range: 参数 t 的积分区间 (t1, t2) + :param function: 被积函数,用参数 t 表示 + :return: 第一类曲线积分的值 + """ + t = sp.symbols('t') + if len(r) == 2: + x, y = r + dx_dt = sp.diff(x, t) + dy_dt = sp.diff(y, t) + ds = sp.sqrt(dx_dt**2 + dy_dt**2) + elif len(r) == 3: + x, y, z = r + dx_dt = sp.diff(x, t) + dy_dt = sp.diff(y, t) + dz_dt = sp.diff(z, t) + ds = sp.sqrt(dx_dt**2 + dy_dt**2 + dz_dt**2) + function = sp.sympify(convert_formula(function)) + lower, upper = t_range + result = sp.integrate(function * ds, (t, lower, upper)) + return result + +def second_kind_curve_integration(P, Q, r, t_range): + """ + 计算第二类曲线积分 + :param P: 向量场的 P 分量,用参数 t 表示 + :param Q: 向量场的 Q 分量,用参数 t 表示 + :param r: 曲线的参数方程,以元组形式给出 (x(t), y(t)) + :param t_range: 参数 t 的积分区间 (t1, t2) + :return: 第二类曲线积分的值 + """ + t = sp.symbols('t') + x, y = r + dx_dt = sp.diff(x, t) + dy_dt = sp.diff(y, t) + P = sp.sympify(convert_formula(P)) + Q = sp.sympify(convert_formula(Q)) + lower, upper = t_range + result = sp.integrate(P * dx_dt + Q * dy_dt, (t, lower, upper)) + return result + +def first_kind_surface_integration(r, u_range, v_range, function): + """ + 计算第一类曲面积分 + :param r: 曲面的参数方程,以元组形式给出 (x(u, v), y(u, v), z(u, v)) + :param u_range: 参数 u 的积分区间 (u1, u2) + :param v_range: 参数 v 的积分区间 (v1, v2) + :param function: 被积函数,用参数 u 和 v 表示 + :return: 第一类曲面积分的值 + """ + u, v = sp.symbols('u v') + x, y, z = r + ru = (sp.diff(x, u), sp.diff(y, u), sp.diff(z, u)) + rv = (sp.diff(x, v), sp.diff(y, v), sp.diff(z, v)) + cross_product = (ru[1]*rv[2] - ru[2]*rv[1], ru[2]*rv[0] - ru[0]*rv[2], ru[0]*rv[1] - ru[1]*rv[0]) + E = sp.sqrt(cross_product[0]**2 + cross_product[1]**2 + cross_product[2]**2) + function = sp.sympify(convert_formula(function)) + u1, u2 = u_range + v1, v2 = v_range + integration1 = sp.integrate(function * E, (u, u1, u2)) + result = sp.integrate(integration1, (v, v1, v2)) + return result + +# 第二类曲面积分 +def second_kind_surface_integration(P, Q, R, r, u_range, v_range): + """ + 计算第二类曲面积分 + :param P: 向量场的 P 分量,用参数 u 和 v 表示 + :param Q: 向量场的 Q 分量,用参数 u 和 v 表示 + :param R: 向量场的 R 分量,用参数 u 和 v 表示 + :param r: 曲面的参数方程,以元组形式给出 (x(u, v), y(u, v), z(u, v)) + :param u_range: 参数 u 的积分区间 (u1, u2) + :param v_range: 参数 v 的积分区间 (v1, v2) + :return: 第二类曲面积分的值 + """ + u, v = sp.symbols('u v') + x, y, z = r + ru = (sp.diff(x, u), sp.diff(y, u), sp.diff(z, u)) + rv = (sp.diff(x, v), sp.diff(y, v), sp.diff(z, v)) + cross_product = (ru[1]*rv[2] - ru[2]*rv[1], ru[2]*rv[0] - ru[0]*rv[2], ru[0]*rv[1] - ru[1]*rv[0]) + P = sp.sympify(convert_formula(P)) + Q = sp.sympify(convert_formula(Q)) + R = sp.sympify(convert_formula(R)) + u1, u2 = u_range + v1, v2 = v_range + integration1 = sp.integrate(P * cross_product[0] + Q * cross_product[1] + R * cross_product[2], (u, u1, u2)) + result = sp.integrate(integration1, (v, v1, v2)) + return result + +# 环量计算 +def circulation(P, Q, r, t_range): + """ + 计算环量 + :param P: 向量场的 P 分量,用参数 t 表示 + :param Q: 向量场的 Q 分量,用参数 t 表示 + :param r: 曲线的参数方程,以元组形式给出 (x(t), y(t)) + :param t_range: 参数 t 的积分区间 (t1, t2) + :return: 环量的值 + """ + return second_kind_curve_integration(P, Q, r, t_range) + +# 流量计算 +def flux(P, Q, R, r, u_range, v_range): + """ + 计算流量 + :param P: 向量场的 P 分量,用参数 u 和 v 表示 + :param Q: 向量场的 Q 分量,用参数 u 和 v 表示 + :param R: 向量场的 R 分量,用参数 u 和 v 表示 + :param r: 曲面的参数方程,以元组形式给出 (x(u, v), y(u, v), z(u, v)) + :param u_range: 参数 u 的积分区间 (u1, u2) + :param v_range: 参数 v 的积分区间 (v1, v2) + :return: 流量的值 + """ + return second_kind_surface_integration(P, Q, R, r, u_range, v_range) + +# 散度计算 +def divergence(P, Q, R): + """ + 计算散度 + :param P: 向量场的 P 分量,用 x, y, z 表示 + :param Q: 向量场的 Q 分量,用 x, y, z 表示 + :param R: 向量场的 R 分量,用 x, y, z 表示 + :return: 散度的值 + """ + x, y, z = sp.symbols('x y z') + P = sp.sympify(convert_formula(P)) + Q = sp.sympify(convert_formula(Q)) + R = sp.sympify(convert_formula(R)) + div = sp.diff(P, x) + sp.diff(Q, y) + sp.diff(R, z) + return div + +# 旋度计算 +def curl(P, Q, R): + """ + 计算旋度 + :param P: 向量场的 P 分量,用 x, y, z 表示 + :param Q: 向量场的 Q 分量,用 x, y, z 表示 + :param R: 向量场的 R 分量,用 x, y, z 表示 + :return: 旋度的值,以元组形式给出 (curl_x, curl_y, curl_z) + """ + x, y, z = sp.symbols('x y z') + P = sp.sympify(convert_formula(P)) + Q = sp.sympify(convert_formula(Q)) + R = sp.sympify(convert_formula(R)) + curl_x = sp.diff(R, y) - sp.diff(Q, z) + curl_y = sp.diff(P, z) - sp.diff(R, x) + curl_z = sp.diff(Q, x) - sp.diff(P, y) + return curl_x, curl_y, curl_z + + +class TestIntegralFunctions(unittest.TestCase): + + def test_first_kind_curve_integration(self): + """测试第一类曲线积分""" + t = sp.symbols('t') + r = (sp.cos(t), sp.sin(t)) + t_range = (0, 2*math.pi) + function = "1" + result = first_kind_curve_integration(r, t_range, function) + # 圆的周长为 2π + self.assertEqual(result, 2*math.pi) + + def test_second_kind_curve_integration(self): + """测试第二类曲线积分""" + t = sp.symbols('t') + r = (sp.cos(t), sp.sin(t)) + t_range = (0, 2*math.pi) + P = "-y" + Q = "x" + result = second_kind_curve_integration(P, Q, r, t_range) + # 环量为 2π + self.assertEqual(result, 2*math.pi) + + def test_first_kind_surface_integration(self): + """测试第一类曲面积分""" + u, v = sp.symbols('u v') + r = (sp.cos(u)*sp.sin(v), sp.sin(u)*sp.sin(v), sp.cos(v)) + u_range = (0, 2*math.pi) + v_range = (0, math.pi) + function = "1" + result = first_kind_surface_integration(r, u_range, v_range, function) + # 单位球的表面积为 4π + self.assertEqual(result, 4*math.pi) + + def test_second_kind_surface_integration(self): + """测试第二类曲面积分""" + u, v = sp.symbols('u v') + r = (sp.cos(u)*sp.sin(v), sp.sin(u)*sp.sin(v), sp.cos(v)) + u_range = (0, 2*math.pi) + v_range = (0, math.pi) + P = "x" + Q = "y" + R = "z" + result = second_kind_surface_integration(P, Q, R, r, u_range, v_range) + # 流量为 4π + self.assertEqual(result, 4*math.pi) + + def test_divergence(self): + """测试散度计算""" + P = "x" + Q = "y" + R = "z" + result = divergence(P, Q, R) + self.assertEqual(result, 3) + + def test_curl(self): + """测试旋度计算""" + P = "x" + Q = "y" + R = "z" + result = curl(P, Q, R) + self.assertEqual(result, (0, 0, 0)) + + +def UI(): + while True: + command = input('请选择将要进行的操作:1.第一类曲线积分 2.第二类曲线积分 3.第一类曲面积分 4.第二类曲面积分 5.环量计算 6.流量计算 7.散度计算 8.旋度计算 9.退出(所有积分区间请以元组或列表的形式给出)') + if command == '1': + r = input('请输入曲线的参数方程,以元组形式给出 (x(t), y(t)) 或 (x(t), y(t), z(t)):') + t_range = input('请输入参数 t 的积分区间:') + function = input('请输入被积函数,用参数 t 表示:') + print('第一类曲线积分值:', first_kind_curve_integration(r, t_range, function)) + if command == '2': + P = input('请输入向量场的 P 分量,用参数 t 表示:') + Q = input('请输入向量场的 Q 分量,用参数 t 表示:') + r = input('请输入曲线的参数方程,以元组形式给出 (x(t), y(t)):') + t_range = input('请输入参数 t 的积分区间:') + print('第二类曲线积分值:', second_kind_curve_integration(P, Q, r, t_range)) + if command == '3': + r = input('请输入曲面的参数方程,以元组形式给出 (x(u, v), y(u, v), z(u, v)):') + u_range = input('请输入参数 u 的积分区间:') + v_range = input('请输入参数 v 的积分区间:') + function = input('请输入被积函数,用参数 u 和 v 表示:') + print('第一类曲面积分:', first_kind_surface_integration(r, u_range, v_range, function)) + if command == '4': + P = input('请输入向量场的 P 分量,用参数 u 和 v 表示:') + Q = input('请输入向量场的 Q 分量,用参数 u 和 v 表示:') + R = input('请输入向量场的 R 分量,用参数 u 和 v 表示:') + r = input('请输入曲面的参数方程,以元组形式给出 (x(u, v), y(u, v), z(u, v)):') + u_range = input('请输入参数 u 的积分区间:') + v_range = input('请输入参数 v 的积分区间:') + print('第二类曲面积分:', second_kind_surface_integration(P, Q, R, r, u_range, v_range)) + if command == '5': + P = input('请输入向量场的 P 分量,用参数 t 表示:') + Q = input('请输入向量场的 Q 分量,用参数 t 表示:') + r = input('请输入曲线的参数方程,以元组形式给出 (x(t), y(t)):') + t_range = input('请输入参数 t 的积分区间:') + print('环量:', circulation(P, Q, r, t_range)) + if command == '6': + P = input('请输入向量场的 P 分量,用参数 u 和 v 表示:') + Q = input('请输入向量场的 Q 分量,用参数 u 和 v 表示:') + R = input('请输入向量场的 R 分量,用参数 u 和 v 表示:') + r = input('请输入曲面的参数方程,以元组形式给出 (x(u, v), y(u, v), z(u, v)):') + u_range = input('请输入参数 u 的积分区间:') + v_range = input('请输入参数 v 的积分区间:') + print('流量:', flux(P, Q, R, r, u_range, v_range)) + if command == '7': + P = input('请输入向量场的 P 分量,用 x, y, z 表示:') + Q = input('请输入向量场的 Q 分量,用 x, y, z 表示:') + R = input('请输入向量场的 R 分量,用 x, y, z 表示:') + print('散度:', divergence(P, Q, R)) + if command == '8': + P = input('请输入向量场的 P 分量,用 x, y, z 表示:') + Q = input('请输入向量场的 Q 分量,用 x, y, z 表示:') + R = input('请输入向量场的 R 分量,用 x, y, z 表示:') + print('旋度:', curl(P, Q, R)) + if command == '9': + break + + +if __name__ == "__main__": + UI() + unittest.main(exit=False) + + + + + + +