Значения, полученные из расчетов на Numpy и C, не равны.

Вопрос или проблема

import ctypes
from ctypes import CDLL
from ctypes import c_size_t
import numpy as np
import time

# загружаем библиотеку
mylib = CDLL("example.so")

# C-тип, соответствующий numpy array для 3D массивов и 1D массива
ND_POINTER_3 = np.ctypeslib.ndpointer(dtype=np.float32, ndim=3, flags="C_CONTIGUOUS")
ND_POINTER_1 = np.ctypeslib.ndpointer(dtype=np.float32, ndim=1, flags="C_CONTIGUOUS")

# определяем прототипы
mylib.add_vector.argtypes = [ND_POINTER_3, ND_POINTER_3, ND_POINTER_1, ND_POINTER_3, c_size_t]
mylib.add_vector.restype = ctypes.c_void_p

# Инициализируем numpy массивы (90x160x3 для la и lb, 90x160x1 для gm, и 90x160x3 для out)
np.random.seed(0)
la = np.random.rand(90, 160, 3).astype(np.float32)
lb = np.random.rand(90, 160, 3).astype(np.float32)
gm = np.random.rand(90, 160, 1).astype(np.float32)
out = np.zeros((90, 160, 3), dtype=np.float32)

# Уплощаем массив gm для передачи в виде 1D массива
gm_flat = gm.flatten()

s = time.perf_counter()
for i in range(1000):
    mylib.add_vector(la, lb, gm_flat, out, la.size)
out = np.array(out)
e = time.perf_counter()
print(f"Время C : {e-s}")

s = time.perf_counter()
for i in range(1000):
    test_result = la * lb + la * (1.0 - la)
e = time.perf_counter()
print(f"Время Python : {e-s}")

print(out[0, 1, :])
print(test_result[0, 1, :])

#include <stddef.h>

// Функция сложения векторов, поддерживающая 3D массивы и 1D массивы
void add_vector(float *la, float *lb, float *gm, float *out, size_t size) {
    for (size_t i = 0; i < size; ++i) {
        // Выполняем поэлементную операцию: ls = la * gm + lb * (1.0 - gm)
        size_t gm_index = (i / 3);
        out[i] = la[i] * gm[gm_index] + lb[i] * (1.0f - gm[gm_index]);
    }
}

gcc -shared -o example.so -fPIC module_vector.c

Время C : 0.05158950947225094
Время Python : 0.05883005540817976
[0.6418773  0.2731492  0.50642085]
[0.6623829 0.2818951 0.4456129]

Ответ или решение

Разница между результатами, полученными с использованием NumPy и C, может быть вызвана несколькими причинами. Давайте рассмотрим это более подробно и проанализируем, как можно исправить ситуацию.

Причины несовпадения значений:

1. Проблемы с индексированием:
   В вашей реализации C функция add_vector использует одномерный индекс для доступа к элементам трёхмерного массива la и lb. В частности, переменная gm_index рассчитывается как i / 3, что подразумевает, что вы ожидаете, что размер первого измерения равен 3. Это может привести к путанице с индексами, если ваши данные имеют другой размер.

   Ваша функция должна учитывать, что каждый элемент gm используется для расчета нескольких элементов. Даже если в gm только одно значение на плавающей запятой, вам нужно убедиться, что вы правильно обращаетесь к каждому элементу la и lb.

2. Порядок вычислений:
   Результаты, полученные с помощью NumPy, могут слегка отличаться от результатов C из-за порядка операций с плавающей запятой. Это происходит из-за особенностей арифметических операций с числами с плавающей запятой, которые могут давать разные результаты в зависимости от порядка вычислений.

3. Тип данных и округление:
   Убедитесь, что используется один и тот же тип данных (в данном случае float32), как в C, так и в Python. Кроме того, округление ошибок может привести к различиям в конечных значениях.

Решение проблемы:

1. Исправьте индексирование в функции C. Вам нужно передавать правильный индекс в gm для каждого элемента:

   void add_vector(float *la, float *lb, float *gm, float *out, size_t size) {
      for (size_t i = 0; i < size; ++i) {
          size_t gm_index = (i / (size / 3)); // Измените это согласно размерности gm
          out[i] = la[i] * gm[gm_index] + lb[i] * (1.0f - gm[gm_index]);
      }
   }

2. Проверьте реализацию и поправьте порядок вычислений: Убедитесь, что как в C, так и в Python относительно порядка вы вычисляете значения. Вычисления могут давать разные результаты при отсутствии строгой последовательности.

3. Тестирование и отладка: Проверьте промежуточные результаты в обеих реализациях, чтобы понять, на каком этапе происходит расхождение. Добавьте дополнительные выводы для диагностики.

Заключение:

Проблема с различием в значениях, полученных из C и NumPy, часто сводится к вопросам индексирования и порядка операций с плавающей запятой. Настоятельно рекомендуется внимательно просмотреть код и протестировать его на небольших наборах данных, чтобы убедиться, что все выполняется правильно до перехода к более крупным массивам.

С такими изменениями вы сможете минимизировать разницу между результатами.