timprove porosity estimation at top BC - sphere - GPU-based 3D discrete element method algorithm with optional fluid coupling
(HTM) git clone git://src.adamsgaard.dk/sphere
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) LICENSE
---
(DIR) commit 9916a978ac5afdf76edb17d6c3e154a064cffadd
(DIR) parent de5a49343d0628cabc393a22841501d2666073af
(HTM) Author: Anders Damsgaard Christensen <adc@geo.au.dk>
Date: Tue, 27 Sep 2016 14:09:42 -0700
improve porosity estimation at top BC
Diffstat:
M src/darcy.cuh | 12 ++++++++----
M src/device.cu | 6 +++---
2 files changed, 11 insertions(+), 7 deletions(-)
---
(DIR) diff --git a/src/darcy.cuh b/src/darcy.cuh
t@@ -380,8 +380,8 @@ __device__ Float weightDist(
return 0.0;
}
-// Find the porosity in each cell on the base of a sphere, centered at the cell
-// center.
+// Find the porosity in each cell on the base of a bilinear weighing scheme,
+// centered at the cell center.
__global__ void findDarcyPorositiesLinear(
const unsigned int* __restrict__ dev_cellStart, // in
const unsigned int* __restrict__ dev_cellEnd, // in
t@@ -610,13 +610,17 @@ __global__ void findDarcyPorositiesLinear(
}
}
+ Float cell_volume = dx*dy*dz;
+ if (z == nz - 1)
+ cell_volume *= 0.75;
+
// Make sure that the porosity is in the interval [0.0;1.0]
//phi = fmin(0.9, fmax(0.1, void_volume/(dx*dy*dz)));
//phi = fmin(1.0, fmax(0.01, 1.0 - solid_volume/(dx*dy*dz)));
- phi = fmin(0.9, fmax(0.1, 1.0 - solid_volume/(dx*dy*dz)));
+ phi = fmin(0.9, fmax(0.1, 1.0 - solid_volume/cell_volume));
//Float phi_new = fmin(1.0, fmax(0.01,
Float phi_new = fmin(0.9, fmax(0.1,
- 1.0 - solid_volume_new/(dx*dy*dz)));
+ 1.0 - solid_volume_new/cell_volume));
Float dphi;
Float3 vp_avg;
(DIR) diff --git a/src/device.cu b/src/device.cu
t@@ -2037,8 +2037,8 @@ __host__ void DEM::startTime()
cudaThreadSynchronize();
}
- // copy porosities to the frictionless lower Z boundary
- if (grid.adaptive == 1) {
+ // copy porosities to the upper Z boundary
+ /*if (grid.adaptive == 1) {
copyDarcyPorositiesToTop<<<dimGridFluid,
dimBlockFluid>>>(
dev_darcy_phi,
t@@ -2046,7 +2046,7 @@ __host__ void DEM::startTime()
dev_darcy_div_v_p,
dev_darcy_vp_avg);
cudaThreadSynchronize();
- }
+ }*/
// Modulate the pressures at the upper boundary cells
if ((darcy.p_mod_A > 1.0e-5 || darcy.p_mod_A < -1.0e-5) &&