Updated writeblock also

This commit is contained in:
André R. Brodtkorb
2018-08-24 15:44:51 +02:00
parent 92fba5dc4f
commit a1902a1330
7 changed files with 69 additions and 136 deletions

View File

@@ -120,11 +120,10 @@ __global__ void FORCEKernel(
//Read into shared memory //Read into shared memory
readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(h0_ptr_, h0_pitch_, Q[0], nx_+1, ny_+1); readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>( h0_ptr_, h0_pitch_, Q[0], nx_+2, ny_+2);
readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hu0_ptr_, hu0_pitch_, Q[1], nx_+1, ny_+1); readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hu0_ptr_, hu0_pitch_, Q[1], nx_+2, ny_+2);
readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hv0_ptr_, hv0_pitch_, Q[2], nx_+1, ny_+1); readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hv0_ptr_, hv0_pitch_, Q[2], nx_+2, ny_+2);
__syncthreads(); __syncthreads();
//Set boundary conditions //Set boundary conditions
noFlowBoundary1(Q, nx_, ny_); noFlowBoundary1(Q, nx_, ny_);
@@ -147,10 +146,9 @@ __global__ void FORCEKernel(
__syncthreads(); __syncthreads();
//Write to main memory //Write to main memory
writeBlock1(h1_ptr_, h1_pitch_, writeBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2, 1, 1>( h1_ptr_, h1_pitch_, Q[0], nx_, ny_);
hu1_ptr_, hu1_pitch_, writeBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2, 1, 1>(hu1_ptr_, hu1_pitch_, Q[1], nx_, ny_);
hv1_ptr_, hv1_pitch_, writeBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2, 1, 1>(hv1_ptr_, hv1_pitch_, Q[2], nx_, ny_);
Q, nx_, ny_);
} }
} // extern "C" } // extern "C"

View File

@@ -166,9 +166,9 @@ __global__ void HLL2Kernel(
//Read into shared memory //Read into shared memory
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(h0_ptr_, h0_pitch_, Q[0], nx_+3, ny_+3); readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>( h0_ptr_, h0_pitch_, Q[0], nx_+4, ny_+4);
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hu0_ptr_, hu0_pitch_, Q[1], nx_+3, ny_+3); readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hu0_ptr_, hu0_pitch_, Q[1], nx_+4, ny_+4);
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hv0_ptr_, hv0_pitch_, Q[2], nx_+3, ny_+3); readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hv0_ptr_, hv0_pitch_, Q[2], nx_+4, ny_+4);
__syncthreads(); __syncthreads();
//Set boundary conditions //Set boundary conditions

View File

@@ -121,21 +121,19 @@ __global__ void HLLKernel(
float* h1_ptr_, int h1_pitch_, float* h1_ptr_, int h1_pitch_,
float* hu1_ptr_, int hu1_pitch_, float* hu1_ptr_, int hu1_pitch_,
float* hv1_ptr_, int hv1_pitch_) { float* hv1_ptr_, int hv1_pitch_) {
const int block_width = BLOCK_WIDTH;
const int block_height = BLOCK_HEIGHT;
//Shared memory variables //Shared memory variables
__shared__ float Q[3][block_height+2][block_width+2]; __shared__ float Q[3][BLOCK_HEIGHT+2][BLOCK_WIDTH+2];
__shared__ float F[3][block_height+1][block_width+1]; __shared__ float F[3][BLOCK_HEIGHT+1][BLOCK_WIDTH+1];
//Read into shared memory //Read into shared memory
readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(h0_ptr_, h0_pitch_, Q[0], nx_+1, ny_+1); readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>( h0_ptr_, h0_pitch_, Q[0], nx_+2, ny_+2);
readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hu0_ptr_, hu0_pitch_, Q[1], nx_+1, ny_+1); readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hu0_ptr_, hu0_pitch_, Q[1], nx_+2, ny_+2);
readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hv0_ptr_, hv0_pitch_, Q[2], nx_+1, ny_+1); readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hv0_ptr_, hv0_pitch_, Q[2], nx_+2, ny_+2);
__syncthreads(); __syncthreads();
//Set boundary conditions
noFlowBoundary1(Q, nx_, ny_); noFlowBoundary1(Q, nx_, ny_);
__syncthreads(); __syncthreads();
@@ -155,16 +153,10 @@ __global__ void HLLKernel(
evolveG1(Q, F, nx_, ny_, dy_, dt_); evolveG1(Q, F, nx_, ny_, dy_, dt_);
__syncthreads(); __syncthreads();
//Q[0][threadIdx.y + 1][threadIdx.x + 1] += 0.1;
// Write to main memory for all internal cells // Write to main memory for all internal cells
writeBlock1(h1_ptr_, h1_pitch_, writeBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2, 1, 1>( h1_ptr_, h1_pitch_, Q[0], nx_, ny_);
hu1_ptr_, hu1_pitch_, writeBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2, 1, 1>(hu1_ptr_, hu1_pitch_, Q[1], nx_, ny_);
hv1_ptr_, hv1_pitch_, writeBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2, 1, 1>(hv1_ptr_, hv1_pitch_, Q[2], nx_, ny_);
Q, nx_, ny_);
} }
} // extern "C" } // extern "C"

