// Copyright 2019 Google LLC
//
// This source code is licensed under the BSD-style license found in the
// LICENSE file in the root directory of this source tree.

$assert DATATYPE in ["F32", "QC4", "QC8"]
$assert ACTIVATION in ["LINEAR", "RELU", "MINMAX"]
$assert ACTIVATION != "LINEAR" or not WASM
$if DATATYPE == "QC4":
  $assert ACTIVATION == "MINMAX"
  $assert NR % 2 == 0
#include <assert.h>

#include "xnnpack/gemm.h"
#include "xnnpack/math.h"
$if DATATYPE in ["QC4", "QC8"] and NR % 4 != 0:
  #include "xnnpack/unaligned.h"


$RANGE_MR = list(reversed(range(MR))) if INC else list(range(MR))
$DATATYPE_SPEC = {"F32": "f32", "QC8": "f32_qc8w", "QC4": "f32_qc4w"}[DATATYPE]
$MIN_F32 = "__builtin_wasm_min_f32" if WASM else "math_min_f32"
$MAX_F32 = "__builtin_wasm_max_f32" if WASM else "math_max_f32"
$KERNEL = "gemminc" if INC else "gemm"
$SUFFIX = {"LINEAR": "", "RELU": "_relu", "MINMAX": "_minmax"}[ACTIVATION]
$PARAMS = {"LINEAR": "struct xnn_f32_default_params", "RELU": "struct xnn_f32_relu_params", "MINMAX": "union xnn_f32_minmax_params"}[ACTIVATION]
void xnn_${DATATYPE_SPEC}_${KERNEL}${SUFFIX}_ukernel_${MR}x${NR}__${"wasm" if WASM else "scalar"}(
    size_t mr,
    size_t nc,
    size_t kc,
    const float* restrict a,
    size_t a_stride,
    $if DATATYPE == "F32":
      const float* restrict w,
    $else:
      const void* restrict w,
    float* restrict c,
    size_t cm_stride,
    size_t cn_stride,
    $if INC:
      const float* restrict acc,
    $if DATATYPE == "QC4":
      const struct xnn_f32_qc4w_minmax_params params[restrict XNN_MIN_ELEMENTS(1)])
    $else:
      const ${PARAMS} params[restrict XNN_MIN_ELEMENTS(1)])
{
  assert(mr != 0);
  assert(mr <= ${MR});
  assert(nc != 0);
  assert(kc != 0);
  assert(kc % sizeof(float) == 0);
  assert(a != NULL);
  assert(w != NULL);
  assert(c != NULL);
  $if INC:
    assert(acc != NULL);

  const float* a0 = a;
  float* c0 = c;
  $for M in range(1, MR):
    const float* a${M} = (const float*) ((uintptr_t) a${M-1} + a_stride);
    float* c${M} = (float*) ((uintptr_t) c${M-1} + cm_stride);
    $if M % 2 == 0:
      if XNN_UNPREDICTABLE(mr <= ${M}) {
        a${M} = a${M-1};
        c${M} = c${M-1};
      }
    $elif M + 1 == MR:
      if XNN_UNPREDICTABLE(mr != ${M+1}) {
        a${M} = a${M-1};
        c${M} = c${M-1};
      }
    $else:
      if XNN_UNPREDICTABLE(mr < ${M+1}) {
        a${M} = a${M-1};
        c${M} = c${M-1};
      }

  $if ACTIVATION == "MINMAX":
    const float vmin = params->scalar.min;
    const float vmax = params->scalar.max;
  $if DATATYPE == "QC4":
    const int32_t vminus_kernel_zero_point = -params->scalar.kernel_zero_point;
  do {
    $if INC:
      $for M in range(MR):
        $for N in range(NR):
          float vacc${M}${N} = acc[${M*NR+N}];
      acc += ${MR*NR};
    $else:
      $if DATATYPE == "F32":
        $for N in range(NR):
          float vacc0${N} = w[${N}];
        w += ${NR};
      $else:
        $for N in range(NR):
          $if NR % 4 == 0:
            float vacc0${N} = ((const float*)w)[${N}];
          $else:
            float vacc0${N} = unaligned_indexed_load_f32(w, ${N});
        w = (const float*) w + ${NR};
      $for M in range(1, MR):
        $for N in range(NR):
          float vacc${M}${N} = vacc0${N};

    size_t k = kc;
    $if DATATYPE == "QC4":
      for (; k >= 2 * sizeof(float); k -= 2 * sizeof(float)) {
        $for M in range(MR):
          const float va${M}0 = *a${M}++;
          const float va${M}1 = *a${M}++;

        $for N in range(NR):
          const uint8_t vbi${N} = ((const uint8_t*) w)[${N}];
        $for N in range(NR):
          const float vb${N}0 = (float) ((int32_t) (vbi${N} & 0xF) + vminus_kernel_zero_point);
        $for N in range(NR):
          const float vb${N}1 = (float) ((int32_t) (vbi${N} >> 4) + vminus_kernel_zero_point);
        w = (const int8_t*) w + ${NR};

        $for M in range(MR):
          $for N in range(NR):
            vacc${M}${N} = math_muladd_f32(va${M}0, vb${N}0, vacc${M}${N});
        $for M in range(MR):
          $for N in range(NR):
            vacc${M}${N} = math_muladd_f32(va${M}1, vb${N}1, vacc${M}${N});
      }
      if XNN_UNLIKELY(k != 0) {
        $for M in range(MR):
          const float va${M} = *a${M}++;

        $for N in range(NR):
          const uint8_t vbi${N} = ((const uint8_t*) w)[${N}];
        $for N in range(NR):
          const float vb${N} = (float) ((int32_t) vbi${N} + vminus_kernel_zero_point);
        w = (const int8_t*) w + ${NR};

        $for M in range(MR):
          $for N in range(NR):
            vacc${M}${N} = math_muladd_f32(va${M}, vb${N}, vacc${M}${N});
      }
    $else:
      do {
        $for M in range(MR):
          const float va${M} = *a${M}++;

        $if DATATYPE == "F32":
          $for N in range(NR):
            const float vb${N} = w[${N}];
          w += ${NR};
        $elif DATATYPE == "QC8":
          $for N in range(NR):
            const float vb${N} = (float) ((const int8_t*) w)[${N}];
          w = (const int8_t*) w + ${NR};

        $for M in range(MR):
          $for N in range(NR):
            vacc${M}${N} = math_muladd_f32(va${M}, vb${N}, vacc${M}${N});

        k -= sizeof(float);
      } while (k != 0);

    $if DATATYPE in ["QC8", "QC4"]:
      $for N in range(NR):
        $if NR % 4 == 0:
          const float vscale${N} = ((const float*)w)[${N}];
        $else:
          const float vscale${N} = unaligned_indexed_load_f32(w, ${N});
      w = (const float*) w + ${NR};
      $for N in range(NR):
        $for M in range(MR):
          vacc${M}${N} *= vscale${N};
    $if ACTIVATION == "MINMAX":
      $for M in range(MR):
        $for N in range(NR):
          vacc${M}${N} = ${MAX_F32}(vacc${M}${N}, vmin);

      $for M in range(MR):
        $for N in range(NR):
          vacc${M}${N} = ${MIN_F32}(vacc${M}${N}, vmax);
    $elif ACTIVATION == "RELU":
      $for M in range(MR):
        $for N in range(NR):
          vacc${M}${N} = ${MAX_F32}(vacc${M}${N}, 0.0f);

    if XNN_LIKELY(nc >= ${NR}) {
      $for M in RANGE_MR:
        $for N in range(NR):
          c${M}[${N}] = vacc${M}${N};
        c${M} = (float*) ((uintptr_t) c${M} + cn_stride);

      $for M in RANGE_MR:
        a${M} = (const void*) ((uintptr_t) a${M} - kc);

      nc -= ${NR};
    } else {
      $for LOG2N in reversed(range(NR.bit_length() - 1)):
        if (nc & ${1 << LOG2N}) {
          $for M in RANGE_MR:
            $for N in range(1 << LOG2N):
              c${M}[${N}] = vacc${M}${N};
            $if LOG2N != 0:
              $for N in range(NR - (1 << LOG2N) - 1):
                vacc${M}${N} = vacc${M}${N + (1 << LOG2N)};
              c${M} += ${1 << LOG2N};
        }

      nc = 0;
    }
  } while (nc != 0);
}
