2D convolution on 28x28 matrix using 5x5 filter. Stride assumed to be 1, padding 0 and input matrix as well as kernel are assumed to be square. Not very fast or usable, but illustrates the core algorithm/what convolution is doing.

conv.c

#include <stdio.h>

/* 
2D convolution as used in machine learning, i.e. actually cross-correlation, not convolution.
But since the kernel would be learned either way, flipping it is not necessary. 
So, let's just call it convolution.
*/

int main() {
    int STRIDE = 1; // Don't let this fool you; stride other than 1 will not work currently
    int H = 28, W = 28, K = 5;
    int RESULT_SIZE = (H-K)/STRIDE+1;
    printf("%dx%d\n", RESULT_SIZE, RESULT_SIZE);

    float input_mat[H][W], kernel[K][K];
    float res[RESULT_SIZE][RESULT_SIZE]; // Result is assumed to be square so implicitly input_mat is as well..

    // Initialize input_mat and kernel with values
    for (int i=0; i<K; i++) {
        for (int j=0; j<K; j++) {
            kernel[i][j] = 2.0;
        }
    }
    for (int i=0; i<H; i++) {
        for (int j=0; j<W; j++) {
            input_mat[i][j] = 5.0;
        }
    }

    // Do regular 2D Convolution using loops
    for (int h=0; h<(H-K+1); h++) {
        for (int w=0; w<(W-K+1); w++) {
            float sum = 0;
            for (int x=0; x<K; x++) {
                for (int y=0; y<K; y++) {
                    sum += input_mat[h+x][w+y] * kernel[x][y];
                }
            }
            res[w][h] = sum;
        }
    }

    // Print result array
    for (int i=0; i<RESULT_SIZE; i++) {
        for (int j=0; j<RESULT_SIZE; j++) {
            printf("%2.f ", res[i][j]);
            if (j == RESULT_SIZE-1) {
                printf("\n");
            }
        }
    }
    return 0;
}