企业官网建设流程全解析

虚拟器官插件开发教程（7）：从细胞到组织——传导、折返与 svFSI 的 EP 算例

版本声明块

• 工具/软件：svFSI / svMultiPhysics（BSD 系许可，EP 内建）；openCARP v19.0（2026-04，本体 Academic Public License v1.1）；Chaste 2026.1（BSD-3-Clause）
• 语言/环境：Python 3.11 + numpy（行波演示）；Fortran/C++ 求解器（官方算例）
• 本文目标：把第 6 篇的"单细胞 qNet"推进到组织尺度——理解传导方程、摸清 svFSI 内建电生理的能与不能、建立波长-折返的基质判据

一句话结论：单域方程dV/dt = D·∂²V/∂x² − I_ion + I_stim的最小 numpy 实现即可复现传导行波，实测 D 每翻倍 CV 约乘 √2（本文实跑 71.5→106.7→154.3→219.9 cm/s）；svFSI/svMultiPhysics 内建cepModel_AP/FN/BO/TTP电生理注册表（C++ 类CepModTtp含I_Na、I_bNa），官方svFSI-Tests的08-cep/03-benchmark_tTP对应 Niederer 2011 N-version benchmark（Phil Trans R Soc A 369:4331，3×7×20 mm 组织块）；openCARP 本体是非商业Academic PL v1.1（商用找 NumeriCor），而"Chaste 是 GPL"属于过时误传（2026.1 为 BSD-3）。

〇、本篇要解决的认知问题

1. 单域（monodomain）与双域（bidomain）模型到底差在哪，插件开发者该选哪个？
2. SimVascular 的 svFSI/svMultiPhysics 究竟有没有心脏电生理？内建哪些模型、官方给了哪些算例？
3. openCARP、carputils、meshalyzer、Chaste 各自的许可证是什么？商用插件哪些能直接链接？
4. Niederer 2011 N-version benchmark 对电生理插件验证有什么用？
5. 动作电位波速怎么变成"传导速度→激动波长→折返基质"的心律失常风险语言？

一、机制解析

1.1 从细胞到组织：只差一项扩散，语义全变

第 6 篇里一个细胞的状态是 ODE 系统；铺成组织后膜电位多了空间耦合，成为反应-扩散（reaction-diffusion）方程：

• 单域模型（monodomain）：只解一个跨膜电位场 V，∂V/∂t = D∇²V − I_ion/C_m + I_stim。各向异性写进张量 D。计算量小，适合插件迭代与批量筛选。
• 双域模型（bidomain）：胞内、胞外两个电位场 coupled，能给出胞外电位/伪 ECG 与刺激阈值场。文献与代码中统称"单域/双域"，中文首次出现给英文原文即 monodomain/bidomain。

两者的工程分水岭：要不要细胞外空间。做传导速度与折返机制，单域够用；要往第 8 篇的伪 ECG 或临床电极对比走，才需要双域（或单域+离线积分的近似）。

尺度链（插件视角）： 通道阻滞(第5篇) → 单细胞AP/qNet(第6篇) → 组织传导/波长(本篇) → 器官ECG(第8篇) Hill/Markov ODE 积分 反应-扩散 PDE 容积导体叠加

1.2 波速为什么是风险语言

行波速度即传导速度（conduction velocity，CV）。两个直接推论：

• 激动波长 λ = CV × ERP（有效不应期）：一个波长内组织来不及恢复兴奋。
• 折返基质：折返环周长 ≥ λ 时，波前追上的组织已脱不应期，折返可持续——这是室速/室颤维持的经典判据；电生理标测（EP study）中的程序刺激（EPS）本质就是在测不应期与传导，判断有没有这种"可被点燃的环路"。

药物插件在这里的输出语义就变了：第 6 篇给"APD90/qNet 延长了"，本篇要回答"传导慢了没、不应期变了没、λ 小到能让哪块组织折返"。智源"虚拟生理心脏"官方叙事中的"组织传导速度/波长/兴奋易损区/折返波"层级（见第 8 篇讨论其器官层）说的就是这一段。

