Альтернативный способ заполнения «спиральной матрицы»

В процессе изучения основ алгоритмизации и программирования в качестве студента еще в середине 2000х мне попалась довольно известная всем задача по заполнению «спиральной» матрицы. Суть состоит в том, чтобы начиная с позиции [1, 1], продвигаясь по часовой стрелке, заполнить квадратную матрицу заданной величины числами в возрастающем порядке. На ее решение было потрачено около двух часов.

Собственно алгоритм заполнения был тривиален, циклы, в общей сложности состоящие из N² итерации, предполагали прохождение по всем элементам массива в требуемом порядке, при этом увеличивая значение итератора на 1 и заполняя им текущий элемент матрицы. Маршрут начинался с элемента [1, 1], далее продвигается по горизонтали до правого верхнего элемента [1, N], после – вниз до нижнего правого угла [N, N], затем до левого нижнего угла [N, 1] и завершал первый круг на столбец ниже отправной точки [2, 1]. В дальнейшем, такое же круговое движение происходило уже в следующем внутреннем круге, и так далее до центра матрицы.

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

Через некоторое время после успешного решения задачи, на волне переоценки собственных способностей, я задумался: а возможно ли разработать универсальную формулу, позволяющей на основании размера матрицы и координат ячейки вычислить ее значение согласно условию задачи – то есть в конечном итоге заполнить массив, перебирая элементы «традиционным» путем двумя вложенными циклами от 1 до N без использования счетчика.

Но, как вы могли догадаться, ничего путного у меня не вышло, и данная затея была заброшена в пользу других, более насущных проблем в изучении программирования.

Спустя эти самые 18 лет, перебирая наиболее интересные задачи и пути их решения для обучения уже следующего поколения представителей неординарной профессии «программист», меня заинтересовала одна статья на ресурсе «Хабр» (https://habr.com/ru/post/261773), описывающая процесс создания формулы для вычисления количества дней в заданном месяце без использования каких-либо условных операторов, циклов и заранее подготовленных данных.

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

Именно после этого я вспомнил ту самую спиральную матрицу и решился на очередную попытку выведения универсальной формулы.

(Следует отметить, что на одном из сайтов я все же обнаружил достаточно короткое решение в виде функции, которое, однако, содержало два условных перехода).

Итак, задача: разработать алгоритм для вычисления значения элемента спиральной матрицы на основании его координат [i, j] и размера самой матрицы, пользуясь только простыми арифметическими вычислениями. Исключение составляет модуль - абсолютное значение числа.

Условие: никаких условных (прошу прощения за тавтологию) переходов или заготовленных данных в виде массивов, словарей и т.д.

Практическая полезность: никакая. Традиционным способом данная задача решается гораздо быстрее, при этом мозг программиста в процессе разработки функции не травмируется.

Очевидно, такую непростую задачу необходимо разделить на несколько этапов, или, если угодно, выделить основные элементы.

1)       Найти закономерность между «координатами» ячейки и ее порядковым номером в первом внешнем «кольце», и вывести соответствующую форму.

2)       Найти связь между, опять же, координатами ячейки и номером кольца, в котором она находится. Составить формулу.

3)       Связать между собой два алгоритма для выведения общей формулы вычисления значения элемента массива.

Входные данные: N – размер квадратной матрицы (массива).

1 этап. Заполнение внешнего кольца. Для начала попытаемся вывести формулу вычисления значения ячейки для внешнего кольца массива, в котором отсчет начинается с ячейки [1, 1] и двигается по часовой стрелке поворачивая на углах [1, N], [N, N]. Заканчивается на строку ниже начальной точки, т.е. [2, 1].

Математический аппарат будет разрабатываться параллельно с написанием кода на языке программирования, в качестве которого выбираем C++, как наиболее распространенный и удобный язык программирования (здесь конечно читатели могут поспорить, но классика есть классика). Но, в принципе, выведенная формула не должна иметь привязки языке.

Для наглядного изучения процесса напишем соответствующий скрипт, объявив в нем массив размером 5x5 (назовем его “a”), элементы которого будем перебирать традиционным способом двумя вложенными циклами от 1 до N. На данном этапе будем работать только с внешним кольцом.

