/*
 * Single Precision Floating Point
 *
 * Copyright (C) 2016  Krzysztof Mazur <krzysiek@podlesie.net>
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 */

#ifdef HAVE_CONFIG_H
#include <config.h>
#endif
#include <svl.h>

#ifdef C6X_ASM
#include "svl-c6x.h"
#endif

#ifndef __svl_norm
static float __svl_norm(const float *a, unsigned int n)
{
	unsigned int i;
	float sum = 0;

	for (i = 0; i < n; i++)
		sum += a[i] * a[i];
	return sum;
}
#endif

float svl_norm(const float *a, unsigned int n)
{
	return __svl_norm(a, n);
}

float svl_norm_w(const float *a, unsigned int size, unsigned int i,
		unsigned int n)
{
	float sum;
	unsigned int chunk;

	chunk = size - i;
	if (chunk > n)
		chunk = n;

	sum = svl_norm(a + i, chunk);
	sum += svl_norm(a, n - chunk);
	return sum;
}

#ifndef __svl_prod
static float __svl_prod(const float *a, const float *b, unsigned int n)
{
	unsigned int i;
	float sum = 0;

	for (i = 0; i < n; i++)
		sum += a[i] * b[i];
	return sum;
}
#endif

float svl_prod(const float *a, const float *b, unsigned int n)
{
	return __svl_prod(a, b, n);
}

float svl_prod_w(const float *a, const float *b,
		unsigned int b_size, unsigned int b_i,
		unsigned int n)
{
	float sum;
	unsigned int chunk;

	chunk = b_size - b_i;
	if (chunk > n)
		chunk = n;

	sum = __svl_prod(a, b + b_i, chunk);
	sum += __svl_prod(a + chunk, b, n - chunk);
	return sum;
}

#ifndef __svl_mpy
static void __svl_mpy(const float *a, float alpha,
		const float *restrict b, float beta,
		unsigned int n, float *c)
{
	unsigned int i;

	for (i = 0; i < n; i++)
		c[i] = alpha * a[i] + beta * b[i];
}
#endif

void svl_mpy(const float *a, float alpha,
		const float *restrict b, float beta,
		unsigned int n, float *c)
{
	return __svl_mpy(a, alpha, b, beta, n, c);
}

void svl_mpy_w(const float *a, float alpha,
		const float *b, float beta,
		unsigned int b_size, unsigned int b_i,
		unsigned int n, float *c)
{
	unsigned int chunk;

	chunk = b_size - b_i;
	if (chunk > n)
		chunk = n;

	__svl_mpy(a, alpha, b + b_i, beta, chunk, c);
	__svl_mpy(a + chunk, alpha, b, beta, n - chunk, c + chunk);
}