• API
• FAQ
• Tools
• Archive
daily pastebin goal
25%
SHARE
TWEET

Untitled

a guest Feb 22nd, 2019 67 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
2. //  /*
3. //    Notes
4. //    -----
5. //    This is basically the implementation of Eq.(28) in the following paper:
6. //    "a generalized wall boundary condition for SPH"
7. //    by : adami et al. 2012
8. //   */
9. //
10. //  singleton_physics  *physics   = get_singleton_physics();
11. //  singleton_geometry *geometry  = get_singleton_geometry();
12. //
13. //  double X = physics->X;
14. //  double gamma_gas = physics->gamma_gas;
15. //  double gamma_liq = physics->gamma_liq;
16. //  double epsilon = 1e-12;
17. //
18. //  for (unsigned int b = 0; b < geometry->num_cell; b++) {
19. //      // find neighbors
20. //      int gi = 0;
21. //      int gj = 0;
22. //      grid_unhash(b,&gi,&gj);
23. //
24. //      int low_i  = gi-1 < 0 ? 0 : gi-1;
25. //      int low_j  = gj-1 < 0 ? 0 : gj-1;
26. //      int high_i = gi+2 > geometry->nx ? geometry->nx : gi+2;
27. //      int high_j = gj+2 > geometry->ny ? geometry->ny : gj+2;
28. //
29. //      for (unsigned int i = geometry->cell_start[b]; i < geometry->cell_end[b]; i++) {
30. //          if (!particles[i].bnd) continue;
31. //          double P0_fluid = particles[i].p0;
32. //          double rho0_fluid = particles[i].rho0;
33. //          double gamma = (particles[i].tag == AIR) ? gamma_gas : gamma_liq;
34. //
35. //          // loop over dummy neighbors
36. //          for (unsigned int ni = low_i; ni < high_i; ni++) {
37. //              for (unsigned int nj = low_j; nj < high_j; nj++) {
38. //                  for (unsigned int j = geometry->cell_start[ni*geometry->ny+nj]; j < geometry->cell_end[ni*geometry->ny+nj]; j++) {
39. //                      if (particles[j].bnd) continue;
40. //                      double Pwall = particles[j].p;
41. //
42. //                      // note: particle "i" from fluid phase is now interacting with particle "j" from dummy wall phase
43. //
44. //                      assert(fabs(P0_fluid) > epsilon);
45. //                      double frac = ((Pwall - X)/P0_fluid + 1.);
46. //                      assert(frac > epsilon);
47. //                      double rho_wall = rho0_fluid * pow(frac , (1./gamma));
48. //                      assert(rho_wall > epsilon);
49. //
50. //                      particles[i].rho = rho_wall;
51. //                  }
52. //              }
53. //          }
54. //      }
55. //  }
56. //}
57.
59.     /*
60.     Notes
61.     -----
62.     This is basically the implementation of Eq.(28) in the following paper:
63.     "a generalized wall boundary condition for SPH"
64.     by : adami et al. 2012
65.      */
66.
67.     singleton_physics  *physics   = get_singleton_physics();
68.     singleton_geometry *geometry  = get_singleton_geometry();
69.
70.     double X = physics->X;
71.     double gamma_gas = physics->gamma_gas;
72.     double gamma_liq = physics->gamma_liq;
73.     double epsilon = 1e-12;
74.
75.     // loop over dummy particles at wall only
76.     for (unsigned int i = 0; i < geometry->n; i++) {
77.         if (!particles[i].bnd) continue;
78.
79.         // assign density/pressure of the interacting actual fluid particle
80.         double rho0_fluid = (particles[i].rho0==physics->rho_wat) ? physics->rho_wat : physics->rho_air;
81.         double P0_fluid = physics->p0L;
82.         double Pwall = particles[i].p;
83.         double gamma = (particles[i].tag == AIR) ? gamma_gas : gamma_liq;
84.         double frac = ((Pwall - X)/P0_fluid + 1.);
85.         assert(frac > epsilon);
86.         double rho_wall = rho0_fluid * pow(frac , (1./gamma));
87.         assert(rho_wall > epsilon);
88.
89.         particles[i].rho = rho_wall;
90.     }
91. }
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.

Top