main
py 1,702 lines 47.7 KB
Raw
1 # SPDX-License-Identifier: Apache-2.0
2 # SPDX-FileCopyrightText: Copyright contributors to the vLLM project
3 # SPDX-FileCopyrightText: Songlin Yang, Yu Zhang
4 #
5 # This file contains code copied from the flash-linear-attention project.
6 # The original source code was licensed under the MIT license and included
7 # the following copyright notice:
8 # Copyright (c) 2023-2025, Songlin Yang, Yu Zhang
9 # ruff: noqa: E501
10
11
12 import torch
13 import torch.nn as nn
14
15 from vllm.model_executor.custom_op import CustomOp
16 from vllm.triton_utils import tl, triton
17 from vllm.utils.math_utils import RCP_LN2, cdiv, next_power_of_2
18
19 from .chunk_delta_h import chunk_gated_delta_rule_fwd_h
20 from .cumsum import chunk_local_cumsum
21 from .fused_recurrent import fused_recurrent_gated_delta_rule_fwd_kernel
22 from .index import prepare_chunk_indices
23 from .l2norm import l2norm_fwd
24 from .op import exp2, log
25 from .solve_tril import solve_tril
26 from .utils import FLA_CHUNK_SIZE, is_amd
27
28 BT_LIST_AUTOTUNE = [32, 64, 128]
29 NUM_WARPS_AUTOTUNE = [2, 4, 8, 16] if is_amd else [4, 8, 16, 32]
30
31
32 def fused_recurrent_kda_fwd(
33 q: torch.Tensor,
34 k: torch.Tensor,
35 v: torch.Tensor,
36 g: torch.Tensor,
37 beta: torch.Tensor,
38 scale: float,
39 initial_state: torch.Tensor,
40 inplace_final_state: bool = True,
41 cu_seqlens: torch.Tensor | None = None,
42 ssm_state_indices: torch.Tensor | None = None,
43 num_accepted_tokens: torch.Tensor | None = None,
44 use_qk_l2norm_in_kernel: bool = False,
45 ) -> tuple[torch.Tensor, torch.Tensor]:
46 B, T, H, K, V = *k.shape, v.shape[-1]
47 HV = v.shape[2]
48 N = B if cu_seqlens is None else len(cu_seqlens) - 1
49 BK, BV = next_power_of_2(K), min(next_power_of_2(V), 8)
50 NK, NV = cdiv(K, BK), cdiv(V, BV)
51 assert NK == 1, "NK > 1 is not supported yet"
52 num_stages = 3
53 num_warps = 1
54
55 o = torch.empty_like(k)
56 if inplace_final_state:
57 final_state = initial_state
58 else:
59 final_state = q.new_empty(T, HV, V, K, dtype=initial_state.dtype)
60
61 stride_init_state_token = initial_state.stride(0)
62 stride_final_state_token = final_state.stride(0)
63
64 if ssm_state_indices is None:
65 stride_indices_seq, stride_indices_tok = 1, 1
66 elif ssm_state_indices.ndim == 1:
67 stride_indices_seq, stride_indices_tok = ssm_state_indices.stride(0), 1
68 else:
69 stride_indices_seq, stride_indices_tok = ssm_state_indices.stride()
70
71 grid = (NK, NV, N * HV)
72 fused_recurrent_gated_delta_rule_fwd_kernel[grid](
73 q=q,
74 k=k,
75 v=v,
76 g=g,
77 beta=beta,
78 o=o,
79 h0=initial_state,
80 ht=final_state,
81 cu_seqlens=cu_seqlens,
82 ssm_state_indices=ssm_state_indices,
83 num_accepted_tokens=num_accepted_tokens,
84 scale=scale,
85 N=N,
86 T=T,
87 B=B,
88 H=H,
89 HV=HV,
90 K=K,
91 V=V,
92 BK=BK,
93 BV=BV,
94 stride_init_state_token=stride_init_state_token,
95 stride_final_state_token=stride_final_state_token,
96 stride_indices_seq=stride_indices_seq,
97 stride_indices_tok=stride_indices_tok,
98 IS_BETA_HEADWISE=beta.ndim == v.ndim,
99 USE_QK_L2NORM_IN_KERNEL=use_qk_l2norm_in_kernel,
100 INPLACE_FINAL_STATE=inplace_final_state,
101 IS_KDA=True,
102 num_warps=num_warps,
103 num_stages=num_stages,
104 )
105
106 return o, final_state
107
108
109 def fused_recurrent_kda(
110 q: torch.Tensor,
111 k: torch.Tensor,
112 v: torch.Tensor,
113 g: torch.Tensor,
114 beta: torch.Tensor = None,
115 scale: float = None,
116 initial_state: torch.Tensor = None,
117 inplace_final_state: bool = True,
118 use_qk_l2norm_in_kernel: bool = True,
119 cu_seqlens: torch.Tensor | None = None,
120 ssm_state_indices: torch.LongTensor | None = None,
121 **kwargs,
122 ) -> tuple[torch.Tensor, torch.Tensor]:
123 if cu_seqlens is not None and q.shape[0] != 1:
124 raise ValueError(
125 f"The batch size is expected to be 1 rather than {q.shape[0]} when using `cu_seqlens`."
126 f"Please flatten variable-length inputs before processing."
127 )
128 if scale is None:
129 scale = k.shape[-1] ** -0.5
130
131 o, final_state = fused_recurrent_kda_fwd(
132 q=q.contiguous(),
133 k=k.contiguous(),
134 v=v.contiguous(),
135 g=g.contiguous(),
136 beta=beta.contiguous(),
137 scale=scale,
138 initial_state=initial_state,
139 inplace_final_state=inplace_final_state,
140 cu_seqlens=cu_seqlens,
141 ssm_state_indices=ssm_state_indices,
142 num_accepted_tokens=None,
143 use_qk_l2norm_in_kernel=use_qk_l2norm_in_kernel,
144 )
145 return o, final_state
146
147
148 @triton.heuristics(
149 {
150 "STORE_RESIDUAL_OUT": lambda args: args["residual_out"] is not None,
151 "HAS_RESIDUAL": lambda args: args["residual"] is not None,
152 "HAS_WEIGHT": lambda args: args["w"] is not None,
153 "HAS_BIAS": lambda args: args["b"] is not None,
154 }
155 )
156 @triton.jit
157 def layer_norm_gated_fwd_kernel(
158 x, # pointer to the input
159 g, # pointer to the gate
160 y, # pointer to the output
161 w, # pointer to the weights
162 b, # pointer to the biases
163 residual, # pointer to the residual
164 residual_out, # pointer to the residual
165 mean, # pointer to the mean
166 rstd, # pointer to the 1/std
167 eps, # epsilon to avoid division by zero
168 T, # number of rows in x
169 H: tl.constexpr, # number of heads
170 g_stride_n: tl.constexpr,
171 D: tl.constexpr, # number of columns in x
172 BT: tl.constexpr,
173 BD: tl.constexpr,
174 ACTIVATION: tl.constexpr,
175 IS_RMS_NORM: tl.constexpr,
176 STORE_RESIDUAL_OUT: tl.constexpr,
177 HAS_RESIDUAL: tl.constexpr,
178 HAS_WEIGHT: tl.constexpr,
179 HAS_BIAS: tl.constexpr,
180 ):
181 i_t = tl.program_id(0)
182
183 o_d = tl.arange(0, BD)
184 m_d = o_d < D
185
186 p_x = tl.make_block_ptr(x, (T, D), (D, 1), (i_t * BT, 0), (BT, BD), (1, 0))
187 b_x = tl.load(p_x, boundary_check=(0, 1)).to(tl.float32)
188 if HAS_RESIDUAL:
189 p_res = tl.make_block_ptr(
190 residual, (T, D), (D, 1), (i_t * BT, 0), (BT, BD), (1, 0)
191 )
192 b_x += tl.load(p_res, boundary_check=(0, 1)).to(tl.float32)
193 if STORE_RESIDUAL_OUT:
194 p_res_out = tl.make_block_ptr(
195 residual_out, (T, D), (D, 1), (i_t * BT, 0), (BT, BD), (1, 0)
196 )
197 tl.store(p_res_out, b_x.to(p_res_out.dtype.element_ty), boundary_check=(0, 1))
198 if not IS_RMS_NORM:
199 b_mean = tl.sum(b_x, axis=1) / D
200 p_mean = tl.make_block_ptr(mean, (T,), (1,), (i_t * BT,), (BT,), (0,))
201 tl.store(p_mean, b_mean.to(p_mean.dtype.element_ty), boundary_check=(0,))
202 b_xbar = tl.where(m_d[None, :], b_x - b_mean[:, None], 0.0)
203 b_var = tl.sum(b_xbar * b_xbar, axis=1) / D
204 else:
205 b_xbar = tl.where(m_d[None, :], b_x, 0.0)
206 b_var = tl.sum(b_xbar * b_xbar, axis=1) / D
207 b_rstd = 1 / tl.sqrt(b_var + eps)
208
209 p_rstd = tl.make_block_ptr(rstd, (T,), (1,), (i_t * BT,), (BT,), (0,))
210 tl.store(p_rstd, b_rstd.to(p_rstd.dtype.element_ty), boundary_check=(0,))
211
212 if HAS_WEIGHT:
213 b_w = tl.load(w + o_d, mask=m_d).to(tl.float32)
214 if HAS_BIAS:
215 b_b = tl.load(b + o_d, mask=m_d).to(tl.float32)
216 b_x_hat = (
217 (b_x - b_mean[:, None]) * b_rstd[:, None]
218 if not IS_RMS_NORM
219 else b_x * b_rstd[:, None]
220 )
221 b_y = b_x_hat * b_w[None, :] if HAS_WEIGHT else b_x_hat
222 if HAS_BIAS:
223 b_y = b_y + b_b[None, :]
224
225 # swish/sigmoid output gate
226 o_t = i_t * BT + tl.arange(0, BT)
227 o_g = (o_t // H) * g_stride_n + (o_t % H) * D
228 b_g = tl.load(
229 g + o_g[:, None] + o_d[None, :],
230 mask=(o_t[:, None] < T) & m_d[None, :],
231 other=0.0,
232 ).to(tl.float32)
233 if ACTIVATION == "swish" or ACTIVATION == "silu":
234 b_y = b_y * b_g * tl.sigmoid(b_g)
235 elif ACTIVATION == "sigmoid":
236 b_y = b_y * tl.sigmoid(b_g)
237
238 # Write output
239 p_y = tl.make_block_ptr(y, (T, D), (D, 1), (i_t * BT, 0), (BT, BD), (1, 0))
240 tl.store(p_y, b_y.to(p_y.dtype.element_ty), boundary_check=(0, 1))
241
242
243 @triton.heuristics(
244 {
245 "STORE_RESIDUAL_OUT": lambda args: args["residual_out"] is not None,
246 "HAS_RESIDUAL": lambda args: args["residual"] is not None,
247 "HAS_WEIGHT": lambda args: args["w"] is not None,
248 "HAS_BIAS": lambda args: args["b"] is not None,
249 }
250 )
251 @triton.jit
252 def layer_norm_gated_fwd_kernel1(
253 x, # pointer to the input
254 g, # pointer to the gate
255 y, # pointer to the output
256 w, # pointer to the weights
257 b, # pointer to the biases
258 residual, # pointer to the residual
259 residual_out, # pointer to the residual
260 mean, # pointer to the mean
261 rstd, # pointer to the 1/std
262 eps, # epsilon to avoid division by zero
263 D: tl.constexpr, # number of columns in x
264 BD: tl.constexpr,
265 ACTIVATION: tl.constexpr,
266 IS_RMS_NORM: tl.constexpr,
267 STORE_RESIDUAL_OUT: tl.constexpr,
268 HAS_RESIDUAL: tl.constexpr,
269 HAS_WEIGHT: tl.constexpr,
270 HAS_BIAS: tl.constexpr,
271 ):
272 i_t = tl.program_id(0)
273 x += i_t * D
274 y += i_t * D
275 g += i_t * D
276 if HAS_RESIDUAL:
277 residual += i_t * D
278 if STORE_RESIDUAL_OUT:
279 residual_out += i_t * D
280
281 o_d = tl.arange(0, BD)
282 m_d = o_d < D
283 b_x = tl.load(x + o_d, mask=m_d, other=0.0).to(tl.float32)
284 if HAS_RESIDUAL:
285 b_x += tl.load(residual + o_d, mask=m_d, other=0.0).to(tl.float32)
286 if STORE_RESIDUAL_OUT:
287 tl.store(residual_out + o_d, b_x, mask=m_d)
288 if not IS_RMS_NORM:
289 b_mean = tl.sum(b_x, axis=0) / D
290 tl.store(mean + i_t, b_mean)
291 b_xbar = tl.where(m_d, b_x - b_mean, 0.0)
292 b_var = tl.sum(b_xbar * b_xbar, axis=0) / D
293 else:
294 b_xbar = tl.where(m_d, b_x, 0.0)
295 b_var = tl.sum(b_xbar * b_xbar, axis=0) / D
296 b_rstd = 1 / tl.sqrt(b_var + eps)
297 tl.store(rstd + i_t, b_rstd)
298
299 if HAS_WEIGHT:
300 b_w = tl.load(w + o_d, mask=m_d).to(tl.float32)
301 if HAS_BIAS:
302 b_b = tl.load(b + o_d, mask=m_d).to(tl.float32)
303 b_x_hat = (b_x - b_mean) * b_rstd if not IS_RMS_NORM else b_x * b_rstd
304 b_y = b_x_hat * b_w if HAS_WEIGHT else b_x_hat
305 if HAS_BIAS:
306 b_y = b_y + b_b
307
308 # swish/sigmoid output gate
309 b_g = tl.load(g + o_d, mask=m_d, other=0.0).to(tl.float32)
310 if ACTIVATION == "swish" or ACTIVATION == "silu":
311 b_y = b_y * b_g * tl.sigmoid(b_g)
312 elif ACTIVATION == "sigmoid":
313 b_y = b_y * tl.sigmoid(b_g)
314
315 # Write output
316 tl.store(y + o_d, b_y, mask=m_d)
317
318
319 def layer_norm_gated_fwd(
320 x: torch.Tensor,
321 g: torch.Tensor,
322 weight: torch.Tensor,
323 bias: torch.Tensor,
324 activation: str = "swish",
325 eps: float = 1e-5,
326 residual: torch.Tensor = None,
327 out_dtype: torch.dtype = None,
328 residual_dtype: torch.dtype = None,
329 is_rms_norm: bool = False,
330 H: int = 1,
331 g_stride_n: int | None = None,
332 ):
333 if residual is not None:
334 residual_dtype = residual.dtype
335 T, D = x.shape
336 if g_stride_n is None:
337 g_stride_n = D
338 assert T % H == 0
339 if residual is not None:
340 assert residual.shape == (T, D)
341 if weight is not None:
342 assert weight.shape == (D,)
343 if bias is not None:
344 assert bias.shape == (D,)
345 # allocate output
346 y = x if out_dtype is None else torch.empty_like(x, dtype=out_dtype)
347 if residual is not None or (
348 residual_dtype is not None and residual_dtype != x.dtype
349 ):
350 residual_out = torch.empty(T, D, device=x.device, dtype=residual_dtype)
351 else:
352 residual_out = None
353 mean = (
354 torch.empty((T,), dtype=torch.float, device=x.device)
355 if not is_rms_norm
356 else None
357 )
358 rstd = torch.empty((T,), dtype=torch.float, device=x.device)
359 # Less than 64KB per feature: enqueue fused kernel
360 MAX_FUSED_SIZE = 65536 // x.element_size()
361 BD = min(MAX_FUSED_SIZE, next_power_of_2(D))
362 if D > BD:
363 raise RuntimeError("This layer norm doesn't support feature dim >= 64KB.")
364 if D <= 512:
365 BT = 16
366 layer_norm_gated_fwd_kernel[(cdiv(T, BT),)](
367 x=x,
368 g=g,
369 y=y,
370 w=weight,
371 b=bias,
372 residual=residual,
373 residual_out=residual_out,
374 mean=mean,
375 rstd=rstd,
376 eps=eps,
377 T=T,
378 H=H,
379 g_stride_n=g_stride_n,
380 D=D,
381 BD=BD,
382 BT=BT,
383 ACTIVATION=activation,
384 IS_RMS_NORM=is_rms_norm,
385 num_warps=8,
386 )
387 else:
388 layer_norm_gated_fwd_kernel1[(T,)](
389 x=x,
390 g=g,
391 y=y,
392 w=weight,
393 b=bias,
394 residual=residual,
395 residual_out=residual_out,
396 mean=mean,
397 rstd=rstd,
398 eps=eps,
399 D=D,
400 BD=BD,
401 ACTIVATION=activation,
402 IS_RMS_NORM=is_rms_norm,
403 num_warps=4,
404 )
405 # residual_out is None if residual is None and residual_dtype == input_dtype
406 return y, mean, rstd, residual_out if residual_out is not None else x
407
408
409 def rms_norm_gated(
410 x: torch.Tensor,
411 g: torch.Tensor,
412 weight: torch.Tensor,
413 bias: torch.Tensor,
414 activation: str = "swish",
415 residual: torch.Tensor | None = None,
416 prenorm: bool = False,
417 residual_in_fp32: bool = False,
418 eps: float = 1e-6,
419 ):
420 x_shape_og = x.shape
421 # reshape input data into 2D tensor
422 x = x.contiguous().reshape(-1, x.shape[-1])
423 D = x.shape[-1]
424 # The tiled kernel supports row-strided gates; kernel1 does not.
425 if D <= 512:
426 H = 1 if g.ndim == 2 else g.shape[-2]
427 g = g.view(-1, H, D)
428 g_stride_n = g.stride(0)
429 else:
430 g = g.contiguous()
431 H = 1
432 g_stride_n = D
433 if residual is not None:
434 assert residual.shape == x_shape_og
435 residual = residual.contiguous().reshape(-1, residual.shape[-1])
436 residual_dtype = (
437 residual.dtype
438 if residual is not None
439 else (torch.float if residual_in_fp32 else None)
440 )
441 y, _, _, residual_out = layer_norm_gated_fwd(
442 x=x,
443 g=g,
444 weight=weight,
445 bias=bias,
446 activation=activation,
447 eps=eps,
448 residual=residual,
449 residual_dtype=residual_dtype,
450 is_rms_norm=True,
451 H=H,
452 g_stride_n=g_stride_n,
453 )
454 y = y.reshape(x_shape_og)
455 return y if not prenorm else (y, residual_out.reshape(x_shape_og))
456
457
458 @CustomOp.register("fused_rms_norm_gated")
459 class FusedRMSNormGated(CustomOp):
460 def __init__(
461 self,
462 hidden_size: int,
463 elementwise_affine: bool = True,
464 eps: float = 1e-5,
465 activation: str = "swish",
466 device: torch.device | None = None,
467 dtype: torch.dtype | None = None,
468 ) -> None:
469 factory_kwargs = {"device": device, "dtype": dtype}
470 super().__init__()
471
472 self.hidden_size = hidden_size
473 self.elementwise_affine = elementwise_affine
474 self.eps = eps
475 self.activation = activation
476
477 if self.activation not in ["swish", "silu", "sigmoid"]:
478 raise ValueError(f"Unsupported activation: {self.activation}")
479
480 if elementwise_affine:
481 self.weight = nn.Parameter(torch.empty(hidden_size, **factory_kwargs))
482 else:
483 self.register_parameter("weight", None)
484 self.register_parameter("bias", None)
485
486 def forward_native(
487 self,
488 x: torch.Tensor,
489 g: torch.Tensor,
490 residual: torch.Tensor | None = None,
491 prenorm: bool = False,
492 residual_in_fp32: bool = False,
493 ) -> torch.Tensor:
494 """Decomposed PyTorch ops for torch.compile/inductor fusion."""
495 # TODO(https://github.com/vllm-project/vllm/issues/36175): implement
496 # native residual/prenorm path and unify with RMSNormGated.
497 # For now, fall back to the triton kernel.
498 if residual is not None or prenorm:
499 return self.forward_cuda(x, g, residual, prenorm, residual_in_fp32)
500 x_float = x.float()
501 variance = x_float.pow(2).mean(dim=-1, keepdim=True)
502 x_normed = x_float * torch.rsqrt(variance + self.eps)
503 if self.weight is not None:
504 x_normed = x_normed * self.weight.float()
505 g_float = g.float()
506 if self.activation in ("swish", "silu"):
507 out = x_normed * g_float * torch.sigmoid(g_float)
508 else: # sigmoid
509 out = x_normed * torch.sigmoid(g_float)
510 return out.to(x.dtype)
511
512 def forward_cuda(
513 self,
514 x: torch.Tensor,
515 g: torch.Tensor,
516 residual: torch.Tensor | None = None,
517 prenorm: bool = False,
518 residual_in_fp32: bool = False,
519 ) -> torch.Tensor:
520 return rms_norm_gated(
521 x,
522 g,
523 self.weight,
524 self.bias,
525 self.activation,
526 residual=residual,
527 eps=self.eps,
528 prenorm=prenorm,
529 residual_in_fp32=residual_in_fp32,
530 )
531
532 def forward_xpu(
533 self,
534 x: torch.Tensor,
535 g: torch.Tensor,
536 residual: torch.Tensor | None = None,
537 prenorm: bool = False,
538 residual_in_fp32: bool = False,
539 ) -> torch.Tensor:
540 return self.forward_cuda(
541 x,
542 g,
543 residual=residual,
544 prenorm=prenorm,
545 residual_in_fp32=residual_in_fp32,
546 )
547
548
549 @triton.heuristics({"IS_VARLEN": lambda args: args["cu_seqlens"] is not None})
550 @triton.autotune(
551 configs=[
552 triton.Config({"BK": BK}, num_warps=num_warps, num_stages=num_stages)
553 for BK in [32, 64]
554 for num_warps in [1, 2, 4, 8]
555 for num_stages in [2, 3, 4]
556 ],
557 key=["BC"],
558 )
559 @triton.jit(do_not_specialize=["T"])
560 def chunk_kda_scaled_dot_kkt_fwd_kernel_intra_sub_inter(
561 q,
562 k,
563 g,
564 beta,
565 A,
566 Aqk,
567 scale,
568 cu_seqlens,
569 chunk_indices,
570 T,
571 H: tl.constexpr,
572 K: tl.constexpr,
573 BT: tl.constexpr,
574 BC: tl.constexpr,
575 BK: tl.constexpr,
576 NC: tl.constexpr,
577 IS_VARLEN: tl.constexpr,
578 ):
579 i_t, i_c, i_bh = tl.program_id(0), tl.program_id(1), tl.program_id(2)
580 i_b, i_h = i_bh // H, i_bh % H
581 i_i, i_j = i_c // NC, i_c % NC
582 if IS_VARLEN:
583 i_n, i_t = (
584 tl.load(chunk_indices + i_t * 2).to(tl.int32),
585 tl.load(chunk_indices + i_t * 2 + 1).to(tl.int32),
586 )
587 bos, eos = (
588 tl.load(cu_seqlens + i_n).to(tl.int32),
589 tl.load(cu_seqlens + i_n + 1).to(tl.int32),
590 )
591 T = eos - bos
592 else:
593 bos, eos = i_b * T, i_b * T + T
594
595 if i_t * BT + i_i * BC >= T:
596 return
597 if i_i <= i_j:
598 return
599
600 q += (bos * H + i_h) * K
601 k += (bos * H + i_h) * K
602 g += (bos * H + i_h) * K
603 A += (bos * H + i_h) * BT
604 Aqk += (bos * H + i_h) * BT
605
606 p_b = tl.make_block_ptr(
607 beta + bos * H + i_h, (T,), (H,), (i_t * BT + i_i * BC,), (BC,), (0,)
608 )
609 b_b = tl.load(p_b, boundary_check=(0,))
610
611 b_A = tl.zeros([BC, BC], dtype=tl.float32)
612 b_Aqk = tl.zeros([BC, BC], dtype=tl.float32)
613 for i_k in range(tl.cdiv(K, BK)):
614 p_q = tl.make_block_ptr(
615 q, (T, K), (H * K, 1), (i_t * BT + i_i * BC, i_k * BK), (BC, BK), (1, 0)
616 )
617 p_k = tl.make_block_ptr(
618 k, (T, K), (H * K, 1), (i_t * BT + i_i * BC, i_k * BK), (BC, BK), (1, 0)
619 )
620 p_g = tl.make_block_ptr(
621 g, (T, K), (H * K, 1), (i_t * BT + i_i * BC, i_k * BK), (BC, BK), (1, 0)
622 )
623 b_kt = tl.make_block_ptr(
624 k, (K, T), (1, H * K), (i_k * BK, i_t * BT + i_j * BC), (BK, BC), (0, 1)
625 )
626 p_gk = tl.make_block_ptr(
627 g, (K, T), (1, H * K), (i_k * BK, i_t * BT + i_j * BC), (BK, BC), (0, 1)
628 )
629
630 o_k = i_k * BK + tl.arange(0, BK)
631 m_k = o_k < K
632 # [BK,]
633 b_gn = tl.load(g + (i_t * BT + i_i * BC) * H * K + o_k, mask=m_k, other=0)
634 # [BC, BK]
635 b_g = tl.load(p_g, boundary_check=(0, 1))
636 b_k = tl.load(p_k, boundary_check=(0, 1)) * exp2(b_g - b_gn[None, :])
637 # [BK, BC]
638 b_gk = tl.load(p_gk, boundary_check=(0, 1))
639 b_kt = tl.load(b_kt, boundary_check=(0, 1))
640 # [BC, BC]
641 b_ktg = b_kt * exp2(b_gn[:, None] - b_gk)
642 b_A += tl.dot(b_k, b_ktg)
643
644 b_q = tl.load(p_q, boundary_check=(0, 1))
645 b_qg = b_q * exp2(b_g - b_gn[None, :]) * scale
646 b_Aqk += tl.dot(b_qg, b_ktg)
647
648 b_A *= b_b[:, None]
649
650 p_A = tl.make_block_ptr(
651 A, (T, BT), (H * BT, 1), (i_t * BT + i_i * BC, i_j * BC), (BC, BC), (1, 0)
652 )
653 tl.store(p_A, b_A.to(A.dtype.element_ty), boundary_check=(0, 1))
654 p_Aqk = tl.make_block_ptr(
655 Aqk, (T, BT), (H * BT, 1), (i_t * BT + i_i * BC, i_j * BC), (BC, BC), (1, 0)
656 )
657 tl.store(p_Aqk, b_Aqk.to(Aqk.dtype.element_ty), boundary_check=(0, 1))
658
659
660 @triton.heuristics({"IS_VARLEN": lambda args: args["cu_seqlens"] is not None})
661 @triton.autotune(
662 configs=[triton.Config({}, num_warps=num_warps) for num_warps in [1, 2, 4, 8]],
663 key=["BK", "BT"],
664 )
665 @triton.jit(do_not_specialize=["T"])
666 def chunk_kda_scaled_dot_kkt_fwd_kernel_intra_sub_intra(
667 q,
668 k,
669 g,
670 beta,
671 A,
672 Aqk,
673 scale,
674 cu_seqlens,
675 chunk_indices,
676 T,
677 H: tl.constexpr,
678 K: tl.constexpr,
679 BT: tl.constexpr,
680 BC: tl.constexpr,
681 BK: tl.constexpr,
682 IS_VARLEN: tl.constexpr,
683 ):
684 i_t, i_i, i_bh = tl.program_id(0), tl.program_id(1), tl.program_id(2)
685 i_b, i_h = i_bh // H, i_bh % H
686 if IS_VARLEN:
687 i_n, i_t = (
688 tl.load(chunk_indices + i_t * 2).to(tl.int32),
689 tl.load(chunk_indices + i_t * 2 + 1).to(tl.int32),
690 )
691 bos, eos = (
692 tl.load(cu_seqlens + i_n).to(tl.int32),
693 tl.load(cu_seqlens + i_n + 1).to(tl.int32),
694 )
695 T = eos - bos
696 else:
697 bos, eos = i_b * T, i_b * T + T
698
699 if i_t * BT + i_i * BC >= T:
700 return
701
702 o_i = tl.arange(0, BC)
703 o_k = tl.arange(0, BK)
704 m_k = o_k < K
705 m_A = (i_t * BT + i_i * BC + o_i) < T
706 o_A = (bos + i_t * BT + i_i * BC + o_i) * H * BT + i_h * BT + i_i * BC
707
708 p_q = tl.make_block_ptr(
709 q + (bos * H + i_h) * K,
710 (T, K),
711 (H * K, 1),
712 (i_t * BT + i_i * BC, 0),
713 (BC, BK),
714 (1, 0),
715 )
716 p_k = tl.make_block_ptr(
717 k + (bos * H + i_h) * K,
718 (T, K),
719 (H * K, 1),
720 (i_t * BT + i_i * BC, 0),
721 (BC, BK),
722 (1, 0),
723 )
724 p_g = tl.make_block_ptr(
725 g + (bos * H + i_h) * K,
726 (T, K),
727 (H * K, 1),
728 (i_t * BT + i_i * BC, 0),
729 (BC, BK),
730 (1, 0),
731 )
732 b_q = tl.load(p_q, boundary_check=(0, 1))
733 b_k = tl.load(p_k, boundary_check=(0, 1))
734 b_g = tl.load(p_g, boundary_check=(0, 1))
735
736 p_b = beta + (bos + i_t * BT + i_i * BC + o_i) * H + i_h
737 b_k = b_k * tl.load(p_b, mask=m_A, other=0)[:, None]
738
739 p_kt = k + (bos + i_t * BT + i_i * BC) * H * K + i_h * K + o_k
740 p_gk = g + (bos + i_t * BT + i_i * BC) * H * K + i_h * K + o_k
741
742 for j in range(0, min(BC, T - i_t * BT - i_i * BC)):
743 b_kt = tl.load(p_kt, mask=m_k, other=0).to(tl.float32)
744 b_gk = tl.load(p_gk, mask=m_k, other=0).to(tl.float32)
745 b_ktg = b_kt[None, :] * exp2(b_g - b_gk[None, :])
746 b_A = tl.sum(b_k * b_ktg, 1)
747 b_A = tl.where(o_i > j, b_A, 0.0)
748 b_Aqk = tl.sum(b_q * b_ktg, 1)
749 b_Aqk = tl.where(o_i >= j, b_Aqk * scale, 0.0)
750 tl.store(A + o_A + j, b_A, mask=m_A)
751 tl.store(Aqk + o_A + j, b_Aqk, mask=m_A)
752 p_kt += H * K
753 p_gk += H * K
754
755
756 def chunk_kda_scaled_dot_kkt_fwd(
757 q: torch.Tensor,
758 k: torch.Tensor,
759 gk: torch.Tensor | None = None,
760 beta: torch.Tensor | None = None,
761 scale: float | None = None,
762 cu_seqlens: torch.Tensor | None = None,
763 chunk_indices: torch.Tensor | None = None,
764 chunk_size: int = FLA_CHUNK_SIZE,
765 output_dtype: torch.dtype = torch.float32,
766 ) -> tuple[torch.Tensor, torch.Tensor]:
767 r"""
768 Compute beta * K * K^T.
769
770 Args:
771 k (torch.Tensor):
772 The key tensor of shape `[B, T, H, K]`.
773 beta (torch.Tensor):
774 The beta tensor of shape `[B, T, H]`.
775 gk (torch.Tensor):
776 The cumulative sum of the gate tensor of shape `[B, T, H, K]` applied to the key tensor. Default: `None`.
777 cu_seqlens (torch.Tensor):
778 The cumulative sequence lengths of the input tensor.
779 Default: None
780 chunk_size (int):
781 The chunk size. Default: 64.
782 output_dtype (torch.dtype):
783 The dtype of the output tensor. Default: `torch.float32`
784
785 Returns:
786 beta * K * K^T of shape `[B, T, H, BT]` where `BT` is the chunk size.
787 """
788 B, T, H, K = k.shape
789 assert K <= 256
790 BT = chunk_size
791 if chunk_indices is None and cu_seqlens is not None:
792 chunk_indices = prepare_chunk_indices(cu_seqlens, BT)
793 NT = cdiv(T, BT) if cu_seqlens is None else len(chunk_indices)
794
795 BC = min(16, BT)
796 NC = cdiv(BT, BC)
797 BK = max(next_power_of_2(K), 16)
798 A = torch.zeros(B, T, H, BT, device=k.device, dtype=output_dtype)
799 Aqk = torch.zeros(B, T, H, BT, device=k.device, dtype=output_dtype)
800 grid = (NT, NC * NC, B * H)
801 chunk_kda_scaled_dot_kkt_fwd_kernel_intra_sub_inter[grid](
802 q=q,
803 k=k,
804 g=gk,
805 beta=beta,
806 A=A,
807 Aqk=Aqk,
808 scale=scale,
809 cu_seqlens=cu_seqlens,
810 chunk_indices=chunk_indices,
811 T=T,
812 H=H,
813 K=K,
814 BT=BT,
815 BC=BC,
816 NC=NC,
817 )
818
819 grid = (NT, NC, B * H)
820 chunk_kda_scaled_dot_kkt_fwd_kernel_intra_sub_intra[grid](
821 q=q,
822 k=k,
823 g=gk,
824 beta=beta,
825 A=A,
826 Aqk=Aqk,
827 scale=scale,
828 cu_seqlens=cu_seqlens,
829 chunk_indices=chunk_indices,
830 T=T,
831 H=H,
832 K=K,
833 BT=BT,
834 BC=BC,
835 BK=BK,
836 )
837 return A, Aqk
838
839
840 @triton.heuristics(
841 {
842 "STORE_QG": lambda args: args["qg"] is not None,
843 "STORE_KG": lambda args: args["kg"] is not None,
844 "IS_VARLEN": lambda args: args["cu_seqlens"] is not None,
845 }
846 )
847 @triton.autotune(
848 configs=[
849 triton.Config({}, num_warps=num_warps, num_stages=num_stages)
850 for num_warps in [2, 4, 8]
851 for num_stages in [2, 3, 4]
852 ],
853 key=["H", "K", "V", "BT", "BK", "BV", "IS_VARLEN"],
854 )
855 @triton.jit(do_not_specialize=["T"])
856 def recompute_w_u_fwd_kernel(
857 q,
858 k,
859 qg,
860 kg,
861 v,
862 beta,
863 w,
864 u,
865 A,
866 gk,
867 cu_seqlens,
868 chunk_indices,
869 T,
870 H: tl.constexpr,
871 K: tl.constexpr,
872 V: tl.constexpr,
873 BT: tl.constexpr,
874 BK: tl.constexpr,
875 BV: tl.constexpr,
876 STORE_QG: tl.constexpr,
877 STORE_KG: tl.constexpr,
878 IS_VARLEN: tl.constexpr,
879 DOT_PRECISION: tl.constexpr,
880 ):
881 i_t, i_bh = tl.program_id(0), tl.program_id(1)
882 i_b, i_h = i_bh // H, i_bh % H
883 if IS_VARLEN:
884 i_n, i_t = (
885 tl.load(chunk_indices + i_t * 2).to(tl.int32),
886 tl.load(chunk_indices + i_t * 2 + 1).to(tl.int32),
887 )
888 bos, eos = (
889 tl.load(cu_seqlens + i_n).to(tl.int32),
890 tl.load(cu_seqlens + i_n + 1).to(tl.int32),
891 )
892 T = eos - bos
893 else:
894 bos, eos = i_b * T, i_b * T + T
895 p_b = tl.make_block_ptr(beta + bos * H + i_h, (T,), (H,), (i_t * BT,), (BT,), (0,))
896 b_b = tl.load(p_b, boundary_check=(0,))
897
898 p_A = tl.make_block_ptr(
899 A + (bos * H + i_h) * BT, (T, BT), (H * BT, 1), (i_t * BT, 0), (BT, BT), (1, 0)
900 )
901 b_A = tl.load(p_A, boundary_check=(0, 1))
902
903 for i_v in range(tl.cdiv(V, BV)):
904 p_v = tl.make_block_ptr(
905 v + (bos * H + i_h) * V,
906 (T, V),
907 (H * V, 1),
908 (i_t * BT, i_v * BV),
909 (BT, BV),
910 (1, 0),
911 )
912 p_u = tl.make_block_ptr(
913 u + (bos * H + i_h) * V,
914 (T, V),
915 (H * V, 1),
916 (i_t * BT, i_v * BV),
917 (BT, BV),
918 (1, 0),
919 )
920 b_v = tl.load(p_v, boundary_check=(0, 1))
921 b_vb = (b_v * b_b[:, None]).to(b_v.dtype)
922 b_u = tl.dot(b_A, b_vb, input_precision=DOT_PRECISION)
923 tl.store(p_u, b_u.to(p_u.dtype.element_ty), boundary_check=(0, 1))
924
925 for i_k in range(tl.cdiv(K, BK)):
926 p_w = tl.make_block_ptr(
927 w + (bos * H + i_h) * K,
928 (T, K),
929 (H * K, 1),
930 (i_t * BT, i_k * BK),
931 (BT, BK),
932 (1, 0),
933 )
934 p_k = tl.make_block_ptr(
935 k + (bos * H + i_h) * K,
936 (T, K),
937 (H * K, 1),
938 (i_t * BT, i_k * BK),
939 (BT, BK),
940 (1, 0),
941 )
942 b_k = tl.load(p_k, boundary_check=(0, 1))
943 b_kb = b_k * b_b[:, None]
944
945 p_gk = tl.make_block_ptr(
946 gk + (bos * H + i_h) * K,
947 (T, K),
948 (H * K, 1),
949 (i_t * BT, i_k * BK),
950 (BT, BK),
951 (1, 0),
952 )
953 b_gk = tl.load(p_gk, boundary_check=(0, 1))
954 b_kb *= exp2(b_gk)
955 if STORE_QG:
956 p_q = tl.make_block_ptr(
957 q + (bos * H + i_h) * K,
958 (T, K),
959 (H * K, 1),
960 (i_t * BT, i_k * BK),
961 (BT, BK),
962 (1, 0),
963 )
964 p_qg = tl.make_block_ptr(
965 qg + (bos * H + i_h) * K,
966 (T, K),
967 (H * K, 1),
968 (i_t * BT, i_k * BK),
969 (BT, BK),
970 (1, 0),
971 )
972 b_q = tl.load(p_q, boundary_check=(0, 1))
973 b_qg = b_q * exp2(b_gk)
974 tl.store(p_qg, b_qg.to(p_qg.dtype.element_ty), boundary_check=(0, 1))
975 if STORE_KG:
976 last_idx = min(i_t * BT + BT, T) - 1
977
978 o_k = i_k * BK + tl.arange(0, BK)
979 m_k = o_k < K
980 b_gn = tl.load(
981 gk + ((bos + last_idx) * H + i_h) * K + o_k, mask=m_k, other=0.0
982 )
983 b_kg = b_k * exp2(b_gn - b_gk)
984
985 p_kg = tl.make_block_ptr(
986 kg + (bos * H + i_h) * K,
987 (T, K),
988 (H * K, 1),
989 (i_t * BT, i_k * BK),
990 (BT, BK),
991 (1, 0),
992 )
993 tl.store(p_kg, b_kg.to(p_kg.dtype.element_ty), boundary_check=(0, 1))
994
995 b_w = tl.dot(b_A, b_kb.to(b_k.dtype))
996 tl.store(p_w, b_w.to(p_w.dtype.element_ty), boundary_check=(0, 1))
997
998
999 def recompute_w_u_fwd(
1000 k: torch.Tensor,
1001 v: torch.Tensor,
1002 beta: torch.Tensor,
1003 A: torch.Tensor,
1004 q: torch.Tensor | None = None,
1005 gk: torch.Tensor | None = None,
1006 cu_seqlens: torch.Tensor | None = None,
1007 chunk_indices: torch.Tensor | None = None,
1008 ) -> tuple[torch.Tensor, torch.Tensor]:
1009 B, T, H, K, V = *k.shape, v.shape[-1]
1010 BT = A.shape[-1]
1011 BK = 64
1012 BV = 64
1013
1014 if chunk_indices is None and cu_seqlens is not None:
1015 chunk_indices = prepare_chunk_indices(cu_seqlens, BT)
1016 NT = cdiv(T, BT) if cu_seqlens is None else len(chunk_indices)
1017
1018 w = torch.empty_like(k)
1019 u = torch.empty_like(v)
1020 kg = torch.empty_like(k) if gk is not None else None
1021 recompute_w_u_fwd_kernel[(NT, B * H)](
1022 q=q,
1023 k=k,
1024 qg=None,
1025 kg=kg,
1026 v=v,
1027 beta=beta,
1028 w=w,
1029 u=u,
1030 A=A,
1031 gk=gk,
1032 cu_seqlens=cu_seqlens,
1033 chunk_indices=chunk_indices,
1034 T=T,
1035 H=H,
1036 K=K,
1037 V=V,
1038 BT=BT,
1039 BK=BK,
1040 BV=BV,
1041 DOT_PRECISION="ieee",
1042 )
1043 return w, u, None, kg
1044
1045
1046 @triton.heuristics({"IS_VARLEN": lambda args: args["cu_seqlens"] is not None})
1047 @triton.autotune(
1048 configs=[
1049 triton.Config({"BK": BK, "BV": BV}, num_warps=num_warps, num_stages=num_stages)
1050 for BK in [32, 64]
1051 for BV in [64, 128]
1052 for num_warps in [2, 4, 8]
1053 for num_stages in [2, 3, 4]
1054 ],
1055 key=["BT"],
1056 )
1057 @triton.jit(do_not_specialize=["T"])
1058 def chunk_gla_fwd_kernel_o(
1059 q,
1060 v,
1061 g,
1062 h,
1063 o,
1064 A,
1065 cu_seqlens,
1066 chunk_indices,
1067 scale,
1068 T,
1069 H: tl.constexpr,
1070 K: tl.constexpr,
1071 V: tl.constexpr,
1072 BT: tl.constexpr,
1073 BK: tl.constexpr,
1074 BV: tl.constexpr,
1075 IS_VARLEN: tl.constexpr,
1076 ):
1077 i_v, i_t, i_bh = tl.program_id(0), tl.program_id(1), tl.program_id(2)
1078 i_b, i_h = i_bh // H, i_bh % H
1079 if IS_VARLEN:
1080 i_tg = i_t
1081 i_n, i_t = (
1082 tl.load(chunk_indices + i_t * 2).to(tl.int32),
1083 tl.load(chunk_indices + i_t * 2 + 1).to(tl.int32),
1084 )
1085 bos, eos = (
1086 tl.load(cu_seqlens + i_n).to(tl.int32),
1087 tl.load(cu_seqlens + i_n + 1).to(tl.int32),
1088 )
1089 T = eos - bos
1090 NT = tl.cdiv(T, BT)
1091 else:
1092 NT = tl.cdiv(T, BT)
1093 i_tg = i_b * NT + i_t
1094 bos, eos = i_b * T, i_b * T + T
1095
1096 m_s = tl.arange(0, BT)[:, None] >= tl.arange(0, BT)[None, :]
1097
1098 b_o = tl.zeros([BT, BV], dtype=tl.float32)
1099 for i_k in range(tl.cdiv(K, BK)):
1100 p_q = tl.make_block_ptr(
1101 q + (bos * H + i_h) * K,
1102 (T, K),
1103 (H * K, 1),
1104 (i_t * BT, i_k * BK),
1105 (BT, BK),
1106 (1, 0),
1107 )
1108 p_g = tl.make_block_ptr(
1109 g + (bos * H + i_h) * K,
1110 (T, K),
1111 (H * K, 1),
1112 (i_t * BT, i_k * BK),
1113 (BT, BK),
1114 (1, 0),
1115 )
1116 p_h = tl.make_block_ptr(
1117 h + (i_tg * H + i_h) * K * V,
1118 (V, K),
1119 (K, 1),
1120 (i_v * BV, i_k * BK),
1121 (BV, BK),
1122 (1, 0),
1123 )
1124
1125 # [BT, BK]
1126 b_q = tl.load(p_q, boundary_check=(0, 1))
1127 b_q = (b_q * scale).to(b_q.dtype)
1128 # [BT, BK]
1129 b_g = tl.load(p_g, boundary_check=(0, 1))
1130 # [BT, BK]
1131 b_qg = (b_q * exp2(b_g)).to(b_q.dtype)
1132 # [BV, BK]
1133 b_h = tl.load(p_h, boundary_check=(0, 1))
1134 # [BT, BV]
1135 if i_k >= 0:
1136 b_o += tl.dot(b_qg, tl.trans(b_h).to(b_qg.dtype))
1137 p_v = tl.make_block_ptr(
1138 v + (bos * H + i_h) * V,
1139 (T, V),
1140 (H * V, 1),
1141 (i_t * BT, i_v * BV),
1142 (BT, BV),
1143 (1, 0),
1144 )
1145 p_o = tl.make_block_ptr(
1146 o + (bos * H + i_h) * V,
1147 (T, V),
1148 (H * V, 1),
1149 (i_t * BT, i_v * BV),
1150 (BT, BV),
1151 (1, 0),
1152 )
1153 p_A = tl.make_block_ptr(
1154 A + (bos * H + i_h) * BT, (T, BT), (H * BT, 1), (i_t * BT, 0), (BT, BT), (1, 0)
1155 )
1156 # [BT, BV]
1157 b_v = tl.load(p_v, boundary_check=(0, 1))
1158 # [BT, BT]
1159 b_A = tl.load(p_A, boundary_check=(0, 1))
1160 b_A = tl.where(m_s, b_A, 0.0).to(b_v.dtype)
1161 b_o += tl.dot(b_A, b_v, allow_tf32=False)
1162 tl.store(p_o, b_o.to(p_o.dtype.element_ty), boundary_check=(0, 1))
1163
1164
1165 def chunk_gla_fwd_o_gk(
1166 q: torch.Tensor,
1167 v: torch.Tensor,
1168 g: torch.Tensor,
1169 A: torch.Tensor,
1170 h: torch.Tensor,
1171 o: torch.Tensor,
1172 scale: float,
1173 cu_seqlens: torch.Tensor | None = None,
1174 chunk_indices: torch.Tensor | None = None,
1175 chunk_size: int = FLA_CHUNK_SIZE,
1176 ):
1177 B, T, H, K, V = *q.shape, v.shape[-1]
1178 BT = chunk_size
1179
1180 if chunk_indices is None and cu_seqlens is not None:
1181 chunk_indices = prepare_chunk_indices(cu_seqlens, chunk_size)
1182 NT = cdiv(T, BT) if cu_seqlens is None else len(chunk_indices)
1183
1184 def grid(meta):
1185 return (cdiv(V, meta["BV"]), NT, B * H)
1186
1187 chunk_gla_fwd_kernel_o[grid](
1188 q=q,
1189 v=v,
1190 g=g,
1191 h=h,
1192 o=o,
1193 A=A,
1194 cu_seqlens=cu_seqlens,
1195 chunk_indices=chunk_indices,
1196 scale=scale,
1197 T=T,
1198 H=H,
1199 K=K,
1200 V=V,
1201 BT=BT,
1202 )
1203 return o
1204
1205
1206 @triton.heuristics(
1207 {
1208 "HAS_BIAS": lambda args: args["g_bias"] is not None,
1209 "IS_VARLEN": lambda args: args["cu_seqlens"] is not None,
1210 }
1211 )
1212 @triton.autotune(
1213 configs=[
1214 triton.Config({"BD": BD}, num_warps=num_warps)
1215 for BD in [32, 64]
1216 for num_warps in [2, 4, 8]
1217 ],
1218 key=["H", "D", "BT", "IS_VARLEN"],
1219 )
1220 @triton.jit(do_not_specialize=["T"])
1221 def kda_gate_cumsum_fwd_kernel(
1222 g,
1223 A,
1224 y,
1225 g_bias,
1226 cu_seqlens,
1227 chunk_indices,
1228 cumsum_scale,
1229 lower_bound,
1230 beta,
1231 threshold,
1232 T,
1233 H: tl.constexpr,
1234 D: tl.constexpr,
1235 BT: tl.constexpr,
1236 BD: tl.constexpr,
1237 HAS_BIAS: tl.constexpr,
1238 IS_VARLEN: tl.constexpr,
1239 USE_LOWER_BOUND: tl.constexpr,
1240 ):
1241 i_d, i_t, i_bh = tl.program_id(0), tl.program_id(1), tl.program_id(2)
1242 i_b, i_h = i_bh // H, i_bh % H
1243 if IS_VARLEN:
1244 i_n, i_t = (
1245 tl.load(chunk_indices + i_t * 2).to(tl.int32),
1246 tl.load(chunk_indices + i_t * 2 + 1).to(tl.int32),
1247 )
1248 bos, eos = (
1249 tl.load(cu_seqlens + i_n).to(tl.int32),
1250 tl.load(cu_seqlens + i_n + 1).to(tl.int32),
1251 )
1252 T = eos - bos
1253 else:
1254 bos = i_b * T
1255
1256 p_g = tl.make_block_ptr(
1257 g + (bos * H + i_h) * D,
1258 (T, D),
1259 (H * D, 1),
1260 (i_t * BT, i_d * BD),
1261 (BT, BD),
1262 (1, 0),
1263 )
1264 p_y = tl.make_block_ptr(
1265 y + (bos * H + i_h) * D,
1266 (T, D),
1267 (H * D, 1),
1268 (i_t * BT, i_d * BD),
1269 (BT, BD),
1270 (1, 0),
1271 )
1272
1273 b_g = tl.load(p_g, boundary_check=(0, 1)).to(tl.float32)
1274 if HAS_BIAS:
1275 o_d = i_d * BD + tl.arange(0, BD)
1276 b_bias = tl.load(g_bias + i_h * D + o_d, mask=o_d < D, other=0.0).to(tl.float32)
1277 b_g = b_g + b_bias[None, :]
1278
1279 b_a = tl.exp(tl.load(A + i_h).to(tl.float32))
1280 if USE_LOWER_BOUND:
1281 b_gate = lower_bound * tl.sigmoid(b_a * b_g)
1282 else:
1283 b_g_scaled = b_g * beta
1284 b_softplus = tl.where(
1285 b_g_scaled > threshold,
1286 b_g,
1287 (1.0 / beta) * log(1.0 + tl.exp(b_g_scaled)),
1288 )
1289 b_gate = -b_a * b_softplus
1290
1291 # Out-of-bounds rows (load returns 0, but softplus/bias can still make
1292 # b_gate non-zero) participate in the dot product. They only contribute to
1293 # out-of-bounds output rows, which are masked away by `boundary_check` on
1294 # the store, so visible output matches unfused gate + chunk-local cumsum.
1295 o_t = tl.arange(0, BT)
1296 m_cumsum = tl.where(o_t[:, None] >= o_t[None, :], 1.0, 0.0)
1297 b_y = tl.dot(m_cumsum, b_gate, allow_tf32=False) * cumsum_scale
1298 tl.store(p_y, b_y.to(p_y.dtype.element_ty), boundary_check=(0, 1))
1299
1300
1301 def fused_kda_gate_chunk_cumsum(
1302 raw_g: torch.Tensor,
1303 A_log: torch.Tensor,
1304 g_bias: torch.Tensor | None = None,
1305 beta: float = 1.0,
1306 threshold: float = 20.0,
1307 lower_bound: float | None = None,
1308 cu_seqlens: torch.Tensor | None = None,
1309 chunk_indices: torch.Tensor | None = None,
1310 chunk_size: int = FLA_CHUNK_SIZE,
1311 output_dtype: torch.dtype | None = torch.float,
1312 ) -> torch.Tensor:
1313 if cu_seqlens is not None:
1314 assert raw_g.shape[0] == 1, (
1315 "Only batch size 1 is supported when cu_seqlens are provided"
1316 )
1317 B, T, H, D = raw_g.shape
1318 if chunk_indices is None and cu_seqlens is not None:
1319 chunk_indices = prepare_chunk_indices(cu_seqlens, chunk_size)
1320 NT = cdiv(T, chunk_size) if cu_seqlens is None else len(chunk_indices)
1321
1322 A_log = A_log.reshape(-1)
1323 if g_bias is not None:
1324 g_bias = g_bias.reshape(-1)
1325 y = torch.empty_like(raw_g, dtype=output_dtype or raw_g.dtype)
1326
1327 def grid(meta):
1328 return (cdiv(meta["D"], meta["BD"]), NT, B * H)
1329
1330 kda_gate_cumsum_fwd_kernel[grid](
1331 g=raw_g,
1332 A=A_log,
1333 y=y,
1334 g_bias=g_bias,
1335 cu_seqlens=cu_seqlens,
1336 chunk_indices=chunk_indices,
1337 # RCP_LN2 folds in the natural-log -> log2 conversion so downstream
1338 # exp2-based kernels reproduce exp(g). Keep this in sync with the
1339 # `use_exp2=True` path in `_chunk_kda_fwd_with_cumulative_g`.
1340 cumsum_scale=RCP_LN2,
1341 lower_bound=lower_bound or 0.0,
1342 beta=beta,
1343 threshold=threshold,
1344 T=T,
1345 H=H,
1346 D=D,
1347 BT=chunk_size,
1348 USE_LOWER_BOUND=lower_bound is not None,
1349 )
1350 return y
1351
1352
1353 def _chunk_kda_fwd_with_cumulative_g(
1354 q: torch.Tensor,
1355 k: torch.Tensor,
1356 v: torch.Tensor,
1357 g: torch.Tensor,
1358 beta: torch.Tensor,
1359 scale: float,
1360 initial_state: torch.Tensor,
1361 output_final_state: bool,
1362 cu_seqlens: torch.Tensor | None = None,
1363 chunk_indices: torch.Tensor | None = None,
1364 chunk_size: int = FLA_CHUNK_SIZE,
1365 ):
1366 # `g` must already be chunk-local cumulatively-summed AND scaled by
1367 # RCP_LN2 (so the downstream exp2-based kernels reproduce exp(g)).
1368 # Use `chunk_kda_fwd` or `chunk_kda_with_fused_gate_fwd` instead of
1369 # calling this helper directly unless that invariant is upheld.
1370 # the intra Aqk is kept in fp32
1371 # the computation has very marginal effect on the entire throughput
1372 A, Aqk = chunk_kda_scaled_dot_kkt_fwd(
1373 q=q,
1374 k=k,
1375 gk=g,
1376 beta=beta,
1377 scale=scale,
1378 cu_seqlens=cu_seqlens,
1379 chunk_indices=chunk_indices,
1380 chunk_size=chunk_size,
1381 output_dtype=torch.float32,
1382 )
1383 A = solve_tril(A=A, cu_seqlens=cu_seqlens, output_dtype=k.dtype)
1384 w, u, _, kg = recompute_w_u_fwd(
1385 k=k,
1386 v=v,
1387 beta=beta,
1388 A=A,
1389 gk=g,
1390 cu_seqlens=cu_seqlens,
1391 chunk_indices=chunk_indices,
1392 )
1393 del A
1394 h, v_new, final_state = chunk_gated_delta_rule_fwd_h(
1395 k=kg,
1396 w=w,
1397 u=u,
1398 gk=g,
1399 initial_state=initial_state,
1400 output_final_state=output_final_state,
1401 cu_seqlens=cu_seqlens,
1402 chunk_indices=chunk_indices,
1403 use_exp2=True,
1404 )
1405 del w, u, kg
1406 o = chunk_gla_fwd_o_gk(
1407 q=q,
1408 v=v_new,
1409 g=g,
1410 A=Aqk,
1411 h=h,
1412 o=v,
1413 scale=scale,
1414 cu_seqlens=cu_seqlens,
1415 chunk_indices=chunk_indices,
1416 chunk_size=chunk_size,
1417 )
1418 del Aqk, v_new, h
1419 return o, final_state
1420
1421
1422 def chunk_kda_fwd(
1423 q: torch.Tensor,
1424 k: torch.Tensor,
1425 v: torch.Tensor,
1426 g: torch.Tensor,
1427 beta: torch.Tensor,
1428 scale: float,
1429 initial_state: torch.Tensor,
1430 output_final_state: bool,
1431 cu_seqlens: torch.Tensor | None = None,
1432 ):
1433 chunk_size = FLA_CHUNK_SIZE
1434 chunk_indices = (
1435 prepare_chunk_indices(cu_seqlens, chunk_size)
1436 if cu_seqlens is not None
1437 else None
1438 )
1439 g = chunk_local_cumsum(
1440 g,
1441 chunk_size=chunk_size,
1442 cu_seqlens=cu_seqlens,
1443 chunk_indices=chunk_indices,
1444 )
1445 # KDA evaluates cumulative gate decays with exp2. Convert from natural-log
1446 # space so exp(x) is preserved as exp2(x / ln(2)).
1447 g = g * RCP_LN2
1448 return _chunk_kda_fwd_with_cumulative_g(
1449 q=q,
1450 k=k,
1451 v=v,
1452 g=g,
1453 beta=beta,
1454 scale=scale,
1455 initial_state=initial_state,
1456 output_final_state=output_final_state,
1457 cu_seqlens=cu_seqlens,
1458 chunk_indices=chunk_indices,
1459 chunk_size=chunk_size,
1460 )
1461
1462
1463 def chunk_kda_with_fused_gate_fwd(
1464 q: torch.Tensor,
1465 k: torch.Tensor,
1466 v: torch.Tensor,
1467 raw_g: torch.Tensor,
1468 beta: torch.Tensor,
1469 A_log: torch.Tensor,
1470 g_bias: torch.Tensor | None,
1471 scale: float,
1472 initial_state: torch.Tensor,
1473 output_final_state: bool,
1474 lower_bound: float | None = None,
1475 cu_seqlens: torch.Tensor | None = None,
1476 ):
1477 chunk_size = FLA_CHUNK_SIZE
1478 chunk_indices = (
1479 prepare_chunk_indices(cu_seqlens, chunk_size)
1480 if cu_seqlens is not None
1481 else None
1482 )
1483 g = fused_kda_gate_chunk_cumsum(
1484 raw_g,
1485 A_log=A_log,
1486 g_bias=g_bias,
1487 cu_seqlens=cu_seqlens,
1488 chunk_indices=chunk_indices,
1489 chunk_size=chunk_size,
1490 lower_bound=lower_bound,
1491 )
1492 return _chunk_kda_fwd_with_cumulative_g(
1493 q=q,
1494 k=k,
1495 v=v,
1496 g=g,
1497 beta=beta,
1498 scale=scale,
1499 initial_state=initial_state,
1500 output_final_state=output_final_state,
1501 cu_seqlens=cu_seqlens,
1502 chunk_indices=chunk_indices,
1503 chunk_size=chunk_size,
1504 )
1505
1506
1507 def chunk_kda(
1508 q: torch.Tensor,
1509 k: torch.Tensor,
1510 v: torch.Tensor,
1511 g: torch.Tensor,
1512 beta: torch.Tensor,
1513 scale: float = None,
1514 initial_state: torch.Tensor = None,
1515 output_final_state: bool = False,
1516 use_qk_l2norm_in_kernel: bool = False,
1517 cu_seqlens: torch.Tensor | None = None,
1518 **kwargs,
1519 ):
1520 if scale is None:
1521 scale = k.shape[-1] ** -0.5
1522
1523 if use_qk_l2norm_in_kernel:
1524 q = l2norm_fwd(q.contiguous())
1525 k = l2norm_fwd(k.contiguous())
1526
1527 o, final_state = chunk_kda_fwd(
1528 q=q,
1529 k=k,
1530 v=v.contiguous(),
1531 g=g.contiguous(),
1532 beta=beta.contiguous(),
1533 scale=scale,
1534 initial_state=initial_state.contiguous(),
1535 output_final_state=output_final_state,
1536 cu_seqlens=cu_seqlens,
1537 )
1538 return o, final_state
1539
1540
1541 def chunk_kda_with_fused_gate(
1542 q: torch.Tensor,
1543 k: torch.Tensor,
1544 v: torch.Tensor,
1545 raw_g: torch.Tensor,
1546 beta: torch.Tensor,
1547 A_log: torch.Tensor,
1548 g_bias: torch.Tensor | None,
1549 scale: float | None = None,
1550 initial_state: torch.Tensor | None = None,
1551 output_final_state: bool = False,
1552 lower_bound: float | None = None,
1553 use_qk_l2norm_in_kernel: bool = False,
1554 cu_seqlens: torch.Tensor | None = None,
1555 **kwargs,
1556 ):
1557 """Run chunk KDA from raw gate projection using fused gate+cumsum."""
1558 if scale is None:
1559 scale = k.shape[-1] ** -0.5
1560
1561 if use_qk_l2norm_in_kernel:
1562 q = l2norm_fwd(q.contiguous())
1563 k = l2norm_fwd(k.contiguous())
1564
1565 o, final_state = chunk_kda_with_fused_gate_fwd(
1566 q=q,
1567 k=k,
1568 v=v.contiguous(),
1569 raw_g=raw_g.contiguous(),
1570 beta=beta.contiguous(),
1571 A_log=A_log,
1572 g_bias=g_bias,
1573 scale=scale,
1574 initial_state=initial_state.contiguous() if initial_state is not None else None,
1575 output_final_state=output_final_state,
1576 lower_bound=lower_bound,
1577 cu_seqlens=cu_seqlens,
1578 )
1579 return o, final_state
1580
1581
1582 @triton.autotune(
1583 configs=[
1584 triton.Config({"BT": bt}, num_warps=nw, num_stages=ns)
1585 for bt in BT_LIST_AUTOTUNE
1586 for nw in NUM_WARPS_AUTOTUNE
1587 for ns in [2, 3]
1588 ],
1589 key=["H", "D"],
1590 )
1591 @triton.jit
1592 def kda_gate_fwd_kernel(
1593 g,
1594 A,
1595 y,
1596 g_bias,
1597 lower_bound,
1598 beta: tl.constexpr,
1599 threshold: tl.constexpr,
1600 T,
1601 H,
1602 D: tl.constexpr,
1603 BT: tl.constexpr,
1604 BD: tl.constexpr,
1605 HAS_BIAS: tl.constexpr,
1606 USE_LOWER_BOUND: tl.constexpr,
1607 ):
1608 i_t, i_h = tl.program_id(0), tl.program_id(1)
1609 n_t = i_t * BT
1610
1611 b_a = tl.exp(tl.load(A + i_h).to(tl.float32))
1612
1613 stride_row = H * D
1614 stride_col = 1
1615
1616 g_ptr = tl.make_block_ptr(
1617 base=g + i_h * D,
1618 shape=(T, D),
1619 strides=(stride_row, stride_col),
1620 offsets=(n_t, 0),
1621 block_shape=(BT, BD),
1622 order=(1, 0),
1623 )
1624
1625 y_ptr = tl.make_block_ptr(
1626 base=y + i_h * D,
1627 shape=(T, D),
1628 strides=(stride_row, stride_col),
1629 offsets=(n_t, 0),
1630 block_shape=(BT, BD),
1631 order=(1, 0),
1632 )
1633
1634 b_g = tl.load(g_ptr, boundary_check=(0, 1)).to(tl.float32)
1635
1636 if HAS_BIAS:
1637 n_d = tl.arange(0, BD)
1638 bias_mask = n_d < D
1639 b_bias = tl.load(g_bias + i_h * D + n_d, mask=bias_mask, other=0.0).to(
1640 tl.float32
1641 )
1642 b_g = b_g + b_bias[None, :]
1643
1644 if USE_LOWER_BOUND:
1645 b_y = lower_bound * tl.sigmoid(b_a * b_g)
1646 else:
1647 g_scaled = b_g * beta
1648 use_linear = g_scaled > threshold
1649 sp = tl.where(use_linear, b_g, (1.0 / beta) * log(1.0 + tl.exp(g_scaled)))
1650 b_y = -b_a * sp
1651
1652 tl.store(y_ptr, b_y.to(y.dtype.element_ty), boundary_check=(0, 1))
1653
1654
1655 def fused_kda_gate(
1656 g: torch.Tensor,
1657 A: torch.Tensor,
1658 head_k_dim: int,
1659 g_bias: torch.Tensor | None = None,
1660 beta: float = 1.0,
1661 threshold: float = 20.0,
1662 lower_bound: float | None = None,
1663 ) -> torch.Tensor:
1664 """
1665 Forward pass for KDA gate:
1666 input g: [..., H*D]
1667 param A: [H] or [1, 1, H, 1]
1668 beta: softplus beta parameter
1669 threshold: softplus threshold parameter
1670 return : [..., H, D]
1671 """
1672 orig_shape = g.shape[:-1]
1673
1674 g = g.view(-1, g.shape[-1])
1675 T = g.shape[0]
1676 HD = g.shape[1]
1677 H = A.numel()
1678 assert H * head_k_dim == HD
1679
1680 y = torch.empty_like(g, dtype=torch.float32)
1681
1682 def grid(meta):
1683 return (cdiv(T, meta["BT"]), H)
1684
1685 kda_gate_fwd_kernel[grid](
1686 g,
1687 A,
1688 y,
1689 g_bias,
1690 lower_bound or 0.0,
1691 beta,
1692 threshold,
1693 T,
1694 H,
1695 head_k_dim,
1696 BD=next_power_of_2(head_k_dim),
1697 HAS_BIAS=g_bias is not None,
1698 USE_LOWER_BOUND=lower_bound is not None,
1699 )
1700
1701 y = y.view(*orig_shape, H, head_k_dim)
1702 return y