严格耦合波分析与平面波展开法的MATLAB实现：从光栅衍射到光子晶体能带

简介：针对周期性介质结构中的光波传播与衍射问题，这份MATLAB实现包提供了严格耦合波分析（RCWA）与平面波展开法的完整代码框架，适用于光学工程、光栅设计、光子晶体等方向的研究者与学生。包内共24个文件，以13个m脚本为主体，涵盖RCWA主程序、网格生成、PML吸收边界、反射/透射计算等模块；另有8个dat数据文件、1个txt说明、1个asv备份及1个md文档，便于快速理解代码结构与运行流程。压缩包仅17KB，小巧精炼。目前已有85人学习下载，是快速上手两种数值方法的实用参考。通过源码可直观掌握傅里叶级数展开、耦合波方程组构建及平面波叠加的实现细节，并可直接改造用于自定义周期性结构分析，显著缩短理论到仿真的转化成本。

1. 严格耦合波分析和平面波展开法的MATLAB实现：先从两条倒空间思路说起

严格耦合波分析（RCWA）和平面波展开法（PWE）的MATLAB实现，是我在处理周期性光栅和光子晶体时最先选择自己动手写的一对频域工具。它们不依赖商业电磁仿真包，也不需要调FDTD网格，只要会用矩阵分解，就能让衍射效率和能带图从你自己的数据里长出来。RCWA 擅长回答“一个周期光栅把入射光分到每个衍射级几次”，PWE 擅长回答“这族周期结构在哪些频率存在传播模式”。适合刚接触光子模拟的人建立坐标系，也适合熟手在收敛性问题上快速拆解算法。

很多人以为这两个方法已经很老，但在超构表面、硅光栅耦合器和光子晶体波导的设计里，标准做法仍然是它们。MATLAB 里实现两者的最大价值不是性能，而是把截断阶数、介电常数傅里叶系数这些容易当成“黑匣子”的参数全部暴露在眼前。后面各章会从原理一直走到可运行的代码和参数边界，你可以把模型从单层光栅一路扩展到二维光子晶体。

2. 先把倒格子语言备好：严格耦合波分析和平面波展开法共同的数学根基

2.1 光栅方程和布洛赫条件：两种方法都在倒空间里截断场

周期性介质里的电磁场不能再写成单一平面波。RCWA 处理的问题是：入射波打在周期为 d 的界面后，散射场的横向波数被周期性调制，产生一系列衍射级，满足光栅方程 (k_{x,m}=k_0 \sin\theta - 2\pi m/d)。这个关系决定了每一级的横向波数，也决定了哪些级次能真实传播、哪些是倏逝波。

PWE 处理的则是无限周期介质中的本征模式，场写成布洛赫展开 (\mathbf{H}(\mathbf{r})=\sum_{\mathbf{G}}\mathbf{H}_{\mathbf{G}} e^{i(\mathbf{k}+\mathbf{G})\cdot\mathbf{r}})，其中 (\mathbf{G}) 是倒格矢，(\mathbf{k}) 被限制在第一布里渊区。两种方法的共同点在于：都把连续空间中的电磁场问题转换成倒空间中的矩阵问题。矩阵的大小完全取决于你保留多少个衍射级或倒格矢，这就是截断阶数。

从代码框架看，RCWA 和 PWE 的步骤高度相似：先构建周期结构介电常数的傅里叶系数，再用这些系数组装某个矩阵，最后做特征值分解或求解线性系统。理解了这一层，后面所有参数调整都有了方向。

2.2 方法选型：算衍射效率用RCWA，算能带用PWE，别用反

选型最大的错误是拿 PWE 去推衍射效率。PWE 给你的只是本征频率，没有外部入射波的耦合强度，算不出“百分之几的光透过了第 m 级”。反过来拿 RCWA 去算无限周期光子晶体的能带也存在问题，你需要引入复频率扫描或者加周期边界，代价明显高于 PWE。我一般的判断标准是：问题里有没有对结构的“外部入射”？有，选 RCWA；没有，只是问模式在哪些频率存在，选 PWE。

2.3 截断阶数与收敛性：决定你矩阵规模的核心参数

RCWA 的一维截断阶数 (N) 决定矩阵宽度为 (2N+1)。PWE 在二维方格子中取倒格矢截断半径 (N_p) 时，基函数数是 ((2N_p+1)^2)。它们的共同规律是：截断越大，矩阵越稠密，计算越慢，结果越逼近真解；但到某个值之后，再增大只是把小数点后第三位推稳而已。

