Slit Validation Runs for SPH_CCA Code
This page presents a series of SPH simulations of flow between parallel plates used to test and validate the SPH_CCA code. The series of simulations was performed within the SALSSA workflow environment, and run on sixteen processors on the EMSL supercomputer Chinook.
Parameters that are held constant for all simulations shown below:
- Particle spacing = 0.333 (Density = 3x3x3=27.0)
- Particle-Particle repulsion: 0.00001
- Interaction cutoff (h): 1.0
- Padding on cutoff: 0.3
The image below shows the series of simulations presented in a graphical format created by the SALSSA workflow environment. Superimposed on the graph are several colored ellipses that highlight important groups of simulation runs. From the bottom up, these groups are (generally in the order that they were performed):
- Green: Initial run of example problem set up by Bruce Palmer (SPH_CCA code author) - to learn how to execute the code and post-process output within SALSSA
- Pink: Runs of slit problem to find combination of parameters (particularly time step and speed of sound) that gives a stable solution
- Purple: Once a basic stable solution was found, a preliminary sequence of runs spanning the parameter space of interest (varying gravity and viscosity simultaneously while maintaining the Reynolds number at 1.0) was performed to determine the necessary time step for stability across the parameter range.
- Blue: This is the main body of simulation runs, performed to evaluate simulation performance for low Reynolds number (Re=1.0) with variable gravity and viscosity.
- Yellow: Based on a stable solution selected from the Blue group, a series of runs was made where only the speed of sound parameter (c) was varied.
- Orange: A few runs were made to evalute the effect of changing the aperture width
Blue Set: Constant Re=1.0; Vary gravity and viscosity
A series of thirteen simulation runs was performed in this suite to examine the impact of varying gravity (body force) and viscosity simultaneously in such a manner so as to maintain a constant Reynolds number of Re=1.0. For all these runs, the speed of sound was set to a value of c=1.0, and the aperture width was set to 20. Starting with the base case of viscosity = 15.5885 and gravity = 0.0005, the viscosity was increased by 20% in each subsequent run, and gravity was increased so as to maintain Re=1.0. The Reynolds number is calculated based on the formula: Re = (rho^2)*g*(L^3)/(12*mu^2) where mu is the viscosity and g is gravity force, rho is density (27.0) and L is aperture width (20). Therefore, for constant Re, g varies in proportion to the viscosity squared.
Selected results are summarized below. The complete output of these simulations is provided in Excel spreadsheet format here (Large File: 75 Mbytes). Note that as the viscosity and gravity were increased, the time step had to be decreased to maintain solution stability (from an initial value of 0.02 to final value of 0.005).
Base Case: Simulation results are summarized below (click on links to see plots):
- Mean and max. velocity over time: Average velocity converges to an assymptotic value of 0.0302. Compared with the analytical solution of V(avg)=0.0289, this is 4.6% too high.
- Min. and max. density over time: Density fluctuations are reasonably small and decrease as the simulation progresses.
- Velocity profile at final time: Simulated velocities are slightly larger than predicted by the analytical solution, and also appear to be somewhat "flattened" near the center of the profile.
As the viscosity and gravity were increased, the simulated velocity profile became more parabolic in nature (less "flattening" in the center), but the mean velocity remained too high relative to the analytical solution. Results from simulation runs #4 and #12 are provided here as representative examples:
Simulation #4: Viscosity = 32.324; gravity = 0.00215. Simulation results are summarized below (click on links to see plots):
- Mean and max. velocity over time: Average velocity converges to an assymptotic value of 0.0632. Compared with the analytical solution of V(avg)=0.0599, this is 5.6% too high.
- Min. and max. density over time: Density fluctuations are reasonably small and decrease as the simulation progresses.
- Velocity profile at final time: Simulated velocities are slightly larger than predicted by the analytical solution but show a nice parabolic profile shape, perhaps with a bit of a "shoulder" at x=5 and x=15.
Simulation #12: Viscosity = 138.99; gravity = 0.03975. Simulation results are summarized below (click on links to see plots):
- Mean and max. velocity over time: Average velocity converges to an assymptotic value of 0.2719. Compared with the analytical solution of V(avg)=0.2574, this is 5.7% too high.
- Min. and max. density over time: Density fluctuations are reasonably small and decrease as the simulation progresses.
- Velocity profile at final time: Simulated velocities are slightly larger than predicted by the analytical solution but show a nice parabolic profile shape, perhaps with a bit of a "shoulder" at x=5 and x=15.
To look at the long-term stability of the solutions, Simulation #4 was run out for a longer time period. Results are shown below:
Simulation #4 (long time): Simulation results are summarized below (click on links to see plots):
- Mean and max. velocity over time: Average velocity converges to an assymptotic value of 0.0633. Compared with the analytical solution of V(avg)=0.0599, this is 5.7% too high.
- Mean and max. velocity over time - expanded view of later time: There is a bit of a "blip" in the velocities around time step 1700, with waves starting around time step 1500. The "blip" corresponds to a short-lived instability in the densities around the same time (see below).
- Min. and max. density over time: Density fluctuations are reasonably small and decrease smoothly up until about 400 time steps, after which they show some instability but continue to remain fairly small overall.
- Velocity profile at final time: Simulated velocities are slightly larger than predicted by the analytical solution but show a nice parabolic profile shape, perhaps with a bit of a "shoulder" at x=5 and x=15.
Yellow Set: Constant gravity and viscosity; vary speed of sound
The speed of sound parameter (celerity or "C") controls the degree of compressibility of the simulated fluid (higher "C" is less compressible). In the set of runs above (Blue set), C=1.0 in all simulations. Here, we select simulation #4 parameters (viscosity = 32.3243 and gravity = 0.00215) and hold those constant while varying the speed of sound parameter to evaluate solution stability and impacts on solution results. In all simulations here, the analytical solution mean velocity is 0.0599 and the Reynolds number is 1.0. A time step of 0.01 was used for all simulations.
Selected results are summarized below. The complete output of these simulations is provided in Excel spreadsheet format here.
Simulation #C1: The speed of sound parameter was doubled to C=2.0. Density and velocity plots (not shown) were very similar to those of Simulation #4 above. The average velocity converged to a value of 0.0629, slightly better than that of Simulation #4 (5.2% error compared to 5.6%). However, the reason for this "improvement" is that the velocity profile was slightly flattened in the center.
Simulation #C2: The speed of sound parameter was again doubled to C=4.0. Simulation results are summarized below (click on links to see plots):
- Mean and max. velocity over time: Average velocity converges to an assymptotic value of 0.0587. Compared with the analytical solution of V(avg)=0.0599, this is 2.0% too low. Note that the convergence here has some waviness and is not smooth as in the runs with lower C values.
- Min. and max. density over time: Density fluctuations are still reasonably small but do not seem to be decreasing over time (particularly the minimum density which is nearly constant).
- Velocity profile at final time: Simulated velocity profile shows a dramatic flattening in the center.
Simulations #C3 through #C8: Since significant errors were introduced by increasing C, we subsequently performed a series of runs where C was decreased by a factor of two and then by one, two, three and four orders of magnitude. In each of these runs the results were essentially the same with the exception of the time evolution of density fluctuations. In all cases, the average velocity converged smoothly (without waviness) to a value of 0.0632, 5.6% too large relative to the analytical solution. In all cases, the velocity profiles had nice parabolic shapes (with perhaps a slight "shoulder" at x=5 and x=15), but with velocities overall too large relative to the analytical solution. The density fluctuations over time became progressively smoother, and the last three simulations (C=0.01, 0.001, and 0.0001) showed very similar density time series. Links to selected plots are shown below:
- C=0.01 Mean and max. velocity over time: Average velocity smoothly converges to an assymptotic value of 0.0632. Compared with the analytical solution of V(avg)=0.0599, this is 5.6% too high.
- C=0.01 Min. and max. density over time: Density fluctuations are small and smooth, with maximum density varying more than minimum density.
- C=0.01 Velocity profile at final time: Simulated velocities are slightly larger than predicted by the analytical solution but show a nice parabolic profile shape, perhaps with a bit of a "shoulder" at x=5 and x=15.
- C=0.001 Min. and max. density over time: Density fluctuations are small and smooth, with maximum density varying more than minimum density.
- C=0.0001 Min. and max. density over time: Density fluctuations are small and smooth, with maximum density varying more than minimum density.
Orange Set: Vary aperture width
It was hypothesized that the overprediction of velocity in the previous sets of simulations might be associated with imperfect enforcement of the no-slip boundary condition at the aperture walls. To test this hypothesis, we increased the size of the aperture from 20 to 30 units, which theoretically should decrease the impact of small wall effects.
Selected results are summarized below. The complete output of these simulations is provided in Excel spreadsheet format here.
Simulation #W1: We used simulation #4 parameters (viscosity = 32.3243, gravity = 0.00215, C=1.0). With increased aperture width, the analytical solution gives an average velocity of 0.1347. Note that the Reynolds number is no longer held at 1.0, but increased in this simulation to a value of 3.37. Simulation results are summarized below (click on links to see plots):
- Mean and max. velocity over time: Average velocity converges smoothly to an assymptotic value of 0.1368. Compared with the analytical solution of V(avg)=0.1347, this is only 1.6% too high.
- Min. and max. density over time: Density fluctuations are reasonably small and decrease as the simulation progresses.
- Velocity profile at final time: Simulated velocity profile matches the analytical solution quite well.
Simulation #W2: To check whether the improvement is due in part to the increase in Reynolds number, we decreased the gravity to 0.000637 in order to make the Re=1.0 again. For these parameters, the analytical solution gives an average velocity of 0.0399. Simulation results are summarized below (click on links to see plots):
- Mean and max. velocity over time: Average velocity converges smoothly to an assymptotic value of 0.0404. Compared with the analytical solution of V(avg)=0.0399, this is only 1.2% too high.
- Min. and max. density over time: Density fluctuations are reasonably small and decrease as the simulation progresses.
- Velocity profile at final time: Simulated velocity profile matches the analytical solution quite well.
Simulation #W3: In the simulation above (W2), the mean velocity appears to still be increasing at the last time step simulated (t=300). So, we ran this simulation out longer to see if the velocity continues to increase. Simulation results are summarized below (click on links to see plots):
- Mean and max. velocity over time: Average velocity converges smoothly to an assymptotic value of 0.0414. Compared with the analytical solution of V(avg)=0.0399, this is 3.6% too high.
- Velocity profile at final time: Simulated velocities are slightly higher than the analytical solution, but much better than for the case with narrower slot width.
Discussion:
Results of initial simulations were a bit perplexing in that the simulated velocities were too large relative to the analytical solution. It is known that the SPH method can introduce "artificial viscosity" into the solution, which would tend to lead to *lower* velocities, not higher. However, the simulation with increased aperture width gave significantly improved results, supporting the hypothesis that the overprediction of velocity is associated with edge effects (imperfect enforcement of no-slip boundary condition). Some methods for improving the no-slip boundary implementation have been proposed but have not yet been implemented in this code.
The results of the simulations are insensitive to the speed of sound parameter as long as it is sufficiently small. However, flow between two parallel plates does not create any "barriers" to flow that would tend to lead to larger density fluctuations in compressible flow. Therefore, in general flow through porous media may be more sensitive to the speed of sound parameter than these simple test simulations.
For the right set of parameters and sufficiently large aperture that boundary effects are minimal, the SPH_CCA code can accurately reproduce the analytical solution for low Reynolds number flow between two parallel plates.