#include <iostream>
using namespace std;
int main()
{
    int N = 5;          // задаем размер матрицы
    int a[N][N];        // и инициализируем ее
    for (int ik = 0; ik < N; ik++)
        for (int jk = 0; jk < N; jk++)
            a[ik][jk] = 0;          // заполнив для удобства нулями
                                      
    for (int ik = 0; ik < N; ik++){     //назовем его "Основной цикл"
        for (int jk = 0; jk < N; jk++){
            if (!(ik == 0 || ik == N - 1 || jk == 0 || jk == N - 1)) 
                continue;      // Временное условие для фильтрации элементов внесшего "кольца"
            int i = ik + 1;     // Номера строк и столбцов приводим в удобный
            int j = jk + 1;     // в математическом плане вид (от 1 до N)  
            	//  ... здесь будем вставлять основной код вычислений
        }   
    }     

   for (int ik = 0; ik < N; ik++){          //Блок "Вывод массива"
        for (int jk = 0; jk < N; jk++){
           printf( "%02d ", a[ik][jk]);	// дополняем число ведущими нулями
        }
        cout << endl;
    }  
    return 0;
}

Очевидно, что, по крайней мере, до правого нижнего угла внешнего кольца суммарное значение координат (i + j) увеличивается ровно на 1 с каждым шагом. Однако первый элемент в таком случае равняется 2, E_1,2 = 3 и т.д. Поэтому необходимо уменьшить значение суммы (i + j) на один. В результате введем переменную Xs = i + j - 1, которая часто будет использоваться в дальнейших вычислениях. Пишем код:

…
int Xs = (i + j – 1);
a[ik][jk] = Xs;
…

В результате запуска первого скрипта получаем массив:

Как видим заполнение от начала до правого нижнего угла произведено корректно. Однако далее у нас происходит уменьшение шага, вместо увеличения.

Очевидно a[ik][jk] = Xsнуждается в дополнении, при котором все его значения до [5, 5] останутся неизменными, но после этой точки начнут увеличиваться на 1.

Но для начала постараемся привести в норму вторую часть кольца, которая заполнялась бы с ячейки (i = 5, j = 4) начиная со значения 10. В данном случае это легко сделать, лишь вычитая от общего количества элементов первого кольца увеличенного на два (равняется периметру первого кольца N * 4 - 4 = 16 плюс 2) значение Xs.

То есть a_1,2 = 4N – 4 + 2 – Xs = 4N – Xs - 2.

…
int Xs = i + j - 1;     
a[ik][jk] = 4 * N - Xs - 2;
…

Однако теперь у нас правильно заполнилась только левая и нижняя стороны кольца, а противоположенные элементы – нет.

На текущем этапе у нас имеются две формулы, пригодные только для определенных участков массива.

1.    a_i,j = Xs = i + j - 1; действует от [1, 1] до [N, N].

2.    a_i,j = 4*N – 2 - X; действует от [N, N] до [2, 1].

Самым простым решением было бы использование условного перехода, однако это не соответствует нашей начальной установке. Здесь необходимо дополнительно отметить, что из всех стандартных кусочно-заданных функций, как отмечено выше, мы будем использовать только модуль (y = |x|) и формулы собственной разработки.

В этой связи, необходимо привести наше уравнение в вид:

Здесь функция   F₁ принимает значение 1 при i, j между  [1, 1] … [N, N] , в остальных случаях – 0. В свою очередь, F₂ , наоборот, принимает значение 1 когда ячейка находится в диапазоне [N, N - 1] … [2, 1], и 0 между [1, 1] … [N, N].

В качестве аргумента функции используется переменная switcher, вычисляемая на основе значении i и j, и опять же, принимающая различные значения в вышеуказанных диапазонах.

Неплохим вариантом выглядит идея получения значений -1, 0 и/или 1 от манипуляций с координатами элемента. Тогда F₁ и F₂ содержали бы простые арифметические операции.

Итак,  мы уже использовали сложение i и j, но оно всегда дает положительное число. Почему бы сейчас не попробовать вычитание?

Заменим в нашем предыдущем листинге значение a[ik][jk].

…
int Xs = i + j - 1;     
a[ik][jk] = j – i;
…

Наиболее легким вариантом видится целочисленное деление на самого себя, но в таком случае мы столкнемся с ошибкой деления на 0, поскольку а[1][1] и а[N][N] уже содержат 0. В данном случае можно прибавить ко всем значениям N и осуществить целочисленное деление на N. Введенную переменную назовем switcher.

…
int switcher =  (j - i + N) / N;
a[ik][jk] =  switcher; // временно подставим в ячейки switcher
…

