第1章:计算模拟的数学基础
说实话,很多人一听到「数学基础」四个字就想打退堂鼓。我当年也是这样。但做了十几年模拟项目后,我慢慢发现——真正在代码里用到的数学,其实就那么几块。今天咱们就把它们捋清楚。
核心观点:计算模拟的本质,就是用数学语言描述物理世界。你不需要成为数学家,但得知道工具箱里有什么工具。
1.1 微积分基础:变化率的语言
微积分这东西,说白了就是研究「变化」的学问。模拟里最常见的就是跟踪某个量随时间怎么变。
导数:瞬间的变化率
导数描述的是「这一刻的变化速度」。比如你开车,速度表显示的就是位移对时间的导数。
# 数值求导的简单实现
def derivative(f, x, h=1e-6):
return (f(x + h) - f(x - h)) / (2 * h)
# 例子:求 x^2 在 x=3 处的导数
f = lambda x: x**2
print(derivative(f, 3)) # 输出接近 6.0
我的经验:实际项目中,我很少手算导数。大部分时候用自动微分库(比如PyTorch的autograd)。但理解导数的几何意义——切线的斜率——能帮你调试很多奇怪的问题。
积分:累积的总量
积分是导数的逆运算。如果说导数看「速度」,积分就算「总路程」。模拟里经常用积分来累加某个量,比如总热量、总流量。
# 数值积分:梯形法则
def integrate(f, a, b, n=1000):
h = (b - a) / n
total = 0.5 * (f(a) + f(b))
for i in range(1, n):
total += f(a + i * h)
return total * h
# 计算 sin(x) 从 0 到 pi 的积分
import math
print(integrate(math.sin, 0, math.pi)) # 接近 2.0
避坑指南:我曾经在模拟流体时,积分步长选得太大,结果总质量不守恒。后来才意识到——数值积分不是精确的,步长越小越准,但计算量也越大。这是个经典的权衡问题。
1.2 线性代数基础:高维数据的骨架
你想想看,模拟里处理的数据很少是单个数字。温度场有几千个点,每个点有坐标和温度值——这就是向量和矩阵的用武之地。
向量:数据的自然表示
向量就是一串有序的数字。在模拟里,一个粒子的位置是三维向量,一个网格点的物理量也是向量。
# 向量运算示例
import numpy as np
v1 = np.array([1.0, 2.0, 3.0])
v2 = np.array([4.0, 5.0, 6.0])
# 点积:衡量两个向量的相似度
dot_product = np.dot(v1, v2)
print(f"点积: {dot_product}")
# 叉积:得到垂直于两个向量的新向量
cross_product = np.cross(v1, v2)
print(f"叉积: {cross_product}")
矩阵:线性变换的载体
矩阵可以看作是对向量做「操作」的工具。旋转、缩放、投影——这些在图形学和物理模拟里天天用。
# 矩阵乘法:旋转二维向量
import numpy as np
theta = np.pi / 4 # 45度
rotation_matrix = np.array([
[np.cos(theta), -np.sin(theta)],
[np.sin(theta), np.cos(theta)]
])
vector = np.array([1.0, 0.0])
rotated = rotation_matrix @ vector
print(f"旋转后的向量: {rotated}")
特征值:系统的「固有频率」
特征值和特征向量听起来高大上,其实它们描述的是:一个矩阵作用在某个方向上,只是拉伸或压缩,不改变方向。这在结构力学里特别重要——建筑物的固有频率就是特征值问题。
关键理解:特征值告诉你系统「最敏感」的方向。我在做地震模拟时,就是靠特征值分析找出建筑最容易损坏的位置。
1.3 概率论与数理统计:不确定性的量化
模拟不是算命。真实世界充满随机性——材料有缺陷、测量有误差、边界条件不确定。概率统计就是处理这些「不确定性」的工具。
基本概念:随机变量与分布
随机变量描述一个不确定的量。分布告诉我们这个量取不同值的可能性有多大。
# 生成正态分布随机数
import numpy as np
import matplotlib.pyplot as plt
# 均值0,标准差1的正态分布
samples = np.random.normal(0, 1, 10000)
# 计算统计量
mean = np.mean(samples)
std = np.std(samples)
print(f"均值: {mean:.3f}, 标准差: {std:.3f}")
蒙特卡洛方法:用随机数解决确定性问题
这个方法我特别喜欢。它的思路很简单:如果一个问题太难直接求解,那就随机采样,用统计结果逼近真实值。
# 用蒙特卡洛方法估算圆周率
import random
def estimate_pi(num_points=100000):
inside = 0
for _ in range(num_points):
x = random.uniform(-1, 1)
y = random.uniform(-1, 1)
if x*x + y*y <= 1:
inside += 1
return 4 * inside / num_points
print(f"估算的圆周率: {estimate_pi()}")
我的习惯:做蒙特卡洛模拟时,我总会先跑一个小样本看看收敛情况。如果结果波动太大,就增加采样数。别一上来就跑百万级样本——那是浪费计算资源。
1.4 数值计算方法简介
理论公式很美,但计算机只能做加减乘除。数值计算方法就是「怎么把连续的问题离散化,让计算机能算」。
离散化:从连续到离散
现实世界是连续的——温度场处处有值。但计算机只能存有限个点。离散化就是把连续空间切成网格,每个网格存一个值。
# 一维热传导方程的离散化
# 理论: dT/dt = alpha * d²T/dx²
# 离散: T_new[i] = T[i] + dt * alpha * (T[i+1] - 2*T[i] + T[i-1]) / dx²
def heat_1d(T, alpha, dt, dx):
n = len(T)
T_new = T.copy()
for i in range(1, n-1):
T_new[i] = T[i] + dt * alpha * (T[i+1] - 2*T[i] + T[i-1]) / (dx*dx)
return T_new
迭代法:逐步逼近真解
很多方程没有解析解,只能一步步逼近。迭代法的核心思想是:猜一个初始值,然后不断修正,直到结果足够精确。
# 牛顿法求方程的根
def newton(f, df, x0, tol=1e-6, max_iter=100):
x = x0
for i in range(max_iter):
fx = f(x)
if abs(fx) < tol:
return x
x = x - fx / df(x)
return x
# 求 x² - 2 = 0 的根(即 sqrt(2))
f = lambda x: x**2 - 2
df = lambda x: 2*x
root = newton(f, df, 1.0)
print(f"sqrt(2) ≈ {root:.6f}")
避坑指南:我曾经用牛顿法求解一个非线性方程组,结果迭代发散。后来发现初始值选得太离谱。记住——迭代法不保证收敛,初始猜测很关键。如果发散,换个初始值或者用更稳健的二分法。
知识体系总览
下面这张图把本章的核心内容串起来了。你可以把它当作一个「数学工具箱」的索引——遇到什么问题,就去翻对应的工具。
嗯,以上就是计算模拟需要的数学基础。别被这些名词吓到——实际用起来,每个工具都有固定的套路。后面几章我们会把这些工具一个个用起来,到时候你就知道它们有多实用了。
公众号:蓝海资料掘金营,微信deep3321