Verification/fda nozzle

[YuVirtonomy]

Case Description: FDA Nozzle Benchmark Model for Medical Device Evaluation

This verification test leverages the FDA's nozzle benchmark model, a simplified and idealized construct of a medical device. The model is designed to simulate critical fluid dynamics that may influence blood damage, such as accelerating and decelerating flows, varying shear stresses, velocities, and recirculating flows.

Model Specifications:

• Geometry:
  The model features an axisymmetric nozzle with a stenotic throat measuring 0.04 meters in length. It includes a 20° connecting cone at one end and a step change in diameter at the other end, facilitating flow from left to right.
• Physical Properties:
  The fluid within the model has a density of 1056 kg/m³ and a dynamic viscosity of 3.5 mPa·s.
• Flow Conditions:
  The Reynolds number, a key parameter in the study, is defined at the throat of the nozzle.

Reference data

Axial velocity along nozzle centerline for throat Reynolds number (a) 500, (b) 2000

Resources and Further Reading:

• Detailed specifications and results for this model can be found on the FDA's dedicated here
• A comparative study utilizing Lattice Boltzmann Methods (LBM) is documented here paper

Numerical Model setup:

Geometry:

Geometry Specification:

The illustration above details the geometry of the FDA nozzle, integral to this experimental setup.

Geometry Specification for the FDA Nozzle Model

In the current model configuration, the FDA nozzle has been deliberately extended by 0.012 meters at both the inlet and outlet to minimize potential numerical artifacts arising from boundary conditions. This modification is crucial for ensuring that the simulation results accurately reflect the fluid dynamics within the device without being unduly influenced by the boundaries.

Orientation and Dimensions:

• Axial Direction: The x-axis has been designated as the axial direction in this revised model setup.
• Extended Geometry: Each end of the nozzle now includes an additional 0.012 meters, providing a buffer that helps stabilize flow before it enters and after it exits the measurement area.

Illustration of the Extended Nozzle Geometry:

Boundary Condition:

This case will firstly focus on Reynolds number = 500, in which the flowrate Q_f is set to 5e-6

• Inlet Boundary Condition: Parabolic velocity inlet where U_max inlet = 2 * Q_f / Area_inlet
• Outlet Boundary Condition: Pressure 0
• U_max: 2 * Q_f / Area_throat
• Dynamic viscosity (μ): 0.0035 [N s/m2]
• Density (ρ): 1056 [kg/m3]
• Resolution (dp): D_throat/ N, where N is the number of particles
• pressure_relaxation: Integration1stHalfWithWallRiemann
• density_relaxation: Integration2ndHalfWithWallRiemann

Simulation Results: Initial Setup

The duration for all simulations is set to 5 seconds, which exceeds the time necessary to reach a convergent state, ensuring the accuracy and reliability of the results. Following the example, the buffer zone width is set to 5*dp.

Particle Count N = 10

The following images display the simulation results for a particle count of N = 10, showcasing different aspects of fluid dynamics within the model.

Axial Velocity Profiles

Velocity:

Pressure:

Velocity Profile Along the Axial Direction:

Velocity Profile Along the Radial Direction at x = -0.064m:

Particle Count N = 15

The following images display the simulation results for a particle count of N = 15, showcasing different aspects of fluid dynamics within the model.

Axial Velocity Profiles

Velocity:

Pressure:

Velocity Profile Along the Axial Direction:

Velocity Profile Along the Radial Direction at x = -0.064m:

Summary

The peak velocity occurs right at the sudden expansion. For N=15, the velocity at the throat region aligns well with experimental data. However, the fluid velocity between the inlet and throat does not exhibit the expected parabolic flow, but more like turbulent flow, which contrasts with FDA experimental data suggesting a parabolic profile in this region. In the sudden expansion region, the fluid velocity appears slower than experimental data, resembling results from simulations using the k-epsilon turbulent model rather than those from CFD simulations. Notably, the flow near the outlet is less symmetric than expected.

Focus on the velocity between inlet and throat region. (N=10)

Further investigating of the non-parabolic profile between inlet and throat region.

Changing Inflow Buffer Zone Width

[Xiangyu-Hu]

@YuVirtonomy First the definition of Reynold number. Second, do you have time average of the simulation? Third, how is the new results?

[Xiangyu-Hu]

@YuVirtonomy Another comment. If is the results is very asymmetric near entrance, there is must be something wrong.

[Xiangyu-Hu]

One more thing. It is suggested for channel flow simulation, one always accelerate to entire with gravity first if the inflow condition prescribes velocity profile.

[Xiangyu-Hu]

@YuVirtonomy the fdA case is too large so it triggers timeout. You may do not run full time length or using less resolution.

[Xiangyu-Hu]

@YuVirtonomy or we can make as special user case in a separate repo.

[YuVirtonomy]

@FengWang3119 can you review this and suggest me some solution?

