Создание дважды стохастической матрицы

[Alexashek]

Думаю все знают что это такое. Это матрица, сумма значений элементов по строкам и столбцам которой равна 1.
Вот я решил написать некий генератор таких матриц, вот пока без идей
Подкиньте кто нибудь мыслишек, буду очень благодарен.

[Mandrivnyk]

Хм...
Ну, как начальный вариант, разве что. Причем для квадратной матрицы. Кстати, я не совсем уверен, может ли дважды стохастическая матрица быть _не_квадратной...

Так вот.
В начале "засеваем" всю матрицу одинаковыми значениями 1/n (n - размерность матрицы).
Затем выбираем случайный элемент Pij и изменяем его на некоторую случайную величину, например, (0 < Δ < 1/2n). Соответственно, меняем один случайный элемент из этого же столбца и один случайный элемент из этой же строки. То есть, если к элементу Pij мы прибавили Δ, то из случайных элементов в строке и столбце мы эту Δ вычитаем. И наоборот.
Затем делаешь коррекцию для элемента, находящегося на пересечении строки и столбца этих двух случайно выбранных элементов.
То есть, по идее, ты меняешь четыре элемента, образующие прямоугольник внутри матрицы. При этом, если я не ошибаюсь, матрица должна остаться дважды стохастической
В общем -- муторное дело получается -)
Этот пункт повторить некоторое, сравнительно большое количество раз.

Можно просто пройтись циклом по всем элементам по очереди, изменяя их на случайную Δ и коррегировать матрицу указанным способом. Так, наверное, будет лучше.

Или вообще -- совместить оба способа -) Сначала циклом, потом случайно.

А вообще, задача интересная -)
Если будут более рациональные идеи -- отпишись.

[Alexashek]

Делюсь мыслями:
Сначала я решил попробовать генерировать числа, в сумме дающие скажем 100 и делать матрицу из них(потом просто нормировать их), но этот вариант провалился, т.к. задача свелась к проблеме перечисления латинских квадратов, которая пока что насколько я знаю является нерешенной.

С дважды стохастичностью все оказалось проще. Известны следующие факты:
-дважды стохастическая матрица квадратная;
-произведение дважды стохастических матриц есть дважды стохастическая матрица;
-справедлива теорема Биркгофа-фон Неймана(точнее следствие из нее): матрица М является дважды стохастической тогда и только тогда когда она представляется в виде суммы матриц перестановок умноженных на коэффициенты c, в сумме дающие единицу.
Отсюда следует следующий алгоритм:
1)сгенерировать случайные подстановки(это не сложно, если надо могу выложить код функции, делающей это);
2)поставить в соответствие сгенерированным подстановком матрицы(матрица подстановок-матрица в каждой строке и столбце которой стоит ровно одна единица, а остальные нули, причем единицы соответствуют переходам в подстановке);
3)сгенерировать коэффициенты, т.е. числа, в сумме дающие 1, и каждое умножить на одну из матриц;
4)просуммировать полученные матрицы и получим искомую дважды стохастическую.

Алгоритм собственно такой, может это поможет, если у кого нибудь возникнет потребность в чем то подобном. К сожалению сейчас пока нет времени написать программу и выложить код, но возможно через пару недель сяду напишу.

[Mandrivnyk]

Реализовал свой алгоритм.
Вроде, работает -)

#include <iostream>
#include <iomanip>
#include <ctime>
#include <cstdlib>

using namespace std;

void showArray(float **Array, int dim);
void randDelta(float **Array, int rind, int cind, int dim);
bool falseDelta(float element, float delta);

int main()
{
    int dim;
    cout << "Введите размерность матрицы:\n";
    cin >> dim;
    float **Array = new float *[dim];
    for (int i = 0; i < dim; i++)
        Array[i] = new float[dim];
    
    for (int i = 0; i < dim; i++)
        for (int j = 0; j < dim; j++)
            Array[i][j] = (float)(1.0 / dim);
    
    srand(time(NULL));

    for (int i = 0; i < dim; i++)
        for (int j = 0; j < dim; j++)
            randDelta(Array, i, j, dim);

    // Полученная матрица
    showArray(Array, dim);

    // Проверка по строкам...
    float check;
    for (int i = 0; i < dim; i++)
    {
        check = 0.0;
        for (int j = 0; j < dim; j++)
            check += Array[i][j];
        cout << check << endl;
    }
    cout << endl;

    // ... и по столбцам
    for (int i = 0; i < dim; i++)
    {
        check = 0.0;
        for (int j = 0; j < dim; j++)
            check += Array[j][i];
        cout << check << endl;
    }

    // Освобождение памяти
    for(int i = 0; i < dim; i++)
        delete[] Array[i];
    delete[] Array;

    return 0;
 }

void randDelta(float **Array, int rind, int cind, int dim)
{
    int randrow, randcolumn;
    float delta;

    // Берем случайный элемент матрицы, отличный от переданного в функцию
    do
    {
        randrow = rand() % dim;
        randcolumn = rand() % dim;
    } while (randrow == rind && randcolumn == cind);

    do
        delta = (float)(rand() % 100) / 277.0 - 50.0 / 277.0;   // Числа взяты "от фонаря"
                                                                // тут получается случайное число от -50/277 до 50/277

    // Проверка выхода за допустимые границы (0..1)
    while (falseDelta(Array[rind][cind], delta)
            || falseDelta(Array[rind][randcolumn], delta)
            || falseDelta(Array[randrow][cind], delta)
            || falseDelta(Array[randrow][randcolumn], delta));

        Array[rind][cind] += delta;
        Array[rind][randcolumn] -= delta;
        Array[randrow][cind] -= delta;
        Array[randrow][randcolumn] += delta;
}

void showArray(float **Array, int dim)
{
    cout << fixed;
    for (int i = 0; i < dim; i++)
    {
        for (int j = 0; j < dim; j++)
            cout << setprecision(4) << Array[i][j] << "\t";
        cout << endl;
    }
    cout << endl;
}

bool falseDelta(float element, float delta)
{
    return element - delta < 0 || element + delta > 1 || element + delta < 0 || element - delta > 1;
}