1.3 工具箱：三条组织级仿真路线与许可证分层（铁律 3）

openCARP 论文口径：Comput Methods Programs Biomed 2021;208:106223（doi:10.1016/j.cmpb.2021.106223）。它继承 CARP/acCELLerate 脉络，组织-器官级功能最全——但立项前先过许可证矩阵：把 openCARP 当"BSD 同名软件"商用是本系列见过的高危误判。一个诚实的替代姿势：把 openCARP 用于研发期对照，商用交付切到 svMultiPhysics/Chaste 路线，插件层用统一契约隔离后端。

1.4 svFSI/svMultiPhysics 内建 EP：有什么、缺什么

官方口径先纠正一个流传很广的错误：“SimVascular 只管血流”——不成立。svFSI 简介原文含 “blood flow simulation including fluid-structure interaction and cardiac electrophysiology”；svMultiPhysics 覆盖 “solid and fluid mechanics, diffusion, and electrophysiology … whole heart dynamics”。

内建电生理模型注册表（svFSI README Features 表，Doxygen 可查实现类）：

TP06 出处 Am J Physiol Heart Circ Physiol 291(5):H2396–H2411，doi:10.1152/ajpheart.00109.2006。注意CepModTtp已经显式拆出I_Na/I_bNa（晚钠）——这给第 9 节"药物阻滞层自建"留了接口：把浓度→阻滞率映射接到电流项缩放，就是插件要写的代码（铁律 4：官方没有 hERG 阻滞→I_Kr 缩放的药理学模块，GUI 也没有 EP 面板，EP 全部在求解器输入文件层配置）。

官方算例库svFSI-Tests（EP 相关部分）：

svFSI 的 JOSS 论文（2022）明确口径 “simulating the complex excitation-contraction coupling…”。N-version benchmark 的价值：多个独立团队在同一几何/同一模型上的互测基准——你的 EP 求解器接进插件后，先复跑它对齐官方数值，再谈药物场景（第 15 篇验证方法学会展开）。浦肯野刺激通路用 SimTK 官方分发的Purkinje Plugin（跨平台安装包 2022-07-20，“used to create a Purkinje network on a surface model of the heart”）——SimVascular 体系里唯一官方电生理刺激插件；患者特异全管线看 SimCardio（docsSimCardio.html）。

二、完整代码与逐行剖析

2.1 工程动作：定位官方算例并跑通（bash）

# 获取算例库与求解器（命令形态按各仓库 README，以官方文档为准）gitclone https://github.com/SimVascular/svFSI-Tests# 官方算例库gitclone https://github.com/SimVascular/svFSI&&cdsvFSI# svFSI 为 Fortran 多物理求解器；官方调用形态（调研核实）：mpiexec-np<N><build>/svFSI-build/bin/svFSI<input>.inp# 输入 .inp，模板为 svFSI_master.inpcd../svFSI-Tests/08-cep/03-benchmark_tTP# Niederer 2011 N-version 基准# 该算例几何为 3×7×20 mm 组织块——正是 1.2 节"环路周长 vs 波长"判据能直接套用的尺寸

要点：svFSI 的 EP 场景全部在.inp输入文件层配置，GUI 不面板化；svMultiPhysics（现行默认后端）走 XML 输入，且无原生 Windows，需 WSL。

2.2 svFSI.inp 关键段解读（示意结构，字段以官方仓库为准）

