CC BY
0
О методах решения проблем компьютерного моделирования при ограниченных и расширенных библиотеках встроенных
функций
Т. Р. Абдульмянов, С.А. Соловьев, О.В. Соловьева, В.К. Епифанов
Аннотация— Рассматривается возможность применения результатов сложных аналитических расчетов, полученных при помощи компьютерной системы MAPLE, в моделировании динамики газа в системе ANSYS. Показано, что, несмотря на ограниченную возможность использования библиотеки встроенных специальных функций в системе ANSYS, результаты сложных аналитических расчетов могут быть подключены в эту систему как пользовательские функции (UDF). Пользовательская функция была разработана на основе аналитической теории либрационных движений коорбитальных астероидов Юпитера Б. Гарфинкеля. На основе этой же теории в данной работе получены формулы динамической вязкости и плотности кольца для случая долгопериодических возмущений. Пользовательская функция (UDF) была получена для динамической вязкости при помощи САВ MAPLE в виде ряда Бесселевых функций нулевого и первого порядков. Затем, Бесселевы функции были преобразованы в тригонометрические функции. Представлены основные требования при разработке и подключении UDF для системы ANSYS, а так же программный код UDF динамической вязкости. Результаты данной работы могут быть использованы для моделирования движения газа внутри кругового диска в системе ANSYS и в исследовании динамики газа и пыли в газопылевых дисках молодых одиночных звезд.
Ключевые слова—Модель либраций, модель турбулентности, динамическая вязкость, компьютерные системы моделирования, пользовательская функция (UDF).
I. Введение
В современные компьютерные средства и методы моделирования физических процессов интегрированы
Статья получена 30 октября 2023.
Т. Р. Абдульмянов, канд. физ.-мат. наук, доцент, Казанский государственный энергетический университет, Казань, Россия (abdulmyanov .tagir@yandex.ru).
С. А. Соловьев, канд. физ.-мат. наук, зав. кафедрой ИСИС, Казанский государственный энергетический университет, Казань, Россия (solovev.sa@kgeu.ru).
О. В. Соловьева, канд. физ.-мат. наук, доцент Казанский государственный энергетический университет, Казань, Россия (solovyeva.ov@yandex.ru).
В. К. Епифанов, ИП-индивидуальный предприниматель, руководитель, г. Йошкар-Ола, Россия (rendy183@gmail.com).
результаты теоретических исследований всех наук, которые были получены ранее в течение нескольких предыдущих столетий. По этой причине можно предполагать, что применение этих средств позволит решить многие задачи, решение которые не представлялось возможным до настоящего времени. Для применения аналитических методов в настоящее время применяются компьютерные системы аналитических вычислений (САВ) [1], особенностью которых является то, что в библиотеки встроенных функций этих систем загружены все основные специальные функции. В результате этого становится возможным получить аналитическое решение многих задач в виде сходящихся рядов. Однако, САВ предназначены для решения вполне определенного круга задач, и не следует в них пытаться найти решение сложных систем уравнений в частных производных. В этом отношении более предпочтительным представляются компьютерные системы, в основу которых положены численные методы решения систем уравнений в частных производных. Но и в этом случае не все проблемы оказываются разрешимыми. Даже в таких универсальных средствах компьютерного
моделирования как ANSYS невозможно вносить изменения в базовую систему уравнений в частных производных. Следовательно, как и несколько десятилетий назад, в настоящее время необходимо искать компромиссы между аналитическими и численными методами моделирования физических процессов, но уже для систем компьютерных вычислений.
Ограничивают область применения системы ANSYS не только жестко привязанная с ней система уравнений в частных производных. К таким ограничениям относятся также размеры моделируемого объекта, единицы измерения величин, формы определения начальных и граничных условий краевых задач и др. Но сложнее всего решать задачи в тех случаях, когда в систему необходимо вносить принципиальные изменения, которые определяются спецификой физического процесса. Например, динамическая вязкость может быть определена экспериментально и внесена в систему компьютерного моделирования [2]-[6]. Если же вязкость экспериментально не определяется, то проблема замкнутости системы основных
дифференциальных уравнений остается открытой. В таком случае можно вязкость попытаться определить из общих уравнений и законов гидродинамики [7]-[10]. Подключение полученной вязкости в систему ANSYS можно выполнить через определение пользовательской функции в соответствующем разделе системы. Для этого UDF должна быть представлена в специальной форме как подпрограмма, написанная на языке C++. В связи с этим в данной работе рассматриваются проблемы, связанные с преобразованием формулы для динамической вязкости в UDF. Формула вязкости была получена ранее в монографии [11]. Эта формула представлена в виде ряда, содержащего функции Бесселя первого и второго порядка. Применение функций Бесселя недопустимо в определении UDF. Поэтому они заменяются их аналогами. Оценивается погрешность такой аппроксимации при помощи графического представления исходной и преобразованной формулы вязкости. Окончательный вид формулы динамической включен в программный код на языке Си с расширением c и представлен в Приложении. Программный код записан как UDF согласно требованиям системы ANSYS.
II. РЕЗУЛЬТАТЫ ОПРЕДЕЛЕНИЯ ДИНАМИЧЕСКОЙ ВЯЗКОСТИ ПРИ ПОМОЩИ САВ MAPLE
Применение компьютерных средств аналитического вычисления позволило, при решении ряда задач, упростить процесс моделирования как с точки зрения возможных ошибок аналитических вычислений, так и сократить трудоемкость вычислений. Однако САВ являются вспомогательными средствами
моделирования. В таких САВ не прикреплены базовые системы уравнений в частных производных. Можно ли считать это преимуществом таких систем или их недостатком? Если САВ применяется как вспомогательное средство моделирования, то это будет ее преимуществом. Но САВ невозможно применить в тех случаях, когда необходимо численное интегрирование сложных систем дифференциальных уравнений. Что же собой тогда представляют современные аналитические методы моделирования? Основу современных аналитических методов моделирования составляют классические аналитические методы моделирования. Для обоснования этого тезиса коротко рассмотрим алгоритм вывода формулы динамической вязкости из уравнений Навье - Стокса. При этом предполагается, что динамическая вязкость заранее неизвестна или ее невозможно определить экспериментально. Например, если газопылевой диск планетной системы образовался 4.5 млрд. лет назад, то неизвестно какова была динамическая вязкость такого диска. Система уравнений гидродинамики в таком случае оказывается незамкнутой и не решается как аналитическими, так и численными методами. Для таких случаев в монографии [11] предложено разделить решение уравнений гидродинамики на два этапа. На первом этапе определяются классы с характерными скоростями движений внутри газового диска. Плотность диска определяется на основе теории плотностных волн [7] в результате решения краевой задачи для волнового
уравнения [11]. После этого, для определенных скоростей движения и плотности, из уравнений Навье -Стокса выводится формула динамической вязкости. На этом этапе были используются вычислительные и графические ресурсы системы MAPLE. На втором этапе динамическая вязкость, определенная на первом этапе, применяется в моделировании движения газа в круговом диске при помощи системы ANSYS. Определенная таким образом формула динамической вязкости представлена в виде пользовательской функции (UDF). Затем, для моделирования движения газа в диске с полученной динамической вязкостью, UDF загружается в соответствующую картотеку вспомогательных файлов системы ANSYS. Следовательно, система уравнений гидродинамики становится замкнутой. Однако при этом необходимо учитывать условия, при выполнении которых была получена формула динамической вязкости. Эти условия специальные для каждого физического объекта и не могут быть универсальными.
Рассмотрим эти условия для конкретных объектов из области небесной механики, а именно, для газопылевых дисков, внутри которых формируются планетные системы. Для того, чтобы определить классы объектов с характерными скоростями движения необходимо, чтобы для описания движения этих объектов были получены общие решения в виде функций времени. Классическими задачами небесной механики, решения которых хорошо известны, являются: задача двух тел и ограниченная задача трех тел. При помощи решений этих задач ранее были определены классы K и Kr орбит с характерным движением [11]. Для этих классов были получены [11] соответствующие им формулы динамической вязкости. Для определения вязкости в Kr применялось аналитическое решение ограниченной задачи трех тел, полученное Гарфинкелем и др. в ряде статей [12]-[16]. Решение Гарфинкеля дает полное представление о топологии либрационных орбит и по этой причине больше других решений подходит для определения динамической вязкости.
Определение двух или трех газопылевых дисков с характерными скоростями движения не может решить все проблемы, которые связаны с моделированием начальных этапов эволюции протопланетных дисков. Например, решение Гарфинкеля содержит резонансный параметр, который является константой и, в результате, либрационные орбиты не будут эволюционирующими. То есть, изменения резонансного параметра должны быть получены при помощи данных наблюдений. Кроме того, в газопылевом диске не существует второго тела из задачи двух или трех тел. Тогда каким образом при помощи уравнений задачи трех или n тел можно представить движение газа и пыли в газопылевых дисках? Такая проблема рассматривалась и была решена еще Гауссом [17] (кольцо Гаусса) при определении вековых возмущений. Плотность распределения и соответствующее кольцо можно определить и для долгопериодических возмущений. Однако неясно как такой переход сделать при моделировании эволюции кольца? В интеграторы численного решения задачи n тел загружены уравнения задачи n тел. Они не могут заменить уравнения гидродинамики. Если же наоборот, исходить из
уравнении гидродинамики, то с некоторого момента необходимо будет рассматривать гравитационное взаимодействие тел, которые формируются в газопылевых дисках. Проблема может показаться вовсе неразрешимой. Однако в этом случае результаты моделирования, полученные в небесной механике, можно применить в определении динамической вязкости для уравнений гидродинамики и, тем самым, перенести наблюдаемые свойства движения небесных тел в газопылевые диски. Для такого перехода от моделирования движения п тел к моделированию движения газа и пыли в газопылевых дисках, при помощи ANSYS, определяется динамическая вязкость.
III. ПРЕОБРАЗОВАНИЕ ФОРМУЛЫ
ДИНАМИЧЕСКОЙ ВЯЗКОСТИ ДЛЯ ЕЕ ПРИМЕНЕНИЯ В СИСТЕМЕ ANSYS
Динамическая вязкость, учитывающая либрации, была получена при помощи уравнения Навье - Стокса в монографии [11] в следующем виде:
n(r, z, t) = —
4
gMp - rI/4 dP r dr
2mG370 + e • G2 sin ф
(G + Г0)3
(I)
G = G0 + G02Aj + 2G(3Al2 /3 + ..., A1 =- 1/3>/6m[a2 -70(l*)]1/2sgn(dll/dt),
At = t -1 = -
"J
J6mG0 l* [a2 - 70]I/2 9
d ^ТШ - 4g0(1 -ll)
(2)
_ (s + q)n _
где g - гравитационная постоянная, М - масса Солнца (протозвезды), m - масса Юпитера, v, e - истинная аномалия и эксцентриситет Кеплеровой орбиты астероида (пылевых частиц), sin ф - является функцией времени t. Функция f0(X ) характеризует либрационное движение [16]. Аргумент этой функции X = X - щ (s + q)/s является функцией времени t и выражен в виде ряда по функциям Якоби [11], X, ф, щ - средняя долгота, истинная аномалия пылевых частиц и среднее движение Юпитера, a - резонансный параметр, s и q - целые числа, определяющие соизмеримости средних движений Юпитера и пылевой частицы.
По определению Ньютона динамическая вязкость в гидродинамике - это трение между слоями жидкости (или газа). Тогда следствием такого трения будет образование вихрей в движении жидкости или газа. В таком смысле динамическая вязкость, которая представлена в формуле (1), влияет на динамику частиц газа и пыли как трение: там, где вязкость п(г) максимальная, относительная скорость движения частиц будет минимальной и наоборот. В результате такого трения будут уменьшаться взаимные расстояния между пылевыми частицами, что приведет к аккумуляции
пылевых частиц и формированию плотных тел разных размеров.
Плотность распределения по закону пропорциональности Ада = тДУТ для случая долгопериодических либраций определим при помощи решения задачи идеального резонанса [12], [16]. Применяя величину Д/, определенную в формуле (2), получим:
-A t--
dl
T -Jë^G^ [a2 - //
'-Pr (r, l ), r
(3)
Т = 7] (а,т) У Г1(0,т) = 0(«2 - /а)^*,
где Т1(0, т) = 4пУ(27т)12 - период малых колебаний вблизи точек Лагранжа Ь4 и Ь5 , ас - критическое значение резонансного параметра а, то есть, при а > ас астероид выходит из резонанса 1/1, Т - либрационный период астероида. Разрешая уравнение (3) относительно р„ получим следующее выражение плотности кольца для случая долгопериодических либраций:
Pr (r,l ) =
I
dl
T^G02[a2 - Л]1'
-/, r
(4)
dl 1
Согласно формуле плотности pr, полученной из равенства (3), масса m Юпитера распределяется вдоль круговой орбиты так, что гравитационные системы Солнце - Юпитер - астероид и Солнце - кольцо с плотностью pr оказываются одинаковыми.
При помощи системы MAPLE формула (1) динамической вязкости для Kr дисков была представлена в виде ряда [11] (Приложение Б). Выражение вязкости в виде этого ряда содержит функции синус, косинус и функции Бесселя нулевого и первого порядка. Функции синус и косинус можно сохранить в выражении вязкости при создании UDF, а функции Бесселя необходимо выразить через тригонометрические функции. Для больших значений r функции Бесселя Jn (r) порядка n c точностью можно представить в виде следующей формулы [18]:
Jn (r ) C0S! r -?-4 J+ O(r-3/2).
(5)
Затем, в выражении динамической вязкости (1) вместо функции Бесселя нулевого и первого порядка при помощи оператора подстановки системы MAPLE и формулы (5) была выполнена замена функций Бесселя на косинусы. Результат такой замены представлен в Приложении Б монографии [11]. Резонансный параметр а, при помощи данных о взаимных расстояний планет Солнечной системы, представлен в виде линейной функции времени: a(t) = a0 ± ai (t - t0), где ао -константа, которая определяется в модели идеального резонанса, а постоянный коэффициент a1 характеризует интенсивность притока пылевых частиц на экваториальный диск. На Рис. 1 изображены графики погрешности: функций Jn(r) при n = 0 и n =1 минус аппроксимирующая функция (4). Согласно Рис.1 погрешность аппроксимации функций Бесселя при r > 1 будет порядка 10-2 и будет уменьшаться с увеличением r.
/
Рис. 1. Графики погрешности аппроксимации функций Бесселя Jo(r) (черный), J\(r) (красный) функциями (4), г (а. е.) полярный радиус.
На Рис. 2 представлены графики изменения динамической вязкости в зависимости от г для исходного ряда (Рис. 2а) и для этого же ряда после замены функций Бесселя (Рис. 2 б) на их аппроксимирующие функции (4).
б
Рис. 2. Графики изменения динамической вязкости n(r, t) (н-с/м2): а) для исходного ряда; б) для этого же ряда после замены функций Бесселя. Время t имеет фиксированное значение: t = 3 млн. лет.
Вертикальные оси на Рис. 2 - это динамическая вязкость t) (н-с/м2) для фиксированного значения времени t = 3 млн. лет, а горизонтальные - полярные радиусы r (a. е.). Сравнение этих двух графиков показывает, что функция вязкости после преобразования в целом сохраняет общую структуру и характер изменения вязкости. Различие графиков на Рис. 2а и Рис. 2б связано с тем, что при суммировании слагаемых ряда погрешность аппроксимации суммируется и дает определенное искажение амплитуды изменения вязкости. Это следует иметь в виду при интерпретации
результатов моделирования в системе ANSYS.
IV. О ВЫЧИСЛИТЕЛЬНЫХ РЕСУРСАХ И ГРАФИЧЕСКИХ ВОЗМОЖНОСТЯХ САВ MAPLE И СИСТЕМЫ ANSYS
Система аналитических вычислений MAPLE разработана для аналитического решения дифференциальных уравнений и простейших систем дифференциальных уравнений. Она имеет расширенную библиотеку специальных функций. Как следствие, аналитические решения во многих случаях получаются либо в виде некоторых сложных функций, в которые входят специальные функции, либо в виде сходящихся рядов специальных функций. Существование в системе MAPLE расширенной библиотеки специальных функций, возможности выполнять сложение, умножение и подстановку рядов и др. делают эту систему незаменимым вспомогательным средством в моделировании многих физических процессов. Поиск решения задач в аналитической форме необходим в тех случаях, когда требуется прогнозирования тех или иных явлений на больших промежутках времени. К недостаткам системы MAPLE можно отнести
а
ограниченные графические возможности этой системы. Эта система не является полноценным средством моделирования физических процессов. В ней невозможно учитывать специфику материалов, геометрию объекта исследования и др.
Компьютерная система ANSYS предназначена для решения большого круга задач моделирования физических процессов. Эта система является полноценным средством моделирования процессов при помощи численных методов интегрирования тех дифференциальных уравнений, которые загружены в эту систему. Кроме того, система ANSYS имеет большие графические возможности в выборе и конструировании геометрии объектов моделирования. Ограничения в применении системы ANSYS связаны с тем, что пользователи этой системы не могут, в случае необходимости, вносить изменения в основные уравнения системы. Разработчики таких систем не могут обеспечить выбор основных дифференциальных уравнений по объективным причинам. Например, неясно как разрешить проблему выбора основных уравнений в небесной механике, астрофизике или гидродинамике. Может быть, не следует соединять несоединимые разделы наук? Однако это будет означать, что моделирование непрерывной эволюции пылевых дисков и этапа формирования тел невозможно. Но в газопылевых дисках небесные тела непрерывно формируются из газа и пыли. Тогда исследователями должны быть созданы модели и для таких случаев. В определенной степени решить такую задачу моделирования можно при помощи подключения в систему ANSYS пользовательских функций (UDF). Приведенная в Приложении к данной работе пользовательская функция (UDF) разработана для моделирования динамической вязкости и исследования процесса формирования тел в газопылевых дисках [19], [20].
V. ОСНОВНЫЕ ТРЕБОВАНИЯ ПРИ РАЗРАБОТКЕ UDF ДЛЯ СИСТЕМЫ ANSYS
Рассмотрим некоторые детали разработки и подключения UDF для системы ANSYS. Пользовательская функция (User Defined Function) должна быть представлена как программа, написанная на языке C или C++. Особенностями оформления загружаемой UDF являются входящие в нее специальные операторы, которые связывают UDF с геометрией моделируемого объекта и сеткой (Mesh):
#include "udf.h"
DEFINE_PROPERTY(cell_viscosity,c,t) {
real mu_current;
real T = CURRENT_TIME;
real xc[ND_ND];
real x0 = 0.0;
real y0 = 0.0;
real r;
C_CENTROID(xc,c,t);
r = sqrt(pow((xc[1] - x0),2.) + pow((xc[2] - y0),2.)); mu_current = (см. Приложение)
return mu_current;
}
Перед загрузкой UDF вязкости была создана геометрия объекта (диск) с определением его размеров, а также расчетная сетка (Рис. 3а). Кроме того, заранее был подготовлен файл с программой UDF (Приложение) и загружен в рабочую зону системы ANSYS компьютера. Только после этого, на третьем этапе моделирования (SETUP), выбирается основная модель (например, k-e модель): в верхней строке меню выбираем User Defined. Далее открывается окно выбора вида UDF: интерпретатор или компилятор (Рис. 2а). Выбираются параметры моделирования и через кнопку Show Fewer Options открываем Environment, где и определяется маршрут загрузки UDF из каталога загрузки (Рис.2а). Затем, закрытся окно выбора UDF (клавиша Close). Следующий этап расчета (Solution): выбор шага интегрирования, числа шагов интегрирования и др. Если все параметры загружены корректно, то выполняется расчет с параллельным выводом невязки для каждого параметра модели (Рис. 2б). При успешном завершении расчета выводится сообщение Done ниже графиков невязки.
Основные уравнения, которые применяются в расчетах и моделировании физических процессов при помощи системы ANSYS, рассматриваются в одном из двух основных руководств, а именно ANSYS Fluent Theory Guide. В разделе 4 этого руководства о турбулентном движении показано, что турбулентная вязкость ^ определяется при помощи усреднения уравнения Навье - Стокса. В отличие от нее динамическая вязкость д, которая рассматривается в данной работе и представлена в виде UDF в приложении, получена из уравнения Навье - Стокса без усреднения. То есть, подключение такого UDF необходимо для решения конкретной задачи о движении газа внутри диска. Такая вязкость не является универсальной для всех задач. По этой причине часто эффективная вязкость рассматривается как сумма нескольких составляющих, некоторой постоянной составляющей д0, динамической вязкости д, молекулярной вязкости дт: д0 + д + дт. Постоянная составляющая до определяется как абсолютная величина минимума вязкости д. То есть, согласно Рис. 2 д0 = п0 = 10 (н-с/м2). Применение конкретной составляющей вязкости означает, что влияние других составляющих вязкости существенно не оказывают влияние на процесс. В противном случае, необходимо найти способы определения этих составляющих. Для каждого конкретного случая необходимо учитывать специфику решаемой задачи и исследуемого объекта. Основные модели турбулентного движения, включенные в ANSYS (k-e модель, k-œ модель и др.), рассматриваются в статьях [21] - [27]. Второе руководство для пользователей системы ANSYS - это ANSYS Fluent Users Guide. Примеры разработки и подключения UDF можно найти в третьем руководстве для пользователей
этой системы: ANSYS Fluent UDF Manual.
VI. ЗАКЛЮЧЕНИЕ
Моделирование физических процессов и явлений представляет собой общую задачу, как для разработчиков компьютерных систем моделирования, так и пользователей таких систем. Пользователи систем [9], [10] в процессе моделирования сталкиваются c такими проблемами, которые относятся к задачам разработчиков. Их изучение может способствовать разработке новых средств компьютерного моделирования. В этом отношении в настоящей работе предлагается алгоритм решения одной из важных задач эволюции газопылевых дисков. Из этого алгоритма вытекает следующая задача: необходимо разработать такую систему моделирования, в которой уравнение Навье - Стокса не усредняется, а динамическая вязкость выводится из общих уравнений гидродинамики.
В результате анализа, проблем применения средств компьютерного моделирования с ограниченными и расширенными библиотеками специальных функций, были решены следующие задачи:
1) при помощи метода Рольфса [7] и теории волн плотности определена плотность сжимающейся оболочки протозвезды [19];
2) для определенной таким образом плотности и для характерных скоростей движения из уравнения Навье -Стокса аналитическими методами получены динамические вязкости для трех классов орбит [11]. Динамические вязкости представлены в виде рядов, содержащих функции Бесселя;
3) для случая Kr дисков получено выражение вязкости, не содержащее функции Бесселя и аппроксимирующее исходную функцию вязкости. Исходная вязкость получена на основе теории либраций Гарфинкеля [13] и определения плотности (4) аналога кольца Гаусса для долгопериодических либраций.
4) преобразованное выражение динамической вязкости представлено в виде UDF для подключения в систему ANSYS;
5) перечислены основные детали разработки и подключения UDF для системы ANSYS
а
б
Рис. 3. Рабочее окно системы ANSYS: а) на третьем этапе моделирования (SETUP); б)
на этапе расчета (Solution).
БЛАГОДАРНОСТИ
Авторы настоящей работы выражают благодарность Организационному комитету Международного конгресса «Современные проблемы компьютерных и информационных наук» ВМК МГУ за приглашение к участию в работе конгресса в г. Москва (Россия) с 30 ноября по 2 декабря 2023г.
Библиография
[1] Т.Р. Абдульмянов, Моделирование вязкой динамики в пылевых дисках вокруг молодых звезд при помощи системы MAPLE. в сб.: Информационные технологии в образовании и науке (ИТОН-2021). Казань, 2021, с. 7 - 13. EDN: TJMIPW
[2] J. O. Hinze, Turbulence. McGraw-Hill Publishing Co., New York. 1975.
[3] Kim S. E., Choudhury D., Patel B. Computations of Complex Turbulent Flows Using the Commercial Code Fluent. In: Salas M.D., Hefner J.N., Sakell L. (eds) Modeling Complex Turbulent Flows. ICASE/LaRC Interdisciplinary Series in Science and Engineering, Vol 7. Springer, Dordrecht, 1999. P. 259-276. https://doi.org/10.1007/978-94-011-4724-8_15
[4] S. Sarkar, M. Y. Hussaini, Computation of the sound generated by isotropic turbulence. NASA Contract Report 93-74, NASA Langley Research Center, Hampton, VA. 1993.
[5] P. Spalart, S. Allmaras, A one-equation turbulence model for aerodynamic flows, Technical Report AIAA-92-0439, American Institute of Aeronautics and Astronautics. 1992.
[6] D. C. Wilcox, Turbulence Modeling for CFD. DCW Industries, Inc. La Canada, California. 1998.
[7] К. Рольфс, Лекции по теории волн плотности. Москва, Наука,1980.
[8] T. Abdulmyanov, On the forms of accretion, of interstellar gas and dust, during the formation of single stars and their planetary systems // Open Astronomy. 2021, Vol. 30, No 3, p. 83 - 90. https://doi.org/10.1515/astro-2021-0010
[9] Т.Р. Абдульмянов, Моделирование квазипериодического движения газа в дисках при помощи системы ANSYS. в сб.: Информационные технологии в образовании и науке (ИТОН-
2022). Казань, 2022, с. 10 - 15. EDN: AXFCDU
[10] Т.Р. Абдульмянов, Влияние динамической вязкости на квазипериодическое движение тел в газопылевых дисках. в сб.: Информационные технологии в образовании и науке (ИТОН-
2023). Казань, 2023, с. 6 - 11. EDN: MGMDGU
[11] Т.Р. Абдульмянов, Моделирование процесса формирования вихревых движений в газопылевых дисках при помощи систем компьютерных вычислений. Казань: Изд. КГЭУ, 2022, 391 c. EDN: JSZNHK
[12] B. Garfinkel, A. Jupp A., C. Williams, A Recursive von Zeipel algorithm for the Ideal Resonance Problem // Astronomical Journal. 76, 157 (1971). https://doi.org/10.1086/111099
[13] B. Garfinkel, Theory of libration // Celestial Mechanics. Vol. 13, issue 2, 229-246 (1976). https://doi.org/10.1007/BF01232726
[14] B. Garfinkel, Theory of the Trojan asteroids: Part I // Astronomical Journal. 82, 368 (1977). https://doi.org/10.1086/112149
[15] B. Garfinkel, Theory of the Trojan asteroids // Celestial Mechanics: Part II. Vol. 18, issue 3, p. 259-275 (1980). https://doi .org/10.1007/BF01230167
[16] Т.Р. Абдульмянов, Определение характеристик либрационных движений астероидов вблизи соизмеримостей 1/1, 4/3, 3/2, 7/3, 5/2, 3/1 средних движений астероидов и Юпитера // Астрономический вестник. Исследования Солнечной системы. 35. № 5, 449 (2001).
[17] K.F. Gauss, Werke . 1866. - Bd. III. 655 p.
[18] Е. П. Аксенов, Специальные функции в небесной механике. Москва, Наука, 1986, 320 с.
[19] Т.Р. Абдульмянов, Общая модель образования небесных тел от начальной конденсации газопылевых частиц до формирования «зародышей» планет // Вестник Московского университета. Серия 3: Физика. Астрономия. 74. 3-15 (2019). EDN: UTZRGV
[20] Т.Р. Абдульмянов, Крупномасштабная структура газопылевых дисков и устойчивость газодинамического равновесия пылевых оболочек молодых звезд // Астрофизический бюллетень. 75, N° 2, с. 132-139 (2020). EDN: YQCJIZ
[21] B. Aupoix, P. R. Spalart. Extensions of the Spalart-Allmaras turbulence model to account for wall roughness // International Journal of Heat and Fluid Flow. 2003. Vol. 24, issue 3, p. 454-462. https://doi .org/10.1016/S0142-727X(03)00043-2
[22] J. Dacles-Mariani, G. G. Zilliac, J. S. Chow, and P. Bradshaw. Numerical/Experimental Study of a Wingtip Vortex in the Near Field // AIAA Journal. 1995. Vol. 33, issue 9, p. 1561-1568. https://doi .org/10.2514/3.12826
[23] F. R. Menter. Two-Equation Eddy-Viscosity Turbulence Models for Engineering Applications // AIAA Journal. 1994. Vol. 32, issue 8, p. 1598-1605. https://doi .org/10.2514/3.12149
[24] F. R. Menter. Review of the SST Turbulence Model Experience from an Industrial Perspective // International Journal of Computational Fluid Dynamics. 2009. Vol. 23, issue 4, p. 305-316. https://doi.org/10.1080/10618560902773387
[25] Orszag S.A., Yakhot V., Flannery W.S., Boysan F., Choudhury D., Maruzewski J., Patel B. Renormalization-group modeling and turbulence simulations. In: So R.M.C., Speziale G., Launder B. E. (eds) International Conference on Near-wall turbulent flows. Tempe, Arizona: Elsevier; 1993. p. 1031-1046.
[26] T.-H. Shih, W.W. Liou, A. Shabbir, Z. Yang, and J. Zhu. A New -Eddy-Viscosity Model for High Reynolds Number Turbulent Flows -Model Development and Validation // Computers & Fluids. 1995. Vol. 24, issue 3. p. 227-238. https://doi.org/10.1016/0045-7930(94)00032-T
[27] V. Yakhot, S. A. Orszag. Renormalization Group Analysis of Turbulence I Basic Theory // Journal of Scientific Computing. 1986. Vol. 1, No. 1, p. 1-51.
On methods for solving computer modeling problems with limited and extended libraries of
built-in functions
Tagir Abdulmyanov, Sergey Soloviev, Olga Solovyova, Vyacheslav Epifanov
Abstract— The possibility of using the results of complex analytical calculations obtained using the MAPLE computer system in modeling gas dynamics in the ANSYS system is considered. It is shown that, despite the limited ability to use the library of built-in special functions in the ANSYS system, the results of complex analytical calculations can be connected to this system as user-defined functions (UDFs). The custom function was developed based on B. Garfinkel's analytical theory of librational motions of Jupiter co-orbital asteroids. Based on the same theory, in this work, formulas for the dynamic viscosity and density of the ring are obtained for the case of long-period disturbances. A user-defined function (UDF) was obtained for dynamic viscosity using MAPLE CAB as a series of zero- and first-order Bessel functions. Then, Bessel functions were converted to trigonometric functions. The basic requirements for developing and connecting a UDF for the ANSYS system and the program code for a dynamic viscosity UDF are presented. The results of this work can be used to simulate the movement of gas inside a circular disk in the ANSYS system, as well as to study the dynamics of gas and dust in gas-dust disks of young single stars.
Keywords - Libration model, turbulence model, dynamic viscosity, computer simulation systems, user-defined function (UDF).
References
[1] Abdulmyanov T. Modeling viscous dynamics in dust discs around young stars using the maple system. In: Proceedings on Digital Technologies for Teaching and Learning (DTTL-2021). Kazan; 2021, p. 7-13. (In Russ.) EDN: TJMIPW
[2] Hinze J. O. Turbulence. McGraw-Hill Publishing Co., New York. 1975.
[3] Kim S.E., Choudhury D., Patel B. Computations of Complex Turbulent Flows Using the Commercial Code Fluent. In: Salas M.D., Hefner J.N., Sakell L. (eds) Modeling Complex Turbulent Flows. ICASE/LaRC Interdisciplinary Series in Science and Engineering, vol 7. Springer, Dordrecht; 1999. p. 259-276. https://doi .org/10.1007/978-94-011 -4724-8_15
[4] Sarkar S., Hussaini M. Y. Computation of the sound generated by isotropic turbulence. NASA Contract Report 93-74, NASA Langley Research Center, Hampton, VA. 1993.
[5] Spalart P., Allmaras S. A one-equation turbulence model for aerodynamic flows, Technical Report AIAA-92-0439, American Institute of Aeronautics and Astronautics. 1992.
[6] Wilcox D.C. Turbulence Modeling for CFD. DCW Industries, Inc. La Canada, California. 1998.
[7] Rohlfs K. Lectures on Density Wave Theory. Lecture Notes in Physics, Vol. 69. Springer, Berlin, Heidelberg; 1977. 184 p. https://doi.org/10.1007/3-540-08448-7
[8] Abdulmyanov T. On the forms of accretion, of interstellar gas and dust, during the formation of single stars and their planetary systems. Open Astronomy. 2021;30(3):83-90. https://doi.org/10.1515/astro-2021-0010
[9] Abdulmyanov T. Simulation of quasi-periodic gas motion in disks using the ansys system. In: Proceedings on Digital Technologies for
Teaching and Learning (DTTL-2022). Kazan; 2022, p. 10-15. (In Russ.) EDN: AXFCDU
[10] Abdulmyanov T. Influence of dynamic viscosity on quasi-periodic motion of bodies in gas and dust disks. In: Proceedings on Digital Technologies for Teaching and Learning (DTTL-2023). Kazan; 2023, p. 6-11. (In Russ.) EDN: MGMDGU.
[11] Abdulmyanov T. Modelirovanie processa formirovanija vihrevyh dvizhenij v gazopylevyh diskah pri pomoshhi sistem kompjuternyh vychislenij. Kazan; KGEU Publ.; 2022. 391 p. (In Russ.) EDN: JSZNHK
[12] Garfinkel B., Jupp A., Williams C. A Recursive von Zeipel Algorithm for the Ideal Resonance Problem. Astronomical Journal. 1971;76:157. https://doi.org/10.1086/111099
[13] Garfinkel B. Theory of libration. Celestial Mechanics. 1976;13(2):229-246. https://doi.org/10.1007/BF01232726
[14] Garfinkel B. Theory of the Trojan asteroids: Part I. Astronomical Journal. 1977;82:368. https://doi.org/10.1086/112149
[15] Garfinkel B. Theory of the Trojan asteroids: Part II. Celestial Mechanics .1980;18(3):259-275. https://doi.org/10.1007/BF01230167
[16] Abdulmyanov T. Opredelenie harakteristik libracionnyh dvizhenij asteroidov vblizi soizmerimostej 1/1, 4/3, 3/2, 7/3, 5/2, 3/1 srednih dvizhenij asteroidov i Jupitera. Solar System Research (Astronomicheskii Vestnik). 2001;35(5):449. (In Russ.)
[17] Gauss K.F. Werke. 1866. Bd. III. 655 p.
[18] Aksenov E.P. Spetsial'nye funktsii v nebesnoj mekhanike [Special functions in celestial mechanics]. Moscow: Nauka; 1986. 320 p. Zbl 0667.70002 (In Russ.)
[19] Abdulmyanov T.R. A General Model of Celestial Body Formation from Initial Condensation of Gas and Dust Particles to the "Embryos" of Planets. Moscow University Physics Bulletin. 2019;74(4):309-322. https://doi.org/10.3103/S0027134919040039
[20] Abdulmyanov T.R. Large-scale structure of gas-and-dust disks and stability of gas-dynamic equilibrium of the dust envelopes of young stars. Astrophysical Bulletin. 2020;75(2):132-139 (2020). (In Russ., abstract in Eng.) EDN: YQCJIZ
[21] Aupoix B., Spalart P.R. Extensions of the Spalart-Allmaras turbulence model to account for wall roughness. International Journal of Heat and Fluid Flow. 2003;24(3):454-462. https://doi.org/10.1016/S0142-727X(03)00043-2
[22] Dacles-Mariani J, Zilliac G.G., Chow J.S., Bradshaw P. Numerical/Experimental Study of a Wingtip Vortex in the Near Field. AIAA Journal. 1995;33(9):1561-1568. https://doi.org/10.2514Z3.12826
[23] Menter F.R. Two-Equation Eddy-Viscosity Turbulence Models for Engineering Applications. AIAA Journal. 1994;32(8):1598-1605. https://doi .org/10.2514/3.12149
[24] Menter F.R. Review of the SST Turbulence Model Experience from an Industrial Perspective. International Journal of Computational Fluid Dynamics. 2009;23(4):305-316. https://doi.org/10.1080/10618560902773387
[25] Orszag S.A., Yakhot V., Flannery W.S., Boysan F., Choudhury D., Maruzewski J., Patel B. Renormalization-group modeling and turbulence simulations. In: So R.M.C., Speziale G., Launder B. E. (eds) International Conference on Near-wall turbulent flows. Tempe, Arizona: Elsevier; 1993. p. 1031-1046.
[26] Orszag S.A., Yakhot V., Flannery W.S., Boysan F., Choudhury D., Maruzewski J., Patel D. Renormalization Group Modeling and Turbulence Simulations. In: International Conference on Near-Wall Turbulent Flows. Tempe, Arizona; 1993. p. 1031-1046.
[27] Shih T.-H., Liou W.W., Shabbir A., Yang Z., Zhu J. A New - Eddy-Viscosity Model for High Reynolds Number Turbulent Flows -Model Development and Validation. Computers & Fluids. 1995;24(3):227-238. https://doi.org/10.1016/0045-7930(94)00032-T
[28] Yakhot V., Orszag S.A. Renormalization Group Analysis of Turbulence I Basic Theory // Journal of Scientific Computing. 1986; 1(1):1-51.
Приложение
#include "udf.h"
DEFINE_PROPERTY(cell_viscosity,c,t)
{
real mu_current;
real T = CURRENT_TIME;
real xc[ND_ND];
real x0 = 0.0;
real y0 = 0.0;
real r;
C_CENTROID(xc,c,t);
r = sqrt(pow((xc[1] - x0),2.) + pow((xc[2] - y0),2.));
mu_current = (3.718*((0.0129*(0.0791*cos(0.263*T)*sqrt(1/r)*cos(0.480*r-0.785) +0.0521*cos(0.603*T)*sqrt(1/r)*cos(1.104*r-0.785)+0.0416*cos(0.946*T)* sqrt(1/r)*cos(1.730*r-0.785)+0.0357*cos(1.289*T)*sqrt(1/r)*cos(2.358*r-0.785) +0.0317*cos(1.633*T)*sqrt(1/r)*cos(2.986*r-0.785)+0.0288*cos(1.976*T)*sqrt(1/r)* cos(3.614*r-0.785)+0.0266*cos(2.320*T)*sqrt(1/r)*cos(4.242*r-0.785)+0.0248*cos(2.663*T)* sqrt(1/r)*cos(4.870*r-0.785)+0.0233*cos(3.007*T)*sqrt(1/r)*cos(5.498*r-0.785) +0.0221*cos(3.350*T)*sqrt(1/r)*cos(6.126*r-0.785)+0.0210*cos(3.694*T)*sqrt(1/r)* cos(6.755*r-0.785)+0.0201*cos(4.037*T)*sqrt(1/r)*cos(7.383*r-0.785)+0.0193*cos(4.381*T)* sqrt(1/r)*cos(8.011*r-0.785)+0.0186*cos(4.724*T)*sqrt(1/r)*cos(8.639*r-0.785)+ 0.0180*cos(5.068*T)*sqrt(1/r)*cos(9.268*r-0.785+0.0174*cos(5.412*T)*sqrt(1/r)* cos(9.896*r-0.785)+0.0169*cos(5.755*T)*sqrt(1/r)*cos(10.524*r-0.785)+0.0164*cos(6.099*T)* sqrt(1/r)*cos(11.153*r-0.785)+0.0159*cos(6.442*T)*sqrt(1/r)*cos(11.781*r-0.785) +0.0155*cos(6.786*T)*sqrt(1/r)*cos(12.409*r-0.785)))/r-0.750*pow(r,0.25)* (0.161*cos(4.724*T)* sqrt(1/r)*cos(8.639*r-2.355)-0.166*cos(5.068*T)*sqrt(1/r)*cos(9.268*r-2.355)-0.172*cos(5.412*T)*sqrt(1/r)*cos(9.896*r-2.355)-0.181*cos(5.755*T)*sqrt(1/r)* cos(10.154*r-2.355)-0.1830*cos(6.099*T)*sqrt(1/r)*cos(11.153*r-2.355)-0.188*cos(6.442*T)*sqrt(1/r)*cos(11.781*r-2.355)-0.193*cos(6.786*T)*sqrt(1/r)* cos(12.409*r-2.355)-0.0380*cos(0.263*T)*sqrt(1/r)*cos(0.480*r-2.355)-0.0575*cos(0.603*T)* sqrt(1/r)*cos(1.104*r-2.355)-0.0715*cos(0.946*T)*sqrt(1/r)*cos(1.730*r-2.355)-0.0841*cos(1.289*T)*sqrt(1/r)*cos(2.358*r-2.355)-0.0946*cos(1.633*T)*sqrt(1/r)* cos(2.986*r-2.355)-0.104*cos(1.976*T)*sqrt(1/r)*cos(3.614*r-2.355)-0.112*cos(2.320*T)* sqrt(1/r)*cos(4.242*r-2.355)-0.120*cos(2.663*T)*sqrt(1/r)*cos(4.870*r-2.355)-0.128*cos(3.007*T)* sqrt(1/r)*cos(5.498*r-2.355)-0.135*cos(3.350*T)*sqrt(1/r)*cos(6.126*r-2.355)-0.1424*cos(3.694*T)*sqrt(1/r)*cos(6.755*r-2.355)-0.1489*cos(4.037*T)*sqrt(1/r)* cos(7.383*r-2.355)-0.155*cos(4.381*T)*sqrt(1/r)*cos(8.011*r-2.355))))/((-1*cos(T)* sin(1.+sin(T))/pow(2-2*cos(1+sin(T))),3./2.)+cos(T)*sin(1+sin(T)))*pow(1.0008+ 0.816*sqrt(0.00191-0.000418*T-0.000955/sqrt(2.-2.*cos(1.+sin(T)))+0.000955*cos(1.+sin(T)))-0.000185*T-0.000424/sqrt(2.-2*cos(1.+sin(T)))+0.000424*cos(1.+sin(T))),3.); return mu_current;
}
