第32章科学计算

学习目标

完成本章学习后，你将能够：

掌握NumPy数组操作：数组创建、索引、切片、形状变换
理解广播机制：数组运算的自动扩展规则
使用NumPy进行数值计算：线性代数、统计运算、随机数生成
应用SciPy科学计算：优化、插值、积分、信号处理
实现矩阵运算：矩阵分解、特征值、线性方程组求解
进行数值优化：最小化、曲线拟合、约束优化
处理科学数据：数据插值、数值积分、傅里叶变换
构建科学计算应用：物理模拟、工程计算实例

32.1 NumPy基础

32.1.1 数组创建

python

import numpy as np
from typing import List, Tuple, Optional, Union
from dataclasses import dataclass


class ArrayFactory:
    @staticmethod
    def from_list(data: List, dtype: Optional[np.dtype] = None) -> np.ndarray:
        return np.array(data, dtype=dtype)

    @staticmethod
    def zeros(shape: Tuple[int, ...], dtype: np.dtype = np.float64) -> np.ndarray:
        return np.zeros(shape, dtype=dtype)

    @staticmethod
    def ones(shape: Tuple[int, ...], dtype: np.dtype = np.float64) -> np.ndarray:
        return np.ones(shape, dtype=dtype)

    @staticmethod
    def empty(shape: Tuple[int, ...], dtype: np.dtype = np.float64) -> np.ndarray:
        return np.empty(shape, dtype=dtype)

    @staticmethod
    def full(shape: Tuple[int, ...], fill_value: Union[int, float]) -> np.ndarray:
        return np.full(shape, fill_value)

    @staticmethod
    def arange(start: int, stop: Optional[int] = None, step: int = 1) -> np.ndarray:
        return np.arange(start, stop, step)

    @staticmethod
    def linspace(start: float, stop: float, num: int = 50) -> np.ndarray:
        return np.linspace(start, stop, num)

    @staticmethod
    def logspace(start: float, stop: float, num: int = 50, base: float = 10.0) -> np.ndarray:
        return np.logspace(start, stop, num, base=base)

    @staticmethod
    def identity(n: int, dtype: np.dtype = np.float64) -> np.ndarray:
        return np.identity(n, dtype=dtype)

    @staticmethod
    def eye(n: int, m: Optional[int] = None, k: int = 0) -> np.ndarray:
        return np.eye(n, m, k)

    @staticmethod
    def diag(values: Union[List, np.ndarray], k: int = 0) -> np.ndarray:
        return np.diag(values, k)

    @staticmethod
    def random_uniform(low: float = 0.0, high: float = 1.0, size: Tuple[int, ...] = None) -> np.ndarray:
        return np.random.uniform(low, high, size)

    @staticmethod
    def random_normal(loc: float = 0.0, scale: float = 1.0, size: Tuple[int, ...] = None) -> np.ndarray:
        return np.random.normal(loc, scale, size)

    @staticmethod
    def random_integers(low: int, high: int, size: Tuple[int, ...] = None) -> np.ndarray:
        return np.random.randint(low, high, size)


@dataclass
class ArrayInfo:
    shape: Tuple[int, ...]
    ndim: int
    size: int
    dtype: np.dtype
    itemsize: int
    nbytes: int

    @classmethod
    def from_array(cls, arr: np.ndarray) -> "ArrayInfo":
        return cls(
            shape=arr.shape,
            ndim=arr.ndim,
            size=arr.size,
            dtype=arr.dtype,
            itemsize=arr.itemsize,
            nbytes=arr.nbytes
        )

    def __repr__(self) -> str:
        return (
            f"ArrayInfo(shape={self.shape}, ndim={self.ndim}, "
            f"size={self.size}, dtype={self.dtype})"
        )


def demonstrate_array_creation():
    arr1d = ArrayFactory.arange(0, 10, 2)
    arr2d = ArrayFactory.zeros((3, 4))
    arr3d = ArrayFactory.random_normal(size=(2, 3, 4))

    print(f"1D数组: {arr1d}")
    print(f"2D数组:\n{arr2d}")
    print(f"数组信息: {ArrayInfo.from_array(arr3d)}")

32.1.2 数组索引与切片

python

class ArrayIndexer:
    def __init__(self, array: np.ndarray):
        self.array = array

    def get_element(self, *indices) -> Union[np.ndarray, np.number]:
        return self.array[indices]

    def set_element(self, value: Union[int, float], *indices) -> None:
        self.array[indices] = value

    def slice_1d(self, start: Optional[int] = None, stop: Optional[int] = None, step: Optional[int] = None) -> np.ndarray:
        return self.array[start:stop:step]

    def slice_2d(
        self,
        row_start: Optional[int] = None,
        row_stop: Optional[int] = None,
        col_start: Optional[int] = None,
        col_stop: Optional[int] = None
    ) -> np.ndarray:
        return self.array[row_start:row_stop, col_start:col_stop]

    def get_rows(self, indices: List[int]) -> np.ndarray:
        return self.array[indices, :]

    def get_columns(self, indices: List[int]) -> np.ndarray:
        return self.array[:, indices]

    def boolean_index(self, condition: np.ndarray) -> np.ndarray:
        return self.array[condition]

    def fancy_index(self, indices: Union[List, np.ndarray]) -> np.ndarray:
        return self.array[indices]

    def where(self, condition: np.ndarray, x: Optional[np.ndarray] = None, y: Optional[np.ndarray] = None) -> np.ndarray:
        if x is not None and y is not None:
            return np.where(condition, x, y)
        return np.where(condition)


class ArraySlicer:
    @staticmethod
    def get_diagonal(array: np.ndarray, offset: int = 0) -> np.ndarray:
        return np.diagonal(array, offset)

    @staticmethod
    def get_upper_triangle(array: np.ndarray, k: int = 0) -> np.ndarray:
        return np.triu(array, k)

    @staticmethod
    def get_lower_triangle(array: np.ndarray, k: int = 0) -> np.ndarray:
        return np.tril(array, k)

    @staticmethod
    def take(array: np.ndarray, indices: Union[List, np.ndarray], axis: Optional[int] = None) -> np.ndarray:
        return np.take(array, indices, axis)

    @staticmethod
    def put(array: np.ndarray, indices: Union[List, np.ndarray], values: Union[List, np.ndarray]) -> None:
        np.put(array, indices, values)

    @staticmethod
    def extract(condition: np.ndarray, array: np.ndarray) -> np.ndarray:
        return np.extract(condition, array)

