tmore work on FD Darcy - 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 72813c78d9d2c75a835a58382f4219361bb1c41b
(DIR) parent b9a5fcbe8fbd1f8109ecb6a0c2a1e98aedc00bc4
(HTM) Author: Anders Damsgaard <adc@geo.au.dk>
Date: Thu, 6 Jun 2013 13:24:51 +0200
more work on FD Darcy
Diffstat:
M src/darcy.cpp | 87 +++++++++++++++++++++++--------
M src/device.cu | 17 ++++++++++++++++-
M src/sphere.h | 18 ++++++++++++++++++
3 files changed, 99 insertions(+), 23 deletions(-)
---
(DIR) diff --git a/src/darcy.cpp b/src/darcy.cpp
t@@ -103,6 +103,7 @@ void DEM::findDarcyGradients()
H = d_H[cellidx]; // cell hydraulic pressure
+ // Second order central difference
d_dH[cellidx].x
= (d_H[idx(ix+1,iy,iz)] - 2.0*H + d_H[idx(ix-1,iy,iz)])/dx2;
t@@ -263,7 +264,9 @@ void DEM::findDarcyVelocities()
// Approximate cell porosity
Float n = cellPorosity(ix, iy, iz);
- // Calculate flux (if 0,0,0 is lower left corner)
+ // Calculate flux
+ // The sign might need to be reversed, depending on the
+ // grid orientation
q.x = -d_K[cellidx]/nu * dH.x;
q.y = -d_K[cellidx]/nu * dH.y;
q.z = -d_K[cellidx]/nu * dH.z;
t@@ -278,26 +281,31 @@ void DEM::findDarcyVelocities()
}
}
-// Find particles with centres inside a spatial interval
-// NOTE: This function is untested and unused
-std::vector<unsigned int> DEM::particlesInCell(
- const Float3 min, const Float3 max)
+// Return the lower corner coordinates of a cell
+Float3 DEM::cellMinBoundaryDarcy(
+ const unsigned int x,
+ const unsigned int y,
+ const unsigned int z)
{
- // Particles radii inside cell will be stored in this vector
- std::vector<unsigned int> pidx;
-
- unsigned int i;
- Float4 x;
- for (i=0; i<np; ++i) {
+ const Float3 x_min = {x*d_dx, y*d_dy, z*d_dz};
+ return x_min;
+}
- // Read the position
- x = k.x[i];
+// Return the lower corner coordinates of a cell
+Float3 DEM::cellMaxBoundaryDarcy(
+ const unsigned int x,
+ const unsigned int y,
+ const unsigned int z)
+{
+ const Float3 x_max = {(x+1)*d_dx, (y+1)*d_dy, (z+1)*d_dz};
+ return x_max;
+}
- if (x.x >= min.x && x.y >= min.y && x.z >= min.z
- && x.x < max.x && x.y < max.y && x.z < max.z) {
- pidx.push_back(i);
- }
- }
+// Return the volume of a cell
+Float DEM::cellVolumeDarcy()
+{
+ const Float cell_volume = d_dx*d_dy*d_dz;
+ return cell_volume;
}
// Find the porosity of a target cell
t@@ -306,10 +314,9 @@ Float DEM::cellPorosity(
const unsigned int y,
const unsigned int z)
{
- const Float3 x_min = {x*d_dx, y*d_dy, z*d_dz};
- const Float3 x_max = {x_min.x+d_dx, x_min.y+d_dy, x_min.z+d_dz};
- const Float cell_volume = d_dx*d_dy*d_dz;
-
+ const Float3 x_min = cellMinBoundaryDarcy(x,y,z);
+ const Float3 x_max = cellMaxBoundaryDarcy(x,y,z);
+ Float cell_volume = cellVolumeDarcy();
Float void_volume = cell_volume;
unsigned int i;
t@@ -331,6 +338,42 @@ Float DEM::cellPorosity(
return n;
}
+// Find particles with centres inside a spatial interval
+// NOTE: This function is untested and unused
+std::vector<unsigned int> DEM::particlesInCell(
+ const Float3 min, const Float3 max)
+{
+ // Particles radii inside cell will be stored in this vector
+ std::vector<unsigned int> pidx;
+
+ unsigned int i;
+ Float4 x;
+ for (i=0; i<np; ++i) {
+
+ // Read the position
+ x = k.x[i];
+
+ if (x.x >= min.x && x.y >= min.y && x.z >= min.z
+ && x.x < max.x && x.y < max.y && x.z < max.z) {
+ pidx.push_back(i);
+ }
+ }
+}
+
+// Add fluid drag to the particles inside each cell
+void DEM::fluidDragDarcy()
+{
+ unsigned int ix, iy, iz, cellidx;
+ for (ix=0; ix<d_nx; ++ix) {
+ for (iy=0; iy<d_ny; ++iy) {
+ for (iz=0; iz<d_nz; ++iz) {
+
+
+ }
+ }
+ }
+}
+
// Solve Darcy flow on a regular, cubic grid
void DEM::initDarcy(const Float cellsizemultiplier)
{
(DIR) diff --git a/src/device.cu b/src/device.cu
t@@ -857,9 +857,24 @@ __host__ void DEM::startTime()
}
// Solve darcy flow through grid
- if (darcy == 1)
+ if (darcy == 1) {
+
+ // Copy device data to host memory
+ transferFromGlobalDeviceMemory();
+
+ // Pause the CPU thread until all CUDA calls previously issued are completed
+ cudaThreadSynchronize();
+
+ // Perform explicit Darcy time step
explDarcyStep();
+ // Transfer data from host to device memory
+ transferToGlobalDeviceMemory();
+
+ // Pause the CPU thread until all CUDA calls previously issued are completed
+ cudaThreadSynchronize();
+ }
+
// Update particle kinematics
if (PROFILING == 1)
startTimer(&kernel_tic);
(DIR) diff --git a/src/sphere.h b/src/sphere.h
t@@ -177,6 +177,24 @@ class DEM {
std::vector<unsigned int> particlesInCell(
const Float3 min, const Float3 max);
+ // Return the lower corner coordinates of a cell
+ Float3 cellMinBoundaryDarcy(
+ const unsigned int x,
+ const unsigned int y,
+ const unsigned int z);
+
+ // Return the upper corner coordinates of a cell
+ Float3 cellMaxBoundaryDarcy(
+ const unsigned int x,
+ const unsigned int y,
+ const unsigned int z);
+
+ // Returns the cell volume
+ Float cellVolumeDarcy();
+
+ // Add fluid drag to the particles inside each cell
+ void fluidDragDarcy();
+
// Find porosity of cell
Float cellPorosity(
const unsigned int x,