收敛性的通用判据是做一组截断扫描。比如 RCWA 中固定波长和角度，比较 (N=8,12,16,20) 时的零级反射效率，变化小于 0.1% 就可以认为稳定。PWE 中固定高对称点，看第一布里渊区边界频率随 (N_p) 的变化，通常 (N_p=6\sim10) 对大多数介质柱结构足够。注意高折射率对比和尖锐边缘会让傅里叶系数衰减变慢，需要更高截断。TM 偏振比 TE 收敛慢，这是由介电常数傅里叶展开的性质决定的，后面避坑章会展开。

3. 严格耦合波分析的MATLAB实现：一维光栅的最小可运行框架

3.1 定义光栅参数和入射条件：波长、周期、占空比和偏振

一个能立即跑起来的 RCWA 模型，参数越少越好。下面这套定义对应一维矩形光栅，只在 x 方向周期变化，入射面是 x-z 平面，层内 y 方向均匀。这也是严格耦合波分析最经典的教学场景。

lambda = 1.55e-6; % 入射自由空间波长，单位 m d = 0.8e-6; % 光栅周期，单位 m f = 0.5; % 占空比：一个周期内高介电材料占比 eps_g = 12.25; % 光栅材料的相对介电常数，比如硅 eps_b = 1.0; % 背景材料，空气 theta = 30 * pi / 180; % 入射角，相对表面法线 h = 0.4e-6; % 光栅层厚度，单位 m N = 12; % 保留 -N..N，共 2*N+1 个衍射级 m = (-N:N).'; % 级次索引列向量

这里的关键参数是周期和波长的比值。当 (d/\lambda) 较大时，可传播的衍射级数多，矩阵即使不涨，也需要更多级次来描述能量分配。占空比 (f) 直接影响介电常数傅里叶系数，定义为高介电材料在一个周期内占的比例。如果你模拟的是槽结构，把背景和光栅材料对调即可。厚度 (h) 在这个单层模型里决定层内相位积累，深沟槽光栅的 RCWA 对层内数值稳定性更敏感，后面会说明。(N) 的初始值取 12 是经验做法，先跑通再扫描收敛。

3.2 组装介电常数Toeplitz矩阵与状态矩阵A

RCWA 的核心是把空间周期的介电常数 (\varepsilon(x)) 展开成傅里叶级数。对于矩形占空比结构，这个级数有解析形式：零阶项是平均介电常数，高阶项正比于 (\mathrm{sinc}) 函数。

deltaE = eps_g - eps_b; eps_m = eps_b + deltaE * f * sinc(m * f); % 介电常数傅里叶谐波 E = toeplitz(eps_m); % 卷积运算对应的 Toeplitz 矩阵 kx = sin(theta) - m * lambda / d; % 归一化横向波数 Kx = diag(kx); A_TE = Kx^2 - E; % TE 偏振状态矩阵

代码中sinc(m*f)使用的是 MATLAB 的归一化 sinc 函数，(\mathrm{sinc}(x)=\sin(\pi x)/(\pi x))。当 (m=0) 时它等于 1，所以eps_m的第一个元素就是 ( \varepsilon_b + \Delta\varepsilon f)，对应空间平均介电常数。toeplitz(eps_m)构造的矩阵第 (i,j) 个元素是eps_m(i-j)，恰好把傅里叶级数在倒空间的卷积变成了矩阵乘法。这是整个 RCWA 最容易写错的地方，很多翻车都是从傅里叶系数正负号或索引错位开始。

横向波数kx直接来自光栅方程。注意这里用了 (\sin\theta - m\lambda/d) 的符号约定，它决定了倒空间索引的方向；无论用正号还是负号，最终衍射效率相同，但和后续边界匹配的相位定义必须保持一致。A_TE = Kx^2 - E来自 TE 偏振的波动方程：电场只有 y 分量，对 z 的二阶导数由横向波数平方和介电常数共同决定。TM 偏振需要换成逆介电常数矩阵，见避坑章。

3.3 求解层内模式：特征值分解与传播常数分支选择

