Skip to content

GitLab

Commit a60ed4d7 authored by lishen's avatar lishen
Browse files

warpctc from github

parent 949fcc19
Hide whitespace changes
Inline Side-by-side
Showing with 2347 additions and 0 deletions
.gitignore 0 → 100644
View file @ a60ed4d7
*~
Makefile
build
\ No newline at end of file
.idea/.gitignore 0 → 100644
View file @ a60ed4d7
# Default ignored files
/shelf/
/workspace.xml
# Editor-based HTTP Client requests
/httpRequests/
# Datasource local storage ignored files
/dataSources/
/dataSources.local.xml
.idea/deployment.xml 0 → 100644
View file @ a60ed4d7
<?xml version="1.0" encoding="UTF-8"?>
<project version="4">
<component name="PublishConfigData" autoUpload="On explicit save action" serverName="10.6.10.62" preserveTimestamps="false" deleteMissingItems="true" createEmptyFolders="true" filePermissions="420" folderPermissions="493" confirmBeforeUploading="false" confirmBeforeDeletion="false" autoUploadExternalChanges="true">
<option name="confirmBeforeDeletion" value="false" />
<option name="confirmBeforeUploading" value="false" />
<serverData>
<paths name="10.6.10.61">
<serverdata>
<mappings>
<mapping local="$PROJECT_DIR$" web="/" />
</mappings>
</serverdata>
</paths>
<paths name="10.6.10.62">
<serverdata>
<mappings>
<mapping deploy="/public/home/lishen/warpctc/warpctc_nv" local="$PROJECT_DIR$" web="/" />
</mappings>
</serverdata>
</paths>
<paths name="10.6.10.69">
<serverdata>
<mappings>
<mapping local="$PROJECT_DIR$" web="/" />
</mappings>
</serverdata>
</paths>
<paths name="10.6.6.220">
<serverdata>
<mappings>
<mapping local="$PROJECT_DIR$" web="/" />
</mappings>
</serverdata>
</paths>
<paths name="lishen_89_docker (c96f91f6-1c12-4fb0-b3d8-cc68ccf3f77c)">
<serverdata>
<mappings>
<mapping local="$PROJECT_DIR$" web="/" />
</mappings>
</serverdata>
</paths>
<paths name="lishen_90">
<serverdata>
<mappings>
<mapping local="$PROJECT_DIR$" web="/" />
</mappings>
</serverdata>
</paths>
<paths name="lishen_93">
<serverdata>
<mappings>
<mapping local="$PROJECT_DIR$" web="/" />
</mappings>
</serverdata>
</paths>
<paths name="lishen_95">
<serverdata>
<mappings>
<mapping local="$PROJECT_DIR$" web="/" />
</mappings>
</serverdata>
</paths>
</serverData>
<option name="myAutoUpload" value="ON_EXPLICIT_SAVE" />
</component>
</project>
\ No newline at end of file
.idea/modules.xml 0 → 100644
View file @ a60ed4d7
<?xml version="1.0" encoding="UTF-8"?>
<project version="4">
<component name="ProjectModuleManager">
<modules>
<module fileurl="file://$PROJECT_DIR$/.idea/warpctc_nv.iml" filepath="$PROJECT_DIR$/.idea/warpctc_nv.iml" />
</modules>
</component>
</project>
\ No newline at end of file
.idea/vcs.xml 0 → 100644
View file @ a60ed4d7
<?xml version="1.0" encoding="UTF-8"?>
<project version="4">
<component name="VcsDirectoryMappings">
<mapping directory="$PROJECT_DIR$" vcs="Git" />
</component>
</project>
\ No newline at end of file
.idea/warpctc_nv.iml 0 → 100644
View file @ a60ed4d7
<?xml version="1.0" encoding="UTF-8"?>
<module type="CPP_MODULE" version="4">
<component name="NewModuleRootManager">
<content url="file://$MODULE_DIR$" />
<orderEntry type="inheritedJdk" />
<orderEntry type="sourceFolder" forTests="false" />
</component>
</module>
\ No newline at end of file
.travis.yml 0 → 100644
View file @ a60ed4d7
sudo: required
dist: xenial
language: python
python:
- "3.6"
install:
## Start conda installation
# https://conda.io/docs/user-guide/tasks/use-conda-with-travis-ci.html
- sudo apt-get update
- wget https://repo.continuum.io/miniconda/Miniconda3-latest-Linux-x86_64.sh -O miniconda.sh;
- bash miniconda.sh -b -p $HOME/miniconda
- export PATH="$HOME/miniconda/bin:$PATH"
- hash -r
- conda config --set always_yes yes --set changeps1 no
- conda update -q conda
# Useful for debugging any issues with conda
- conda info -a
## End conda installation
- conda create -q -n test-environment python=$TRAVIS_PYTHON_VERSION pytest pytest-flakes
- source activate test-environment
- conda install pytorch torchvision -c pytorch
- mkdir build && cd build && cmake .. && make && sudo make install && cd .. && sudo ldconfig
- cd pytorch_binding && python setup.py install && cd .. && sudo ldconfig
script:
- cd pytorch_binding && py.test tests && cd ..
- cd pytorch_binding && py.test --flakes --ignore=warpctc_pytorch/_warp_ctc && cd ..
CMakeLists.txt 0 → 100644
View file @ a60ed4d7
IF (APPLE)
cmake_minimum_required(VERSION 3.4)
ELSE()
cmake_minimum_required(VERSION 2.8)
ENDIF()
project(ctc_release)
IF (NOT APPLE)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -O2")
ENDIF()
IF (APPLE)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++14 -O2")
add_definitions(-DAPPLE)
ENDIF()
include_directories(include)
FIND_PACKAGE(CUDA 6.5)
MESSAGE(STATUS "cuda found ${CUDA_FOUND}")
option(WITH_GPU "compile warp-ctc with cuda." ${CUDA_FOUND})
option(WITH_OMP "compile warp-ctc with openmp." ON)
if(NOT WITH_OMP)
add_definitions(-DCTC_DISABLE_OMP)
endif()
if (WITH_OMP)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -fopenmp")
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -Xcompiler -fopenmp")
endif()
# need to be at least 30 or __shfl_down in reduce wont compile
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -gencode arch=compute_30,code=sm_30 -O2")
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -gencode arch=compute_35,code=sm_35")
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -gencode arch=compute_50,code=sm_50")
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -gencode arch=compute_52,code=sm_52")
IF(CMAKE_CXX_COMPILER_VERSION VERSION_GREATER 5)
SET(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -D_MWAITXINTRIN_H_INCLUDED -D_FORCE_INLINES")
ENDIF()
IF (CUDA_VERSION GREATER 7.6)
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -gencode arch=compute_60,code=sm_60")
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -gencode arch=compute_61,code=sm_61")
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -gencode arch=compute_62,code=sm_62")
ENDIF()
IF (CUDA_VERSION GREATER 8.9)
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} -gencode arch=compute_70,code=sm_70")
ENDIF()
if (NOT APPLE)
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS} --std=c++14")
set(CUDA_NVCC_FLAGS "${CUDA_NVCC_FLAGS}")
ENDIF()
IF (APPLE)
EXEC_PROGRAM(uname ARGS -v OUTPUT_VARIABLE DARWIN_VERSION)
STRING(REGEX MATCH "[0-9]+" DARWIN_VERSION ${DARWIN_VERSION})
MESSAGE(STATUS "DARWIN_VERSION=${DARWIN_VERSION}")
#for el capitain have to use rpath
IF (DARWIN_VERSION LESS 15)
set(CMAKE_SKIP_RPATH TRUE)
ENDIF ()
ELSE()
#always skip for linux
set(CMAKE_SKIP_RPATH TRUE)
ENDIF()
IF (WITH_GPU)
MESSAGE(STATUS "Building shared library with GPU support")
CUDA_ADD_LIBRARY(warpctc SHARED src/ctc_entrypoint.cu src/reduce.cu)
IF (!Torch_FOUND)
TARGET_LINK_LIBRARIES(warpctc ${CUDA_curand_LIBRARY})
ENDIF()
add_executable(test_cpu tests/test_cpu.cpp tests/random.cpp )
TARGET_LINK_LIBRARIES(test_cpu warpctc)
SET_TARGET_PROPERTIES(test_cpu PROPERTIES COMPILE_FLAGS "${CMAKE_CXX_FLAGS} --std=c++14")
cuda_add_executable(test_gpu tests/test_gpu.cu tests/random.cpp )
TARGET_LINK_LIBRARIES(test_gpu warpctc ${CUDA_curand_LIBRARY})
SET_TARGET_PROPERTIES(test_gpu PROPERTIES COMPILE_FLAGS "${CMAKE_CXX_FLAGS} --std=c++14")
INSTALL(TARGETS warpctc
RUNTIME DESTINATION "bin"
LIBRARY DESTINATION "lib"
ARCHIVE DESTINATION "lib")
INSTALL(FILES include/ctc.h DESTINATION "include")
ELSE()
MESSAGE(STATUS "Building shared library with no GPU support")
if (NOT APPLE)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++14 -O2")
ENDIF()
ADD_LIBRARY(warpctc SHARED src/ctc_entrypoint.cpp)
add_executable(test_cpu tests/test_cpu.cpp tests/random.cpp )
TARGET_LINK_LIBRARIES(test_cpu warpctc)
SET_TARGET_PROPERTIES(test_cpu PROPERTIES COMPILE_FLAGS "${CMAKE_CXX_FLAGS} --std=c++14")
INSTALL(TARGETS warpctc
RUNTIME DESTINATION "bin"
LIBRARY DESTINATION "lib"
ARCHIVE DESTINATION "lib")
INSTALL(FILES include/ctc.h DESTINATION "include")
ENDIF()
README.md 0 → 100644
View file @ a60ed4d7
# PyTorch bindings for Warp-ctc
[![Build Status](https://travis-ci.org/SeanNaren/warp-ctc.svg?branch=pytorch_bindings)](https://travis-ci.org/SeanNaren/warp-ctc)
This is an extension onto the original repo found [here](https://github.com/baidu-research/warp-ctc).
## Installation
Install [PyTorch](https://github.com/pytorch/pytorch#installation) v0.4.
`WARP_CTC_PATH` should be set to the location of a built WarpCTC
(i.e. `libwarpctc.so`). This defaults to `../build`, so from within a
new warp-ctc clone you could build WarpCTC like this:
```bash
git clone https://github.com/SeanNaren/warp-ctc.git
cd warp-ctc
mkdir build; cd build
cmake ..
make
```
Now install the bindings:
```bash
cd pytorch_binding
python setup.py install
```
If you try the above and get a dlopen error on OSX with anaconda3 (as recommended by pytorch):
```bash
cd ../pytorch_binding
python setup.py install
cd ../build
cp libwarpctc.dylib /Users/$WHOAMI/anaconda3/lib
```
This will resolve the library not loaded error. This can be easily modified to work with other python installs if needed.
Example to use the bindings below.
```python
import torch
from warpctc_pytorch import CTCLoss
ctc_loss = CTCLoss()
# expected shape of seqLength x batchSize x alphabet_size
probs = torch.FloatTensor([[[0.1, 0.6, 0.1, 0.1, 0.1], [0.1, 0.1, 0.6, 0.1, 0.1]]]).transpose(0, 1).contiguous()
labels = torch.IntTensor([1, 2])
label_sizes = torch.IntTensor([2])
probs_sizes = torch.IntTensor([2])
probs.requires_grad_(True) # tells autograd to compute gradients for probs
cost = ctc_loss(probs, labels, probs_sizes, label_sizes)
cost.backward()
```
## Documentation
```
CTCLoss(size_average=False, length_average=False)
# size_average (bool): normalize the loss by the batch size (default: False)
# length_average (bool): normalize the loss by the total number of frames in the batch. If True, supersedes size_average (default: False)
forward(acts, labels, act_lens, label_lens)
# acts: Tensor of (seqLength x batch x outputDim) containing output activations from network (before softmax)
# labels: 1 dimensional Tensor containing all the targets of the batch in one large sequence
# act_lens: Tensor of size (batch) containing size of each output sequence from the network
# label_lens: Tensor of (batch) containing label length of each example
```
\ No newline at end of file
include/contrib/moderngpu/LICENSE 0 → 100644
View file @ a60ed4d7
/******************************************************************************
* Copyright (c) 2013, NVIDIA CORPORATION. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
* * Neither the name of the NVIDIA CORPORATION nor the
* names of its contributors may be used to endorse or promote products
* derived from this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
******************************************************************************/
include/contrib/moderngpu/include/device/ctaloadbalance.cuh 0 → 100644
View file @ a60ed4d7
/******************************************************************************
* Copyright (c) 2013, NVIDIA CORPORATION. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
* * Neither the name of the NVIDIA CORPORATION nor the
* names of its contributors may be used to endorse or promote products
* derived from this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
******************************************************************************/
/******************************************************************************
*
* Code and text by Sean Baxter, NVIDIA Research
* See http://nvlabs.github.io/moderngpu for repository and documentation.
*
******************************************************************************/
#pragma once
#include "ctasearch.cuh"
#include "loadstore.cuh"
namespace mgpu {
////////////////////////////////////////////////////////////////////////////////
// DeviceLoadBalancingSearch
// Upper Bound search from A (needles) into B (haystack). The A values are
// natural numbers from aBegin to aEnd. bFirst is the index of the B value at
// bBegin in shared memory.
template<int VT, bool RangeCheck>
MGPU_DEVICE void DeviceSerialLoadBalanceSearch(const int* b_shared, int aBegin,
int aEnd, int bFirst, int bBegin, int bEnd, int* a_shared) {
int bKey = b_shared[bBegin];
#pragma unroll
for(int i = 0; i < VT; ++i) {
bool p;
if(RangeCheck)
p = (aBegin < aEnd) && ((bBegin >= bEnd) || (aBegin < bKey));
else
p = aBegin < bKey;
if(p)
// Advance A (the needle).
a_shared[aBegin++] = bFirst + bBegin;
else
// Advance B (the haystack).
bKey = b_shared[++bBegin];
}
}
////////////////////////////////////////////////////////////////////////////////
// CTALoadBalance
// Computes upper_bound(counting_iterator<int>(first), b_global) - 1.
// Unlike most other CTA* functions, CTALoadBalance loads from global memory.
// This returns the loaded B elements at the beginning or end of shared memory
// depending on the aFirst argument.
// CTALoadBalance requires NT * VT + 2 slots of shared memory.
template<int NT, int VT, typename InputIt>
MGPU_DEVICE int4 CTALoadBalance(int destCount, InputIt b_global,
int sourceCount, int block, int tid, const int* mp_global,
int* indices_shared, bool loadPrecedingB) {
int4 range = ComputeMergeRange(destCount, sourceCount, block, 0, NT * VT,
mp_global);
int a0 = range.x;
int a1 = range.y;
int b0 = range.z;
int b1 = range.w;
if(!b0) loadPrecedingB = false;
// Load one trailing term from B. If we're already at the end, fill the
// end of the buffer with destCount.
int aCount = a1 - a0;
int bCount = b1 - b0;
int extended = b1 < sourceCount;
int loadCount = bCount + extended;
int fillCount = NT * VT + 1 - loadCount - aCount;
int* a_shared = indices_shared;
int* b_shared = indices_shared + aCount + (int)loadPrecedingB;
// Load the B values.
// DeviceMemToMemLoop<NT>(bCount + extended + (int)loadPrecedingB,
// b_global + b0 - (int)loadPrecedingB, tid,
// b_shared - (int)loadPrecedingB);
for(int i = tid - (int)loadPrecedingB; i < bCount + extended; i += NT)
b_shared[i] = b_global[b0 + i];
// Fill the end of the array with destCount.
for(int i = tid + extended; i < fillCount; i += NT)
b_shared[bCount + i] = destCount;
__syncthreads();
// Run a merge path to find the start of the serial merge for each thread.
int diag = VT * tid;
int mp = MergePath<MgpuBoundsUpper>(mgpu::counting_iterator<int>(a0),
aCount, b_shared, bCount, diag, mgpu::less<int>());
int a0tid = a0 + mp;
int b0tid = diag - mp;
// Subtract 1 from b0 because we want to return upper_bound - 1.
DeviceSerialLoadBalanceSearch<VT, false>(b_shared, a0tid, a1, b0 - 1,
b0tid, bCount, a_shared - a0);
__syncthreads();
b0 -= (int)loadPrecedingB;
return make_int4(a0, a1, b0, b1);
}
} // namespace mgpu
include/contrib/moderngpu/include/device/ctamerge.cuh 0 → 100644
View file @ a60ed4d7
/******************************************************************************
* Copyright (c) 2013, NVIDIA CORPORATION. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
* * Neither the name of the NVIDIA CORPORATION nor the
* names of its contributors may be used to endorse or promote products
* derived from this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
******************************************************************************/
/******************************************************************************
*
* Code and text by Sean Baxter, NVIDIA Research
* See http://nvlabs.github.io/moderngpu for repository and documentation.
*
******************************************************************************/
#pragma once
#include "ctasearch.cuh"
#include "loadstore.cuh"
#include "sortnetwork.cuh"
namespace mgpu {
////////////////////////////////////////////////////////////////////////////////
// SerialMerge
template<int VT, bool RangeCheck, typename T, typename Comp>
MGPU_DEVICE void SerialMerge(const T* keys_shared, int aBegin, int aEnd,
int bBegin, int bEnd, T* results, int* indices, Comp comp) {
T aKey = keys_shared[aBegin];
T bKey = keys_shared[bBegin];
#pragma unroll
for(int i = 0; i < VT; ++i) {
bool p;
if(RangeCheck)
p = (bBegin >= bEnd) || ((aBegin < aEnd) && !comp(bKey, aKey));
else
p = !comp(bKey, aKey);
results[i] = p ? aKey : bKey;
indices[i] = p ? aBegin : bBegin - !RangeCheck;
if(p) aKey = keys_shared[++aBegin];
else bKey = keys_shared[++bBegin];
}
__syncthreads();
}
////////////////////////////////////////////////////////////////////////////////
// FindMergeFrame and FindMergesortInterval help mergesort (both CTA and global
// merge pass levels) locate lists within the single source array.
// Returns (offset of a, offset of b, length of list).
MGPU_HOST_DEVICE int3 FindMergesortFrame(int coop, int block, int nv) {
// coop is the number of CTAs or threads cooperating to merge two lists into
// one. We round block down to the first CTA's ID that is working on this
// merge.
int start = ~(coop - 1) & block;
int size = nv * (coop>> 1);
return make_int3(nv * start, nv * start + size, size);
}
// Returns (a0, a1, b0, b1) into mergesort input lists between mp0 and mp1.
MGPU_HOST_DEVICE int4 FindMergesortInterval(int3 frame, int coop, int block,
int nv, int count, int mp0, int mp1) {
// Locate diag from the start of the A sublist.
int diag = nv * block - frame.x;
int a0 = frame.x + mp0;
int a1 = min(count, frame.x + mp1);
int b0 = min(count, frame.y + diag - mp0);
int b1 = min(count, frame.y + diag + nv - mp1);
// The end partition of the last block for each merge operation is computed
// and stored as the begin partition for the subsequent merge. i.e. it is
// the same partition but in the wrong coordinate system, so its 0 when it
// should be listSize. Correct that by checking if this is the last block
// in this merge operation.
if(coop - 1 == ((coop - 1) & block)) {
a1 = min(count, frame.x + frame.z);
b1 = min(count, frame.y + frame.z);
}
return make_int4(a0, a1, b0, b1);
}
////////////////////////////////////////////////////////////////////////////////
// ComputeMergeRange
MGPU_HOST_DEVICE int4 ComputeMergeRange(int aCount, int bCount, int block,
int coop, int NV, const int* mp_global) {
// Load the merge paths computed by the partitioning kernel.
int mp0 = mp_global[block];
int mp1 = mp_global[block + 1];
int gid = NV * block;
// Compute the ranges of the sources in global memory.
int4 range;
if(coop) {
int3 frame = FindMergesortFrame(coop, block, NV);
range = FindMergesortInterval(frame, coop, block, NV, aCount, mp0,
mp1);
} else {
range.x = mp0; // a0
range.y = mp1; // a1
range.z = gid - range.x; // b0
range.w = min(aCount + bCount, gid + NV) - range.y; // b1
}
return range;
}
////////////////////////////////////////////////////////////////////////////////
// CTA mergesort support
template<int NT, int VT, typename T, typename Comp>
MGPU_DEVICE void CTABlocksortPass(T* keys_shared, int tid, int count,
int coop, T* keys, int* indices, Comp comp) {
int list = ~(coop - 1) & tid;
int diag = min(count, VT * ((coop - 1) & tid));
int start = VT * list;
int a0 = min(count, start);
int b0 = min(count, start + VT * (coop / 2));
int b1 = min(count, start + VT * coop);
int p = MergePath<MgpuBoundsLower>(keys_shared + a0, b0 - a0,
keys_shared + b0, b1 - b0, diag, comp);
SerialMerge<VT, true>(keys_shared, a0 + p, b0, b0 + diag - p, b1, keys,
indices, comp);
}
template<int NT, int VT, bool HasValues, typename KeyType, typename ValType,
typename Comp>
MGPU_DEVICE void CTABlocksortLoop(ValType threadValues[VT],
KeyType* keys_shared, ValType* values_shared, int tid, int count,
Comp comp) {
#pragma unroll
for(int coop = 2; coop <= NT; coop *= 2) {
int indices[VT];
KeyType keys[VT];
CTABlocksortPass<NT, VT>(keys_shared, tid, count, coop, keys,
indices, comp);
if(HasValues) {
// Exchange the values through shared memory.
DeviceThreadToShared<VT>(threadValues, tid, values_shared);
DeviceGather<NT, VT>(NT * VT, values_shared, indices, tid,
threadValues);
}
// Store results in shared memory in sorted order.
DeviceThreadToShared<VT>(keys, tid, keys_shared);
}
}
////////////////////////////////////////////////////////////////////////////////
// CTAMergesort
// Caller provides the keys in shared memory. This functions sorts the first
// count elements.
template<int NT, int VT, bool Stable, bool HasValues, typename KeyType,
typename ValType, typename Comp>
MGPU_DEVICE void CTAMergesort(KeyType threadKeys[VT], ValType threadValues[VT],
KeyType* keys_shared, ValType* values_shared, int count, int tid,
Comp comp) {
// Stable sort the keys in the thread.
if(VT * tid < count) {
if(Stable)
OddEvenTransposeSort<VT>(threadKeys, threadValues, comp);
else
OddEvenMergesort<VT>(threadKeys, threadValues, comp);
}
// Store the locally sorted keys into shared memory.
DeviceThreadToShared<VT>(threadKeys, tid, keys_shared);
// Recursively merge lists until the entire CTA is sorted.
CTABlocksortLoop<NT, VT, HasValues>(threadValues, keys_shared,
values_shared, tid, count, comp);
}
template<int NT, int VT, bool Stable, typename KeyType, typename Comp>
MGPU_DEVICE void CTAMergesortKeys(KeyType threadKeys[VT],
KeyType* keys_shared, int count, int tid, Comp comp) {
int valuesTemp[VT];
CTAMergesort<NT, VT, Stable, false>(threadKeys, valuesTemp, keys_shared,
(int*)keys_shared, count, tid, comp);
}
template<int NT, int VT, bool Stable, typename KeyType, typename ValType,
typename Comp>
MGPU_DEVICE void CTAMergesortPairs(KeyType threadKeys[VT],
ValType threadValues[VT], KeyType* keys_shared, ValType* values_shared,
int count, int tid, Comp comp) {
CTAMergesort<NT, VT, Stable, true>(threadKeys, threadValues, keys_shared,
values_shared, count, tid, comp);
}
////////////////////////////////////////////////////////////////////////////////
// DeviceMergeKeysIndices
template<int NT, int VT, bool LoadExtended, typename It1, typename It2,
typename T, typename Comp>
MGPU_DEVICE void DeviceMergeKeysIndices(It1 a_global, int aCount, It2 b_global,
int bCount, int4 range, int tid, T* keys_shared, T* results, int* indices,
Comp comp) {
int a0 = range.x;
int a1 = range.y;
int b0 = range.z;
int b1 = range.w;
if(LoadExtended) {
bool extended = (a1 < aCount) && (b1 < bCount);
aCount = a1 - a0;
bCount = b1 - b0;
int aCount2 = aCount + (int)extended;
int bCount2 = bCount + (int)extended;
// Load one element past the end of each input to avoid having to use
// range checking in the merge loop.
DeviceLoad2ToShared<NT, VT, VT + 1>(a_global + a0, aCount2,
b_global + b0, bCount2, tid, keys_shared);
// Run a Merge Path search for each thread's starting point.
int diag = VT * tid;
int mp = MergePath<MgpuBoundsLower>(keys_shared, aCount,
keys_shared + aCount2, bCount, diag, comp);
// Compute the ranges of the sources in shared memory.
int a0tid = mp;
int b0tid = aCount2 + diag - mp;
if(extended) {
SerialMerge<VT, false>(keys_shared, a0tid, 0, b0tid, 0, results,
indices, comp);
} else {
int a1tid = aCount;
int b1tid = aCount2 + bCount;
SerialMerge<VT, true>(keys_shared, a0tid, a1tid, b0tid, b1tid,
results, indices, comp);
}
} else {
// Use the input intervals from the ranges between the merge path
// intersections.
aCount = a1 - a0;
bCount = b1 - b0;
// Load the data into shared memory.
DeviceLoad2ToShared<NT, VT, VT>(a_global + a0, aCount, b_global + b0,
bCount, tid, keys_shared);
// Run a merge path to find the start of the serial merge for each
// thread.
int diag = VT * tid;
int mp = MergePath<MgpuBoundsLower>(keys_shared, aCount,
keys_shared + aCount, bCount, diag, comp);
// Compute the ranges of the sources in shared memory.
int a0tid = mp;
int a1tid = aCount;
int b0tid = aCount + diag - mp;
int b1tid = aCount + bCount;
// Serial merge into register.
SerialMerge<VT, true>(keys_shared, a0tid, a1tid, b0tid, b1tid, results,
indices, comp);
}
}
////////////////////////////////////////////////////////////////////////////////
// DeviceMerge
// Merge pairs from global memory into global memory. Useful factorization to
// enable calling from merge, mergesort, and locality sort.
template<int NT, int VT, bool HasValues, bool LoadExtended, typename KeysIt1,
typename KeysIt2, typename KeysIt3, typename ValsIt1, typename ValsIt2,
typename KeyType, typename ValsIt3, typename Comp>
MGPU_DEVICE void DeviceMerge(KeysIt1 aKeys_global, ValsIt1 aVals_global,
int aCount, KeysIt2 bKeys_global, ValsIt2 bVals_global, int bCount,
int tid, int block, int4 range, KeyType* keys_shared, int* indices_shared,
KeysIt3 keys_global, ValsIt3 vals_global, Comp comp) {
KeyType results[VT];
int indices[VT];
DeviceMergeKeysIndices<NT, VT, LoadExtended>(aKeys_global, aCount,
bKeys_global, bCount, range, tid, keys_shared, results, indices, comp);
// Store merge results back to shared memory.
DeviceThreadToShared<VT>(results, tid, keys_shared);
// Store merged keys to global memory.
aCount = range.y - range.x;
bCount = range.w - range.z;
DeviceSharedToGlobal<NT, VT>(aCount + bCount, keys_shared, tid,
keys_global + NT * VT * block);
// Copy the values.
if(HasValues) {
DeviceThreadToShared<VT>(indices, tid, indices_shared);
DeviceTransferMergeValuesShared<NT, VT>(aCount + bCount,
aVals_global + range.x, bVals_global + range.z, aCount,
indices_shared, tid, vals_global + NT * VT * block);
}
}
} // namespace mgpu
include/contrib/moderngpu/include/device/ctascan.cuh 0 → 100644
View file @ a60ed4d7
/******************************************************************************
* Copyright (c) 2013, NVIDIA CORPORATION. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
* * Neither the name of the NVIDIA CORPORATION nor the
* names of its contributors may be used to endorse or promote products
* derived from this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
******************************************************************************/
/******************************************************************************
*
* Code and text by Sean Baxter, NVIDIA Research
* See http://nvlabs.github.io/moderngpu for repository and documentation.
*
******************************************************************************/
#pragma once
#include "../mgpuenums.h"
#include "deviceutil.cuh"
#include "intrinsics.cuh"
namespace mgpu {
////////////////////////////////////////////////////////////////////////////////
// CTAReduce
template<int NT, typename Op = mgpu::plus<int> >
struct CTAReduce {
typedef typename Op::first_argument_type T;
enum { Size = NT, Capacity = NT };
struct Storage { T shared[Capacity]; };
MGPU_DEVICE static T Reduce(int tid, T x, Storage& storage, Op op = Op()) {
storage.shared[tid] = x;
__syncthreads();
// Fold the data in half with each pass.
#pragma unroll
for(int destCount = NT / 2; destCount >= 1; destCount /= 2) {
if(tid < destCount) {
// Read from the right half and store to the left half.
x = op(x, storage.shared[destCount + tid]);
storage.shared[tid] = x;
}
__syncthreads();
}
T total = storage.shared[0];
__syncthreads();
return total;
}
};
#if __CUDA_ARCH__ >= 300
template<int NT>
struct CTAReduce<NT, mgpu::plus<int> > {
typedef mgpu::plus<int> Op;
typedef int T;
enum { Size = NT, Capacity = WARP_SIZE };
struct Storage { int shared[Capacity]; };
MGPU_DEVICE static int Reduce(int tid, int x, Storage& storage,
Op op = Op()) {
const int NumSections = WARP_SIZE;
const int SecSize = NT / NumSections;
int lane = (SecSize - 1) & tid;
int sec = tid / SecSize;
// In the first phase, threads cooperatively find the reduction within
// their segment. The segments are SecSize threads (NT / WARP_SIZE)
// wide.
#pragma unroll
for(int offset = 1; offset < SecSize; offset *= 2)
x = shfl_add(x, offset, SecSize);
// The last thread in each segment stores the local reduction to shared
// memory.
if(SecSize - 1 == lane) storage.shared[sec] = x;
__syncthreads();
// Reduce the totals of each input segment. The spine is WARP_SIZE
// threads wide.
if(tid < NumSections) {
x = storage.shared[tid];
#pragma unroll
for(int offset = 1; offset < NumSections; offset *= 2)
x = shfl_add(x, offset, NumSections);
storage.shared[tid] = x;
}
__syncthreads();
int reduction = storage.shared[NumSections - 1];
__syncthreads();
return reduction;
}
};
template<int NT>
struct CTAReduce<NT, mgpu::maximum<int> > {
typedef mgpu::maximum<int> Op;
enum { Size = NT, Capacity = WARP_SIZE };
struct Storage { int shared[Capacity]; };
MGPU_DEVICE static int Reduce(int tid, int x, Storage& storage,
Op op = Op()) {
const int NumSections = WARP_SIZE;
const int SecSize = NT / NumSections;
int lane = (SecSize - 1) & tid;
int sec = tid / SecSize;
#pragma unroll
for(int offset = 1; offset < SecSize; offset *= 2)
x = shfl_max(x, offset, SecSize);
if(SecSize - 1 == lane) storage.shared[sec] = x;
__syncthreads();
if(tid < NumSections) {
x = storage.shared[tid];
#pragma unroll
for(int offset = 1; offset < NumSections; offset *= 2)
x = shfl_max(x, offset, NumSections);
storage.shared[tid] = x;
}
__syncthreads();
int reduction = storage.shared[NumSections - 1];
__syncthreads();
return reduction;
}
};
#endif // __CUDA_ARCH__ >= 300
////////////////////////////////////////////////////////////////////////////////
// CTAScan
template<int NT, typename Op = mgpu::plus<int> >
struct CTAScan {
typedef typename Op::result_type T;
enum { Size = NT, Capacity = 2 * NT + 1 };
struct Storage { T shared[Capacity]; };
MGPU_DEVICE static T Scan(int tid, T x, Storage& storage, T* total,
MgpuScanType type = MgpuScanTypeExc, T identity = (T)0, Op op = Op()) {
storage.shared[tid] = x;
int first = 0;
__syncthreads();
#pragma unroll
for(int offset = 1; offset < NT; offset += offset) {
if(tid >= offset)
x = op(storage.shared[first + tid - offset], x);
first = NT - first;
storage.shared[first + tid] = x;
__syncthreads();
}
*total = storage.shared[first + NT - 1];
if(MgpuScanTypeExc == type)
x = tid ? storage.shared[first + tid - 1] : identity;
__syncthreads();
return x;
}
MGPU_DEVICE static T Scan(int tid, T x, Storage& storage) {
T total;
return Scan(tid, x, storage, &total, MgpuScanTypeExc, (T)0, Op());
}
};
////////////////////////////////////////////////////////////////////////////////
// Special partial specialization for CTAScan<NT, ScanOpAdd> on Kepler.
// This uses the shfl intrinsic to reduce scan latency.
#if __CUDA_ARCH__ >= 300
template<int NT>
struct CTAScan<NT, mgpu::plus<int> > {
typedef mgpu::plus<int> Op;
enum { Size = NT, NumSegments = WARP_SIZE, SegSize = NT / NumSegments };
enum { Capacity = NumSegments + 1 };
struct Storage { int shared[Capacity + 1]; };
MGPU_DEVICE static int Scan(int tid, int x, Storage& storage, int* total,
MgpuScanType type = MgpuScanTypeExc, int identity = 0, Op op = Op()) {
// Define WARP_SIZE segments that are NT / WARP_SIZE large.
// Each warp makes log(SegSize) shfl_add calls.
// The spine makes log(WARP_SIZE) shfl_add calls.
int lane = (SegSize - 1) & tid;
int segment = tid / SegSize;
// Scan each segment using shfl_add.
int scan = x;
#pragma unroll
for(int offset = 1; offset < SegSize; offset *= 2)
scan = shfl_add(scan, offset, SegSize);
// Store the reduction (last element) of each segment into storage.
if(SegSize - 1 == lane) storage.shared[segment] = scan;
__syncthreads();
// Warp 0 does a full shfl warp scan on the partials. The total is
// stored to shared[NumSegments]. (NumSegments = WARP_SIZE)
if(tid < NumSegments) {
int y = storage.shared[tid];
int scan = y;
#pragma unroll
for(int offset = 1; offset < NumSegments; offset *= 2)
scan = shfl_add(scan, offset, NumSegments);
storage.shared[tid] = scan - y;
if(NumSegments - 1 == tid) storage.shared[NumSegments] = scan;
}
__syncthreads();
// Add the scanned partials back in and convert to exclusive scan.
scan += storage.shared[segment];
if(MgpuScanTypeExc == type) {
scan -= x;
if(identity && !tid) scan = identity;
}
*total = storage.shared[NumSegments];
__syncthreads();
return scan;
}
MGPU_DEVICE static int Scan(int tid, int x, Storage& storage) {
int total;
return Scan(tid, x, storage, &total, MgpuScanTypeExc, 0);
}
};
#endif // __CUDA_ARCH__ >= 300
////////////////////////////////////////////////////////////////////////////////
// CTABinaryScan
template<int NT>
MGPU_DEVICE int CTABinaryScan(int tid, bool x, int* shared, int* total) {
const int NumWarps = NT / WARP_SIZE;
int warp = tid / WARP_SIZE;
int lane = (WARP_SIZE - 1);
// Store the bit totals for each warp.
uint bits = __ballot(x);
shared[warp] = popc(bits);
__syncthreads();
#if __CUDA_ARCH__ >= 300
if(tid < NumWarps) {
int x = shared[tid];
int scan = x;
#pragma unroll
for(int offset = 1; offset < NumWarps; offset *= 2)
scan = shfl_add(scan, offset, NumWarps);
shared[tid] = scan - x;
}
__syncthreads();
#else
// Thread 0 scans warp totals.
if(!tid) {
int scan = 0;
#pragma unroll
for(int i = 0; i < NumWarps; ++i) {
int y = shared[i];
shared[i] = scan;
scan += y;
}
shared[NumWarps] = scan;
}
__syncthreads();
#endif // __CUDA_ARCH__ >= 300
// Add the warp scan back into the partials.
int scan = shared[warp] + __popc(bfe(bits, 0, lane));
*total = shared[NumWarps];
__syncthreads();
return scan;
}
} // namespace mgpu
include/contrib/moderngpu/include/device/ctasearch.cuh 0 → 100644
View file @ a60ed4d7
/******************************************************************************
* Copyright (c) 2013, NVIDIA CORPORATION. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
* * Neither the name of the NVIDIA CORPORATION nor the
* names of its contributors may be used to endorse or promote products
* derived from this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
******************************************************************************/
/******************************************************************************
*
* Code and text by Sean Baxter, NVIDIA Research
* See http://nvlabs.github.io/moderngpu for repository and documentation.
*
******************************************************************************/
#pragma once
#include "deviceutil.cuh"
#include "../mgpudevice.cuh"
namespace mgpu {
template<MgpuBounds Bounds, typename IntT, typename It, typename T,
typename Comp>
MGPU_HOST_DEVICE void BinarySearchIt(It data, int& begin, int& end, T key,
int shift, Comp comp) {
IntT scale = (1<< shift) - 1;
int mid = (int)((begin + scale * end)>> shift);
T key2 = data[mid];
bool pred = (MgpuBoundsUpper == Bounds) ?
!comp(key, key2) :
comp(key2, key);
if(pred) begin = mid + 1;
else end = mid;
}
template<MgpuBounds Bounds, typename IntT, typename T, typename It,
typename Comp>
MGPU_HOST_DEVICE int BiasedBinarySearch(It data, int count, T key, int levels,
Comp comp) {
int begin = 0;
int end = count;
if(levels >= 4 && begin < end)
BinarySearchIt<Bounds, IntT>(data, begin, end, key, 9, comp);
if(levels >= 3 && begin < end)
BinarySearchIt<Bounds, IntT>(data, begin, end, key, 7, comp);
if(levels >= 2 && begin < end)
BinarySearchIt<Bounds, IntT>(data, begin, end, key, 5, comp);
if(levels >= 1 && begin < end)
BinarySearchIt<Bounds, IntT>(data, begin, end, key, 4, comp);
while(begin < end)
BinarySearchIt<Bounds, int>(data, begin, end, key, 1, comp);
return begin;
}
template<MgpuBounds Bounds, typename T, typename It, typename Comp>
MGPU_HOST_DEVICE int BinarySearch(It data, int count, T key, Comp comp) {
int begin = 0;
int end = count;
while(begin < end)
BinarySearchIt<Bounds, int>(data, begin, end, key, 1, comp);
return begin;
}
////////////////////////////////////////////////////////////////////////////////
// MergePath search
template<MgpuBounds Bounds, typename It1, typename It2, typename Comp>
MGPU_HOST_DEVICE int MergePath(It1 a, int aCount, It2 b, int bCount, int diag,
Comp comp) {
typedef typename std::iterator_traits<It1>::value_type T;
int begin = max(0, diag - bCount);
int end = min(diag, aCount);
while(begin < end) {
int mid = (begin + end)>> 1;
T aKey = a[mid];
T bKey = b[diag - 1 - mid];
bool pred = (MgpuBoundsUpper == Bounds) ?
comp(aKey, bKey) :
!comp(bKey, aKey);
if(pred) begin = mid + 1;
else end = mid;
}
return begin;
}
////////////////////////////////////////////////////////////////////////////////
// SegmentedMergePath search
template<typename InputIt, typename Comp>
MGPU_HOST_DEVICE int SegmentedMergePath(InputIt keys, int aOffset, int aCount,
int bOffset, int bCount, int leftEnd, int rightStart, int diag, Comp comp) {
// leftEnd and rightStart are defined from the origin, and diag is defined
// from aOffset.
// We only need to run a Merge Path search if the diagonal intersects the
// segment that strides the left and right halves (i.e. is between leftEnd
// and rightStart).
if(aOffset + diag <= leftEnd) return diag;
if(aOffset + diag >= rightStart) return aCount;
bCount = min(bCount, rightStart - bOffset);
int begin = max(max(leftEnd - aOffset, 0), diag - bCount);
int end = min(diag, aCount);
while(begin < end) {
int mid = (begin + end)>> 1;
int ai = aOffset + mid;
int bi = bOffset + diag - 1 - mid;
bool pred = !comp(keys[bi], keys[ai]);
if(pred) begin = mid + 1;
else end = mid;
}
return begin;
}
////////////////////////////////////////////////////////////////////////////////
// BalancedPath search
template<bool Duplicates, typename IntT, typename InputIt1, typename InputIt2,
typename Comp>
MGPU_HOST_DEVICE int2 BalancedPath(InputIt1 a, int aCount, InputIt2 b,
int bCount, int diag, int levels, Comp comp) {
typedef typename std::iterator_traits<InputIt1>::value_type T;
int p = MergePath<MgpuBoundsLower>(a, aCount, b, bCount, diag, comp);
int aIndex = p;
int bIndex = diag - p;
bool star = false;
if(bIndex < bCount) {
if(Duplicates) {
T x = b[bIndex];
// Search for the beginning of the duplicate run in both A and B.
// Because
int aStart = BiasedBinarySearch<MgpuBoundsLower, IntT>(a, aIndex, x,
levels, comp);
int bStart = BiasedBinarySearch<MgpuBoundsLower, IntT>(b, bIndex, x,
levels, comp);
// The distance between the merge path and the lower_bound is the
// 'run'. We add up the a- and b- runs and evenly distribute them to
// get a stairstep path.
int aRun = aIndex - aStart;
int bRun = bIndex - bStart;
int xCount = aRun + bRun;
// Attempt to advance b and regress a.
int bAdvance = max(xCount>> 1, bRun);
int bEnd = min(bCount, bStart + bAdvance + 1);
int bRunEnd = BinarySearch<MgpuBoundsUpper>(b + bIndex,
bEnd - bIndex, x, comp) + bIndex;
bRun = bRunEnd - bStart;
bAdvance = min(bAdvance, bRun);
int aAdvance = xCount - bAdvance;
bool roundUp = (aAdvance == bAdvance + 1) && (bAdvance < bRun);
aIndex = aStart + aAdvance;
if(roundUp) star = true;
} else {
if(aIndex && aCount) {
T aKey = a[aIndex - 1];
T bKey = b[bIndex];
// If the last consumed element in A (aIndex - 1) is the same as
// the next element in B (bIndex), we're sitting at a starred
// partition.
if(!comp(aKey, bKey)) star = true;
}
}
}
return make_int2(aIndex, star);
}
} // namespace mgpu
include/contrib/moderngpu/include/device/ctasegreduce.cuh 0 → 100644
View file @ a60ed4d7
/******************************************************************************
* Copyright (c) 2013, NVIDIA CORPORATION. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
* * Neither the name of the NVIDIA CORPORATION nor the
* names of its contributors may be used to endorse or promote products
* derived from this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
******************************************************************************/
/******************************************************************************
*
* Code and text by Sean Baxter, NVIDIA Research
* See http://nvlabs.github.io/moderngpu for repository and documentation.
*
******************************************************************************/
#pragma once
#include "ctasegscan.cuh"
#include "ctasearch.cuh"
namespace mgpu {
////////////////////////////////////////////////////////////////////////////////
// Segmented reduce utility functions.
// Extract the upper-bound indices from the coded ranges. Decrement to include
// the first addressed row/segment.
struct SegReduceRange {
int begin;
int end;
int total;
bool flushLast;
};
MGPU_DEVICE SegReduceRange DeviceShiftRange(int limit0, int limit1) {
SegReduceRange range;
range.begin = 0x7fffffff & limit0;
range.end = 0x7fffffff & limit1;
range.total = range.end - range.begin;
range.flushLast = 0 == (0x80000000 & limit1);
range.end += !range.flushLast;
return range;
}
// Reconstitute row/segment indices from a starting row index and packed end
// flags. Used for pre-processed versions of interval reduce and interval Spmv.
template<int VT>
MGPU_DEVICE void DeviceExpandFlagsToRows(int first, int endFlags,
int rows[VT + 1]) {
rows[0] = first;
#pragma unroll
for(int i = 0; i < VT; ++i) {
if((1<< i) & endFlags) ++first;
rows[i + 1] = first;
}
}
////////////////////////////////////////////////////////////////////////////////
// After loading CSR terms into shared memory, each thread binary searches
// (upper-bound) to find its starting point. Each thread then walks forward,
// emitting the csr0-relative row indices to register.
template<int NT, int VT>
MGPU_DEVICE int DeviceExpandCsrRows(int tidOffset, int* csr_shared,
int numRows, int end, int rows[VT + 1], int rowStarts[VT]) {
// Each thread binary searches for its starting row.
int row = BinarySearch<MgpuBoundsUpper>(csr_shared, numRows, tidOffset,
mgpu::less<int>()) - 1;
// Each thread starts at row and scans forward, emitting row IDs into
// register. Store the CTA-local row index (starts at 0) to rows and the
// start of the row (globally) to rowStarts.
int curOffset = csr_shared[row];
int nextOffset = (row + 1 < numRows) ? csr_shared[row + 1] : end;
rows[0] = row;
rowStarts[0] = curOffset;
int endFlags = 0;
#pragma unroll
for(int i = 1; i <= VT; ++i) {
// Advance the row cursor when the iterator hits the next row offset.
if(tidOffset + i == nextOffset) {
// Set an end flag when the cursor advances to the next row.
endFlags |= 1<< (i - 1);
// Advance the cursor and load the next row offset.
++row;
curOffset = nextOffset;
nextOffset = (row + 1 < numRows) ? csr_shared[row + 1] : end;
}
rows[i] = row;
if(i < VT) rowStarts[i] = curOffset;
}
__syncthreads();
return endFlags;
}
////////////////////////////////////////////////////////////////////////////////
// DeviceSegReducePrepare
// Expand non-empty interval of CSR elements into row indices. Compute end-flags
// by comparing adjacent row IDs.
// DeviceSegReducePrepare may be called either by a pre-processing kernel or by
// the kernel that actually evaluates the segmented reduction if no preprocesing
// is desired.
struct SegReduceTerms {
int endFlags;
int tidDelta;
};
template<int NT, int VT>
MGPU_DEVICE SegReduceTerms DeviceSegReducePrepare(int* csr_shared, int numRows,
int tid, int gid, bool flushLast, int rows[VT + 1], int rowStarts[VT]) {
// Pass a sentinel (end) to point to the next segment start. If we flush,
// this is the end of this tile. Otherwise it is INT_MAX
int endFlags = DeviceExpandCsrRows<NT, VT>(gid + VT * tid, csr_shared,
numRows, flushLast ? (gid + NT * VT) : INT_MAX, rows, rowStarts);
// Find the distance to to scan to compute carry-in for each thread. Use the
// existance of an end flag anywhere in the thread to determine if carry-out
// values from the left should propagate through to the right.
int tidDelta = DeviceFindSegScanDelta<NT>(tid, rows[0] != rows[VT],
csr_shared);
SegReduceTerms terms = { endFlags, tidDelta };
return terms;
}
////////////////////////////////////////////////////////////////////////////////
// CTASegReduce
// Core segmented reduction code. Supports fast-path and slow-path for intra-CTA
// segmented reduction. Stores partials to global memory.
// Callers feed CTASegReduce::ReduceToGlobal values in thread order.
template<int NT, int VT, bool HalfCapacity, typename T, typename Op>
struct CTASegReduce {
typedef CTASegScan<NT, Op> SegScan;
enum {
NV = NT * VT,
Capacity = HalfCapacity ? (NV / 2) : NV
};
union Storage {
typename SegScan::Storage segScanStorage;
T values[Capacity];
};
template<typename DestIt>
MGPU_DEVICE static void ReduceToGlobal(const int rows[VT + 1], int total,
int tidDelta, int startRow, int block, int tid, T data[VT],
DestIt dest_global, T* carryOut_global, T identity, Op op,
Storage& storage) {
// Run a segmented scan within the thread.
T x, localScan[VT];
#pragma unroll
for(int i = 0; i < VT; ++i) {
x = i ? op(x, data[i]) : data[i];
localScan[i] = x;
if(rows[i] != rows[i + 1]) x = identity;
}
// Run a parallel segmented scan over the carry-out values to compute
// carry-in.
T carryOut;
T carryIn = SegScan::SegScanDelta(tid, tidDelta, x,
storage.segScanStorage, &carryOut, identity, op);
// Store the carry-out for the entire CTA to global memory.
if(!tid) carryOut_global[block] = carryOut;
dest_global += startRow;
if(HalfCapacity && total > Capacity) {
// Add carry-in to each thread-local scan value. Store directly
// to global.
#pragma unroll
for(int i = 0; i < VT; ++i) {
// Add the carry-in to the local scan.
T x2 = op(carryIn, localScan[i]);
// Store on the end flag and clear the carry-in.
if(rows[i] != rows[i + 1]) {
carryIn = identity;
dest_global[rows[i]] = x2;
}
}
} else {
// All partials fit in shared memory. Add carry-in to each thread-
// local scan value.
#pragma unroll
for(int i = 0; i < VT; ++i) {
// Add the carry-in to the local scan.
T x2 = op(carryIn, localScan[i]);
// Store reduction when the segment changes and clear the
// carry-in.
if(rows[i] != rows[i + 1]) {
storage.values[rows[i]] = x2;
carryIn = identity;
}
}
__syncthreads();
// Cooperatively store reductions to global memory.
for(int index = tid; index < total; index += NT)
dest_global[index] = storage.values[index];
__syncthreads();
}
}
};
} // namespace mgpu
include/contrib/moderngpu/include/device/ctasegscan.cuh 0 → 100644
View file @ a60ed4d7
/******************************************************************************
* Copyright (c) 2013, NVIDIA CORPORATION. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
* * Neither the name of the NVIDIA CORPORATION nor the
* names of its contributors may be used to endorse or promote products
* derived from this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
******************************************************************************/
/******************************************************************************
*
* Code and text by Sean Baxter, NVIDIA Research
* See http://nvlabs.github.io/moderngpu for repository and documentation.
*
******************************************************************************/
#pragma once
#include "ctascan.cuh"
namespace mgpu {
////////////////////////////////////////////////////////////////////////////////
// DeviceFindSegScanDelta
// Runs an inclusive max-index scan over binary inputs.
template<int NT>
MGPU_DEVICE int DeviceFindSegScanDelta(int tid, bool flag, int* delta_shared) {
const int NumWarps = NT / 32;
int warp = tid / 32;
int lane = 31 & tid;
uint warpMask = 0xffffffff>> (31 - lane); // inclusive search
uint ctaMask = 0x7fffffff>> (31 - lane); // exclusive search
uint warpBits = __ballot(flag);
delta_shared[warp] = warpBits;
__syncthreads();
if(tid < NumWarps) {
uint ctaBits = __ballot(0 != delta_shared[tid]);
int warpSegment = 31 - clz(ctaMask & ctaBits);
int start = (-1 != warpSegment) ?
(31 - clz(delta_shared[warpSegment]) + 32 * warpSegment) : 0;
delta_shared[NumWarps + tid] = start;
}
__syncthreads();
// Find the closest flag to the left of this thread within the warp.
// Include the flag for this thread.
int start = 31 - clz(warpMask & warpBits);
if(-1 != start) start += ~31 & tid;
else start = delta_shared[NumWarps + warp];
__syncthreads();
return tid - start;
}
////////////////////////////////////////////////////////////////////////////////
// CTASegScan
template<int NT, typename _Op = mgpu::plus<int> >
struct CTASegScan {
typedef _Op Op;
typedef typename Op::result_type T;
enum { NumWarps = NT / 32, Size = NT, Capacity = 2 * NT };
union Storage {
int delta[NumWarps];
T values[Capacity];
};
// Each thread passes the reduction of the LAST SEGMENT that it covers.
// flag is set to true if there's at least one segment flag in the thread.
// SegScan returns the reduction of values for the first segment in this
// thread over the preceding threads.
// Return the value init for the first thread.
// When scanning single elements per thread, interpret the flag as a BEGIN
// FLAG. If tid's flag is set, its value belongs to thread tid + 1, not
// thread tid.
// The function returns the reduction of the last segment in the CTA.
MGPU_DEVICE static T SegScanDelta(int tid, int tidDelta, T x,
Storage& storage, T* carryOut, T identity = (T)0, Op op = Op()) {
// Run an inclusive scan
int first = 0;
storage.values[first + tid] = x;
__syncthreads();
#pragma unroll
for(int offset = 1; offset < NT; offset += offset) {
if(tidDelta >= offset)
x = op(storage.values[first + tid - offset], x);
first = NT - first;
storage.values[first + tid] = x;
__syncthreads();
}
// Get the exclusive scan.
x = tid ? storage.values[first + tid - 1] : identity;
*carryOut = storage.values[first + NT - 1];
__syncthreads();
return x;
}
MGPU_DEVICE static T SegScan(int tid, T x, bool flag, Storage& storage,
T* carryOut, T identity = (T)0, Op op = Op()) {
// Find the left-most thread that covers the first segment of this
// thread.
int tidDelta = DeviceFindSegScanDelta<NT>(tid, flag, storage.delta);
return SegScanDelta(tid, tidDelta, x, storage, carryOut, identity, op);
}
};
} // namespace mgpu
include/contrib/moderngpu/include/device/ctasegsort.cuh 0 → 100644
View file @ a60ed4d7
/******************************************************************************
* Copyright (c) 2013, NVIDIA CORPORATION. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
* * Neither the name of the NVIDIA CORPORATION nor the
* names of its contributors may be used to endorse or promote products
* derived from this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
******************************************************************************/
/******************************************************************************
*
* Code and text by Sean Baxter, NVIDIA Research
* See http://nvlabs.github.io/moderngpu for repository and documentation.
*
******************************************************************************/
#pragma once
#include "ctascan.cuh"
#include "ctasearch.cuh"
#include "loadstore.cuh"
#include "sortnetwork.cuh"
namespace mgpu {
template<int VT, typename T, typename Comp>
MGPU_DEVICE void SegmentedSerialMerge(const T* keys_shared, int aBegin,
int aEnd, int bBegin, int bEnd, T results[VT], int indices[VT],
int leftEnd, int rightStart, Comp comp, bool sync = true) {
bEnd = min(rightStart, bEnd);
T aKey = keys_shared[aBegin];
T bKey = keys_shared[bBegin];
#pragma unroll
for(int i = 0; i < VT; ++i) {
bool p;
// If A has run out of inputs, emit B.
if(aBegin >= aEnd)
p = false;
else if(bBegin >= bEnd || aBegin < leftEnd)
// B has hit the end of the middle segment.
// Emit A if A has inputs remaining in the middle segment.
p = true;
else
// Emit the smaller element in the middle segment.
p = !comp(bKey, aKey);
results[i] = p ? aKey : bKey;
indices[i] = p ? aBegin : bBegin;
if(p) aKey = keys_shared[++aBegin];
else bKey = keys_shared[++bBegin];
}
if(sync) { __syncthreads(); }
}
////////////////////////////////////////////////////////////////////////////////
// CTASegsortPass
template<int NT, int VT, typename T, typename Comp>
MGPU_DEVICE void CTASegsortPass(T* keys_shared, int* ranges_shared, int tid,
int pass, T results[VT], int indices[VT], int2& activeRange, Comp comp) {
// Locate the intervals of the input lists.
int3 frame = FindMergesortFrame(2<< pass, tid, VT);
int a0 = frame.x;
int b0 = frame.y;
int listLen = frame.z;
int list = tid>> pass;
int listParity = 1 & list;
int diag = VT * tid - frame.x;
// Fetch the active range for the list this thread's list is merging with.
int siblingRange = ranges_shared[1 ^ list];
int siblingStart = 0x0000ffff & siblingRange;
int siblingEnd = siblingRange>> 16;
// Create a new active range for the merge.
int leftEnd = listParity ? siblingEnd : activeRange.y;
int rightStart = listParity ? activeRange.x : siblingStart;
activeRange.x = min(activeRange.x, siblingStart);
activeRange.y = max(activeRange.y, siblingEnd);
int p = SegmentedMergePath(keys_shared, a0, listLen, b0, listLen, leftEnd,
rightStart, diag, comp);
int a0tid = a0 + p;
int b0tid = b0 + diag - p;
SegmentedSerialMerge<VT>(keys_shared, a0tid, b0, b0tid, b0 + listLen,
results, indices, leftEnd, rightStart, comp);
// Store the ranges to shared memory.
if(0 == diag)
ranges_shared[list>> 1] =
(int)bfi(activeRange.y, activeRange.x, 16, 16);
}
////////////////////////////////////////////////////////////////////////////////
// CTASegsortLoop
template<int NT, int VT, bool HasValues, typename KeyType, typename ValType,
typename Comp>
MGPU_DEVICE int2 CTASegsortLoop(KeyType threadKeys[VT],
ValType threadValues[VT], KeyType* keys_shared, ValType* values_shared,
int* ranges_shared, int tid, int2 activeRange, Comp comp) {
const int NumPasses = sLogPow2<NT>::value;
#pragma unroll
for(int pass = 0; pass < NumPasses; ++pass) {
int indices[VT];
CTASegsortPass<NT, VT>(keys_shared, ranges_shared, tid, pass,
threadKeys, indices, activeRange, comp);
if(HasValues) {
// Exchange values through shared memory.
DeviceThreadToShared<VT>(threadValues, tid, values_shared);
DeviceGather<NT, VT>(NT * VT, values_shared, indices, tid,
threadValues);
}
// Store results in shared memory in sorted order.
DeviceThreadToShared<VT>(threadKeys, tid, keys_shared);
}
return activeRange;
}
////////////////////////////////////////////////////////////////////////////////
// CTASegsort
// Pass keys and values in register. On return, values are returned in register
// and keys returned in shared memory.
template<int NT, int VT, bool Stable, bool HasValues, typename KeyType,
typename ValType, typename Comp>
MGPU_DEVICE int2 CTASegsort(KeyType threadKeys[VT], ValType threadValues[VT],
int tid, int headFlags, KeyType* keys_shared, ValType* values_shared,
int* ranges_shared, Comp comp) {
if(Stable)
// Odd-even transpose sort.
OddEvenTransposeSortFlags<VT>(threadKeys, threadValues, headFlags,
comp);
else
// Batcher's odd-even mergesort.
OddEvenMergesortFlags<VT>(threadKeys, threadValues, headFlags, comp);
// Record the first and last occurrence of head flags in this segment.
int blockEnd = 31 - clz(headFlags);
if(-1 != blockEnd) blockEnd += VT * tid;
int blockStart = ffs(headFlags);
blockStart = blockStart ? (VT * tid - 1 + blockStart) : (NT * VT);
ranges_shared[tid] = (int)bfi(blockEnd, blockStart, 16, 16);
// Store back to shared mem. The values are in VT-length sorted lists.
// These are merged recursively.
DeviceThreadToShared<VT>(threadKeys, tid, keys_shared);
int2 activeRange = CTASegsortLoop<NT, VT, HasValues>(threadKeys,
threadValues, keys_shared, values_shared, ranges_shared, tid,
make_int2(blockStart, blockEnd), comp);
return activeRange;
}
template<int NT, int VT, bool Stable, typename KeyType, typename Comp>
MGPU_DEVICE int2 CTASegsortKeys(KeyType threadKeys[VT], int tid, int headFlags,
KeyType* keys_shared, int* ranges_shared, Comp comp) {
int valuesTemp[VT];
return CTASegsort<NT, VT, Stable, false>(threadKeys, valuesTemp, tid,
headFlags, keys_shared, (int*)keys_shared, ranges_shared, comp);
}
template<int NT, int VT, bool Stable, typename KeyType, typename ValType,
typename Comp>
MGPU_DEVICE int2 CTASegsortPairs(KeyType threadKeys[VT],
ValType threadValues[VT], int tid, int headFlags, KeyType* keys_shared,
ValType* values_shared, int* ranges_shared, Comp comp) {
return CTASegsort<NT, VT, Stable, true>(threadKeys, threadValues, tid,
headFlags, keys_shared, values_shared, ranges_shared, comp);
}
////////////////////////////////////////////////////////////////////////////////
// DeviceSegBlocksort
// Load keys and values from global memory, sort in shared memory, and store
// back to global memory. Store the left-most and right-most encountered
// headflag locations to ranges_global to prepare for the next pass.
// This function is factored out of the blocksort kernel to allow easier
// customization of that kernel - we have two implementations currently:
// sort over indices and sort over bitfield.
template<int NT, int VT, bool Stable, bool HasValues, typename InputIt1,
typename InputIt2, typename KeyType, typename ValType, typename OutputIt1,
typename OutputIt2, typename Comp>
MGPU_DEVICE void DeviceSegBlocksort(InputIt1 keys_global,
InputIt2 values_global, int count2, KeyType* keys_shared,
ValType* values_shared, int* ranges_shared, int headFlags, int tid,
int block, OutputIt1 keysDest_global, OutputIt2 valsDest_global,
int* ranges_global, Comp comp) {
// Load keys into register in thread order.
int gid = NT * VT * block;
KeyType threadKeys[VT];
DeviceGlobalToShared<NT, VT>(count2, keys_global + gid, tid, keys_shared);
DeviceSharedToThread<VT>(keys_shared, tid, threadKeys);
// Load the values from global memory and into register in thread order.
ValType threadValues[VT];
if(HasValues) {
DeviceGlobalToShared<NT, VT>(count2, values_global + gid, tid,
values_shared);
DeviceSharedToThread<VT>(values_shared, tid, threadValues);
}
// Run the CTA segmented blocksort.
int2 activeRange = CTASegsort<NT, VT, Stable, HasValues>(threadKeys,
threadValues, tid, headFlags, keys_shared, values_shared, ranges_shared,
comp);
// Store the keys to global memory.
DeviceSharedToGlobal<NT, VT>(count2, keys_shared, tid,
keysDest_global + gid);
if(HasValues) {
// Store the values to global memory.xk b
DeviceThreadToShared<VT>(threadValues, tid, values_shared);
DeviceSharedToGlobal<NT, VT>(count2, values_shared, tid,
valsDest_global + gid, false);
}
// Store the 16-bit packed ranges. These are used by all merge kernels and
// the first level of global segmented merge path partitioning.
if(!tid)
ranges_global[block] = bfi(activeRange.y, activeRange.x, 16, 16);
}
////////////////////////////////////////////////////////////////////////////////
// DeviceIndicesToHeadFlags
// Load indices from an array and cooperatively turn into a head flag bitfield
// for each thread.
template<int NT, int VT>
MGPU_DEVICE int DeviceIndicesToHeadFlags(const int* indices_global,
const int* partitions_global, int tid, int block, int count2,
int* words_shared, byte* flags_shared) {
const int FlagWordsPerThread = MGPU_DIV_UP(VT, 4);
int gid = NT * VT * block;
int p0 = partitions_global[block];
int p1 = partitions_global[block + 1];
int headFlags = 0;
if(p1 > p0 || count2 < NT * VT) {
// Clear the flag bytes, then loop through the indices and poke in flag
// values.
#pragma unroll
for(int i = 0; i < FlagWordsPerThread; ++i)
words_shared[NT * i + tid] = 0;
__syncthreads();
for(int index = p0 + tid; index < p1; index += NT) {
int headFlag = indices_global[index];
flags_shared[headFlag - gid] = 1;
}
__syncthreads();
// Combine all the head flags for this thread.
int first = VT * tid;
int offset = first / 4;
int prev = words_shared[offset];
int mask = 0x3210 + 0x1111 * (3 & first);
#pragma unroll
for(int i = 0; i < FlagWordsPerThread; ++i) {
// Gather the next four flags.
int next = words_shared[offset + 1 + i];
int x = prmt(prev, next, mask);
prev = next;
// Set the head flag bits.
if(0x00000001 & x) headFlags |= 1<< (4 * i);
if(0x00000100 & x) headFlags |= 1<< (4 * i + 1);
if(0x00010000 & x) headFlags |= 1<< (4 * i + 2);
if(0x01000000 & x) headFlags |= 1<< (4 * i + 3);
}
__syncthreads();
// Set head flags for out-of-range keys.
int outOfRange = min(VT, first + VT - count2);
if(outOfRange > 0)
headFlags = bfi(0xffffffff, headFlags, VT - outOfRange, outOfRange);
// Clear head flags above VT.
headFlags &= (1<< VT) - 1;
}
return headFlags;
}
////////////////////////////////////////////////////////////////////////////////
// SegSortSupport
struct SegSortSupport {
int* ranges_global;
int2* ranges2_global;
int4* mergeList_global;
int* copyList_global;
int2* queueCounters_global;
int2* nextCounters_global;
byte* copyStatus_global;
};
////////////////////////////////////////////////////////////////////////////////
// DeviceSegSortMerge
template<int NT, int VT, bool HasValues, typename KeyType, typename ValueType,
typename Comp>
MGPU_DEVICE void DeviceSegSortMerge(const KeyType* keys_global,
const ValueType* values_global, int2 segmentRange, int tid,
int block, int4 range, int pass, KeyType* keys_shared,
int* indices_shared, KeyType* keysDest_global, ValueType* valsDest_global,
Comp comp) {
const int NV = NT * VT;
int gid = NV * block;
// Load the local compressed segment indices.
int a0 = range.x;
int aCount = range.y - range.x;
int b0 = range.z;
int bCount = range.w - range.z;
DeviceLoad2ToShared<NT, VT, VT>(keys_global + a0, aCount, keys_global + b0,
bCount, tid, keys_shared);
////////////////////////////////////////////////////////////////////////////
// Run a merge path to find the starting point for each thread to merge.
// If the entire warp fits into the already-sorted segments, we can skip
// sorting it and leave its keys in shared memory. Doing this on the warp
// level rather than thread level (also legal) gives slightly better
// performance.
int segStart = segmentRange.x;
int segEnd = segmentRange.y;
int listParity = 1 & (block>> pass);
int warpOffset = VT * (~31 & tid);
bool sortWarp = listParity ?
// The spliced segment is to the left (segStart).
(warpOffset < segStart) :
// The spliced segment is to the right (segEnd).
(warpOffset + 32 * VT > segEnd);
KeyType threadKeys[VT];
int indices[VT];
if(sortWarp) {
int diag = VT * tid;
int mp = SegmentedMergePath(keys_shared, 0, aCount, aCount, bCount,
listParity ? 0 : segEnd, listParity ? segStart : NV, diag, comp);
int a0tid = mp;
int a1tid = aCount;
int b0tid = aCount + diag - mp;
int b1tid = aCount + bCount;
// Serial merge into register. All threads in the CTA so we hoist the
// check for list parity outside the function call to simplify the
// logic. Unlike in the blocksort, this does not cause warp divergence.
SegmentedSerialMerge<VT>(keys_shared, a0tid, a1tid, b0tid, b1tid,
threadKeys, indices, listParity ? 0 : segEnd,
listParity ? segStart : NV, comp, false);
}
__syncthreads();
// Store sorted data in register back to shared memory. Then copy to global.
if(sortWarp)
DeviceThreadToShared<VT>(threadKeys, tid, keys_shared, false);
__syncthreads();
DeviceSharedToGlobal<NT, VT>(aCount + bCount, keys_shared, tid,
keysDest_global + gid);
////////////////////////////////////////////////////////////////////////////
// Use the merge indices to gather values from global memory. Store directly
// to valsDest_global.
if(HasValues) {
// Transpose the gather indices to help coalesce loads.
if(sortWarp)
DeviceThreadToShared<VT>(indices, tid, indices_shared, false);
else {
#pragma unroll
for(int i = 0; i < VT; ++i)
indices_shared[VT * tid + i] = VT * tid + i;
}
__syncthreads();
DeviceTransferMergeValuesShared<NT, VT>(aCount + bCount,
values_global + a0, values_global + b0, aCount, indices_shared,
tid, valsDest_global + NV * block);
}
}
////////////////////////////////////////////////////////////////////////////////
// DeviceSegSortCopy
template<int NT, int VT, bool HasValues, typename KeyType, typename ValueType>
MGPU_DEVICE void DeviceSegSortCopy(const KeyType* keys_global,
const ValueType* values_global, int tid, int block, int count,
KeyType* keysDest_global, ValueType* valsDest_global) {
int gid = NT * VT * block;
int count2 = min(NT * VT, count - gid);
DeviceGlobalToGlobal<NT, VT>(count2, keys_global + gid, tid,
keysDest_global + gid);
if(HasValues)
DeviceGlobalToGlobal<NT, VT>(count2, values_global + gid, tid,
valsDest_global + gid);
}
} // namespace mgpu
include/contrib/moderngpu/include/device/ctasortedsearch.cuh 0 → 100644
View file @ a60ed4d7
/******************************************************************************
* Copyright (c) 2013, NVIDIA CORPORATION. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
* * Neither the name of the NVIDIA CORPORATION nor the
* names of its contributors may be used to endorse or promote products
* derived from this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
******************************************************************************/
/******************************************************************************
*
* Code and text by Sean Baxter, NVIDIA Research
* See http://nvlabs.github.io/moderngpu for repository and documentation.
*
******************************************************************************/
#pragma once
#include "../mgpudevice.cuh"
#include "ctasearch.cuh"
namespace mgpu {
////////////////////////////////////////////////////////////////////////////////
// DeviceSerialSearch
template<int VT, MgpuBounds Bounds, bool RangeCheck, bool IndexA, bool MatchA,
bool IndexB, bool MatchB, typename T, typename Comp>
MGPU_DEVICE int3 DeviceSerialSearch(const T* keys_shared, int aBegin,
int aEnd, int bBegin, int bEnd, int aOffset, int bOffset, int* indices,
Comp comp) {
const int FlagA = IndexA ? 0x80000000 : 1;
const int FlagB = IndexB ? 0x80000000 : 1;
T aKey = keys_shared[aBegin];
T bKey = keys_shared[bBegin];
T aPrev, bPrev;
if(aBegin > 0) aPrev = keys_shared[aBegin - 1];
if(bBegin > 0) bPrev = keys_shared[bBegin - 1];
int decisions = 0;
int matchCountA = 0;
int matchCountB = 0;
#pragma unroll
for(int i = 0; i < VT; ++i) {
bool p;
if(RangeCheck && aBegin >= aEnd) p = false;
else if(RangeCheck && bBegin >= bEnd) p = true;
else p = (MgpuBoundsUpper == Bounds) ?
comp(aKey, bKey) :
!comp(bKey, aKey);
if(p) {
// aKey is smaller than bKey, so it is inserted before bKey.
// Save bKey's index (bBegin + first) as the result of the search
// and advance to the next needle in A.
bool match = false;
if(MatchA) {
// Test if there is an element in B that matches aKey.
if(MgpuBoundsUpper == Bounds) {
// Upper Bound: We're inserting aKey after bKey. If there
// is a match for aKey it must be bPrev. Check that bPrev
// is in range and equal to aKey.
// The predicate test result !comp(aKey, bPrev) was
// established on the previous A-advancing iteration (it
// failed the comp(aKey, bKey) test to get us to this
// point). Check the other half of the equality condition
// with a second comparison.
bool inRange = !RangeCheck || (bBegin > aEnd);
match = inRange && !comp(bPrev, aKey);
} else {
// Lower Bound: We're inserting aKey before bKey. If there
// is a match for aKey, it must be bKey. Check that bKey
// is in range and equal to aKey.
// The predicate test !comp(bKey, aKey) has established one
// half of the equality condition. We establish the other
// half with a second comparison.
bool inRange = !RangeCheck || (bBegin < bEnd);
match = inRange && !comp(aKey, bKey);
}
}
int index = 0;
if(IndexA) index = bOffset + bBegin;
if(match) index |= FlagA;
if(IndexA || MatchA) indices[i] = index;
matchCountA += match;
// Mark the decision bit to indicate that this iteration has
// progressed A (the needles).
decisions |= 1<< i;
aPrev = aKey;
aKey = keys_shared[++aBegin];
} else {
// aKey is larger than bKey, so it is inserted after bKey (but we
// don't know where yet). Advance the B index to the next element in
// the haystack to continue the search for the current needle.
bool match = false;
if(MatchB) {
if(MgpuBoundsUpper == Bounds) {
// Upper Bound: aKey is not smaller than bKey. We advance to
// the next haystack element in B. If there is a match in A
// for bKey it must be aKey. By entering this branch we've
// verified that !comp(aKey, bKey). Making the reciprocal
// comparison !comp(bKey, aKey) establishes aKey == bKey.
bool inRange = !RangeCheck ||
((bBegin < bEnd) && (aBegin < aEnd));
match = inRange && !comp(bKey, aKey);
} else {
// Lower Bound: bKey is smaller than aKey. We advance to the
// next element in B. If there is a match for bKey, it must
// be aPrev. The previous A-advancing iteration proved that
// !comp(bKey, aPrev). We test !comp(aPrev, bKey) for the
// other half of the equality condition.
bool inRange = !RangeCheck ||
((bBegin < bEnd) && (aBegin > 0));
match = inRange && !comp(aPrev, bKey);
}
}
int index = 0;
if(IndexB) index = aOffset + aBegin;
if(match) index |= FlagB;
if(IndexB || MatchB) indices[i] = index;
matchCountB += match;
// Keep the decision bit cleared to indicate that this iteration
// has progressed B (the haystack).
bPrev = bKey;
bKey = keys_shared[++bBegin];
}
}
return make_int3(decisions, matchCountA, matchCountB);
}
////////////////////////////////////////////////////////////////////////////////
// CTASortedSearch
// Take keys in shared memory and return indices and b-match flags in shared
// memory.
// NOTE: This function doesn't do any strided-to-thread order transposes so
// using an even number of values per thread will incur no additional bank
// conflicts.
template<int NT, int VT, MgpuBounds Bounds, bool IndexA, bool MatchA,
bool IndexB, bool MatchB, typename T, typename Comp>
MGPU_DEVICE int2 CTASortedSearch(T* keys_shared, int aStart, int aCount,
int aEnd, int a0, int bStart, int bCount, int bEnd, int b0, bool extended,
int tid, int* indices_shared, Comp comp) {
// Run a merge path to find the start of the serial search for each thread.
int diag = VT * tid;
int mp = MergePath<Bounds>(keys_shared + aStart, aCount,
keys_shared + bStart, bCount, diag, comp);
int a0tid = mp;
int b0tid = diag - mp;
// Serial search into register.
int3 results;
int indices[VT];
if(extended)
results = DeviceSerialSearch<VT, Bounds, false, IndexA, MatchA, IndexB,
MatchB>(keys_shared, a0tid + aStart, aEnd, b0tid + bStart, bEnd,
a0 - aStart, b0 - bStart, indices, comp);
else
results = DeviceSerialSearch<VT, Bounds, true, IndexA, MatchA, IndexB,
MatchB>(keys_shared, a0tid + aStart, aEnd, b0tid + bStart, bEnd,
a0 - aStart, b0 - bStart, indices, comp);
__syncthreads();
// Compact the indices into shared memory. Use the decision bits (set is A,
// cleared is B) to select the destination.
int decisions = results.x;
b0tid += aCount;
#pragma unroll
for(int i = 0; i < VT; ++i) {
if((1<< i) & decisions) {
if(IndexA || MatchA) indices_shared[a0tid++] = indices[i];
} else {
if(IndexB || MatchB) indices_shared[b0tid++] = indices[i];
}
}
__syncthreads();
// Return the match counts for A and B keys.
return make_int2(results.y, results.z);
}
} // namespace mgpu
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment