简介北京理工大学2009级《数值分析》期末试题及答案PDF适合计算机学院本科生及备考数值分析课程的学生使用。内容覆盖数值解法、插值法、迭代法、线性方程组求解等核心章节包含填空题、判断题和计算题多种题型答案附有详细求解过程便于对照自测与查漏补缺。资源共1个PDF文件压缩包大小约91KB体量小巧、便于下载后随时翻阅。目前已有259人学习浏览对期末复习和知识点整理具有一定参考价值。通过这份试题读者可以系统回顾有效数字、差商公式、样条插值、牛顿下山法、线性方程组迭代收敛条件等高频考点并结合列主元高斯消元、牛顿-埃尔米特插值、雅可比与高斯-赛德尔迭代等典型计算题检验掌握程度适合考前冲刺和专题强化。1. 数值分析不是背答案的学科试题是校准算法判断的标尺“数值分析”是个有点奇特的领域它在数学系里归计算数学在工科里又常以“数值计算”的面目出现试卷看起来都是手算推导但背后真正的问题永远是“这个算法到计算机上跑起来会偏差多少”。我前几天整理资料时又看到一份标题为“北京理工大学2009级数值分析试题及答案.pdf”的旧试卷第一反应不是去翻答案而是想把它当成一张算法边界检验清单。对搞 IT 或偏工程的读者来说这比记答案更值。试题覆盖插值、拟合、数值积分、常微分方程初值问题、线性方程组求解每一类都对应一个可编程的算法你真正需要做的是把题面公式变成可复现代码再用误差分析验证实现到底有没有写错。2. 复习数值分析先立骨架误差类型、收敛阶与题型映射拿到任何一份数值分析试卷最高效的用法不是从头做到尾而是先扫一眼题型分布。误差分析、插值拟合、数值积分、ODE、线性方程组几乎是固定班底区别只在参数和侧重点。如果你手里同时有教材、课件和往年试题我会建议先画一张题型映射表把“理论概念”和“可运行参数”对应起来再动手写代码。2.1 绝对误差、相对误差与有效数字一道题里的三个坑误差分析通常排在试卷最前面定义简单但最容易在小数位和数量级上翻车。绝对误差是近似值与真实值之差相对误差是绝对误差除以真实值有效数字则从第一位非零数字开始数。一个常见误用是拿绝对误差直接判断结果“好不好”而不看真实值的量级。比如测距时误差 0.1 米对 1000 米量级很小对 1 米量级就是不可接受。编程里这个概念直接对应浮点精度。np.float64大约能保证 15 到 16 位十进制有效数字np.float32只有约 7 位。你可以用一个小脚本体会import numpy as np x_true 3.14159 x_approx 3.14 abs_err abs(x_true - x_approx) rel_err abs_err / x_true digits np.floor(np.log10(abs(x_true) / abs_err)).astype(int) print(fabs_err {abs_err:.6g}) print(frel_err {rel_err:.6g}) print(f粗略有效数字位数 {digits})这段代码里值得留意的不是打印结果而是估算有效数字的方式。log10(真实值 / 绝对误差)只能给出数量级上的估计严格定义还要考虑小数位和四舍五入规则。工程上更常用np.finfo(np.float64).eps来了解机器精度而不是手工判断有效位数。回答试卷题时规范写法是先写科学计数法再数非零位到末位之间的位数两者不要混用。2.2 收敛速度与数值稳定性先于算法存在凡是涉及步长或迭代次数的题目绕不开“收敛阶”。收敛阶描述的是当步长 h 减小时误差下降的速度如果步长减半后误差缩小为原来的约 1/2^p就说该方法是 p 阶收敛。很多教材会把 p 的推导过程用泰勒展开写满一整页但实际验证要简单得多跑一遍 h再跑一遍 h/2看误差比值取对数。这个做法在“数值分析实验”里远比背定理有用。稳定性是另一个更隐蔽的坑。同一个算法在无限精度下收敛不代表在浮点有限字长下稳定。典型例子是中心差分求导公式(f(xh) - f(x-h)) / (2h)理论上误差随 h 减小而减小但 h 太小时两个相近数相减会把舍入误差放大到不可接受。试卷题目如果只让你写公式通常不会暴露这个问题可只要在计算机上做实验就必须寻找一个最优步长通常落在eps^(1/3)附近。这种地方恰恰是我在复习数值分析时觉得最值钱的部分背公式只是第一步知道公式在什么条件下失效才是真正掌握。2.3 用教材目录做题型映射表东南大学数值分析教材电子版与考题对照如果你参考的是常见的《东南大学数值分析教材电子版》目录结构大约包含误差理论、插值法、曲线拟合、数值积分、常微分方程数值解、线性方程组直接解法与迭代解法。对照一份类似“北京理工大学2009级数值分析试题及答案.pdf”的题目我会建立下面这张映射表教材章节典型试题题型常考察参数代码对应操作误差理论有效数字、误差限机器精度 epsnp.finfofloat32/64插值法牛顿插值、拉格朗日插值节点顺序、差商表np.vanderscipy.interpolate曲线拟合最小二乘拟合基函数、条件数np.linalg.lstsq数值积分复合梯形、Simpson 公式分割数 n、奇点位置scipy.integrate.quad常微分方程显式欧拉、RK4步长 h、局部截断误差solve_ivp手写 RK4线性方程组直接法、迭代法条件数、谱半径scipy.linalg.solve幂迭代这张表是我复习时打印出来贴在屏幕边的。每做一道题就先判断它属于哪一行然后挑出该行的参数去写最小用例。比对着答案抄一遍更有效因为你实际是在用代码重放一遍试卷的设计逻辑。所谓“把试卷当成实验报告”说的就是这个流程。3. 从“数值分析试题”到可运行代码插值与拟合的可复现实现试题和真实代码之间的落差往往不在公式本身而在下标和求值顺序。比如牛顿插值的差商表手算时可以很规整地列出来但一到 Python 循环里i和j就特别容易错位。我的处理方法是先保留一个明确的数学定义再写一个最小例子最后用 SciPy 或手工计算交叉验证。3.1 牛顿插值最容易写错的下标牛顿插值的优势在于增加节点时可以复用之前的差商这对纸质试卷很友好因为老师可以按步骤给分。换成代码后核心是差商表的填充顺序import numpy as np def newton_interp(x, y): n len(x) dd np.zeros((n, n)) dd[:, 0] y for j in range(1, n): for i in range(n - j): dd[i, j] (dd[i 1, j - 1] - dd[i, j - 1]) / (x[i j] - x[i]) coeff dd[0, :].copy() def evaluate(t): ans coeff[-1] for k in range(len(coeff) - 1, 0, -1): ans ans * (t - x[k - 1]) coeff[k - 1] return ans return evaluate, dd x np.array([0., 1., 2., 3.]) y np.array([1., 2., 0., 5.]) p, dd newton_interp(x, y) for t in [0.5, 1.5, 2.5]: print(fp({t}) {p(t):.6f})外层循环j表示差商阶数内层循环i表示起点位置。dd[i, j]表示以x[i]为起点、步长为j的 j 阶差商。熟练之后可以将它压缩成只保存两行但考试和调试阶段用二维数组更直观。evaluate里使用的 Horner 算法把多项式求值从多次乘方变成连续乘加避免重复计算。重点关注节点的分布如果x[ij]与x[i]非常接近分母会极小并放大舍入误差。实际工程中应先排序、去重再决定是否换用分段插值。3.2 最小二乘拟合先定基函数再看条件数插值要求曲线穿过所有已知点拟合则相反只要求整体误差最小。试卷里常见的是多项式最小二乘编程时要额外注意基函数的选择和矩阵条件数。直接解正规方程(A^T A) c A^T y在数学上等价但数值上会放大条件数我更推荐直接调用基于 QR 分解或 SVD 的np.linalg.lstsqdef least_squares(x, y, deg): A np.vander(x, deg 1) coeff, _, _, _ np.linalg.lstsq(A, y, rcondNone) cond np.linalg.cond(A) return coeff, cond x np.array([0., 1., 2., 3.]) y np.array([1., 2., 0., 5.]) coeff, cond least_squares(x, y, 2) print(coeff , coeff) print(cond(A) , cond)np.vander生成的是幂基矩阵列方向从高次到低次。对于范围很大的数据比如x在 0 到 1000000 之间cond(A)会迅速超过 1e12此时即使拟合残差很小系数也可能严重失稳。数值分析教材里会用平移和缩放处理实际代码中更简单的做法是把x归一化到[-1, 1]再做拟合最后再映射回原始坐标。试题数据通常比较温和但你如果直接拿这套代码去处理真实生产数据条件数检查不能省。3.3 参数说明与误差观察用残差而不是点到点差拟合质量不能只看某几个点的测试误差我一般会同时看三个指标最大绝对残差、均方根误差、系数向量的模。如果残差很小但系数很大大概率是过拟合也就是多项式次数选得过高。插值问题不需要这个担心因为插值多项式通过所有点但拟合问题必须做。residual y - np.polyval(coeff, x) print(max residual:, np.max(np.abs(residual))) print(rmse:, np.sqrt(np.mean(residual**2)))把次数从 1 逐渐提高到len(x) - 1观察 rmse 的下降速度。如果从二次升到三次之后 rmse 几乎不变说明多项式能力已经饱和再升次数只会让曲线在中间剧烈振荡。这个现象是“数值分析及其MATLAB实验”里最常见的思考题之一也是理解插值与拟合差别的关键场景。4. 数值分析实验里的数值积分与 ODE参数、步长与稳定性数值积分和常微分方程初值问题几乎是同一类思想把连续问题离散成有限步操作再研究误差如何随步长变化。考试题习惯让你手算一个复合梯形或 Simpson 公式然后用几步 RK4 求解一个简单 ODE。代码实现时重要的是把接口统一让不同算法能在同一组参数下直接比较。4.1 复合梯形与 Simpson 公式的精度边界复合梯形公式把所有子区间首尾相加复合 Simpson 公式则要求区间数必须为偶数。实现时最常见的错误是切片位置没对齐导致偶节点和奇节点算反。下面这段代码把两个函数写成相同签名方便直接对比import numpy as np def trapezoid(f, a, b, n): h (b - a) / n xs np.linspace(a, b, n 1) ys f(xs) return h * (0.5 * ys[0] np.sum(ys[1:-1]) 0.5 * ys[-1]) def simpson(f, a, b, n): n n (n % 2) # 保证偶数 h (b - a) / n xs np.linspace(a, b, n 1) ys f(xs) return h / 3 * (ys[0] ys[-1] 4 * np.sum(ys[1:-1:2]) 2 * np.sum(ys[2:-2:2])) f lambda x: np.sin(x) for n in [4, 8, 16]: t trapezoid(f, 0, np.pi, n) s simpson(f, 0, np.pi, n) print(fn{n}: trapezoid{t:.8f}, simpson{s:.8f})对sin(x)在 0 到 π 上积分真实值是 2。可以看到 Simpson 公式在 n16 时的误差远小于梯形公式。ys[1:-1:2]取的是奇数下标点ys[2:-2:2]取的是偶数下标点这个切片关系一旦写错结果会差距很大。更要注意的是端点梯形公式两边系数是 1/2Simpson 是 1/3 和 4/3 交替最后端点只出现一次。题目里如果给出表格形式的函数值建议先用np.asarray把数据读进来再按下标计算别手敲数值。4.2 常微分方程初值问题步长与稳定性MATLAB ode45 对应的参数ODE 数值解在试卷上常考显式欧拉和经典四阶 Runge-Kutta。显式欧拉写起来最简单但对步长的限制非常苛刻方程一旦稍带刚性步长稍大就会发散。四阶 RK4 的稳定域更大但也不是万能。MATLAB 里常用ode45它背后是 Dormand-Prince 对采用自适应步长内部参数RelTol默认 1e-3AbsTol默认 1e-6。答题时不能简单写“用 ode45”因为阅卷老师要看的是你能否给出固定步长下的迭代过程和误差估计。我自己调试试卷题时习惯手写一个定步长 RK4方便和答案对比def rk4(f, t0, y0, t_end, n): h (t_end - t0) / n t t0 y y0 hist [(t, y)] for _ in range(n): k1 f(t, y) k2 f(t h / 2, y h * k1 / 2) k3 f(t h / 2, y h * k2 / 2) k4 f(t h, y h * k3) y y h * (k1 2 * k2 2 * k3 k4) / 6 t t h hist.append((t, y)) return np.array(hist) # dy/dt -y真实解为 exp(-t) sol rk4(lambda t, y: -y, 0, 1, 2, 4) exact np.exp(-2) print(sol[-1], exact, abs(sol[-1, 1] - exact))参数含义比较直观t0是初始时间y0是初始状态t_end是终止时间n将区间等分为多少步。k1是起点斜率k2和k3是半步处的中心斜率k4是终点斜率。四阶的含义是局部截断误差为 O(h^5)整体误差为 O(h^4)。如果你把 n 改成 8 或 16会发现最后一步误差按接近 16 倍的比例缩小这是判断 RK4 实现是否正确的快速方法。4.3 一个可复现的小实验把步长减半观察数值分析里的误差数值分析实验的常见做法是把步长减半观察误差是否按理论阶数下降。对上面这个dy/dt -y的问题我通常这样组织for n in [10, 20, 40]: sol rk4(lambda t, y: -y, 0, 1, 2, n) err abs(sol[-1, 1] - np.exp(-2)) print(fn{n}, err{err:.3e})理论预期是 n 翻倍后误差缩小约 16 倍。如果结果只有 4 倍甚至更小说明问题不在公式而在步长已经进入舍入误差主导区或者代码里系数、下标写错。这也是我在做“数值分析实验”时最看重的验证方式用收敛阶来证明实现是对的而不是只看最终值是否像答案。5. 试题里的线性方程组直接法与迭代法的边界线性方程组在数值分析试卷里通常分两块直接法考 LU 分解、选主元、条件数迭代法考 Jacobi、Gauss-Seidel 以及收敛判定。两者不是替代关系而是适用场景不同。直接法适合中小规模稠密矩阵迭代法适合大规模稀疏矩阵。5.1 LU 分解与高斯消元看起来相同数值表现不同手算高斯消元和 LU 分解本质上是同一过程但程序实现时建议直接用scipy.linalg.lu它默认执行部分选主元并返回置换矩阵 P。求解过程分两步先解Ly P b再解Ux y。为了验证自己写的回代逻辑可以再用内置solve做参考解import numpy as np from scipy.linalg import lu, solve A np.array([[4., 3., 2.], [1., 3., 7.], [2., 1., 9.]]) b np.array([1., 2., 3.]) P, L, U lu(A) y solve(L, P b) x_lu solve(U, y) x_ref solve(A, b) print(x_lu) print(x_ref) print(diff , np.max(np.abs(x_lu - x_ref)))lu返回的 L 是单位下三角U 是上三角。如果 A 接近奇异U 的对角元会非常小此时求解结果对 b 的微小扰动异常敏感。判断指标是np.linalg.cond(A)条件数越大误差上界越大。解答试卷时如果题目要求“用 LU 分解求解”必须写出 L、U 和中间向量而不能只贴solve的调用结果。真实编程中我也至少会保留这个最小参考解用来排除自己实现回代时的边界错误。5.2 Jacobi、Gauss-Seidel 与收敛条件迭代法在代码上只比直接法多一个循环但收敛性质完全不同。Jacobi 更新时全部使用旧值Gauss-Seidel 则立刻使用最新值。一个最简单的 Jacobi 实现如下def jacobi(A, b, x0, tol1e-8, max_iter1000): D np.diag(A) L np.tril(A, -1) U np.triu(A, 1) x x0.copy() for it in range(max_iter): x_new (b - (L U) x) / D if np.linalg.norm(x_new - x, ordnp.inf) tol: return x_new, it 1 x x_new return x, max_iter A np.array([[4., 1., 1.], [1., 5., 2.], [1., 2., 4.]]) b np.array([1., 2., 3.]) x, it jacobi(A, b, np.zeros(3)) print(x , x, iterations , it)这里D是对角线向量(b - (LU) x) / D是逐元素除法向量化的写法比较紧凑。tol是迭代终止容差max_iter是防发散上限。试卷上常考的结论是对角占优矩阵保证 Jacobi 和 Gauss-Seidel 收敛但充分必要条件是迭代矩阵的谱半径小于 1。对角占优只是充分条件不是必要条件。遇到不占优但题目说收敛的矩阵不要怀疑题出错了应该直接算谱半径。5.3 用参考解验证实现相对残差比绝对误差更可靠在上一段代码中直接比较x与x_ref的绝对差值会受解的量级影响。更好的指标是相对残差r np.linalg.norm(A x - b) / (np.linalg.norm(A) * np.linalg.norm(x)) print(relative residual , r)这个值描述的是“误差相对于问题规模的占比”不会因为 b 变成 1000 倍而变化太多。实际工程里如果相对残差在 1e-10 量级通常可以认为迭代收敛到满意状态。若相对残差卡在 1e-6 下不去首先查条件数其次查容差设置最后再去怀疑实现。这个判断顺序也是“数值分析及其MATLAB实验”里隐含考察的内容。6. 最后一道“应用题”用残差与收敛阶验证你的数值实现无论试题属于哪一类我最终都会把它收敛到一个统一动作构造一个已知参考解的用例跑几步看误差是否按理论阶数下降。这个方法可以筛掉绝大多数下标记错、系数抄错、切片越界的低级 bug。具体做法是每次写一个新算法先不着急处理真实数据而是用解析解非常明确的小问题做冒烟测试。例如插值就选y x^3取 4 个节点插值结果必须在这些节点上精确等于原值。最小二乘就选一条直线加一点噪声比较系数是否在噪声范围内接近真值。数值积分就选sin(x)或exp(x)。ODE 就选dy/dt y。这些参考解都是临时的成本几乎为零但能立刻暴露循环边界是否多走了一位。真正值得投入精力的地方是收敛阶测试。以复合 Simpson 为例令 n 从 8 增加到 16误差应从某个值下降到约 1/16如果只下降到 1/4说明实现的不是 Simpson 而是某种二阶方法通常是因为4*ys[1:-1:2]和2*ys[2:-2:2]这两个切片的权重被写反了。对 RK4 同理n 翻倍误差下降约 16 倍如果只有 2 倍多半是 k1、k2、k3、k4 的系数写错了。用线性方程组的迭代法时验证方式则是固定同一初值把tol从 1e-6 改成 1e-10观察迭代次数是否按预期增加同时相对残差是否有明显下降。最后还有一个我自己的约定任何算法代码都必须保留一个--check模式运行小型参考用例并自动比较误差。这个习惯直接来自复习“数值分析试题及答案”时的经验答案可以告诉我结果但只有误差行为能告诉我为什么结果可信。如果你在做题或实验时遇到结果与答案不一致先不要改参数停下来跑一遍收敛阶测试通常下一步就清楚了。本文还有配套的精品资源点击获取 数值分析试题解析：从误差分析到线性方程组的可复现代码
简介北京理工大学2009级《数值分析》期末试题及答案PDF适合计算机学院本科生及备考数值分析课程的学生使用。内容覆盖数值解法、插值法、迭代法、线性方程组求解等核心章节包含填空题、判断题和计算题多种题型答案附有详细求解过程便于对照自测与查漏补缺。资源共1个PDF文件压缩包大小约91KB体量小巧、便于下载后随时翻阅。目前已有259人学习浏览对期末复习和知识点整理具有一定参考价值。通过这份试题读者可以系统回顾有效数字、差商公式、样条插值、牛顿下山法、线性方程组迭代收敛条件等高频考点并结合列主元高斯消元、牛顿-埃尔米特插值、雅可比与高斯-赛德尔迭代等典型计算题检验掌握程度适合考前冲刺和专题强化。1. 数值分析不是背答案的学科试题是校准算法判断的标尺“数值分析”是个有点奇特的领域它在数学系里归计算数学在工科里又常以“数值计算”的面目出现试卷看起来都是手算推导但背后真正的问题永远是“这个算法到计算机上跑起来会偏差多少”。我前几天整理资料时又看到一份标题为“北京理工大学2009级数值分析试题及答案.pdf”的旧试卷第一反应不是去翻答案而是想把它当成一张算法边界检验清单。对搞 IT 或偏工程的读者来说这比记答案更值。试题覆盖插值、拟合、数值积分、常微分方程初值问题、线性方程组求解每一类都对应一个可编程的算法你真正需要做的是把题面公式变成可复现代码再用误差分析验证实现到底有没有写错。2. 复习数值分析先立骨架误差类型、收敛阶与题型映射拿到任何一份数值分析试卷最高效的用法不是从头做到尾而是先扫一眼题型分布。误差分析、插值拟合、数值积分、ODE、线性方程组几乎是固定班底区别只在参数和侧重点。如果你手里同时有教材、课件和往年试题我会建议先画一张题型映射表把“理论概念”和“可运行参数”对应起来再动手写代码。2.1 绝对误差、相对误差与有效数字一道题里的三个坑误差分析通常排在试卷最前面定义简单但最容易在小数位和数量级上翻车。绝对误差是近似值与真实值之差相对误差是绝对误差除以真实值有效数字则从第一位非零数字开始数。一个常见误用是拿绝对误差直接判断结果“好不好”而不看真实值的量级。比如测距时误差 0.1 米对 1000 米量级很小对 1 米量级就是不可接受。编程里这个概念直接对应浮点精度。np.float64大约能保证 15 到 16 位十进制有效数字np.float32只有约 7 位。你可以用一个小脚本体会import numpy as np x_true 3.14159 x_approx 3.14 abs_err abs(x_true - x_approx) rel_err abs_err / x_true digits np.floor(np.log10(abs(x_true) / abs_err)).astype(int) print(fabs_err {abs_err:.6g}) print(frel_err {rel_err:.6g}) print(f粗略有效数字位数 {digits})这段代码里值得留意的不是打印结果而是估算有效数字的方式。log10(真实值 / 绝对误差)只能给出数量级上的估计严格定义还要考虑小数位和四舍五入规则。工程上更常用np.finfo(np.float64).eps来了解机器精度而不是手工判断有效位数。回答试卷题时规范写法是先写科学计数法再数非零位到末位之间的位数两者不要混用。2.2 收敛速度与数值稳定性先于算法存在凡是涉及步长或迭代次数的题目绕不开“收敛阶”。收敛阶描述的是当步长 h 减小时误差下降的速度如果步长减半后误差缩小为原来的约 1/2^p就说该方法是 p 阶收敛。很多教材会把 p 的推导过程用泰勒展开写满一整页但实际验证要简单得多跑一遍 h再跑一遍 h/2看误差比值取对数。这个做法在“数值分析实验”里远比背定理有用。稳定性是另一个更隐蔽的坑。同一个算法在无限精度下收敛不代表在浮点有限字长下稳定。典型例子是中心差分求导公式(f(xh) - f(x-h)) / (2h)理论上误差随 h 减小而减小但 h 太小时两个相近数相减会把舍入误差放大到不可接受。试卷题目如果只让你写公式通常不会暴露这个问题可只要在计算机上做实验就必须寻找一个最优步长通常落在eps^(1/3)附近。这种地方恰恰是我在复习数值分析时觉得最值钱的部分背公式只是第一步知道公式在什么条件下失效才是真正掌握。2.3 用教材目录做题型映射表东南大学数值分析教材电子版与考题对照如果你参考的是常见的《东南大学数值分析教材电子版》目录结构大约包含误差理论、插值法、曲线拟合、数值积分、常微分方程数值解、线性方程组直接解法与迭代解法。对照一份类似“北京理工大学2009级数值分析试题及答案.pdf”的题目我会建立下面这张映射表教材章节典型试题题型常考察参数代码对应操作误差理论有效数字、误差限机器精度 epsnp.finfofloat32/64插值法牛顿插值、拉格朗日插值节点顺序、差商表np.vanderscipy.interpolate曲线拟合最小二乘拟合基函数、条件数np.linalg.lstsq数值积分复合梯形、Simpson 公式分割数 n、奇点位置scipy.integrate.quad常微分方程显式欧拉、RK4步长 h、局部截断误差solve_ivp手写 RK4线性方程组直接法、迭代法条件数、谱半径scipy.linalg.solve幂迭代这张表是我复习时打印出来贴在屏幕边的。每做一道题就先判断它属于哪一行然后挑出该行的参数去写最小用例。比对着答案抄一遍更有效因为你实际是在用代码重放一遍试卷的设计逻辑。所谓“把试卷当成实验报告”说的就是这个流程。3. 从“数值分析试题”到可运行代码插值与拟合的可复现实现试题和真实代码之间的落差往往不在公式本身而在下标和求值顺序。比如牛顿插值的差商表手算时可以很规整地列出来但一到 Python 循环里i和j就特别容易错位。我的处理方法是先保留一个明确的数学定义再写一个最小例子最后用 SciPy 或手工计算交叉验证。3.1 牛顿插值最容易写错的下标牛顿插值的优势在于增加节点时可以复用之前的差商这对纸质试卷很友好因为老师可以按步骤给分。换成代码后核心是差商表的填充顺序import numpy as np def newton_interp(x, y): n len(x) dd np.zeros((n, n)) dd[:, 0] y for j in range(1, n): for i in range(n - j): dd[i, j] (dd[i 1, j - 1] - dd[i, j - 1]) / (x[i j] - x[i]) coeff dd[0, :].copy() def evaluate(t): ans coeff[-1] for k in range(len(coeff) - 1, 0, -1): ans ans * (t - x[k - 1]) coeff[k - 1] return ans return evaluate, dd x np.array([0., 1., 2., 3.]) y np.array([1., 2., 0., 5.]) p, dd newton_interp(x, y) for t in [0.5, 1.5, 2.5]: print(fp({t}) {p(t):.6f})外层循环j表示差商阶数内层循环i表示起点位置。dd[i, j]表示以x[i]为起点、步长为j的 j 阶差商。熟练之后可以将它压缩成只保存两行但考试和调试阶段用二维数组更直观。evaluate里使用的 Horner 算法把多项式求值从多次乘方变成连续乘加避免重复计算。重点关注节点的分布如果x[ij]与x[i]非常接近分母会极小并放大舍入误差。实际工程中应先排序、去重再决定是否换用分段插值。3.2 最小二乘拟合先定基函数再看条件数插值要求曲线穿过所有已知点拟合则相反只要求整体误差最小。试卷里常见的是多项式最小二乘编程时要额外注意基函数的选择和矩阵条件数。直接解正规方程(A^T A) c A^T y在数学上等价但数值上会放大条件数我更推荐直接调用基于 QR 分解或 SVD 的np.linalg.lstsqdef least_squares(x, y, deg): A np.vander(x, deg 1) coeff, _, _, _ np.linalg.lstsq(A, y, rcondNone) cond np.linalg.cond(A) return coeff, cond x np.array([0., 1., 2., 3.]) y np.array([1., 2., 0., 5.]) coeff, cond least_squares(x, y, 2) print(coeff , coeff) print(cond(A) , cond)np.vander生成的是幂基矩阵列方向从高次到低次。对于范围很大的数据比如x在 0 到 1000000 之间cond(A)会迅速超过 1e12此时即使拟合残差很小系数也可能严重失稳。数值分析教材里会用平移和缩放处理实际代码中更简单的做法是把x归一化到[-1, 1]再做拟合最后再映射回原始坐标。试题数据通常比较温和但你如果直接拿这套代码去处理真实生产数据条件数检查不能省。3.3 参数说明与误差观察用残差而不是点到点差拟合质量不能只看某几个点的测试误差我一般会同时看三个指标最大绝对残差、均方根误差、系数向量的模。如果残差很小但系数很大大概率是过拟合也就是多项式次数选得过高。插值问题不需要这个担心因为插值多项式通过所有点但拟合问题必须做。residual y - np.polyval(coeff, x) print(max residual:, np.max(np.abs(residual))) print(rmse:, np.sqrt(np.mean(residual**2)))把次数从 1 逐渐提高到len(x) - 1观察 rmse 的下降速度。如果从二次升到三次之后 rmse 几乎不变说明多项式能力已经饱和再升次数只会让曲线在中间剧烈振荡。这个现象是“数值分析及其MATLAB实验”里最常见的思考题之一也是理解插值与拟合差别的关键场景。4. 数值分析实验里的数值积分与 ODE参数、步长与稳定性数值积分和常微分方程初值问题几乎是同一类思想把连续问题离散成有限步操作再研究误差如何随步长变化。考试题习惯让你手算一个复合梯形或 Simpson 公式然后用几步 RK4 求解一个简单 ODE。代码实现时重要的是把接口统一让不同算法能在同一组参数下直接比较。4.1 复合梯形与 Simpson 公式的精度边界复合梯形公式把所有子区间首尾相加复合 Simpson 公式则要求区间数必须为偶数。实现时最常见的错误是切片位置没对齐导致偶节点和奇节点算反。下面这段代码把两个函数写成相同签名方便直接对比import numpy as np def trapezoid(f, a, b, n): h (b - a) / n xs np.linspace(a, b, n 1) ys f(xs) return h * (0.5 * ys[0] np.sum(ys[1:-1]) 0.5 * ys[-1]) def simpson(f, a, b, n): n n (n % 2) # 保证偶数 h (b - a) / n xs np.linspace(a, b, n 1) ys f(xs) return h / 3 * (ys[0] ys[-1] 4 * np.sum(ys[1:-1:2]) 2 * np.sum(ys[2:-2:2])) f lambda x: np.sin(x) for n in [4, 8, 16]: t trapezoid(f, 0, np.pi, n) s simpson(f, 0, np.pi, n) print(fn{n}: trapezoid{t:.8f}, simpson{s:.8f})对sin(x)在 0 到 π 上积分真实值是 2。可以看到 Simpson 公式在 n16 时的误差远小于梯形公式。ys[1:-1:2]取的是奇数下标点ys[2:-2:2]取的是偶数下标点这个切片关系一旦写错结果会差距很大。更要注意的是端点梯形公式两边系数是 1/2Simpson 是 1/3 和 4/3 交替最后端点只出现一次。题目里如果给出表格形式的函数值建议先用np.asarray把数据读进来再按下标计算别手敲数值。4.2 常微分方程初值问题步长与稳定性MATLAB ode45 对应的参数ODE 数值解在试卷上常考显式欧拉和经典四阶 Runge-Kutta。显式欧拉写起来最简单但对步长的限制非常苛刻方程一旦稍带刚性步长稍大就会发散。四阶 RK4 的稳定域更大但也不是万能。MATLAB 里常用ode45它背后是 Dormand-Prince 对采用自适应步长内部参数RelTol默认 1e-3AbsTol默认 1e-6。答题时不能简单写“用 ode45”因为阅卷老师要看的是你能否给出固定步长下的迭代过程和误差估计。我自己调试试卷题时习惯手写一个定步长 RK4方便和答案对比def rk4(f, t0, y0, t_end, n): h (t_end - t0) / n t t0 y y0 hist [(t, y)] for _ in range(n): k1 f(t, y) k2 f(t h / 2, y h * k1 / 2) k3 f(t h / 2, y h * k2 / 2) k4 f(t h, y h * k3) y y h * (k1 2 * k2 2 * k3 k4) / 6 t t h hist.append((t, y)) return np.array(hist) # dy/dt -y真实解为 exp(-t) sol rk4(lambda t, y: -y, 0, 1, 2, 4) exact np.exp(-2) print(sol[-1], exact, abs(sol[-1, 1] - exact))参数含义比较直观t0是初始时间y0是初始状态t_end是终止时间n将区间等分为多少步。k1是起点斜率k2和k3是半步处的中心斜率k4是终点斜率。四阶的含义是局部截断误差为 O(h^5)整体误差为 O(h^4)。如果你把 n 改成 8 或 16会发现最后一步误差按接近 16 倍的比例缩小这是判断 RK4 实现是否正确的快速方法。4.3 一个可复现的小实验把步长减半观察数值分析里的误差数值分析实验的常见做法是把步长减半观察误差是否按理论阶数下降。对上面这个dy/dt -y的问题我通常这样组织for n in [10, 20, 40]: sol rk4(lambda t, y: -y, 0, 1, 2, n) err abs(sol[-1, 1] - np.exp(-2)) print(fn{n}, err{err:.3e})理论预期是 n 翻倍后误差缩小约 16 倍。如果结果只有 4 倍甚至更小说明问题不在公式而在步长已经进入舍入误差主导区或者代码里系数、下标写错。这也是我在做“数值分析实验”时最看重的验证方式用收敛阶来证明实现是对的而不是只看最终值是否像答案。5. 试题里的线性方程组直接法与迭代法的边界线性方程组在数值分析试卷里通常分两块直接法考 LU 分解、选主元、条件数迭代法考 Jacobi、Gauss-Seidel 以及收敛判定。两者不是替代关系而是适用场景不同。直接法适合中小规模稠密矩阵迭代法适合大规模稀疏矩阵。5.1 LU 分解与高斯消元看起来相同数值表现不同手算高斯消元和 LU 分解本质上是同一过程但程序实现时建议直接用scipy.linalg.lu它默认执行部分选主元并返回置换矩阵 P。求解过程分两步先解Ly P b再解Ux y。为了验证自己写的回代逻辑可以再用内置solve做参考解import numpy as np from scipy.linalg import lu, solve A np.array([[4., 3., 2.], [1., 3., 7.], [2., 1., 9.]]) b np.array([1., 2., 3.]) P, L, U lu(A) y solve(L, P b) x_lu solve(U, y) x_ref solve(A, b) print(x_lu) print(x_ref) print(diff , np.max(np.abs(x_lu - x_ref)))lu返回的 L 是单位下三角U 是上三角。如果 A 接近奇异U 的对角元会非常小此时求解结果对 b 的微小扰动异常敏感。判断指标是np.linalg.cond(A)条件数越大误差上界越大。解答试卷时如果题目要求“用 LU 分解求解”必须写出 L、U 和中间向量而不能只贴solve的调用结果。真实编程中我也至少会保留这个最小参考解用来排除自己实现回代时的边界错误。5.2 Jacobi、Gauss-Seidel 与收敛条件迭代法在代码上只比直接法多一个循环但收敛性质完全不同。Jacobi 更新时全部使用旧值Gauss-Seidel 则立刻使用最新值。一个最简单的 Jacobi 实现如下def jacobi(A, b, x0, tol1e-8, max_iter1000): D np.diag(A) L np.tril(A, -1) U np.triu(A, 1) x x0.copy() for it in range(max_iter): x_new (b - (L U) x) / D if np.linalg.norm(x_new - x, ordnp.inf) tol: return x_new, it 1 x x_new return x, max_iter A np.array([[4., 1., 1.], [1., 5., 2.], [1., 2., 4.]]) b np.array([1., 2., 3.]) x, it jacobi(A, b, np.zeros(3)) print(x , x, iterations , it)这里D是对角线向量(b - (LU) x) / D是逐元素除法向量化的写法比较紧凑。tol是迭代终止容差max_iter是防发散上限。试卷上常考的结论是对角占优矩阵保证 Jacobi 和 Gauss-Seidel 收敛但充分必要条件是迭代矩阵的谱半径小于 1。对角占优只是充分条件不是必要条件。遇到不占优但题目说收敛的矩阵不要怀疑题出错了应该直接算谱半径。5.3 用参考解验证实现相对残差比绝对误差更可靠在上一段代码中直接比较x与x_ref的绝对差值会受解的量级影响。更好的指标是相对残差r np.linalg.norm(A x - b) / (np.linalg.norm(A) * np.linalg.norm(x)) print(relative residual , r)这个值描述的是“误差相对于问题规模的占比”不会因为 b 变成 1000 倍而变化太多。实际工程里如果相对残差在 1e-10 量级通常可以认为迭代收敛到满意状态。若相对残差卡在 1e-6 下不去首先查条件数其次查容差设置最后再去怀疑实现。这个判断顺序也是“数值分析及其MATLAB实验”里隐含考察的内容。6. 最后一道“应用题”用残差与收敛阶验证你的数值实现无论试题属于哪一类我最终都会把它收敛到一个统一动作构造一个已知参考解的用例跑几步看误差是否按理论阶数下降。这个方法可以筛掉绝大多数下标记错、系数抄错、切片越界的低级 bug。具体做法是每次写一个新算法先不着急处理真实数据而是用解析解非常明确的小问题做冒烟测试。例如插值就选y x^3取 4 个节点插值结果必须在这些节点上精确等于原值。最小二乘就选一条直线加一点噪声比较系数是否在噪声范围内接近真值。数值积分就选sin(x)或exp(x)。ODE 就选dy/dt y。这些参考解都是临时的成本几乎为零但能立刻暴露循环边界是否多走了一位。真正值得投入精力的地方是收敛阶测试。以复合 Simpson 为例令 n 从 8 增加到 16误差应从某个值下降到约 1/16如果只下降到 1/4说明实现的不是 Simpson 而是某种二阶方法通常是因为4*ys[1:-1:2]和2*ys[2:-2:2]这两个切片的权重被写反了。对 RK4 同理n 翻倍误差下降约 16 倍如果只有 2 倍多半是 k1、k2、k3、k4 的系数写错了。用线性方程组的迭代法时验证方式则是固定同一初值把tol从 1e-6 改成 1e-10观察迭代次数是否按预期增加同时相对残差是否有明显下降。最后还有一个我自己的约定任何算法代码都必须保留一个--check模式运行小型参考用例并自动比较误差。这个习惯直接来自复习“数值分析试题及答案”时的经验答案可以告诉我结果但只有误差行为能告诉我为什么结果可信。如果你在做题或实验时遇到结果与答案不一致先不要改参数停下来跑一遍收敛阶测试通常下一步就清楚了。本文还有配套的精品资源点击获取