状态矩阵A_TE的特征值就是层内归一化传播常数的平方。求解特征值后，还要处理一个物理分支问题：每个特征值开平方得到两个符号，必须选择沿正 z 方向传播或衰减的那一支。

[U, D] = eig(A_TE); q2 = diag(D); q = sqrt(q2); flip = (real(q) < 0) | (abs(real(q)) < 1e-12 & imag(q) < 0); q(flip) = -q(flip); q = q(:); % 确保列向量

这段代码看起来简单，但它是最容易引入玄学错误的地方。sqrt对负实数默认返回正虚部，对复数默认返回实部非负的主分支，但不保证满足沿正 z 方向的辐射条件。我一般把传播模式要求为实部大于 0，倏逝模式要求为虚部大于 0，再统一翻转符号。特征向量U不需要翻转，因为 (q) 和 (-q) 对应同一个特征向量列。

特征向量矩阵U的每一列代表层内一个“耦合波模式”，它的行对应不同衍射级。后面的边界匹配就是把这些模式按边界条件叠加，而不是直接对每个级次单独求解。

3.4 边界匹配与衍射效率：用分块线性系统解出反射和透射系数

单层光栅的边界条件可以写成 4 组未知量：反射系数 (R)、透射系数 (T)、层内正向模式幅度 (c_+)、反向模式幅度 (c_-)。把电场和磁场的切向连续性组装成一个大分块矩阵，一次\求解就能得到所有系数。

k0 = 2*pi/lambda; kIz = sqrt(eps_b*k0^2 - (kx*k0).^2); % 入射介质 z 向波数 kIIz = sqrt(eps_b*k0^2 - (kx*k0).^2); % 出射介质，这里是悬浮膜对称结构 n = 2*N+1; I0 = zeros(n,1); I0(N+1) = 1; % 零级入射振幅为 1 P = exp(1i*k0*q*h); % 正向传播相位 Pm = exp(-1i*k0*q*h); % 反向传播相位 M = [eye(n), zeros(n), -U, -U; -diag(kIz), zeros(n), -k0*U*diag(q), k0*U*diag(q); zeros(n), -eye(n), U*diag(P), U*diag(Pm); zeros(n), -diag(kIIz), k0*U*diag(q)*diag(P), -k0*U*diag(q)*diag(Pm)]; rhs = [-I0; -diag(kIz)*I0; zeros(2*n,1)]; sol = M \ rhs; R = sol(1:n); % 反射系数，对应 m = -N..N T = sol(n+1:2*n); % 透射系数

组装逻辑来自四条关系式：在 (z=0) 面电场和磁场切向连续，在 (z=h) 面电场和磁场切向连续。代码第一行对应电场连续，第二行对应磁场连续，第三、四行是光栅层内部的场传播关系。矩阵中-k0*U*diag(q)前面的k0来自对 z 方向求导后归一化，不能丢掉；丢掉它会得到看起来很合理但能量不守恒的结果。

衍射效率的计算需要按能流归一化。电场系数平方只是振幅比，真实衍射效率要乘上对应级次在 z 方向的波数实部比：

kIz0 = kIz(N+1); ref_eff = abs(R).^2 .* real(kIz) ./ real(kIz0); tra_eff = abs(T).^2 .* real(kIIz) ./ real(kIz0); sum(ref_eff) + sum(tra_eff) % 悬浮膜无吸收时应该等于 1

这个归一化对反射和透射都成立，因为能流的 z 分量正比于 (\mathrm{Re}(k_z)|E|^2)。上下介质不同时，只要分别用对应的入射和出射波数，公式仍然成立。如果总和不等于 1，先别怀疑算法，按第 5 章的分支选择去查代码。

提示：这个线性系统适合单层和较厚的均匀层。层数很多或者有深度光栅时，直接拼接矩阵会出现数值病态，工程上应该改用散射矩阵（S matrix）递推。

4. 平面波展开法的MATLAB实现：二维光子晶体能带的完整流程

4.1 从晶格常数到倒格矢：生成波矢基组

平面波展开法要先把实空间晶格翻译成倒空间基组。以二维正方晶格为例，原胞为边长 (a) 的正方形，倒格矢也是正方形，坐标间隔 (2\pi/a)。先生成波矢基组，后面的矩阵组装都在这个基组里完成。

