# Result

**Some snippets from your converted document:**

Суперкомпьютерные дни в России 2016 // Russian Supercomputing Days 2016 // RussianSCDays.org A parallel algorithm for 3D modeling of monochromatic acoustic field by the method of integral equations ∗ M.S. Malovichko1, N.I. Khokhlov1, N.B. Yavich1, M.S. Zhdanov2 Moscow Institute of Physics and Technology1, The University of Utah2 We present a parallel algorithm for solution of the three-dimensional Helmholtz equa- tion in the frequency domain by the method of volume integral equations. The al- gorithm is applied to seismic forward modeling. The method of integral equations reduces the size of the problem by dividing the geologic model into the anomalous and background parts, but leads to a dense system matrix. A tolerable memory con- sumption and numerical complexity were achieved by applying an iterative solver, accompanied by an effective matrix-vector multiplication operation, based on the fast Fourier transform. We used OpenMP to speed up the matrix-vector multiplication, while MPI was used to speed up the equation system solver, and also for parallelizing across multiple sources. Practical examples and efficiency tests are presented. Keywords: Fredholm integral equations, acoustics, seismics, MPI, OpenMP. 1. Introduction The method of integral equations (IE) is a well-known wavefield modeling technique [2,9]. It has a number of attractive properties. It requires discretization for the anomalous volume only. It drastically reduces the size of a problem in many applications. The IE formulation does not require boundary conditions. It also can be naturally used to compute the Fréchet derivatives and for this reason it is very attractive for inverse problems. The IE method is based on the discretization of the original integral equation and results in a complex-valued dense matrix. In a straightforward implementation, the computational burden becomes prohibitive for large three-dimensional problems, encountered in seismology. Yet some studies, many of them for the two-dimensional solution, have been reported [4–6,8,11]. Recently, the full method of IE, accompanied with an iterative solver, was applied to the visco-acoustic three-dimensional modeling [1]. The authors of [1] have implemented Green’s integral operator for the free space background model via the Fourier transform (this idea itself is quite old). An effective matrix-free implementation allowed them application of the BiCGStab method to resulting system of linear equations with the cost per iteration of O(N logN ), where N is the number of model cells. They showed, that a realistically large acoustic model can be simulated almost as effectively with the IE method, as it is for finite-different schemes. In this study, we design a parallel version of the solver for the half-space host medium, and study the efficiency of parallelization. 2. Problem formulation We consider a 3D model consisting of a half-space host medium and an anomalous volume D, confined within the lower half-space. The host medium has piece-wise constant density and acoustic velocity: ρ0 and c0 for the upper half-space; ρb and cb for the lower half space. The anomalous volume has the same background density ρb , and arbitrary distribution of velocity c = c(r), where r is the position vector. Let ω be the circular frequency. Assuming a lossless medium, we define the anomalous velocity ca = c − cb , wavenumber k = ω/c, background wavenumber kb = ω/cb and parameter χ = 1/c2 − 1/c2b . The pressure response p to a point ∗ This research was supported by the Russian Science Foundation, project No. 16-11-10188. 4 Суперкомпьютерные дни в России 2016 // Russian Supercomputing Days 2016 // RussianSCDays.org source at a given frequency ω satisfies to the following Lippmann-Schwinger equation: p − G[ω 2 χp] = pb , (1) where pb is the background field, G is the integral operator, defined for any scalar field f , Z G[f ] = g(r|r0 )f (r0 )dV 0 , (2) D where g is the background Green’s function. Let the anomalous volume be covered by N = Nx Ny Nz cubical cells. After discretization we receive the following matrix equation: Au = b, (3) where u = (p1 ..pN ) is the vector of N unknown values approximating p in the center of each cell, b = (pb,1 ..pb,N ) is known values of the background field in each cell, A is the scattering matrix, A = (I − GX),

**Recently converted files (publicly available):**