The article presents a parallel algorithm for numerical simulation of two-phase flows in digital rock images. The algorithm is based on the phase-field method and the finite-difference approach, which includes solving the Poisson equation using conjugate gradients with a spectral preconditioner. The domain decomposition method and CUDA and MPI technologies are used to parallelize the algorithm. The algorithm is applied to calculate the relative phase permeabilities using the steady-state method for oil and water, oil and CO2 in a supercritical form. Based on the calculated two-phase permeabilities, ternary relative permeability diagrams are constructed for a three-phase system using the Stone I, Stone II, and Baker models. The justification of the choice of liquid parameters for two-phase modeling is carried out. В статье представлен параллельный алгоритм для численного моделирования двухфазных потоков в цифровых изображениях горных пород. Основу алгоритма составляют метод фазового поля и конечно-разностный подход, включающий решение уравнения Пуассона методом сопряженных градиентов со спектральным предобуславливателем. Для параллелизации алгоритма применяются метод декомпозиции расчетной области и технологии CUDA и MPI. Представленный алгоритм используется для расчета относительных фазовых проницаемостей методом установившегося режима фильтрации для нефти и воды, нефти и CO2 в сверхкритической форме. На основе рассчитанных двухфазных проницаемостей проводится построение тернарных диаграмм относительной проницаемости для трехфазной системы с использованием моделей Стоуна I, Стоуна II и Бейкера. Проводится обоснование выбора параметров жидкостей для двухфазного моделирования.

This article focuses on the development and investigation of the matrix-based modifications of the Ant Colony Optimization (ACO) method for solving parametric optimization problems using CUDA-based Graphics Processing Units (GPUs). Unlike most existing studies focusing on the Traveling Salesman Problem (TSP), the proposed method is specialized for finding optimal values of independent parameters, which leads to a different graph structure and requires the development of novel parallelization schemes. Original parallelization strategies are proposed, including a three-stage decomposition of the algorithm and introduction of “parameter-level” parallelism in the problem dimensions. Various strategies for parallelization and memory management were developed and tested: from the basic MatrixACO-Basic modification to the optimized versions of MatrixACO-Optimization (featuring warp reduction and optimized hash-table handling) and MatrixACO-Transported (utilizing transposed data storage). Experiments were conducted on NVIDIA Tesla V100 16GB SXM2 GPU, the Jetson ORIN platform, and various personal computers using multi-dimensional test functions. Results demonstrate that for a problem with 128 parameters the optimized version of MatrixACO-Optimization provides an execution time of 0.13 ms per graph layer, which corresponds to a speed-up of more than 30 times compared to the optimized CPU implementation based on OpenMP, and in the case of a very large-scale problem with 65536 parameters it speeds up by almost 45 times. For a problem with 512 parameters, the GPU realization shows 76% higher energy efficiency compared to the CPU version, as well as an almost eightfold (7.9 times) reduction in Energy-Delay Product. The study confirms the necessity of developing specialized GPU implementations for parametric tasks and reveals the excellents calability of the proposed algorithms. Статья посвящена разработке и исследованию матричных модификаций метода муравьиных колоний (ACO) для решения параметрических задач оптимизации с использованием графических процессоров (GPU) на архитектуре CUDA. В отличие от большинства существующих исследований, фокусирующихся на задаче коммивояжера (TSP), предложенный метод специализирован для поиска оптимальных значений независимых параметров, что приводит к иной структуре графа и обуславливает необходимость разработки новых схем распараллеливания. Представлены оригинальные схемы параллелизма, включающие разделение алгоритма на три этапа и введение “параметрического” параллелизма по измерениям задачи. Разработаны и протестированы различные стратегии распараллеливания и управления памятью: от базовой модификации MatrixACO-Basic до оптимизированных версий MatrixACO-Optimization (с warp-редукцией и оптимизацией хэш-таблицы) и MatrixACO-Transported (с транспонированным хранением данных). Эксперименты проводились на графических ускорителях NVIDIA Tesla V100 16GB SXM2, платформе Jetson ORIN и различных персональных компьютерах с использованием многомерных тестовых функций. Результаты показывают, что для задачи с 128 параметрами оптимизированная версия MatrixACO-Optimization обеспечивает время обработки одного слоя графа, равное 0.13 мс, что соответствует ускорению более чем в 30 раз по сравнению с оптимизированной CPU-реализацией на основе OpenMP, а в случае задачи с 65536 параметрами способствует ускорению почти в 45 раз. Для задачи с 512 параметрами GPU-реализация демонстрирует на 76% более высокую энергоэффективность по сравнению с CPU-версией, а также почти восьмикратное (в 7.9 раза) снижение Energy-Delay Product. Исследование подтверждает необходимость создания специализированных GPU-реализаций для параметрических задач и выявляет отличную масштабируемость предложенных алгоритмов.