a = 1; % 晶格常数，后面归一化会用到 R = 0.2*a; % 圆孔半径 eps_b = 12.0; % 背景介质介电常数 eps_a = 1.0; % 孔内介质，空气 Np = 8; % 倒格矢截断半径，以 2*pi/a 为单位 f = pi*R^2/a^2; % 孔占空比 g = 2*pi/a; [gx, gy] = meshgrid((-Np:Np)*g, (-Np:Np)*g); G = [gx(:).'; gy(:).']; % 2 x nG 的倒格矢矩阵

这里 (N_p=8) 表示每维保留 (2N_p+1=17) 个倒格矢，总基函数数为 289。对于圆孔光子晶体，这个数值已经能把前几条带算到工程可用精度。如果你做的是三角晶格或复杂原胞，网格生成方式要换成晶格对应的倒格矢公式，但后续矩阵组装逻辑完全一致。G的列顺序由meshgrid决定，第一列是倒空间原点，这正好对应后面的 Toeplitz 索引。

4.2 计算结构因子和傅里叶系数：圆孔的解析形式

圆孔的傅里叶系数有解析解，不需要用实空间网格做 FFT。关键在于把“孔内为空气、背景为高介电”的介电函数写成背景加修正项。修正项的傅里叶系数就是圆孔结构因子。

gn = sqrt(sum(G.^2, 1)); S = zeros(size(gn)); S(gn < 1e-12) = f; % 零倒格矢处的结构因子等于填充率 nz = gn >= 1e-12; S(nz) = 2 * f * besselj(1, gn(nz)*R) ./ (gn(nz)*R); epsC = eps_b * double(gn < 1e-12) - (eps_b - eps_a) * S; invC = (1/eps_b) * double(gn < 1e-12) - (1/eps_b - 1/eps_a) * S;

代码里besselj(1, x)是一阶贝塞尔函数，这个公式来自圆盘在倒空间的填充因子，(S(\mathbf{G})=f\frac{2J_1(GR)}{GR})。当 (G=0) 时结构因子等于填充率 (f)，于是epsC(0)正好是空间平均介电常数 (\varepsilon_a f + \varepsilon_b(1-f))。invC不是简单地对epsC取倒数，而是对倒数介电函数 (1/\varepsilon(\mathbf{r})) 做同样的展开。这一步对 PWE 的 TM 模式是必须的，跳过它就会出现本征频率偏移。

4.3 组装本征方程并扫描高对称路径

TM 极化（电场沿 z 方向）的本征方程是 [ \sum_{\mathbf{G}'}\varepsilon^{-1}(\mathbf{G}-\mathbf{G}')(\mathbf{k}+\mathbf{G})\cdot(\mathbf{k}+\mathbf{G}')E_{\mathbf{G}'} = \frac{\omega^2}{c^2}E_{\mathbf{G}}. ] 写成矩阵形式时，介电常数卷积矩阵和点积矩阵按元素相乘，而不是普通矩阵乘法。

