constrainPressure主要是更新压力第二类边界条件，公式如下：
\begin{equation}
\left( \nabla p_{rgh} \right)_f \cdot\bfn_f=
\frac{\left(\mathbf{HbyA}_f^{*}- \frac{1}{{{A^n_{\mathrm{P},f}}}}(\bfg\cdot\bfh\nabla\rho)_f - \mathbf{U}_f \right)\cdot\bfS_f}
{
|\bfS_f|
\frac{1}{{{A^n_{\mathrm{P},f}}}}
}
\end{equation}
上述公式与代码并不一致。在OpenFOAM中，constrainPressure为
多乘了一个密度。应该改为：
forAll(pBf, patchi) { if (isA<fixedFluxPressureFvPatchScalarField>(pBf[patchi])) { refCast<fixedFluxPressureFvPatchScalarField> ( pBf[patchi] ).updateCoeffs ( ( phiHbyABf[patchi] - MRF.relative(SfBf[patchi] & UBf[patchi], patchi) ) /(magSfBf[patchi]*rhorAUBf[patchi]/rho.boundaryField()[patchi]) ); } }在非常老的OpenFOAM版本中，看起来是正确的，与公式一致
setSnGrad<fixedFluxPressureFvPatchScalarField> ( p_rgh.boundaryField(), ( phiHbyA.boundaryField() - fvOptions.relative(mesh.Sf().boundaryField() & U.boundaryField()) )/(mesh.magSf().boundaryField()*rAUf.boundaryField()) );