# ---- svFSI_master.inp 电生理相关关键段（教学示意，逐段语义解读） ---- [EP] cepModel cepModel_TTP ; 内建注册表名：AP/FN/BO/TTP 四选一 ; 选 TTP = ten Tusscher-Panfilov 全离子模型（CepModTtp） ; 药物插件的挂载点：把 Hill 阻滞率接到 I_Na/I_bNa/IKr 电流缩放 ;（官方无药理学层，这段"接到"的代码就是第 17 篇 cardiotox 插件本体） cellType 2 ; 透壁细胞类型 endo/M/epi 选档 → 第 8 篇伪ECG的梯度来源 stimNodesFile stim.dat ; 起搏点文件；浦肯野方案改用 Purkinje Plugin 生成的网络文件 [Diffusion] D tensor 1 ; 各向异性扩散张量：需要纤维方向数据（无纤维=各向同性退化） ; 纤维方向缺失时 CV 会被系统性高估——传导结论先过此检查

读这段的正确姿势：不求逐字背诵，记住三件事——模型名在注册表里选、起搏靠输入文件、各向异性靠纤维数据。

2.3 1D 单域行波：CV 标度 + ERP + 波长（Python，已实跑）

# -*- coding: utf-8 -*-"""第7篇实跑验证：1D 单域(monodomain)反应-扩散行波 + CV 标度 + ERP + 波长/折返基质 反应项 = 第6篇同款"最小平台期细胞"（玩具演示模型，非 TP04/TP06/ORd 本体； 真实传导速度以 svFSI-Tests 08-cep/03-benchmark_tTP（Niederer 2011 3×7×20mm）等官方算例为准）。 """importnumpyasnp E_NA,E_K,E_CA=60.0,-85.0,60.0defxinf(V,Vh,k):return1.0/(1.0+np.exp(np.clip(-(V-Vh)/k,-50,50)))# 激活门defhinf(V,Vh,k):return1.0/(1.0+np.exp(np.clip((V-Vh)/k,-50,50)))# 失活门defionic(V,h,s,f,w,gNa=4.0,gKr=0.06,gCa=0.05,gNal=0.002,gK1=0.09):m=xinf(V,-40.,4.)f1=1.0/(1.0+np.exp(np.clip((V+40.)/10.,-50,50)))# IK1 内向整流因子return(gNa*m*h*(V-E_NA)+gCa*s*f*(V-E_CA)+gNal*(V-E_NA)+gKr*w*(V-E_K)+gK1*f1*(V-E_K))# 合电流（外向为正）defgate_step(V,h,s,f,w,dt,tau_h=2.,tau_s=30.,tau_f=150.,tau_w=80.):h+=dt*((hinf(V,-60.,3.)-h)/tau_h)s+=dt*((xinf(V,-25.,7.)-s)/tau_s)f+=dt*((hinf(V,-45.,5.)-f)/tau_f)w+=dt*((xinf(V,-15.,8.)-w)/tau_w)returnh,s,f,wdefmonodomain_1d(D=0.002,L=8.0,dx=0.01,T=120.0,dt=None):"""dV/dt = D·∂²V/∂x² − I_ion + I_stim；显式 Euler。 dt 自动满足扩散稳定条件 dt ≤ 0.4·dx²/D——单域方程的第一号数值坑。"""ifdtisNone:dt=min(0.05,0.4*dx*dx/D)nx=int(L/dx);nt=int(T/dt)V=np.full(nx,-85.);h=hinf(V,-60.,3.);s=xinf(V,-25.,7.)f=hinf(V,-45.,5.);w=xinf(V,-15.,8.)t_up=np.full(nx,np.nan);x=np.arange(nx)*dxforitinrange(nt):t=it*dt lap=np.zeros(nx)lap[1:-1]=(V[2:]-2*V[1:-1]+V[:-2])/dx**2# 中心差分；两端保持 0=Neumann 零流I=np.where(x<dx*2,40.,0.)ift<1.0elsenp.zeros(nx)# 左端 1ms 局域刺激V+=dt*(D*lap-ionic(V,h,s,f,w)+I)# 扩散项+反应项：行波的全部来源h,s,f,w=gate_step(V,h,s,f,w,dt)fresh=np.isnan(t_up)&(V>-10.)# 阈值上时：记录各点激活时刻t_up[fresh]=t fin=np.isfinite(t_up);xs,ts=x[fin],t_up[fin]keep=xs>dx*20# 去掉刺激附近非线性启动段p=np.polyfit(xs[keep],ts[keep],1)# t_up≈x/CV+b → 斜率倒数=CVresid=np.max(np.abs(ts[keep]-np.polyval(p,xs[keep])))return1.0/p[0],resid,int(keep.sum()),dtprint("== CV 随扩散系数（=细胞间耦合强度/纤维化程度代理）的标度 ==")prev=NoneforDin[0.001,0.002,0.004,0.008]:cv,resid,n,dt=monodomain_1d(D=D)line=f" D={D:g}dt={dt:.3g}→ CV={cv*1000:6.1f}cm/s (拟合残差{resid:.2f}ms, n={n})"ifprev:line+=f" CV 比={cv/prev:.2f}vs √2={np.sqrt(2):.2f}"prev=cv;print(line)deferp_cardiac_like(dt=0.01,SI_list=None):"""向量化 S1-S2：所有间期一次积分。返回能再次激发的最短 S1-S2 间期。"""ifSI_listisNone:SI_list=np.arange(40,520,10).astype(float)N=len(SI_list)V=np.full(N,-85.);h=hinf(V,-60.,3.);s=xinf(V,-25.,7.)f=hinf(V,-45.,5.);w=xinf(V,-15.,8.)fired=np.zeros(N,bool);NT=int((SI_list[-1]+200)/dt)foritinrange(NT):t=it*dt I=np.where(t<1.,40.,0.)+np.where((t>=SI_list)&(t<SI_list+1.),40.,0.)V+=dt*(-ionic(V,h,s,f,w)+I)h,s,f,w=gate_step(V,h,s,f,w,dt)fired|=(t>=SI_list+1.)&~fired&(V>-30.)# S2 后出现可传播兴奋=已脱敏idx=np.where(fired)[0]returnfloat(SI_list[idx[0]])ifidx.sizeelsefloat('nan')erp=erp_cardiac_like()cv,*_=monodomain_1d(D=0.002)lam=(cv*10)*erp# cm/ms→mm/ms(×10)，再×ERP(ms)=波长(mm)print(f"\nERP(再激发最短间期) ≈{erp:.0f}ms；D=0.002 时 CV={cv*1000:.0f}cm/s")print(f"激动波长 λ = CV×ERP ={lam:.1f}mm")print("== 折返基质判据（环周长 ≥ 波长 → 折返可维持），以 3×7×20 mm 组织块为例 ==")for(a,b)in[(3,7),(7,20),(3,20)]:loop=2*(a+b)print(f"{a}×{b}mm 环路 周长={loop:3d}mm vs λ={lam:.1f}mm → "+("可维持折返"ifloop>=lamelse"折返自灭（波前撞上不应期波尾）"))print("== 敏感性：ERP 缩短（缺血/晚钠抑制类比）让小块组织也能折返 ==")forkin[1.0,0.5,0.25]:l2=lam*kprint(f" ERP×{k:4.2f}→ λ={l2:5.1f}mm : 3×7环(20mm) "+("可折返"if20>=l2else"不能")+f"; 7×20环(54mm) "+("可折返"if54>=l2else"不能"))

