pinn求解三维固体力学问题（pytorch）

pinn求解三维固体力学问题（pytorch）

论文名称

《An introduction to programming Physics-Informed Neural Network-based computational solid mechanics》

问题描述

• 几何模型：

• 完整立方体边长       。
• 利用对称性，仅建模立方体的八分之一（如右图），以降低计算成本。

• 载荷条件：

• 在立方体顶部和底部表面 (       和       ) 施加沿z方向的分布力       。
• 载荷函数为：      

• 边界条件：

• 位移边界条件 (在对称面上施加)：      

• 材料属性：

• 杨氏模量       
• 泊松比       
• 为平面应力问题，Lamé常数        由        和        通过公式计算得到。

神经网络架构设计

该PINN模型采用多个独立的前馈神经网络来近似位移场。

物理信息嵌入与损失函数

• 主要物理约束：
1. 平衡方程 (在计算域        内)：       
2. 力边界条件 (在力边界        上)：       
3. 位移边界条件 (在位移边界        上)：       
• 损失函数类型：
• 在此3D案例中，论文仅使用了配点损失函数 (Collocation Loss Function) 进行求解。
• 能量基损失函数 (Energy-based Loss Function) 可以通过修改程序实现，但3D示例中未使用。
• **损失函数构建 (Collocation)**：
• 损失函数是平衡方程残差和力边界条件残差的均方误差 (MSE) 之和：      
• 其中，       是域内样本点总数，       是力边界上的样本点数。
• 边界条件精确施加：
• 对于固定的位移边界条件（如        ），采用了硬边界的方法来精确满足，而非将其作为损失项。这简化了损失函数，降低了训练复杂度。具体方法为：      ，其中        是原始FNN输出，       是固定边界位置。

数值实现流程一些关键问题

由于在三维问题下，采样点会指数增加，例如采用    ，不可避免有内存超限问题。

训练阶段

def closure(self):
    self.optimizer.zero_grad()
    x, x1u, x1b, x2u, x2b, x3u, x3b = self.x_train

    ### 第1段：平衡方程残差（内域点）
    Gex, Gey, Gez = self.pinn.forward_ge(x)
    l1 = Loss_GE(Gex, Gey, Gez)
    l1.backward()                    
    l1_val = l1.item()
    del Gex, Gey, Gez               

    ### 第2段：x面边界条件
    s11u, s121u, s131u, s121b, s131b = self.pinn.forward_bc_x(x1u, x1b)
    l_x = Loss_BC_x(s11u, s121u, s131u, s121b, s131b)
    l_x.backward()                  
    l_x_val = l_x.item()
    del s11u, s121u, s131u, s121b, s131b

    ### 第3段：y面边界条件
    s22u, s122u, s232u, s122b, s232b = self.pinn.forward_bc_y(x2u, x2b)
    l_y = Loss_BC_y(s22u, s122u, s232u, s122b, s232b)
    l_y.backward()
    l_y_val = l_y.item()
    del s22u, s122u, s232u, s122b, s232b

    ### 第4段：z面边界条件
    s33u, s133u, s233u, s133b, s233b = self.pinn.forward_bc_z(x3u, x3b)
    l_z = Loss_BC_z(s33u, s133u, s233u, s133b, s233b, self.y_train)
    l_z.backward()
    l_z_val = l_z.item()
    del s33u, s133u, s233u, s133b, s233b

    ### 汇总 loss 数值
    loss_val = l1_val + l_x_val + l_y_val + l_z_val
    return torch.tensor(loss_val)

• 将一次大的 forward+backward 拆分为4段独立计算，每段算完立即 backward 释放计算图。由于 zero_grad() 只在开头调用一次，4次 backward() 的梯度自动累积，数学上完全等价。

后处理阶段

n = xyz.shape[0]
batch_size = 2000
s1_list, s2_list, s3_list = [], [], []
s12_list, s13_list, s23_list = [], [], []

for i in range(0, n, batch_size):
    xyz_b = torch.tensor(xyz[i:i+batch_size], dtype=torch.float32,
                         requires_grad=True).to(device)

    # 只需一组微分 + 一次 Material 调用
    U_x, U_y, U_z, U_xx, U_xy, U_xz, U_yy, U_yz, U_zz = pinn.dif_x(xyz_b)
    V_x, V_y, V_z, V_xx, V_xy, V_xz, V_yy, V_yz, V_zz = pinn.dif_y(xyz_b)
    W_x, W_y, W_z, W_xx, W_xy, W_xz, W_yy, W_yz, W_zz = pinn.dif_z(xyz_b)
    _, _, _, _, _, _, s1_t, s2_t, s3_t, s12_t, s23_t, s13_t, _, _, _ = Material(
        U_x, U_y, U_z, V_x, V_y, V_z, W_x, W_y, W_z,
        U_xx, U_xy, U_xz, U_yy, U_yz, U_zz,
        V_xx, V_xy, V_xz, V_yy, V_yz, V_zz,
        W_xx, W_xy, W_xz, W_yy, W_yz, W_zz, pinn.E, pinn.mu)

    # detach() 断开计算图，只保留数值
    s1_list.append(s1_t.detach().cpu().numpy())
    s2_list.append(s2_t.detach().cpu().numpy())

s1 = np.concatenate(s1_list)
# ... 拼接其余分量

每个批次中：

• detach() 将结果从计算图断开，只保留纯数值
• 下一轮循环中变量（xyz_b, U_x, s1_t 等）被重新赋值，旧对象的引用计数归零
• Python 引用计数机制立即回收旧对象，关联的计算图跟随释放

二维物理问题描述

• 模型：

• 几何与建模：
• 方板边长       
• 利用几何对称性，仅建模方板的四分之一（如图12(b)），以降低计算成本
• 载荷条件：
• 在板的右侧（      ）施加沿x方向的分布力      
• 载荷函数为：      
• 边界条件：
• 位移边界条件 (在对称面上施加)：      
• 力（牵引力）边界条件：已由        在右侧边界施加
• 材料属性：
• 杨氏模量       
• 泊松比       
• 为平面应力问题，Lamé常数        由以下公式计算：