We propose algorithms for efficient time integration of large systems of oscillatory second order ordinary differential equations (ODEs) whose solution can be expressed in terms of trigonometric matrix functions. Our algorithms are based on a residual notion for second order ODEs, which allows to extend the ``residual-time restarting'' Krylov subspace framework -- which was recently introduced for exponential and $\varphi$-functions occurring in time integration of first order ODEs -- to our setting. We then show that the computational cost can be further reduced in many cases by using our restarting in the Gautschi cosine scheme. We analyze residual convergence in terms of Faber and Chebyshev series and supplement these theoretical results by numerical experiments illustrating the efficiency of the proposed methods.
The main purpose of the paper is to present some powerful data on the advantage of the rational approximation procedure based on Hermite-Pad\'e polynomials over the Pad\'e approximation procedure. The first part of the paper is devoted to some numerical examples in this direction. The second part will be devoted to some theoretical results. In particular, we demonstrate our ideas about the advantage of rational Hermite-Pad\'e approximants over Pad\'e approximants analyzing the analytical structure of the frequency function $\nu$ of the free Van der Pol equation.
Rational approximation recently emerged as an efficient numerical tool for the solution of exterior wave propagation problems. Currently, this technique is limited to wave media which are invariant along the main propagation direction. We propose a new model order reduction-based approach for compressing unbounded waveguides with layered inclusions. It is based on the solution of a nonlinear rational least squares problem using the RKFIT method. We show that approximants can be converted into an accurate finite difference representation within a rational Krylov framework. Numerical experiments indicate that RKFIT computes more accurate grids than previous analytic approaches and even works in the presence of pronounced scattering resonances. Spectral adaptation effects allow for finite difference grids with dimensions near or even below the Nyquist limit.
An efficient Krylov subspace algorithm for computing actions of the \varphi matrix function for large matrices is proposed. This matrix function is widely used in exponential time integration, Markov chains, and network analysis and many other applications. Our algorithm is based on a reliable residual based stopping criterion and a new efficient restarting procedure. We analyze residual convergence and prove, for matrices with numerical range in the stable complex half-plane, that the restarted method is guaranteed to converge for any Krylov subspace dimension. Numerical tests demonstrate efficiency of our approach for solving large scale evolution problems resulting from discretized in space time-dependent PDEs, in particular, diffusion and convection-diffusion problems.
In this paper a new restarting method for Krylov subspace matrix exponential evaluations is proposed. Since our restarting technique essentially employs the residual, some convergence results for the residual are given. We also discuss how the restart length can be adjusted after each restart cycle, which leads to an adaptive restarting procedure. Numerical tests arepresented to compare our restarting with three other restarting methods. Some of the algorithms described in this paper are a part of the Octave/Matlab package expmARPACK available at http://team.kiam.ru/botchev/expm/ .
An efficient Krylov subspace algorithm for computing actions of the φ matrix function for large matrices is proposed. This matrix function is widely used in exponential time integration, Markov chains and network analysis and many other applications. Our algorithm is based on a reliable residual based stopping criterion and a new efficient restarting procedure. For matrices with numerical range in the stable complex half plane, we analyze residual convergence and prove that the restarted method is guaranteed to converge for any Krylov subspace dimension. Numerical tests demonstrate efficiency of our approach for solving large scale evolution problems resulting from discretized in space time-dependent PDEs, in particular, diffusion and convection-diffusion problems.
PreviousNext No AccessSEG Technical Program Expanded Abstracts 2020Quality control of ultra-deep resistivity imaging using fast 3D electromagnetic modelingAuthors: Sofia DavydychevaVladimir DruskinLeonid KnizhnermanMichael RabinovichSofia Davydycheva3DEM HoldingSearch for more papers by this author, Vladimir DruskinWorcester Polytechnic InstituteSearch for more papers by this author, Leonid Knizhnerman3DEM HoldingSearch for more papers by this author, and Michael RabinovichBPSearch for more papers by this authorhttps://doi.org/10.1190/segam2020-3425809.1 SectionsSupplemental MaterialAboutPDF/ePub ToolsAdd to favoritesDownload CitationsTrack CitationsPermissions ShareFacebookTwitterLinked InRedditEmail AbstractFull 3D electromagnetic modeling based on optimal multiscale Lebedev’s discretization of Maxwell’s equations with general anisotropy is used for reliable modeling and quality control of ultra-deep electromagnetic (EM) measurements. With this fast and robust full 3D anisotropic modeling software we calculate the tool responses in provided 2D/3D models and compare them to the respective 1D modeling results or to the real data whenever available. The observed differences are analyzed for the resistivity and for the directional curves. We often observe significant differences between 3D and local 1D responses approximating real data since the ultra-deep tool length is typically comparable with size of detected 3D anomalies and cannot be neglected especially when transmitters and receivers are in different structures. Using the fast and accurate 3D modeling we also study various effects as 3D features, anisotropy, changing dips, etc. The accurate full 3D modeling with general anisotropy is critically important for quality control of the standard deep resistivity image and distance-to-bed answer products based on the local 1D inversion. We illustrate our modeling and QC approach on the ultra-deep resistivity data acquired in a horizontal well in the North Sea. 3D modeling & 2D/3D inversion can significantly improve answer products for the ultra-deep resistivity tools.Tuesday, October 13, 2020Session Start Time: 1:50 PMPresentation Time: 2:40 PMLocation: 351DPresentation Type: OralKeywords: 3D, electromagnetics, finite difference, borehole measurements, logging while drillingPermalink: https://doi.org/10.1190/segam2020-3425809.1FiguresReferencesRelatedDetails SEG Technical Program Expanded Abstracts 2020ISSN (print):1052-3812 ISSN (online):1949-4645Copyright: 2020 Pages: 3887 publication data© 2020 Published in electronic format with permission by the Society of Exploration GeophysicistsPublisher:Society of Exploration Geophysicists HistoryPublished Online: 30 Sep 2020 CITATION INFORMATION Sofia Davydycheva, Vladimir Druskin, Leonid Knizhnerman, and Michael Rabinovich, (2020), "Quality control of ultra-deep resistivity imaging using fast 3D electromagnetic modeling," SEG Technical Program Expanded Abstracts : 380-384. https://doi.org/10.1190/segam2020-3425809.1 Plain-Language Summary Keywords3Delectromagneticsfinite differenceborehole measurementslogging while drillingPDF DownloadLoading ...
A new construction of an absorbing boundary condition for indefinite Helmholtz problems on unbounded domains is presented. This construction is based on a near-best uniform rational interpolant of the inverse square root function on the union of a negative and positive real interval, designed with the help of a classical result by Zolotarev. Using Krein's interpretation of a Stieltjes continued fraction, this interpolant can be converted into a three-term finite difference discretization of a perfectly matched layer (PML) which converges exponentially fast in the number of grid points. The convergence rate is asymptotically optimal for both propagative and evanescent wave modes. Several numerical experiments and illustrations are included.
Rational Arnoldi is a powerful method for approximating functions of large sparse matrices times a vector. The selection of asymptotically optimal parameters for this method is crucial for its fast convergence. We present and investigate a novel strategy for the automated parameter selection when the function to be approximated is of Cauchy–Stieltjes (or Markov) type, such as the matrix square root or the logarithm. The performance of this approach is demonstrated by numerical examples involving symmetric and nonsymmetric matrices. These examples suggest that our black-box method performs at least as well, and typically better, as the standard rational Arnoldi method with parameters being manually optimized for a given matrix.
PreviousNext No AccessSEG Technical Program Expanded Abstracts 2012Schemes for improving efficiency of pixel-based inversion algorithms for electromagnetic logging-while-drilling measurementsAuthors: Y. LinA. AbubakarT. M. HabashyG. PanM. LiV. DruskinL. KnizhnermanY. LinSchlumberger-Doll Research, USASearch for more papers by this author, A. AbubakarSchlumberger-Doll Research, USASearch for more papers by this author, T. M. HabashySchlumberger-Doll Research, USASearch for more papers by this author, G. PanSchlumberger-Doll Research, USASearch for more papers by this author, M. LiSchlumberger-Doll Research, USASearch for more papers by this author, V. DruskinSchlumberger-Doll Research, USASearch for more papers by this author, and L. KnizhnermanCentral Geophysical Expedition, RussiaSearch for more papers by this authorhttps://doi.org/10.1190/segam2012-0208.1 SectionsAboutPDF/ePub ToolsAdd to favoritesDownload CitationsTrack CitationsPermissions ShareFacebookTwitterLinked InRedditEmail Abstract We present an application of the two-and-half dimensional (2.5D) multiplicative-regularized Gauss-Newton inversion method for the interpretation of the directional resistivity measurements [electromagnetic logging-while-drilling (LWD) measurements]. In this work, an inversion algorithm is employed for obtaining a two-dimensional (2D) resistivity distribution (pixel-based) of the subsurface. Several modifications on both forward and inversion algorithms have been implemented to minimize the computational time and memory usage as well as to improve the accuracy of the algorithms. Numerical examples demonstrate the feasibility of using this 2.5D pixel-based inversion algorithm for electromagnetic LWD data. Permalink: https://doi.org/10.1190/segam2012-0208.1FiguresReferencesRelatedDetailsCited byTreatment of singularity in the method of boundary integral equations for 2.5D electromagnetic modelingGleb Dyatlov, Dmitry Kushnir, and Yuliy Dashevsky16 February 2017 | GEOPHYSICS, Vol. 82, No. 2Efficient 2.5D electromagnetic modeling using boundary integral equationsGleb Dyatlov, Elizaveta Onegova, and Yuliy Dashevsky27 March 2015 | GEOPHYSICS, Vol. 80, No. 3 SEG Technical Program Expanded Abstracts 2012ISSN (print):1052-3812 ISSN (online):1949-4645Copyright: 2012 Pages: 4609 Publisher:Society of Exploration Geophysicists HistoryPublished Online: 25 Oct 2012 CITATION INFORMATION Y. Lin, A. Abubakar, T. M. Habashy, G. Pan, M. Li, V. Druskin, and L. Knizhnerman, (2012), "Schemes for improving efficiency of pixel-based inversion algorithms for electromagnetic logging-while-drilling measurements," SEG Technical Program Expanded Abstracts : 1-5. https://doi.org/10.1190/segam2012-0208.1 Plain-Language Summary PDF DownloadLoading ...
For large scale problems, an effective approach for solving the algebraic Lyapunov equation consists of projecting the problem onto a significantly smaller space and then solving the reduced order matrix equation. Although Krylov subspaces have been used for a long time, only more recent developments have shown that rational Krylov subspaces can be a competitive alternative to the classical and very popular alternating direction implicit (ADI) recurrence. In this paper we develop a convergence analysis of the rational Krylov subspace method (RKSM) based on the Kronecker product formulation and on potential theory. Moreover, we propose new enlightening relations between this approach and the ADI method. Our results provide solid theoretical ground for recent numerical evidence of the superiority of RKSM over ADI when the involved parameters cannot be computed optimally, as is the case in many practical application problems.
The extended Krylov subspace method has recently arisen as a competitive method for solving large-scale Lyapunov equations. Using the theoretical framework of orthogonal rational functions, in this paper we provide a general a priori error estimate when the known term has rank-one. Special cases, such as symmetric coefficient matrix, are also treated. Numerical experiments confirm the proved theoretical assertions.
Time-domain problems for controlled-source electromagnetic exploration require accurate discretization of the solution for multiple spacial and temporal scales. Therefore, forward simulation using conventional computational methods becomes computationally expensive, even without accounting for induced-polarization (IP) effects. These effects create another complication caused by the presence of a convolution integral in the time-domain Maxwell system. We suggested a novel, fast, and robust algorithm to solve the 3D time-domain electromagnetic (EM) problems that can be considered as a generalization of the spectral Lanczos decomposition method. The new method also allowed us to incorporate the IP effects without significant cost increase. The discretized large-scale Maxwell system was projected onto a small subspace consisting of the Laplace-domain solutions (the so-called parameter-dependent Krylov subspace) for an optimally chosen set of Laplace parameters. The projected system preserved stability and passivity of the original problem. Moreover, our approach (even without the IP effects) yielded an optimal solution within a wide class of computational algorithms that included the conventional time-domain finite-difference, discrete Fourier transform and spectral Lanczos decomposition methods. Numerical examples for the controlled-source EM problem showed that the new algorithm produces accurate solutions on time intervals spanning from milliseconds to hundreds of seconds with the cost of (at most) 60 time steps of the implicit finite-difference time domain scheme. This showed significant improvement even compared with results for nonpolarized media reported in recent literature. Additionally, the new algorithm had the unique capability to accurately handle large-scale 3D models, including the IP effects.
We developed a 2.5D finite‐difference (FD) code for modeling EM tool responses for 2D formation conductivity distributions and 3D well trajectories. The code is primarily developed for well placement and for 3D formation evaluation applications to interpret responses of the new generation LWD deep directional EM tools and tensor induction tools in high angle and horizontal (HA/HZ) wells. Frequency‐domain Maxwell equations are solved using electric field formulation. The problem is discretized using a combination of staggered uniform FD grids in the tool region based on skin depth and tool orientation, optimal grids outside the tool region in the x‐z plane, and a spatial Fourier transform in the y‐direction using Zoltarov discretization. Material averaging is applied in conjunction with optimal grids, resulting in a small problem size, allowing use of a direct LUD solver. We present the directional tool responses in realistic 2D scenarios such as an unconformity structure with internal layering. In addition, we also model the coupling of faults and nearby boundaries and present how that affect the 1D real‐time inversion ignoring the fault.
The modeling of the controlled-source electromagnetic (CSEM) and single-well and crosswell electromagnetic (EM) configurations requires fine gridding to take into account the 3D nature of the geometries encountered in these applications that include geological structures with complicated shapes and exhibiting large variations in conductivities such as the sea-floor bathymetry, the land topography, and targets with complex geometries and large contrasts in conductivities. Such problems significantly increase the computational cost of the conventional finite-difference (FD) approaches mainly due to the large condition numbers of the corresponding linear systems. To handle these problems, we employ a volume integral equation (IE) approach to arrive at an effective preconditioning operator for our FD solver. We refer to this new hybrid algorithm as the finite-difference integral equation method (FDIE). This FDIE preconditioning operator is divergence free and is based on a magnetic field formulation. Similar to the Lippman-Schwinger IE method, this scheme allows us to use a background elimination approach to reduce the computational domain, resulting in a smaller size stiffness matrix. Furthermore, it yields a linear system whose condition number is close to that of the conventional Lippman-Schwinger IE approach, significantly reducing the condition number of the stiffness matrix of the FD solver. Moreover, the FD framework allows us to substitute convolution operations by the inversion of banded matrices, which significantly reduces the computational cost per iteration of the hybrid method compared to the standard IE approaches. Also, well-established FD homogenization and optimal gridding algorithms make the FDIE more appropriate for the discretization of strongly inhomogeneous media. Some numerical studies are presented to illustrate the accuracy and effectiveness of the presented solver for CSEM, single-well, and crosswell EM applications.
Rational Arnoldi is a powerful method for approximating functions of large sparse matrices times a vector. The selection of asymptotically optimal parameters for this method is crucial for its fast convergence. We present a heuristic for the automated pole selection when the function to be approximated is of Markov type, such as the matrix square root. The performance of this approach is demonstrated at several numerical examples. (© 2011 Wiley‐VCH Verlag GmbH & Co. KGaA, Weinheim)