实跑输出（Python 3.10 + numpy 2.2.6 复验一致）：

D=0.001 dt=0.04 → CV= 71.5 cm/s (拟合残差 0.26 ms, n=777) D=0.002 dt=0.02 → CV= 106.7 cm/s (拟合残差 0.24 ms, n=776) CV 比=1.49 vs √2=1.41 D=0.004 dt=0.01 → CV= 154.3 cm/s (拟合残差 0.20 ms, n=776) CV 比=1.45 vs √2=1.41 D=0.008 dt=0.005 → CV= 219.9 cm/s (拟合残差 0.14 ms, n=774) CV 比=1.42 vs √2=1.41 ERP(再激发最短间期) ≈ 40 ms；D=0.002 时 CV=107 cm/s 激动波长 λ = CV×ERP = 42.7 mm 3×7 mm 环路 周长= 20 mm vs λ=42.7 mm → 折返自灭 7×20 mm 环路 周长= 54 mm vs λ=42.7 mm → 可维持折返 ERP×0.25 → λ=10.7 mm : 3×7环(20mm) 可折返

逐行要害：

1. dt = min(0.05, 0.4*dx*dx/D)：显式 Euler 解扩散的稳定性条件是dt ≤ dx²/(2D)，留 0.4 安全系数。CV 比接近理论标度 √2 本身就说明网格/步长进入了"波速已收敛区"——这行是插件里最便宜的自检。
2. t_up记录"每点首次越过 −10 mV"的时刻表，polyfit一次线性拟合出 1/CV；残差即行波质量指标（本文 0.14–0.26 ms），残差大=波还没进入稳定传播段或被边界反射污染。
3. ERP 用向量化 S1-S2 程序刺激（所有间期一次积分），与电生理实验室的期前刺激协议同构——把实验范式写成代码，比背公式更重要。
4. 玩具参数下 CV 数量级（几十~两百 cm/s）不代表人心室肌真实值；把它当"CV 对 D 的标度关系与测量方法"的验证。真实定量结论交给cepModel_TTP官方算例（传导各向异性还需要纤维方向数据，缺失时结果系统性偏快——svFSI 官方边界）。

