C++：有限差分求解非粘性时变汉堡方程和拉克斯-温德罗夫方法（附带源码）-CSDN博客

本文链接：https://blog.csdn.net/m0_61840987/article/details/145081279

项目名称：C++ 使用有限差分法求解非粘性时变汉堡方程与拉克斯-温德罗夫方法

项目背景

非粘性时变汉堡方程（Non-Viscous Time-Dependent Burgers' Equation）是一个经典的偏微分方程，广泛应用于流体动力学、气体动力学、气象学等领域。它的形式为：

其中，u(x,t) 是流体的速度，x 是空间坐标，t 是时间。边界条件为：

时变汉堡方程中包含非线性对流项，数值求解通常需要采用合适的差分方法。拉克斯-温德罗夫（Lax-Wendroff）方法是求解此类方程的一种常见数值方法，尤其适用于处理非线性方程，具有较好的稳定性和精度。

项目目标

问题描述：
- 给定非粘性时变汉堡方程：
- 通过有限差分法数值求解该方程，给定初始条件 u(x, 0) = f(x) 和边界条件 u(0, t) = u_0，u(L,t) = u_L。
数值方法：
- 使用有限差分法对方程进行离散化，采用拉克斯-温德罗夫（Lax-Wendroff）方法来处理非线性项。
- 使用显式差分法进行空间和时间离散化。
项目功能：
- 输入问题的参数，包括空间区间、时间区间、初始条件、边界条件等。
- 使用有限差分法求解汉堡方程，并使用拉克斯-温德罗夫方法进行非线性项的稳定处理。
- 输出解的数值结果。

拉克斯-温德罗夫方法（Lax-Wendroff Method）

拉克斯-温德罗夫方法是一种高精度的二阶差分方法，用于求解一类具有对流项的偏微分方程。在拉克斯-温德罗夫方法中，时间和空间的离散化通过二阶精度的差分格式完成，通常用于处理非线性方程。

拉克斯-温德罗夫方法的离散化形式如下：

其中，u_i^n 是在第 n 时间步和第 i 空间节点的解，f(u) = u^2 / 2 是对流方程的非线性项。

C++代码实现

以下是 C++ 实现的有限差分法求解非粘性时变汉堡方程，采用拉克斯-温德罗夫方法进行稳定处理。

完整代码

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

class BurgersSolver {
private:
    int Nx, Nt;            // 空间和时间步数
    double L, T;            // 空间区间和时间区间
    double u0, uL;          // 边界条件
    double (*f)(double);    // 初始条件函数
    double dx, dt;          // 空间和时间步长
    vector<double> u_old, u_current, u_new;  // 解向量

    // 初始条件 f(x)
    double f_initial(double x) {
        return sin(M_PI * x / L);  // 假设初始条件为正弦函数
    }

    // 边界条件处理
    void applyBoundaryConditions() {
        u_current[0] = u0;          // 左边界条件
        u_current[Nx - 1] = uL;     // 右边界条件
        u_old[0] = u0;              // 左边界条件
        u_old[Nx - 1] = uL;         // 右边界条件
    }

public:
    // 构造函数，初始化变量
    BurgersSolver(double L, double T, int Nx, int Nt, double u0, double uL, double (*f)(double))
        : L(L), T(T), Nx(Nx), Nt(Nt), u0(u0), uL(uL), f(f) {
        dx = L / (Nx - 1);
        dt = T / (Nt - 1);

        u_old.resize(Nx);
        u_current.resize(Nx);
        u_new.resize(Nx);
    }

    // 使用拉克斯-温德罗夫方法求解非粘性时变汉堡方程
    void solve() {
        // 初始化 u(x,0)
        for (int i = 0; i < Nx; ++i) {
            double x = i * dx;
            u_old[i] = f_initial(x);  // 初始条件
            u_current[i] = f_initial(x);  // 初始条件
        }

        // 应用边界条件
        applyBoundaryConditions();

        // 时间步进
        for (int n = 0; n < Nt - 1; ++n) {
            // 使用拉克斯-温德罗夫方法更新解
            for (int i = 1; i < Nx - 1; ++i) {
                u_new[i] = u_current[i] - (dt / (2 * dx)) * (u_current[i+1]*u_current[i+1] - u_current[i-1]*u_current[i-1])
                           + (dt * dt / (2 * dx * dx)) * (u_current[i+1] - 2 * u_current[i] + u_current[i-1]);
            }

            // 应用边界条件
            applyBoundaryConditions();

            // 更新解
            u_old = u_current;
            u_current = u_new;

            // 输出当前时间步的解
            if (n % 50 == 0) {
                cout << "时间步 " << n << " 的温度分布：" << endl;
                for (int i = 0; i < Nx; ++i) {
                    cout << setw(10) << u_current[i];
                }
                cout << endl;
            }
        }
    }
};

// 初始条件函数 f(x)
double f(double x) {
    return sin(M_PI * x);  // 假设初始条件为正弦函数
}

int main() {
    double L = 1.0;       // 空间区间长度
    double T = 0.5;       // 时间区间长度
    int Nx = 100;         // 空间网格数
    int Nt = 200;         // 时间步数
    double u0 = 0.0;      // 左边界条件
    double uL = 0.0;      // 右边界条件

    // 创建求解器对象
    BurgersSolver solver(L, T, Nx, Nt, u0, uL, f);

    // 求解非粘性时变汉堡方程
    solver.solve();

    return 0;
}

代码解读

BurgersSolver 类：
- 构造函数：初始化空间长度 L、时间总长度 T、空间步数 Nx、时间步数 Nt、边界条件 u_0、u_L，以及初始条件函数f(x)。
- solve()：使用拉克斯-温德罗夫方法进行时间推进，依次计算每个时间步的解，并更新解向量。每隔 50 个时间步输出一次当前的温度分布。
初始条件函数：
- f_initial(x)：假设初始条件为正弦波函数，用户可以根据需求修改。
边界条件处理：
- 在每个时间步中，左边界 u(0,t) = u_0 和右边界 Lu(L,t) = u_L 被设置为常数。
拉克斯-温德罗夫方法：
- 拉克斯-温德罗夫方法通过组合一阶和二阶差分来进行空间和时间的离散化，提升了数值解的精度。