Simulation set up

Re = 500
Inflow condition:
Parabolic velocity profile is specified at the entrance, with the maximum velocity at the centerline approximately equal to 0.088.
Outflow condition:
Pressure outflow condition with pressure set to 1.0
Resolution:
Resolution is defined by throat diameter / Number of Particles (NP)
For NP = 15, total fluid particles is about 0.9 million

Simulation result with NP = 15

Velocity profile of NP = 15, t= 5s

axial velocity near by inlet, nozzle (sudden expansion side), and near by outlet. The convergence is calculated based on the axial velocity at nozzle

axial velocity at x = -0.088
axial velocity at x = -0.064
axial velocity at x = -0.048
axial velocity at x = -0.02
axial velocity at x = -0.008
axial velocity at x = 0
axial velocity at x = 0.008
axial velocity at x = 0.016
axial velocity at x = 0.024
axial velocity at x = 0.032
axial velocity at x = 0.06
axial velocity at x = 0.08

numbers of particles vs time

Issue

Asymmetry Post-Sudden Expansion: The simulation results exhibit unexpected asymmetrical behavior following the sudden expansion area.
Velocity Dissipation: There is noticeable velocity dissipation after the sudden expansion, which should not be as pronounced according to theoretical predictions.

[Xiangyu-Hu]

Is this a 3d simulation? Also, do you have obtained results by time average?

[YuVirtonomy]

Is this a 3d simulation? Also, do you have obtained results by time average?
it is 3D case.

No, this is not the time-averaged data. However, the figure Velocity and Convergence Rate shows the velocity evolution over time. In the graph, the orange line represents the velocity of an observer placed in the axial direction near the inlet. The blue line indicates the velocity of an observer near the throat and sudden expansion region. The green line depicts the velocity of an observer near the outlet. The red lines is the converge rate based on the blue line. As time progresses, we observe that the velocity of observer near the outlet decreases over time.

[Xiangyu-Hu]

Is this a 3d simulation? Also, do you have obtained results by time average?
it is 3D case.

No, this is not the time-averaged data. However, the figure Velocity and Convergence Rate shows the velocity evolution over time. In the graph, the orange line represents the velocity of an observer placed in the axial direction near the inlet. The blue line indicates the velocity of an observer near the throat and sudden expansion region. The green line depicts the velocity of an observer near the outlet. The red lines is the converge rate based on the blue line. As time progresses, we observe that the velocity of observer near the outlet decreases over time.

Do you mean that the asymmetrical result is converged solution?

[YuVirtonomy]

Is this a 3d simulation? Also, do you have obtained results by time average?
it is 3D case.

No, this is not the time-averaged data. However, the figure Velocity and Convergence Rate shows the velocity evolution over time. In the graph, the orange line represents the velocity of an observer placed in the axial direction near the inlet. The blue line indicates the velocity of an observer near the throat and sudden expansion region. The green line depicts the velocity of an observer near the outlet. The red lines is the converge rate based on the blue line. As time progresses, we observe that the velocity of observer near the outlet decreases over time.

Do you mean that the asymmetrical result is converged solution?

I have updated the simulation result using average velocity data.

As this figure indicates, the simulation results are converged at the inlet, nozzle, and outlet. The velocity measurements are taken along the x-coordinate (specify the x-value), with the y and z coordinates fixed at zero, ensuring measurements in the axial centerline.

this is the velocity profile at the outlet, obviously asymmetric.

The velocity profile at the outlet, as shown, is clearly asymmetric.
Overall, the observed velocity dissipation may be attributed to the flow not being symmetric.

I will provide the results with a resolution of Np = 10 and 20 at a later time.

[Xiangyu-Hu]

It is very strange that time average results of a symmetrical geometry is not symmetrical.

[YuVirtonomy]

It is very strange that time average results of a symmetrical geometry is not symmetrical.

Yes, especially at low resolution.
As you can see in the attached gif, the asymmetry becomes more pronounced as the resolution decreases.

Resolution Np = 10 for end_time =5s

Do you think without using wall relaxation will improve this problem?

[Xiangyu-Hu]

It is very strange that time average results of a symmetrical geometry is not symmetrical.

Yes, especially at low resolution. As you can see in the attached gif, the asymmetry becomes more pronounced as the resolution decreases. Resolution Np = 10 for end_time =5s

Do you think without using wall relaxation will improve this problem?

The results is very interesting. Do you have other runs which its result gives the same asymmetric bias direction?
@Bo-Zhang1995 Its seems this is quite good case for your paper on Eulerian SPH.

[Xiangyu-Hu]

You may use dissipative Riemann solver to have a look.

[Xiangyu-Hu]

@Xiangyu-Hu

You may use dissipative Riemann solver to have a look.

I assume you are referring to using Integration1stHalfWithWall with the DissipativeRiemannSolver?

I think that is the 2nd half for the momentum equation.

