              gpu_particleCharge[i];
 
     rT11 += factor * (gpu_kxRef[imageID] * gpu_comDx[i]);
-    rT12 += factor * 0.5 * (gpu_kxRef[imageID] * gpu_comDy[i] + gpu_kyRef[imageID] * gpu_comDx[i]);
-    rT13 += factor * 0.5 * (gpu_kxRef[imageID] * gpu_comDz[i] + gpu_kzRef[imageID] * gpu_comDx[i]);
-    rT22 += factor * (gpu_kyRef[imageID] * gpu_comDy[i]);
-    rT13 += factor * 0.5 * (gpu_kyRef[imageID] * gpu_comDz[i] + gpu_kzRef[imageID] * gpu_comDy[i]);
-    rT33 += factor * (gpu_kzRef[imageID] * gpu_comDz[i]);
+    rT12 += factor * 0.5 * (gpu_kxRef[imageID] * shared_coords[particleID * 6 + 4] + gpu_kyRef[imageID] * shared_coords[particleID * 6 + 3]);
+    rT13 += factor * 0.5 * (gpu_kxRef[imageID] * shared_coords[particleID * 6 + 5] + gpu_kzRef[imageID] * shared_coords[particleID * 6 + 3]);
+    rT22 += factor * (gpu_kyRef[imageID] * shared_coords[particleID * 6 + 4]);
+    rT13 += factor * 0.5 * (gpu_kyRef[imageID] * shared_coords[particleID * 6 + 5] + gpu_kzRef[imageID] * shared_coords[particleID * 6 + 4]);
+    rT33 += factor * (gpu_kzRef[imageID] * shared_coords[particleID * 6 + 5]);
   }
 
   atomicAdd(&gpu_rT11[imageID], rT11);

  //Intramolecular part
  for(int particleID = 0; particleID < PARTICLE_PER_BLOCK; particleID++) {
    dot = DotProductGPU(gpu_kxRef[imageID], gpu_kyRef[imageID],
                        gpu_kzRef[imageID], shared_coords[particleID * 3],
                        shared_coords[particleID * 3 + 1],
                        shared_coords[particleID * 3 + 2]);

    double dotsin, dotcos;
    sincos(dot, &dotsin, &dotcos);

                                    uint imageSize,
                                    uint atomNumber)
{
  __shared__ double shared_coords[PARTICLE_PER_BLOCK*6];
  int imageID = blockIdx.x * blockDim.x + threadIdx.x;
  int offset_coordinates_index = blockIdx.y * PARTICLE_PER_BLOCK;
  if(threadIdx.x < PARTICLE_PER_BLOCK) {
    shared_coords[threadIdx.x * 6    ] = gpu_x[offset_coordinates_index + threadIdx.x];
    shared_coords[threadIdx.x * 6 + 1] = gpu_y[offset_coordinates_index + threadIdx.x];
    shared_coords[threadIdx.x * 6 + 2] = gpu_z[offset_coordinates_index + threadIdx.x];
    shared_coords[threadIdx.x * 6 + 3] = gpu_comDx[offset_coordinates_index + threadIdx.x];
    shared_coords[threadIdx.x * 6 + 4] = gpu_comDy[offset_coordinates_index + threadIdx.x];
    shared_coords[threadIdx.x * 6 + 5] = gpu_comDz[offset_coordinates_index + threadIdx.x];
  }

  if(imageID >= imageSize)
