论文阅读:Failure-informed adaptive sampling for PINNs, Part II: combining with re-sampling and subset simulation
Failure-informed adaptive sampling for PINNs, Part II: combining with re-sampling and subset simulation
Annealing failure-informed PINNs
上回说到,FI-PINN会在训练过程中不断往训练集中增加训练点,并不能人为限制训练点的数量,这就使得训练过程变得不可控,同时难以预计训练所要花费的时间。作者受到R3 Sampling的影响,提出了AFI-PINN。
AFI-PINN 框架分为两个阶段。首先是确定何时更新训练集,这称为重启阶段。第二个是具体如何对点进行重采样。对于第一阶段,一个简单的解决方案是定期交换点,即在整个训练过程中每经过给定数量的epoch就重新采样。作者在本文中采取的做法是,当训练损失在几个时期内不再减少时,就对搭配点进行重新采样。
比较特殊的是使用退火方式的自适应采样过程。当需要在第 T restart T_\text {restart} Trestart 个epoch重新采样时,首先使用采样方法(例如上回使用的 SAIS 方法或下文中使用的子集模拟)估计失败概率 P ^ F \hat P_F P^F ,并同时从故障区域生成 N F N_F NF 个搭配点组成的候选点集 D adaptive D_\text{adaptive} Dadaptive。如果失败概率 P ^ F \hat P_F P^F 小于预定义的容差 ρ \rho ρ,则可以提前终止整个训练。否则就要使用余弦退火方法更新训练点。具体来说,新的训练集将由三部分组成:自适应搭配数据集 D adaptive D_\text{adaptive} Dadaptive、从初始数据集 D c D_c Dc 生成的大小为 N c ( 1 − η ) N_c(1 − \eta) Nc(1−η) 的预选数据集 D ^ c \hat D_c D^c ,以及从先验分布中生成大小为 η N c − N F \eta N_c − N_F ηNc−NF 的新数据集 D prior D_\text{prior} Dprior。
其中, 0 < η < 1 0 < \eta < 1 0<η<1 是重采样比例。通过这样的方法,就可以使更新后数据点的数量保持不变。其中,包含预选点的原因是网络需要保留之前训练中学到的特征。而自适应样本和先验样本有助于支持网络分别从故障区域和问题域的其他部分学习新特征。这种类型的组合将极大地提高训练效率,同时还降低计算成本并加快收敛速度。
为了确定重采样比例 η,作者使用了基于重新启动和最大训练轮数的余弦退火方式。通过这个方法可以自适应地选择数据集的比例,以保持新数据集不同部分之间的平衡。具体来说,假设第
k
k
k 个重新采样时期为
T
restart
(
k
)
T^{(k)}_\text{restart}
Trestart(k),总训练轮数为
T
max
T_\text{max}
Tmax。当重新采样时,比例可以计算为:
η
=
a
(
1
+
b
cos
(
π
(
T
r
e
s
t
a
r
t
(
k
)
−
T
s
(
k
)
)
T
m
a
x
−
T
s
(
k
)
)
)
\eta=a\left(1+b\cos\left(\frac{\pi(T_{restart}^{(k)}-T_{s}^{(k)})}{T_{max}-T_{s}^{(k)}}\right)\right)
η=a(1+bcos(Tmax−Ts(k)π(Trestart(k)−Ts(k))))
其中,
a
a
a 和
b
b
b 是两个正数,用来控制采样比例的上下界,
T
s
(
k
)
T^{(k)}_s
Ts(k) 为上一次重新采样时的轮数,即
T
s
(
k
)
=
T
restart
(
k
−
1
)
T^{(k)}_s = T^{(k-1)}_\text{restart}
Ts(k)=Trestart(k−1)
上图展示了训练过程中比例的趋势,其中 η 1 = a ( 1 + b ) ∈ ( 0 , 1 ) \eta 1 = a(1 + b) \in (0, 1) η1=a(1+b)∈(0,1), η 2 = a ( 1 − b ) ∈ ( 0 , 1 ) \eta 2 = a(1 − b) \in (0, 1) η2=a(1−b)∈(0,1)。当 N F < N c η ≔ N s N_F < N_c\eta \coloneqq N_s NF<Ncη:=Ns 时,可以从先验中抽取 N s − N F N_s − N_F Ns−NF 个样本,记为 D prior D_\text{prior} Dprior。否则,可以通过从候选自适应数据集中抽取 N s N_s Ns 个样本来构建自适应样本。在这种情况下, D prior = ⊘ D_\text{prior} = \oslash Dprior=⊘。
完整算法如下:
子集模拟
作者在文章中采用基于MCMC方法的子集模拟(SS)来完成自适应采样。 SS 是最常用的方差减少技术之一,作者在本文中用其估计失效概率
P
^
F
\hat P_F
P^F 。子集模拟的主要思想是定义一系列嵌套的失败区域,其中包括目标失败区域
Ω
F
\Omega_F
ΩF
Ω
F
⊂
Ω
F
m
⊂
Ω
F
m
−
1
⊂
…
Ω
F
1
⊂
Ω
F
0
=
Ω
\Omega_F \subset \Omega_{F_m} \subset \Omega_{F_{m-1}} \subset \dots \Omega_{F_1} \subset \Omega_{F_0} = \Omega
ΩF⊂ΩFm⊂ΩFm−1⊂…ΩF1⊂ΩF0=Ω
于是,失败概率可以被分解如下:
P
F
=
P
(
Ω
F
)
=
P
(
Ω
F
1
∣
Ω
F
0
)
P
(
Ω
F
2
∣
Ω
F
1
)
…
P
(
Ω
F
m
∣
Ω
F
m
−
1
)
P
(
Ω
F
∣
Ω
F
m
)
P_F = P(\Omega_F)=P(\Omega_{F_1}|\Omega_{F_0})P(\Omega_{F_2}|\Omega_{F_1})\dots P(\Omega_{F_m}|\Omega_{F_{m-1}})P(\Omega_{F}|\Omega_{F_m})
PF=P(ΩF)=P(ΩF1∣ΩF0)P(ΩF2∣ΩF1)…P(ΩFm∣ΩFm−1)P(ΩF∣ΩFm)
其中,
P
(
Ω
F
k
+
1
∣
Ω
F
k
)
=
P
(
Ω
F
k
+
1
)
P
(
Ω
F
k
)
,
0
≤
k
≤
m
−
1
P(\Omega_{F_{k+1}}|\Omega_{F_k})= \frac{P(\Omega_{F_{k+1}})}{P(\Omega_{F_{k}})}, 0 \le k \le m-1
P(ΩFk+1∣ΩFk)=P(ΩFk)P(ΩFk+1),0≤k≤m−1 表示中间区域
Ω
F
m
\Omega_{F_m}
ΩFm 中失败概率的条件概率,
P
(
Ω
F
∣
Ω
F
m
)
P(\Omega_{F}|\Omega_{F_m})
P(ΩF∣ΩFm) 表示最终的条件概率。
在子集模拟的过程中,中间区域失败的条件概率被设定为一个固定值
ρ
∈
(
0
,
1
)
\rho \in (0,1)
ρ∈(0,1) 。假设
N
p
=
N
s
ρ
N_p = N_s \rho
Np=Nsρ 与
ρ
−
1
\rho ^{-1}
ρ−1 均为正数,那么,
N
p
N_p
Np 就是马尔科夫链的数量,
ρ
−
1
\rho ^{-1}
ρ−1 就是每个链的采样数目。作者根据其他SS文献,将
ρ
\rho
ρ 设置为了
0.1
0.1
0.1 。采样数目
N
s
N_s
Ns 则控制着结果的精度。具体来说,子集模拟首先从无条件概率开始生成
N
s
N_s
Ns 。使用直接蒙特卡罗方法从先验分布
ω
(
x
)
\omega (x)
ω(x) 获取
x
x
x 个样本。计算相应的
Q
Q
Q 值并按降序排列,得到一个有序列表
S
0
=
{
x
0
(
i
)
}
i
=
1
N
s
S_0 = \{x^{(i)}_0 \}^{N_s}_{i=1}
S0={x0(i)}i=1Ns ,此时,中间失效区域为 $ \Omega_{F_1} \coloneqq {x : Q(x) >\epsilon ^{(1)}_r }$,其中
ϵ
r
(
1
)
\epsilon ^{(1)}_r
ϵr(1) 为
x
x
x 在第 0 次中第
(
N
p
+
1
)
(N_p + 1)
(Np+1) 大样本值,即
ϵ
r
(
1
)
=
x
0
(
N
p
+
1
)
\epsilon ^{(1)}_r = x^{(N_p+1)}_0
ϵr(1)=x0(Np+1)。以此类推,可以得到第
k
+
1
k+1
k+1 次时中间失败区域表示如下:
Ω
F
k
+
1
≔
{
x
:
Q
(
x
)
>
ϵ
r
(
k
+
1
)
}
\Omega_{F_{k+1}} \coloneqq \{x : Q(x) >\epsilon ^{(k+1)}_r \}
ΩFk+1:={x:Q(x)>ϵr(k+1)}
于是,条件概率
P
(
Ω
F
k
+
1
∣
Ω
F
k
)
P(\Omega_{F_{k+1}}|\Omega_{F_k})
P(ΩFk+1∣ΩFk) 可以表示如下:
P
(
Ω
F
k
+
1
∣
Ω
F
k
)
=
∫
Ω
F
k
+
1
ρ
(
x
∣
Ω
F
k
)
d
x
≈
1
N
s
∑
i
=
1
N
s
I
Ω
F
k
+
1
(
x
k
(
i
)
)
=
p
.
\mathcal{P}(\Omega_{\mathcal{F}_{k+1}}|\Omega_{\mathcal{F}_{k}})=\int_{\Omega_{\mathcal{F}_{k+1}}}\rho(\mathbf{x}|\Omega_{\mathcal{F}_{k}})d\mathbf{x}\approx\frac{1}{N_{s}}\sum_{i=1}^{N_{s}}\mathcal{I}_{\Omega_{\mathcal{F}_{k+1}}}(\mathbf{x}_{k}^{(i)})=p.
P(ΩFk+1∣ΩFk)=∫ΩFk+1ρ(x∣ΩFk)dx≈Ns1i=1∑NsIΩFk+1(xk(i))=p.
值得注意的是,样本
{
x
k
(
i
)
}
i
=
1
N
p
\{x^{(i)}_k \}^{N_p}_{i=1}
{xk(i)}i=1Np 也遵循条件密度
ρ
(
⋅
∣
Ω
F
k
+
1
)
\rho(\cdot|\Omega_{\mathcal{F}_{k+1}})
ρ(⋅∣ΩFk+1) 。因此,为了保持样本大小恒定,需要从
ρ
(
⋅
∣
Ω
F
k
+
1
)
\rho(\cdot|\Omega_{\mathcal{F}_{k+1}})
ρ(⋅∣ΩFk+1) 生成最终的
(
N
s
−
N
p
)
(N_s − N_p)
(Ns−Np) 样本。在实践中,修改后的Metropolis Hasting算法(MMA)可用于生成初始种子
{
x
k
(
i
)
}
i
=
1
N
p
\{x^{(i)}_k \}^{N_p}_{i=1}
{xk(i)}i=1Np的
N
p
N_p
Np 条链。在采样过程中,每个链将接受 $ \frac 1 p − 1$ 个新样本。因此,新的样本集
S
k
+
1
S_{k+1}
Sk+1可以表示为
{
x
k
+
1
(
i
)
}
i
=
1
N
s
=
{
x
k
+
1
(
i
)
}
i
=
1
N
s
−
N
p
∪
{
x
k
(
i
)
}
i
=
1
N
p
\{\mathbf{x}_{k+1}^{(i)}\}_{i=1}^{N_{s}}=\{\mathbf{x}_{k+1}^{(i)}\}_{i=1}^{N_{s}-N_{p}}\cup\{\mathbf{x}_{k}^{(i)}\}_{i=1}^{N_{p}}
{xk+1(i)}i=1Ns={xk+1(i)}i=1Ns−Np∪{xk(i)}i=1Np
不断重复上述过程,第
k
k
k 次时失败区域的样本数量可以表示为:
N
F
k
=
∑
i
=
1
N
s
I
Ω
F
(
x
k
(
i
)
)
.
N_{\mathcal{F}_k}=\sum_{i=1}^{N_s}\mathcal{I}_{\Omega_{\mathcal{F}}}(\mathbf{x}_k^{(i)}).
NFk=i=1∑NsIΩF(xk(i)).
如果
N
p
>
N
F
k
N_p > N_{F_k}
Np>NFk,因为中间失效区域接近真实失效区域,则可以终止仿真。如果模拟在
k
=
m
k = m
k=m 处结束,则条件概率
P
(
Ω
F
∣
Ω
F
m
)
P(\Omega_{F}|\Omega_{F_m})
P(ΩF∣ΩFm) 可以估计为:
P
(
Ω
F
∣
Ω
F
m
)
=
∫
Ω
F
ρ
(
x
∣
Ω
F
m
)
d
x
≈
1
N
s
∑
i
=
1
N
s
I
Ω
F
(
x
m
(
i
)
)
=
N
F
m
N
s
:
=
q
\mathcal{P}(\Omega_{\mathcal{F}}|\Omega_{\mathcal{F}_{m}})=\int_{\Omega_{\mathcal{F}}}\rho(\mathbf{x}|\Omega_{\mathcal{F}_{m}})d\mathbf{x}\approx\frac{1}{N_{s}}\sum_{i=1}^{N_{s}}\mathcal{I}_{\Omega_{\mathcal{F}}}(\mathbf{x}_{m}^{(i)})=\frac{N_{\mathcal{F}_{m}}}{N_{s}}:=q
P(ΩF∣ΩFm)=∫ΩFρ(x∣ΩFm)dx≈Ns1i=1∑NsIΩF(xm(i))=NsNFm:=q
此时,失败概率可作如下估计:
P
F
=
P
(
Ω
F
1
∣
Ω
F
0
)
P
(
Ω
F
2
∣
Ω
F
1
)
⋯
P
(
Ω
F
m
∣
Ω
F
m
−
1
)
P
(
Ω
F
∣
Ω
F
m
)
≈
p
m
q
=
P
^
F
S
S
.
\begin{aligned}P_{\mathcal{F}}&=\mathcal{P}(\Omega_{\mathcal{F}_1}|\Omega_{\mathcal{F}_0})\mathcal{P}(\Omega_{\mathcal{F}_2}|\Omega_{\mathcal{F}_1})\cdots\mathcal{P}(\Omega_{\mathcal{F}_m}|\Omega_{\mathcal{F}_{m-1}})\mathcal{P}(\Omega_{\mathcal{F}}|\Omega_{\mathcal{F}_m})\\&\approx p^mq=\hat{P}_{\mathcal{F}}^{SS}.\end{aligned}
PF=P(ΩF1∣ΩF0)P(ΩF2∣ΩF1)⋯P(ΩFm∣ΩFm−1)P(ΩF∣ΩFm)≈pmq=P^FSS.
完整算法如下:
实验结果
作者在具有两个和四个峰值的二维问题、一个瞬态波动方程和一个十维泊松方程共四个样例上测试了R-FIPINN(基于残差的)、G-FIPINN(基于残差梯度的)。
The time-dependent wave equation
∂ 2 u ∂ t 2 − 3 ∂ 2 u ∂ x 2 = 0 , ( t , x ) ∈ [ 0 , 6 ] × [ − 5 , 5 ] , u ( 0 , x ) = 1 cosh ( 2 x ) − 0.5 cosh ( 2 ( x − 10 ) ) − 0.5 cosh ( 2 ( x + 10 ) ) , ∂ u ∂ t ( 0 , x ) = 0 , , u ( t , − 5 ) = u ( t , 5 ) = 0 , \begin{aligned} &\frac{\partial^{2}u}{\partial t^{2}}-3\frac{\partial^{2}u}{\partial x}^{2}=0,\quad(t,x)\in[0,6]\times[-5,5], \\ &u(0,x)=\frac{1}{\cosh(2x)}-\frac{0.5}{\cosh(2(x-10))}-\frac{0.5}{\cosh(2(x+10))}, \\ &\frac{\partial u}{\partial t}(0,x)=0,, \\ &u(t,-5)=u(t,5)=0, \end{aligned} ∂t2∂2u−3∂x∂2u2=0,(t,x)∈[0,6]×[−5,5],u(0,x)=cosh(2x)1−cosh(2(x−10))0.5−cosh(2(x+10))0.5,∂t∂u(0,x)=0,,u(t,−5)=u(t,5)=0,
其精确解为:
u
(
t
,
x
)
=
0.5
cosh
(
2
(
x
−
3
t
)
)
−
0.5
cosh
(
2
(
x
−
10
+
3
t
)
)
+
0.5
cosh
(
2
(
x
+
3
t
)
)
−
0.5
cosh
(
2
(
x
+
10
−
3
t
)
)
.
\begin{aligned}u(t,x)&=\frac{0.5}{\cosh(2(x-\sqrt3t))}-\frac{0.5}{\cosh(2(x-10+\sqrt3t))}\\&+\frac{0.5}{\cosh(2(x+\sqrt3t))}-\frac{0.5}{\cosh(2(x+10-\sqrt3t))}.\end{aligned}
u(t,x)=cosh(2(x−3t))0.5−cosh(2(x−10+3t))0.5+cosh(2(x+3t))0.5−cosh(2(x+10−3t))0.5.
上图为四种算法在不同采样点数量下达到的 L 2 L^2 L2 误差。
上图为三种方法训练结果与绝对误差可视化
上图为三种算法最终的训练点分布图
上图为三种方法在训练过程中损失函数的变化
上图为两种算法在训练过程中失败概率的变化。
更多样例见原文
总结
这篇文章弥补了之前FI-PINN方法在采样点数量上无法控制的不足,同时,引入子集模拟方法来估计失败概率。通过实验对比也可以看出效果确实比原来更好了。但作者还没有提供代码,因此需要自行复现。
相关链接: