openfoam 不可压缩计算中的能量守恒
-
各位老师好，
我最近在用openfoam计算一个简单的强迫振动问题，想要验证能量守恒的关系(Loose Leaf for Modern Compressible Flow With Historical Perspective (John Anderson) )。
因为我是计算不可压缩的问题，所以能量方程有一定简化。这里我等号左侧第一二项为0（没有热量变化和旋转），左侧最后一项也为0（没有外力）
模型如下:
后面通过输出每个时刻的速度压力场处理公式里的各项，但是结果相差很大，可以说完全不守恒。我以为不可压缩问题天然满足能量守恒公式，结果和算例贴在下面。
c2case.7z这个结果相差一个数量级
我分析有两个原因，一个是选的能量公式 有误，在openfoam里输出每项的时候有误，所以在comsol里验证了一次，结果和openfoam的相差不大
（后面整理了再贴出来）。如何计算质量守恒和能量守恒 | COMSOL 博客感觉很奇怪，想请教各位老师。
-
我是输出各个边界上的速度压力场自己积分的这一项。代码里是B3,按照这样计算是有非零值的，不知道哪里出了问题。```
// 计算速度场的梯度 scalar rho = 1000; // 流体密度 const scalar mu = 1e-3; // 粘性系数 tmp<volTensorField> gradUTmp = fvc::grad(U); const volTensorField& gradU = gradUTmp(); // 初始化总的B3和剪切通量 scalar totalB3 = 0.0; scalar totalViscousFlux = 0.0; // 遍历所有边界patch forAll(mesh.boundaryMesh(), patchI) { const polyPatch& patch = mesh.boundaryMesh()[patchI]; word patchName = patch.name(); // 如果对特定的patch感兴趣，可以添加条件检查 if ((patchName == "front")||(patchName == "back")) { continue; // 跳过名为"front"或"back"的边界 } // 初始化当前patch的B3和剪切通量 scalar B3Patch = 0.0; scalar ViscousFluxPatch = 0.0; // 循环遍历patch上的所有面 forAll(patch, faceId) { const vector& Sf = patch.faceNormals()[faceId]; // 从patch对象获取归一化法向向量 scalar faceArea = mag(patch.faceAreas()[faceId]); // 计算面积向量的模长作为实际面积 const vector& U_boundary = U.boundaryField()[patchI][faceId]; const Tensor<double>& faceGradU = gradU.boundaryField()[patchI][faceId]; // 计算应力张量（包含法向和切向应力） Tensor<double> stressTensor = mu * (faceGradU + faceGradU.T()); // 计算壁面剪切应力（只保留切向分量） vector wallShearStress = stressTensor & Sf; wallShearStress -= (wallShearStress & Sf) * Sf / (magSqr(Sf) + VSMALL); // 壁面剪切应力与速度的点积 scalar wallShearStressDotU = wallShearStress & U_boundary; const volScalarField& pressure = mesh.lookupObject<volScalarField>("p"); scalar p_boundary = pressure.boundaryField()[patchI][faceId]; scalar pressureFlux = p_boundary * (U_boundary & Sf); // 根据面积加权的壁面剪切应力与速度的点积以及压力与速度的通量 ViscousFluxPatch += faceArea * wallShearStressDotU; B3Patch += faceArea * pressureFlux; } // 输出当前边界上的B3和剪切通量 Info << "Patch: " << patchName << ", Integrated B3: " << B3Patch << ", Integrated Viscous Flux: " << ViscousFluxPatch << endl; // 将当前patch的B3和剪切通量累加到总和中 totalB3 += B3Patch; totalViscousFlux += ViscousFluxPatch; } // 输出总的B3和剪切通量 Info << "Total Integrated B3 on all patches: " << totalB3 << endl; Info << "Total Integrated Viscous Flux on all patches: " << totalViscousFlux << endl;
但是案例结果并不满足。