Теперь осталось немногое для завершения данного этапа, а именно создать функции  F₁  и F_2. F₁ должна возвращать такое же значение, какое ей передали в качестве аргумента, т.е. F₁ (switcher) = switcher. В таком случае F₁ (switcher) * Xs работает только для диапазона от [1, 1] до [N, N], в остальных случаях равняется нулю. Вторая часть уравнения, должна действовать наоборот. В таком случае она должна возвращать значение по модулю switcher – 1, т.е. F_2.(switcher) = |switcher – 1|.

Итак, пишем:

…
a[ik][jk] =  switcher * Xs + abs(switcher - 1) * (4 * N - 2 - Xs);
…

Проверяем:

Отлично, этот этап завершен успешно, мы получили требуемые значения для внешнего кольца массива.

2 этап. Альтернативная система координат.

Нам удалось заполнить внешнее первое кольцо требуемыми данными.   Однако, что произойдет, если мы снимем фильтр, и попытаемся вычислить данные для остальных элементов массива? Для этого необходимо удалить участок кода if (!(ik == 0 || ik == N - 1 || jk == 0 || jk == N - 1)) continue;  

Чего и следовало ожидать, следующее кольцо матрицы представлено неверными значениями, что, однако, не скрывает некую  системность в числах. Так, по крайней мере, сохранен возрастающий порядок. Более того, прирост значений составляет 1, за исключением стыка между нижним правым элементом и следующей ячейкой. Данный признак является благоприятным знаком того, что мы находимся на верном пути.

Очевидно, в ячейке [2, 2] вместо 03 должно находиться число 17, следующее значение после конечного элемента внешнего круга. Небольшое размышление показывает, что необходимо ввести определенный коэффициент, выправляющий значение ячейки в зависимости от «глубины залегания» кольца в котором он находится. Т.е. требуется вычислить степень удаленности элемента массива от центра (или наоборот) на основании номера строки и номера столбца.

Для этого, введем альтернативную систему координат. Так, обычная нумерация ячеек в матрице начинается с левого верхнего угла и аналогична двухмерной системе координат, с центром в указанном месте. Для тех, кто еще помнит древний Turbo Pascal с синим экраном – графики в нем чертились именно таким образом (да и в большинстве современных графических сред разработки ПО принята подобная система координат при работе с визуальными компонентами).

Мы же, перенесем начало координат в центр нашей матрицы следующим образом:

…
a[ik][jk] =  abs(N / 2 + 1 -  i);
…

Поскольку нам требуется только расстояние ячейки от центра, мы берем только абсолютные значения.

Сейчас элементы массива заполнены номерами строк относительно его центра. Если мы проделаем ту же самую операцию со стоблцами, мы получим такой же массив, только с вертикальным распределением значений.

Введем две новые переменные Ic, Jc (c обозначает center).

…
Ic = abs(i -  N / 2  - 1);
Jc = abs(j -  N / 2  - 1);
…

Теперь необходимо определить расстояние ячейки от центра, независимо в каком направлении она удалена.

Опять же, самый простой способ что-то узнать – попробовать вывести суммы номеров строк и столбцов.

…
a[ik][jk] =  Ic + Jc;
…

Пока ничего путного не выходит, но проявляется закономерность – чем дальше от центра, тем больше сумма. В дальнейшем это пригодится. Теперь рассмотрим, что нам даст разница между Ic и Jc.

…
a[ik][jk] =  Ic - Jc;
…

Пока также неопределенно. Но если изучить расположение чисел, и сложить поэлементно их абсолютные значения с предыдущим массивом, то получим матрицу с уникальными числами для каждого кольца.

…
a[ik][jk] =  abs(Ic – Jc) + Ic + Jc;
…

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

…
Ring = N / 2 -  (abs(Ic – Jc) + Ic + Jc) / 2; 
a[ik][jk] =  Ring;
…

Замечательно. Однако, при размере матрицы N = 6 данная формула работает не совсем корректно, так она в качестве центра считает только один элемент (что является справедливым для матриц с нечетной размерностью, как в предыдущих примерах).

N = 6

При четном размере центральный квадрат из четырех элементов должен считаться центром матрицы. Возникает вопрос: как это реализовать? Здесь к нам опять на помощь приходит целочисленное деление и кусочно-заданная функция.

Но для этого вернемся немного назад, к вычислению Ic и Jc. Попробуем запустить наш скрипт при N = 6, и посмотрим значения Ic = |i - N / 2  - 1|.

