PIV_2:流体中的POD分解深度解析与Python实践

0.写在前面

• 参考材料:斯坦福的教程 CME 345: MODEL REDUCTION Proper Orthogonal Decomposition (POD) & 百度网盘copy (提取码：ktdv)
• 我的中文翻译版本：百度网盘（提取码：2gvm）

1.基于Python实现的POD对流场分析应用

• 建立仿真流场
• POD学习（SVD分解）
• 基本模态展示，重构结果展示

1.1 建立仿真流场时间序列

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
sns.set()
%matplotlib inline

#- 合成100个时刻的流场
N = 32
x,y = np.meshgrid(np.linspace(0,1,N),np.linspace(0,1,N))
u,v = [],[] 
for t in range(100):
    T = 1+t/100.0
    u.append(np.cos(2*np.pi*x/T+np.pi/2)*np.cos(2*np.pi*y/T))
    v.append(np.sin(2*np.pi*x/T+np.pi/2)*np.sin(2*np.pi*y/T))

plt.figure(figsize=(12,12))
plt.quiver(u[0],v[0])
plt.title("T= 1")
plt.figure(figsize=(12,12))
plt.quiver(u[99],v[99])
plt.title("T=2")

<matplotlib.text.Text at 0x7f41902d8be0>

output_2_1.png output_2_2.png

1.2 POD分解：得到正交基函数（basic models）

#- 做POD分解
snapshot_u =  np.swapaxes(np.array(u).reshape(100,-1),0,1)
snapshot_v =  np.swapaxes(np.array(v).reshape(100,-1),0,1)
snapshot = np.concatenate((snapshot_u,snapshot_v),axis = 0)
print(snapshot.shape)
U,s,V =np.linalg.svd(snapshot) # SVD就实现了POD分解
print(U.shape,V.shape)
print(s[:10]) #可见，前10个模态就包含了足够的信息了。

(2048, 100)
(2048, 2048) (100, 100)
[  1.87380136e+02   1.12594810e+02   5.27326645e+01   1.42560358e+01
   3.27259977e+00   5.29841962e-01   3.17332170e-02   1.23968446e-03
   3.60398843e-05   8.13446663e-07]

#- 展示前5个基本模态
for i in range(4):
    model_u,model_v = U[:N**2,i].reshape(N,N),U[N**2:,i].reshape(N,N)
    plt.figure(figsize=(4,4))
    plt.quiver(model_u,model_v)
    plt.title("model:"+str(i))

output_4_0.png output_4_1.png output_4_2.png output_4_3.png

1.3 POD投影重构（reconstruction）

#- 利用前S个模态进行重构
S = 5
R = []
for i in range(100):
    snapshot_i = snapshot[:,i]
    reconstruction = np.zeros_like(snapshot_i)
    for s in range(S):
#         print(np.dot(U[:,s],snapshot_i))
        reconstruction = reconstruction + np.dot(U[:,s],snapshot_i)*U[:,s]
    R.append(reconstruction)
R = np.array(R)
R = np.swapaxes(R,0,1)

#- 重构的结果展示
for i in range(0,100,20):
    plt.figure(figsize=(16,8))
    model_u,model_v = snapshot[:N**2,i].reshape(N,N),snapshot[N**2:,i].reshape(N,N)
    plt.subplot(121)
    plt.quiver(model_u,model_v)
    plt.title("original")
    
    model_u,model_v = R[:N**2,i].reshape(N,N),R[N**2:,i].reshape(N,N)
    plt.subplot(122)
    plt.quiver(model_u,model_v)
    plt.title("Reconstruction with %d basic moda"%S)
    

output_6_0.png output_6_1.png output_6_2.png output_6_3.png output_6_4.png

2.总结

• POD分解就是一个SVD操作，和PCA人脸相同
• POD分解对降维操作，然后再进行仿真分析等，叫做Reduced Order Method (ROM)在流场分析中广泛使用
• 结合POD和光流法，参考文章"Discrete Orthogonal Decomposition and Variational Fluid Flow Estimation"
• POD分解的缺点，只能很好处理time- invariant的数据，能否处理time-variant的数据呢？DMD（Dynamic Model Decomposition），这个有空再看看，争取搞清楚。