# Planar Poiseuille Flow Tutorial This tutorial documents the analytical solution implemented for incompressible planar Poiseuille flow with either constant-viscosity or power-law rheology. The solution gives the streamwise velocity and temperature profiles across a planar channel. The implementation assumes the exact pressure is zero because periodic boundary conditions are used in the streamwise direction. --- ## 1. Problem setup Consider fully developed planar Poiseuille flow between two parallel plates. The cross-channel coordinate is $y$, and the channel walls are located at $$ y = -H, \qquad y = H, $$ where $$ H = \frac{h_{\max} - h_{\min}}{2}. $$ The flow is driven by a uniform streamwise momentum source $S_u$. The velocity field is $$ \mathbf{u}(y) = \begin{bmatrix} u(y) \\ 0 \end{bmatrix} $$ in two dimensions, or $$ \mathbf{u}(y) = \begin{bmatrix} u(y) \\ 0 \\ 0 \end{bmatrix} $$ in three dimensions. The wall temperatures are prescribed as $$ T(-H) = T_l, \qquad T(H) = T_u. $$ The exact Lagrange pressure is $$ p(y) = 0. $$ --- ## 2. Parameters | Symbol | Code parameter | Description | |---|---|---| | $h_{\min}$ | `H min` | Lower geometric bound used to define the channel height | | $h_{\max}$ | `H max` | Upper geometric bound used to define the channel height | | $H$ | computed | Channel half-height, $H=(h_{\max}-h_{\min})/2$ | | $T_l$ | `Lower wall temperature` | Lower-wall temperature | | $T_u$ | `Upper wall temperature` | Upper-wall temperature | | $S_u$ | `Momentum source` | Streamwise momentum source | | $\nu$ | `Kinematic viscosity` or `Flow consistency index` | Constant kinematic viscosity or power-law consistency index | | $n$ | `Flow behavior index` | Power-law index; $n=1$ for a Newtonian fluid | | $k$ | `Thermal conductivity` | Thermal conductivity | For a constant-property Newtonian fluid, the implementation sets $$ n = 1. $$ For a power-law fluid, $n$ is read from the input as the flow behavior index. --- ## 3. Analytical velocity solution The streamwise velocity profile is $$ u(y) = \frac{n}{n+1} \left(\frac{S_u}{\nu}\right)^{1/n} H^{(n+1)/n} \left[ 1 - \left(\frac{|y|}{H}\right)^{(n+1)/n} \right]. $$ The transverse velocity components are zero: $$ v(y) = 0, $$ and, in three dimensions, $$ w(y) = 0. $$ ### Newtonian limit For $n=1$, the velocity reduces to the standard parabolic Poiseuille profile: $$ u(y) = \frac{1}{2} \frac{S_u}{\nu} \left( H^2 - y^2 \right). $$ --- ## 4. Analytical temperature solution The analytical temperature profile is $$ T(y) = \frac{n^2}{(3n+1)(2n+1)} \frac{\nu}{k} \left(\frac{S_u}{\nu}\right)^{(n+1)/n} \left[ H^{(3n+1)/n} - |y|^{(3n+1)/n} \right] \\ + \frac{1}{2} \frac{y}{H} \left(T_u - T_l\right) + \frac{1}{2} \left(T_u + T_l\right). $$ This expression contains two contributions: 1. A viscous-heating contribution caused by the imposed momentum source. 2. A linear conductive contribution satisfying the imposed wall temperatures. The boundary values are recovered directly: $$ T(-H) = T_l, \qquad T(H) = T_u. $$ ### Newtonian limit For $n=1$, $$ T(y) = \frac{1}{20} \frac{\nu}{k} \left(\frac{S_u}{\nu}\right)^2 \left( H^4 - y^4 \right) + \frac{1}{2} \frac{y}{H} \left(T_u - T_l\right) + \frac{1}{2} \left(T_u + T_l\right). $$ --- ## 5. Python script for plotting the analytical solution Save the following as `plot_planar_poiseuille_analytical_solution.py` and run ```bash python plot_planar_poiseuille_analytical_solution.py ``` The script generates velocity and temperature plots and saves them as PNG files. ```python import numpy as np import matplotlib.pyplot as plt def planar_poiseuille_solution( y, h_min=-1.0, h_max=1.0, T_l=300.0, T_u=310.0, S_u=1.0, nu=1.0, k=1.0, n=1.0, ): """ Analytical solution for planar Poiseuille flow with viscous heating. Parameters ---------- y : array_like Cross-channel coordinate. h_min, h_max : float Geometric bounds used to define the channel half-height, H = (h_max - h_min) / 2. T_l, T_u : float Lower- and upper-wall temperatures. S_u : float Streamwise momentum source. nu : float Kinematic viscosity for n = 1, or power-law consistency index for n != 1. k : float Thermal conductivity. n : float Flow behavior index. Use n = 1 for a Newtonian fluid. Returns ------- u : ndarray Streamwise velocity. T : ndarray Temperature. p : ndarray Exact Lagrange pressure, equal to zero. """ y = np.asarray(y, dtype=float) H = 0.5 * (h_max - h_min) if H <= 0.0: raise ValueError("The channel half-height H must be positive.") if nu <= 0.0: raise ValueError("nu must be positive.") if k <= 0.0: raise ValueError("k must be positive.") if n <= 0.0: raise ValueError("n must be positive.") if np.any(np.abs(y) > H): raise ValueError("All y values must satisfy |y| <= H.") abs_y = np.abs(y) u = ( n / (n + 1.0) * (S_u / nu) ** (1.0 / n) * H ** ((n + 1.0) / n) * (1.0 - (abs_y / H) ** ((n + 1.0) / n)) ) T = ( n**2 / ((3.0 * n + 1.0) * (2.0 * n + 1.0)) * nu / k * (S_u / nu) ** ((n + 1.0) / n) * (H ** ((3.0 * n + 1.0) / n) - abs_y ** ((3.0 * n + 1.0) / n)) ) T += 0.5 * y * (T_u - T_l) / H + 0.5 * (T_u + T_l) p = np.zeros_like(y) return u, T, p def main(): # Example parameters. Modify these to match your input deck. h_min = -1.0 h_max = 1.0 T_l = 300.0 T_u = 310.0 S_u = 1.0 nu = 1.0 k = 1.0 n = 1.0 H = 0.5 * (h_max - h_min) y = np.linspace(-H, H, 500) u, T, p = planar_poiseuille_solution( y, h_min=h_min, h_max=h_max, T_l=T_l, T_u=T_u, S_u=S_u, nu=nu, k=k, n=n, ) plt.figure() plt.plot(u, y) plt.xlabel("Streamwise velocity, u(y)") plt.ylabel("Cross-channel coordinate, y") plt.title("Analytical velocity profile") plt.grid(True) plt.tight_layout() plt.savefig("analytical_velocity_profile.png", dpi=300) plt.figure() plt.plot(T, y) plt.xlabel("Temperature, T(y)") plt.ylabel("Cross-channel coordinate, y") plt.title("Analytical temperature profile") plt.grid(True) plt.tight_layout() plt.savefig("analytical_temperature_profile.png", dpi=300) plt.show() if __name__ == "__main__": main() ``` --- ## 6. Notes for verification Useful checks are: $$ u(-H) = u(H) = 0, $$ $$ T(-H) = T_l, \qquad T(H) = T_u, $$ $$ p(y) = 0. $$ For $n=1$, the velocity profile should be parabolic and the viscous-heating contribution to temperature should be quartic in $y$. --- ## 7. Running this case with VERTEX-CFD on the Helios CPU platform The VERTEX-CFD input file for this analytical solution is located at ```bash vertex-cfd/examples/inputs/incompressible/incompressible_2d_planar_poiseuille.xml ``` These instructions assume that VERTEX-CFD has already been built on the Helios CPU platform and that the Helios setup environment is available at ```bash /software/build/vertex-cfd/setup-env.sh ``` Runs may be performed from your home directory. ### 7.1. Stage the input file Create a run directory and copy the planar Poiseuille input file into it: ```bash mkdir -p ~/planar_poiseuille cd ~/planar_poiseuille cp /path/to/vertex-cfd/examples/inputs/incompressible/incompressible_2d_planar_poiseuille.xml . ``` If you are working from a clone of the VERTEX-CFD repository, replace `/path/to/vertex-cfd` with the path to your source tree. The copied input file should be named ```bash incompressible_2d_planar_poiseuille.xml ``` The submission script below sets ```bash BASEFILE=incompressible_2d_planar_poiseuille ``` so the input file is resolved as ```bash INPUTFILE=${BASEFILE}.xml ``` ### 7.2. Create a SLURM submission script Create a file named ```bash touch run-planar-poiseuille-on-helios-cpu.sh ``` and copy the following contents into it: ```bash #!/bin/bash #SBATCH -N 1 #SBATCH --ntasks-per-node 2 #SBATCH --cpus-per-task 1 #SBATCH --output output%j.txt #SBATCH --error error%j.txt #SBATCH -J planar-poiseuille #SBATCH -t 1-00:00 #SBATCH --mail-type=ALL #SBATCH --mail-user=YOUR-EMAIL@ornl.gov # Load the VERTEX-CFD Helios CPU environment. source /software/build/vertex-cfd/setup-env.sh module list export OMP_PROC_BIND=true export OMP_PLACES=threads # Input file for this verification case. # This assumes the job is submitted from the directory containing the XML file. BASEFILE=incompressible_2d_planar_poiseuille INPUTFILE=${BASEFILE}.xml # Update this path to point to your VERTEX-CFD executable. # For example, this may be in your build or install directory. VERTEXCFD=/PATH/TO/BUILD/DIRECTORY/src/vertexcfd mpirun $VERTEXCFD --i=$INPUTFILE # Recombine decomposed Exodus output files. # Each MPI rank writes one file during the simulation. EPU=/software/build/vertex-cfd/trilinos/16.0.0/bin/epu $EPU -auto ${BASEFILE}_solution.exo.${SLURM_NTASKS}.000 # Remove decomposed files after successful recombination. rm ${BASEFILE}_solution.exo.${SLURM_NTASKS}.* ``` Before submitting, edit the following entries as needed: - `#SBATCH -N 2`: number of nodes. - `#SBATCH --ntasks-per-node 185`: MPI ranks per node. - `#SBATCH --mail-user=YOUR-EMAIL@ornl.gov`: your email address. - `VERTEXCFD=/PATH/TO/BUILD/DIRECTORY/src/vertexcfd`: path to your VERTEX-CFD executable. For a small verification run, you may want to reduce the node count and MPI rank count to match the size of the mesh and the resources available on Helios. ### 7.3. Submit the job Submit the job using ```bash sbatch run-planar-poiseuille-on-helios-cpu.sh ``` During the run, VERTEX-CFD writes decomposed Exodus solution files, one per MPI rank. After the solver finishes, the script recombines them using Trilinos `epu`: ```bash /software/build/vertex-cfd/trilinos/16.0.0/bin/epu -auto incompressible_2d_planar_poiseuille_solution.exo.${SLURM_NTASKS}.000 ``` The script then removes the decomposed per-rank files: ```bash rm incompressible_2d_planar_poiseuille_solution.exo.${SLURM_NTASKS}.* ``` If the job exits before recombination, you can manually rerun the `epu` command after checking how many decomposed files were written. ### 7.4. Expected output For this case, the main recombined solution file should have a name similar to ```bash incompressible_2d_planar_poiseuille_solution.exo ``` depending on the output settings in the XML input file. You can inspect the solution in ParaView, for example on a Helios visualization node, and compare the numerical fields against the analytical profiles documented above: $$ u(y) = \frac{n}{n+1} \left(\frac{S_u}{\nu}\right)^{1/n} H^{(n+1)/n} \left[ 1 - \left(\frac{|y|}{H}\right)^{(n+1)/n} \right], $$ $$ T(y) = \frac{n^2}{(3n+1)(2n+1)} \frac{\nu}{k} \left(\frac{S_u}{\nu}\right)^{(n+1)/n} \left[ H^{(3n+1)/n} - |y|^{(3n+1)/n} \right] \\ + \frac{1}{2} \frac{y}{H} \left(T_u - T_l\right) + \frac{1}{2} \left(T_u + T_l\right). $$ The exact Lagrange pressure for this periodic setup is $$ p(y)=0. $$