View File

@@ -157,9 +157,9 @@ __global__ void KP07DimsplitKernel(
//Read into shared memory //Read into shared memory
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(h0_ptr_, h0_pitch_, Q[0], nx_+3, ny_+3); readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>( h0_ptr_, h0_pitch_, Q[0], nx_+4, ny_+4);
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hu0_ptr_, hu0_pitch_, Q[1], nx_+3, ny_+3); readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hu0_ptr_, hu0_pitch_, Q[1], nx_+4, ny_+4);
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hv0_ptr_, hv0_pitch_, Q[2], nx_+3, ny_+3); readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hv0_ptr_, hv0_pitch_, Q[2], nx_+4, ny_+4);
__syncthreads(); __syncthreads();

View File

@@ -141,9 +141,9 @@ __global__ void KP07Kernel(
//Read into shared memory //Read into shared memory
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(h0_ptr_, h0_pitch_, Q[0], nx_+3, ny_+3); readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>( h0_ptr_, h0_pitch_, Q[0], nx_+4, ny_+4);
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hu0_ptr_, hu0_pitch_, Q[1], nx_+3, ny_+3); readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hu0_ptr_, hu0_pitch_, Q[1], nx_+4, ny_+4);
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hv0_ptr_, hv0_pitch_, Q[2], nx_+3, ny_+3); readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hv0_ptr_, hv0_pitch_, Q[2], nx_+4, ny_+4);
__syncthreads(); __syncthreads();

View File

