0

I'd like to write functions that do common matrix operations. This can be done by 2-dim arrays or pointer arithmetic. I prefer a pointer version. Now with pointers I'd write a function like this:

void matmult(double *a, double *b, double *c, int m, int n, int k);

The problem is that I have to use a cast when I pass 2-dim arrays to the function. Is there a good solution to avoid this problem?

Works without cast (of course), but I want to avoid compiler warnings.

Update: the arrays are defined as 2-dim array and the calling function looks like this:

// M, N, K are constants
double a[M][N];
double b[N][K];
double c[M][K];
matmult((double *)a, (double *)b, (double *)c, M, N, K);

The function matmult is a straight forward implementation of matrix multiplication (three nested for loops using pointers)

*(c + i*k + j) += *(a + i*n + p) * *(b + p*k + j);

I just want to get rid of the cast.

bernhard67
  • 3
  • 2
  • well it depends on what your 2d array looks like, show calling code and show the code in matmult – pm100 Nov 20 '22 at 22:22
  • Why don't you store the matrix dimensions in the matrix itself as well as an array of numbers; it's what math does. – Neil Nov 21 '22 at 06:55
  • @Neil What exactly would you do? Note, that the matrix dimensions are not the problem as the caller must ensure them to have proper dimensions. The problem is the type mismatch between double * and (*doube)[]. – bernhard67 Nov 22 '22 at 07:09
  • https://stackoverflow.com/q/42094465/16835308 – Zakk Nov 22 '22 at 10:39

2 Answers2

0
static void matmult(double a[][N], double (*b)[K], double c[M][K], int m, int n, int k);

This prototype works for passing 2-dimensional arrays. See this C FAQ entry or Array to pointer decay and passing multidimensional arrays to functions. The M in c[M][K] is totally superfluous in C, but may be useful as self-documentation, (or may confuse the reader more.)

However, this is not very encapsulated and I would hesitate to programme a general matrix algorithm with this. The size is an integral part of the matrix itself. I might use C99's flexible array member to store the size and the data together, as long as it's not part of other structures.

#include <stdlib.h>
#include <stdio.h>
#include <errno.h>
#include <stdarg.h>
#include <limits.h>

struct matrix {
    unsigned x, y;
    double data[]; /* flexible array member */
};

/** Constructor; data is uninitialized. */
static struct matrix *matrix(const unsigned x, const unsigned y) {
    if(!x || !y || x >= UINT_MAX / y) { errno = ERANGE; return 0; }
    struct matrix *const m = malloc(offsetof(struct matrix, data)
        + sizeof *m->data * x * y);
    if(!m) return 0;
    m->x = x, m->y = y;
    return m;
}

/** Constructor; IBM has a useful extension that allows stack construction. */
static struct matrix *matrix_init(const unsigned x, const unsigned y, ...) {
    struct matrix *const m = matrix(x, y);
    if(!m) return 0;
    va_list argp;
    va_start(argp, y);
    for(unsigned i = 0, size = x * y; i < size; i++)
        m->data[i] = va_arg(argp, double);
    va_end(argp);
    return m;
}

static void matrix_print(const struct matrix *const m) {
    if(!m) { printf("null\n"); return; }
    for(unsigned y = 0; y < m->y; y++) {
        printf("[");
        for(unsigned x = 0; x < m->x; x++)
            printf("%s%4.f", x ? ", " : "", m->data[y * m->x + x]);
        printf("]\n");
    }
}

static struct matrix *matmult(const struct matrix *a,
    const struct matrix *b) {
    if(!a || !b || a->y != b->x) { errno = EDOM; return 0; }
    struct matrix *const c = matrix(b->x, a->y);
    if(!c) return 0;
    for(unsigned y = 0; y < a->y; y++)
        for(unsigned x = 0; x < b->x; x++)
            c->data[y * b->x + x] = 0.0;
    /* implement:
     *(c + i*k + j) += *(a + i*n + p) * *(b + p*k + j); */
    return c;
}

int main(void) {
    struct matrix *a = 0, *b = 0, *c = 0;
    int success = EXIT_SUCCESS;
    if(!(a = matrix_init(2, 2, 1.0, 2.0, 3.0, 4.0))
        || !(b = matrix_init(2, 2, 3.0, 20.0, 1.0, 0.0))) goto catch;
    matrix_print(a), printf("*\n"), matrix_print(b);
    if(!(c = matmult(a, b))) goto catch;
    printf("=\n"), matrix_print(c);
    goto finally;
catch:
    success = EXIT_FAILURE;
    perror("matrix");
finally:
    free(a), free(b), free(c);
    return success;
}
Neil
  • 1,767
  • 2
  • 16
  • 22
0

This is how I finally resolved my problem thanks to the hint of @Neil:

void matrixmult(double (*a)[], double (*b)[], double (*c)[], int m, int n, int k)
{
   for (...)
       for (..)
          for (...)
              *((double *)c + i*k + j) = ...
}

Now, the function can be called without cast and I am avoiding the restrictions of the C99 VLA feature (dimensions must be passed first).

bernhard67
  • 3
  • 2