相似 API 对比：svFSI 的 CEP 是"求解器内建、输入文件配置"；openCARP 的对应物是.inp + carputils编排并可经CellML 导入外部细胞模型（许可证提醒：openCARP 本体非商业！）；Chaste 则是"代码生成"路线（CellML→C++，chaste_codegen），Myokit 能导出 C/CUDA/OpenCL 代码嵌入任意 PDE 框架——插件层选择多，契约层不变：进来的是阻滞率数组，出去的是 APD/CV/λ。

三、常见报错与排查

1. ValueError: The truth value of an array with more than one element is ambiguous。现象：第 4 篇能跑的标量 HH 速率函数搬进 1D 网格立刻崩。根因：def f(u): return 10.0 if abs(u)<1e-7 else ...是标量写法，numpy 数组不能进if。解法：np.where向量化奇点替换（本文调试实录：改为np.where(np.abs(u)<1e-6, 10.0, u/(1.0-np.exp(-u/10.0)))并用np.errstate屏蔽 0/0 警告）。
2. 全域瞬间一起兴奋（CV 测出来是天文数字）。现象：激活时刻表平坦、拟合斜率≈0、CV≈1e15 cm/s。根因两种：扩散稳定条件被违反导致数值爆掉后被截断，或反应项激活门"瞬时化"+强耦合让传播退化为非再生性扩散。解法：先dt ≤ 0.4dx²/D，再让激活门有有限时间常数；本文最终用与第 6 篇同源的最小平台期细胞作反应项（HH 瞬时 m 门不适合粗网格）。
3. RuntimeWarning: overflow encountered in exp/invalid value in add随后 NaN 一片。根因：exp参数没有 clip（V 一旦过冲 ±500 mV 即溢出污染全场）。解法：所有 Boltzmann 写np.clip(..., -50, 50)；批量作业加 NaN 自检行（铁律 8/10）。
4. “Chaste 是 GPL 所以商用要慎”。这是旧网页/旧课件的过时说法：Chaste 2026.1 为 3-clause BSD，可商用链接。反例要小心：FDA/CiPA 的 R 仓库确是GPL-3.0，闭源分发插件须做进程级隔离（第 6 篇细讲）——别把两件事混着记。
5. 把 openCARP 直接接进商用插件。许可证是Academic Public License v1.1（非商业），商用须向 NumeriCor 购买；其子组件再分层：carputils 是 Apache-2.0、meshalyzer 是 GPL-3.0——链接策略逐个核（铁律 3）。

四、动手练习

