Адаптивные сетки, представляя дискретную модель физической области, нашли важное применение при решении сложных задач численного моделирования во многих областях науки, таких как гидродинамика, динамика твердого тела, теория упругости, термодинамика и т.д. [2, 3, 6, 7, 13].
Интенсивное развитие вычислительной техники значительно расширило возможности численных алгоритмов решения прикладных задач. Известны многочисленные теоретические и экспериментальные результаты, показывающие преимущества использования адаптивных сеток при решении сложных многомерных задач [12] и др.
Традиционными методами построения адаптивных сеток являются метод эквираспределения [6], Томпсона [13], эллиптический метод [3], алгебраические методы [7], конформных отображений [2] и т.д. В основе большинства применяемых методов лежат теории дифференциальных уравнений, вариационного исчисления и дифференциальной геометрии. При этом построение адаптивных сеток требует решения сложных систем нелинейных дифференциальных уравнений с частными производными, что накладывает ряд известных ограничений.
В настоящее время высокую эффективность построения адаптивных сеток с заданной плотностью на сложной физической области демонстрируют нейросетевые алгоритмы [4, 5, 8, 11] и т.д. Развитие современных нейросетевых моделей обязано классической теории самоорганизующихся карт Кохонена (SOM – Self-Organizing Maps) (например, [9, 10]). Этот метод обладает рядом преимуществ по сравнению с классическими методами, он прост в реализации и может использовать произвольные начальные данные, например равномерные треугольные, прямоугольные или шестиугольные сетки. Данный метод не зависит от размерностей вычислительного или физического пространств. Это значит, что с помощью этого метода можно строить адаптивные сетки любой размерности на одно-, двух- и трехмерных областях.
Наиболее важными свойствами модели SOM являются, во-первых, ее способность строить отображения произвольной размерности на основе самоорганизации и методов теории вероятностей и, во-вторых, наличие у нее внутреннего параллелизма, как у представителя нейронных сетей. Однако при использовании базовой модели SOM возникает ряд проблем, упоминающихся во многих источниках, например [5], к числу которых относятся граничный эффект, наличие «мертвых» нейронов и нарушение сохранения топологии при отображении. Эти проблемы ограничивают возможности применения модели SOM, в том числе и для построения адаптивных сеток. Было установлено, что граничный эффект приводит к неплотному прилеганию сетки к границе физической области, «мертвые» нейроны появляются при выходе узлов сетки за границу области, а нарушение сохранения топологии приводит к вырожденности сетки.
Для решения указанных проблем разработаны модифицированные методы, в основе которых лежат идеи чередования базового алгоритма для внутренних и внешних узлов сетки [11], а также усовершенствования функции соседства нейронов с усложнением параметров функций расстояния и скорости обучения относительно дискретного времени [4, 5].
Настоящая работа посвящена совершенствованию нейросетевых алгоритмов построения адаптивных сеток на сложных пространственных областях, характеризующихся пористой структурой, программой реализации и апробации предложенных алгоритмов, выработке рекомендаций по выбору параметров обучения. В основе применяемых алгоритмов лежит идея чередования базового алгоритма для граничных и внутренних узлов сетки, а также использования упрощенных формул функции соседства нейронов. Вычислительные эксперименты показали эффективность рекомендаций по выбору параметров обучения нейросети на односвязной плоской области [1] и в случае трехмерной полой области.
Основные идеи
Как известно, нейросетевая модель SOM состоит из карты нейронов, определенным образом расположенных друг относительно друга в пространстве, которое называется выходным. Наиболее распространенным расположением нейронов в карте является их расположение в узлах прямоугольной сетки. С каждым нейроном связан весовой вектор, являющийся элементом пространства, которое называется входным пространством. Весовой вектор корректируется в процессе итерационного алгоритма обучения модели SOM, используя элементы обучающего множества, состоящего из конечного числа векторов входного пространства. Считается, что обучающее множество – это выборка из вероятностного распределения, которое требуется аппроксимировать весовыми векторами нейронов в результате обучения.
Рассмотрим алгоритм модифицированной модели SOM для пространственной адаптивной сетки на многосвязной пространственной области. Под многосвязной областью будем подразумевать область G, для которой существует замкнутая поверхность, не стягиваемая в пределах области непрерывной деформацией в точку.
Пусть в области G выделены внутренние границы G1, G2, ..., Gn–1, внешняя граница Gn и внутренность Gn+1. Такое разбиение позволит использовать базовую модель SOM на всей области G, на ее внешней Gn и внутренних G1, G2, ..., Gn–1 границах, а также внутренности Gn+1. Пусть Wijk = N×N×N – это равномерная кубическая сетка, зафиксированная на области G. Подмножество индексов граничных узлов Wn описывает узлы Gn, т.е. узлы, расположенные на границе области G, Wn+1 описывает Gn+1 – внутренние узлы G, а W1, W2, ..., Wn–1 – узлы, расположенные на соответствующих внутренних границах G1, G2, ..., Gn–1.
В этих обозначениях модифицированный алгоритм можно описать следующей последовательностью шагов.
1. Устанавливаются начальные веса wijk всех нейронов в вершинах трехмерной сетки.
2. На первой макроитерации (s = 1), соответствующей дискретному времени t ∈ [1, t0], выполняется:
2.1. Генерируется случайная точка x(t) во всей области G;
2.2. Определяется нейрон победитель в евклидовой метрике d. Фиксируются сеточные координаты нейрона победителя iBMU, jBMU, kBMU (BMU – best matching unit) для образца x(t);
2.3. Настраиваются новые весовые значения нейронов сети по формуле
(1)
где
(2)
– функция соседства нейронов; t – номер итерации; δ(t) – функция скорости обучения; h(d,t) – функция расстояния.
3. Осуществление вырезов:
3.1. Удаление нейронов, попавших в полости области G.
3.2. Определение внутренних границ сетки W1, W2, ..., Wn–1.
4. На каждой макроитерации s > 1 выполняется:
4.1. В течение T1(s) итераций:
4.1.1. Генерируется случайная точка x(t) на одной из границ (внешней/внутренней) области G, т.е. точка принадлежащая одной из областей G1, G2, ..., Gn;
4.1.2. Определяется нейрон победитель BMU из граничных узлов сетки, индексы которых входят в Wi, где i – индекс области Gi, из которой был выбран случайный объект;
4.1.3. Настраиваются новые весовые значения граничных узлов сетки, индексы которых входят в Wi, по формуле (1).
4.2. В течение T2(s) итераций:
4.2.1. Генерируется случайная точка x(t) во всей области G;
4.2.2. Определяется нейрон победитель BMU среди всех узлов сетки;
4.2.3. Если в шаге 4.2.2 победил нейрон сетки, индексы которого входят в Wi, 
, то случайно сгенерированная точка x(t) заменяется на этот нейрон;
4.2.4. Настраиваются новые весовые значения внутренних узлов сетки, индексы которых входят в Wn+1, по формуле (1).
5. Повторяются макроитерации до тех пор, пока изменения положений узлов не станут достаточно малыми.
Первая макроитерация в классической терминологии соответствует этапу упорядочивания. Здесь важным показателем является правильность предварительного расположения сетки в области. Последующие макроитерации уточняют расположение сетки относительно границы и внутренней области G.
В работе рассмотрена упрощенная функция соседства (2). Так за скорость обучения выбрана убывающая степенная функция δ(t) = t–0,2, за функцию расстояния выбрана функция Гаусса
где 
– радиус обучения, параметр a:
(3)
Следует отметить, что использование формулы (3) в качестве параметра a(t) продиктовано работой [1].
В представленном алгоритме наибольшее время занимает процесс нахождения ближайшего узла сетки к заданной точке пространства. Поскольку поиск BMU применяется не один раз, то минимизация времени и числа операций является важной подзадачей. Чтобы ускорить вычисления, производится дихотомия над узлами сетки – разбиение области сетки на одинаковые элементарные ее подобласти – кубики. При этом сложность рекурсивного поиска можно оценить как O(log), когда как сложность простого перебора узлов сетки N×N×N равна O(N3).
Еще одним относительно трудоемким процессом при построении сетки является процесс корректировки узлов. Так как сдвиг узла определяется по формуле, которая не зависит от других узлов сетки, то эти операции можно выполнять параллельно. Поэтому при реализации алгоритма операции корректировки узлов целесообразно распараллелить на N3 процессов.
Вычислительные эксперименты
Вычислительные эксперименты проводились на симметричной выпуклой области G – куб с полостью в форме сферы. Выбор симметричной формы для вычислительных экспериментов обусловлен очевидностью выбора правильного расположения адаптивной сетки в процессе вариаций параметров обучения.
При анализе качества построения сетки в модифицированном алгоритме целесообразно разделить этапы упорядочивания и уточнения. В следующих вычислительных экспериментах исследовалось качество построения сети при вариациях параметров a и t0. За основу на первом этапе обучения принята [1] линейная зависимость параметра a от дискретного времени t.
На рис. 1 приведены итоги стадии упорядочивания в зависимости от числа итераций t0. Как видно (рис. 1, в), с увеличением числа итераций сеть, приспосабливаясь к особенностям области, на 100000 итерации приняла оптимальное положение. Под оптимальным положением понимается правильная ориентация сетки – для квадратной области в соответствии с расположением вершин квадрата.
Уменьшение нижнего предела параметра a на этапе упорядочивания при неизменных t0 и a(1) приводит к ухудшению качества предварительного построения сетки (рис. 2). Это объясняется тем, что с уменьшением a(t0) уменьшается и радиус обучения σ(t), а следовательно, основное смещение узлов сетки происходит в окрестности сгенерированных точек. В результате нарушается не только гладкость сетки, но и увеличивается ее размер, что значительно ухудшает качество адаптации на втором этапе обучения.
а 
б 
в
Рис. 1. Предварительное расположение сетки в области при a(1) = 100, a(t0) = 25 по истечении t0 итераций: а – t0 = 1000; б – t0 = 5000; и в – t0 = 100000
а 
б 
в
Рис. 2. Предварительное расположение сетки в области при a(1) = 100, t0 = 100000 и различных значений параметра a: а – a(t0) = 25; б – a(t0) = 10; в – a(t0) = 5
а 
б 
в
Рис. 3. Результат применения алгоритма для полой трехмерной области: а – вид сбоку; б – вид снизу; в – в разрезе
Исследовалась зависимость качества построения адаптивной сетки от верхнего предела параметра a при неизменных t0 и a(t0). Выявлено, что при уменьшении верхнего предела параметра a качество адаптации сетки значительно снижается. Вычисления показывают, что при достаточно большом радиусе начального обучения, например при a(1) = 100, в первых итерациях сеть «сжимается» к случайно сгенерированным точкам, уменьшаясь в размерах. После чего медленно разворачивается, адаптируясь к особенностям области. При уменьшении начального радиуса сеть минует стадию «сжатия», и начальные размеры сети уменьшают возможности правильной адаптации.
На стадии уточнения эксперименты подтвердили [1] целесообразность использования экспоненциально убывающей функции (3), расположенной достаточно близко к оси дискретного времени в последних итерациях. Следует отметить, что использование линейно убывающей функции сохраняет качество адаптации, но значительно увеличивает время расчетов. Пример построения адаптивной сетки для функции (3) приведен на рис. 3.
Выводы
Представленная реализация нейросетевого алгоритма доказывает возможность построения адаптивных сеток на сложных физических областях, характеризующихся пористой структурой. В основе модифицированного алгоритма лежит идея чередования базового алгоритма для граничных и внутренних узлов сетки в сочетании с упрощенной функцией соседства. Для получения лучшего результата в модифицированном алгоритме SOM рекомендуется в качестве параметра a Гауссовой функции расстояния использовать линейно убывающую функцию на первой макроитерации и экспоненциально убывающую в последующей. Вычислительные эксперименты по адаптации симметричной трехмерной полой области подтвердили рекомендации по выбору параметров обучения нейросети на односвязной плоской области [1].