32.1.3 数组形状操作

python

class ArrayReshaper:
    def __init__(self, array: np.ndarray):
        self.array = array

    def reshape(self, new_shape: Tuple[int, ...]) -> np.ndarray:
        return self.array.reshape(new_shape)

    def flatten(self, order: str = "C") -> np.ndarray:
        return self.array.flatten(order)

    def ravel(self, order: str = "C") -> np.ndarray:
        return np.ravel(self.array, order)

    def transpose(self, axes: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return np.transpose(self.array, axes)

    def swapaxes(self, axis1: int, axis2: int) -> np.ndarray:
        return np.swapaxes(self.array, axis1, axis2)

    def squeeze(self, axis: Optional[Union[int, Tuple[int, ...]]] = None) -> np.ndarray:
        return np.squeeze(self.array, axis)

    def expand_dims(self, axis: int) -> np.ndarray:
        return np.expand_dims(self.array, axis)

    def tile(self, reps: Tuple[int, ...]) -> np.ndarray:
        return np.tile(self.array, reps)

    def repeat(self, repeats: Union[int, List[int]], axis: Optional[int] = None) -> np.ndarray:
        return np.repeat(self.array, repeats, axis)

    def split(self, indices_or_sections: Union[int, List[int]], axis: int = 0) -> List[np.ndarray]:
        return np.split(self.array, indices_or_sections, axis)

    def hsplit(self, indices_or_sections: Union[int, List[int]]) -> List[np.ndarray]:
        return np.hsplit(self.array, indices_or_sections)

    def vsplit(self, indices_or_sections: Union[int, List[int]]) -> List[np.ndarray]:
        return np.vsplit(self.array, indices_or_sections)


class ArrayCombiner:
    @staticmethod
    def concatenate(arrays: List[np.ndarray], axis: int = 0) -> np.ndarray:
        return np.concatenate(arrays, axis)

    @staticmethod
    def stack(arrays: List[np.ndarray], axis: int = 0) -> np.ndarray:
        return np.stack(arrays, axis)

    @staticmethod
    def hstack(arrays: List[np.ndarray]) -> np.ndarray:
        return np.hstack(arrays)

    @staticmethod
    def vstack(arrays: List[np.ndarray]) -> np.ndarray:
        return np.vstack(arrays)

    @staticmethod
    def dstack(arrays: List[np.ndarray]) -> np.ndarray:
        return np.dstack(arrays)

    @staticmethod
    def column_stack(arrays: List[np.ndarray]) -> np.ndarray:
        return np.column_stack(arrays)

    @staticmethod
    def row_stack(arrays: List[np.ndarray]) -> np.ndarray:
        return np.row_stack(arrays)

32.2 数组运算

32.2.1 基本运算

python

class ArrayOperations:
    def __init__(self, array: np.ndarray):
        self.array = array

    def add(self, other: Union[np.ndarray, int, float]) -> np.ndarray:
        return self.array + other

    def subtract(self, other: Union[np.ndarray, int, float]) -> np.ndarray:
        return self.array - other

    def multiply(self, other: Union[np.ndarray, int, float]) -> np.ndarray:
        return self.array * other

    def divide(self, other: Union[np.ndarray, int, float]) -> np.ndarray:
        return self.array / other

    def floor_divide(self, other: Union[np.ndarray, int, float]) -> np.ndarray:
        return self.array // other

    def power(self, other: Union[np.ndarray, int, float]) -> np.ndarray:
        return self.array ** other

    def mod(self, other: Union[np.ndarray, int, float]) -> np.ndarray:
        return self.array % other

    def negative(self) -> np.ndarray:
        return -self.array

    def reciprocal(self) -> np.ndarray:
        return 1.0 / self.array

    def abs(self) -> np.ndarray:
        return np.abs(self.array)

    def sign(self) -> np.ndarray:
        return np.sign(self.array)

    def sqrt(self) -> np.ndarray:
        return np.sqrt(self.array)

    def square(self) -> np.ndarray:
        return np.square(self.array)

    def exp(self) -> np.ndarray:
        return np.exp(self.array)

    def log(self) -> np.ndarray:
        return np.log(self.array)

    def log10(self) -> np.ndarray:
        return np.log10(self.array)

    def log2(self) -> np.ndarray:
        return np.log2(self.array)


class TrigonometricOperations:
    def __init__(self, array: np.ndarray):
        self.array = array

    def sin(self) -> np.ndarray:
        return np.sin(self.array)

    def cos(self) -> np.ndarray:
        return np.cos(self.array)

    def tan(self) -> np.ndarray:
        return np.tan(self.array)

    def arcsin(self) -> np.ndarray:
        return np.arcsin(self.array)

    def arccos(self) -> np.ndarray:
        return np.arccos(self.array)

    def arctan(self) -> np.ndarray:
        return np.arctan(self.array)

    def sinh(self) -> np.ndarray:
        return np.sinh(self.array)

    def cosh(self) -> np.ndarray:
        return np.cosh(self.array)

    def tanh(self) -> np.ndarray:
        return np.tanh(self.array)

    def deg2rad(self) -> np.ndarray:
        return np.deg2rad(self.array)

    def rad2deg(self) -> np.ndarray:
        return np.rad2deg(self.array)


class RoundingOperations:
    @staticmethod
    def around(array: np.ndarray, decimals: int = 0) -> np.ndarray:
        return np.around(array, decimals)

    @staticmethod
    def floor(array: np.ndarray) -> np.ndarray:
        return np.floor(array)

    @staticmethod
    def ceil(array: np.ndarray) -> np.ndarray:
        return np.ceil(array)

    @staticmethod
    def trunc(array: np.ndarray) -> np.ndarray:
        return np.trunc(array)

    @staticmethod
    def rint(array: np.ndarray) -> np.ndarray:
        return np.rint(array)

    @staticmethod
    def fix(array: np.ndarray) -> np.ndarray:
        return np.fix(array)

32.2.2 广播机制

python

class BroadcastingDemo:
    @staticmethod
    def explain_broadcasting(shape1: Tuple[int, ...], shape2: Tuple[int, ...]) -> Tuple[Tuple[int, ...], str]:
        shape1 = np.array(shape1)
        shape2 = np.array(shape2)

        max_len = max(len(shape1), len(shape2))
        shape1 = np.pad(shape1, (max_len - len(shape1), 0), constant_values=1)
        shape2 = np.pad(shape2, (max_len - len(shape2), 0), constant_values=1)

        result_shape = []
        explanation = []

        for i, (s1, s2) in enumerate(zip(shape1, shape2)):
            if s1 == s2:
                result_shape.append(s1)
                explanation.append(f"维度{i}: {s1} == {s2}, 保持 {s1}")
            elif s1 == 1:
                result_shape.append(s2)
                explanation.append(f"维度{i}: {s1} -> {s2} (广播)")
            elif s2 == 1:
                result_shape.append(s1)
                explanation.append(f"维度{i}: {s2} -> {s1} (广播)")
            else:
                raise ValueError(f"无法广播: 维度{i} 形状 {s1} 和 {s2} 不兼容")

        return tuple(result_shape), "\n".join(explanation)

    @staticmethod
    def broadcast_to(array: np.ndarray, shape: Tuple[int, ...]) -> np.ndarray:
        return np.broadcast_to(array, shape)

    @staticmethod
    def broadcast_arrays(*arrays: np.ndarray) -> List[np.ndarray]:
        return np.broadcast_arrays(*arrays)


class UniversalFunctions:
    @staticmethod
    def add(arrays: List[np.ndarray]) -> np.ndarray:
        return np.add.reduce(arrays)

    @staticmethod
    def multiply(arrays: List[np.ndarray]) -> np.ndarray:
        return np.multiply.reduce(arrays)

    @staticmethod
    def accumulate_add(array: np.ndarray, axis: int = 0) -> np.ndarray:
        return np.add.accumulate(array, axis)

    @staticmethod
    def accumulate_multiply(array: np.ndarray, axis: int = 0) -> np.ndarray:
        return np.multiply.accumulate(array, axis)

    @staticmethod
    def outer_add(a: np.ndarray, b: np.ndarray) -> np.ndarray:
        return np.add.outer(a, b)

    @staticmethod
    def outer_multiply(a: np.ndarray, b: np.ndarray) -> np.ndarray:
        return np.multiply.outer(a, b)

    @staticmethod
    def apply_along_axis(func: callable, array: np.ndarray, axis: int = 0) -> np.ndarray:
        return np.apply_along_axis(func, axis, array)

    @staticmethod
    def apply_over_axes(func: callable, array: np.ndarray, axes: Tuple[int, ...]) -> np.ndarray:
        return np.apply_over_axes(func, array, axes)

    @staticmethod
    def vectorize(func: callable) -> callable:
        return np.vectorize(func)

32.2.3 统计运算

python

class StatisticalOperations:
    def __init__(self, array: np.ndarray):
        self.array = array

    def sum(self, axis: Optional[int] = None, keepdims: bool = False) -> Union[np.ndarray, np.number]:
        return np.sum(self.array, axis=axis, keepdims=keepdims)

    def mean(self, axis: Optional[int] = None, keepdims: bool = False) -> Union[np.ndarray, np.number]:
        return np.mean(self.array, axis=axis, keepdims=keepdims)

    def median(self, axis: Optional[int] = None, keepdims: bool = False) -> Union[np.ndarray, np.number]:
        return np.median(self.array, axis=axis, keepdims=keepdims)

    def std(self, axis: Optional[int] = None, ddof: int = 0, keepdims: bool = False) -> Union[np.ndarray, np.number]:
        return np.std(self.array, axis=axis, ddof=ddof, keepdims=keepdims)

    def var(self, axis: Optional[int] = None, ddof: int = 0, keepdims: bool = False) -> Union[np.ndarray, np.number]:
        return np.var(self.array, axis=axis, ddof=ddof, keepdims=keepdims)

    def min(self, axis: Optional[int] = None, keepdims: bool = False) -> Union[np.ndarray, np.number]:
        return np.min(self.array, axis=axis, keepdims=keepdims)

    def max(self, axis: Optional[int] = None, keepdims: bool = False) -> Union[np.ndarray, np.number]:
        return np.max(self.array, axis=axis, keepdims=keepdims)

    def argmin(self, axis: Optional[int] = None) -> Union[np.ndarray, int]:
        return np.argmin(self.array, axis=axis)

    def argmax(self, axis: Optional[int] = None) -> Union[np.ndarray, int]:
        return np.argmax(self.array, axis=axis)

    def percentile(self, q: Union[float, List[float]], axis: Optional[int] = None) -> Union[np.ndarray, float]:
        return np.percentile(self.array, q, axis=axis)

    def quantile(self, q: Union[float, List[float]], axis: Optional[int] = None) -> Union[np.ndarray, float]:
        return np.quantile(self.array, q, axis=axis)

    def cumsum(self, axis: Optional[int] = None) -> np.ndarray:
        return np.cumsum(self.array, axis=axis)

    def cumprod(self, axis: Optional[int] = None) -> np.ndarray:
        return np.cumprod(self.array, axis=axis)

    def histogram(self, bins: int = 10, range: Optional[Tuple[float, float]] = None) -> Tuple[np.ndarray, np.ndarray]:
        return np.histogram(self.array, bins=bins, range=range)

    def unique(self, return_counts: bool = False) -> Union[np.ndarray, Tuple[np.ndarray, np.ndarray]]:
        if return_counts:
            return np.unique(self.array, return_counts=True)
        return np.unique(self.array)


class CorrelationOperations:
    @staticmethod
    def corrcoef(x: np.ndarray, y: Optional[np.ndarray] = None) -> np.ndarray:
        return np.corrcoef(x, y)

    @staticmethod
    def cov(x: np.ndarray, y: Optional[np.ndarray] = None) -> np.ndarray:
        return np.cov(x, y)

    @staticmethod
    def correlate(a: np.ndarray, v: np.ndarray, mode: str = "valid") -> np.ndarray:
        return np.correlate(a, v, mode)

    @staticmethod
    def convolve(a: np.ndarray, v: np.ndarray, mode: str = "full") -> np.ndarray:
        return np.convolve(a, v, mode)

32.3 线性代数

32.3.1 矩阵运算

python

class LinearAlgebra:
    @staticmethod
    def dot(a: np.ndarray, b: np.ndarray) -> np.ndarray:
        return np.dot(a, b)

    @staticmethod
    def matmul(a: np.ndarray, b: np.ndarray) -> np.ndarray:
        return np.matmul(a, b)

    @staticmethod
    def inner(a: np.ndarray, b: np.ndarray) -> np.ndarray:
        return np.inner(a, b)

    @staticmethod
    def outer(a: np.ndarray, b: np.ndarray) -> np.ndarray:
        return np.outer(a, b)

    @staticmethod
    def tensordot(a: np.ndarray, b: np.ndarray, axes: Union[int, Tuple[List[int], List[int]]] = 2) -> np.ndarray:
        return np.tensordot(a, b, axes)

    @staticmethod
    def einsum(subscripts: str, *operands: np.ndarray) -> np.ndarray:
        return np.einsum(subscripts, *operands)

    @staticmethod
    def trace(a: np.ndarray, offset: int = 0) -> Union[int, float]:
        return np.trace(a, offset)

    @staticmethod
    def det(a: np.ndarray) -> float:
        return np.linalg.det(a)

    @staticmethod
    def inv(a: np.ndarray) -> np.ndarray:
        return np.linalg.inv(a)

    @staticmethod
    def pinv(a: np.ndarray, rcond: float = 1e-15) -> np.ndarray:
        return np.linalg.pinv(a, rcond)

    @staticmethod
    def matrix_rank(a: np.ndarray, tol: Optional[float] = None) -> int:
        return np.linalg.matrix_rank(a, tol)

    @staticmethod
    def solve(a: np.ndarray, b: np.ndarray) -> np.ndarray:
        return np.linalg.solve(a, b)

    @staticmethod
    def lstsq(a: np.ndarray, b: np.ndarray, rcond: Optional[float] = None) -> Tuple[np.ndarray, np.ndarray, int, np.ndarray]:
        return np.linalg.lstsq(a, b, rcond)


class MatrixDecomposition:
    @staticmethod
    def eig(a: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        return np.linalg.eig(a)

    @staticmethod
    def eigvals(a: np.ndarray) -> np.ndarray:
        return np.linalg.eigvals(a)

    @staticmethod
    def eigh(a: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        return np.linalg.eigh(a)

    @staticmethod
    def eigvalsh(a: np.ndarray) -> np.ndarray:
        return np.linalg.eigvalsh(a)

    @staticmethod
    def svd(a: np.ndarray, full_matrices: bool = True) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
        return np.linalg.svd(a, full_matrices)

    @staticmethod
    def qr(a: np.ndarray, mode: str = "reduced") -> Tuple[np.ndarray, np.ndarray]:
        return np.linalg.qr(a, mode)

    @staticmethod
    def cholesky(a: np.ndarray) -> np.ndarray:
        return np.linalg.cholesky(a)

    @staticmethod
    def norm(a: np.ndarray, ord: Optional[Union[int, str]] = None, axis: Optional[int] = None) -> Union[float, np.ndarray]:
        return np.linalg.norm(a, ord, axis)

    @staticmethod
    def cond(a: np.ndarray, p: Optional[Union[int, str]] = None) -> float:
        return np.linalg.cond(a, p)


class LinearSystemSolver:
    @staticmethod
    def solve_triangular(a: np.ndarray, b: np.ndarray, lower: bool = False) -> np.ndarray:
        from scipy.linalg import solve_triangular
        return solve_triangular(a, b, lower=lower)

    @staticmethod
    def solve_banded(l_and_u: Tuple[int, int], ab: np.ndarray, b: np.ndarray) -> np.ndarray:
        from scipy.linalg import solve_banded
        return solve_banded(l_and_u, ab, b)

    @staticmethod
    def solve_circulant(c: np.ndarray, b: np.ndarray) -> np.ndarray:
        from scipy.linalg import solve_circulant
        return solve_circulant(c, b)

    @staticmethod
    def solve_toeplitz(c: np.ndarray, b: np.ndarray) -> np.ndarray:
        from scipy.linalg import solve_toeplitz
        return solve_toeplitz(c, b)

32.3.2 矩阵应用示例

python

class MatrixApplications:
    @staticmethod
    def least_squares_fit(x: np.ndarray, y: np.ndarray, degree: int = 1) -> np.ndarray:
        A = np.vander(x, degree + 1)
        coefficients, _, _, _ = np.linalg.lstsq(A, y, rcond=None)
        return coefficients

    @staticmethod
    def pca(X: np.ndarray, n_components: int = 2) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
        X_centered = X - np.mean(X, axis=0)
        cov_matrix = np.cov(X_centered, rowvar=False)
        eigenvalues, eigenvectors = np.linalg.eigh(cov_matrix)
        idx = np.argsort(eigenvalues)[::-1]
        eigenvalues = eigenvalues[idx]
        eigenvectors = eigenvectors[:, idx]
        components = eigenvectors[:, :n_components]
        X_transformed = np.dot(X_centered, components)
        return X_transformed, components, eigenvalues[:n_components]

    @staticmethod
    def svd_compression(image: np.ndarray, k: int) -> np.ndarray:
        if len(image.shape) == 3:
            compressed_channels = []
            for channel in range(image.shape[2]):
                U, S, Vt = np.linalg.svd(image[:, :, channel])
                compressed = np.dot(U[:, :k], np.dot(np.diag(S[:k]), Vt[:k, :]))
                compressed_channels.append(compressed)
            return np.stack(compressed_channels, axis=2)
        else:
            U, S, Vt = np.linalg.svd(image)
            return np.dot(U[:, :k], np.dot(np.diag(S[:k]), Vt[:k, :]))

    @staticmethod
    def solve_linear_equations(A: np.ndarray, b: np.ndarray) -> np.ndarray:
        return np.linalg.solve(A, b)

    @staticmethod
    def matrix_power(A: np.ndarray, n: int) -> np.ndarray:
        return np.linalg.matrix_power(A, n)

    @staticmethod
    def exponential(A: np.ndarray) -> np.ndarray:
        from scipy.linalg import expm
        return expm(A)

32.4 SciPy科学计算

32.4.1 优化与拟合

python

from scipy.optimize import minimize, curve_fit, root_scalar, brentq
from scipy.optimize import linprog, minimize_scalar
from typing import Callable


class OptimizationSolver:
    @staticmethod
    def minimize_scalar(
        func: Callable[[float], float],
        bounds: Optional[Tuple[float, float]] = None,
        method: str = "brent"
    ) -> dict:
        if bounds:
            result = minimize_scalar(func, bounds=bounds, method="bounded")
        else:
            result = minimize_scalar(func, method=method)
        return {"x": result.x, "fun": result.fun, "success": result.success}

    @staticmethod
    def minimize_multivariate(
        func: Callable[[np.ndarray], float],
        x0: np.ndarray,
        method: str = "BFGS",
        bounds: Optional[List[Tuple[float, float]]] = None,
        constraints: Optional[List[dict]] = None
    ) -> dict:
        result = minimize(func, x0, method=method, bounds=bounds, constraints=constraints)
        return {
            "x": result.x,
            "fun": result.fun,
            "success": result.success,
            "nit": result.nit
        }

    @staticmethod
    def curve_fit_wrapper(
        func: Callable,
        xdata: np.ndarray,
        ydata: np.ndarray,
        p0: Optional[np.ndarray] = None,
        bounds: Optional[Tuple[np.ndarray, np.ndarray]] = None
    ) -> Tuple[np.ndarray, np.ndarray]:
        popt, pcov = curve_fit(func, xdata, ydata, p0=p0, bounds=bounds)
        return popt, np.sqrt(np.diag(pcov))

    @staticmethod
    def find_root(
        func: Callable[[float], float],
        bracket: Tuple[float, float],
        method: str = "brentq"
    ) -> float:
        if method == "brentq":
            return brentq(func, bracket[0], bracket[1])
        else:
            result = root_scalar(func, bracket=bracket, method=method)
            return result.root

    @staticmethod
    def linear_programming(
        c: np.ndarray,
        A_ub: Optional[np.ndarray] = None,
        b_ub: Optional[np.ndarray] = None,
        A_eq: Optional[np.ndarray] = None,
        b_eq: Optional[np.ndarray] = None,
        bounds: Optional[List[Tuple[float, float]]] = None
    ) -> dict:
        result = linprog(c, A_ub=A_ub, b_ub=b_ub, A_eq=A_eq, b_eq=b_eq, bounds=bounds)
        return {
            "x": result.x,
            "fun": result.fun,
            "success": result.success
        }


class CurveFitting:
    @staticmethod
    def polynomial_fit(x: np.ndarray, y: np.ndarray, degree: int) -> Tuple[np.ndarray, np.ndarray]:
        coefficients = np.polyfit(x, y, degree)
        polynomial = np.poly1d(coefficients)
        return coefficients, polynomial

    @staticmethod
    def exponential_fit(x: np.ndarray, y: np.ndarray) -> Tuple[Tuple[float, float], Callable]:
        def exp_func(x, a, b):
            return a * np.exp(b * x)

        log_y = np.log(y[y > 0])
        x_valid = x[y > 0]
        b, log_a = np.polyfit(x_valid, log_y, 1)
        a = np.exp(log_a)

        return (a, b), lambda x: exp_func(x, a, b)

    @staticmethod
    def power_law_fit(x: np.ndarray, y: np.ndarray) -> Tuple[Tuple[float, float], Callable]:
        def power_func(x, a, b):
            return a * np.power(x, b)

        mask = (x > 0) & (y > 0)
        log_x = np.log(x[mask])
        log_y = np.log(y[mask])
        b, log_a = np.polyfit(log_x, log_y, 1)
        a = np.exp(log_a)

        return (a, b), lambda x: power_func(x, a, b)

    @staticmethod
    def gaussian_fit(x: np.ndarray, y: np.ndarray) -> Tuple[Tuple[float, float, float], Callable]:
        def gaussian(x, amplitude, mean, sigma):
            return amplitude * np.exp(-((x - mean) ** 2) / (2 * sigma ** 2))

        mean_guess = x[np.argmax(y)]
        sigma_guess = (x.max() - x.min()) / 4
        amplitude_guess = y.max()

        popt, _ = curve_fit(
            gaussian, x, y,
            p0=[amplitude_guess, mean_guess, sigma_guess]
        )

        return tuple(popt), lambda x: gaussian(x, *popt)

32.4.2 插值与积分

python

from scipy.interpolate import interp1d, interp2d, CubicSpline, UnivariateSpline
from scipy.integrate import quad, dblquad, odeint, solve_ivp


class Interpolation:
    @staticmethod
    def linear_interpolation(x: np.ndarray, y: np.ndarray, x_new: np.ndarray) -> np.ndarray:
        f = interp1d(x, y, kind="linear")
        return f(x_new)

    @staticmethod
    def cubic_spline(x: np.ndarray, y: np.ndarray, x_new: np.ndarray) -> np.ndarray:
        cs = CubicSpline(x, y)
        return cs(x_new)

    @staticmethod
    def univariate_spline(
        x: np.ndarray,
        y: np.ndarray,
        s: float = 0,
        x_new: Optional[np.ndarray] = None
    ) -> Tuple[Callable, Optional[np.ndarray]]:
        spl = UnivariateSpline(x, y, s=s)
        if x_new is not None:
            return spl, spl(x_new)
        return spl, None

    @staticmethod
    def bilinear_interpolation(
        x: np.ndarray,
        y: np.ndarray,
        z: np.ndarray,
        x_new: np.ndarray,
        y_new: np.ndarray
    ) -> np.ndarray:
        f = interp2d(x, y, z, kind="linear")
        return f(x_new, y_new)


class NumericalIntegration:
    @staticmethod
    def integrate_1d(
        func: Callable[[float], float],
        a: float,
        b: float,
        **kwargs
    ) -> Tuple[float, float]:
        result, error = quad(func, a, b, **kwargs)
        return result, error

    @staticmethod
    def integrate_2d(
        func: Callable[[float, float], float],
        a: float,
        b: float,
        gfun: Callable[[float], float],
        hfun: Callable[[float], float],
        **kwargs
    ) -> Tuple[float, float]:
        result, error = dblquad(func, a, b, gfun, hfun, **kwargs)
        return result, error

    @staticmethod
    def trapezoidal(y: np.ndarray, x: Optional[np.ndarray] = None, dx: float = 1.0) -> float:
        return np.trapz(y, x=x, dx=dx)

    @staticmethod
    def simpson(y: np.ndarray, x: Optional[np.ndarray] = None, dx: float = 1.0) -> float:
        from scipy.integrate import simpson
        return simpson(y, x=x, dx=dx)

    @staticmethod
    def cumulative_trapezoid(y: np.ndarray, x: Optional[np.ndarray] = None, dx: float = 1.0) -> np.ndarray:
        from scipy.integrate import cumulative_trapezoid
        return cumulative_trapezoid(y, x=x, dx=dx, initial=0)


class DifferentialEquations:
    @staticmethod
    def solve_ode(
        func: Callable[[float, np.ndarray], np.ndarray],
        y0: np.ndarray,
        t_span: Tuple[float, float],
        t_eval: Optional[np.ndarray] = None,
        method: str = "RK45"
    ) -> dict:
        result = solve_ivp(func, t_span, y0, method=method, t_eval=t_eval)
        return {
            "t": result.t,
            "y": result.y,
            "success": result.success
        }

    @staticmethod
    def solve_ode_system(
        func: Callable[[float, np.ndarray], np.ndarray],
        y0: np.ndarray,
        t: np.ndarray
    ) -> np.ndarray:
        return odeint(func, y0, t)

    @staticmethod
    def harmonic_oscillator(omega: float, damping: float = 0.0):
        def equations(t, y):
            x, v = y
            return [v, -omega**2 * x - damping * v]
        return equations

    @staticmethod
    def lotka_volterra(alpha: float, beta: float, gamma: float, delta: float):
        def equations(t, y):
            prey, predator = y
            dprey = alpha * prey - beta * prey * predator
            dpredator = delta * prey * predator - gamma * predator
            return [dprey, dpredator]
        return equations

32.4.3 信号处理

python

from scipy import signal
from scipy.fft import fft, ifft, fft2, ifft2, fftfreq


class SignalProcessing:
    @staticmethod
    def fft_1d(data: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        n = len(data)
        freq = fftfreq(n)
        spectrum = fft(data)
        return freq, spectrum

    @staticmethod
    def ifft_1d(spectrum: np.ndarray) -> np.ndarray:
        return ifft(spectrum)

    @staticmethod
    def fft_2d(data: np.ndarray) -> np.ndarray:
        return fft2(data)

    @staticmethod
    def ifft_2d(spectrum: np.ndarray) -> np.ndarray:
        return ifft2(spectrum)

    @staticmethod
    def lowpass_filter(data: np.ndarray, cutoff: float, fs: float, order: int = 5) -> np.ndarray:
        nyquist = 0.5 * fs
        normal_cutoff = cutoff / nyquist
        b, a = signal.butter(order, normal_cutoff, btype="low", analog=False)
        return signal.filtfilt(b, a, data)

    @staticmethod
    def highpass_filter(data: np.ndarray, cutoff: float, fs: float, order: int = 5) -> np.ndarray:
        nyquist = 0.5 * fs
        normal_cutoff = cutoff / nyquist
        b, a = signal.butter(order, normal_cutoff, btype="high", analog=False)
        return signal.filtfilt(b, a, data)

    @staticmethod
    def bandpass_filter(
        data: np.ndarray,
        lowcut: float,
        highcut: float,
        fs: float,
        order: int = 5
    ) -> np.ndarray:
        nyquist = 0.5 * fs
        low = lowcut / nyquist
        high = highcut / nyquist
        b, a = signal.butter(order, [low, high], btype="band", analog=False)
        return signal.filtfilt(b, a, data)

    @staticmethod
    def spectrogram(data: np.ndarray, fs: float, nperseg: int = 256) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
        freqs, times, Sxx = signal.spectrogram(data, fs, nperseg=nperseg)
        return freqs, times, Sxx

    @staticmethod
    def find_peaks(data: np.ndarray, height: Optional[float] = None, distance: Optional[int] = None) -> np.ndarray:
        peaks, _ = signal.find_peaks(data, height=height, distance=distance)
        return peaks

    @staticmethod
    def resample(data: np.ndarray, num: int) -> np.ndarray:
        return signal.resample(data, num)

    @staticmethod
    def convolve(signal1: np.ndarray, signal2: np.ndarray, mode: str = "full") -> np.ndarray:
        return np.convolve(signal1, signal2, mode)

    @staticmethod
    def correlate(signal1: np.ndarray, signal2: np.ndarray, mode: str = "full") -> np.ndarray:
        return np.correlate(signal1, signal2, mode)

32.5 随机数与概率

32.5.1 随机数生成

python

class RandomGenerator:
    def __init__(self, seed: Optional[int] = None):
        self.rng = np.random.default_rng(seed)

    def uniform(self, low: float = 0.0, high: float = 1.0, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.uniform(low, high, size)

    def normal(self, loc: float = 0.0, scale: float = 1.0, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.normal(loc, scale, size)

    def integers(self, low: int, high: Optional[int] = None, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.integers(low, high, size)

    def choice(self, a: Union[int, np.ndarray], size: Optional[Tuple[int, ...]] = None, replace: bool = True, p: Optional[np.ndarray] = None) -> np.ndarray:
        return self.rng.choice(a, size=size, replace=replace, p=p)

    def permutation(self, x: Union[int, np.ndarray]) -> np.ndarray:
        return self.rng.permutation(x)

    def shuffle(self, x: np.ndarray) -> None:
        self.rng.shuffle(x)

    def poisson(self, lam: float = 1.0, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.poisson(lam, size)

    def exponential(self, scale: float = 1.0, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.exponential(scale, size)

    def binomial(self, n: int, p: float, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.binomial(n, p, size)

    def beta(self, a: float, b: float, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.beta(a, b, size)

    def gamma(self, shape: float, scale: float = 1.0, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.gamma(shape, scale, size)

    def chisquare(self, df: float, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.chisquare(df, size)

    def standard_t(self, df: float, size: Optional[Tuple[int, ...]] = None) -> np.ndarray:
        return self.rng.standard_t(df, size)


class ProbabilityDistributions:
    @staticmethod
    def pdf_normal(x: np.ndarray, mean: float = 0.0, std: float = 1.0) -> np.ndarray:
        from scipy.stats import norm
        return norm.pdf(x, mean, std)

    @staticmethod
    def cdf_normal(x: np.ndarray, mean: float = 0.0, std: float = 1.0) -> np.ndarray:
        from scipy.stats import norm
        return norm.cdf(x, mean, std)

    @staticmethod
    def pdf_binomial(k: np.ndarray, n: int, p: float) -> np.ndarray:
        from scipy.stats import binom
        return binom.pmf(k, n, p)

    @staticmethod
    def pdf_poisson(k: np.ndarray, mu: float) -> np.ndarray:
        from scipy.stats import poisson
        return poisson.pmf(k, mu)

    @staticmethod
    def pdf_exponential(x: np.ndarray, scale: float = 1.0) -> np.ndarray:
        from scipy.stats import expon
        return expon.pdf(x, scale=scale)

    @staticmethod
    def pdf_uniform(x: np.ndarray, loc: float = 0.0, scale: float = 1.0) -> np.ndarray:
        from scipy.stats import uniform
        return uniform.pdf(x, loc=loc, scale=scale)

32.6 科学计算应用实例

32.6.1 物理模拟

python

class PhysicsSimulation:
    @staticmethod
    def projectile_motion(
        v0: float,
        angle: float,
        g: float = 9.81,
        dt: float = 0.01
    ) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
        angle_rad = np.radians(angle)
        vx = v0 * np.cos(angle_rad)
        vy = v0 * np.sin(angle_rad)

        t = [0]
        x = [0]
        y = [0]

        while y[-1] >= 0:
            t.append(t[-1] + dt)
            x.append(x[-1] + vx * dt)
            y.append(y[-1] + vy * dt)
            vy -= g * dt

        return np.array(t), np.array(x), np.array(y)

    @staticmethod
    def pendulum_simulation(
        length: float,
        initial_angle: float,
        initial_velocity: float = 0,
        g: float = 9.81,
        t_max: float = 10,
        dt: float = 0.01
    ) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
        def derivatives(t, state):
            theta, omega = state
            return [omega, -g / length * np.sin(theta)]

        t = np.arange(0, t_max, dt)
        initial_state = [np.radians(initial_angle), initial_velocity]

        solution = odeint(lambda y, t: derivatives(t, y), initial_state, t)

        theta = solution[:, 0]
        omega = solution[:, 1]
        x = length * np.sin(theta)
        y = -length * np.cos(theta)

        return t, x, y

    @staticmethod
    def heat_equation_1d(
        L: float,
        T: float,
        nx: int = 100,
        nt: int = 1000,
        alpha: float = 0.01,
        initial_condition: Optional[Callable[[np.ndarray], np.ndarray]] = None
    ) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
        dx = L / (nx - 1)
        dt = T / nt

        x = np.linspace(0, L, nx)
        u = np.zeros((nt + 1, nx))

        if initial_condition:
            u[0, :] = initial_condition(x)
        else:
            u[0, :] = np.sin(np.pi * x / L)

        r = alpha * dt / dx**2

        for n in range(nt):
            for i in range(1, nx - 1):
                u[n + 1, i] = u[n, i] + r * (u[n, i + 1] - 2 * u[n, i] + u[n, i - 1])
            u[n + 1, 0] = 0
            u[n + 1, -1] = 0

        return x, np.linspace(0, T, nt + 1), u

    @staticmethod
    def wave_equation_1d(
        L: float,
        T: float,
        c: float = 1.0,
        nx: int = 100,
        nt: int = 1000,
        initial_displacement: Optional[Callable[[np.ndarray], np.ndarray]] = None,
        initial_velocity: Optional[Callable[[np.ndarray], np.ndarray]] = None
    ) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
        dx = L / (nx - 1)
        dt = T / nt

        x = np.linspace(0, L, nx)
        u = np.zeros((nt + 1, nx))

        if initial_displacement:
            u[0, :] = initial_displacement(x)
        else:
            u[0, :] = np.sin(np.pi * x / L)

        if initial_velocity:
            v0 = initial_velocity(x)
            u[1, 1:-1] = u[0, 1:-1] + dt * v0[1:-1]
        else:
            u[1, 1:-1] = u[0, 1:-1]

        r = (c * dt / dx) ** 2

        for n in range(1, nt):
            for i in range(1, nx - 1):
                u[n + 1, i] = 2 * u[n, i] - u[n - 1, i] + r * (u[n, i + 1] - 2 * u[n, i] + u[n, i - 1])
            u[n + 1, 0] = 0
            u[n + 1, -1] = 0

        return x, np.linspace(0, T, nt + 1), u

32.6.2 数值方法

python

class NumericalMethods:
    @staticmethod
    def numerical_derivative(func: Callable[[float], float], x: float, h: float = 1e-5) -> float:
        return (func(x + h) - func(x - h)) / (2 * h)

    @staticmethod
    def numerical_gradient(func: Callable[[np.ndarray], float], x: np.ndarray, h: float = 1e-5) -> np.ndarray:
        grad = np.zeros_like(x)
        for i in range(len(x)):
            x_plus = x.copy()
            x_minus = x.copy()
            x_plus[i] += h
            x_minus[i] -= h
            grad[i] = (func(x_plus) - func(x_minus)) / (2 * h)
        return grad

    @staticmethod
    def numerical_hessian(func: Callable[[np.ndarray], float], x: np.ndarray, h: float = 1e-5) -> np.ndarray:
        n = len(x)
        hessian = np.zeros((n, n))
        for i in range(n):
            for j in range(n):
                x_pp = x.copy()
                x_pm = x.copy()
                x_mp = x.copy()
                x_mm = x.copy()

                x_pp[i] += h
                x_pp[j] += h
                x_pm[i] += h
                x_pm[j] -= h
                x_mp[i] -= h
                x_mp[j] += h
                x_mm[i] -= h
                x_mm[j] -= h

                hessian[i, j] = (func(x_pp) - func(x_pm) - func(x_mp) + func(x_mm)) / (4 * h * h)
        return hessian

    @staticmethod
    def gradient_descent(
        func: Callable[[np.ndarray], float],
        grad_func: Callable[[np.ndarray], np.ndarray],
        x0: np.ndarray,
        learning_rate: float = 0.01,
        max_iter: int = 1000,
        tol: float = 1e-6
    ) -> Tuple[np.ndarray, List[float]]:
        x = x0.copy()
        history = [func(x)]

        for _ in range(max_iter):
            grad = grad_func(x)
            x = x - learning_rate * grad
            history.append(func(x))

            if np.linalg.norm(grad) < tol:
                break

        return x, history

    @staticmethod
    def newton_method(
        func: Callable[[np.ndarray], float],
        grad_func: Callable[[np.ndarray], np.ndarray],
        hess_func: Callable[[np.ndarray], np.ndarray],
        x0: np.ndarray,
        max_iter: int = 100,
        tol: float = 1e-6
    ) -> Tuple[np.ndarray, List[float]]:
        x = x0.copy()
        history = [func(x)]

        for _ in range(max_iter):
            grad = grad_func(x)
            hess = hess_func(x)

            try:
                delta = np.linalg.solve(hess, grad)
            except np.linalg.LinAlgError:
                break

            x = x - delta
            history.append(func(x))

            if np.linalg.norm(delta) < tol:
                break

        return x, history

32.7 知识图谱

32.7.1 NumPy数组架构

NumPy数组内存布局

┌─────────────────────────────────────────────────────────────┐
│                    ndarray对象                              │
├─────────────────────────────────────────────────────────────┤
│ data      数据指针                                         │
│ dtype     数据类型                                         │
│ shape     形状维度                                         │
│ strides   步长                                             │
│ flags     内存布局标志                                     │
└─────────────────────────────────────────────────────────────┘

数据类型:
┌─────────────────────────────────────────┐
│ int8/16/32/64    整数类型              │
│ float16/32/64    浮点类型              │
│ complex64/128    复数类型              │
│ bool             布尔类型              │
│ object           Python对象            │
└─────────────────────────────────────────┘

广播规则:
┌─────────────────────────────────────────┐
│ 1. 维度从右向左对齐                     │
│ 2. 缺失维度视为1                        │
│ 3. 维度为1可广播                        │
│ 4. 其他维度必须相等                     │
└─────────────────────────────────────────┘

32.7.2 SciPy模块结构

SciPy模块层次

┌─────────────────────────────────────────┐
│ scipy.linalg     线性代数              │
│ scipy.optimize   优化算法              │
│ scipy.integrate  数值积分              │
│ scipy.interpolate 插值                 │
│ scipy.fft        傅里叶变换            │
│ scipy.signal     信号处理              │
│ scipy.stats      统计函数              │
│ scipy.sparse     稀疏矩阵              │
│ scipy.ndimage    多维图像              │
│ scipy.io         数据输入输出          │
└─────────────────────────────────────────┘

32.8 技术选型指南

32.8.1 数值计算库选型

场景	推荐库	原因
通用数组计算	NumPy	标准选择
科学计算	SciPy	算法丰富
大规模数据	Dask	分布式
GPU加速	CuPy	CUDA支持

32.8.2 线性代数选型

操作规模	推荐方案	说明
小规模	numpy.linalg	简单直接
中规模	scipy.linalg	功能更多
大规模	scipy.sparse	稀疏矩阵
超大规模	PETSc/SLEPc	并行计算

32.9 常见问题与解决方案

32.9.1 内存不足

python

# 问题：大数组内存溢出
# 解决方案：使用内存映射

import numpy as np

# 创建内存映射文件
arr = np.memmap('large_array.dat', dtype='float64', mode='w+', shape=(1000000, 1000))

# 或使用分块处理
def process_in_chunks(filename, chunk_size=10000):
    for i in range(0, total_size, chunk_size):
        chunk = np.load(filename, mmap_mode='r')[i:i+chunk_size]
        process(chunk)

32.9.2 广播错误

python

# 问题：形状不匹配
a = np.array([[1, 2, 3]])  # (1, 3)
b = np.array([[1], [2]])   # (2, 1)
# a + b  # 可以广播

# 解决方案：理解广播规则
# (1, 3) + (2, 1) → (2, 3)

32.9.3 性能优化

python

# 问题：循环太慢
# 解决方案：向量化操作

# 慢
result = []
for i in range(1000000):
    result.append(i ** 2)

# 快
result = np.arange(1000000) ** 2

32.10 本章小结

本章详细介绍了Python科学计算的核心概念和实践：

NumPy数组：创建、索引、切片、形状变换
数组运算：基本运算、广播机制、统计运算
线性代数：矩阵运算、分解、线性方程组求解
SciPy优化：最小化、曲线拟合、约束优化
插值与积分：数值插值、数值积分、微分方程
信号处理：FFT、滤波、频谱分析
随机数生成：各种概率分布的随机数
应用实例：物理模拟、数值方法

练习题

实现一个矩阵运算库，支持加减乘除、转置、行列式计算
使用最小二乘法拟合多项式曲线
实现PCA降维算法并应用于数据可视化
求解常微分方程组，模拟捕食者-猎物系统
实现快速傅里叶变换并进行频谱分析

第32章 科学计算 ​

学习目标 ​

32.1 NumPy基础 ​

32.1.1 数组创建 ​

32.1.2 数组索引与切片 ​

32.1.3 数组形状操作 ​

32.2 数组运算 ​

32.2.1 基本运算 ​

32.2.2 广播机制 ​

32.2.3 统计运算 ​

32.3 线性代数 ​

32.3.1 矩阵运算 ​

32.3.2 矩阵应用示例 ​

32.4 SciPy科学计算 ​

32.4.1 优化与拟合 ​

32.4.2 插值与积分 ​

32.4.3 信号处理 ​

32.5 随机数与概率 ​

32.5.1 随机数生成 ​

32.6 科学计算应用实例 ​

32.6.1 物理模拟 ​

32.6.2 数值方法 ​

32.7 知识图谱 ​

32.7.1 NumPy数组架构 ​

32.7.2 SciPy模块结构 ​

32.8 技术选型指南 ​

32.8.1 数值计算库选型 ​

32.8.2 线性代数选型 ​

32.9 常见问题与解决方案 ​

32.9.1 内存不足 ​

32.9.2 广播错误 ​

32.9.3 性能优化 ​

32.10 本章小结 ​

练习题 ​

扩展阅读 ​

第32章科学计算

学习目标

32.1 NumPy基础

32.1.1 数组创建

32.1.2 数组索引与切片

32.1.3 数组形状操作

32.2 数组运算

32.2.1 基本运算

32.2.2 广播机制

32.2.3 统计运算

32.3 线性代数

32.3.1 矩阵运算

32.3.2 矩阵应用示例

32.4 SciPy科学计算

32.4.1 优化与拟合

32.4.2 插值与积分

32.4.3 信号处理

32.5 随机数与概率

32.5.1 随机数生成

32.6 科学计算应用实例

32.6.1 物理模拟

32.6.2 数值方法

32.7 知识图谱

32.7.1 NumPy数组架构

32.7.2 SciPy模块结构

32.8 技术选型指南

32.8.1 数值计算库选型

32.8.2 线性代数选型

32.9 常见问题与解决方案

32.9.1 内存不足

32.9.2 广播错误

32.9.3 性能优化

32.10 本章小结

练习题

扩展阅读