Creating a Stochastic matrix using dynamic C style arrays on C++

Question

So im trying to create an function that creates a Stochastic matrix (a square matrix where each column summing is 1), it tried dividing every element of the column to the summming of it ,but sometimes it seems to throw unexpected values. The values of the matrix need to have 2 decimal places so i cannot change that. What am i doing wrong?

Some of the tests i did: enter image description here enter image description here In the second one, the summing of some of the columns isn't 1.

#include <iostream>
#include <cmath>   //fmod
#include <math.h> //round()

using namespace std;

template<typename T>
T** crear_matriz(int filas, int columnas) {
    T** matriz = new T * [filas];
    for (size_t i = 0; i < filas; i++) {
        matriz[i] = new T[columnas];
    }
    return matriz;
}

template<typename T>
void borra_matriz(T** matrix, int filas) {
    for (size_t i = 0; i < filas; i++) {
        delete[] matrix[i];
    }
    delete[] matrix;
}

double** CrearMatrizEstocastica(double ** matrix,int rows, int column) {
    for (int j = 0; j < column; j++) {
        double columnsum = 0;
        for (int i = 0; i < rows; i++) {
            double random_num = static_cast<double>(rand()) / RAND_MAX;
            matrix[i][j] = random_num;
            columnsum += matrix[i][j];
        }
        for (int i = 0; i < rows; i++) {
            matrix[i][j] /= columnsum;
            matrix[i][j] = round(matrix[i][j] * 10) / 10;         
        }
        
    }
    cout << "Schocastic Matrix created!" << endl;
    return matrix;
}
int main()
{
    srand(time(NULL));
    double** MatrizEstocastica=crear_matriz<double>(4,4) ;
MatrizEstocastica = CrearMatrizEstocastica(MatrizEstocastica,4,4);
MostrarMatrizEstocastica(MatrizEstocastica, 4, 4);
borra_matriz<double>(MatrizEstocastica , 4);

    return 0;
}

https://stackoverflow.com/questions/25385173/what-is-a-debugger-and-how-can-it-help-me-diagnose-problems?r=Saves_AllUserSaves — πάντα ῥεῖ, Apr 22 '23 at 08:05
Do not post plain text, code or error messages as images here please! — πάντα ῥεῖ, Apr 22 '23 at 08:51
Does this answer your question? [Is floating point math broken?](https://stackoverflow.com/questions/588004/is-floating-point-math-broken) — jabaa, Apr 22 '23 at 11:31

score 0 · Answer 1 · answered Apr 22 '23 at 10:05

Of course, the sum isn't 1. Example with 3 numbers : { 0.33, 0.33, 0.34 } sum is about 1. After rounding it becomes { 0.3, 0.3, 0.3 } sum is about 0.9.

You can propagate the errors to next value, so do not lose decimals

    for (int j = 0; j < column; j++) {
        double columnsum = 0;
        for (int i = 0; i < rows; i++) {
            double  random_num = rand();
            matrix[i][j] = random_num;
            columnsum += matrix[i][j];
        }
        double  propagate = 0.;
        for (int i = 0; i < rows; i++) {
            double  new_val = matrix[i][j] / columnsum + propagate;
            double  rounded_val = round(new_val * 10) / 10;
            matrix[i][j] = rounded_val;
            propagate = new_val - rounded_val;
        }
    }

Note that the sum may still not be 1. Because floating numbers are stored internally in binary and decimal numbers may not be exact values (here 0, 0.5 and 1 are exact, but 0.1, 0.2, 0.3, 0.4, 0.6, 0.7, 0.8, 0.9 are approximations).

Creating a Stochastic matrix using dynamic C style arrays on C++

1 Answers1