The paper presents a numerical algorithm of estimating the formation factor for digital rock samples constructed using microtomographic images. The algorithm is based on the numerical solution of the three-dimensional Poisson equation with rapidly changing high-contrast coefficients. The arising system of linear algebraic equations turns out to be ill-conditioned, so to speed up the convergence it is necessary to use a preconditioner. The conjugate gradient method is applied to solve the system and the preconditioner is constructed as the inverse Laplace operator corresponding to the homogeneous model. In its turn, the spectral decomposition of two tridiagonal matrices appearing in the approximation of one-dimensional derivatives is used for inversion. The resulting series of one dimensional problems is solved by the Thomas algorithm. The used preconditioner can be effectively applied both to the original problem with rapidly changing high-contrast coefficients and to the problem, in which the solution is calculated only in the pore space. In the latter case, the method provides the same accuracy of results as the original one for inhomogeneous models, but it converges almost twice as fast. Implementation using graphics processors allows solving problems up to 109 unknowns with a single GPU. В работе представлен численный алгоритм оценки параметра пористости для цифровых образцов горных пород, построенных с использованием микротомографических изображений. В основе алгоритма лежит численное решение трехмерного уравнения Пуассо на с быстро изменяющимися высококонтрастными коэффициентами. Возникающая при этом система линейных алгебраических уравнений оказывается плохо обусловленной, поэтому для ускорения сходимости требуется использование предобуславливателя. Для решения применяется метод сопряженных градиентов, а предобуславливатель строится как обратный оператор Лапласа, соответствующий однородной модели. В свою очередь, для обращения используется спектральное разложение двух трехдиагональных матриц, соответствующих аппроксимации одномерных производных. Полученная серия одномерных задач решается методом прогонки. Используемый предобуславливатель может быть эффективно применен как к исходной задаче с быстро изменяющимися высококонтрастными коэффициентами, так и к задаче, в которой решение вычисляется только в поровом пространстве. В этом случае данный метод обеспечивает такую же точность результатов, как и исходный подход для неоднородных моделей, однако сходится почти вдвое быстрее исходного. Реализация с использованием графических процессоров позволяет решать задачи размером до 109 неизвестных с помощью одного GPU.

При решении современных задач математического моделирования широко применяется динамическая адаптация расчетной сетки. Особенности реализации численных методов при динамической адаптации требуют соблюдения определенных соотношений между линейными размерами смежных элементов, которые называются структурным критерием. Структурный критерий позволяет формировать буферный слой ячеек между областями, чей уровень измельчения отличается более чем вдвое. В работе рассматривается алгоритм формирования буферного слоя ячеек в листовой модели адаптации расчетной сетки, реализованный на распределенной вычислительной системе. Получены оценки алгоритмической сложности и эффективности алгоритма на распределенной вычислительной системе. Эксперименты показали возможность алгоритма формировать буферный слой в областях сложной формы, однако при этом его слабую масштабируемость с ростом числа процессов. Adaptive mesh refinement (AMR) is widely used in solving modern mathematical modeling problems. The specifics of the implementation of numerical methods for AMR require compliance with certain ratios between the linear sizes of adjacent elements, which are called the structural criterion. Structural criterion allows creation of buffer layer of cells between regions, which cell division level difference is greater than two. In this paper we discuss the algorithm for buffer layer generation for tree-based adaptive mesh refinement on distributed computational system. Estimates of the algorithmic complexity and efficiency of the algorithm on distributed computational system are obtained. Experiments have shown the possibility of algorithm for forming a buffer layer in areas of complex shape, but at the same time its poor scalability with an increasing number of processes.

Многие задачи физики, механики и других наук связаны с решением краевых задач для дифференциальных уравнений дробного порядка. Найти точные решения данных задач чрезвычайно сложно, и в этом случае приходится искать приближенные решения. В настоящей работе предлагается математический метод приближенного решения краевой задачи для дифференциальных уравнений дробного порядка. Для производных дробного порядка мы используем определение подчиненной производной дробного порядка. Мы используем модель нейронной сети прямого распространения с одним скрытым слоем. Обучение модели осуществляется в режиме с учителем, при этом для оптимизации функции ошибки и коррекции параметров нейронной сети применяется алгоритм обратного распространения ошибки. Для иллюстрации нашего метода была разработана компьютерная программа для проведения экспериментов, в которых полученные результаты сравниваются с аналитическими расчетами. Many problems in physics, mechanics and other sciences are related to solving boundary value problems for fractional differential equations. Finding exact solutions to these problems is very difficult, and in this case, we have to look for approximate solutions. This paper proposes a mathematical method for approximate solving of a boundary value problems for fractional differential equations. For fractional derivatives we use the definition of a conformable fractional derivative. We use a feedforward neural network model with one hidden layer. The model is trained in a supervised learning mode using the backpropagation algorithm to optimize the error function and update the neural network parameters. To illustrate our method, a computer program was developed to conduct experiments in which the obtained results are compared with analytical solutions.