I configured Integration2ndHalfWithWall to use the DissipativeRiemannSolver as follows: using Integration2ndHalfWithWallDissipativeRiemann = Integration2ndHalfWithWall<DissipativeRiemannSolver>; The results with resolution = throat_diameter/10 are shown below:

Interesting, it seems that there is a parameter between dissipative to limited riemann solver may be leads to symmetrical solution but without too much damping.

[YuVirtonomy]

@Xiangyu-Hu

Interesting, it seems that there is a parameter between dissipative to limited riemann solver may be leads to symmetrical solution but without too much damping.

One thing I noticed is that with increasing resolution, the simulation tends to become more symmetric. The following results were obtained using the following configuration:

**Integration1stHalfCorrectionWithWallRiemann**
**Integration2ndHalfWithWall**
The results are snapshots of about 1 second, as NP = 20 experiences some issues after 1 second. However, it seems 1 second is sufficient to consider reaching a steady state.

NP = 10 (0.25 million fluid particles)

NP = 15 (1 million fluid particles)

NP = 20 (2 million fluid particles)

Axial Velocity Profiles at Different Resolutions

Figure of the axial velocity profiles measured along the x-axis (axial direction), at three different resolutions (NP10, NP15, and NP20).

For NP = 20, the velocity at x = 0.08 (near the outlet) can reach 0.35 m/s, while in experimental data, it is about 0.42 to 0.5 m/s. I will provide more details with NP = 20 later, as this simulation terminates before I can summarize the results.

@Xiangyu-Hu any suggestion of what to try next?

[Xiangyu-Hu]

@Xiangyu-Hu

Interesting, it seems that there is a parameter between dissipative to limited riemann solver may be leads to symmetrical solution but without too much damping.

One thing I noticed is that with increasing resolution, the simulation tends to become more symmetric. The following results were obtained using the following configuration:

**Integration1stHalfCorrectionWithWallRiemann** **Integration2ndHalfWithWall** The results are snapshots of about 1 second, as NP = 20 experiences some issues after 1 second. However, it seems 1 second is sufficient to consider reaching a steady state.

NP = 10 (0.25 million fluid particles)

NP = 15 (1 million fluid particles)

NP = 20 (2 million fluid particles)

Axial Velocity Profiles at Different Resolutions

Figure of the axial velocity profiles measured along the x-axis (axial direction), at three different resolutions (NP10, NP15, and NP20).

For NP = 20, the velocity at x = 0.08 (near the outlet) can reach 0.35 m/s, while in experimental data, it is about 0.42 to 0.5 m/s. I will provide more details with NP = 20 later, as this simulation terminates before I can summarize the results.

@Xiangyu-Hu any suggestion of what to try next?

It is very strange that time average results of a symmetrical geometry is not symmetrical.

