/*
 * matutil.c
 */

#include <stdio.h>
#include "matutil.h"

vector new_vector(int n)
{
    return (vector)malloc(sizeof(scalar) * n);
}

void free_vector(vector v)
{
    free(v);
}

matrix new_matrix(int nrow, int ncol)
{
    int i;
    vector ap;
    matrix a;

    if ((a = (matrix)malloc(nrow * sizeof(vector))) == NULL)
        return NULL;
    if ((ap = (vector)malloc(nrow * ncol * sizeof(scalar))) == NULL) {
        free(a);
        return NULL;
    }
    for (i = 0; i < nrow; i++)
        a[i] = ap + (i * ncol);
    return a;
}

void free_matrix(matrix a)
{
    free(a[0]);
    free(a);
}

double dotprod(int n, vector u, vector v)
{
    int i, n4;
    double s;

    s = 0;
    n4 = n & 3; /* n % 4 */
    for (i = 0; i < n4; i++)
        s += u[i] * v[i];
    for (i = n4; i < n; i += 4)
        s += u[i] * v[i] + u[i + 1] * v[i + 1]
            + u[i + 2] * v[i + 2] + u[i + 3] * v[i + 3];
    return s;
}