| 1 | #include "libtorch_stable/torch_utils.h" |
| 2 | |
| 3 | #include <cmath> |
| 4 | |
| 5 | #include "libtorch_stable/core/math.hpp" |
| 6 | #include "cuda_compat.h" |
| 7 | #include "libtorch_stable/dispatch_utils.h" |
| 8 | #include "quantization/w8a8/fp8/common.cuh" |
| 9 | |
| 10 | #ifndef USE_ROCM |
| 11 | #include <cuda_bf16.h> |
| 12 | #include <cuda_fp16.h> |
| 13 | #include <cuda_fp8.h> |
| 14 | #else |
| 15 | #include <hip/hip_bf16.h> |
| 16 | #include <hip/hip_fp16.h> |
| 17 | #include <hip/hip_fp8.h> |
| 18 | |
| 19 | typedef __hip_bfloat162 __nv_bfloat162; |
| 20 | typedef __hip_bfloat16 __nv_bfloat16; |
| 21 | typedef __hip_bfloat16_raw __nv_bfloat16_raw; |
| 22 | #if defined(HIP_FP8_TYPE_OCP) |
| 23 | typedef __hip_fp8_e4m3 __nv_fp8_e4m3; |
| 24 | typedef __hip_fp8x4_e4m3 __nv_fp8x4_e4m3; |
| 25 | #else |
| 26 | // ROCm 6.2 fallback: only *_fnuz types exist |
| 27 | typedef __hip_fp8_e4m3_fnuz __nv_fp8_e4m3; |
| 28 | typedef __hip_fp8x4_e4m3_fnuz __nv_fp8x4_e4m3; |
| 29 | #endif |
| 30 | #endif |
| 31 | |
| 32 | namespace vllm { |
| 33 | |
| 34 | template <typename T> |
| 35 | __device__ __forceinline__ T silu_kernel(const T& x) { |
| 36 | // x * sigmoid(x) |
| 37 | return (T)(((float)x) / (1.0f + expf((float)-x))); |
| 38 | } |
| 39 | |
| 40 | // Activation and gating kernel template. |
| 41 | template <typename scalar_t, scalar_t (*ACT_FN)(const scalar_t&), |
| 42 | typename fp8_type> |
| 43 | __global__ void act_and_mul_quant_kernel( |
| 44 | fp8_type* __restrict__ out, // [..., d] |
| 45 | const scalar_t* __restrict__ input, // [..., 2, d] |
| 46 | const float* scale, const int d) { |
| 47 | const int32_t blocks_per_token = gridDim.y; |
| 48 | |
| 49 | const int32_t elems_per_128bit_load = (128 / 8) / sizeof(scalar_t); |
| 50 | |
| 51 | // We don't expect the hidden dimension to exceed 32 bits so int32 should |
| 52 | // be safe here. |
| 53 | const int32_t tgt_elems_per_block = div_ceil(d, blocks_per_token); |
| 54 | const int32_t elems_per_block = |
| 55 | round_to_next_multiple_of(tgt_elems_per_block, elems_per_128bit_load); |
| 56 | const int32_t block_start = blockIdx.y * elems_per_block; |
| 57 | int32_t block_end = block_start + elems_per_block; |
| 58 | block_end = block_end > d ? d : block_end; |
| 59 | |
| 60 | // token_idx is 64 bit to prevent 32 bit overflow when the number of tokens |
| 61 | // is very large |
| 62 | const int64_t token_idx = blockIdx.x; |
| 63 | const scalar_t* __restrict__ x_ptr = input + token_idx * 2 * d; |
| 64 | const scalar_t* __restrict__ y_ptr = input + token_idx * 2 * d + d; |
| 65 | fp8_type* __restrict__ out_ptr = out + token_idx * d; |
| 66 | |
| 67 | // 128-bit vectorized code |
| 68 | const int32_t vec_loop_end = |
| 69 | round_to_previous_multiple_of(elems_per_128bit_load, block_end); |
| 70 | const int32_t vec_end_idx = vec_loop_end / elems_per_128bit_load; |
| 71 | const int32_t vec_start_idx = block_start / elems_per_128bit_load; |
| 72 | |
| 73 | const int4* __restrict__ x_128bit_ptr = reinterpret_cast<const int4*>(x_ptr); |
| 74 | const int4* __restrict__ y_128bit_ptr = reinterpret_cast<const int4*>(y_ptr); |
| 75 | int2* __restrict__ out_128bit_ptr = reinterpret_cast<int2*>(out_ptr); |
| 76 | |
| 77 | float inverted_scale = 1 / *scale; |
| 78 | #pragma unroll |
| 79 | for (int32_t vec_idx = vec_start_idx + threadIdx.x; vec_idx < vec_end_idx; |
| 80 | vec_idx += blockDim.x) { |
| 81 | const int4 x_128bit = VLLM_LDG(&x_128bit_ptr[vec_idx]); |
| 82 | const int4 y_128bit = VLLM_LDG(&y_128bit_ptr[vec_idx]); |
| 83 | using scalar_128bit_vec_t = std::array<scalar_t, elems_per_128bit_load>; |
| 84 | using scalar_64bit_vec_t = std::array<fp8_type, elems_per_128bit_load>; |
| 85 | |
| 86 | scalar_64bit_vec_t out_vec; |
| 87 | const auto x_vec = reinterpret_cast<scalar_128bit_vec_t const&>(x_128bit); |
| 88 | const auto y_vec = reinterpret_cast<scalar_128bit_vec_t const&>(y_128bit); |
| 89 | |
| 90 | #pragma unroll |
| 91 | for (int i = 0; i < elems_per_128bit_load; i++) { |
| 92 | out_vec[i] = scaled_fp8_conversion<true, fp8_type>( |
| 93 | ACT_FN(x_vec[i]) * y_vec[i], inverted_scale); |
| 94 | } |
| 95 | |
| 96 | out_128bit_ptr[vec_idx] = reinterpret_cast<const int2&>(out_vec); |
| 97 | } |
| 98 | |
| 99 | // Scalar cleanup code |
| 100 | if (block_end > vec_loop_end) { |
| 101 | for (int64_t idx = vec_loop_end + threadIdx.x; idx < block_end; |
| 102 | idx += blockDim.x) { |
| 103 | const scalar_t x = VLLM_LDG(&x_ptr[idx]); |
| 104 | const scalar_t y = VLLM_LDG(&y_ptr[idx]); |
| 105 | out_ptr[idx] = |
| 106 | scaled_fp8_conversion<true, fp8_type>(ACT_FN(x) * y, inverted_scale); |
| 107 | } |
| 108 | } |
| 109 | } |
| 110 | |
| 111 | __device__ __forceinline__ float silu(float x) { |
| 112 | return __fdividef(x, (1.f + expf(-x))); |
| 113 | } |
| 114 | |
| 115 | __device__ __forceinline__ float2 silu2(float2 x) { |
| 116 | return make_float2(silu(x.x), silu(x.y)); |
| 117 | } |
| 118 | |
| 119 | __device__ __forceinline__ __nv_bfloat162 silu2_v2(float2 x) { |
| 120 | #ifndef USE_ROCM |
| 121 | return make_bfloat162(__float2bfloat16_rn(silu(x.x)), |
| 122 | __float2bfloat16_rn(silu(x.y))); |
| 123 | #else |
| 124 | return __float22bfloat162_rn(make_float2(silu(x.x), silu(x.y))); |
| 125 | #endif |
| 126 | } |
| 127 | |
| 128 | #ifndef USE_ROCM |
| 129 | __device__ __forceinline__ float warp_max(float v) { |
| 130 | static constexpr unsigned FULL_MASK = 0xffffffffu; |
| 131 | for (int offset = 1; offset < WARP_SIZE; offset *= 2) { |
| 132 | v = fmaxf(v, __shfl_xor_sync(FULL_MASK, v, offset)); |
| 133 | } |
| 134 | return v; |
| 135 | } |
| 136 | |
| 137 | __device__ __forceinline__ __nv_bfloat16 warp_max(__nv_bfloat16 v) { |
| 138 | static constexpr unsigned FULL_MASK = 0xffffffffu; |
| 139 | for (int offset = 1; offset < WARP_SIZE; offset *= 2) { |
| 140 | v = __hmax(v, __shfl_xor_sync(FULL_MASK, v, offset)); |
| 141 | } |
| 142 | return v; |
| 143 | } |
| 144 | #endif |
| 145 | |
| 146 | template <typename T, typename U> |
| 147 | __device__ __forceinline__ void cp_async4(T* _smem_ptr, const U* _glob_ptr) { |
| 148 | #if __CUDACC_VER_MAJOR__ >= 11 && __CUDA_ARCH__ >= 800 |
| 149 | auto smem_ptr = reinterpret_cast<void*>(_smem_ptr); |
| 150 | auto glob_ptr = reinterpret_cast<const void*>(_glob_ptr); |
| 151 | const int BYTES = 16; |
| 152 | uint32_t smem = static_cast<uint32_t>(__cvta_generic_to_shared(smem_ptr)); |
| 153 | asm volatile( |
| 154 | "{\n" |
| 155 | " cp.async.cg.shared.global [%0], [%1], %2;\n" |
| 156 | "}\n" ::"r"(smem), |
| 157 | "l"(glob_ptr), "n"(BYTES)); |
| 158 | #else |
| 159 | _smem_ptr[0] = _glob_ptr[0]; |
| 160 | #endif |
| 161 | } |
| 162 | |
| 163 | __device__ __forceinline__ void cp_async_fence() { |
| 164 | #if __CUDACC_VER_MAJOR__ >= 11 && __CUDA_ARCH__ >= 800 |
| 165 | asm volatile("cp.async.commit_group;\n" ::); |
| 166 | #else |
| 167 | #endif |
| 168 | } |
| 169 | |
| 170 | template <int N> |
| 171 | __device__ __forceinline__ void cp_async_wait() { |
| 172 | #if __CUDACC_VER_MAJOR__ >= 11 && __CUDA_ARCH__ >= 800 |
| 173 | asm volatile("cp.async.wait_group %0;\n" ::"n"(N)); |
| 174 | #else |
| 175 | #endif |
| 176 | } |
| 177 | |
| 178 | template <> |
| 179 | __device__ __forceinline__ void cp_async_wait<0>() { |
| 180 | #if __CUDACC_VER_MAJOR__ >= 11 && __CUDA_ARCH__ >= 800 |
| 181 | asm volatile("cp.async.wait_all;\n" ::); |
| 182 | #else |
| 183 | #endif |
| 184 | } |
| 185 | |
| 186 | __device__ __forceinline__ float clip(float v, float mmin, float mmax) { |
| 187 | return fminf(mmax, fmaxf(v, mmin)); |
| 188 | } |
| 189 | |
| 190 | __device__ __forceinline__ __nv_bfloat16 clip(__nv_bfloat16 v, |
| 191 | __nv_bfloat16 mmin, |
| 192 | __nv_bfloat16 mmax) { |
| 193 | return __hmin(mmax, __hmax(v, mmin)); |
| 194 | } |
| 195 | |
| 196 | __device__ __forceinline__ __nv_bfloat162 clip(__nv_bfloat162 v, |
| 197 | __nv_bfloat162 mmin, |
| 198 | __nv_bfloat162 mmax) { |
| 199 | return __hmin2(mmax, __hmax2(v, mmin)); |
| 200 | } |
| 201 | |
| 202 | // We use the following values for fp8 min/max: |
| 203 | // __nv_fp8_e4m3 = (-448, +448) |
| 204 | // __nv_fp8_e4m3uz = (-240.0, +240.0) |
| 205 | // It is currently assumed that only |
| 206 | template <class T> |
| 207 | constexpr __nv_bfloat16 get_fp8_max() { |
| 208 | static_assert(std::is_same_v<T, c10::Float8_e4m3fn> || |
| 209 | std::is_same_v<T, c10::Float8_e4m3fnuz>); |
| 210 | if constexpr (std::is_same_v<T, c10::Float8_e4m3fn>) { |
| 211 | return __nv_bfloat16(__nv_bfloat16_raw{.x = 17376}); |
| 212 | } else { |
| 213 | return __nv_bfloat16(__nv_bfloat16_raw{.x = 17264}); |
| 214 | } |
| 215 | } |
| 216 | |
| 217 | template <class T> |
| 218 | constexpr __nv_bfloat16 get_fp8_min() { |
| 219 | static_assert(std::is_same_v<T, c10::Float8_e4m3fn> || |
| 220 | std::is_same_v<T, c10::Float8_e4m3fnuz>); |
| 221 | if constexpr (std::is_same_v<T, c10::Float8_e4m3fn>) { |
| 222 | return __nv_bfloat16(__nv_bfloat16_raw{.x = 50144}); |
| 223 | } else { |
| 224 | return __nv_bfloat16(__nv_bfloat16_raw{.x = 50032}); |
| 225 | } |
| 226 | } |
| 227 | |
| 228 | template <typename Idx_t> |
| 229 | __device__ __forceinline__ int warp_expert_search( |
| 230 | int idx, int n, const Idx_t* __restrict__ input, Idx_t val) { |
| 231 | const Idx_t* input_ptr = input + idx; |
| 232 | int base_offset = 0; |
| 233 | |
| 234 | for (;;) { |
| 235 | bool move_on = (idx < n && *input_ptr <= val); |
| 236 | |
| 237 | unsigned mask = __ballot_sync(0xffffffff, move_on); |
| 238 | |
| 239 | if (mask != 0xffffffffu) { |
| 240 | int last_lane = 31 - __clz(mask); |
| 241 | return base_offset + last_lane; |
| 242 | } |
| 243 | |
| 244 | input_ptr += 32; |
| 245 | base_offset += 32; |
| 246 | idx += 32; |
| 247 | } |
| 248 | } |
| 249 | |
| 250 | template <int num_parallel_tokens> |
| 251 | __device__ __forceinline__ void token_bounds(int32_t n_tokens, |
| 252 | int32_t worker_id, |
| 253 | int32_t& n_tokens_lower, |
| 254 | int32_t& n_tokens_upper) { |
| 255 | if (n_tokens < num_parallel_tokens && worker_id < n_tokens) { |
| 256 | if (worker_id >= num_parallel_tokens) return; |
| 257 | n_tokens_lower = worker_id; |
| 258 | n_tokens_upper = worker_id + 1; |
| 259 | } else { |
| 260 | int32_t chunk_size = n_tokens / num_parallel_tokens; |
| 261 | int32_t residual = n_tokens - chunk_size * num_parallel_tokens; |
| 262 | auto calc_id = [&](int32_t id) { |
| 263 | if (id < residual) |
| 264 | return min(n_tokens, id * (chunk_size + 1)); |
| 265 | else |
| 266 | return min(n_tokens, id * chunk_size + residual); |
| 267 | }; |
| 268 | n_tokens_lower = calc_id(worker_id); |
| 269 | n_tokens_upper = calc_id(worker_id + 1); |
| 270 | } |
| 271 | } |
| 272 | |
| 273 | template <int BLOCK_COUNT, int SMEM_SIZE_BYTES_Y, typename fp8_type, |
| 274 | typename scale_t, int THREADS, typename Idx_t, bool CEIL_UE8M0, |
| 275 | int GROUP_SIZE = 128, int NUM_STAGES = 3> |
| 276 | __global__ void silu_mul_fp8_quant_deep_gemm_kernel( |
| 277 | const __nv_bfloat16* __restrict__ _input, fp8_type* __restrict__ _y_q, |
| 278 | scale_t* __restrict__ _y_s, const int32_t* __restrict__ tokens_per_expert, |
| 279 | // sizes |
| 280 | Idx_t E, Idx_t T, Idx_t H, |
| 281 | // strides (in elements) |
| 282 | Idx_t stride_i_e, Idx_t stride_i_t, Idx_t stride_i_h, Idx_t stride_yq_e, |
| 283 | Idx_t stride_yq_t, Idx_t stride_yq_h, Idx_t stride_ys_e, Idx_t stride_ys_t, |
| 284 | Idx_t stride_ys_g, Idx_t stride_ys_p, Idx_t stride_counts_e) { |
| 285 | #ifndef USE_ROCM |
| 286 | static constexpr int NUM_WARPS = THREADS / WARP_SIZE; |
| 287 | |
| 288 | static constexpr int LOAD_STAGE_SIZE = 2 * GROUP_SIZE / 8; |
| 289 | static constexpr int LOAD_STAGE_MOD = NUM_STAGES * LOAD_STAGE_SIZE; |
| 290 | |
| 291 | static constexpr int COMPUTE_STAGE_SIZE = 2 * GROUP_SIZE / 4; |
| 292 | static constexpr int COMPUTE_STAGE_MOD = COMPUTE_STAGE_SIZE * NUM_STAGES; |
| 293 | |
| 294 | extern __shared__ __align__(16) __int128_t smem_128[]; |
| 295 | |
| 296 | int* s_expert_offsets = |
| 297 | reinterpret_cast<int*>(smem_128 + (SMEM_SIZE_BYTES_Y / 16)); |
| 298 | |
| 299 | static constexpr __nv_bfloat16 fp8_min = get_fp8_min<fp8_type>(); |
| 300 | static constexpr __nv_bfloat16 fp8_max = get_fp8_max<fp8_type>(); |
| 301 | // We assign EPS with it's 16-bit unsigned counterpart to allow constexpr. |
| 302 | static constexpr __nv_bfloat16 EPS = (__nv_bfloat16_raw{.x = 11996}); |
| 303 | int tid = threadIdx.x; |
| 304 | int warp_id = tid >> 5; |
| 305 | int lane_id = tid & 0x1f; |
| 306 | |
| 307 | int running_sum{}; |
| 308 | if (!warp_id) { |
| 309 | for (int i = 0; i < E; i += WARP_SIZE) { |
| 310 | bool valid = (i + threadIdx.x) < E; |
| 311 | int value = |
| 312 | (valid ? tokens_per_expert[i + threadIdx.x * stride_counts_e] : 0) + |
| 313 | (!lane_id ? running_sum : 0); |
| 314 | |
| 315 | for (int offset = 1; offset < 32; offset *= 2) { |
| 316 | int n = __shfl_up_sync(0xFFFFFFFFu, value, offset); |
| 317 | if (lane_id >= offset) value += n; |
| 318 | } |
| 319 | |
| 320 | if (valid) { |
| 321 | s_expert_offsets[i + threadIdx.x + 1] = value; |
| 322 | } |
| 323 | |
| 324 | running_sum = __shfl_sync(0xFFFFFFFFu, value, WARP_SIZE - 1); |
| 325 | } |
| 326 | |
| 327 | if (!lane_id) { |
| 328 | s_expert_offsets[0] = 0; |
| 329 | } |
| 330 | } |
| 331 | |
| 332 | __syncthreads(); |
| 333 | |
| 334 | int32_t total_tokens = s_expert_offsets[E]; |
| 335 | |
| 336 | const int warp_position_yq = warp_id * (H / NUM_WARPS); |
| 337 | const int warp_position_scales = warp_id * (H / (GROUP_SIZE * NUM_WARPS)); |
| 338 | |
| 339 | // A single block will handle tokens_per_block tokens. |
| 340 | // Each block i iterates over tokens of a slice of n_tokens = |
| 341 | // expert_counts[i], with the size of chunk being |
| 342 | // (n_tokens / NUM_PARALLEL_TOKENS) + residual, instead of |
| 343 | // updiv(n_tokens, NUM_PARALLEL_TOKENS) for better scheduling. |
| 344 | |
| 345 | // Each warp will get space to store its hidden dim for gate and up. |
| 346 | __int128_t* s_hidden_load = smem_128 + warp_id * ((2 * 128 / 8) * NUM_STAGES); |
| 347 | __int128_t* smem_load_ptr = s_hidden_load + lane_id; |
| 348 | |
| 349 | const __nv_bfloat16 fp8_inv = __hdiv(__float2bfloat16(1.f), fp8_max); |
| 350 | |
| 351 | int32_t compute_pipeline_offset_64 = 0; |
| 352 | int32_t load_stage_offset{}; |
| 353 | const __nv_bfloat16 one_bf16 = __float2bfloat16_rn(1.f); |
| 354 | |
| 355 | __int64_t* smem_compute_ptr = reinterpret_cast<__int64_t*>(smem_128) + |
| 356 | warp_id * (2 * (GROUP_SIZE / 4) * NUM_STAGES) + |
| 357 | lane_id; |
| 358 | __int64_t* s_gate64_ptr = smem_compute_ptr; |
| 359 | __int64_t* s_up64_ptr = smem_compute_ptr + GROUP_SIZE / 4; |
| 360 | |
| 361 | int tokens_lower, tokens_upper; |
| 362 | |
| 363 | token_bounds<BLOCK_COUNT>(total_tokens, blockIdx.x, tokens_lower, |
| 364 | tokens_upper); |
| 365 | |
| 366 | Idx_t expert_id{}, expert_offset{}, next_expert_offset{}; |
| 367 | int token_id = tokens_lower; |
| 368 | int32_t t_load{}; |
| 369 | |
| 370 | if (token_id < tokens_upper) { |
| 371 | expert_id = warp_expert_search<int>(lane_id, E, s_expert_offsets, token_id); |
| 372 | expert_offset = s_expert_offsets[expert_id]; |
| 373 | next_expert_offset = s_expert_offsets[expert_id + 1]; |
| 374 | } else { |
| 375 | // This thread block has no work to do. |
| 376 | return; |
| 377 | } |
| 378 | |
| 379 | int t_load_bound = H / (GROUP_SIZE * NUM_WARPS); |
| 380 | |
| 381 | Idx_t base_i = ((expert_id * stride_i_e) / 8) + |
| 382 | (token_id - expert_offset) * stride_i_t / 8; |
| 383 | const Idx_t gate_warp_offset = |
| 384 | warp_id * ((stride_i_h * H) / (8 * NUM_WARPS)) + (lane_id & 0b1111); |
| 385 | |
| 386 | const __int128_t* input_128_ptr = |
| 387 | reinterpret_cast<const __int128_t*>(_input) + gate_warp_offset + |
| 388 | ((lane_id < 16) ? 0 : ((H * stride_i_h) / 8)); |
| 389 | __int128_t* load_ptr = const_cast<__int128_t*>(input_128_ptr + base_i); |
| 390 | |
| 391 | auto token_offset = token_id - expert_offset; |
| 392 | |
| 393 | auto load_and_advance_y_pred = [&] { |
| 394 | if (t_load < t_load_bound) { |
| 395 | // Here we are simply continuing to load data |
| 396 | // from the current token. |
| 397 | auto smem_load_ptr_staged = smem_load_ptr + load_stage_offset; |
| 398 | |
| 399 | // It is very important that LOAD_STAGE_SIZE is constexpr to avoid |
| 400 | // unnecessary ALU ops. |
| 401 | load_stage_offset += LOAD_STAGE_SIZE; |
| 402 | load_stage_offset %= LOAD_STAGE_MOD; |
| 403 | |
| 404 | cp_async4(smem_load_ptr_staged, load_ptr); |
| 405 | load_ptr += GROUP_SIZE / 8; |
| 406 | ++t_load; |
| 407 | } else if (token_id + 1 < tokens_upper) { |
| 408 | // We loaded everything from the current token, let's move on |
| 409 | // to the next one, and we checked that we have more tokens to load. |
| 410 | ++token_id; |
| 411 | t_load = 0; |
| 412 | if (token_id >= next_expert_offset) { |
| 413 | // We need to find the next expert. |
| 414 | do { |
| 415 | // This is a loop because it's possible |
| 416 | // that some experts are assigned 0 tokens. |
| 417 | // NOTE: We are guaranteed that there's at least |
| 418 | // one more token left so we don't have to check for |
| 419 | // expert_id bounds. |
| 420 | ++expert_id; |
| 421 | // This skips 1 memory read. |
| 422 | expert_offset = next_expert_offset; |
| 423 | next_expert_offset = s_expert_offsets[expert_id + 1]; |
| 424 | } while (next_expert_offset == expert_offset); |
| 425 | |
| 426 | base_i = expert_id * (stride_i_e / 8); |
| 427 | token_offset = 0; |
| 428 | load_ptr = const_cast<__int128_t*>(input_128_ptr + base_i); |
| 429 | } else { |
| 430 | // We remain within the same expert, so just |
| 431 | // move by H/4 __int128_t (2 * H/8). |
| 432 | base_i += stride_yq_t / 4; |
| 433 | token_offset++; |
| 434 | } |
| 435 | |
| 436 | load_ptr = const_cast<__int128_t*>(input_128_ptr + base_i); |
| 437 | |
| 438 | auto smem_load_ptr_staged = smem_load_ptr + load_stage_offset; |
| 439 | |
| 440 | // It is very important that LOAD_STAGE_SIZE is constexpr to avoid |
| 441 | // unnecessary ALU ops. |
| 442 | load_stage_offset += LOAD_STAGE_SIZE; |
| 443 | load_stage_offset %= LOAD_STAGE_MOD; |
| 444 | |
| 445 | cp_async4(smem_load_ptr_staged, load_ptr); |
| 446 | load_ptr += GROUP_SIZE / 8; |
| 447 | ++t_load; |
| 448 | } |
| 449 | // We fence even if there is nothing to load to simplify pipelining. |
| 450 | cp_async_fence(); |
| 451 | }; |
| 452 | |
| 453 | // We need to warm-up the pipeline. |
| 454 | #pragma unroll |
| 455 | for (int i = 0; i < NUM_STAGES - 1; i++) { |
| 456 | load_and_advance_y_pred(); |
| 457 | } |
| 458 | |
| 459 | __nv_fp8x4_e4m3* y_q_base_ptr = |
| 460 | reinterpret_cast<__nv_fp8x4_e4m3*>(_y_q) + lane_id; |
| 461 | |
| 462 | Idx_t scale_group_offset = 0; |
| 463 | if constexpr (std::is_same<scale_t, uint8_t>::value) { |
| 464 | // packed int32_t format |
| 465 | int pack_id = warp_position_scales / 4; |
| 466 | int scale_in_pack = warp_position_scales % 4; |
| 467 | scale_group_offset = pack_id * stride_ys_p + scale_in_pack * stride_ys_g; |
| 468 | } else { |
| 469 | scale_group_offset = warp_position_scales * stride_ys_g; |
| 470 | } |
| 471 | |
| 472 | scale_t* const y_scale_base_ptr = _y_s + scale_group_offset; |
| 473 | |
| 474 | for (auto j = tokens_lower; j < tokens_upper; j++) { |
| 475 | int current_group_id = warp_position_scales; // Running count of which |
| 476 | // group is being processed |
| 477 | const Idx_t base_ys = expert_id * stride_ys_e; |
| 478 | auto y_s_ptr = y_scale_base_ptr + base_ys + token_offset * stride_ys_t; |
| 479 | __nv_fp8x4_e4m3* y_q_ptr = |
| 480 | y_q_base_ptr + (expert_id * stride_yq_e + token_offset * stride_yq_t + |
| 481 | warp_position_yq * stride_yq_h) / |
| 482 | 4; |
| 483 | const int COMPUTE_LIMIT = H / (GROUP_SIZE * NUM_WARPS); |
| 484 | |
| 485 | for (int i = 0; i < COMPUTE_LIMIT; i++) { |
| 486 | cp_async_wait<NUM_STAGES - 2>(); |
| 487 | __syncthreads(); |
| 488 | load_and_advance_y_pred(); |
| 489 | |
| 490 | __int64_t* gate64_ptr = s_gate64_ptr + compute_pipeline_offset_64; |
| 491 | __int64_t* up64_ptr = s_up64_ptr + compute_pipeline_offset_64; |
| 492 | |
| 493 | // COMPUTE_STAGE_SIZE/MOD must also be constexpr! |
| 494 | compute_pipeline_offset_64 += COMPUTE_STAGE_SIZE; |
| 495 | compute_pipeline_offset_64 %= COMPUTE_STAGE_MOD; |
| 496 | |
| 497 | __int64_t gate64 = *gate64_ptr; |
| 498 | __int64_t up64 = *up64_ptr; |
| 499 | |
| 500 | // Compute |
| 501 | __nv_bfloat162 res[2]; |
| 502 | __nv_bfloat162* s_up_comp = reinterpret_cast<__nv_bfloat162*>(&up64); |
| 503 | __nv_bfloat162* s_gate_comp = reinterpret_cast<__nv_bfloat162*>(&gate64); |
| 504 | |
| 505 | #pragma unroll |
| 506 | for (int32_t k = 0; k < 2; ++k) { |
| 507 | __nv_bfloat162 gate = silu2_v2(__bfloat1622float2(s_gate_comp[k])); |
| 508 | res[k] = __hmul2(gate, s_up_comp[k]); |
| 509 | } |
| 510 | |
| 511 | auto _y_max2 = __hmax2(__habs2(res[0]), __habs2(res[1])); |
| 512 | |
| 513 | _y_max2.x = __hmax(__hmax(_y_max2.x, _y_max2.y), EPS); |
| 514 | |
| 515 | __nv_bfloat16 y_s = __hmul(warp_max(_y_max2.x), fp8_inv); |
| 516 | |
| 517 | if constexpr (CEIL_UE8M0) { |
| 518 | y_s = hexp2(hceil(hlog2(y_s))); |
| 519 | } |
| 520 | |
| 521 | __nv_bfloat16 inv_y = __hdiv(one_bf16, y_s); |
| 522 | |
| 523 | auto y_s2 = make_bfloat162(inv_y, inv_y); |
| 524 | |
| 525 | #pragma unroll |
| 526 | for (int32_t k = 0; k < 2; ++k) { |
| 527 | res[k] = clip(__hmul2(res[k], y_s2), __bfloat162bfloat162(fp8_min), |
| 528 | __bfloat162bfloat162(fp8_max)); |
| 529 | } |
| 530 | |
| 531 | *y_q_ptr = __nv_fp8x4_e4m3(res[0], res[1]); |
| 532 | y_q_ptr += WARP_SIZE * stride_yq_h; |
| 533 | |
| 534 | if (!lane_id) { |
| 535 | // Store scales. |
| 536 | if constexpr (std::is_same<scale_t, uint8_t>::value) { |
| 537 | // Packed UE8M0 format. Remove Mantissa. |
| 538 | *y_s_ptr = reinterpret_cast<int16_t&>(y_s) >> 7; |
| 539 | |
| 540 | bool const jump_pack = (current_group_id + 1) % 4 == 0; |
| 541 | // Minus 3 because we need to get to the first group in the |
| 542 | // next pack. |
| 543 | y_s_ptr += jump_pack ? (stride_ys_p - 3) : stride_ys_g; |
| 544 | |
| 545 | } else { |
| 546 | // float32 format |
| 547 | static_assert(std::is_same<scale_t, float>::value); |
| 548 | *y_s_ptr = y_s; |
| 549 | y_s_ptr += stride_ys_g; |
| 550 | } |
| 551 | |
| 552 | current_group_id += 1; |
| 553 | } |
| 554 | } |
| 555 | } |
| 556 | #endif |
| 557 | } |
| 558 | |
| 559 | } // namespace vllm |
| 560 | |
| 561 | // Launch activation, gating, and quantize kernel. |
| 562 | #define LAUNCH_ACTIVATION_GATE_KERNEL(KERNEL) \ |
| 563 | int d = input.size(-1) / 2; \ |
| 564 | int64_t num_tokens = input.numel() / input.size(-1); \ |
| 565 | dim3 grid(num_tokens, num_tokens > 16 ? num_tokens > 32 ? 1 : 2 : 4); \ |
| 566 | dim3 block(std::min(d, 512)); \ |
| 567 | const torch::stable::accelerator::DeviceGuard device_guard( \ |
| 568 | input.get_device_index()); \ |
| 569 | const cudaStream_t stream = \ |
| 570 | get_current_cuda_stream(input.get_device_index()); \ |
| 571 | VLLM_STABLE_DISPATCH_FLOATING_TYPES( \ |
| 572 | input.scalar_type(), "act_and_mul_kernel", [&] { \ |
| 573 | VLLM_STABLE_DISPATCH_FP8_TYPES( \ |
| 574 | out.scalar_type(), "act_and_mul_quant_kernel_fp8_type", [&] { \ |
| 575 | vllm::act_and_mul_quant_kernel<scalar_t, KERNEL<scalar_t>, \ |
| 576 | fp8_t> \ |
| 577 | <<<grid, block, 0, stream>>>( \ |
| 578 | out.mutable_data_ptr<fp8_t>(), \ |
| 579 | input.const_data_ptr<scalar_t>(), \ |
| 580 | scale.const_data_ptr<float>(), d); \ |
| 581 | }); \ |
| 582 | }); |
| 583 | |
| 584 | void silu_and_mul_quant(torch::stable::Tensor& out, // [..., d] |
| 585 | torch::stable::Tensor& input, // [..., 2 * d] |
| 586 | torch::stable::Tensor& scale) { |
| 587 | STD_TORCH_CHECK( |
| 588 | out.scalar_type() == torch::headeronly::ScalarType::Float8_e4m3fn || |
| 589 | out.scalar_type() == torch::headeronly::ScalarType::Float8_e4m3fnuz); |
| 590 | STD_TORCH_CHECK( |
| 591 | input.scalar_type() == torch::headeronly::ScalarType::Half || |
| 592 | input.scalar_type() == torch::headeronly::ScalarType::BFloat16, |
| 593 | "Input must be FP16 or BF16"); |
| 594 | STD_TORCH_CHECK(input.size(-1) % 2 == 0); |
| 595 | LAUNCH_ACTIVATION_GATE_KERNEL(vllm::silu_kernel); |
| 596 | } |
| 597 | |
| 598 | void persistent_masked_m_silu_mul_quant( |
| 599 | const torch::stable::Tensor& input, // (E, T, 2*H) |
| 600 | const torch::stable::Tensor& tokens_per_expert, // (E) |
| 601 | torch::stable::Tensor& y_q, // (E, T, H) [OUT] |
| 602 | torch::stable::Tensor& y_s, // (E, T, H//group_size) [OUT] |
| 603 | bool cast_scale_ue8m0) { |
| 604 | #ifndef USE_ROCM |
| 605 | |
| 606 | // This kernel currently only supports H % 128 == 0 and assumes a |
| 607 | // fixed GROUP_SIZE of 128. |
| 608 | static constexpr int GROUP_SIZE = 128; |
| 609 | |
| 610 | STD_TORCH_CHECK(input.scalar_type() == |
| 611 | torch::headeronly::ScalarType::BFloat16); |
| 612 | STD_TORCH_CHECK( |
| 613 | y_q.scalar_type() == torch::headeronly::ScalarType::Float8_e4m3fn || |
| 614 | y_q.scalar_type() == torch::headeronly::ScalarType::Float8_e4m3fnuz); |
| 615 | STD_TORCH_CHECK(input.size(-1) % (GROUP_SIZE * 2) == 0); |
| 616 | |
| 617 | bool const is_packed_ue8m0 = |
| 618 | (y_s.scalar_type() == torch::headeronly::ScalarType::Int && |
| 619 | cast_scale_ue8m0); |
| 620 | STD_TORCH_CHECK(y_s.scalar_type() == torch::headeronly::ScalarType::Float || |
| 621 | is_packed_ue8m0); |
| 622 | |
| 623 | using Idx_t = int64_t; |
| 624 | |
| 625 | Idx_t E = input.size(0); |
| 626 | Idx_t T = input.size(1); |
| 627 | Idx_t H = input.size(2) / 2; |
| 628 | Idx_t stride_i_e = input.stride(0); |
| 629 | Idx_t stride_i_t = input.stride(1); |
| 630 | Idx_t stride_i_h = input.stride(2); |
| 631 | Idx_t stride_yq_e = y_q.stride(0); |
| 632 | Idx_t stride_yq_t = y_q.stride(1); |
| 633 | Idx_t stride_yq_h = y_q.stride(2); |
| 634 | |
| 635 | Idx_t stride_counts_e = tokens_per_expert.stride(0); |
| 636 | |
| 637 | int const NUM_GROUPS = H / GROUP_SIZE; |
| 638 | |
| 639 | const cudaStream_t stream = get_current_cuda_stream(input.get_device_index()); |
| 640 | |
| 641 | // TODO: Get this from cuda_arch ? |
| 642 | static constexpr int SILU_V2_BLOCK_COUNT = 132 * 32; |
| 643 | |
| 644 | #define KERNEL(BLOCK_COUNT, scale_t, STRIDE_YS_E, STRIDE_YS_T, STRIDE_YS_G, \ |
| 645 | STRIDE_YS_P, CEIL_UE8M0, THREAD_COUNT, STAGES) \ |
| 646 | static constexpr int NUM_WARPS = THREAD_COUNT / WARP_SIZE; \ |
| 647 | int sms = SILU_V2_BLOCK_COUNT; \ |
| 648 | static constexpr int max_shared_mem_bytes = \ |
| 649 | GROUP_SIZE * 2 * STAGES * NUM_WARPS * 2; \ |
| 650 | dim3 grid(sms), block(THREAD_COUNT); \ |
| 651 | const torch::stable::accelerator::DeviceGuard device_guard( \ |
| 652 | input.get_device_index()); \ |
| 653 | VLLM_STABLE_DISPATCH_FP8_TYPES( \ |
| 654 | y_q.scalar_type(), "silu_mul_fp8_quant_deep_gemm_kernel", [&] { \ |
| 655 | vllm::silu_mul_fp8_quant_deep_gemm_kernel< \ |
| 656 | BLOCK_COUNT, max_shared_mem_bytes, fp8_t, scale_t, THREAD_COUNT, \ |
| 657 | Idx_t, CEIL_UE8M0, GROUP_SIZE, STAGES> \ |
| 658 | <<<grid, block, max_shared_mem_bytes + (E + 1) * 16, stream>>>( \ |
| 659 | reinterpret_cast<const __nv_bfloat16*>( \ |
| 660 | input.const_data_ptr()), \ |
| 661 | y_q.mutable_data_ptr<fp8_t>(), \ |
| 662 | reinterpret_cast<scale_t*>(y_s.mutable_data_ptr()), \ |
| 663 | reinterpret_cast<const int32_t*>( \ |
| 664 | tokens_per_expert.const_data_ptr()), \ |
| 665 | E, T, H, stride_i_e, stride_i_t, stride_i_h, stride_yq_e, \ |
| 666 | stride_yq_t, stride_yq_h, STRIDE_YS_E, STRIDE_YS_T, \ |
| 667 | STRIDE_YS_G, STRIDE_YS_P, stride_counts_e); \ |
| 668 | }); |
| 669 | |
| 670 | #define LAUNCH_ON_H(scale_t, STRIDE_YS_E, STRIDE_YS_T, STRIDE_YS_G, \ |
| 671 | STRIDE_YS_P, CEIL_UE8M0) \ |
| 672 | if (H >= 4096 && (NUM_GROUPS % 8) == 0) { \ |
| 673 | /* 8 warp config */ \ |
| 674 | static constexpr int NUM_STAGES = 4; \ |
| 675 | static constexpr int THREAD_COUNT = 256; \ |
| 676 | KERNEL(SILU_V2_BLOCK_COUNT, scale_t, STRIDE_YS_E, STRIDE_YS_T, \ |
| 677 | STRIDE_YS_G, STRIDE_YS_P, CEIL_UE8M0, THREAD_COUNT, NUM_STAGES); \ |
| 678 | } else { \ |
| 679 | /* 1 warp config */ \ |
| 680 | static constexpr int THREAD_COUNT = 32; \ |
| 681 | KERNEL(SILU_V2_BLOCK_COUNT, scale_t, STRIDE_YS_E, STRIDE_YS_T, \ |
| 682 | STRIDE_YS_G, STRIDE_YS_P, CEIL_UE8M0, THREAD_COUNT, 2); \ |
| 683 | } |
| 684 | |
| 685 | Idx_t stride_ys_e = y_s.stride(0); |
| 686 | Idx_t stride_ys_t = y_s.stride(1); |
| 687 | Idx_t stride_ys_g = y_s.stride(2); |
| 688 | Idx_t stride_ys_p = 0; |
| 689 | if (!cast_scale_ue8m0) { |
| 690 | STD_TORCH_CHECK(!is_packed_ue8m0); |
| 691 | LAUNCH_ON_H(float, stride_ys_e, stride_ys_t, stride_ys_g, stride_ys_p, |
| 692 | false); |
| 693 | return; |
| 694 | } |
| 695 | |
| 696 | if (!is_packed_ue8m0) { |
| 697 | // UE8M0 but not packed |
| 698 | LAUNCH_ON_H(float, stride_ys_e, stride_ys_t, stride_ys_g, stride_ys_p, |
| 699 | true); |
| 700 | return; |
| 701 | } |
| 702 | |
| 703 | STD_TORCH_CHECK(cast_scale_ue8m0 && is_packed_ue8m0); |
| 704 | STD_TORCH_CHECK(y_s.scalar_type() == torch::headeronly::ScalarType::Int); |
| 705 | |
| 706 | // Int32 packed ue8m0 scales tensor. |
| 707 | // Let E, T, G be the number to experts, number of tokens and number of groups |
| 708 | // respectively. Let, E = 2, T = 4, G = 6, in this case the int32 scales |
| 709 | // tensor are of shape [1, 4, 2] and stride [8, 1, 4]. The scales are expected |
| 710 | // to be arranged as follows, |
| 711 | // [[T0G0-T0G1-T0G2-T0G3, T0G4-T0G5-X-X,], |
| 712 | // [T1G0-T1G1-T1G2-T1G3, T1G4-T1G5-X-X,] |
| 713 | // [T2G0-T2G1-T2G2-T2G3, T2G4-T2G5-X-X,] |
| 714 | // [T3G0-T3G1-T3G2-T3G3, T3G4-T3G5-X-X,]] |
| 715 | // where, TxGy is the scale ue8m0 scale value of Token x, Group y. |
| 716 | // |
| 717 | // In memory (in bytes) the scale values are arranged as, |
| 718 | // [T0G0, T0G1, T0G2, T0G3, T1G0, T1G2, T1G3, T1G4, T2G0, T2G1, T2G3, T2G4, |
| 719 | // T3G0, T3G1, T3G2, T3G3, T0G4, T0G5, X, X, T1G4, T1G5, X, X, T2G4, T2G5, |
| 720 | // X, X, T3G4, T3G5, X, X] |
| 721 | // |
| 722 | // An Int32 tensor of size [1, 4, 2] and stride [8, 1, 4] can be represented |
| 723 | // as an uint8 tensor of shape [1, 2, 4, 4] and stride [32, 16, 4, 1]. In |
| 724 | // english, ignoring the Experts dimension, the original int32 tensor is |
| 725 | // simply treated as two packed [4, 4] uint8 tensor (or two [4, 1] int32 |
| 726 | // tensor). The following strides setting reflects this change. Caveat: This |
| 727 | // means that the G dimension is no longer contiguous. i.e. Note that to move |
| 728 | // from G3 to G4, we need to jump along the packing dimension. The kernel |
| 729 | // handles this case. |
| 730 | |
| 731 | stride_ys_e *= sizeof(int32_t); |
| 732 | stride_ys_p = T * sizeof(int32_t); // Packing dimension |
| 733 | stride_ys_t = sizeof(int32_t); |
| 734 | stride_ys_g = 1; |
| 735 | |
| 736 | LAUNCH_ON_H(uint8_t, stride_ys_e, stride_ys_t, stride_ys_g, stride_ys_p, |
| 737 | true); |
| 738 | |
| 739 | #endif |
| 740 | } |