基于TensorFlow的PINN求解热传导方程时物理损失计算报错

我尝试用TensorFlow实现物理信息神经网络（PINN）求解热传导方程 C·T'+K·T=Q'，其中T是杆体温度，C是51×51对角热容矩阵，K是51×51对角热导矩阵，Q'是15193×51的热载荷。我有15193组温度数据（每60秒采集一次51个节点温度），神经网络输入为时间，输出为51个节点温度，但计算物理损失时出现错误。

原代码 #

import os 
os.environ['TF_CPP_MIN_LOG_LEVEL']='2'
import numpy as np
import tensorflow as tf

# Constants
t0 = 0
dt = 60
t_end = 15193 * dt

# Load data
C = np.genfromtxt("C_th.txt", dtype=float, delimiter=",")
K = np.genfromtxt("K_th.txt", dtype=float, delimiter=",")
Q_prime = np.genfromtxt("thermal_loads.txt", dtype=float, skip_header=1, delimiter=",")[:15193]
temperatures = np.genfromtxt("temperature_rod.txt", dtype=float, skip_header=1, delimiter=",")[:15193]

# Convert to TensorFlow tensors and appropriate dtypes
C = tf.cast(C, tf.float32)
K = tf.cast(K, tf.float32)
Q_prime = tf.cast(Q_prime, tf.float32)
temperatures = tf.cast(temperatures, tf.float32)

# Neural Network Model
class PINNModel(tf.keras.Model):
    def __init__(self):
        super(PINNModel, self).__init__()
        self.dense_1 = tf.keras.layers.Dense(20, activation=tf.nn.tanh)
        self.dense_2 = tf.keras.layers.Dense(20, activation=tf.nn.tanh)
        self.dense_3 = tf.keras.layers.Dense(20, activation=tf.nn.tanh)
        self.output_layer = tf.keras.layers.Dense(51, activation=None)

    def call(self, inputs):
        x = self.dense_1(inputs)
        x = self.dense_2(x)
        x = self.dense_3(x)
        return self.output_layer(x)

# Preparing the data
time_inputs = np.arange(t0, t_end, dt)
time_inputs = tf.reshape(tf.cast(time_inputs, dtype=tf.float32), (-1, 1))  # Convert to column vector
dataset = tf.data.Dataset.from_tensor_slices((time_inputs, temperatures)).batch(32)

# Define the model and optimizer
model = PINNModel()
optimizer = tf.keras.optimizers.Adam(learning_rate=0.001)

# Loss functions
def data_loss(predictions, true_temperatures):
    return tf.reduce_mean(tf.square(predictions - true_temperatures))

def physics_loss(time_inputs, predictions, C, K, Q_prime):
    with tf.GradientTape(persistent=True) as tape:
        tape.watch(time_inputs)
        dT_dt = tape.gradient(predictions, time_inputs)
    C_dT_dt = C @ dT_dt
    KT = K @ predictions
    residuals = C_dT_dt + KT - Q_prime
    return tf.reduce_mean(tf.square(residuals))


def combined_loss(time_inputs, true_temperatures, predictions, C, K, Q_prime):
    return data_loss(predictions, true_temperatures) + physics_loss(time_inputs, predictions, C, K, Q_prime)

predictions = model(time_inputs)
print("predictions====================", predictions, tf.shape(predictions))
print("time===============", time_inputs, tf.shape(time_inputs))

print(physics_loss(time_inputs, predictions, C, K, Q_prime))

报错信息 #

ValueError: Attempt to convert a value (None) with an unsupported type (<class 'NoneType'>) to a Tensor.

问题原因 #

1. 梯度追踪失效：physics_loss中传入的predictions是预先计算好的模型输出，和当前tape.watch的time_inputs没有梯度关联，导致tape.gradient返回None。
2. 维度不匹配：普通矩阵乘法C @ dT_dt无法适配批量数据的维度（C是51×51，dT_dt初始形状为(15193,1)），会引发维度错误。

修复后的完整代码 #

import os 
os.environ['TF_CPP_MIN_LOG_LEVEL']='2'
import numpy as np
import tensorflow as tf

# Constants
t0 = 0
dt = 60
t_end = 15193 * dt

# Load data
C = np.genfromtxt("C_th.txt", dtype=float, delimiter=",")
K = np.genfromtxt("K_th.txt", dtype=float, delimiter=",")
Q_prime = np.genfromtxt("thermal_loads.txt", dtype=float, skip_header=1, delimiter=",")[:15193]
temperatures = np.genfromtxt("temperature_rod.txt", dtype=float, skip_header=1, delimiter=",")[:15193]

# Convert to TensorFlow tensors and appropriate dtypes
C = tf.cast(C, tf.float32)
K = tf.cast(K, tf.float32)
Q_prime = tf.cast(Q_prime, tf.float32)
temperatures = tf.cast(temperatures, tf.float32)

# Neural Network Model
class PINNModel(tf.keras.Model):
    def __init__(self):
        super(PINNModel, self).__init__()
        self.dense_1 = tf.keras.layers.Dense(20, activation=tf.nn.tanh)
        self.dense_2 = tf.keras.layers.Dense(20, activation=tf.nn.tanh)
        self.dense_3 = tf.keras.layers.Dense(20, activation=tf.nn.tanh)
        self.output_layer = tf.keras.layers.Dense(51, activation=None)

    def call(self, inputs):
        x = self.dense_1(inputs)
        x = self.dense_2(x)
        x = self.dense_3(x)
        return self.output_layer(x)

# Preparing the data
time_inputs = np.arange(t0, t_end, dt)
time_inputs = tf.reshape(tf.cast(time_inputs, dtype=tf.float32), (-1, 1))  # Convert to column vector
dataset = tf.data.Dataset.from_tensor_slices((time_inputs, temperatures)).batch(32)

# Define the model and optimizer
model = PINNModel()
optimizer = tf.keras.optimizers.Adam(learning_rate=0.001)

# Loss functions
def data_loss(predictions, true_temperatures):
    return tf.reduce_mean(tf.square(predictions - true_temperatures))

def physics_loss(time_inputs, model, C, K, Q_prime_batch):
    with tf.GradientTape(persistent=True) as tape:
        tape.watch(time_inputs)
        # 在tape内重新计算模型输出，确保梯度可追踪
        predictions = model(time_inputs)
        # 对每个节点的预测值求时间导数，得到(批量大小,51)的导数矩阵
        dT_dt = tf.stack([tape.gradient(predictions[:, i], time_inputs) for i in range(51)], axis=1)
        dT_dt = tf.squeeze(dT_dt, axis=2)  # 调整为(批量大小,51)的形状
    
    # 用matvec处理批量矩阵-向量乘法，适配维度
    C_dT_dt = tf.linalg.matvec(C, dT_dt)
    KT = tf.linalg.matvec(K, predictions)
    residuals = C_dT_dt + KT - Q_prime_batch
    return tf.reduce_mean(tf.square(residuals))

def combined_loss(time_inputs, true_temperatures, model, C, K, Q_prime_batch):
    predictions = model(time_inputs)
    return data_loss(predictions, true_temperatures) + physics_loss(time_inputs, model, C, K, Q_prime_batch)

# 测试物理损失计算（取第一个batch验证）
for batch_time, batch_temp in dataset.take(1):
    batch_q = Q_prime[:batch_time.shape[0]]
    print("物理损失值：", physics_loss(batch_time, model, C, K, batch_q))

关键修复点 #

• 梯度追踪修正：在physics_loss内部调用model(time_inputs)生成预测值，确保GradientTape能捕捉输入到输出的梯度流，避免返回None。
• 导数维度适配：通过循环对每个节点的预测值求导，得到和输出维度一致的导数矩阵，适配后续运算。
• 批量矩阵运算：使用tf.linalg.matvec替代普通矩阵乘法，自动处理批量数据的维度匹配问题。
• 热载荷批量匹配：训练时传入对应批量的Q'数据，确保时间步与热载荷一一对应。

内容的提问来源于stack exchange，提问作者Mayssa