3
$\begingroup$

I am having difficulties applying the boundary condition

$$\frac{\partial \vec{V}}{\partial t} + u\frac{\partial \vec{V}}{\partial x} = \frac{1}{\operatorname{Re}}\frac{\partial ^2 \vec{V}}{\partial y^2},$$

where $\vec{V}=(u,v)$ to the momentum equation at the exit (right side of the domain) of the Incompressible Navier-Stokes Equation. The problem lies in the diffusion and advection operators.

Consider incompressible NSE in non-dimensional form on square $1\times 1$ domain

\begin{equation}\label{eqn:nse-matrix} \begin{bmatrix} \mathbf{I} && 0 \\ 0 && 0 \end{bmatrix} \frac{\partial }{\partial t} \begin{pmatrix} \vec{V} \\ P \end{pmatrix} = \begin{bmatrix} \hat{L} && - \hat{G} \\ -\hat{D} && 0 \end{bmatrix} \begin{bmatrix} \vec{V} \\ P \end{bmatrix} + \begin{pmatrix} -\mathbf{\hat{H}}(\vec{V})\\ 0 \end{pmatrix} + \text{bc}_{\vec{V},P}, \end{equation} where letters stand for operators: $L$ - Laplacian, $G$ - gradient, $D$ - divergence, $H$ - non-linear advective terms.

Apply the following discretization schemes (trying to keep them second order in time except pressure):

  1. Viscous - Implicit trapezoidal - Crank Nicholson scheme (second-order method in time)
    \begin{equation}\label{eqn:viscous-crank-nicholson} \hat{L}\vec{V}=\frac{1}{2}\left(\hat{L}\vec{V}^{n+1}+\hat{L}\vec{V}^n\right). \end{equation}

  2. Nonlinear - Explicit Adams-Bashforth (second-order method in time) \begin{equation}\label{eqn:nonlinear-adams-bashforth} \mathbf{\hat{H}}(\vec{V}) = \frac{3}{2}\mathbf{\hat{H}}(\vec{V}^{n}) - \frac{1}{2}\mathbf{\hat{H}}(\vec{V})^{n-1}. \end{equation}

  3. Pressure - Implicit Euler \begin{equation}\label{eqn:pressure-implicit-euler} \hat{G}P = \hat{G}P^{n+1}. \end{equation}

To make things easier, consider an issue for Laplacian first. Let the system be 1-dimensional on $M$ equidistant nodes on which we evaluate the equation. The boundary condition then simplifies to

$$\frac{\partial u}{\partial t} + u\frac{\partial u}{\partial x} = 0.$$

Let us write one equation at the right-most point $M$ using the second-order schemes in space. Our flow is allowed to go both directions and, hence, I think that the one-sided second derivative will not be applicable, that's why we are restricted to using central differences

$$\frac{\partial u}{\partial t} -\frac{1}{2}\frac{u^{n+1}_{M+1}-2u^{n+1}_{M}+u^{n+1}_{M-1}}{\Delta x^2}+\hat{G}P^{n+1}=\frac{1}{2}\frac{u^{n}_{M+1}-2u^{n}_{M}+u^{n}_{M-1}}{\Delta x^2}-\left(\frac{3}{2}\mathbf{\hat{H}}(\vec{V}^{n}) - \frac{1}{2}\mathbf{\hat{H}}(\vec{V})^{n-1} \right).$$

The problem here is that the velocity components $u_{M+1}^{n+1},u_{M+1}^{n}$ lie outside the domain.

The same issue goes to the advective terms, I am using $$\frac{\partial uu}{\partial x}=\frac{\frac{u_{M+1}+u_{M}}{2}-\frac{u_{M}+u_{M-1}}{2}}{\Delta x}$$ in conservative form for approximation (did not include in the 1D case because the conservative form only works together with continuity equation).

I am also unsure how to properly apply the $v$ component boundary condition, in 2D case. A staggered grid is considered, hence, velocity components lie at the centres of each interval and, therefore, I do not think it is appropriate to directly apply the boundary condition at the corresponding rightmost point of $v$, which lies not at $x=1$, but at $x=1-\Delta x/2$ coordinate.