Yes, especially at low resolution. As you can see in the attached gif, the asymmetry becomes more pronounced as the resolution decreases.

    [
      
    
        ![output10_0717](https://private-user-images.githubusercontent.com/111049554/349850134-9fdbf7a0-1f7f-4d7b-ad71-c31f23db2c47.gif?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEyODkxMzIsIm5iZiI6MTcyMTI4ODgzMiwicGF0aCI6Ii8xMTEwNDk1NTQvMzQ5ODUwMTM0LTlmZGJmN2EwLTFmN2YtNGQ3Yi1hZDcxLWMzMWYyM2RiMmM0Ny5naWY_WC1BbXotQWxnb3JpdGhtPUFXUzQtSE1BQy1TSEEyNTYmWC1BbXotQ3JlZGVudGlhbD1BS0lBVkNPRFlMU0E1M1BRSzRaQSUyRjIwMjQwNzE4JTJGdXMtZWFzdC0xJTJGczMlMkZhd3M0X3JlcXVlc3QmWC1BbXotRGF0ZT0yMDI0MDcxOFQwNzQ3MTJaJlgtQW16LUV4cGlyZXM9MzAwJlgtQW16LVNpZ25hdHVyZT05YjVhODVmZWE0OWJkNjdiMzRmYjFhNDczZjVjYTRmYTRkZmY2MWRhOTFkMTIwMTg5YTNlZGM1YTE4ZTgzOWI0JlgtQW16LVNpZ25lZEhlYWRlcnM9aG9zdCZhY3Rvcl9pZD0wJmtleV9pZD0wJnJlcG9faWQ9MCJ9.QI5MH0qCCaGgbhM_siw_JBF-jNjlGPl4mUysMcwjf3A)
      ](https://private-user-images.githubusercontent.com/111049554/349850134-9fdbf7a0-1f7f-4d7b-ad71-c31f23db2c47.gif?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEyODkxMzIsIm5iZiI6MTcyMTI4ODgzMiwicGF0aCI6Ii8xMTEwNDk1NTQvMzQ5ODUwMTM0LTlmZGJmN2EwLTFmN2YtNGQ3Yi1hZDcxLWMzMWYyM2RiMmM0Ny5naWY_WC1BbXotQWxnb3JpdGhtPUFXUzQtSE1BQy1TSEEyNTYmWC1BbXotQ3JlZGVudGlhbD1BS0lBVkNPRFlMU0E1M1BRSzRaQSUyRjIwMjQwNzE4JTJGdXMtZWFzdC0xJTJGczMlMkZhd3M0X3JlcXVlc3QmWC1BbXotRGF0ZT0yMDI0MDcxOFQwNzQ3MTJaJlgtQW16LUV4cGlyZXM9MzAwJlgtQW16LVNpZ25hdHVyZT05YjVhODVmZWE0OWJkNjdiMzRmYjFhNDczZjVjYTRmYTRkZmY2MWRhOTFkMTIwMTg5YTNlZGM1YTE4ZTgzOWI0JlgtQW16LVNpZ25lZEhlYWRlcnM9aG9zdCZhY3Rvcl9pZD0wJmtleV9pZD0wJnJlcG9faWQ9MCJ9.QI5MH0qCCaGgbhM_siw_JBF-jNjlGPl4mUysMcwjf3A)
    
    
      
        
          
        
        
          
          
        
      
      [
        
          
        
      ](https://private-user-images.githubusercontent.com/111049554/349850134-9fdbf7a0-1f7f-4d7b-ad71-c31f23db2c47.gif?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEyODkxMzIsIm5iZiI6MTcyMTI4ODgzMiwicGF0aCI6Ii8xMTEwNDk1NTQvMzQ5ODUwMTM0LTlmZGJmN2EwLTFmN2YtNGQ3Yi1hZDcxLWMzMWYyM2RiMmM0Ny5naWY_WC1BbXotQWxnb3JpdGhtPUFXUzQtSE1BQy1TSEEyNTYmWC1BbXotQ3JlZGVudGlhbD1BS0lBVkNPRFlMU0E1M1BRSzRaQSUyRjIwMjQwNzE4JTJGdXMtZWFzdC0xJTJGczMlMkZhd3M0X3JlcXVlc3QmWC1BbXotRGF0ZT0yMDI0MDcxOFQwNzQ3MTJaJlgtQW16LUV4cGlyZXM9MzAwJlgtQW16LVNpZ25hdHVyZT05YjVhODVmZWE0OWJkNjdiMzRmYjFhNDczZjVjYTRmYTRkZmY2MWRhOTFkMTIwMTg5YTNlZGM1YTE4ZTgzOWI0JlgtQW16LVNpZ25lZEhlYWRlcnM9aG9zdCZhY3Rvcl9pZD0wJmtleV9pZD0wJnJlcG9faWQ9MCJ9.QI5MH0qCCaGgbhM_siw_JBF-jNjlGPl4mUysMcwjf3A)
    
   [ ![output10_0717](https://private-user-images.githubusercontent.com/111049554/349850134-9fdbf7a0-1f7f-4d7b-ad71-c31f23db2c47.gif?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEyODkxMzIsIm5iZiI6MTcyMTI4ODgzMiwicGF0aCI6Ii8xMTEwNDk1NTQvMzQ5ODUwMTM0LTlmZGJmN2EwLTFmN2YtNGQ3Yi1hZDcxLWMzMWYyM2RiMmM0Ny5naWY_WC1BbXotQWxnb3JpdGhtPUFXUzQtSE1BQy1TSEEyNTYmWC1BbXotQ3JlZGVudGlhbD1BS0lBVkNPRFlMU0E1M1BRSzRaQSUyRjIwMjQwNzE4JTJGdXMtZWFzdC0xJTJGczMlMkZhd3M0X3JlcXVlc3QmWC1BbXotRGF0ZT0yMDI0MDcxOFQwNzQ3MTJaJlgtQW16LUV4cGlyZXM9MzAwJlgtQW16LVNpZ25hdHVyZT05YjVhODVmZWE0OWJkNjdiMzRmYjFhNDczZjVjYTRmYTRkZmY2MWRhOTFkMTIwMTg5YTNlZGM1YTE4ZTgzOWI0JlgtQW16LVNpZ25lZEhlYWRlcnM9aG9zdCZhY3Rvcl9pZD0wJmtleV9pZD0wJnJlcG9faWQ9MCJ9.QI5MH0qCCaGgbhM_siw_JBF-jNjlGPl4mUysMcwjf3A) ](https://private-user-images.githubusercontent.com/111049554/349850134-9fdbf7a0-1f7f-4d7b-ad71-c31f23db2c47.gif?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEyODkxMzIsIm5iZiI6MTcyMTI4ODgzMiwicGF0aCI6Ii8xMTEwNDk1NTQvMzQ5ODUwMTM0LTlmZGJmN2EwLTFmN2YtNGQ3Yi1hZDcxLWMzMWYyM2RiMmM0Ny5naWY_WC1BbXotQWxnb3JpdGhtPUFXUzQtSE1BQy1TSEEyNTYmWC1BbXotQ3JlZGVudGlhbD1BS0lBVkNPRFlMU0E1M1BRSzRaQSUyRjIwMjQwNzE4JTJGdXMtZWFzdC0xJTJGczMlMkZhd3M0X3JlcXVlc3QmWC1BbXotRGF0ZT0yMDI0MDcxOFQwNzQ3MTJaJlgtQW16LUV4cGlyZXM9MzAwJlgtQW16LVNpZ25hdHVyZT05YjVhODVmZWE0OWJkNjdiMzRmYjFhNDczZjVjYTRmYTRkZmY2MWRhOTFkMTIwMTg5YTNlZGM1YTE4ZTgzOWI0JlgtQW16LVNpZ25lZEhlYWRlcnM9aG9zdCZhY3Rvcl9pZD0wJmtleV9pZD0wJnJlcG9faWQ9MCJ9.QI5MH0qCCaGgbhM_siw_JBF-jNjlGPl4mUysMcwjf3A)
  
    [
      
    
        ![output10_0717](https://private-user-images.githubusercontent.com/111049554/349850134-9fdbf7a0-1f7f-4d7b-ad71-c31f23db2c47.gif?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEyODkxMzIsIm5iZiI6MTcyMTI4ODgzMiwicGF0aCI6Ii8xMTEwNDk1NTQvMzQ5ODUwMTM0LTlmZGJmN2EwLTFmN2YtNGQ3Yi1hZDcxLWMzMWYyM2RiMmM0Ny5naWY_WC1BbXotQWxnb3JpdGhtPUFXUzQtSE1BQy1TSEEyNTYmWC1BbXotQ3JlZGVudGlhbD1BS0lBVkNPRFlMU0E1M1BRSzRaQSUyRjIwMjQwNzE4JTJGdXMtZWFzdC0xJTJGczMlMkZhd3M0X3JlcXVlc3QmWC1BbXotRGF0ZT0yMDI0MDcxOFQwNzQ3MTJaJlgtQW16LUV4cGlyZXM9MzAwJlgtQW16LVNpZ25hdHVyZT05YjVhODVmZWE0OWJkNjdiMzRmYjFhNDczZjVjYTRmYTRkZmY2MWRhOTFkMTIwMTg5YTNlZGM1YTE4ZTgzOWI0JlgtQW16LVNpZ25lZEhlYWRlcnM9aG9zdCZhY3Rvcl9pZD0wJmtleV9pZD0wJnJlcG9faWQ9MCJ9.QI5MH0qCCaGgbhM_siw_JBF-jNjlGPl4mUysMcwjf3A)
      ](https://private-user-images.githubusercontent.com/111049554/349850134-9fdbf7a0-1f7f-4d7b-ad71-c31f23db2c47.gif?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEyODkxMzIsIm5iZiI6MTcyMTI4ODgzMiwicGF0aCI6Ii8xMTEwNDk1NTQvMzQ5ODUwMTM0LTlmZGJmN2EwLTFmN2YtNGQ3Yi1hZDcxLWMzMWYyM2RiMmM0Ny5naWY_WC1BbXotQWxnb3JpdGhtPUFXUzQtSE1BQy1TSEEyNTYmWC1BbXotQ3JlZGVudGlhbD1BS0lBVkNPRFlMU0E1M1BRSzRaQSUyRjIwMjQwNzE4JTJGdXMtZWFzdC0xJTJGczMlMkZhd3M0X3JlcXVlc3QmWC1BbXotRGF0ZT0yMDI0MDcxOFQwNzQ3MTJaJlgtQW16LUV4cGlyZXM9MzAwJlgtQW16LVNpZ25hdHVyZT05YjVhODVmZWE0OWJkNjdiMzRmYjFhNDczZjVjYTRmYTRkZmY2MWRhOTFkMTIwMTg5YTNlZGM1YTE4ZTgzOWI0JlgtQW16LVNpZ25lZEhlYWRlcnM9aG9zdCZhY3Rvcl9pZD0wJmtleV9pZD0wJnJlcG9faWQ9MCJ9.QI5MH0qCCaGgbhM_siw_JBF-jNjlGPl4mUysMcwjf3A)
    
    
      
        
          
        
        
          
          
        
      
      [
        
          
        
      ](https://private-user-images.githubusercontent.com/111049554/349850134-9fdbf7a0-1f7f-4d7b-ad71-c31f23db2c47.gif?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEyODkxMzIsIm5iZiI6MTcyMTI4ODgzMiwicGF0aCI6Ii8xMTEwNDk1NTQvMzQ5ODUwMTM0LTlmZGJmN2EwLTFmN2YtNGQ3Yi1hZDcxLWMzMWYyM2RiMmM0Ny5naWY_WC1BbXotQWxnb3JpdGhtPUFXUzQtSE1BQy1TSEEyNTYmWC1BbXotQ3JlZGVudGlhbD1BS0lBVkNPRFlMU0E1M1BRSzRaQSUyRjIwMjQwNzE4JTJGdXMtZWFzdC0xJTJGczMlMkZhd3M0X3JlcXVlc3QmWC1BbXotRGF0ZT0yMDI0MDcxOFQwNzQ3MTJaJlgtQW16LUV4cGlyZXM9MzAwJlgtQW16LVNpZ25hdHVyZT05YjVhODVmZWE0OWJkNjdiMzRmYjFhNDczZjVjYTRmYTRkZmY2MWRhOTFkMTIwMTg5YTNlZGM1YTE4ZTgzOWI0JlgtQW16LVNpZ25lZEhlYWRlcnM9aG9zdCZhY3Rvcl9pZD0wJmtleV9pZD0wJnJlcG9faWQ9MCJ9.QI5MH0qCCaGgbhM_siw_JBF-jNjlGPl4mUysMcwjf3A)
    
   [ ](https://private-user-images.githubusercontent.com/111049554/349850134-9fdbf7a0-1f7f-4d7b-ad71-c31f23db2c47.gif?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEyODkxMzIsIm5iZiI6MTcyMTI4ODgzMiwicGF0aCI6Ii8xMTEwNDk1NTQvMzQ5ODUwMTM0LTlmZGJmN2EwLTFmN2YtNGQ3Yi1hZDcxLWMzMWYyM2RiMmM0Ny5naWY_WC1BbXotQWxnb3JpdGhtPUFXUzQtSE1BQy1TSEEyNTYmWC1BbXotQ3JlZGVudGlhbD1BS0lBVkNPRFlMU0E1M1BRSzRaQSUyRjIwMjQwNzE4JTJGdXMtZWFzdC0xJTJGczMlMkZhd3M0X3JlcXVlc3QmWC1BbXotRGF0ZT0yMDI0MDcxOFQwNzQ3MTJaJlgtQW16LUV4cGlyZXM9MzAwJlgtQW16LVNpZ25hdHVyZT05YjVhODVmZWE0OWJkNjdiMzRmYjFhNDczZjVjYTRmYTRkZmY2MWRhOTFkMTIwMTg5YTNlZGM1YTE4ZTgzOWI0JlgtQW16LVNpZ25lZEhlYWRlcnM9aG9zdCZhY3Rvcl9pZD0wJmtleV9pZD0wJnJlcG9faWQ9MCJ9.QI5MH0qCCaGgbhM_siw_JBF-jNjlGPl4mUysMcwjf3A) Resolution Np = 10 for end_time =5s

Do you think without using wall relaxation will improve this problem?

The results is very interesting. Do you have other runs which its result gives the same asymmetric bias direction? @Bo-Zhang1995 Its seems this is quite good case for your paper on Eulerian SPH.

Please have a check whether the biased results is repeatable in different runs.

[Xiangyu-Hu]

@YuVirtonomy could you also check geometry, i.e. the level set? May be there is some defects. I suddenly thought about this when I use new geometry in a recent pull request #626.

[YuVirtonomy]

@YuVirtonomy could you also check geometry, i.e. the level set? May be there is some defects. I suddenly thought about this when I use new geometry in a recent pull request #626.

Here are the results showing the overlap of the wall level set shape with the wall mesh:

Resolution Np 5

Resolution Np 10

Resolution Np 15

Resolution Np 20

Resolution Np 25

Overall, aside from the resolution, the 90-degree corner becomes sharper. Otherwise, the level set and mesh almost align perfectly.

[YuVirtonomy]

Update with Np = 25 result

Simulation set up

Re = 500
Inflow condition:
Parabolic velocity profile is specified at the entrance, with the maximum velocity at the centerline approximately equal to 0.088.
Outflow condition:
Pressure outflow condition with pressure set to 1.0
Resolution:
Resolution is defined by throat diameter / Number of Particles (NP)
For NP = 25, total fluid particles is about 4 million, Total run time = 1.5s

Simulation result with NP = 25

Velocity profile of NP = 25, t= 1.5s

axial velocity at x = -0.088
axial velocity at x = -0.064

axial velocity at x = -0.048

axial velocity at x = -0.02

axial velocity at x = -0.008

axial velocity at x = 0

axial velocity at x = 0.008

axial velocity at x = 0.016

axial velocity at x = 0.024

axial velocity at x = 0.032

axial velocity at x = 0.06

axial velocity at x = 0.08

Compare with Np = 15, we can see the flow profile getting closer to the experimental data

[Xiangyu-Hu]

@YuVirtonomy could you have careful check whether there are small solid pieces with the fluid domain (which may pollute the fluid domain). You can check it by cut the level set zero contour.

[YuVirtonomy]

@YuVirtonomy could you have careful check whether there are small solid pieces with the fluid domain (which may pollute the fluid domain). You can check it by cut the level set zero contour.

I have carefully checked the fluid domain for any small solid pieces that might contaminate it. I inspected the zero level set contour, it appears that there are no solid particles polluting the fluid domain and the fluid and solid level set are perfectly align

[Xiangyu-Hu]

@YuVirtonomy I am more interested on the lower resolution results and wonder why the the simulation results are asymmetric. Basically, I am thinking the possibility of polluted fluid domain if the the asymmetric bias is reproducible. Another thing is that we can use the linear reproducing wall boundary condition when computing viscous force. This may increase the resolution near the wall.

[Xiangyu-Hu]

Please also have a check on the new pull request #629 in which a linear corrected viscous force is implemented.

[YuVirtonomy]

Issue Focus: Sudden Expansion Part to Outlet

To address the issue between the sudden expansion part and the outlet, the concentrate part has been removed. Below are the details of the simulations and the observed results. The simulation ran for 2 seconds, which is sufficient to consider as convergence. The inlet is set so that the velocity profile at x = 0 aligns with experimental data.

Simulation Results

Velocity Measurements

The graphs below show the comparison between the simulation results using ViscousForceWithWall and ViscousForceWithWallCorrection.

• Position Z vs. avg_velocity X: Represents the simulation result using ViscousForceWithWall.
• Position Z vs. avg_velocity X (ViscCorrection): Represents the simulation result using ViscousForceWithWallCorrection.

Velocity Measured at x = 0

Velocity Measured at x = 0.08

Additional Algorithms

Details of other algorithms used in the simulations:

• LinearGradientCorrectionMatrixComplex
• Integration1stHalfCorrectionWithWallRiemann
• Integration2ndHalfWithWallRiemann
• DensitySummationPressureComplex
• TransportVelocityCorrectionCorrectedComplex

Summary

With the use of ViscousForceWithWallCorrection, no apparent improvement is observed. However, in some previous tests, it was found that increasing the wall resolution (wall_boundary.defineAdaptationRatios(1.15, 2.0);) can improve the result, though it is still far from the experimental data. A case with wall_boundary.defineAdaptationRatios(1.15, 2.0); is still running, and updates will be provided later.

[Xiangyu-Hu]

@YuVirtonomy Another important test you should have a look is using gravity to replace pressure boundary condition and have a look on the results.

[YuVirtonomy]

Visualization of Cross Section View at x = 0.08

For this setup, I used a resolution of NP = 10 with the second half without a Riemann solver. I created a clip to better visualize the cross-section view at x = 0.08, as shown below:

Original Mesh vs. Rotated Mesh

The cross-section view where A to C is using the original mesh, and D to F use the mesh rotated 90 degrees along the axial side:

Summary

The results do not exhibit symmetry across different simulations.

Details of algorithms

• Wall has run 1000 relaxation step
• LinearGradientCorrectionMatrixComplex
• Integration1stHalfCorrectionWithWallRiemann
• Integration2ndHalfWithWall**NO**Riemann
• DensitySummationPressureComplex
• TransportVelocityCorrectionCorrectedComplex

[YuVirtonomy]

@YuVirtonomy Another important test you should have a look is using gravity to replace pressure boundary condition and have a look on the results.

Do you meant directly apply acceleration? but over here i only set pressure zero at outlet

[Xiangyu-Hu]

@YuVirtonomy Good. The results show no preference on any asymmetric direction, it just gives unsteady solution. This quite reasonable,like the experimental results. Therefore, you need use temporal average to obtained the flow velocity profiles.

[Xiangyu-Hu]

@YuVirtonomy Another important test you should have a look is using gravity to replace pressure boundary condition and have a look on the results.

Do you meant directly apply acceleration? but over here i only set pressure zero at outlet

I mean do not use pressure boundary condition to check whether pressure gradient computing has any serious problem.
If you use gravity force, you do not need pressure condition to impose flow.

[YuVirtonomy]

@YuVirtonomy Good. The results show no preference on any asymmetric direction, it just gives unsteady solution. This quite reasonable,like the experimental results. Therefore, you need use temporal average to obtained the flow velocity profiles.

Simulation Results

In this simulation, ( N_p = 8 ) was used to verify the asymmetric issue. The simulation runs for a total of 5 seconds, with the inflow velocity increasing from 0 to the target velocity within the first 0.25 seconds. Results are recorded every 0.0005 seconds. The average velocity is calculated over the last 2000 steps, representing the average velocity for the most recent 1-second interval.

Simulation Output

The average velocity of cross section at x = 0.08

The average velocity of cross section at x = 0.08

Set legend at velocity = 0.07 to visualize the peak value

Algorithm Used

• PressureCondition
• TransportVelocityCorrectionCorrectedComplex
• DensitySummationPressureComplex
• Integration2ndHalfWithWallNoRiemann
• Integration1stHalfCorrectionWithWallRiemann
• LinearGradientCorrectionMatrixComplex
• ViscousForceWithWallCorrection

Summary

As shown above, although the result is smoother and exhibits some asymmetry, the last figure still does not show significant asymmetry.

Question:
Why doesn't the flow exhibit as much unsteadiness with better resolution?

Hypothesis:
I believe one reason for the asymmetry is the arrangement of wall particles near the sudden expansion region. The arrangement of these particles significantly influences the jet direction near the wall. This makes sense because, with increased resolution, the jet appears more asymmetric.

To verify this, I set the entire throat channel, plus a bit of extension, with a parabolic inflow velocity and ran exactly the same test as above.

The left side is the original simulation from above, and the right side is the result setting the entire throat channel with parabolic velocity to avoid the influence of wall particles near the sudden expansion part. It can clearly be seen that the asymmetry improves significantly. However, the velocity dissipation is still high.

[YuVirtonomy]

@YuVirtonomy Another important test you should have a look is using gravity to replace pressure boundary condition and have a look on the results.

Do you meant directly apply acceleration? but over here i only set pressure zero at outlet

I mean do not use pressure boundary condition to check whether pressure gradient computing has any serious problem. If you use gravity force, you do not need pressure condition to impose flow.

@YuVirtonomy Another important test you should have a look is using gravity to replace pressure boundary condition and have a look on the results.

Do you meant directly apply acceleration? but over here i only set pressure zero at outlet

I mean do not use pressure boundary condition to check whether pressure gradient computing has any serious problem. If you use gravity force, you do not need pressure condition to impose flow.

@YuVirtonomy Another important test you should have a look is using gravity to replace pressure boundary condition and have a look on the results.

Do you meant directly apply acceleration? but over here i only set pressure zero at outlet

I mean do not use pressure boundary condition to check whether pressure gradient computing has any serious problem. If you use gravity force, you do not need pressure condition to impose flow.

I have tried to avoid using the pressure condition by using DensitySummationFreeStreamComplex directly, as done in the T-shaped pipe example.

With ( N_p = 10 ), the simulation failed due to an issue at the outlet. I was unable to capture the last step before the simulation failed.

With ( N_p = 15 ), the simulation seemed to work.

I set up two cases:

1. Using the regular inflow condition with a 20-resolution inflow boundary length.
2. Assigning the entire throat channel with an inflow condition to avoid the influence of wall particles at the sudden expansion region.

The boundary conditions and algorithms used follow the T-shaped pipe example but utilize Integration2ndHalfWithWallRiemann.

For the version where the entire throat channel is assigned with a velocity condition, the simulation failed at 2.9 seconds, likely due to reverse flow issues at the outlet.

For the version with a 20-resolution inflow boundary length, the simulation reached 4.6 seconds without issues.

The upper set of figures shows the entire throat channel assigned with the velocity condition, and the lower half shows the 20-resolution inflow boundary length version.

Velocity profile at x = 0.08, with entire throat channel assigned with velocity condition.

Velocity profile at x = 0.08, with 20-resolution inflow boundary length version.

Summary

It seems that prescribing the entire throat channel with a velocity condition can generally increase asymmetry. It may be worth to repeat the case with Integration2ndHalfWithWallNoRiemann.

[Xiangyu-Hu]

@YuVirtonomy
Your summary is quite confusing. My questions are (1) does the asymmetry has direction preference, or does the asymmetric decreases with increase longer time average. (2) are the last two figures time average results? if they are, how long the average is carried out?
Furthermore, as I mentioned before, we will not use Riemann solver in Integration2ndHalf for viscous flow.

[Xiangyu-Hu]

@YuVirtonomy Also, how long time the experiment has used for time average?

[YuVirtonomy]

@YuVirtonomy Also, how long time the experiment has used for time average?

In their study. for the Re =500 condition using Particle Image Velocimetry (PIV), approximately 500 image pairs were captured with intervals between frames varying by region: 1.4 ms for the entrance, 0.395 ms for the throat, and 1.5 ms for the recirculation areas. Thus, the total experiment duration for capturing all images was estimated to be around 750 ms, assuming sequential image capture without overlap.

[YuVirtonomy]

@YuVirtonomy Your summary is quite confusing. My questions are (1) does the asymmetry has direction preference, or does the asymmetric decreases with increase longer time average. (2) are the last two figures time average results? if they are, how long the average is carried out? Furthermore, as I mentioned before, we will not use Riemann solver in Integration2ndHalf for viscous flow.

(1) The asymmetry does not have a direction preference, and it does not decrease with a longer time average. In fact, the flow always remains unstable at low resolution.

(2) Yes, the last two figures are time-averaged results. The average was carried out over 1 second, using 2000 frames in this 1-second period.