Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-04-02 03:10:29

0001 #include "CUDADataFormats/Track/interface/TrajectoryStateSoAT.h"
0002 
0003 using Vector5d = Eigen::Matrix<double, 5, 1>;
0004 using Matrix5d = Eigen::Matrix<double, 5, 5>;
0005 
0006 __host__ __device__ Matrix5d loadCov(Vector5d const& e) {
0007   Matrix5d cov;
0008   for (int i = 0; i < 5; ++i)
0009     cov(i, i) = e(i) * e(i);
0010   for (int i = 0; i < 5; ++i) {
0011     for (int j = 0; j < i; ++j) {
0012       double v = 0.3 * std::sqrt(cov(i, i) * cov(j, j));  // this makes the matrix pos defined
0013       cov(i, j) = (i + j) % 2 ? -0.4 * v : 0.1 * v;
0014       cov(j, i) = cov(i, j);
0015     }
0016   }
0017   return cov;
0018 }
0019 
0020 using TS = TrajectoryStateSoAT<128>;
0021 
0022 __global__ void testTSSoA(TS* pts, int n) {
0023   assert(n <= 128);
0024 
0025   Vector5d par0;
0026   par0 << 0.2, 0.1, 3.5, 0.8, 0.1;
0027   Vector5d e0;
0028   e0 << 0.01, 0.01, 0.035, -0.03, -0.01;
0029   auto cov0 = loadCov(e0);
0030 
0031   TS& ts = *pts;
0032 
0033   int first = threadIdx.x + blockIdx.x * blockDim.x;
0034 
0035   for (int i = first; i < n; i += blockDim.x * gridDim.x) {
0036     ts.copyFromDense(par0, cov0, i);
0037     Vector5d par1;
0038     Matrix5d cov1;
0039     ts.copyToDense(par1, cov1, i);
0040     Vector5d delV = par1 - par0;
0041     Matrix5d delM = cov1 - cov0;
0042     for (int j = 0; j < 5; ++j) {
0043       assert(std::abs(delV(j)) < 1.e-5);
0044       for (auto k = j; k < 5; ++k) {
0045         assert(cov0(k, j) == cov0(j, k));
0046         assert(cov1(k, j) == cov1(j, k));
0047         assert(std::abs(delM(k, j)) < 1.e-5);
0048       }
0049     }
0050   }
0051 }
0052 
0053 #ifdef __CUDACC__
0054 #include "HeterogeneousCore/CUDAUtilities/interface/requireDevices.h"
0055 #include "HeterogeneousCore/CUDAUtilities/interface/cudaCheck.h"
0056 #endif
0057 
0058 int main() {
0059 #ifdef __CUDACC__
0060   cms::cudatest::requireDevices();
0061 #endif
0062 
0063   TS ts;
0064 
0065 #ifdef __CUDACC__
0066   TS* ts_d;
0067   cudaCheck(cudaMalloc(&ts_d, sizeof(TS)));
0068   testTSSoA<<<1, 64>>>(ts_d, 128);
0069   cudaCheck(cudaGetLastError());
0070   cudaCheck(cudaMemcpy(&ts, ts_d, sizeof(TS), cudaMemcpyDefault));
0071   cudaCheck(cudaDeviceSynchronize());
0072 #else
0073   testTSSoA(&ts, 128);
0074 #endif
0075 }