符号学习初学代码——从开普勒第三定律到万有引力定律

文章探讨了SymbolicRegression(SR)和Physics-InformedNeuralNetworks(PINN)在解决热传导方程中的应用。首先,通过SR生成了一个数学模型,然后利用PINN进行训练,实现了数值解的预测。展示了如何在Python中实现这两者的结合以及训练过程中的损失函数和优化策略。
摘要由CSDN通过智能技术生成
备注

PINN——physics informed neural network   SR——symbolic regression

代码详细分析见评论区链接

一、SR_test

import numpy as np

T = np.array([0.241,0.615,1,1.881,11.862]).reshape(-1, 1)
R = np.array([0.381,0.723,1,1.524,5.023]).reshape(-1, 1)
M = np.array([0.056,0.8364,1,0.1089,322.27]).reshape(-1, 1)
X = np.column_stack((M.ravel(), R.ravel())) 
Y = M * R / T**2
from pysr import PySRRegressor

model = PySRRegressor(
    niterations=100,  
    binary_operators=["+", "-", "*", "/", "^"],
    unary_operators=[
        #"cos","sin","tan",
        "exp","log",
        #"sqrt",
        #"inv(x) = 1/x",
    ],
    maxsize=7,
    maxdepth=3,
    #extra_sympy_mappings={"inv": lambda x: 1 / x},
    elementwise_loss="loss(prediction, target) = (prediction - target)^2",
)
equation = model.fit(X, Y.ravel())
print(equation)
from sympy import symbols, sympify

def predict(x, equation):
    # Parse the equation string into a SymPy expression
    x_symbols = symbols('x:' + str(len(x)))
    expr = sympify(equation)
    
    # Substitute the values of x into the expression
    expr_subs = expr.subs(list(zip(x_symbols, x)))
    
    # Evaluate the expression
    return expr_subs.evalf()

equation = '(x0 ^ 0.9818996) / (x1 * x1)'
t = 29.457
r = 9.539
m = 96.5
f = m*r/t**2
print(f)

new_data = [m,r]
predictions = predict(new_data, equation)
print(predictions)

二、SR+PINN

import numpy as np 
import matplotlib.pyplot as plt 
import torch
import torch.nn as nn
import torch.optim as optim

device = "cuda" if torch.cuda.is_available() else "cpu" #条件表达式

# same_seed 确定结果一致,保证实验的可重复性
def same_seed(seed):
    np.random.seed(seed) #确保NumPy生成的随机数也是可重复的
    torch.manual_seed(seed) #PyTorch中使用的一些随机数生成操作也会是可重复的
    if torch.cuda.is_available():
        torch.cuda.manual_seed(seed) #设置当前GPU的随机数种子
        torch.cuda.manual_seed_all(seed) #设置所有可用的GPU的随机数种子
    torch.backends.cudnn.benchmark = False #关闭benchmark模式意味着每次运行时都使用相同的卷积算法,从而提高可重复性
    torch.backends.cudnn.deterministic = True #确保每次运行时相同输入产生相同的输出,增加实验的可控性
same_seed(9258)
'''seed_value = np.random.randint(1, 10000)  
same_seed(seed_value)
current_state_tuple = np.random.get_state()
print(f"当前的随机数生成器种子是: {current_state_tuple[1][0]}")'''

# 搭建网络
class PINN(nn.Module):
    def __init__(self):
        super(PINN, self).__init__()
        self.fc1 = nn.Linear(2, 50)  
        self.fc2 = nn.Linear(50, 50)
        self.fc3 = nn.Linear(50, 1)

    def forward(self, x):
        input_data = x
        h1 = torch.tanh(self.fc1(input_data))
        h2 = torch.tanh(self.fc2(h1))
        output = self.fc3(h2)
        return output


# 定义训练函数
def train_PINN(model, epochs, lr, hist):
    optimizer = optim.Adam(model.parameters(), lr=lr)
    mse_loss = nn.MSELoss()

    for step in range(1, epochs+1): #固定迭代次数
        T = np.array([0.241,0.615,1,1.881,11.862]).reshape(-1, 1)
        R = np.array([0.381,0.723,1,1.524,5.023]).reshape(-1, 1)
        M = np.array([0.056,0.8364,1,0.1089,322.27]).reshape(-1, 1) 
        X = np.column_stack((M.ravel(), R.ravel())) 
        Y = M * R / T**2  
        T = torch.tensor(T, dtype=torch.float32)
        R = torch.tensor(R, dtype=torch.float32)
        M = torch.tensor(M, dtype=torch.float32)
        X = torch.tensor(X, dtype=torch.float32)
        Y = torch.tensor(Y, dtype=torch.float32)
        y_hat = model(X) 
        data_loss = mse_loss(y_hat, Y)
        #equation_loss = mse_loss(y_hat*(R**1.919), M**0.95927) # PDE微分函数
 
        optimizer.zero_grad() #梯度不会累积,因为 PyTorch 默认会在进行反向传播时累积梯度
        loss = data_loss #+ equation_loss
        loss.backward()
        optimizer.step() #更新模型参数
        hist.append(loss.item())
        if step % 1000 == 0:
            print(f'{step}/{epochs} loss {loss.item():.5f}')
       

#参数设置
lr = 1e-3 #学习率
epochs=8000 #迭代次数
hist = []


# 训练 PINN 模型
model=PINN() #初始化
train_PINN(model, epochs, lr, hist) 

plt.plot(hist)
plt.xlabel("update step")
plt.ylabel("loss")
#plt.yscale("log")
plt.grid(True)
plt.show()

model.eval() #将模型切换到评估模式
T = np.array([0.615,1,1.881,11.862]).reshape(-1, 1)
R = np.array([0.723,1,1.524,5.023]).reshape(-1, 1)
M = np.array([0.8364,1,0.1089,322.27]).reshape(-1, 1) 
X = np.column_stack((M.ravel(), R.ravel())) 
Y = M * R / T**2  
T = torch.tensor(T, dtype=torch.float32)
R = torch.tensor(R, dtype=torch.float32)
M = torch.tensor(M, dtype=torch.float32)
X = torch.tensor(X, dtype=torch.float32)
Y = torch.tensor(Y, dtype=torch.float32)
with torch.no_grad(): #关闭 PyTorch 的梯度追踪
    F_predict = model(X)
F_true = Y
l2 = np.sqrt(np.sum((F_predict.numpy()-F_true.numpy())**2))
print(l2)

评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值