В работе приводится развитая модель связанных физико-химических процессов в реагирующих средах, в соответствии с которой цементный композит рассматривается как гетерогенная смесь реагентов и наполнителя с заданными концентрациями и пористостью. Фильтрация воды за счет капиллярных и термокапиллярных сил обеспечивает конвективный тепло- и массоперенос. Задача тепло- и массопереноса представляется в виде двухтемпературной нестационарной нелинейной краевой задачи теплопроводности каркаса и жидкости с источниками тепла химической природы и стоками, вызванными испарением воды. Проверка адекватности разработанной схемы решения краевой задачи теплопроводности реализуется путем решения тестовой одномерной задачи с постоянными коэффициентами. Сходимость схемы решения задачи подтверждена в процессе моделирования процесса гидратации цементного композита при различных шагах по времени и пространству. The paper presents a developed model of related physical and chemical processes in reacting medium, according to which a cement composite is considered as a heterogeneous mixture of reagents and filler with specified concentrations and porosity. Water filtration due to capillary and thermocapillary forces ensures convective heat and mass transfer. The heat and mass transfer problem is presented as a two-temperature non-stationary nonlinear boundary value problem of heat conductivity of the frame and liquid phase with heat sources of chemical nature and sinks caused by water evaporation. Verification of the adequacy of the developed scheme for solving the boundary value problem of thermal conductivity is realized by solving a test one-dimensional problem with constant coefficients. The convergence of the problem solution scheme is confirmed in the process of modeling the hydration process of a cement composite at different steps in time and space.

В работе представлен подход к построению предобуславливателя для численного решения уравнения Пуассона для существенно неоднородной среды в приложении к задачам вычислительной физики горных пород. В качестве предобуславливателя используется оператор, обратный к дискретному оператору Лапласа, но для упрощенной — слоистой — модели среды. Для обращения оператора Лапласа в этом случае используется спектральное разложение по одному из пространственных направлений и метод прогонки для серии одномерных задач по второму направлению. Такой подход к построению предобуславливателя обеспечивает независимость числа итераций от размера решаемой задачи, что подтверждается серией численных экспериментов. Важной особенностью предложенного подхода является именно использование слоистых моделей среды для построения предобуславливателя, что увеличивает скорость сходимости метода сопряженных градиентов на 10–40% в сравнении с использованием предобуславливателя, основанного на обращении оператора Лапласа для однородной среды. При этом ускорение зависит от контраста коэффициентов исходной модели, с ростом контраста повышается и эффективность предложенного подхода. The paper presents an approach to constructing a preconditioner for numerically solving the Poisson equation for a substantially inhomogeneous medium in application to problems of computational rock physics. The preconditioner is an operator inverse to the discrete Laplace operator, but for a simplified, layered model of the medium. In this case, the Laplace operator is inverted using spectral decomposition in one of the spatial directions and a sweep method for a series of one-dimensional problems in the second direction. This approach to constructing a preconditioner ensures that the number of iterations does not depend on the size of the problem being solved, which is confirmed by a series of numerical experiments. An important feature of the proposed approach is the use of layered models of the medium to construct the preconditioner, which increases the convergence rate of the conjugate gradient method by 10–40% compared to using a preconditioner based on inverting the Laplace operator for a homogeneous medium. In this case, the acceleration depends on the contrast of the coefficients of the original model; with increasing contrast, the efficiency of the proposed approach also increases.