В данный момент требуется привести массив в симметричную форму. Легче всего это сделать, увеличив значение на 1 всех строк после 3-й (то есть вторую половину).

…
Ic = abs(i - N / 2  - 1) + (i - 1) / (N /2) * ((N-1) % 2);
…

Вторая часть данного выражения работает только при четном N, и приводит номера строк в симметричный вид по горизонтали. 

Тоже самое проделываем и Jc.

…
Jc = abs(j - N / 2  - 1) + (j - 1)/(N /2) * ((N-1) % 2);
…

В итоге получаем симметричную матрицу уже по вертикали. Теперь проверяем значение Ring для матрицы с размером 6.

… 
a[ik][jk] =  Ring;
…

Все работает нормально. Второй этап завершен.

3 этап. Соединение. На данном этапе мы объединим полученные данные и выведем искомую матрицу (через формулу вычисления значения заданной ячейки).

Но для начала, вернемся к первому этапу и выведем на экран:

…
a[ik][jk] =  switcher * Xs + abs(switcher - 1) * (4 * N - 2 - Xs);
…

Дальнейший порядок действий представляется следующим образом: привести значения элементов внутренних колец к «спиральному» виду (т.е. заполнить их начинания с верхнего левого угла по спирали возрастающими значениями от 1), далее вычислить коэффициент прироста, которой обеспечит нормальный переход от одного кольца к нижестоящему (в нашем примере от ячейки [1,2] к ячейке [2,2], при этом меняя значение с 16 на 17)

Описывать естественным языком или даже академическим стилем математические вычисления крайне тяжело, но еще сложнее понимать чьи-то не всегда ясные наработки изложенные подобным образом. Поэтому если вы столкнулись с определенными трудностями в понимании материала, просто читайте дальше и следите за вычислениями – скоро все будет ясно.

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

Xs = (i – Ring) + (j – Ring) – 1.

Теперь попробуем вывести

…
a[ik][jk] =  switcher * Xs + abs(switcher - 1) * (4 * N - 2 - Xs);
…

Как вы могли заметить верхние и правые элементы внутренних колец пришли в требуемое значение. Однако нижние и левые стороны приняли гораздо большие значения, чем ожидалось. Это связано с тем, что выражение 4 * N - 2 - Xs во второй части функции вычисляет значения исходя из размера внешнего кольца, которое нужно уменьшить, заменив на N – 2 * Ring. То есть формула будет работать в соответствии с размером текущего кольца.

Итак:

…
a[ik][jk] =  switcher * (Xs) + abs(switcher - 1) * (4 * (N - 2 * Ring) - 2 - Xs); 
…

Все правильно, с первой задачей данного этапа мы справились . Остается последний шаг – вычислить коэффициент прироста, о котором говорилось выше, и связать между собой все кольца.

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

Coef = N² – (N – 2Ring)²

Воспользовавшись правилами разложения квадратов разницы элементов ((a−b)²=a²−2ab+b²), можно сократить до 4Ring(N - Ring).

Теперь этот коэффициент нужно добавить к нашей основной формуле.

…
a[ik][jk] =  Coef + switcher * (Xs) + abs(switcher - 1) * (4 * (N - 2 * Ring) - 2 - Xs);
…

Требуемый результат!

Полный код участка вычисления значений представлены следующим образом:

…
int switcher =  (j - i + N) / N;
int Ic = abs(i - N / 2  - 1) + (i - 1)/(N /2) * ((N-1) % 2);
int Jc = abs(j - N / 2  - 1) + (j - 1)/(N /2) * ((N-1) % 2);
int Ring = N / 2 - (abs(Ic - Jc) + Ic + Jc) / 2;
int Xs = i - Ring + j - Ring - 1;    
int Coef =  4 * Ring * (N - Ring);
a[ik][jk] =  Coef + switcher * Xs + abs(switcher - 1) * (4 * (N - 2 * Ring) - 2 - Xs);
…

Можно конечно попробовать развернуть единую формулу, заменив дополнительные переменные выражениями, основанными только на i, j и N, и попытаться сократить ее математическими методами. Но поверьте мне, я пытался, получилось нечто такое невообразимое (на половину страницы), что я решил оставить все как есть.

(PS. Не работает только при N = 1, так как возникает ошибка деления на ноль. Но как говорится, «чуть-чуть – не считается»).