# Comenzamos importando las librerías a utilizar:

# Para crear gráficos y visualizaciones:
import matplotlib.pyplot as plt
# Para operaciones numéricas eficientes, especialmente con arrays:
import numpy as np
# Para manejar y analizar datos en estructuras tipo DataFrame:
import pandas as pd
# Para medir tiempos de ejecución o hacer pausas:
import time
# Para interactuar con el sistema (manejo de rutas y argumentos):
import sys

# Se define la función 'mydot' para calcular el producto punto entre dos
# vectores:

# Calcula el producto punto entre dos vectores 1D.
# Parámetros:
# u : array-like de dimensión (d,)
# x : array-like de dimensión (d,)
# Retorna s : El valor del producto punto u · x

def mydot(u, x):
  s = 0.0
  for ui, xi in zip(u, x):
    s += ui * xi
  return s

# Se define la función 'getmaxdot' para encontrar el vector columna en una matriz
# que produce el máximo producto punto con un vector dado

# Se define la función 'getmaxdot' donde:
# u: es el vector columna (d,)
# X: es la matriz de tamaño (d, n)
# Regresa (índice, valor) del máximo producto punto

def getmaxdot(u,X):
  maxpos = 0
  maxdot = mydot(u, X[:, 0])

  _, ncols = X.shape
  for i in range(1, ncols):
    d = mydot(u, X[:, i])
    if d > maxdot:
      maxpos = i
      maxdot = d
  return maxpos, maxdot

# Ejemplo de prueba de las funciones definidas

# Creamos un vector 'u' de tamaño 4 con valores aleatorios entre 0 y 1
u = np.random.rand(4)

# Creamos una matriz 'X' de tamaño 4x1000 con valores aleatorios entre 0 y 1
# Cada columna representa un vector con el que compararemos 'u'
X = np.random.rand(4, 1000)

# Usamos la función 'getmaxdot' para encontrar la columna de X que tiene
# el producto punto máximo con el vector u
idx, val = getmaxdot(u, X)

# Mostramos el índice de la columna con el producto punto máximo y su valor
print(f"Mejor índice: {idx}, Valor del producto punto: {val}")

Mejor índice: 991, Valor del producto punto: 1.6486541749755446

# Se genera una matriz (d,n) con vectores columna en la esfera unitaria

# Parámetros:
# d : Dimensión de cada vector
# n : Número de vectores columna
# dtype : Tipo de datos de la matriz generada
# Retorna:
# M : Matriz con columnas normalizadas a longitud 1

def generate_unit_sphere_matrix(d, n, dtype=np.float32):
  # Matriz normal (0, 1)
  M = np.random.randn(d, n).astype(dtype)
  # Normalizamos las columnas
  M /= np.linalg.norm(M, axis=0, keepdims=True)
  return M

# Generación del algoritmo en versión vectorizada con producto matricial

# Encuentra, para cada vector columna en Q, la columna de X que tiene
# el producto punto máximo, utilizando operaciones vectorizadas para eficiencia
# Parámetros: Q - Matriz cuyas columnas son los vectores de consulta.
# X - Matriz cuyas columnas se comparan con cada columna de Q.
# Retorna: max_indices : indices de las columnas de X que producen el producto
# punto máximo para cada vector columna de Q.

def top_maxdot(Q, X):
  # Q.T @ X produce una matriz de tamaño (m, n) donde cada elemento (i,j)
  # es el producto punto entre la columna i de Q y la columna j de X
  dot_products = Q.T @ X
  max_indices = np.argmax(dot_products, axis=1)
  return max_indices

# Generación del algoritmo en versión 'getmaxdot'

# Se retoma la versión vectorizada para encontrar el índice de la columna con
# producto punto máximo
def top1_maxdot(Q, X):
    # Algoritmo vectorizado con producto matricial:
    dot_products = Q.T @ X
    max_indices = np.argmax(dot_products, axis=1)
    return max_indices

# Establecemos una versión por bloques para manejar matrices muy grandes y
# ahorrar memoria

# Parámetros: Q - vectores de consulta. X - Matriz a comparar
# block_size : Tamaño de cada bloque de columnas de X para procesamiento
# Retorna: best_idx - índices de las columnas de X con producto punto máximo
def top1_maxdot_blocked(Q, X, block_size=5000):
    m = Q.shape[1]
    n = X.shape[1]
    best_idx = np.zeros(m, dtype=int)
    best_val = np.full(m, -np.inf, dtype=np.float32)

    # Procesamos X en bloques
    for start in range(0, n, block_size):
      end = min(start + block_size, n)
      dots = Q.T @ X[:, start:end]
      local_idx = np.argmax(dots, axis=1)
      local_val = np.max(dots, axis=1)

      # Actualizamos solo donde encontramos un valor mayor al anterior
      mask = local_val > best_val
      best_val[mask] = local_val[mask]
      best_idx[mask] = start + local_idx[mask]

    return best_idx

