mirror of
https://github.com/snowykami/mbcp.git
synced 2024-11-22 06:07:37 +08:00
🐛 fix ε accuracy
This commit is contained in:
parent
a74811291f
commit
5a0e2f189c
143
main.py
143
main.py
@ -8,3 +8,146 @@ Copyright (C) 2020-2024 LiteyukiStudio. All Rights Reserved
|
||||
@File : main.py
|
||||
@Software: PyCharm
|
||||
"""
|
||||
import logging
|
||||
|
||||
from mbcp.mp_math.line import Line3
|
||||
from mbcp.mp_math.plane import Plane3
|
||||
from mbcp.mp_math.point import Point3
|
||||
|
||||
# def ac8s4e4():
|
||||
# """
|
||||
# 第八章第四节例4
|
||||
# 问题:求与两平面x-4z-3=0和2x-y-5z-1=0的交线平行且过点(-3, 2, 5)的直线方程。
|
||||
# """
|
||||
# correct_ans = Line3(4, 3, 1, 1)
|
||||
#
|
||||
# pl1 = Plane3(1, 0, -4, -3)
|
||||
# pl2 = Plane3(2, -1, -5, -1)
|
||||
# p = Point3(-3, 2, 5)
|
||||
# """解法1"""
|
||||
# # 求直线方向向量s
|
||||
# s = pl1.normal @ pl2.normal
|
||||
# actual_ans = Line3.from_point_and_direction(p, s)
|
||||
#
|
||||
# logging.info(f"正确答案:{correct_ans} 实际答案:{actual_ans}")
|
||||
# assert actual_ans == correct_ans
|
||||
#
|
||||
# """解法2"""
|
||||
# # 过点p且与pl1平行的平面pl3
|
||||
# pl3 = pl1.cal_parallel_plane3(p)
|
||||
# # 过点p且与pl2平行的平面pl4
|
||||
# pl4 = pl2.cal_parallel_plane3(p)
|
||||
# # 求pl3和pl4的交线
|
||||
# actual_ans = pl3.cal_intersection_line3(pl4)
|
||||
# print(pl3, pl4, actual_ans)
|
||||
#
|
||||
# logging.info(f"正确答案:{correct_ans} 实际答案:{actual_ans}")
|
||||
# assert actual_ans == correct_ans
|
||||
#
|
||||
#
|
||||
# ac8s4e4()
|
||||
import logging
|
||||
|
||||
from mbcp.mp_math.mp_math_typing import RealNumber
|
||||
from mbcp.mp_math.utils import Approx
|
||||
|
||||
|
||||
def three_var_func(x: RealNumber, y: RealNumber) -> RealNumber:
|
||||
return x ** 3 * y ** 2 - 3 * x * y ** 3 - x * y + 1
|
||||
|
||||
|
||||
class TestPartialDerivative:
|
||||
# 样例来源:同济大学《高等数学》第八版下册 第九章第二节 例6
|
||||
def test_2v_1o_1v(self):
|
||||
"""测试二元函数关于第一个变量(x)的一阶偏导 df/dx"""
|
||||
|
||||
from mbcp.mp_math.utils import Approx
|
||||
from mbcp.mp_math.equation import get_partial_derivative_func
|
||||
|
||||
partial_derivative_func = get_partial_derivative_func(three_var_func, 0)
|
||||
|
||||
# assert partial_derivative_func(1, 2, 3) == 4.0
|
||||
def df_dx(x, y):
|
||||
"""原函数关于x的偏导"""
|
||||
return 3 * (x ** 2) * (y ** 2) - 3 * (y ** 3) - y
|
||||
|
||||
logging.info(f"Expected: {df_dx(1, 2)}, Actual: {partial_derivative_func(1, 2)}")
|
||||
assert Approx(partial_derivative_func(1, 2)) == df_dx(1, 2)
|
||||
|
||||
def test_2v_1o_2v(self):
|
||||
"""测试二元函数关于第二个变量(y)的一阶偏导 df/dy"""
|
||||
|
||||
from mbcp.mp_math.utils import Approx
|
||||
from mbcp.mp_math.equation import get_partial_derivative_func
|
||||
|
||||
partial_derivative_func = get_partial_derivative_func(three_var_func, 1)
|
||||
|
||||
def df_dy(x, y):
|
||||
"""原函数关于y的偏导"""
|
||||
return 2 * (x ** 3) * y - 9 * x * (y ** 2) - x
|
||||
|
||||
logging.info(f"Expected: {df_dy(1, 2)}, Actual: {partial_derivative_func(1, 2)}")
|
||||
assert Approx(partial_derivative_func(1, 2)) == df_dy(1, 2)
|
||||
|
||||
def test_2v_2o_12v(self):
|
||||
"""高阶偏导d^2f/(dxdy)"""
|
||||
|
||||
from mbcp.mp_math.utils import Approx
|
||||
from mbcp.mp_math.equation import get_partial_derivative_func
|
||||
|
||||
partial_derivative_func = get_partial_derivative_func(three_var_func, (0, 1))
|
||||
|
||||
def df_dxdy(x, y):
|
||||
"""原函数关于y和x的偏导"""
|
||||
return 6 * x ** 2 * y - 9 * y ** 2 - 1
|
||||
|
||||
logging.info(f"Expected: {df_dxdy(1, 2)}, Actual: {partial_derivative_func(1, 2)}")
|
||||
assert Approx(partial_derivative_func(1, 2)) == df_dxdy(1, 2)
|
||||
|
||||
def test_2v_2o_1v2(self):
|
||||
"""二阶偏导d^2f/(dx^2)"""
|
||||
|
||||
from mbcp.mp_math.utils import Approx
|
||||
from mbcp.mp_math.equation import get_partial_derivative_func
|
||||
|
||||
partial_derivative_func = get_partial_derivative_func(three_var_func, (0, 0))
|
||||
|
||||
def df_dydx(x, y):
|
||||
"""原函数关于x和y的偏导"""
|
||||
return 6 * x * y ** 2
|
||||
|
||||
logging.info(f"Expected: {df_dydx(1, 2)}, Actual: {partial_derivative_func(1, 2)}")
|
||||
assert Approx(partial_derivative_func(1, 2)) == df_dydx(1, 2)
|
||||
|
||||
def test_2v_3o_1v3(self):
|
||||
"""高阶偏导d^3f/(dx^3)"""
|
||||
|
||||
from mbcp.mp_math.utils import Approx
|
||||
from mbcp.mp_math.equation import get_partial_derivative_func
|
||||
|
||||
partial_derivative_func = get_partial_derivative_func(three_var_func, (0, 0, 0))
|
||||
|
||||
def d3f_dx3(x, y):
|
||||
"""原函数关于x的三阶偏导"""
|
||||
return 6 * (y ** 2)
|
||||
|
||||
logging.info(f"Expected: {d3f_dx3(1, 2)}, Actual: {partial_derivative_func(1, 2)}")
|
||||
assert Approx(partial_derivative_func(1, 2)) == d3f_dx3(1, 2)
|
||||
|
||||
def test_possible_error(self):
|
||||
from mbcp.mp_math.equation import get_partial_derivative_func
|
||||
def two_vars_func(x: RealNumber, y: RealNumber) -> RealNumber:
|
||||
return x ** 2 * y ** 2
|
||||
|
||||
partial_func = get_partial_derivative_func(two_vars_func, 0)
|
||||
partial_func_2 = get_partial_derivative_func(two_vars_func, (0, 0))
|
||||
assert Approx(partial_func_2(1, 2)) == 8
|
||||
|
||||
|
||||
TestPartialDerivative().test_2v_1o_1v()
|
||||
TestPartialDerivative().test_2v_1o_2v()
|
||||
TestPartialDerivative().test_2v_2o_12v()
|
||||
TestPartialDerivative().test_2v_2o_1v2()
|
||||
TestPartialDerivative().test_2v_3o_1v3()
|
||||
|
||||
TestPartialDerivative().test_possible_error()
|
||||
|
94
mbcp/mp_math/angle.py
Normal file
94
mbcp/mp_math/angle.py
Normal file
@ -0,0 +1,94 @@
|
||||
# -*- coding: utf-8 -*-
|
||||
"""
|
||||
Copyright (C) 2020-2024 LiteyukiStudio. All Rights Reserved
|
||||
|
||||
@Time : 2024/8/26 上午6:29
|
||||
@Author : snowykami
|
||||
@Email : snowykami@outlook.com
|
||||
@File : angle.py
|
||||
@Software: PyCharm
|
||||
"""
|
||||
from typing import overload
|
||||
|
||||
from .const import PI
|
||||
|
||||
|
||||
class AnyAngle:
|
||||
def __init__(self, value: float, is_radian: bool = False):
|
||||
"""
|
||||
任意角度。
|
||||
Args:
|
||||
value: 角度或弧度值
|
||||
is_radian: 是否为弧度,默认为否
|
||||
"""
|
||||
if is_radian:
|
||||
self.radian = value
|
||||
else:
|
||||
self.radian = value * PI / 180
|
||||
|
||||
@property
|
||||
def complementary(self) -> "AnyAngle":
|
||||
"""
|
||||
余角:两角的和为90°。
|
||||
Returns:
|
||||
余角
|
||||
"""
|
||||
return AnyAngle(PI / 2 - self.minimum_positive.radian, is_radian=True)
|
||||
|
||||
@property
|
||||
def supplementary(self) -> "AnyAngle":
|
||||
"""
|
||||
补角:两角的和为180°。
|
||||
Returns:
|
||||
补角
|
||||
"""
|
||||
return AnyAngle(PI - self.minimum_positive.radian, is_radian=True)
|
||||
|
||||
@property
|
||||
def degree(self) -> float:
|
||||
"""
|
||||
角度。
|
||||
Returns:
|
||||
弧度
|
||||
"""
|
||||
return self.radian * 180 / PI
|
||||
|
||||
@property
|
||||
def minimum_positive(self) -> "AnyAngle":
|
||||
"""
|
||||
最小正角。
|
||||
Returns:
|
||||
最小正角度
|
||||
"""
|
||||
return AnyAngle(self.radian % (2 * PI))
|
||||
|
||||
@property
|
||||
def maximum_negative(self) -> "AnyAngle":
|
||||
"""
|
||||
最大负角。
|
||||
Returns:
|
||||
最大负角度
|
||||
"""
|
||||
return AnyAngle(-self.radian % (2 * PI), is_radian=True)
|
||||
|
||||
def __add__(self, other: "AnyAngle") -> "AnyAngle":
|
||||
return AnyAngle(self.radian + other.radian, is_radian=True)
|
||||
|
||||
def __sub__(self, other: "AnyAngle") -> "AnyAngle":
|
||||
return AnyAngle(self.radian - other.radian, is_radian=True)
|
||||
|
||||
def __mul__(self, other: float) -> "AnyAngle":
|
||||
return AnyAngle(self.radian * other, is_radian=True)
|
||||
|
||||
@overload
|
||||
def __truediv__(self, other: float) -> "AnyAngle":
|
||||
...
|
||||
|
||||
@overload
|
||||
def __truediv__(self, other: "AnyAngle") -> float:
|
||||
...
|
||||
|
||||
def __truediv__(self, other):
|
||||
if isinstance(other, AnyAngle):
|
||||
return self.radian / other.radian
|
||||
return AnyAngle(self.radian / other, is_radian=True)
|
@ -46,13 +46,15 @@ class CurveEquation:
|
||||
|
||||
def get_partial_derivative_func(func: MultiVarsFunc, var: int | tuple[int, ...], epsilon: Number = EPSILON) -> MultiVarsFunc:
|
||||
"""
|
||||
求N元函数一阶偏导函数。
|
||||
求N元函数一阶偏导函数。这玩意不太稳定,慎用。
|
||||
Args:
|
||||
func: 函数
|
||||
var: 变量位置,可为整数(一阶偏导)或整数元组(高阶偏导)
|
||||
epsilon: 偏移量
|
||||
Returns:
|
||||
偏导函数
|
||||
Raises:
|
||||
ValueError: 无效变量类型
|
||||
"""
|
||||
if isinstance(var, int):
|
||||
def partial_derivative_func(*args: Var) -> Var:
|
||||
|
@ -5,49 +5,67 @@ Copyright (C) 2020-2024 LiteyukiStudio. All Rights Reserved
|
||||
@Time : 2024/8/6 下午12:57
|
||||
@Author : snowykami
|
||||
@Email : snowykami@outlook.com
|
||||
@File : line.py
|
||||
@File : other.py
|
||||
@Software: PyCharm
|
||||
"""
|
||||
from typing import TYPE_CHECKING
|
||||
import math
|
||||
from typing import TYPE_CHECKING, overload
|
||||
from .vector import Vector3
|
||||
|
||||
if TYPE_CHECKING:
|
||||
from .point import Point3 # type: ignore
|
||||
from .angle import AnyAngle
|
||||
from .plane import Plane3
|
||||
from .point import Point3
|
||||
|
||||
|
||||
class Line3:
|
||||
def __init__(self, a: float, b: float, c: float, d: float):
|
||||
"""
|
||||
三维空间中的直线。
|
||||
:param a:
|
||||
:param b:
|
||||
:param c:
|
||||
:param d:
|
||||
Args:
|
||||
a: 直线方程的系数a
|
||||
b: 直线方程的系数b
|
||||
c: 直线方程的系数c
|
||||
d: 直线方程的常数项d
|
||||
"""
|
||||
self.a = a
|
||||
self.b = b
|
||||
self.c = c
|
||||
self.d = d
|
||||
|
||||
def __str__(self):
|
||||
return f"Line3({self.a}, {self.b}, {self.c}, {self.d})"
|
||||
def cal_angle(self, other: "Line3") -> "AnyAngle":
|
||||
"""
|
||||
计算直线和直线或面之间的夹角。
|
||||
Args:
|
||||
other: 另一条直线或面
|
||||
Returns:
|
||||
夹角弧度
|
||||
Raises:
|
||||
TypeError: 不支持的类型
|
||||
"""
|
||||
if isinstance(other, Line3):
|
||||
return self.direction.cal_angle(other.direction)
|
||||
elif isinstance(other, Plane3):
|
||||
return self.direction.cal_angle(other.normal).complementary # 方向向量和法向量的夹角的余角
|
||||
else:
|
||||
raise TypeError(f"Unsupported type: {type(other)}")
|
||||
|
||||
def get_perpendicular(self, point: "Point3") -> "Line3":
|
||||
@property
|
||||
def direction(self) -> "Vector3":
|
||||
"""
|
||||
获取直线经过指定点p的垂线。
|
||||
:param point: 指定点p,直线外的点
|
||||
:return: 垂直于self且过点p的直线
|
||||
直线的方向向量。
|
||||
Returns:
|
||||
方向向量
|
||||
"""
|
||||
a = -self.b
|
||||
b = self.a
|
||||
c = 0
|
||||
d = -(a * point.x + b * point.y + self.c * point.z)
|
||||
return Line3(a, b, c, d)
|
||||
return Vector3(self.a, self.b, self.c)
|
||||
|
||||
def get_intersection(self, line: "Line3") -> "Point3":
|
||||
def cal_intersection(self, line: "Line3") -> "Point3":
|
||||
"""
|
||||
获取两条直线的交点。
|
||||
:param line:
|
||||
:return:
|
||||
计算两条直线的交点。
|
||||
Args:
|
||||
line: 另一条直线
|
||||
Returns:
|
||||
交点
|
||||
"""
|
||||
|
||||
if self.is_parallel(line):
|
||||
@ -70,21 +88,39 @@ class Line3:
|
||||
|
||||
return Point3(x, y, z)
|
||||
|
||||
def cal_perpendicular(self, point: "Point3") -> "Line3":
|
||||
"""
|
||||
计算直线经过指定点p的垂线。
|
||||
Args:
|
||||
point: 指定点
|
||||
Returns:
|
||||
垂线
|
||||
"""
|
||||
a = -self.b
|
||||
b = self.a
|
||||
c = 0
|
||||
d = -(a * point.x + b * point.y + self.c * point.z)
|
||||
return Line3(a, b, c, d)
|
||||
|
||||
def is_parallel(self, line: "Line3") -> bool:
|
||||
"""
|
||||
判断两条直线是否平行。
|
||||
直线平行的条件是它们的法向量成比例
|
||||
:param line:
|
||||
:return:
|
||||
Args:
|
||||
line: 另一条直线
|
||||
Returns:
|
||||
是否平行
|
||||
"""
|
||||
return self.a * line.b == self.b * line.a and self.c * line.b == self.d * line.a
|
||||
return self.direction.is_parallel(line.direction)
|
||||
|
||||
def is_collinear(self, line: "Line3") -> bool:
|
||||
"""
|
||||
判断两条直线是否共线。
|
||||
直线共线的条件是它们的法向量成比例且常数项也成比例
|
||||
:param line:
|
||||
:return:
|
||||
Args:
|
||||
line: 另一条直线
|
||||
Returns:
|
||||
是否共线
|
||||
"""
|
||||
return self.is_parallel(line) and (self.d * line.b - self.b * line.d) / (self.a * line.b - self.b * line.a) == 0
|
||||
|
||||
@ -92,10 +128,56 @@ class Line3:
|
||||
"""
|
||||
判断两条直线是否共面。
|
||||
两条直线共面的条件是它们的方向向量和法向量的叉乘为零向量
|
||||
:param line:
|
||||
:return:
|
||||
Args:
|
||||
line: 另一条直线
|
||||
Returns:
|
||||
是否共面
|
||||
"""
|
||||
direction1 = (-self.c, 0, self.a)
|
||||
direction2 = (line.c, -line.b, 0)
|
||||
cross_product = direction1[0] * direction2[1] - direction1[1] * direction2[0]
|
||||
return cross_product == 0
|
||||
|
||||
@classmethod
|
||||
def from_point_and_direction(cls, point: "Point3", direction: "Vector3") -> "Line3":
|
||||
"""
|
||||
工厂函数 由点和方向向量构造直线(点向式构造)。
|
||||
Args:
|
||||
point: 点
|
||||
direction: 方向向量
|
||||
Returns:
|
||||
直线
|
||||
"""
|
||||
a, b, c = direction.x, direction.y, direction.z
|
||||
d = -(a * point.x + b * point.y + c * point.z)
|
||||
return cls(a, b, c, d)
|
||||
|
||||
@classmethod
|
||||
def from_two_points(cls, p1: "Point3", p2: "Point3") -> "Line3":
|
||||
"""
|
||||
工厂函数 由两点构造直线。
|
||||
Args:
|
||||
p1: 点1
|
||||
p2: 点2
|
||||
Returns:
|
||||
直线
|
||||
"""
|
||||
direction = p2 - p1
|
||||
return cls.from_point_and_direction(p1, direction)
|
||||
|
||||
def __eq__(self, other) -> bool:
|
||||
"""
|
||||
判断两条直线是否等价。
|
||||
Args:
|
||||
other:
|
||||
|
||||
Returns:
|
||||
|
||||
"""
|
||||
return self.a / other.a == self.b / other.b == self.c / other.c == self.d / other.d
|
||||
|
||||
def __repr__(self):
|
||||
return f"Line3({self.a}, {self.b}, {self.c}, {self.d})"
|
||||
|
||||
def __str__(self):
|
||||
return f"Line3({self.a}, {self.b}, {self.c}, {self.d})"
|
||||
|
@ -1,10 +1,124 @@
|
||||
# -*- coding: utf-8 -*-
|
||||
"""
|
||||
Copyright (C) 2020-2024 LiteyukiStudio. All Rights Reserved
|
||||
|
||||
@Time : 2024/8/6 下午1:45
|
||||
@Author : snowykami
|
||||
@Email : snowykami@outlook.com
|
||||
@File : plane.py
|
||||
@Software: PyCharm
|
||||
平面模块
|
||||
"""
|
||||
import math
|
||||
from typing import TYPE_CHECKING
|
||||
|
||||
import numpy as np
|
||||
|
||||
from .vector import Vector3
|
||||
from .line import Line3
|
||||
from .point import Point3
|
||||
|
||||
if TYPE_CHECKING:
|
||||
from .angle import AnyAngle
|
||||
|
||||
|
||||
class Plane3:
|
||||
def __init__(self, a: float, b: float, c: float, d: float):
|
||||
"""
|
||||
平面方程:ax + by + cz + d = 0
|
||||
Args:
|
||||
a:
|
||||
b:
|
||||
c:
|
||||
d:
|
||||
"""
|
||||
self.a = a
|
||||
self.b = b
|
||||
self.c = c
|
||||
self.d = d
|
||||
|
||||
def cal_angle(self, other: "Line3 | Plane3") -> "AnyAngle":
|
||||
"""
|
||||
计算平面与平面之间的夹角。
|
||||
Args:
|
||||
other: 另一个平面
|
||||
Returns:
|
||||
夹角弧度
|
||||
Raises:
|
||||
TypeError: 不支持的类型
|
||||
"""
|
||||
if isinstance(other, Line3):
|
||||
return self.normal.cal_angle(other.direction).complementary
|
||||
elif isinstance(other, Plane3):
|
||||
return AnyAngle(math.acos(self.normal * other.normal / (self.normal.length * other.normal.length)), is_radian=True)
|
||||
else:
|
||||
raise TypeError(f"Unsupported type: {type(other)}")
|
||||
|
||||
def cal_distance(self, other: "Plane3 | Point3") -> float:
|
||||
"""
|
||||
计算平面与平面或点之间的距离。
|
||||
Args:
|
||||
other: 另一个平面或点
|
||||
Returns:
|
||||
距离
|
||||
Raises:
|
||||
TypeError: 不支持的类型
|
||||
"""
|
||||
if isinstance(other, Plane3):
|
||||
return 0
|
||||
elif isinstance(other, Point3):
|
||||
return abs(self.a * other.x + self.b * other.y + self.c * other.z + self.d) / (self.a ** 2 + self.b ** 2 + self.c ** 2) ** 0.5
|
||||
else:
|
||||
raise TypeError(f"Unsupported type: {type(other)}")
|
||||
|
||||
def cal_intersection_line3(self, other: "Plane3") -> "Line3":
|
||||
"""
|
||||
计算两平面的交线。该方法有问题,待修复。
|
||||
Args:
|
||||
other: 另一个平面
|
||||
Returns:
|
||||
交线
|
||||
"""
|
||||
# 计算两法向量的叉积作为交线的方向向量
|
||||
s = self.normal.cross(other.normal) # 交线的方向向量
|
||||
# 联立两平面方程求交线的一点
|
||||
# 两平面方程联立得到的方程组
|
||||
# | a1x + b1y + c1z = -d1
|
||||
# | a2x + b2y + c2z = -d2
|
||||
# 用numpy解方程组
|
||||
a = np.array([[self.a, self.b, self.c], [other.a, other.b, other.c]])
|
||||
b = np.array([-self.d, -other.d])
|
||||
p = np.linalg.lstsq(a, b, rcond=None)[0]
|
||||
return Line3.from_point_and_direction(Point3(*p), s)
|
||||
|
||||
def cal_parallel_plane3(self, point: "Point3") -> "Plane3":
|
||||
"""
|
||||
计算平行于该平面且过指定点的平面。
|
||||
Args:
|
||||
point: 指定点
|
||||
Returns:
|
||||
平面
|
||||
"""
|
||||
return Plane3.from_point_and_normal(point, self.normal)
|
||||
|
||||
@property
|
||||
def normal(self) -> "Vector3":
|
||||
"""
|
||||
平面的法向量。
|
||||
Returns:
|
||||
法向量
|
||||
"""
|
||||
return Vector3(self.a, self.b, self.c)
|
||||
|
||||
@classmethod
|
||||
def from_point_and_normal(cls, point: "Point3", normal: "Vector3") -> "Plane3":
|
||||
"""
|
||||
工厂函数 由点和法向量构造平面(点法式构造)。
|
||||
Args:
|
||||
point: 平面上的一点
|
||||
normal: 法向量
|
||||
Returns:
|
||||
平面
|
||||
"""
|
||||
a, b, c = normal.x, normal.y, normal.z
|
||||
d = -a * point.x - b * point.y - c * point.z # d = -ax - by - cz
|
||||
return cls(a, b, c, d)
|
||||
|
||||
def __repr__(self):
|
||||
return f"Plane3({self.a}, {self.b}, {self.c}, {self.d})"
|
||||
|
||||
def __str__(self):
|
||||
return f"{self.a}x + {self.b}y + {self.c}z + {self.d} = 0"
|
||||
|
@ -44,6 +44,7 @@ class Point3:
|
||||
:param other:
|
||||
:return:
|
||||
"""
|
||||
from .vector import Vector3
|
||||
return Vector3(self.x - other.x, self.y - other.y, self.z - other.z)
|
||||
|
||||
def __truediv__(self, other: float) -> "Point3":
|
||||
|
@ -1,2 +0,0 @@
|
||||
# py.typed
|
||||
partial
|
@ -1,9 +1,10 @@
|
||||
import math
|
||||
from typing import overload, TYPE_CHECKING
|
||||
|
||||
import numpy as np
|
||||
from .point import Point3
|
||||
|
||||
if TYPE_CHECKING:
|
||||
from .point import Point3 # type: ignore
|
||||
from .angle import AnyAngle
|
||||
|
||||
|
||||
class Vector3:
|
||||
@ -14,75 +15,59 @@ class Vector3:
|
||||
:param y:
|
||||
:param z:
|
||||
"""
|
||||
self._x = x
|
||||
self._y = y
|
||||
self._z = z
|
||||
self._length = (x ** 2 + y ** 2 + z ** 2) ** 0.5
|
||||
self._normalized = self / self._length
|
||||
self.x = x
|
||||
self.y = y
|
||||
self.z = z
|
||||
|
||||
def __str__(self):
|
||||
return f"Vector3({self._x}, {self._y}, {self._z})"
|
||||
def cal_angle(self, other: 'Vector3') -> 'AnyAngle':
|
||||
"""
|
||||
计算两个向量之间的夹角。
|
||||
Args:
|
||||
other: 另一个向量
|
||||
Returns:
|
||||
夹角
|
||||
"""
|
||||
return AnyAngle(math.acos(self * other / (self.length * other.length)), is_radian=True)
|
||||
|
||||
def _unset_properties(self):
|
||||
self._length = None
|
||||
self._normalized = None
|
||||
def is_parallel(self, other: 'Vector3') -> bool:
|
||||
"""
|
||||
判断两个向量是否平行。
|
||||
Args:
|
||||
other: 另一个向量
|
||||
Returns:
|
||||
是否平行
|
||||
"""
|
||||
return self @ other == Vector3(0, 0, 0)
|
||||
|
||||
@property
|
||||
def x(self):
|
||||
return self._x
|
||||
|
||||
@x.setter
|
||||
def x(self, value):
|
||||
self._x = value
|
||||
self._unset_properties()
|
||||
|
||||
@property
|
||||
def y(self):
|
||||
return self._y
|
||||
|
||||
@y.setter
|
||||
def y(self, value):
|
||||
self._y = value
|
||||
self._unset_properties()
|
||||
|
||||
@property
|
||||
def z(self):
|
||||
return self._z
|
||||
|
||||
@z.setter
|
||||
def z(self, value):
|
||||
self._z = value
|
||||
self._unset_properties()
|
||||
def cross(self, other: 'Vector3') -> 'Vector3':
|
||||
"""
|
||||
向量积 叉乘:V1 @ V2 -> V3
|
||||
Args:
|
||||
other:
|
||||
Returns:
|
||||
叉乘结果,为0向量则两向量平行,否则垂直于两向量
|
||||
"""
|
||||
return Vector3(self.y * other.z - self.z * other.y,
|
||||
self.z * other.x - self.x * other.z,
|
||||
self.x * other.y - self.y * other.x)
|
||||
|
||||
@property
|
||||
def length(self) -> float:
|
||||
"""
|
||||
向量的模。
|
||||
:return:
|
||||
Returns:
|
||||
模
|
||||
"""
|
||||
if self._length is None:
|
||||
self._length = (self._x ** 2 + self._y ** 2 + self._z ** 2) ** 0.5
|
||||
return self._length
|
||||
return math.sqrt(self.x ** 2 + self.y ** 2 + self.z ** 2)
|
||||
|
||||
@property
|
||||
def normalized(self) -> 'Vector3':
|
||||
def unit(self) -> 'Vector3':
|
||||
"""
|
||||
返回该向量的单位向量。
|
||||
:return:
|
||||
获取该向量的单位向量。
|
||||
Returns:
|
||||
单位向量
|
||||
"""
|
||||
if self._normalized is None:
|
||||
self._normalized = self / self.length
|
||||
return self._normalized
|
||||
|
||||
def normalize(self):
|
||||
"""
|
||||
自体归一化。
|
||||
"""
|
||||
self._x /= self.length
|
||||
self._y /= self.length
|
||||
self._z /= self.length
|
||||
self._length = 1
|
||||
self._normalized = self
|
||||
return self / self.length
|
||||
|
||||
@overload
|
||||
def __add__(self, other: 'Vector3') -> 'Vector3':
|
||||
@ -96,16 +81,28 @@ class Vector3:
|
||||
"""
|
||||
V + P -> P\n
|
||||
V + V -> V
|
||||
:param other:
|
||||
:return:
|
||||
Args:
|
||||
other:
|
||||
Returns:
|
||||
|
||||
"""
|
||||
if isinstance(other, Vector3):
|
||||
return Vector3(self._x + other.x, self._y + other.y, self._z + other.z)
|
||||
return Vector3(self.x + other.x, self.y + other.y, self.z + other.z)
|
||||
elif isinstance(other, Point3):
|
||||
return Point3(self._x + other.x, self._y + other.y, self._z + other.z)
|
||||
return Point3(self.x + other.x, self.y + other.y, self.z + other.z)
|
||||
else:
|
||||
raise TypeError(f"unsupported operand type(s) for +: 'Vector3' and '{type(other)}'")
|
||||
|
||||
def __eq__(self, other):
|
||||
"""
|
||||
判断两个向量是否相等。
|
||||
Args:
|
||||
other:
|
||||
Returns:
|
||||
是否相等
|
||||
"""
|
||||
return self.x == other.x and self.y == other.y and self.z == other.z
|
||||
|
||||
def __radd__(self, other: 'Point3') -> 'Point3':
|
||||
"""
|
||||
P + V -> P\n
|
||||
@ -113,7 +110,7 @@ class Vector3:
|
||||
:param other:
|
||||
:return:
|
||||
"""
|
||||
return Point3(self._x + other.x, self._y + other.y, self._z + other.z)
|
||||
return Point3(self.x + other.x, self.y + other.y, self.z + other.z)
|
||||
|
||||
@overload
|
||||
def __sub__(self, other: 'Vector3') -> 'Vector3':
|
||||
@ -127,25 +124,28 @@ class Vector3:
|
||||
"""
|
||||
V - P -> P\n
|
||||
V - V -> V
|
||||
:param other:
|
||||
:return:
|
||||
Args:
|
||||
other:
|
||||
Returns:
|
||||
"""
|
||||
if isinstance(other, Vector3):
|
||||
return Vector3(self._x - other.x, self._y - other.y, self._z - other.z)
|
||||
return Vector3(self.x - other.x, self.y - other.y, self.z - other.z)
|
||||
elif isinstance(other, Point3):
|
||||
return Point3(self._x - other.x, self._y - other.y, self._z - other.z)
|
||||
return Point3(self.x - other.x, self.y - other.y, self.z - other.z)
|
||||
else:
|
||||
raise TypeError(f"unsupported operand type(s) for -: 'Vector3' and '{type(other)}'")
|
||||
|
||||
def __rsub__(self, other: Point3):
|
||||
"""
|
||||
P - V -> P
|
||||
:param other:
|
||||
:return:
|
||||
Args:
|
||||
other:
|
||||
Returns:
|
||||
|
||||
"""
|
||||
|
||||
if isinstance(other, Point3):
|
||||
return Point3(other.x - self._x, other.y - self._y, other.z - self._z)
|
||||
return Point3(other.x - self.x, other.y - self.y, other.z - self.z)
|
||||
else:
|
||||
raise TypeError(f"unsupported operand type(s) for -: '{type(other)}' and 'Vector3'")
|
||||
|
||||
@ -159,34 +159,52 @@ class Vector3:
|
||||
|
||||
def __mul__(self, other):
|
||||
"""
|
||||
乘法。包括点乘和数乘。
|
||||
:param other:
|
||||
:return:
|
||||
点乘法。包括点乘和数乘。
|
||||
V * V -> float\n
|
||||
Args:
|
||||
other:
|
||||
Returns:
|
||||
float
|
||||
Raises:
|
||||
TypeError: 不支持的类型
|
||||
"""
|
||||
if isinstance(other, (int, float)):
|
||||
return Vector3(self._x * other, self._y * other, self._z * other)
|
||||
return Vector3(self.x * other, self.y * other, self.z * other)
|
||||
elif isinstance(other, Vector3):
|
||||
return self._x * other.x + self._y * other.y + self._z * other.z
|
||||
return self.x * other.x + self.y * other.y + self.z * other.z
|
||||
else:
|
||||
raise TypeError(f"unsupported operand type(s) for *: 'Vector3' and '{type(other)}'")
|
||||
|
||||
def __rmul__(self, other: float) -> 'Vector3':
|
||||
"""
|
||||
右乘。
|
||||
:param other:
|
||||
:return:
|
||||
Args:
|
||||
other:
|
||||
Returns:
|
||||
乘积
|
||||
"""
|
||||
return Vector3(self._x * other, self._y * other, self._z * other)
|
||||
return Vector3(self.x * other, self.y * other, self.z * other)
|
||||
|
||||
def __matmul__(self, other: 'Vector3') -> 'Vector3':
|
||||
"""
|
||||
叉乘。
|
||||
:param other: 另一个向量
|
||||
:return: 叉乘结果向量
|
||||
向量积 叉乘:V1 @ V2 -> V3
|
||||
Args:
|
||||
other:
|
||||
Returns:
|
||||
叉乘结果,为0向量则两向量平行,否则垂直于两向量
|
||||
"""
|
||||
return Vector3(self._y * other.z - self._z * other.y,
|
||||
self._z * other.x - self._x * other.z,
|
||||
self._x * other.y - self._y * other.x)
|
||||
return Vector3(self.y * other.z - self.z * other.y,
|
||||
self.z * other.x - self.x * other.z,
|
||||
self.x * other.y - self.y * other.x)
|
||||
|
||||
def __truediv__(self, other: float) -> 'Vector3':
|
||||
return Vector3(self._x / other, self._y / other, self._z / other)
|
||||
return Vector3(self.x / other, self.y / other, self.z / other)
|
||||
|
||||
def __neg__(self):
|
||||
return Vector3(-self.x, -self.y, -self.z)
|
||||
|
||||
def __repr__(self):
|
||||
return f"Vector3({self.x}, {self.y}, {self.z})"
|
||||
|
||||
def __str__(self):
|
||||
return f"Vector3({self.x}, {self.y}, {self.z})"
|
||||
|
@ -1,2 +0,0 @@
|
||||
# py.typed
|
||||
partial
|
33
tests/test_line3.py
Normal file
33
tests/test_line3.py
Normal file
@ -0,0 +1,33 @@
|
||||
# -*- coding: utf-8 -*-
|
||||
"""
|
||||
Copyright (C) 2020-2024 LiteyukiStudio. All Rights Reserved
|
||||
|
||||
@Time : 2024/8/26 上午7:54
|
||||
@Author : snowykami
|
||||
@Email : snowykami@outlook.com
|
||||
@File : test_line3.py
|
||||
@Software: PyCharm
|
||||
"""
|
||||
import logging
|
||||
|
||||
from mbcp.mp_math.point import Point3
|
||||
from mbcp.mp_math.vector import Vector3
|
||||
from mbcp.mp_math.line import Line3
|
||||
|
||||
|
||||
class TestLine3:
|
||||
|
||||
def test_point_and_normal_factory(self):
|
||||
"""
|
||||
测试通过点和法向量构造直线
|
||||
"""
|
||||
correct_ans = Line3(1, -2, 3, -8)
|
||||
|
||||
p = Point3(2, -3, 0)
|
||||
n = Vector3(1, -2, 3)
|
||||
|
||||
actual_ans = Line3.from_point_and_direction(p, n)
|
||||
logging.info(f"正确答案:{correct_ans} 实际答案:{actual_ans}")
|
||||
assert actual_ans == correct_ans
|
||||
|
||||
|
29
tests/test_plane3.py
Normal file
29
tests/test_plane3.py
Normal file
@ -0,0 +1,29 @@
|
||||
# -*- coding: utf-8 -*-
|
||||
"""
|
||||
Copyright (C) 2020-2024 LiteyukiStudio. All Rights Reserved
|
||||
|
||||
@Time : 2024/8/26 上午9:07
|
||||
@Author : snowykami
|
||||
@Email : snowykami@outlook.com
|
||||
@File : test_plane3.py
|
||||
@Software: PyCharm
|
||||
"""
|
||||
import logging
|
||||
|
||||
from mbcp.mp_math.line import Line3
|
||||
from mbcp.mp_math.plane import Plane3
|
||||
|
||||
|
||||
class TestPlane3:
|
||||
|
||||
def test_intersection_line3(self):
|
||||
"""
|
||||
测试平面的交线
|
||||
"""
|
||||
correct_ans = Line3(4, 3, 1, 1)
|
||||
|
||||
pl1 = Plane3(1, 0, -4, 23)
|
||||
pl2 = Plane3(2, -1, -5, 33)
|
||||
actual_ans = pl1.cal_intersection_line3(pl2)
|
||||
logging.info(f"正确答案:{correct_ans} 实际答案:{actual_ans}")
|
||||
assert actual_ans == correct_ans
|
@ -1,22 +0,0 @@
|
||||
# -*- coding: utf-8 -*-
|
||||
"""
|
||||
Copyright (C) 2020-2024 LiteyukiStudio. All Rights Reserved
|
||||
|
||||
@Time : 2024/8/6 下午1:05
|
||||
@Author : snowykami
|
||||
@Email : snowykami@outlook.com
|
||||
@File : test_vector.py
|
||||
@Software: PyCharm
|
||||
"""
|
||||
from mbcp.mp_math.vector import Vector3
|
||||
from mbcp.mp_math.point import Point3
|
||||
|
||||
|
||||
def test_v():
|
||||
v1 = Vector3(1, 2, 3)
|
||||
v2 = Vector3(4, 5, 6)
|
||||
v3 = v1 + v2
|
||||
assert v3._x == 5
|
||||
assert v3._y == 7
|
||||
assert v3._z == 9
|
||||
print("test_v 1111passed")
|
52
tests/test_vector3.py
Normal file
52
tests/test_vector3.py
Normal file
@ -0,0 +1,52 @@
|
||||
# -*- coding: utf-8 -*-
|
||||
"""
|
||||
Copyright (C) 2020-2024 LiteyukiStudio. All Rights Reserved
|
||||
|
||||
@Time : 2024/8/26 上午6:58
|
||||
@Author : snowykami
|
||||
@Email : snowykami@outlook.com
|
||||
@File : test_question_1.py
|
||||
@Software: PyCharm
|
||||
"""
|
||||
import logging
|
||||
|
||||
from mbcp.mp_math.vector import Vector3
|
||||
|
||||
|
||||
class TestVector3:
|
||||
|
||||
"""测试问题集"""
|
||||
def test_vector_cross_product(self):
|
||||
"""
|
||||
测试向量叉乘
|
||||
Returns:
|
||||
"""
|
||||
v1 = Vector3(1, 2, 3)
|
||||
v2 = Vector3(3, 4, 5)
|
||||
actual_ans = v1 @ v2
|
||||
correct_ans = Vector3(-2, 4, -2)
|
||||
logging.info(f"正确答案{correct_ans} 实际答案{v1 @ v2}")
|
||||
|
||||
assert correct_ans == actual_ans
|
||||
|
||||
def test_determine_vector_parallel(self):
|
||||
"""
|
||||
测试判断向量是否平行
|
||||
Returns:
|
||||
"""
|
||||
v1 = Vector3(1, 2, 3)
|
||||
v2 = Vector3(3, 6, 9)
|
||||
actual_ans = v1.is_parallel(v2)
|
||||
correct_ans = True
|
||||
logging.info("v1和v2是否平行:%s", v1.is_parallel(v2))
|
||||
assert correct_ans == actual_ans
|
||||
|
||||
v1 = Vector3(1, 2, 3)
|
||||
v2 = Vector3(3, 6, 8)
|
||||
actual_ans = v1.is_parallel(v2)
|
||||
correct_ans = False
|
||||
logging.info("v1和v2是否平行:%s", v1.is_parallel(v2))
|
||||
assert correct_ans == actual_ans
|
||||
|
||||
|
||||
|
43
tests/test_word_problem.py
Normal file
43
tests/test_word_problem.py
Normal file
@ -0,0 +1,43 @@
|
||||
# -*- coding: utf-8 -*-
|
||||
"""
|
||||
应用题测试集
|
||||
"""
|
||||
import logging
|
||||
|
||||
from mbcp.mp_math.line import Line3
|
||||
from mbcp.mp_math.plane import Plane3
|
||||
from mbcp.mp_math.point import Point3
|
||||
|
||||
|
||||
class TestWordProblem:
|
||||
|
||||
def test_c8s4e4(self):
|
||||
"""
|
||||
同济大学《高等数学》第八版 下册 第八章第四节例4
|
||||
问题:求与两平面x-4z-3=0和2x-y-5z-1=0的交线平行且过点(-3, 2, 5)的直线方程。
|
||||
"""
|
||||
correct_ans = Line3(4, 3, 1, 1)
|
||||
|
||||
pl1 = Plane3(1, 0, -4, -3)
|
||||
pl2 = Plane3(2, -1, -5, -1)
|
||||
p = Point3(-3, 2, 5)
|
||||
"""解法1"""
|
||||
# 求直线方向向量s
|
||||
s = pl1.normal @ pl2.normal
|
||||
actual_ans = Line3.from_point_and_direction(p, s)
|
||||
|
||||
logging.info(f"正确答案:{correct_ans} 实际答案:{actual_ans}")
|
||||
assert actual_ans == correct_ans
|
||||
|
||||
"""解法2"""
|
||||
# 过点p且与pl1平行的平面pl3
|
||||
pl3 = pl1.cal_parallel_plane3(p)
|
||||
# 过点p且与pl2平行的平面pl4
|
||||
pl4 = pl2.cal_parallel_plane3(p)
|
||||
# 求pl3和pl4的交线
|
||||
actual_ans = pl3.cal_intersection_line3(pl4)
|
||||
print(pl3, pl4, actual_ans)
|
||||
|
||||
logging.info(f"正确答案:{correct_ans} 实际答案:{actual_ans}")
|
||||
assert actual_ans == correct_ans
|
||||
|
Loading…
Reference in New Issue
Block a user