Макроуровневые модели динамики газовзвесей часто представляют собой дифференциальные уравнения в частных производных с релаксационными слагаемыми, описывающими передачу импульса и энергии от газа к частицам и наоборот. Для ультрадисперсных частиц время релаксации намного короче, чем время, на котором рассматривается динамика среды. В работе исследуется лагранжев метод моделирования динамики газовзвесей “двухжидкостная гидродинамика сглаженных частиц” (Two-Fluid Smoothed Particle Hydrodynamics, TFSPH). TFSPH подразумевает, что каждая фаза (газ и частицы) моделируется своим набором частиц. В рамках TFSPH известны два подхода к расчету релаксационного взаимодействия (трения), которое определяется разностью скоростей между несущей и дисперсной фазами: частица-частица и частица-сетка. Ранее в вычислительных экспериментах было установлено, что для малых времен релаксации в подходе частица-частица имеет место избыточная диссипация волн, а подход частица-сетка свободен от этого недостатка. В работе впервые дано объяснение этому явлению средствами вычислительной математики. To simulate the dynamics of gas-dust mixtures at the macroscale level, it is necessary to solve numerically a gas-dynamic equations with relaxation terms describing the transfer of momentum and energy from gas to particles and vice versa. For ultrafine particles, the time of relaxation is much shorter than the time at which the dynamics of the medium is considered. We study the Lagrangian method for gas-dust mixture simulation named Two-Fluid Smoothed Particle Hydrodynamics (TFSPH). TFSPH implies that each phase (gas and particles) is modeled by its own set of particles. Within the framework of TFSPH, two approaches to calculating the relaxation interaction (drag) are known: particle-particle and particle-mesh. Previously it was found in numerical experiments that for small relaxation times, the particle-particle approach suffers from waves overdissipation, while the particle-mesh approach is free from this drawback. We provide the first explanation of this phenomenon using computational mathematics.

Рассматривается решение уравнений Максвелла для идеальных проводников в однородной среде с помощью поверхностного интегрального уравнения электрического поля, дополненного уравнением непрерывности заряда. Полученная система уравнений дискретизируется с помощью метода Галеркина. Приводится способ построения численной схемы, позволяющей проводить расчет импеданса системы питания печатной платы совместно с размещенными на ней сосредоточенными электрическими компонентами, а также рассмотрены граничные условия, характерные для микроэлектронных устройств. Полученная в ходе дискретизации система линейных уравнений решается итерационно, при этом для ускорения умножения плотных блоков матрицы на вектор неизвестных применяется быстрый метод мультиполей. Верификация методики проведена на примере имеющей аналитическое решение задачи рассеяния монохроматической волны на сфере, а также на задаче моделирования импеданса печатной платы, где проведено сравнение с результатами расчетов с использованием коммерческого программного обеспечения. The solution of Maxwell’s equations for perfect electric conductors in a homogeneous medium is considered using an augmented electric field integral equation (A-EFIE). The resulting system of equations is discretized using the Galerkin method. The method for constructing a numerical scheme is presented, which allows calculating the impedance of a printed circuit board power system together with the lumped electrical components placed on it. Additionally, boundary conditions typical for microelectronic devices are considered. The system of linear equations obtained through discretization is solved iteratively, with the fast multipole method used to accelerate the matrix-vector product. The methodology was verified using the example of scattering of a monochromatic wave by a sphere, for which an analytical solution exists, as well as through the impedance modeling of a printed circuit board, where the results were compared with those obtained from commercial software.

This paper presents a new version of the terrestrial carbon cycle module for the INM RAS Earth system model family (INM-CM). Its main difference from the previous one is a more detailed account of anthropogenic impacts on land ecosystems. This article describes a new land use database for INM-CM that takes into account changes in the spatial distribution of vegetation from 1850 to 2100. Harvesting in cultivated areas is also accommodated in the new version of the model. Numerical simulations are performed with the original and modified versions of the terrestrial carbon cycle module, covering both the historical period (1850–2014) and the possible future scenario (2015–2100). Estimates of the change in global land carbon stocks compared to the end of the pre-industrial period are obtained. The paper also proposes an effective way to prepare the initial state of carbon pools using the terrestrial carbon cycle module in stand-alone mode. В данной работе представлена новая версия модуля углеродного цикла суши для семейства Моделей Земной системы ИВМ РАН. Ее основное отличие от предыдущей заключается в более подробном учете антропогенного влияния на наземные экосистемы. В статье описывается подготовка новой базы данных землепользования, иллюстрирующей изменения в пространственном распределении растительности с 1850 по 2100 годы. Кроме того, новая версия модели учитывает сбор урожая на возделываемых территориях. С исходной и модифицированной версиями модуля наземного углеродного цикла проведены вычислительные эксперименты, охватывающие как исторический период (1850–2014), так и возможный сценарий будущего (2015–2100). Получены оценки изменения глобальных запасов углерода суши по сравнению с концом прединдустриального периода. Также в работе предложен экономичный способ формирования начального состояния для углеродных пулов с помощью автономной версии модуля наземного углеродного цикла.