I took a look at the corresponding papers (SciHub alert!) One and Two but they do not provide the details, unfortunately. They only state time discretization \begin{gathered} \frac{\mathbf{V}^{n+1}-\mathbf{V}^n}{\Delta t}+\frac{u^n}{2}\left[\left(\frac{\partial \mathbf{V}}{\partial x}\right)^{n+1}+\left(\frac{\partial \mathbf{V}}{\partial x}\right)^n\right] \\ =v\left(\frac{\partial^2 \mathbf{V}}{\partial y^2}\right)^n+O\left(\Delta t^2\right). \end{gathered}

The only solution I came up with so far is expressing outer points $u_{M+1}^{n+1},u_{M+1}^{n}$ in terms of $u_{M}^{n+1},u_{M}^{n},u_{M-1}^{n+1},u_{M-1}^{n}$ using linear extrapolation. But this only works for $x$-momentum. And, moreover, assuming linear interpolation makes whole diffusion term 0…

Any help is appreciated! Thanks for reading!

Edit: My next idea turned out to be wrong as well. I thought one could keep grid unrefined in direction perpendicular to the Bc, and use three point one sided diffusion scheme. It broke the symmetry of diffusion matrix.

$\endgroup$
5
  • 1
    $\begingroup$ This is not any kind of physical boundary condition. What is it supposed to represent? $\endgroup$ Dec 14, 2023 at 6:11
  • $\begingroup$ @WolfgangBangerth hello, this is a version (derivation on page 242 of sci-hub.se/https://doi.org/10.1006/jcph.1993.1140) of absorbing BC of a very popular paper, comes from assuming wave equation is true at the outlet and allowing waves to propagate only outwards, killing any waves that reflect back into the domain. Very popular paper by Engquist and Majda math.mcgill.ca/~gantumur/docs/down/Engquist77.pdf original equations 1.24 1.25. Based on papers this removes inward reflected normal waves only, but is still sufficient looking at the results. $\endgroup$
    – 2Napasa
    Dec 14, 2023 at 6:40
  • $\begingroup$ I think there might be an insight inside of my silly head. Wrote it on paper, will typeset tomorrow to accept criticism :( $\endgroup$
    – 2Napasa
    Dec 14, 2023 at 7:12
  • 1
    $\begingroup$ You are supposed to have outflow at right boundary. So you will discretize the $V$ equation in the last layer of cells adjacent to your right boundary; then for $u V_x$ with upwind scheme only needs points to the left which are available, and for $V_{yy}$ also you have all required stencil values. $\endgroup$
    – cfdlab
    Dec 14, 2023 at 14:02
  • $\begingroup$ @cfdlab Thanks! I tried conditions like these. Upwind isn’t an option at the outlet, unfortunately, my flow is allowed to go both directions. Hence, will need a point outside the domain. That’s why I switched to explicit Adam’s Bashforth scheme for advection, for which I can exe extrapolation to approximate ghost points outside. And try applying this non-standard OpenBC from paper. $\endgroup$
    – 2Napasa
    Dec 14, 2023 at 17:23

1 Answer 1

3
$\begingroup$

I managed to resolve the issue,

Need to compute the values only for the inner part of the domain, then apply BC for the boundary values. The only issue is that there is a backward scheme at the boundary.

Below I will use half indices for the values at nodes, integer indices for the centre interval coordinate.

The last unknown element of the system in $x$ direction is $u^{n+1}_{M-\frac{1}{2},j}$. To apply non uniform Laplacian we need to express $u^{n+1}_{M+\frac{1}{2},j}$ using BC. \begin{equation*} \begin{gathered} \frac{{u}^{n+1}_{M+\frac{1}{2},j}-{u}^n_{M+\frac{1}{2},j}}{\Delta t} +\frac{u^n_{M+\frac{1}{2},j}}{2} \left[ \left( \frac{u^{n+1}_{M+\frac{1}{2},j}-u^{n+1}_{M-\frac{1}{2},j}}{\Delta x_M} \right) +\left( \frac{u^{n}_{M+\frac{1}{2},j}-u^{n}_{M-\frac{1}{2},j}}{\Delta x_M} \right) \right] \\ =\nu\left(\frac{1}{h_c h_s}u^{n}_{M+\frac{1}{2},j-1} - \frac{2}{h_n h_s}u^{n}_{M+\frac{1}{2},j} + \frac{1}{h_c h_n}u^n_{M+\frac{1}{2},j+1}\right) \end{gathered} \end{equation*} leads to \begin{equation} \begin{gathered} {{u}^{n+1}_{M+\frac{1}{2},j}}\left(\frac{1}{\Delta t} + \frac{u^n_{M+\frac{1}{2},j}}{2\Delta x_M}\right) = u^{n+1}_{M-\frac{1}{2},j}\frac{u^n_{M+\frac{1}{2},j}}{2\Delta x_M}\\ -\frac{u^n_{M+\frac{1}{2},j}}{2}\frac{u^{n}_{M+\frac{1}{2},j}-u^{n}_{M-\frac{1}{2},j}}{\Delta x_M} +\frac{{u}^n_{M+\frac{1}{2},j}}{\Delta t} +\nu\left(\frac{1}{h_c h_s}u^{n}_{M+\frac{1}{2},j-1} - \frac{2}{h_n h_s}u^{n}_{M+\frac{1}{2},j} + \frac{1}{h_c h_n}u^n_{M+\frac{1}{2},j+1}\right). \end{gathered} \end{equation} Note that $\left(\frac{1}{\Delta t} + \frac{u^n_{M+\frac{1}{2},j}}{2\Delta x_M}\right)$ is never zero regardless of the flow direction, hence dividing by it is allowed. The second line of fully discretized equation above is treated explicitly since we posses the flow data at time steps $n$ and below. Denote explicit terms as \begin{equation*} bc_e^n=-\frac{u^n_{M+\frac{1}{2},j}}{2}\frac{u^{n}_{M+\frac{1}{2},j}-u^{n}_{M-\frac{1}{2},j}}{\Delta x_M} +\frac{{u}^n_{M+\frac{1}{2},j}}{\Delta t} +\nu\left(\frac{1}{h_c h_s}u^{n}_{M+\frac{1}{2},j-1} - \frac{2}{h_n h_s}u^{n}_{M+\frac{1}{2},j} + \frac{1}{h_c h_n}u^n_{M+\frac{1}{2},j+1}\right), \end{equation*} which simplifies to \begin{equation} \begin{gathered} {{u}^{n+1}_{M+\frac{1}{2},j}}\left(\frac{1}{\Delta t} + \frac{u^n_{M+\frac{1}{2},j}}{2\Delta x_M}\right) = u^{n+1}_{M-\frac{1}{2},j}\frac{u^n_{M+\frac{1}{2},j}}{2\Delta x_M}+bc_e^n. \end{gathered} \end{equation}

Next, we plug above BC into Laplacian, which finally leads us to \begin{equation} \begin{aligned} \left.\frac{\partial ^2 u}{\partial x^2}\right|_{M-\frac{1}{2},j}^{n+1} &=\frac{1}{h_c h_e}u^{n+1}_{M+\frac{1}{2},j} - \frac{2}{h_w h_e}u^{n+1}_{M-\frac{1}{2},j} + \frac{1}{h_c h_w}u^{n+1}_{M-\frac{1}{2}-1,j}\\ &=\frac{1}{h_c h_e}\left(u^{n+1}_{M-\frac{1}{2},j}\frac{u^n_{M+\frac{1}{2},j}}{2h_e}+bc^n_e\right)- \frac{2}{h_w h_e}u^{n+1}_{M-\frac{1}{2},j} + \frac{1}{h_c h_w}u^{n+1}_{M-\frac{1}{2}-1,j}\\ &=\left[ \frac{1}{h_c h_e}\left(\frac{u^n_{M+\frac{1}{2},j}}{2h_e} \right) -\frac{2}{h_w h_e}\right]u^{n+1}_{M-\frac{1}{2},j} + \frac{1}{h_c h_w}u^{n+1}_{M-\frac{1}{2}-1,j} + \frac{bc^n_e}{h_c h_e}. \end{aligned} \end{equation}

$\endgroup$
1
  • $\begingroup$ I checked for $v$ velocity component, it works, uses central difference scheme with ghost velocity outside. I interpolate the boundary value of velocity. $\endgroup$
    – 2Napasa
    Dec 29, 2023 at 4:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.