@@ -115,58 +115,42 @@ void LxFKernel(
float* hu1_ptr_, int hu1_pitch_, float* hu1_ptr_, int hu1_pitch_,
float* hv1_ptr_, int hv1_pitch_) { float* hv1_ptr_, int hv1_pitch_) {
const int block_width = BLOCK_WIDTH; const int tx = threadIdx.x;
const int block_height = BLOCK_HEIGHT; const int ty = threadIdx.y;
//Index of cell within domain
const int ti = blockDim.x*blockIdx.x + threadIdx.x + 1; //Skip global ghost cells, i.e., +1
const int tj = blockDim.y*blockIdx.y + threadIdx.y + 1;
__shared__ float Q[3][block_height+2][block_width+2]; __shared__ float Q[3][BLOCK_HEIGHT+2][BLOCK_WIDTH+2];
__shared__ float F[3][block_height][block_width+1]; __shared__ float F[3][BLOCK_HEIGHT][BLOCK_WIDTH+1];
__shared__ float G[3][block_height+1][block_width]; __shared__ float G[3][BLOCK_HEIGHT+1][BLOCK_WIDTH];
//Read into shared memory //Read into shared memory including ghost cells
readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(h0_ptr_, h0_pitch_, Q[0], nx_+1, ny_+1); readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>( h0_ptr_, h0_pitch_, Q[0], nx_+2, ny_+2);
readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hu0_ptr_, hu0_pitch_, Q[1], nx_+1, ny_+1); readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hu0_ptr_, hu0_pitch_, Q[1], nx_+2, ny_+2);
readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hv0_ptr_, hv0_pitch_, Q[2], nx_+1, ny_+1); readBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2>(hv0_ptr_, hv0_pitch_, Q[2], nx_+2, ny_+2);
__syncthreads(); __syncthreads();
//Set boundary conditions //Set boundary conditions
noFlowBoundary1(Q, nx_, ny_); noFlowBoundary1(Q, nx_, ny_);
__syncthreads(); __syncthreads();
//Compute fluxes along the x and y axis //Compute fluxes along the x and y axis
computeFluxF<block_width, block_height>(Q, F, g_, dx_, dt_); computeFluxF<BLOCK_WIDTH, BLOCK_HEIGHT>(Q, F, g_, dx_, dt_);
computeFluxG<block_width, block_height>(Q, G, g_, dy_, dt_); computeFluxG<BLOCK_WIDTH, BLOCK_HEIGHT>(Q, G, g_, dy_, dt_);
__syncthreads(); __syncthreads();
//Evolve for all cells
//Evolve for all internal cells const int i = tx + 1; //Skip local ghost cells, i.e., +1
if (ti > 0 && ti < nx_+1 && tj > 0 && tj < ny_+1) { const int j = ty + 1;
//Index of thread within block Q[0][j][i] += (F[0][ty][tx] - F[0][ty ][tx+1]) * dt_ / dx_
const int tx = threadIdx.x; + (G[0][ty][tx] - G[0][ty+1][tx ]) * dt_ / dy_;
const int ty = threadIdx.y; Q[1][j][i] += (F[1][ty][tx] - F[1][ty ][tx+1]) * dt_ / dx_
+ (G[1][ty][tx] - G[1][ty+1][tx ]) * dt_ / dy_;
const int i = tx + 1; //Skip local ghost cells, i.e., +1 Q[2][j][i] += (F[2][ty][tx] - F[2][ty ][tx+1]) * dt_ / dx_
const int j = ty + 1; + (G[2][ty][tx] - G[2][ty+1][tx ]) * dt_ / dy_;
const float h1 = Q[0][j][i] + (F[0][ty][tx] - F[0][ty ][tx+1]) * dt_ / dx_
+ (G[0][ty][tx] - G[0][ty+1][tx ]) * dt_ / dy_;
const float hu1 = Q[1][j][i] + (F[1][ty][tx] - F[1][ty ][tx+1]) * dt_ / dx_
+ (G[1][ty][tx] - G[1][ty+1][tx ]) * dt_ / dy_;
const float hv1 = Q[2][j][i] + (F[2][ty][tx] - F[2][ty ][tx+1]) * dt_ / dx_
+ (G[2][ty][tx] - G[2][ty+1][tx ]) * dt_ / dy_;
float* const h_row = (float*) ((char*) h1_ptr_ + h1_pitch_*tj); //Write to main memory
float* const hu_row = (float*) ((char*) hu1_ptr_ + hu1_pitch_*tj); writeBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2, 1, 1>( h1_ptr_, h1_pitch_, Q[0], nx_, ny_);
float* const hv_row = (float*) ((char*) hv1_ptr_ + hv1_pitch_*tj); writeBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2, 1, 1>(hu1_ptr_, hu1_pitch_, Q[1], nx_, ny_);
writeBlock<BLOCK_WIDTH+2, BLOCK_HEIGHT+2, 1, 1>(hv1_ptr_, hv1_pitch_, Q[2], nx_, ny_);
h_row[ti] = h1;
hu_row[ti] = hu1;
hv_row[ti] = hv1;
}
} }
} // extern "C" } // extern "C"

View File