Etoep = toeplitz(epsC); invEtoep = toeplitz(invC); function [freq, evecs] = pwe_eig(kx, ky, G, invEtoep) Kgx = kx + G(1,:).'; Kgy = ky + G(2,:).'; Amat = Kgx * Kgx.' + Kgy * Kgy.'; % (k+G_i) 点乘 (k+G_j) H = invEtoep .* Amat; H = 0.5 * (H + H'); % 显式埃尔米特化，压制数值噪声 [evecs, evals] = eig(H); freq = sqrt(max(real(diag(evals)), 0)); freq = sort(freq); freq = freq(freq > 1e-8); % 丢弃数值零模 end

函数中invEtoep .* Amat是关键：Amat(i,j)是两个波矢点乘，invEtoep(i,j)是倒数介电常数的傅里叶系数，两者逐元素相乘才等价于对 (\mathbf{G}') 求和。显式埃尔米特化对中心对称结构没有物理副作用，却能把eig返回的极小虚部清掉。频率计算取sqrt后要过滤接近 0 的模式，因为平面波基底存在零频率平凡解，它们不代表真实传播模式。

4.4 画出归一化能带：频率排序与实部提取

能带图需要沿着布里渊区边界的高对称路径扫描波矢 (\mathbf{k})。正方晶格常用路径是 (\Gamma \to X \to M \to \Gamma)，坐标以 (2\pi/a) 为单位分别是 ((0,0)\to(0.5,0)\to(0.5,0.5)\to(0,0))。

npts = 30; path1 = [linspace(0, 0.5, npts); zeros(1, npts)]; path2 = [0.5*ones(1, npts); linspace(0, 0.5, npts)]; path3 = [linspace(0.5, 0, npts); linspace(0.5, 0, npts)]; klist = [path1, path2, path3]; xaxis = [linspace(0,1,npts), linspace(1,2,npts), linspace(2,3,npts)]; Nb = 6; bands = zeros(Nb, length(klist)); for i = 1:length(klist) freq = pwe_eig(klist(1,i), klist(2,i), G, invEtoep); bands(:,i) = freq(1:Nb) * a / (2*pi); % 归一化频率 a/λ end plot(xaxis, bands, 'LineWidth', 1.2); set(gca, 'XTick', [0 1 2 3], 'XTickLabel', {'Gamma','X','M','Gamma'});

归一化频率是 (\omega a/(2\pi c))，也就是 (a/\lambda)。扫描时每个 (\mathbf{k}) 点都重新组装Amat，但invEtoep只依赖结构，可以在循环外算好。频率排序后直接取前 (N_b) 条带，就能看到带隙位置。如果某些带在带隙边缘出现小突起，优先怀疑倒格矢截断不够，而不是物理现象。

5. 严格耦合波分析和平面波展开法会遇到的5个坑：从现象到修复

5.1 现象：衍射效率总和不等于 1

单层悬浮膜、无吸收介质的情况下，反射加透射效率应该严格等于 1。我第一次跑通时，总效率只有 0.82，看上去像是能量凭空消失，实际原因在传播常数的分支选择。

原因主要有两个：一是对q取sqrt后没有统一翻转符号，导致层内相位因子exp(1i*k0*q*h)指向了错误的衰减方向；二是kIz和kIIz对倏逝级次取了错误的虚部符号，让边界匹配方程把能量“吸走”了。

解决方式是先做自检：取一个没有光栅或占空比趋近 0/1 的极限，这时多层反射率应该趋近无光栅薄膜的解析解。再对零吸收结构打印效率总和，如果在 (1\pm10^{-6}) 内，分支选择基本没问题。如果差值很大，逐项检查q的实部和虚部符号，特别要关注real(q)<0的翻转逻辑。

5.2 现象：TM 偏振收敛很慢，甚至振荡

TE 偏振用介电常数的 Toeplitz 矩阵E，收敛还算正常。TM 偏振如果照搬同样流程，衍射效率会随截断阶数振荡，很多级次的能量跳来跳去，看起来像数值噪声。

原因是 TM 极化的波动方程里出现的是 (1/\varepsilon) 算符，而 (1/\varepsilon) 的傅里叶系数不等于介电常数傅里叶系数的倒数。直接对每阶系数取倒数，等于改变了材料的空间分布，当然不收敛。

解决方式是在组装状态矩阵时使用逆介电常数的 Toeplitz 矩阵，而不是1/E的元素级逆。也就是说，先对 (1/\varepsilon(\mathbf{r})) 做傅里叶展开，再构造 Toeplitz。这样 TM 的收敛速度会明显改善。更严格的还是遵循傅里叶变换的卷积规则，不同极化方式处理方式不同：TE 用 (\varepsilon) 的 Toeplitz，TM 用 (\varepsilon^{-1}) 的 Toeplitz。

5.3 现象：PWE 能带图出现虚假平带

用 PWE 算能带时，本征值谱里经常出现一条完全平坦的带，频率不随波矢变化。它不是物理上的缺陷态，而是矩阵组装时引入了额外的零空间。

原因通常是逆介电常数矩阵的零倒格矢项和结构因子处理不一致。比如圆孔结构里epsC(0)算错，或者toeplitz(invC)的顺序没有和G的列顺序对齐，导致卷积矩阵不对称，平面波基底不能正确解耦。

解决方式是先把倒空间原点确认清楚：G的第一列必须是 ((0,0))，epsC的第一个元素必须对应原点。然后再看本征值取实部时有没有被sqrt引入负值。过滤频率时，阈值 (10^{-8}) 足够去掉伪零模，但如果H没有埃尔米特化，会留下幅值小于 (10^{-4}) 的虚部分量，推荐代码里显式做0.5*(H+H')。

5.4 现象：RCWA 薄层线性系统病态

单层线性系统在光栅层很薄时能跑通，但层数变多或层厚差异巨大时，M \ rhs开始不稳定，效率总和剧烈波动。我的经验是不要试图扩展这个线性系统到多层，那不是小修小补，而是换算法。

原因是正向相位因子P和反向相位因子Pm在薄层里值非常接近，边界匹配矩阵的列几近线性相关，条件数爆炸。散射矩阵方法把每一层看成独立的散射单元，通过递归合并避免直接求大线性系统，数值稳定性好很多。

常见做法是保留单层代码作为验证工具，但正式的多层光栅全部改用 S 矩阵更新。如果只是自己研究单层光栅，线性系统就够了；如果要搭滤波片或层叠光栅，尽早切换到 S 矩阵，不要在薄层病态问题上花太多时间。

5.5 现象：二维 RCWA 结果和三维仿真对不上

严格耦合波分析这个名字容易让人误以为它是全矢量三维方法，但前面给出的实现只适用于“结构在一维方向周期、另一个横向方向均匀”的问题，也就是入射面内衍射。如果光栅在 x 和 y 两个方向都有周期，或者入射波倾斜到包含 y 分量，衍射级次需要用二维倒格子索引，代码里的横向波数会变成两个分量。

原因很简单：二维光栅的 TE/TM 偏振不再能用一个标量状态矩阵描述，必须写出完整的 3×3 电磁场分量耦合。很多人在这个问题上翻车，是因为拿一维代码去算二维光栅，结果效率复数或者角度不对。

解决方式是检查问题维度：结构只在 x 方向周期、入射面在 x-z 平面，才适用上面的标量公式；否则需要升级到三维 RCWA，矩阵规模随两个方向的截断阶数乘积增长，参数调整原则仍然一致，但代码复杂度不可同日而语。

6. 从能跑通到可发表：验证方法、一个加速技巧和我的自检习惯

6.1 用光栅方程和薄膜极限做交叉验证

我每次写完新的 RCWA 或 PWE 代码，第一件事不是看目标曲线，而是先跑两个廉价验证。RCWA 验证入射角和衍射级次位置：固定波长和周期，把零级透射峰对应的衍射角与解析光栅方程比对。PWE 验证极限结构：把圆孔半径调到很小，能带应该逼近均匀介质的色散线 (\omega = c|\mathbf{k}|/\sqrt{\varepsilon_b})。

这两个验证都不需要额外工具，代码里打印几行数值就能完成。如果极限对不上，不要继续往下调参数，问题一定在傅里叶系数或矩阵组装中。我见过太多同学在复杂结构上调整半天，最后发现是 sinc 函数多了一个 π。

6.2 加速参数扫描：复用 Toeplitz 矩阵

做波长扫描或入射角扫描时，RCWA 的结构矩阵E只依赖光栅几何，不依赖入射波长和角度。PWE 的invEtoep只依赖介质分布，不依赖波矢。把这两项从循环里提出来，可以省掉大量重复构造时间。

更深一层的技巧是在 PWE 中对同一批倒格矢，先算好 (G_x,G_y) 的坐标矩阵，每个 k 点只需要做一次实对称矩阵特征分解。对 MEMS 光栅或温调滤波器这类参数扫描场景，这种复用能把单次扫描从几分钟压到几十秒。特征分解这一步无法避免，但至少不要重复计算不变的部分。

6.3 我的自检习惯和最后的建议

我现在的习惯是在每个脚本开头写一个小节“自检开关”：一个极限结构测试、一个功率守恒断言。RCWA 里断言abs(sum(ref_eff)+sum(tra_eff)-1) < 1e-6，PWE 里断言最低带频率在均匀介质极限下与解析值误差小于 1%。有了这两个断言，后面改动参数时，数值上稍微不干净就会立刻报警。

严格耦合波分析和平面波展开法其实是一对很好的“数值显微镜”，它们让你看到截断阶数如何在倒空间里决定精度。先跑通单层光栅的 TE 模型，再扩展 TM；先算最基础的圆孔正方晶格能带，再加超胞。每个扩展都在已有的矩阵框架上加一个维度，而不是推翻重来。希望这篇实践笔记能帮你少走几个我走过的弯路，也希望你在跑通第一条能带曲线时，能感受到把矩阵变成物理图像的那种踏实感。

本文还有配套的精品资源，点击获取