# Se define la función que calcula el producto punto entre dos vectores 1D
def mydot(u, x):
  s = 0.0
  for ui, xi in zip(u, x):
    s += ui * xi
  return s

# Se calcula el índice de la columna de X con producto punto máximo con u
# usando bucles
def getmaxdot(u, X):
  maxpos = 0
  maxdot = mydot(u, X[:, 0])
  _, ncols = X.shape
  for i in range(1, ncols):
    d = mydot(u, X[:, i])
    if d > maxdot:
      maxpos = i
      maxdot = d
  return maxpos

# Se calculan los índices de columna con producto punto máximo para cada vector
# columna de Q usando la getmaxdot
# Retorna: ndarray de índices máximos por columna de Q
def getmaxdot_for_Q(Q, X):
  results = []
  for j in range(Q.shape[1]):
    results.append(getmaxdot(Q[:, j], X))
  return np.array(results)

# Función que compara el rendimiento de diferentes métodos para encontrar la
# columna de X que produce el producto punto máximo para cada columna de Q
# Parámetros: Q - Matriz de vectores de consulta, X - Matriz de vectores a
# comparar .
# Método a usar: "vectorized" o "getmaxdot"

def comparison(Q, X, method="vectorized", block_threshold=2e8):
  d, m = Q.shape
  _, n = X.shape

  if method == "vectorized":
    # Si el tamaño del problema supera el umbral, usamos la versión por bloques:
    if m * n > block_threshold:
      start = time.time()
      result = top1_maxdot_blocked(Q, X)
      elapsed = time.time() - start
    else:
      start = time.time()
      result = top1_maxdot(Q,X)
      elapsed = time.time() - start
    ops = 2 * m * n * d
  else:
      # Método'getmaxdot' usando bucles
      start = time.time()
      result = getmaxdot_for_Q(Q, X)
      elapsed = time.time() - start
      ops = 2 * m * n * d

  # Estimación de memoria utilizada por las matrices Q y X (en bytes)
  memory = Q.nbytes + X.nbytes

  return {
      "method": method,
      "m": m,
      "n": n,
      "ops": ops,
      "time_sec": elapsed,
      "memory_MB": memory / (1024**2)
  }

# Experimentos para comparar rendimiento de los métodos

# Dimensión de los vectores:
d = 8
# Diferentes tamaños de conjuntos de consulta Q
m_values = [10**3, 10**4, 10**5]
# Diferentes tamaños de conjuntos de vectores X
n_values = [10**3, 10**4, 10**5]

results = []

for m in m_values:
  for n in n_values:
    print(f"Caso m={m}, n={n}")
    Q = generate_unit_sphere_matrix(d, m)
    X = generate_unit_sphere_matrix(d, n)

    # Método vectorizado
    res_vec = comparison(Q, X, method="vectorized")
    results.append(res_vec)
    print("Método vectorizado:", res_vec)

    # Método 'getmaxdot' (se establece el límite a 10^4 para hacer eficiente la
    # comparación)
    if m <= 1e4 and n <= 1e4:
      res_get = comparison(Q, X, method="getmaxdot")
      results.append(res_get)
      print("Método getmaxdot:", res_get)

Caso m=1000, n=1000
Método vectorizado: {'method': 'vectorized', 'm': 1000, 'n': 1000, 'ops': 16000000, 'time_sec': 0.006758451461791992, 'memory_MB': 0.06103515625}
Método getmaxdot: {'method': 'getmaxdot', 'm': 1000, 'n': 1000, 'ops': 16000000, 'time_sec': 5.321566820144653, 'memory_MB': 0.06103515625}
Caso m=1000, n=10000
Método vectorizado: {'method': 'vectorized', 'm': 1000, 'n': 10000, 'ops': 160000000, 'time_sec': 0.04374432563781738, 'memory_MB': 0.335693359375}
Método getmaxdot: {'method': 'getmaxdot', 'm': 1000, 'n': 10000, 'ops': 160000000, 'time_sec': 42.50253486633301, 'memory_MB': 0.335693359375}
Caso m=1000, n=100000
Método vectorizado: {'method': 'vectorized', 'm': 1000, 'n': 100000, 'ops': 1600000000, 'time_sec': 0.3472719192504883, 'memory_MB': 3.082275390625}
Caso m=10000, n=1000
Método vectorizado: {'method': 'vectorized', 'm': 10000, 'n': 1000, 'ops': 160000000, 'time_sec': 0.029251575469970703, 'memory_MB': 0.335693359375}
Método getmaxdot: {'method': 'getmaxdot', 'm': 10000, 'n': 1000, 'ops': 160000000, 'time_sec': 41.9398980140686, 'memory_MB': 0.335693359375}
Caso m=10000, n=10000
Método vectorizado: {'method': 'vectorized', 'm': 10000, 'n': 10000, 'ops': 1600000000, 'time_sec': 0.24840974807739258, 'memory_MB': 0.6103515625}
Método getmaxdot: {'method': 'getmaxdot', 'm': 10000, 'n': 10000, 'ops': 1600000000, 'time_sec': 416.4901270866394, 'memory_MB': 0.6103515625}
Caso m=10000, n=100000
Método vectorizado: {'method': 'vectorized', 'm': 10000, 'n': 100000, 'ops': 16000000000, 'time_sec': 3.0887436866760254, 'memory_MB': 3.35693359375}
Caso m=100000, n=1000
Método vectorizado: {'method': 'vectorized', 'm': 100000, 'n': 1000, 'ops': 1600000000, 'time_sec': 0.2484276294708252, 'memory_MB': 3.082275390625}
Caso m=100000, n=10000
Método vectorizado: {'method': 'vectorized', 'm': 100000, 'n': 10000, 'ops': 16000000000, 'time_sec': 4.3575780391693115, 'memory_MB': 3.35693359375}
Caso m=100000, n=100000
Método vectorizado: {'method': 'vectorized', 'm': 100000, 'n': 100000, 'ops': 160000000000, 'time_sec': 29.75477123260498, 'memory_MB': 6.103515625}