1. 把D=0.002改0.008，确认 CV 比在 1.4±0.1（判定：输出两行 CV 相除）；再把dx加密一倍、D 不变，确认 CV 变化 <10%（网格收敛检查）。
2. 令SI_list步长从 10 ms 改 2 ms，观察 ERP 估值的抖动是否 <5 ms（判定：两次输出差值）；解释为什么 S1-S2 协议是"测"出来的不应期而不是模型参数。
3. 用官方svFSI-Tests复跑08-cep/03-benchmark_tTP（WSL 或 Linux），与你所在小组的插件 CV/APD 输出按 Niederer N-version 方式列表对比（判定：产出一张三行对照表即可交付）。

五、小结与下一篇预告

本篇把第 6 篇的单细胞输出推到组织：单域方程 + 最小平台期细胞就能实测 CV∝√D 与 λ=CV×ERP，并套在 3×7×20 mm 的 benchmark 几何上判折返基质；svFSI/svMultiPhysics 的cepModel_TTP/CepModTtp与官方算例库是插件的现成宿主，openCARP/Chaste 各有分工但许可证先核再用。细胞层参数（APD/ERP）来自第 6 篇的 qNet 体系，组织层的 ERP/CV 又将成为第 8 篇器官级伪 ECG 的输入。下一篇（08）跨到器官尺度：三透壁 AP 组装伪 ECG，提取 J-Tpeak 与 Tpeak-Tend，接通 CiPA 第四工作流的风险评分。

本篇认知问题回显（FAQ）

Q1：心脏组织仿真的单域模型和双域模型区别是什么，插件该选哪个？
A：单域（monodomain）只解跨膜电位 V，方程∂V/∂t=D∇²V−I_ion/C_m+I_stim；双域（bidomain）解胞内胞外两场，可出胞外电位/伪ECG。做传导速度、波长、折返筛机制选单域（计算省一半以上），要与真实体表/电极信号比对才上双域。

Q2：svFSI 和 svMultiPhysics 内建了哪些心脏电生理模型，官方算例在哪里？
A：模型注册表为cepModel_AP（Aliev-Panfilov）、cepModel_FN（Fitzhugh-Nagumo）、cepModel_BO（Bueno-Orovio）、cepModel_TTP（ten Tusscher-Panfilov 全离子，类CepModTtp含I_Na/I_bNa）。算例在 GitHub 组织仓库svFSI-Tests：08-cep/01-2Dsqr_AP、03-benchmark_tTP、04-2Dspiral_BO、05-Purkinje；刺激网络用 SimTK 官方 Purkinje Plugin。

Q3：openCARP 和 Chaste 的许可证分别是什么，商用插件能直接用吗？
A：openCARP 本体是 Academic Public License v1.1（非商业，商用向 NumeriCor 购买），carputils 为 Apache-2.0、meshalyzer 为 GPL-3.0，需分层核查；Chaste 2026.1 为 3-clause BSD 可商用，"Chaste 是 GPL"是过时误传。svMultiPhysics/svZeroDSolver 为 BSD-3、svMorph 为 MIT。

Q4：Niederer 2011 N-version benchmark 在电生理插件验证中的作用是什么？
A：它是多团队在相同 3×7×20 mm 几何、相同 TP06 类模型上的互测基准（svFSI-Tests 提供08-cep/03-benchmark_tTP算例，模板svFSI_master.inp）。插件接入新求解器后端后先复跑该算例对齐参考数值，再进药物场景，能把"实现错误"与"科学分歧"分开，属于模型验证的第一道闸。

Q5：传导速度和有效不应期怎么组合成折返风险判据？
A：激动波长 λ=CV×ERP；折返环周长 ≥ λ 时波前遇到的是已脱不应期的组织，折返可维持。例：本文实跑 CV=107 cm/s、ERP=40 ms → λ=42.7 mm，3×7 mm 小环（20 mm）不足以折返，7×20 mm 环（54 mm）可以；ERP 缩到 1/4 时小环也变危险——这就是插件在组织层要产出的风险语言。