@@ -56,7 +56,7 @@ inline __device__ __host__ float clamp(const float f, const float a, const float
template<int sm_width, int sm_height> template<int sm_width, int sm_height>
__device__ void readBlock(float* ptr_, int pitch_, __device__ void readBlock(float* ptr_, int pitch_,
float shmem[sm_height][sm_width], float shmem[sm_height][sm_width],
const int max_x, const int max_y) { const int max_x_, const int max_y_) {
//Index of block within domain //Index of block within domain
const int bx = blockDim.x * blockIdx.x; const int bx = blockDim.x * blockIdx.x;
@@ -64,13 +64,13 @@ __device__ void readBlock(float* ptr_, int pitch_,
//Read into shared memory //Read into shared memory
for (int j=threadIdx.y; j<sm_height; j+=blockDim.y) { for (int j=threadIdx.y; j<sm_height; j+=blockDim.y) {
const int l = clamp(by + j, 0, max_y); // Clamp out of bounds const int l = clamp(by + j, 0, max_y_-1); // Clamp out of bounds
//Compute the pointer to current row in the arrays //Compute the pointer to current row in the arrays
const float* const row = (float*) ((char*) ptr_ + pitch_*l); const float* const row = (float*) ((char*) ptr_ + pitch_*l);
for (int i=threadIdx.x; i<sm_width; i+=blockDim.x) { for (int i=threadIdx.x; i<sm_width; i+=blockDim.x) {
const int k = clamp(bx + i, 0, max_x); // clamp out of bounds const int k = clamp(bx + i, 0, max_x_-1); // Clamp out of bounds
shmem[j][i] = row[k]; shmem[j][i] = row[k];
} }
@@ -81,52 +81,26 @@ __device__ void readBlock(float* ptr_, int pitch_,
/**
* Reads a block of data with two ghost cells for the shallow water equations
*/
__device__ void readBlock2(float* h_ptr_, int h_pitch_,
float* hu_ptr_, int hu_pitch_,
float* hv_ptr_, int hv_pitch_,
float Q[3][BLOCK_HEIGHT+4][BLOCK_WIDTH+4],
const int nx_, const int ny_) {
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(h_ptr_, h_pitch_, Q[0], nx_+3, ny_+3);
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hu_ptr_, hu_pitch_, Q[1], nx_+3, ny_+3);
readBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4>(hv_ptr_, hv_pitch_, Q[2], nx_+3, ny_+3);
}
/** /**
* Writes a block of data to global memory for the shallow water equations. * Writes a block of data to global memory for the shallow water equations.
*/ */
__device__ void writeBlock1(float* h_ptr_, int h_pitch_, template<int sm_width, int sm_height, int offset_x=0, int offset_y=0>
float* hu_ptr_, int hu_pitch_, __device__ void writeBlock(float* ptr_, int pitch_,
float* hv_ptr_, int hv_pitch_, float shmem[sm_height][sm_width],
float Q[3][BLOCK_HEIGHT+2][BLOCK_WIDTH+2], const int width, const int height) {
const int nx_, const int ny_) {
//Index of thread within block
const int tx = threadIdx.x;
const int ty = threadIdx.y;
//Index of cell within domain //Index of cell within domain
const int ti = blockDim.x*blockIdx.x + threadIdx.x + 1; //Skip global ghost cells, i.e., +1 const int ti = blockDim.x*blockIdx.x + threadIdx.x + offset_x;
const int tj = blockDim.y*blockIdx.y + threadIdx.y + 1; const int tj = blockDim.y*blockIdx.y + threadIdx.y + offset_y;
//Only write internal cells //Only write internal cells
if (ti > 0 && ti < nx_+1 && tj > 0 && tj < ny_+1) { if (ti < width+offset_x && tj < height+offset_y) {
const int i = tx + 1; //Skip local ghost cells, i.e., +1 //Index of thread within block
const int j = ty + 1; const int tx = threadIdx.x + offset_x;
const int ty = threadIdx.y + offset_y;
float* const h_row = (float*) ((char*) h_ptr_ + h_pitch_*tj);
float* const hu_row = (float*) ((char*) hu_ptr_ + hu_pitch_*tj);
float* const hv_row = (float*) ((char*) hv_ptr_ + hv_pitch_*tj);
h_row[ti] = Q[0][j][i]; float* const row = (float*) ((char*) ptr_ + pitch_*tj);
hu_row[ti] = Q[1][j][i]; row[ti] = shmem[ty][tx];
hv_row[ti] = Q[2][j][i];
} }
} }
@@ -134,6 +108,9 @@ __device__ void writeBlock1(float* h_ptr_, int h_pitch_,
/** /**
* Writes a block of data to global memory for the shallow water equations. * Writes a block of data to global memory for the shallow water equations.
*/ */
@@ -142,27 +119,9 @@ __device__ void writeBlock2(float* h_ptr_, int h_pitch_,
float* hv_ptr_, int hv_pitch_, float* hv_ptr_, int hv_pitch_,
float Q[3][BLOCK_HEIGHT+4][BLOCK_WIDTH+4], float Q[3][BLOCK_HEIGHT+4][BLOCK_WIDTH+4],
const int nx_, const int ny_) { const int nx_, const int ny_) {
//Index of thread within block writeBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4, 2, 2>( h_ptr_, h_pitch_, Q[0], nx_, ny_);
const int tx = threadIdx.x; writeBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4, 2, 2>(hu_ptr_, hu_pitch_, Q[1], nx_, ny_);
const int ty = threadIdx.y; writeBlock<BLOCK_WIDTH+4, BLOCK_HEIGHT+4, 2, 2>(hv_ptr_, hv_pitch_, Q[2], nx_, ny_);
//Index of cell within domain
const int ti = blockDim.x*blockIdx.x + threadIdx.x + 2; //Skip global ghost cells, i.e., +2
const int tj = blockDim.y*blockIdx.y + threadIdx.y + 2;
//Only write internal cells
if (ti > 1 && ti < nx_+2 && tj > 1 && tj < ny_+2) {
const int i = tx + 2; //Skip local ghost cells, i.e., +2
const int j = ty + 2;
float* const h_row = (float*) ((char*) h_ptr_ + h_pitch_*tj);
float* const hu_row = (float*) ((char*) hu_ptr_ + hu_pitch_*tj);
float* const hv_row = (float*) ((char*) hv_ptr_ + hv_pitch_*tj);
h_row[ti] = Q[0][j][i];
hu_row[ti] = Q[1][j][i];
hv_row[ti] = Q[2][j][i];
}
} }