# Guardamos los resultados en un DataFrame
df = pd.DataFrame(results)
df["problem_size"] = df["m"] * df["n"]
df["ops_per_sec"] = df["ops"] / df["time_sec"]

print("Resultados Finales:")
print(df)

Resultados Finales:
        method       m       n           ops    time_sec  memory_MB  \
0   vectorized    1000    1000      16000000    0.006758   0.061035   
1    getmaxdot    1000    1000      16000000    5.321567   0.061035   
2   vectorized    1000   10000     160000000    0.043744   0.335693   
3    getmaxdot    1000   10000     160000000   42.502535   0.335693   
4   vectorized    1000  100000    1600000000    0.347272   3.082275   
5   vectorized   10000    1000     160000000    0.029252   0.335693   
6    getmaxdot   10000    1000     160000000   41.939898   0.335693   
7   vectorized   10000   10000    1600000000    0.248410   0.610352   
8    getmaxdot   10000   10000    1600000000  416.490127   0.610352   
9   vectorized   10000  100000   16000000000    3.088744   3.356934   
10  vectorized  100000    1000    1600000000    0.248428   3.082275   
11  vectorized  100000   10000   16000000000    4.357578   3.356934   
12  vectorized  100000  100000  160000000000   29.754771   6.103516   

    problem_size   ops_per_sec  
0        1000000  2.367406e+09  
1        1000000  3.006633e+06  
2       10000000  3.657617e+09  
3       10000000  3.764481e+06  
4      100000000  4.607341e+09  
5       10000000  5.469791e+09  
6       10000000  3.814983e+06  
7      100000000  6.440971e+09  
8      100000000  3.841628e+06  
9     1000000000  5.180100e+09  
10     100000000  6.440507e+09  
11    1000000000  3.671764e+09  
12   10000000000  5.377289e+09

# Tiempo
plt.figure(figsize=(8,5))
for method in df["method"].unique():
    subset = df[df["method"] == method]
    plt.plot(subset["problem_size"], subset["time_sec"], marker="o",
             label=method)

plt.xscale("log")
plt.yscale("log")
plt.xlabel("Tamaño del problema (m × n)")
plt.ylabel("Tiempo (s)")
plt.title("Comparación de tiempos: vectorizado vs getmaxdot", fontweight='bold')
plt.legend()
plt.grid(True, which="both", ls="--", alpha=0.6)
plt.show()

# Memoria
plt.figure(figsize=(8,5))
for method in df["method"].unique():
    subset = df[df["method"] == method]
    plt.plot(subset["problem_size"], subset["memory_MB"], marker="o", label=method)

plt.xscale("log")
plt.xlabel("Tamaño del problema (m × n)")
plt.ylabel("Memoria (MB)")
plt.title("Uso de memoria: vectorizado vs getmaxdot", fontweight='bold')
plt.legend()
plt.grid(True, which="both", ls="--", alpha=0.6)
plt.show()

# Eficiencia (ops/sec)
plt.figure(figsize=(8,5))
for method in df["method"].unique():
    subset = df[df["method"] == method]
    plt.plot(subset["problem_size"], subset["ops_per_sec"], marker="o", label=method)

plt.xscale("log")
plt.yscale("log")
plt.xlabel("Tamaño del problema (m × n)")
plt.ylabel("Operaciones por segundo")
plt.title("Eficiencia: vectorizado vs getmaxdot", fontweight='bold')
plt.legend()
plt.grid(True, which="both", ls="--", alpha=0.6)
plt.show()