CC BY
20
www.volsu.ru
КОМПЬЮТЕРНОЕ МОДЕЛИРОВАНИЕ
DOI: https://doi.oгg/10.15688/j■volsu1.2017.2.4
УДК 517.957+514.752 ББК 32.973.26-018.2
УНИВЕРСАЛЬНЫЙ ПРОГРАММНЫЙ КОМПЛЕКС
ДЛЯ РЕШЕНИЯ МНОГОМЕРНЫХ ВАРИАЦИОННЫХ ЗАДАЧ1
Елена Геннадиевна Григорьева
Кандидат физико-математических наук,
доцент кафедры компьютерных наук и экспериментальной математики, Волгоградский государственный университет e_grigoreva@mail.ru, e_grigoreva@volsu.ru, kiem@volsu.ru просп. Университетский, 100, 400062 г. Волгоград, Российская Федерация
о
см <
<
к к ¡г к
ч
к
Владимир Александрович Клячин
Доктор физико-математических наук,
заведующий кафедрой компьютерных наук и экспериментальной математики, Волгоградский государственный университет klchnv@mail.ru, klyachin.va@volsu.ru, kiem@volsu.ru просп. Университетский, 100, 400062 г. Волгоград, Российская Федерация
Алексей Александрович Клячин
Доктор физико-математических наук, заведующий кафедрой математического анализа и теории функций, ^ Волгоградский государственный университет
^ klyachin-aa@yandex.ru, matf@volsu.ru
ч просп. Университетский, 100, 400062 г. Волгоград, Российская Федерация
Ы
га т
Аннотация. В статье рассматривается задача численного расчета равновесных поверхностей произвольной топологии в классе полигональных моделей. Такого рода поверхности являются решением вариационной задачи на минимум функционала типа площади при наличии ограничений интеграль-@ ного вида. В статье представлено программное решение таких вариационных
задач в двух вариантах. В первом — в виде пакета на языке Python, который реализован на базе пакета NumPy и интегрирован в среду 3D-моделирования Blender. Второй вариант предназначен главным образом для решения непараметрических вариационных задач и выполнен стандартными средствами C++.
Ключевые слова: экстремальная поверхность, кусочно-линейная аппроксимация, триангуляция, пакет NumPy, краевая задача.
Одним из широко распространенных методов решения краевых задач математической физики является вариационный метод, основанный на минимизации некоторого интегрального функционала в том или ином классе функций. Например, при заданных краевых условиях можно свести интегрирование уравнений статической теории упругости к отысканию минимума потенциальной энергии деформации упругого тела. Другим примером может служить задача поиска поверхности, имеющей минимальную площадь и покрывающей фиксированный или максимальный объем. Данная задача возникает в практике проектирования современных сооружений в качестве покрытий или ограждений [7].
Известно, что метод конечных элементов (МКЭ), использующий построение триангуляции расчетной области, является одним из самых применяемых методов приближенного решения вычислительных задач математической физики, в том числе и при решении упомянутых вариационных задач. Стоит назвать такие популярные программные комплексы, как COMSOL, ANSYS, ИСПА, «Диана», COSMOS, MSC/NASTRAIN, SAMSAF и другие, повсеместно используемые в инженерных расчетах, в которых алгоритмы имитационного моделирования существенно используют МКЭ. В настоящей работе мы приводим свою реализацию вариационного метода, основанного на аппроксимации кусочно-линейными функциями и поверхностями. Мы исходим из идеи создания универсальной программной системы, основная часть которой не должна зависеть от используемой триангуляции, интегрального функционала, расчетной области и граничных условий.
1. Вычисление экстремальных поверхностей произвольной топологии
Пусть Ф : R3 ^ R — положительная, однородная степени 1, выпуклая функция такая, что Ф(^) = Ф(-£,). Рассмотрим ограниченную область D С R3 с кусочно-гладкой границей и пусть ф(ж) G С(D) — положительная непрерывная функция. Для гладкой поверхности М С R3 такой, что М С dD рассмотрим величину
где £ — вектор единичной нормали к поверхности М, а(х) — положительная функция. Нас будет интересовать решение вариационной задачи на минимум этого значения при условии, что граница дМ фиксирована, а интеграл
Введение
м
(1)
остается неизменным.
Поверхность М будем аппроксимировать поверхностью М, которая задается конечным множеством точек Р = {Р.\ € М3,г = 1,...,М} и совокупностью треугольников Т = {Т*,] = 1,...,К}. Каждый треугольник Т* однозначно определяется тремя целыми числами к^, ^ ,т* € {1,...,Ж}, такими, что точки Р^. ,Р1} ,Рт} являются вершинами треугольника Тj.
Если М — кусочно-линейная аппроксимация поверхности М, то приближенное значение введенного функционала будет вычисляться согласно формуле
к
Р(М) = £ )|Т*| + а(р*)|Т*|, (2)
j=l
где М представляет собой объединение треугольников Тj; — нормаль к плоскости этих треугольников, а |Т*|,р* — площади и центры треугольников соответственно. Поверхность М можно рассматривать как точку Р в пространстве ЬМ размерности 3N, где N — число вершин на М. Если задать в каждой точке-вершине поверхности вектор к, то этим векторам будет соответствовать вектор в ЬМ. Тогда оказывается, если функционал рассмотреть как числовую функцию : ЬМ ^ М, то имеет место формула
(ЯР \. 1 к 1 к
Ж( Р0 = 2^ (УФ(2^|Т* I)) х 1 *+ 2 С а(Р*)** Х 1 *,к) + / *=1 *=1
1 к
+ -(£ Уа(р*)|Т*|,к>, 3 *=1
где сумма идет по всем треугольникам, имеющим вершину р. € М, £,*•, как и выше, обозначает нормаль к соответствующему треугольнику и I* — вектор, соединяющий вершины треугольника, не совпадающие с . таким образом, что эта вершина при обходе остается слева. Заметим, что если для ]-го треугольника обозначить через а*, Ь* векторы, направленные вдоль сторон треугольника из вершины . , то эта формула может быть переписана в более простом виде, с использованием однородности функции Ф(£,)
/яр У 2 к
(Р^ =2(ЕУф(а* х Ъ*) х I*,к>.
Пример 1. Рассмотрим случай
Ф(4) = III, а(х) = 0, ф(ж) = 0. Тогда, как нетрудно посчитать, получим
(!(р))=2 (£ ^х .*>• (3)
Этот случай соответствует вычислению минимальных кусочно-линейных поверхностей.
Пример 2. Рассмотрим случай
Ф(4) = |£|, а(х) = 0, ф(ж) = 1. Тогда, как нетрудно посчитать, получим
2 Ь х I, ,h). (4)
з=
Этот случай соответствует моделированию кусочно-линейных поверхностей постоянной средней кривизны — равновесных поверхностей раздела сред с постоянным давлением без учета гравитационных сил [10; 12].
Пример 3. Рассмотрим случай
Ф(Ь) = |Ь|, а(х) = д •хз, ф(х) = 1.
Тогда, как нетрудно посчитать, получим
1 1
' X l3,h) + -<>>*,ез)Ь x l3,h) + (5)
/яр V i к i к
ж(Р) = 2<Е Ь X h,h) + 2 (£<Pj, ез)Ь X l3,h) + V / 3 = 1 i=l
1 к
+ 1 <Е ез1Т, l,h)
3 = 1
Этот случай соответствует моделированию капиллярных кусочно-линейных поверхностей — равновесных поверхностей раздела сред с постоянным давлением с учетом гравитационных сил [10; 12].
Для итерационного процесса градиентного спуска мы также должны учесть, что интеграл (1) остается неизменным. Ясно, что из постановки вариационной задачи следует, что это условие будет выполнено при малых деформациях Н(х) поверхности М, если
J ф(х)к(х)ё,8 = 0.
м
Это условие для кусочно-линейной поверхности можно записать в виде
К
£h(# )l Tj | = 0,
=
где р* — как и выше, центр треугольника Т^. Полученное соотношение представляет собой линейное уравнение для значений смещений к в вершинах кусочно-линейной поверхности М. В пространстве ЬМ множество решений этого уравнения задает гиперплоскость С, ортогональную проекцию на которую обозначим через пс : ЬМ ^ С. Тогда формулы итерационного процесса можно записать в виде
рк+1 = рк - (^Р(Рк)),
где величины ек определяют смещения вдоль направления, обратного к направлению градиента VР(Рк).
2. Объектно-ориентированная модель вычислений
Реализацию вычислений было решено выполнить, преследуя две цели. С одной стороны, весь программный комплекс должен быть интегрирован в какую-либо систему визуализации, а с другой стороны, желательно было получить код вычислений, по максимуму не зависящий от такой интеграции. В связи с этим была выбрана объектно-ориентированная модель организации вычислений на базе шаблона проектирования классов «мост» [9;11]. Этот шаблон используется в случаях, когда требуется отделить абстракцию от ее реализации так, чтобы то и другое можно было изменять независимо. В нашей задаче было решено реализацию интегрировать в пакет Blender и выполнить на языке программирования Python. Это, в частности, позволит задавать топологию искомой экстремальной поверхности, а также моделировать искомую поверхность в интерактивном режиме работы в окне ЭЭ-вида программы Blender. Кроме этого, реализация основана на использовании числовых массивов модуля NumPy. Этот модуль представляет собой расширение языка программирования Python, реализующее множество различных операций при работе с многомерными массивами данных. Несмотря на то что этот пакет предназначен для использования его в контексте интерпретируемого языка, тем не менее массивы NumPy являются Си-подобными [8]. Так что эти массивы могут управляться стандартными инструкциями языка Python. Удобство использования массивов NumPy, в частности, состоит и в том, что для этих массивов перегружены арифметические операции и математические функции.
Опишем структуру программного комплекса (см. рис. 1). Весь код организован в пакет MinimalSurface, содержащий несколько модулей с требуемыми классами. Некоторые модули, ранее описанные в работе [2], были выбраны за основу и использованы для решения задачи моделирования экстремальных поверхностей:
ndimvar.py,
geometry.py,
triangulation.py.
Модуль ndimvar.py реализует абстрактный класс NDimVar, инкапсулирующий обобщенный итерационный процесс, в частности итерационный процесс метода градиентного спуска. Приведем отрывок кода, который реализует отдельную итерацию градиентного спуска.
def iteration(s,x,k):
if s.sublin: # если задача на условный экстремум
s.create_subLin() # на каждом шаге
# строится гиперплоскость С return x-s.epsilon(k)*s.proj(s.X(x))
# вычисляем новую точку
# с учетом проекции
# на гиперплоскость C
else:
return x-s.epsilon(k)*s.X(x) # вычисляем новую точку
# без учета проекции
# на гиперплоскость C
Рис. 1. Диаграмма вычислительных классов
В подклассах класса NDimVar требуется реализация метода вычисления градиента минимизируемой функции, который зависит от функционала F (M). Модуль geometry.py реализует классы для вычисления геометрических характеристик граней кусочно-линейной поверхности: длины сторон, площади и углы треугольников. Модуль triangulation.py реализует базовый класс Triangulation для представления в системе вычислений кусочно-линейной поверхности, а также автоматизирует процесс вычисления граничных точек поверхности [1]. Дополнительно к перечисленным модулям реализованы модули minsquare.py, tentsurface.py. В модуле minsquare.py реализованы классы MinimalSurface и CapillarSurface, производные класса NDimVar, в которых переопределены виртуальные методы вычисления по формулам (3)-(5). В модуле tentsurface.py представлен класс MinimalTent для аналогичных вычислений так называемых тентовых поверхностей, как экстремалей функционала, соответствующего функциям
Ф(£) = |£|, а(х) = Aq + А -хз, ф(х) = 0.
Кроме перечисленных выше классов, непосредственно задействованных в вычислениях, реализованы классы интеграции в пакет Blender и предоставляющие взаимодействие среды Blender с указанными классами в соответствии с шаблоном «мост». Такая схема позволяет не только изменять и совершенствовать интерфейс пользователя, но и интегрировать классы вычислений в другие системы визуализации, например, в систему matplotlib, независимо от перечисленных выше классов. Соответствующие классы реализованы в модуле panel_addon.py. В этом модуле содержатся классы MinimalSurfacePanel, FindBoundary, FindMinSurface, FindCMCSurface, FindCapillarSurface, FindTentSurface, Settings.
л
>=
Рис. 2. Примеры результата расчета минимальной поверхности Е.Г. Григорьева, В.А. Клячин, А.А. Клячин. Универсальный программный комплекс
Рис. 3. Пример результата расчета капиллярной поверхности
Примеры визуализации результатов расчета равновесных поверхностей представлены на рисунках 2 и 3.
3. Решение непараметрических вариационных задач
Будем рассматривать следующую вариационную задачу
J G(x,u,Du)dx + j F(x,u)dS ^ min u|aQ\r = ФЫ\г,
Q Г
где G = G(xi ,...,xn, u1 ,...,up ,...,£?), F = F (xi ,...,xn, u1 ), $ = (E\,..., ^ ),i = = 1,...,p, u = (u1, ...,up) : О ^ Rp и О — заданная область в Rn. Приведем некоторые примеры.
Пример 4. Задача нахождения минимальной поверхности, натянутой на заданный контур, формулируется так
yjl + и2 + и2у dxdy ^ min, u\dQ = (p.
В этом случае мы имеем С = у/1 + (£,})2 + для р = 1 и п = 2. Пример 5. Рассмотрим следующую краевую задачу
д ( и \ + д / и | = к для (Ж,,/) € п, (6)
М л/1 + % + Л2/ М л/1 + П + Л
(V f, V
= Ф, (7)
ЯП
ч/1 + V|2
где О — область в R2 с кусочно-гладкой границей $ О, v = (v, v2) — вектор единичной внешней нормали в точках границы дО, где он существует, f £ С2(О) П С1 (О), постоянная к > 0 и ф — заданная непрерывная функция на дО. Известно [12], что решение f описывает капиллярную поверхность, находящуюся в равновесии, над областью О с заданным углом контакта у (в случае, когда ф = cos у) между этой поверхностью и стенкой капиллярной трубки. Данная задача эквивалентна минимизации интеграла
е Ы = + ^ + ~9l dxdy + к// g2 (х> y^dxdy -J дФаз-
О О, ЯП
В этом случае G = л/1 + (£,1)2 + (^2)2 + К(u1 )2, F = и1 ф, р = 1,п = 2 и Г = М.
Пример 6. Задача определения тензора деформации упругого тела сводится к решению системы уравнений
яр
(А + Н) — + НАиг + Хг = 0, г = 1, 2, 3,
ОХг
где А, Н — постоянные, определяемые упругой средой, е = д^ + д^ + д^З, Х = = (Х^Х2,Х3) — вектор внешних сил. Решение различных краевых задач для этой системы может быть сведено к минимизации функционала
з
V (и\и2,и3) = Ц Е ^гиЧв +
дО г=1
+\т (а (± % )2+ 2н £ (яи; )2+ н Е (£+% )2+ I *1**з.
О \ \ г— 1 / г— 1 г^у г— 1
Так же как и в первых двух примерах, не сложно выписать функции С и ^:
/А -Л2 .„-А „ .,\2
/ з \2 3 2 3
с = 2 [м Е ^г +2Н Е(^)2 + Н Е + + ЕигХг
\г=1 / г=1 г=з г =1
р = Е*
з =1
4. Описание входных данных
Мы предполагаем, что расчетная область представляет собой набор треугольников {Тк}^=1, образующих триангуляцию. Обозначим через Р1,...,Рт все вершины этих треугольников. Будем рассматривать класс кусочно-линейных функций и : П ^ И., которые являются непрерывными в П ив каждом треугольнике Тк представляют собой функцию вида и(х, у) = акх + Ьку + ск. Каждая такая функция однозначно определяется набором ее значений и1,...,ит в соответствующих вершинах Р1,...,Рт триангуляции.
Теперь дадим описание входных файлов, содержащих вершины триангуляции, ее треугольники и функции, задающие краевые условия. Итак, узлы триангуляции хранятся в отдельном файле points.pnt, в котором первые 8 байт отводятся под количество точек, а дальше каждая точка занимает 20 байт — две координаты и информация о принадлежности границе:
Номер точки Координата х (8 байт) Координата у (8 байт) Принадлежность границе (0 — внутренняя точка = 0 — граничная точка (4 байта)
0 Хо Уо Ьо
1 Х1 У1 &1
2 Х2 У2 Ь2
Триангуляция этих точек хранится в другом файле triangles.trg, который имеет следующую структуру: вначале 8 байт отводится под число треугольников, а следующие данные можно изобразить в виде таблицы:
Номер треугольника Номер точки, соответствующей 1-й вершине (8 байт) Номер точки соответствующей 2-й вершине (8 байт) Номер точки соответствующей 3-й вершине (8 байт)
0 Р0 Р0 Р0
1 Pi Pi Pi
2 Р2 р'2 pi
Нужно заметить, что номера точек соответствуют их номерам из первого файла, содержащего все вершины триангуляции. Далее, так как кусочно-линейная функция однозначно определяется ее значениями в вершинах Р\,...,Рт, то мы в отдельном файле function.val храним эти значения. В начале этого файла записывается количество точек, а затем идут значения функции в том же порядке, что и точки в первом файле points.pnt. В такого типа файлах мы записываем значения функций, задающих различные краевые условия (Дирихле, Неймана, смешанные условия и др.) и результаты вычислений.
Для случая трехмерных вариационных задач структура входных файлов аналогичная. Например, файл с вершинами триангуляции будет выглядеть так: первые 8 байт отводятся под количество точек, а дальше хранятся декартовы координаты вершин триангуляции.
Так как в случае пространственных вариационных задач триангуляция области состоит из тетраэдров, то соответствующим образом меняется и структура файла triangles.trg.
5. Задание функционала и его минимизация
Минимизация функционала осуществляется с помощью метода градиентного спуска и метода покоординатного спуска с применением алгоритма «золотого сечения» для одномерной минимизации. Отметим, что вершины триангуляции являются объектами класса SQPoint, а треугольники (тетраэдры) — объектами класса SQTreug (SQTetr). Для задания функционала достаточно определить функции G и F. Для функции G создан класс SQIntegrand, имеющий виртуальные функции:
IntG(SQPoint, SQPoint, double), Gu(SQPoint, SQPoint, double), gradG(SQPoint, SQPoint, double),
IntG2(SQPoint point,SQPoint px, SQPoint py, double u, double v), gradGu2(SQPoint point, SQPoint px, SQPoint py, double u, double v), gradGv2(SQPoint point, SQPoint px, SQPoint py, double u, double v), Gu2(SQPoint point, SQPoint px, SQPoint py, double u, double v), Gv2(SQPoint point, SQPoint px, SQPoint py, double u, double v).
Все функции этого класса реализуются в наследуемых классах, которые и определяют конкретный функционал. Например, если в качестве минимизируемого функционала взять площадь поверхности
//утг^««„
П
то мы определяем класс
class SQSurfArea : public SQIntegrand {
public:
double IntG(SQPoint point,SQPoint vect, double u)
{
double value;
value=sqrt(1+vect.x*vect.x+vect.y*vect.y);
return value; }
double Gu(SQPoint point,SQPoint vect, double u)
{
return 0.0;
}
SQPoint gradG(SQPoint point, SQPoint vect, double u) {
SQPoint grd(0.0,0.0);
grd.x=vect.x/sqrt(1+vect.x*vect.x+vect.y*vect.y); grd.y=vect.y/sqrt(1+vect.x*vect.x+vect.y*vect.y); return grd;
}
};
Функция gradG() является градиентом интегранда IntG по переменным ^, а функция Gu() представляет собой производную интегранда по переменным иг. Аналогично для случая, если поиск минимума осуществляется в классе отображений, задаваемых двумя функциями и = и(х, у), v = v(x, у). Для этого предназначены функции IntG2(), gradGu2(), gradGv2(), Gu2(), Gv2().
Для задания функции F создан класс SQIntegrandB
class SQIntegrandB {
public:
virtual double IntF(SQPoint point, double u) {
return 0.0; }
};
Для минимизации функционала используются две функции minimum_func_mix() и minimum_gold_sect(), которые осуществляют один шаг в алгоритме градиентного спуска и метода «золотого сечения» соответственно. Каждая из этих функций имеет одинаковый набор аргументов. Например, рассмотрим функцию
minimum_func_mix(p_Func, p_Func2, p_Func_psi, p_Points, p_Treug, p_Ind, p_PT, Numb, Ntr);
и поясним входящие в нее переменные. Итак, переменная p_Func является указателем на переменные типа double. Этот массив содержит значения нулевого приближения в методе градиентного спуска. Длина массива равна количеству точек триангуляции. Второй аргумент p_Func2 — это указатель на массив, в который записываются результаты вычисления. Переменная p_Func_psi является указателем на массив данных, в которых хранятся значения функции ф (краевое условие Неймана). Четвертая переменная p_Points является указателем на массив, состоящий из объектов класса SQPoint_Mesh. Этот класс является потомком класса SQPoint, в котором добавлено поле bnd, отвечающее за принадлежность точки к границе расчетной области и тип граничного условия. Следующая переменная p_Treug является указателем на массив объектов класса SQTreug треугольников используемой триангуляции. Шестой аргумент p_Ind — это указатель на объекты класса SQIndexT. Этот класс содержит поля IndA, IndB, IndC, которые являются номерами точек вершин треугольника. Индекс треугольника в массиве p_Treug совпадает с соответствующим индексом элемента массива p_Ind. Переменные Numb, Ntr содержат количество точек и треугольников триангуляции. Теперь поясним значение переменной p_PT. Для сравнения значений функционала или вычисления направления спуска мы используем для каждой вершины Pi только те треугольники, которые содержат эту точку в качестве своей вершины. Поэтому в отдельный массив p_PT мы помещаем соответствующую информацию. Выглядит этот массив так:
Номер точки Количество треугольников у данной вершины (4 байта) Номера треугольников
0 щ грО грО гр 0 11 ■ 12 ■ ■■■■ 1п0
1 П\ гр! гр\ гр1 1 1 ■ 12 ■ ■■■■ 1п1
2 П2 гр2 гр2 гр2 11 ■ 12 ■>■■■■> 1П2
Дадим несколько простых примеров решения различных краевых задач, используя описанные выше процедуры.
Пример 7. Предположим, нам нужно решить следующую задачу
JJ(их + иу) dxdy ^ min, п
где П = [0, 2] х [0, 2] и функции и = и(х,у) удовлетворяют краевому условию
и(х, 0) = х2, и(х, 2) = X2 - 4, и(0, у) = -у2, и(2, у) = 4 - у2.
Понятно, что решением этой задачи является гармоническая функция и = х2 — у2. Для приближенного решения мы используем следующую триангуляцию нашей области. Разбиваем по каждой переменной 30-ю прямолинейными отрезками квадрат [0,2] х х [0, 2] на 900 квадратов. Затем каждый из них делим диагональю пополам и получаем разбиение исходного квадрата на 1 800 треугольников. Результаты численного решения
поставленной задачи можно увидеть из следующей таблицы. В ней значения искомого и приближенного решений взяты для фиксированного значения х.
Номер точки Приближенное значение Точное значение Погрешность, %
0 0,266674 0,266667 0,00272361
1 0,413341 0,413333 0,00176132
2 0,568896 0,568889 0,00126476
3 0,73334 0,733333 0,000964057
4 0,906674 0,906667 0,000761486
5 1,0889 1,08889 0,000600262
6 1,28001 1,28 0,000478745
7 1,48001 1,48 0,00037956
8 1,68889 1,68889 0,000303671
9 1,90667 1,90667 0,000243241
10 2,13334 2,13333 0,000187804
11 2,36889 2,36889 0,000133937
12 2,61334 2,61333 8,79675e-005
13 2,86667 2,86667 4,26486e-005
14 3,12889 3,12889 0
15 -1 -1 0
16 -0,995555 -0,995556 9,77761e-005
17 -0,98222 -0,982222 0,000197784
18 -0,959997 -0,96 0,000282251
19 -0,928886 -0,928889 0,000363266
20 -0,888885 -0,888889 0,000438596
21 -0,839996 -0,84 0,000521741
22 -0,782217 -0,782222 0,000629807
23 -0,71555 -0,715556 0,000758859
24 -0,639994 -0,64 0,000911162
25 -0,555549 -0,555556 0,00111466
26 -0,462216 -0,462222 0,00140913
27 -0,359993 -0,36 0,00188718
28 -0,248882 -0,248889 0,00282618
29 -0,128882 -0,128889 0,00556455
Пример 8. Рассмотрим теперь другой пример, в котором имеется краевое условие на нормальную производную. Пусть нужно найти решение уравнения
д2и д2и
— + — = 0, 0 < х < 2, 0 <у< 2, дх2 ду2
удовлетворяющее условиям
и(х, 0) = х2, и(х, 2) = х2 — 4, и(0, у) = — у2, их(2, у) = 4. Ясно, что решением этой краевой задачи является гармоническая функция и = х2 — у2,
на которой достигает своего минимума функционал
2
2 Л(и2х+иУ ^ - и(2, у)(1у'
П 0
Для поиска приближенного решения мы используем рассмотренную выше триангуляцию квадрата [0, 2] х [0, 2]. Ниже приводится таблица полученных результатов с использованием метода «золотого сечения».
Номер точки Приближенное значение Точное значение Погрешность, %
0 -1 -1 0
1 -0,995556 -0,995556 3,69844e-006
2 -0,982222 -0,982222 7,51804e-006
3 -0,96 -0,96 1,13207e-005
4 -0,928889 -0,928889 1,52227e-005
5 -0,888889 -0,888889 1,92366e-005
6 -0,84 -0,84 2,37347e-005
7 -0,782222 -0,782222 2,90668e-005
8 -0,715555 -0,715556 3,55389e-005
9 -0,64 -0,64 4,34191e-005
10 -0,555555 -0,555556 5,45779e-005
11 -0,462222 -0,462222 7,05253e-005
12 -0,36 -0,36 9,64248e-005
13 -0,248889 -0,248889 0,000147995
14 -0,128889 -0,128889 0,000300313
15 4,08992е-007 0 -
16 0,137778 0,137778 0,000312917
17 0,284445 0,284444 0,000158239
18 0,44 0,44 0,000106217
19 0,604445 0,604444 8,07618e-005
20 0,777778 0,777778 6,50983e-005
21 0,960001 0,96 5,45993e-005
22 1,15111 1,15111 4,72088e-005
23 1,35111 1,35111 4,14874e-005
24 1,56 1,56 3,70554e-005
25 1,77778 1,77778 3,33873e-005
26 2,00445 2,00444 3,05363e-005
27 2,24 2,24 2,80369e-005
28 2,48445 2,48444 2,60704e-005
29 2,73778 2,73778 2,44785e-005
Пример 9. Точным решением задачи
и + ^ I ■ ~у I =0, 0 < х < 2, 1 <у< 10,
дх\ y/l + иХ + и1) дУ\ y/l + иХ +
и(х, 1) = 1п(\/х2 + 1 + х), и(х, 10) = 1п(^ж2 + 100 + ^Х2 + 99), 0 < ж < 2, и(0,у) = 1п(у + л/у2—!), и(2, у) = 1п(^ 4 + у2 + ^3 + У2), 1 < У < 10,
является функция и(х,у) = 1п(^/х2 + у2 + у/х^+у2—!), описывающая минимальную поверхность — катеноид. Приводим соответствующие результаты вычислений. Нужно отметить, что в этом примере мы разделили оба отрезка [0, 2] и [1,10] на 20 равных частей. В таблице даны значения для фиксированного у = 5, 5.
Номер точки Приближенное значение Точное значение Погрешность, %
0 2,40606 2,40606 0
1 2,42379 2,42367 0,00478255
2 2,44704 2,4468 0,00996285
3 2,47504 2,47467 0,0146479
4 2,50695 2,50649 0,018244
5 2,54198 2,54146 0,0205082
6 2,57937 2,57882 0,0214784
7 2,61847 2,61791 0,0213603
8 2,65871 2,65816 0,0204247
9 2,6996 2,69909 0,0189375
10 2,74078 2,74031 0,0171238
11 2,78193 2,78151 0,0151552
12 2,82281 2,82243 0,013152
13 2,86323 2,86291 0,0111911
14 2,90307 2,9028 0,00931607
15 2,94223 2,94201 0,00754603
16 2,98062 2,98045 0,00588293
17 3,01821 3,01808 0,00431703
18 3,05497 3,05489 0,00283069
19 3,09089 3,09084 0,00140078
Обоснование используемого вариационного метода для уравнений минимальной поверхности и равновесной капиллярной поверхности дано в работах [5] и [6]. Результаты расчетов и их визуализация для этих уравнений приведены в работах [3] и [4].
ПРИМЕЧАНИЕ
1 Работа выполнена при финансовой поддержке РФФИ и Администрации Волгоградской области (проект № 15-41-02517 р_поволжье_а).
СПИСОК ЛИТЕРАТУРЫ
1. Григорьева, Е. Г. Алгоритмы идентификации граничных точек полигональных 3D моделей / Е. Г. Григорьева, В. А. Клячин // Геометрический анализ и его приложения : материалы III Междунар. науч. шк.-конф. — Волгоград : Изд-во ВолГУ, 2016. — C. 114-116.
2. Григорьева, Е. Г. Численное исследование устойчивости равновесных поверхностей с использованием пакета NumPy / Е. Г. Григорьева, В. А. Клячин // Вестник Волгоградского
государственного университета. Серия 1, Математика. Физика. — 2015. — № 2 (27). — C. 17-30.
3. Клячин, А. А. Визуализация расчета формы поверхностей минимальной площади / А. А. Клячин, В. А. Клячин, Е. Г. Григорьева // Научная визуализация. Электронный журнал. — 2014. — Т. 6, № 2. — C. 34-42.
4. Клячин, А. А. Визуализация устойчивости и расчета формы равновесной капиллярной поверхности / А. А. Клячин, В. А. Клячин, Е. Г. Григорьева // Научная визуализация. Электронный журнал. — 2016. — Т. 8, № 2. — C. 37-52.
5. Клячин, А. А. О равномерной сходимости кусочно-линейных решений уравнения минимальной поверхности / А. А. Клячин, М. А. Гацунаев // Уфимский математический журнал. — 2014. — № 6 (3). — C. 3-16.
6. Клячин, А. А. О равномерной сходимости кусочно-линейных решений уравнения равновесной капиллярной поверхности / А. А. Клячин // Сибирский журнал индустриальной математики. — 2015. — Т. 18, № 2 (62). — C. 52-62.
7. Михайленко, В. Е. Конструирование форм современных архитектурных сооружений / В. Е. Михайленко, С. Н. Ковалев. — Киев : Буд1вельник, 1978. — 138 с.
8. Олифант, Т. Многомерные итераторы NumPy / Т. Олифант // Идеальный код. — СПб. : Питер, 2011. — C. 341-358.
9. Приемы объектно-ориентированного проектирования. Паттерны проектирования / Э. Гамма, Р. Хелм, Р. Джонсон, Д. Влиссидес, Г. Буч. — СПб. : Питер, 2007. — 366 с.
10. Саранин, В. А. Равновесие жидкостей и его устойчивость / В. А. Саранин. — М. : Институт компьютерных исследований, 2002. — 144 с.
11. Тепляков, С. Паттерны проектирования на платформе .NET / С. Тепляков. — СПб. : Питер, 2016. — 320 с.
12. Финн, Р. Равновесные капиллярные поверхности. Математическая теория / Р. Финн. — М. : Наука, 1989. — 312 с.
REFERENCES
1. Grigoryeva E.G., Klyachin V.A. Algoritmy identifikatsii granichnykh tochek poligonalnykh 3D modeley [Identifying Algorithm for Boundary Points of 3D Models]. Geometricheskiy analiz i ego prilozheniya: materialy III Mezhdunar. nauch. shk.-konf. Volgograd, Izd-vo VolGU Publ., 2016, pp. 114-116.
2. Grigoryeva E.G., Klyachin V.A. Chislennoe issledovanie ustoychivocti ravnovesnykh poverkhnostey s ispolzovaniem paketa NumPy [Numerical Study of the Stability of Equilibrium Surfaces Using NumPy Package]. Vestnik Volgogradskogo gosudarstvennogo universiteta. Seriya 1, Matematika. Fizika [Science Journal of Volgograd State University. Mathematics. Physics], 2015, no. 2 (27), pp. 17-30.
3. Klyachin A.A., Klyachin V.A., Grigoryeva E.G. Vizualizatsiya rascheta formy poverkhnostey minimalnoy ploshchadi [Visualization of Calculation of Minimal Area Surfaces]. Nauchnaya vizualizatsiya. Elektronnyy zhurnal, 2014, vol. 6, no. 2, pp. 34-42.
4. Klyachin A.A., Klyachin V.A., Grigoryeva E.G. Vizualizatsiya ustoychivosti i rascheta formy ravnovesnoy kapillyarnoy poverkhnosti [Visualization of Stability and Calculation of the Shape of the Capillary Equilibrium Surface]. Nauchnaya vizualizatsiya. Elektronnyy zhurnal, 2016, vol. 8, no. 2, pp. 37-52.
5. Klyachin A.A., Gatsunaev M.A. O ravnomernoy skhodimosti kusochno-lineynykh resheniy uravneniya minimalnoy poverkhnosti [On Uniform Convergence of Piecewise Linear Solutions of the Minimal Surface Equation]. Ufimskiy matematicheskiy zhurnal [Ufa Mathematical Journal], 2014, no. 6 (3), pp. 3-16.
6. Klyachin A.A. O ravnomernoy skhodimosti kusochno-lineynykh resheniy uravneniya ravnovesnoy kapillyarnoy poverkhnosti [On the Uniform Convergence of Piecewise Linear Solutions of an Equilibrium Capillary Surface Equation]. Sibirskiy zhurnal industrialnoy matematiki [Journal of Applied and Industrial Mathematics], 2015, vol. 18, no. 2 (62), pp. 52-62.
7. Mikhaylenko V.E., Kovalev S.N. Konstruirovanie form sovremennykh arkhitekturnykh sooruzheniy [Designing Forms of Modern Architectural Structures]. Kiev, Budivelnik Publ., 1978. 138 p.
8. Oliphant T. Mnogomernye iteratory NumPy [Multidimensional Iterators in NumPy]. Idealnyy kod [Beautiful Code] Saint Petersburg, Piter Publ., 2011, pp. 341-358.
9. Gamma E., Helm R., Johnson R., Vlissides D., Booch G. Priemy obyektno-orientirovannogo proektirovaniya. Patterny proektirovaniya [Design Patterns: Elements of Reusable Object-Oriented Software]. Saint Petersburg, Piter Publ., 2007. 366 p.
10. Saranin V.A. Ravnovesie zhidkostey i ego ustoychivost [Equilibrium of Fluids and Its Stability]. Moscow, Institut kompyuternykh issledovaniy Publ., 2002. 144 p.
11. Teplyakov S. Patterny proektirovaniya na platforme .NET [Design Patterns for .NET Framework]. Saint Petersburg, Piter Publ., 2016. 320 p.
12. Finn R. Ravnovesnye kapillyarnye poverkhnosti. Matematicheskaya teoriya [Equilibrium Capillar Surfces. Mathematical Theory]. Moscow, Nauka Publ., 1989. 312 p.
Candidate of Physical and Mathematical Sciences,
Associate Professor of Department of Computer Science and Experimental Mathematics, Volgograd State University
e_grigoreva@mail.ru, e_grigoreva@volsu.ru, kiem@volsu.ru Prosp. Universitetsky, 100, 400062 Volgograd, Russian Federation
Doctor of Physical and Mathematical Sciences,
Head of Department of Computer Science and Experimental Mathematics, Volgograd State University
klchnv@mail.ru, klyachin.va@volsu.ru, kiem@volsu.ru
Prosp. Universitetsky, 100, 400062 Volgograd, Russian Federation
Doctor of Physical and Mathematical Sciences,
Head of Department of Mathematical Analysis and Function Theory,
Volgograd State University
klyachin-aa@yandex.ru, matf@volsu.ru
Prosp. Universitetsky, 100, 400062 Volgograd, Russian Federation
Abstract. The article considers the problem of numerical calculation of the piecewise linear surfaces M, wich is extremal for the functional type
To solve this problem, we obtain formulas for the calculation of the gradient of the functional in the space of polygonal surfaces M as a set P of its vertices
UNIVERSAL SOFTWARE FOR SOLVING MULTIDIMENSIONAL VARIATIONAL PROBLEMS
Elena Gennadievna Grigoryeva
Vladimir Aleksandrovich Klyachin
Aleksey Aleksandrovich Klyachin
M
+
1 k
+ - (£ V«(p • )\Tj l,h).
3 i=i
For the organization of calculations on the basis of these formulas, we construct a system of classes in the Python programming language. The system is organized as a class package wich is integrated in the simulation environment 3D Blender. This solution allows the user for the program Blender in interactive mode to perform the desired surface topology modeling, to set the boundary conditions and to visualize the results of calculation.
In this paper, we present our implementation of a variational method based on approximation of piecewise-linear functions and surfaces. We proceed from the idea of creating a universal software system, most of which does not depend on the use of triangulation, integral functional, computational domain and boundary conditions. In this article we present the structure of input data, which is stored in the boundary conditions, and the vertices and the triangles of the triangulation and the description of the main program procedures. We are considering a number of illustrative examples of solving boundary value problems of mathematical physics.
Key words: extremal surface, piecewise linear approximation, triangulation, NumPy package, boundary problem.
