Isaac Harari*
Paul E. Barbone and Joshua M. Montgomery
Department of Solid Mechanics, Materials & Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace & Mechanical Engineering, Boston University
110 Cummington Street, Boston, MA 02215, U.S.A.
We develop formulations for finite element computation of exterior acoustics problems. A prominent feature of the formulations is the lack of integration over the unbounded domain, simplifying the task of discretization and potentially leading to numerous additional benefits. These formulations provide a suitable basis for hybrid asymptotic-numerical methods in scattering, non-reflecting boundary conditions and infinite elements.
Keywords: finite element methods, unbounded domains, acoustics, Lagrange multipliers, absorbing boundary conditions, infinite elements.
*This research was performed in part while the author was visiting the Department of Aerospace & Mechanical Engineering, Boston University.
Isaac Harari
Karl Grosh
Thomas J.R. Hughes, Manish Malhotra, Peter M. Pinsky, James R. Stewart
Lonny L. Thompson
Tel-Aviv University
University of Michigan
Stanford University
Clemson University
The study of structural acoustics involves modeling acoustic radiation and scattering, primarily in exterior regions, coupled with elastic and structural wave propagation. This paper reviews recent progress in finite element analysis that renders computation a practical tool for solving problems of structural acoustics.
The cost-effectiveness of finite element methods is composed of several ingredients. Boundary-value problems in unbounded domains are inappropriate for direct discretization. Employing DtN methodology yields an equivalent problem that is suitable for finite element analysis by posing impedance relations at an artificial exterior boundary. Well-posedness of the resulting continuous formulations is discussed, leading to simple guidelines for practical implementation and verifying that DtN boundary conditions provide a suitable basis for computation.
Approximation by Galerkin finite element methods results in spurious dispersion, degrading with reduced wave resolution. Accuracy is improved by Galerkin/least-squares and related technologies on the basis of detailed examinations of discrete errors in simplified settings, relaxing wave-resolution requirements. This methodology is applied to time-harmonic problems of acoustics and coupled problems of structural acoustics. Space-time finite element methods based on time-discontinuous Galerkin/least-squares are derived for transient problems of structural acoustics. Numerical results validate the superior performance of Galerkin/least-squares finite elements for problems of structural acoustics.
A comparative study of the cost of computation demonstrates that Galerkin/least-squares finite element methods are economically competitive with boundary element methods, the prevailing numerical approach to exterior problems of acoustics. Efficient iterative methods are derived for solving the large-scale matrix problems that arise in structural acoustics computation of realistic configurations at high wavenumbers. An a posteriori error estimator and adaptive strategy are developed for time-harmonic acoustic problems and the role of adaptivity in reducing the cost of computation is addressed.
Isaac Harari and Zion Shohet
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
Problems in unbounded domains can be solved using domain-based computation by introducing an artificial boundary, and then selecting appropriate boundary conditions. The DtN method, which specifies such boundary conditions, is investigated in this work for wave problems in elastic solids. The DtN method defines an exact relation between the displacement field and its normal and tangential tractions on an artificial boundary. This relation is expressed in terms of an infinite series. The DtN boundary conditions are shown to be non-reflective, thus uniqueness of the solution is guaranteed. For practical purposes the full DtN operator is truncated. The truncated DtN operator fails to completely inhibit reflections of higher modes, resulting in loss of uniqueness at characteristic wave numbers of higher harmonics. Guidelines for determining a sufficient number of terms in the truncated operator to retain uniqueness of the solution at any given wave number are derived. The validity of these guidelines is examined and verified by numerical examples. Local DtN boundary conditions are also investigated, and it is shown that local boundary conditions guarantee uniqueness of the solution for all wave numbers, regardless of the number of terms in the operator. This property is used here to modify the truncated DtN operator and to enhance its capability to retain uniqueness of solutions. A modified DtN operator, combining the truncated operator with the local one, is introduced. The modified DtN operator is shown to retain uniqueness of solutions regardless of the number of terms and regardless of the wave number.
Isaac Harari
Paul E. Barbone
Michael Slavutin and Rami Shalom
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace & Mechanical Engineering, Boston University
110 Cummington Street, Boston, MA 02215, U.S.A.
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
A novel approach to the development of infinite element formulations for exterior problems of time-harmonic acoustics is presented. This approach is based on a functional which provides a general framework for domain-based computation of exterior problems. Special cases include non-reflecting boundary conditions (such as the DtN method). A prominent feature of this formulation is the lack of integration over the unbounded domain, simplifying the task of discretization. The original formulation is generalized to account for derivative discontinuities across infinite element boundaries, typical of standard infinite element approximations. Continuity between finite elements and infinite elements is enforced weakly, precluding compatibility requirements. Various infinite element approximations for two-dimensional configurations with circular interfaces are presented. Implementation requirements are relatively simple. Numerical results demonstrate the good performance of this scheme.
Isaac Harari, Rami Shalom
Paul E. Barbone
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace & Mechanical Engineering, Boston University
110 Cummington Street, Boston, MA 02215, U.S.A.
Higher-order approximations are implemented in a novel approach to infinite element formulations for exterior problems of time-harmonic acoustics. This approach is based on a functional which provides a general framework for domain-based computation of exterior problems. Two prominent features of this formulation simplify the task of discretization: the lack of integration over the unbounded domain, and weak enforcement of continuity between finite elements and infinite elements. Consequently, the infinite elements mesh the interface only and need not match the finite elements on the interface. Various infinite element approximations for two-dimensional configurations with circular interfaces are presented. Numerical results demonstrate the good performance of this scheme.
Isaac Harari
Igor Patlashenko
Dan Givoli
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
Asher Space Research Institute,
Technion---Israel Institute of Technology
32000 Haifa, Israel
Department of Aerospace Engineering,
Technion---Israel Institute of Technology
32000 Haifa, Israel
Dirichlet-to-Neumann (DtN) boundary conditions for unbounded wave guides in two and three dimensions are derived and analyzed, defining problems that are suitable for finite element analysis. In the most general cases considered wave numbers may vary in arbitrary cross sections. The full DtN operator, in the form of an infinite series, is exact. Nonunique solutions may occur when this operator is truncated. Simple criteria for the number of terms in the truncated operator that guarantee unique solutions are presented. A simple modification of the truncated operator leads to uniqueness for any number of terms. Numerical results validate the performance of DtN formulations for wave guides and confirm the criteria for uniqueness.
Isaac Harari and Shmuel Haham
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
In this work we develop finite element methods for exterior problems of time-harmonic elastic waves. We employ DtN exterior boundary conditions on an artificial boundary to obtain a well-posed equivalent problem in a bounded region that is suitable for domain-based computation.
Galerkin, Galerkin/least-squares, and Galerkin/gradient least-squares finite element methods are presented. In acoustic wave problems the Galerkin/least-squares method counters the numerical difficulties that result from employing the Galerkin method, but proves insufficient for elastic waves. A specialization of the Galerkin/gradient least-squares method to problems of elastic waves yields good results. Numerical examples validate these conclusions.
Isaac Harari
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
A general framework for domain-based computation in unbounded regions, with particular reference to time-harmonic acoustics, is based on a functional for partitioned problems, weakly enforcing continuity of inner and outer fields. Two prominent features of this formulation simplify the task of discretization: the lack of integration over the unbounded domain, and accommodation of incompatible approximations.
The specific formulation depends on the representation of the outer field. Examples of formulations employing DtN boundary conditions are reproduced. A novel approach to infinite element formulations is presented, leading to various approximations for two-dimensional configurations with circular interfaces. Numerical results demonstrate the good performance of these schemes.
Isaac Harari, Parama Barai
Paul E. Barbone
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace & Mechanical Engineering, Boston University
110 Cummington Street, Boston, MA 02215, U.S.A.
The numerical and spectral performance of novel infinite elements for exterior problems of time-harmonic acoustics are examined. The formulation is based on a functional which provides a general framework for domain-based computation of exterior problems. Two prominent features simplify the task of discretization: the infinite elements mesh the interface only and need not match the finite elements on the interface.
Various infinite element approximations for two-dimensional configurations with circular interfaces are reviewed. Numerical results demonstrate the good performance of these schemes. A simple study points to the proper interpretation of spectral results for the formulation. The spectral properties of these infinite elements are examined with a view to the representation of physics and efficient numerical solution.
Isaac Harari
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
A reference solution for cylindrical waves is obtained by WKB approximation. Finite element representations of both cylindrical and spherical waves satisfy the resolution-dependent plane wave dispersion relation, and approximate the desired decay rates independent of mesh resolution.
Paul E. Barbone
Isaac Harari
Department of Aerospace & Mechanical Engineering, Boston University
110 Cummington Street, Boston, MA 02215, U.S.A.
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
We examine the problem of finding the H1 projection onto a finite element space of an unknown field satisfying a specified boundary value problem. Solving the projection problem typically requires knowing the exact solution. We circumvent this issue and obtain a Petrov-Galerkin formulation which achieves H1 optimality. Requiring weighting functions to be defined locally on the element level permits only approximate H1 optimality in multi-dimensional configurations. We present simple guidelines for obtaining the weight functions, and include a numerical example for the Helmholtz equation.
Isaac Harari and Michael Slavutin
Eli Turkel
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
School of Mathematical Sciences,
Tel-Aviv University
69978 Ramat Aviv, Israel
A symmetric PML formulation that is suitable for finite element computation of time-harmonic acoustic waves in exterior domains is analyzed. Dispersion analysis displays the dependence of the discrete representation of the PML parameters on mesh refinement. Stabilization by modification of the coefficients is employed to improve PML performance, in conjunction with standard stabilized finite elements in the Helmholtz region. Numerical results validate the good performance of this finite element PML approach.
Isaac Harari and Parama Barai
Paul E. Barbone
Michael Slavutin
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace & Mechanical Engineering, Boston University
110 Cummington Street, Boston, MA 02215, U.S.A.
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
Three-dimensional infinite elements for exterior problems of time-harmonic acoustics are developed. The infinite elements mesh only the outer boundary of the finite element domain and need not match the finite elements on the interface. A four-noded infinite element, based on separation of variables in spherical coordinates, is presented. Singular behavior of associated Legendre functions at the poles is circumvented. Numerical results validate the good performance of this approach.
Isaac Harari
Leoplodo P. Franca and Saulo P. Oliveira
Department of Solid Mechanics, Materials and Structures,
Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mathematics, University of Colorado at Denver
Denver, CO 80217-3364, U.S.A.
The dependence of the computation of advective-diffusive transport phenomena on the orientation of the mesh with respect to the flow direction is analyzed. Poor performance of the classical Galerkin finite element method in the convection-dominated regime is alleviated by stabilization. We propose definitions of the stability parameter that rationally incorporate the flow direction. Numerical tests compare the performance of the proposed methods with established techniques.
Charbel Farhat
Isaac Harari
Leoplodo P. Franca
Department of Aerospace Engineering Sciences,
University of Colorado at Boulder
Boulder, CO 80309-0429, U.S.A.
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mathematics, University of Colorado at Denver
Denver, CO 80217-3364, U.S.A.
We propose a finite element based discretization method in which the standard polynomial field is enriched within each element by a nonconforming field that is added to it. The enrichment contains free-space solutions of the homogeneous differential equation that are not represented by the underlying polynomial field. Continuity of the enrichment across element interfaces is enforced weakly by Lagrange multipliers. In this manner, we expect to attain high coarse-mesh accuracy without significant degradation of conditioning, assuring good performance of the computation at any mesh resolution. Examples of application to acoustics and advection-diffusion are presented.
Isaac Harari
Carnot L. Nogueira
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace Engineering Sciences,
University of Colorado at Boulder
Boulder, CO 80309-0429, U.S.A.
Galerkin/least squares (GLS) modification improves the performance of finite element computation of time-harmonic acoustics at high wave numbers. The design of the GLS resolution-dependent method parameter for two-dimensional computation in previous work is based on dispersion analysis of one-dimensional and square bilinear elements. We analyze the dispersion of linear triangular finite elements, and define method parameters that eliminate dispersion on a hexagonal patch. Numerical tests compare the performance of the proposed method with established techniques on structured and unstructured triangular meshes. Based on this work, we propose a method parameter that may be used for computation with both linear triangular and bilinear quadrilateral elements.
Ali Nesliturk
Isaac Harari
Department of Applied Mathematics,
Izmir Institute of Technology
Guzelbahce koyu, Urla / Izmir, TURKIYE
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
The nearly-optimal Petrov-Galerkin (NOPG) method is employed to improve finite element computation of convection-dominated transport phenomena. The design of the NOPG method for convection-diffusion is based on consideration of the advective limit. Nonetheless, the resulting method is applicable to the entire admissible range of problem parameters. An investigation of the stability properties of this method leads to a coercivity inequality. The convergence features of the NOPG method for convection-diffusion are studied in an error analysis that is based on the stability estimates. The proposed method compares favorably to the performance of an established technique on several numerical tests.
Isaac Harari
Sérgio Frey
Leoplodo P. Franca
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Thermal Sciences and Energy Systems Group (GESTE)
Mechanical Engineering Department, UFRGS
Sarmento Leite No. 425
90050-170 Porto Alegre/RS, Brazil
Department of Mathematics, University of Colorado at Denver
Denver, CO 80217-3364, U.S.A.
In a recent paper studying finite element computation of heat transfer processes with dominant sources, for which the classical Galerkin method proves unstable, the authors conclude that Galerkin/least-squares (GLS) stabilization is insufficient while Galerkin-gradient/least-squares (GGLS) stabilization provides good results. It is the intention of this manuscript to correct these conclusions, that are based on a GLS method with a suboptimal parameter and on mislabelling a combined stabilized method as GGLS.
Charbel Farhat
Isaac Harari
Ulrich Hetmaniuk
Department of Aerospace Engineering Sciences,
University of Colorado at Boulder
Boulder, CO 80309-0429, U.S.A.
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace Engineering Sciences,
University of Colorado at Boulder
Boulder, CO 80309-0429, U.S.A.
We present a discontinuous Galerkin method for the solution of the Helmholtz equation in the mid-frequency regime. Our approach is based on the discontinuous enrichment method in which the standard polynomial field is enriched within each finite element by a nonconforming field that contains free-space solutions of the homogeneous partial differential equation to be solved. Hence, for the Helmholtz equation, the enrichment field is chosen here as the superposition of plane waves. We enforce a weak continuity of these plane waves across the element interfaces by suitable Lagrange multipliers. Preliminary results obtained for two-dimensional model problems discretized by uniform meshes reveal that the proposed discontinuous Galerkin method enables the development of elements that are far more competitive than both the standard linear and the standard quadratic Galerkin elements for the discretization of Helmholtz problems.
Isaac Harari
Charbel Farhat and Ulrich Hetmaniuk
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace Engineering Sciences,
University of Colorado at Boulder
Boulder, CO 80309-0429, U.S.A.
We analyze the dispersion properties of elements obtained by a discontinuous Galerkin method with Lagrange multipliers. The dispersion analysis of these elements presents a challenge in that the Lagrange multiplier degrees of freedom are directional, and hence an unbounded mesh is made up of more than one repeating pattern. Two approaches to overcome this difficulty are presented. The similarity in the two sets of results offers mutual validation of the two approaches.
Charbel Farhat
Isaac Harari
Ulrich Hetmaniuk
Department of Aerospace Engineering Sciences,
University of Colorado at Boulder
Boulder, CO 80309-0429, U.S.A.
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace Engineering Sciences,
University of Colorado at Boulder
Boulder, CO 80309-0429, U.S.A.
Computation naturally separates scales of a problem according to the mesh size. A variety of improved numerical methods are described by multiscale considerations, differing in the treatment of the unresolved, fine scales. The discontinuous enrichment method provides a unique multiscale approach to computation by employing fine scales that contain solutions of the homogeneous partial differential equation in a discontinuous framework. The method thus combines relative ease of implementation with improved numerical performance. These properties are demonstrated for both multiscale wave and transport problems, pointing to the potential of considerable savings in computational resources.
Isaac Harari
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Solutions of direct time-integration schemes that converge in time to conventional semidiscrete formulations may be polluted at small time steps by spurious spatial oscillations, along with attendant overshoot in time. This degradation is not an artifact of the time-marching scheme, but rather a property of the solution of the semidiscrete formulation itself. A critical time step for the onset of spatial oscillations is derived by analogy to singularly perturbed elliptic problems. We then propose a simple procedure of spatial stabilization to remove this pathology from implicit time-integration schemes, without affecting unconditional temporal stability. Spatially stabilized implicit time integration is free of spurious spatial oscillations at small time steps, and numerical experience points to improved temporal accuracy as well.
Isaac Harari
Rabia Djellouli
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace Engineering Sciences,
University of Colorado at Boulder
Boulder, CO 80309-0429, U.S.A.
The computation of exterior wave problems at low wave numbers can become prohibitively expensive when higher circumferential modes are significant. An analysis of the effect of wave number on scattering problems, with local absorbing boundary conditions specified on simple shapes as on-surface radiation conditions, provides guidelines for satisfactory performance. Excessive computational cost may be avoided for most practical applications.
Isaac Harari
Frédéric Magoulès
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Institut Elie Cartan de Nancy,
Université Henri Poincaré, BP 239
54506 Vandoeuvre les Nancy Cedex, FRANCE
Least-squares stabilization stands out among the numerous approaches that have been proposed for relaxing resolution requirements of Galerkin computations for acoustics, by combining substantial improvement in performance with extremely simple implementation. The Galerkin/least-squares and Galerkin-gradient/least-squares methods are quite similar for structured meshes of linear finite elements. A series of numerical tests compares the two methods for several configurations with different kinds of boundary conditions employing structured and unstructured meshes. Various definitions of the resolution-dependent stability parameters are considered, along with different definitions of the mesh size upon which they depend.
Isaac Harari
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Many of the current issues and methodologies related to finite element methods for time-harmonic acoustics are reviewed. The effective treatment of unbounded domains is a major challenge. Most prominent among the approaches that have been developed for this purpose are absorbing boundary conditions, infinite elements, and absorbing layers. Standard computational methods are unable to cope with wave phenomena at short wave lengths due to resolutions required to control dispersion and pollution errors, leading to prohibitive computational demands. Since computation naturally separates the scales of a problem according to the mesh size, multiscale considerations provide a useful framework for viewing these difficulties and developing methods to counter them. Other issues addressed are related to the efficient solution of systems of specialized algebraic equations, and inverse problems of acoustics. The tremendous progress that has been made in all of the above areas in recent years will surely continue, leading to many more exciting developments.
Isaac Harari
Radek Tezaur and Charbel Farhat
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Aerospace Engineering Sciences,
University of Colorado at Boulder
Boulder, CO 80309-0429, U.S.A.
One-dimensional analyses provide novel definitions of the Galerkin/least-squares stability parameter for quadratic interpolation. A new approach to the dispersion analysis of the Lagrange multiplier approximation in discontinuous Galerkin methods is presented. A series of computations comparing the performance of Q2 Galerkin and GLS methods with Q-8-2 DGM on large-scale problems shows superior DGM results on analogous meshes, both structured and unstructured. The degradation of the Q2 GLS stabilization on unstructured meshes may be a consequence of inadequate one-dimensional analysis used to derive the stability parameter.
Isaac Harari and Uri Albocher
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
The equation of an absorbing layer for time-harmonic elastic waves, based on the Perfectly Matched Layer (PML) concept, is formulated in a manner that is easily implemented in finite element software. In the proposed approach, the layer is viewed as an anisotropic material with continuously varying complex material properties. The effect of the PML parameters on its discrete representation is investigated through dispersion analysis. Guidelines for proper selection of the PML parameters are presented. The formulation, with the parameters selected according to the guidelines proposed, performs well in computation.
Slava Krylov, Isaac Harari, and Yaron Cohen
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Electrostatically actuated micro structures are inherently nonlinear and can become unstable. Pull-in instability is encountered as a basic static instability mechanism. We demonstrate that the parametric excitation of a micro structure by periodic (AC) voltages may have a stabilizing effect and permits an increase of the steady (DC) component of the actuation voltage beyond the pull-in value. An elastic string as well as a cantilever beam are considered in order to illustrate the influence of fast-scale excitation on the slow-scale behavior. The main conclusions about the stability are drawn using the simplest model of a parametrically excited system described by the Mathieu and Hill's equations. Theoretical results are verified by numerical analysis of micro structure subject to nonlinear electrostatic forces and performed by using Galerkin decomposition with undamped linear modes as base functions. The parametric stabilization of a cantilever beam is demonstrated experimentally.
Robert C. Reiner, Jr. and Rabia Djellouli
Isaac Harari
Department of Mathematics,
California State University Northridge
Northridge, CA 91330-8313, U.S.A.
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
The computation of exterior wave problems at low wave numbers can become prohibitively expensive when higher circumferential modes are significant, yet local absorbing boundary conditions specified on simple shapes provide satisfactory performance for most scattering problems. An analysis of the effect of wave number and eccentricity on scattering problems, with local absorbing boundary conditions specified on an ellipse and a prolate spheroid as on-surface radiation conditions, shows that at moderate values of eccentricity, good performance extends to relatively low wave numbers. Excessive computational cost may be avoided for many practical applications.
Robert C. Reiner, Jr. and Rabia Djellouli
Isaac Harari
Department of Mathematics,
California State University Northridge
Northridge, CA 91330-8313, U.S.A.
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
A mathematical and numerical analysis is performed to assess the performance of second order Bayliss-Gunzburger-Turkel condition (BGT2) when applied for solving low frequency acoustic scattering problems in the case of elliptical-shaped scatterers. This investigation suggests that BGT2 retains an acceptable level of accuracy for relatively low wavenumber. A damping effect is incorporated to the BGT2 condition in order to extend the range of satisfactory performance. This damping procedure consists in adding only a constant imaginary part to the wavenumber. The numerical results indicate that the modified version of BGT2 extends the range of satisfactory performance by improving the level of accuracy by up to two orders of magnitude. Guidelines on the appropriate choice of the damping coefficient are provided.
Slava Krylov, Isaac Harari, and David Gadasi
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
In multiphysics problems, a solid body is in interaction with various three-dimensional fields that generate complex patterns of rapidly varying distributed loading on the solid. Since three-dimensional computation requires excessive resources, methods of reduction to structural models are traditionally exploited in mechanics for the analysis of slender bodies. Although such procedures are well-established, the reduction of loads is often performed in an ad hoc manner which is not sufficient for many coupled problems. In the present work we develop rigorous Structural Reduction (SR) procedures by using a variational framework to consistently convert three-dimensional interface data to the form required by structural representations. The approach is illustrated using the Euler-Bernoulli and Timoshenko beam theories. Some of the loading terms and boundary conditions of the four resulting structural problems (namely, tension, torsion and two bending problems) which are formulated in terms of the original three-dimensional problem could not be derived by ad hoc considerations. Numerical results show that the use of the SR procedures greatly economizes computation and provides insight into the mechanical behavior while preserving a level of accuracy comparable with the fully three-dimensional solution.
Isaac Harari and Kirill Gosteev
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
A comparison of two bubble-enriched methods, derived by different considerations, indicates that the methods are identical in some cases. Thus, series representations of auxiliary functions, derived independently for the two methods, turn out to be equivalent prior to truncation. Three such series for time-harmonic acoustics are considered. Dispersion analysis points to the more efficient series representation and provides guidelines for the number of terms to be retained. Numerical tests confirm the validity of these practical guidelines.
Hashem M. Mourad and John Dolbow
Isaac Harari
Department of Civil and Environmental Engineering, Duke University
Durham, NC 27708-0287, U.S.A.
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
We examine a bubble-stabilized finite element method for enforcing Dirichlet constraints on embedded interfaces. By embedded we refer to problems of general interest wherein the geometry of the interface is assumed independent of some underlying bulk mesh. As such, the robust imposition of Dirichlet constraints with a Lagrange multiplier field is not trivial. To focus issues, we consider a simple one-sided problem that is representative of a wide class of evolving interface problems. The bulk field is decomposed into coarse and fine scales, giving rise to coarse-scale and fine-scale one- sided sub-problems. The fine-scale solution is approximated with bubble functions, permitting static condensation and giving rise to a stabilized form bearing strong analogy with a classical method. Importantly, the method is simple to implement, readily extends to multiple dimensions, obviates the need to specify any free stabilization parameters, and gives rise to a symmetric, positive-definite system of equations. The performance of the method is then examined through several numerical examples. The accuracy of the Lagrange multiplier is compared to results obtained using a local version of the domain integral method. The variational multiscale approach is shown to both stabilize the Lagrange multiplier and improve the accuracy of the post-processed fluxes.
Eran Grosu and Isaac Harari
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Solutions of direct time-integration schemes for elastodynamics that converge in time to conventional semidiscrete formulations may be polluted at small time steps by noncausal oscillations. This degradation is not an artifact of the time-marching scheme, but rather a property of the solution of the semidiscrete formulation itself. An analogy to singularly perturbed elliptic problems provides an upper bound on the time step for the onset of these oscillations. A simple procedure of spatial stabilization is proposed to remove this pathology from implicit time-integration schemes, without affecting unconditional temporal stability. Spatially stabilized implicit time integration methods are free of noncausal oscillations at small time steps.
Isaac Harari
Guillermo Hauke
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
Litec -- Departamento de Mecánica de Fluidos,
Centro Politécnico Superior
C/Maria de Luna 3, Zaragoza 50018, Spain
Solutions of direct time-integration schemes for transient advection-diffusion-reaction problems that converge in time to conventional semidiscrete formulations may be polluted at small time steps by spurious spatial oscillations. This degradation is not an artifact of the time-marching scheme, but rather a property of the solution of the semidiscrete Galerkin approximation itself. An analogy to steady advection-diffusion-reaction problems with a modified reaction coefficient by the Rothe method of discretizing in time prior to spatial discretization provides an upper bound on the time step for the onset of spatial instability. Spatial stabilization removes this pathology, leading to stabilized implicit time integration schemes that are free of spurious oscillations at small time steps.
Assad A. Oberai
Paul E. Barbone
Isaac Harari
Department of Mechanical, Aerospace, and Nuclear Engineering,
Rensselaer Polytechnic Institute
Troy, NY 12180, USA
Department of Aerospace and Mechanical Engineering, Boston University
Boston, MA 02215, USA
Department of Solid Mechanics, Materials and Systems, Tel-Aviv University
69978 Ramat Aviv, Israel
An alternative variational framework suitable for pure advection is obtained by discarding the Galerkin part of stabilized methods. The resulting scheme is similar to the least-squares approach, but with the adjoint operator in the weighting slot. This formulation is not restricted to solenoidal (i.e., divergence free) velocities. Initial numerical results for such configurations show that the method is promising.
Paul E. Barbone
Assad A. Oberai
Isaac Harari
Department of Aerospace and Mechanical Engineering, Boston University
Boston, MA 02215, USA
Department of Mechanical, Aerospace, and Nuclear Engineering,
Rensselaer Polytechnic Institute
Troy, NY 12180, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
We consider the direct (i.e. non-iterative) solution of the inverse problem of heat conduction for which at least two interior temperature fields are available. The strong form of the problem for the single, unknown, thermal conductivity field, is governed by two partial dirential equations of pure advective transport. The given temperature fields must satisfy a compatibility condition for the problem to have a solution. We introduce a novel variational formulation, the Adjoint Weighted Equation (AWE), for solving the two-field problem. In this case, the gradients of two given temperature fields must be linearly independent in the entire domain, a weaker condition than the compatibility required by the strong form. We show that the solution of the AWE formulation is equivalent to that of the strong form when both are well posed. We prove that the Galerkin discretization of the AWE formulation leads to a stable, convergent numerical method that has optimal rates of convergence. We show computational examples that confirm these optimal rates. The AWE formulation shows good numerical performance on problems with both smooth and rough coefficients and solutions.
Eran Grosu and Isaac Harari
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
The discontinuous enrichment method (DEM) for boundary-value problems governed by the Helmholtz equation is reviewed. Quadrilateral and triangular DEM elements for acoustics are considered. Conditioning considerations indicate preferred representations of the oscillatory basis functions. Dispersion properties are used to rate the performance of different element configurations. Numerical results indicate that high-order DEM elements exhibit little pollution. The dispersion and the numerical results are in a good agreement. The proposed configurations of DEM elements become more competitive as the enrichment and Lagrange multipliers are enhanced.
Eran Grosu and Isaac Harari
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Hexahedral and tetrahedral elements are proposed for the discontinuous enrichment method (DEM) for boundary-value problems governed by the Helmholtz equation. A procedure for obtaining sets of approximately uniform spherical enrichment directions, with high flexibility in the choice of their number, is constructed. Conditioning considerations indicate preferred representations of the oscillatory basis functions. Dispersion properties are used to rate the performance of different element configurations. Numerical tests assess accuracy, indicating that high-order DEM elements exhibit little pollution. The dispersion and the numerical results are in a good agreement. The proposed configurations of DEM elements become more competitive as the enrichment and Lagrange multipliers are enhanced.
John Dolbow
Isaac Harari
Department of Civil and Environmental Engineering, Duke University
Durham, North Carolina 27708, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
A stabilized finite element method based on the Nitsche technique for enforcing constraints leads to an efficient computational procedure for embedded interface problems. We consider cases in which the jump of a field across the interface is given, as well as cases in which the primary field on the interface is given. The finite element mesh need not be aligned with the interface geometry. We present closed-form analytical expressions for interfacial stabilization terms, and simple procedures for accurate flux evaluations. Representative numerical examples demonstrate the effectiveness of the proposed methodology.
Isaac Harari
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Standard, low-order, continuous Galerkin finite element methods are unable to cope with wave phenomena at short wave lengths because the computational effort required to resolve the waves and control numerical dispersion errors becomes prohibitive. The failure to adequately represent subgrid scales misses not only the fine-scale part of the solution, but often causes severe pollution of the solution on the resolved scale as well. Since computation naturally separates the scales of a problem according to the mesh size, multiscale considerations provide a useful framework for viewing these difficulties and developing methods to counter them. The Galerkin/least-squares method arises in multiscale settings, and its stability parameter is defined by dispersion considerations. Bubble enriched methods employ auxiliary functions that are usually expressed in the form of infinite series. Dispersion analysis provides guidelines for the implementation of the series representation in practice. In the discontinuous enrichment method, the fine scales are spanned by free-space homogeneous solutions of the governing equations. These auxiliary functions may be discontinuous across element boundaries, and continuity is enforced weakly by Lagrange multipliers.
Isaac Harari and Nir Makmel
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
In this paper we analyze the dispersion properties of plane-strain elements obtained by the discontinuous enrichment method (DEM). In DEM, the standard polynomial field is enriched within each element by a non-conforming field that is added to it. The dispersion analysis is performed on Lagrange multipliers, which weakly enforce the continuity of the solution and thus are associated with the inter-element traction. We use an eigenvalue-based procedure to compute the phase and polarization errors. The dispersion analysis is corroborated by results of numerical experiments.
Uri Albocher
Assad A. Oberai
Paul E. Barbone
Isaac Harari
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mechanical, Aerospace, and Nuclear Engineering,
Rensselaer Polytechnic Institute
Troy, NY 12180, USA
Department of Mechanical Engineering, Boston University
Boston, MA 02215, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
The adjoint-weighted equation (AWE) formulation is applied to incompressible isotropic plane stress elasticity for the solution of inverse problems. For the strong form of the problem to have a solution, a restrictive strain compatibility condition must hold. In the AWE variational formulation on the other hand, a much milder sufficient condition is that uniaxial stress states be precluded, making the method more robust to noisy data. We prove that the formulation leads to a stable converging numerical method, and show through computations that the method performs well.
Paul E. Barbone, Carlos E. Rivas
Isaac Harari, Uri Albocher
Assad A. Oberai, Yixiao Zhang
College of Engineering, Boston University
Boston, MA 02215, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mechanical, Aerospace, and Nuclear Engineering,
Rensselaer Polytechnic Institute
Troy, NY 12180, USA
We describe a novel variational formulation of inverse elasticity problems given interior data. The class of problems considered is rather general and includes, as special cases, plane deformations, compressibility and incompressiblity in isotropic materials, 3D deformations, and anisotropy. The strong form of this problem is governed by equations of pure advective transport. The variational formulation is based on a generalization of the adjoint-weighted variational equation (AWE) formulation, originally developed for flow of a passive scalar. We describe how to apply AWE to various cases, and prove several properties. We prove that the Galerkin discretization of the AWE formulation leads to a stable, convergent numerical method, and prove optimal rates of convergence. The numerical examples demonstrate optimal convergence of the method with mesh renement for multiple unknown material parameters, graceful performance in the presence of noise, and robust behavior of the method when the target solution is C1, C0, or discontinuous.
Isaac Harari
John Dolbow
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Civil and Environmental Engineering, Duke University
Durham, North Carolina 27708, USA
A stabilized finite element method based on the Nitsche technique for enforcing constraints leads to an efficient computational procedure for embedded interface problems, in which the finite element mesh need not be aligned with the interface geometry. We consider cases in which the jump of a field across the interface is given, as well as cases in which the primary field on the interface is given. Optimal rates of convergence hold. Representative numerical examples demonstrate the effectiveness of the proposed methodology.
Anand Embar and John Dolbow
Isaac Harari
Department of Civil and Environmental Engineering, Duke University
Durham, North Carolina 27708, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
A key challenge while employing non-interpolatory basis functions in finite element methods is the robust imposition of Dirichlet boundary conditions. The current work studies the weak enforcement of such conditions for B-spline basis functions, with application to both second and fourth-order problems. This is achieved using concepts borrowed from Nitsche's method, which is a stabilized method for imposing constraints on surfaces. Conditions for the stability of the system of equations are derived for each class of problem. Stability parameters in the Nitsche weak form are then evaluated by solving a local generalized eigenvalue problem at the Dirichlet boundary. The approach is designed to work equally well when the grid used to build the splines conforms to the physical boundary of interest as well as the more general case when it does not. Through several numerical examples, the approach is shown to yield optimal rates of convergence.
Isaac Harari, Igor Sokolov, and Slava Krylov
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Structural models are well-established for the governing operators in solid mechanics, yet the reduction of loads (data) is often performed in an ad hoc manner, which may be inadequate for the complex load distributions that often arise in modern applications. In the present work we consistently convert three-dimensional data to the form required by the Kirchhoff thin plate theory, in a variational framework. We provide formulas for all types of resultant structural loads and boundary conditions in terms of the original three-dimensional data, including proper specification of corner forces, in forms that are readily incorporated into computational tools. In particular, we find that distributed couples, engendered by in-plane components of three-dimensional loads, contribute to an effective distributed transverse force and boundary shear force, the latter generalizing the notion of the celebrated Kirchhoff equivalent force. However, in virtual work we advocate a representation of the twisting moment in a form that involves neither the Kirchhoff equivalent force nor corner forces. An interpretation of the structural deflections as through-the-thickness averages of the continuum displacements, rather than their values on the midplane, yields explicit formulas for the thin plate essential boundary data. The formulation facilitates the solution of problems that would otherwise pose formidable challenges. Numerical results confirm that appropriate use of the thin plate model economizes computation and provides insight into the mechanical behavior, while preserving a level of accuracy comparable with the full three-dimensional solution.
Isaac Harari, Ran Ganel, and Eran Grosu
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Standard finite element methods for wave phenomena entail considerable computational effort at short wavelengths in order to control numerical dispersion and pollution errors. Stabilized methods such as Galerkin/least-squares combine improvement in performance with simple implementation. The mesh-dependent stability parameter is often defined by dispersion considerations. The application to elastic waves must account for polarization errors as well. Balancing dispersion and polarization errors provides the best performance.
Yixiao Zhang, Assad A. Oberai
Paul E. Barbone
Isaac Harari
Department of Mechanical, Aerospace, and Nuclear Engineering,
Rensselaer Polytechnic Institute
Troy, NY 12180, USA
Department of Mechanical Engineering, Boston University
Boston, MA 02215, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
We consider the problem of determining the distribution of the complex-valued shear modulus for an incompressible linear viscoelastic material undergoing infinitesimal time-harmonic deformation, given the knowledge of the displacement field in its interior. In particular we focus on the two-dimensional problems of anti-plane shear and plane stress. These problems are motivated by applications in biomechanical imaging, where the material modulus distributions are used to detect and/or diagnose cancerous tumors. We analyze the well-posedness of the strong form of these problems and conclude that for the solution to exist the measured displacement field is required to satisfy rather restrictive compatibility conditions. We propose a weak, or a variational formulation, and prove the existence and uniqueness of solutions under milder conditions on measured data. This formulation is derived by weighting the original partial differential equation for the shear modulus by the adjoint operator acting on the complex-conjugate of the weighting functions. For this reason we refer to it as the complex adjoint weighted equation (CAWE). We consider a straightforward finite element discretization of these equations with total variation regularization, and test its performance with synthetically generated and experimentally measured data. We find that the CAWE method is in general less diffusive than a corresponding least squares solution, and that the total variation regularization significantly improves its performance in the presence of noise.
Igor Sokolov, Slava Krylov, and Isaac Harari
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
In MEMS devices, solids, often slender in geometry, are in nonlinear interaction with complex three-dimensional electrostatic fields The computational cost of solving these coupled problems can be reduced considerably by the use of structural models. A geometrically exact planar beam model is used for the solid, with particular attention to normal tractions on the interface that arise from electrostatic pressure distribution. The weakly coupled problem is solved with a staggered strategy. The resulting scheme provides accuracy comparable to that obtained by full, three-dimensional representations of the solid, at costs that may be reduced significantly.
Isaac Harari and Evgeny Shavelzon
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
A stabilized variational formulation, based on Nitsche's method for enforcing boundary constraints, leads to an efficient procedure for embedding kinematic boundary conditions in thin plate bending. The absence of kinematic admissibility constraints allows the use of non-conforming meshes with non-interpolatory approximations, thereby providing added flexibility in addressing the C1-continuity requirements typical of these problems. Work-conjugate pairs weakly enforce kinematic boundary conditions. The pointwise enforcement of corner deflections is key to good performance in the presence of corners. Stabilization parameters are determined from local generalized eigenvalue problems, guaranteeing coercivity of the discrete bilinear form. The accuracy of the approach is verified by representative computations with bicubic C2 B-splines, exhibiting optimal rates of convergence and robust performance with respect to values of the stabilization parameters.
Uri Albocher
Paul E. Barbone
Michael S. Richards
Assad A. Oberai
Isaac Harari
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mechanical Engineering, Boston University
Boston, MA 02215, USA
School of Medicine and Dentistry, University of Rochester
Rochester, NY 14642, USA
Department of Mechanical, Aerospace, and Nuclear Engineering,
Rensselaer Polytechnic Institute
Troy, NY 12180, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
We apply the adjoint weighted equation method (AWE) for the direct solution of inverse problems of incompressible plane strain elasticity involving high levels of displacement noise and discontinuous solutions. We show that when the noisy displacements are applied directly to the AWE formulation, the reconstruction of the shear modulus can be very poor. On the other hand we demonstrate that by smoothing the displacements and appending a regularization term to the AWE formulation, a dramatic improvement in the reconstruction can be achieved. With these improvements, the advantages of the AWE method as a direct (non-iterative) solution approach can be extended to a wider range of problems.
Mengyu Lan and Haim Waisman
Isaac Harari
Department of Civil Engineering & Engineering Mechanics,
Columbia University
New York, NY 10027, USA
Department of Mechanical Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
A new analytical approach, within the extended finite element (XFEM) framework, is proposed to compute Strain Energy Release Rates (SERRs) directly from Irwin's integral. Crack tip enrichment functions in XFEM allow for evaluation of integral quantities in closed form (for some crack configurations studied) and therefore results in an accurate and efficient method.
Several benchmark examples on pure and mixed mode problems are studied. In particular we analyze the effects of high order enrichments, mesh refinement and the integration limits of Irwin's integral. The results indicate that high order enrichment functions have significant effect on the convergence, in particular when the integral limits are finite. When the integral limits tend to zero, simpler SERR expressions are obtained and high order terms vanish. Nonetheless, these terms contribute indirectly via coefficients of first order terms.
The numerical results show that high accuracy can be achieved with high order enrichment terms and mesh refinement. However, the effect of the integral limits remains an open question, with finite integration intervals chosen as h=2 tend to give more accurate results.
Mengyu Lan and Haim Waisman
Isaac Harari
Department of Civil Engineering & Engineering Mechanics,
Columbia University
New York, NY 10027, USA
Department of Mechanical Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
An analytical formulation based on Irwin's integral and combined with the extended finite element method (XFEM) is proposed to extract mixed-mode components of Strain Energy Release Rates (SERRs) in linear elastic fracture mechanics. The proposed formulation extends previous work to cracks in arbitrary orientations and is therefore suited for crack propagation problems. In essence, the approach employs high order enrichment functions and evaluates Irwin's integral in closed form once the linear system is solved and the algebraic degrees of freedom are determined.
Several benchmark examples are investigated including off-center cracks, inclined cracks and two crack growth problems. On all these problems, the method is shown to work well, giving accurate results. Moreover, due to its analytical nature, no special postprocessing is required. Thus we conclude that this method may provide a good and simple alternative to the popular J-integral method. Moreover, it may circumvent some of the limitations of the J-integral in 3D modeling and in problems involving branching and coalescence of cracks.
Timo Saksala
Delphine Brancherie
Isaac Harari
Adnan Ibrahimbegovic
Department of Mechanics and Design,
Tampere University of Technology
FIN-33101, Tampere, Finland
Laboratoire Roberval, Université de Compiègne
BP 20529 – 60205, Compiègne, France
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Ecole Normale Supérieure de Cachan, LMT-Cachan
61, avenue du président Wilson, 94235 Cachan, France
In this paper, a novel constitutive model combining continuum damage with embedded discontinuity is developed for explicit dynamic analyses of quasi-brittle failure phenomena. The model is capable of describing the rate-dependent behavior in dynamics and the three phases in failure of quasi-brittle materials. The first phase is always linear elastic, followed by the second phase corresponding to fracture-process zone creation, represented with rate-dependent continuum damage with isotropic hardening formulated by utilizing consistency approach. The third and final phase, involving nonlinear softening, is formulated by using an embedded displacement discontinuity model with constant displacement jumps both in normal and tangential directions. The proposed model is capable of describing the rate-dependent ductile to brittle transition typical of cohesive materials (e.g. rocks, ice). The model is implemented in the finite element setting by using the CST elements. The displacement jump vector is solved for implicitly at the local (finite element) level along with a viscoplastic return mapping algorithm, while the global equations of motion are solved with explicit time-stepping scheme. The model performance is illustrated by several numerical simulations, including both material point and structural tests. The final validation example concerns the dynamic Brazilian disc test on rock material (Kuru granite) under plane stress assumption.
Uri Albocher
Paul E. Barbone
Assad A. Oberai
Isaac Harari
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mechanical Engineering, Boston University
Boston, MA 02215, USA
Department of Mechanical, Aerospace, and Nuclear Engineering,
Rensselaer Polytechnic Institute
Troy, NY 12180, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
The uniqueness of an inverse problem of isotropic incompressible three dimensional elasticity aimed at reconstructing material modulus distributions is considered. We show that given a single strain field and no boundary conditions, the solution for the shear modulus can involve arbitrary functions of space. On the other hand, having two independent strain fields, leads to a favorable solution space of up to five arbitrary constants. We solve inverse problems with two strain fields given, using the adjoint weighted equation method and impose five discrete constraints. The method exhibits good numerical performance with optimal rates of convergence.
Isaac Harari and Eran Grosu
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
An efficient procedure for embedding kinematic boundary conditions in the biharmonic equation, for problems such as the pure streamfunction formulation of the Navier-Stokes equations and thin plate bending, is based on a stabilized variational formulation, obtained by Nitsche's approach for enforcing boundary constraints. The absence of kinematic admissibility constraints allows the use of non-conforming meshes with non-interpolatory approximations, thereby providing added flexibility in addressing the higher continuity requirements typical of these problems. Variationally conjugate pairs weakly enforce kinematic boundary conditions. The use of a scaling factor leads to a formulation with a single stabilization parameter. For plates, the enforcement of tangential derivatives of deflections obviates the need for pointwise enforcement of corner values in the presence of corners. The single stabilization parameter is determined from a local generalized eigenvalue problem, guaranteeing coercivity of the discrete bilinear form. The accuracy of the approach is verified by representative computations with bicubic B-splines, providing guidance to the determination of the scaling, exhibiting optimal rates of convergence, and robust performance with respect to values of the stabilization parameter.
Gan Song, Haim Waisman, and Mengyu Lan
Isaac Harari
Department of Civil Engineering & Engineering Mechanics,
Columbia University
New York, NY 10027, USA
Department of Mechanical Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
The current paper extends our recent work on the extraction of Stress Intensity Factors (SIFs) from Irwin's integral, using high order XFEM formulation [1, 2]. By matching leading r terms in XFEM with analytical expansion of Irwin's Integral, a closed form solution for SIFs is obtained. Hence SIFs can directly be obtained upon solution of the XFEM discrete system, which were shown to be accurate on structured rectangular meshes. However, extension to unstructured quadrilateral elements may not be trivial.
To this end, we extend the formulation to triangular elements, where significant simplification of the closed form expression is obtained. Moreover, the formulation is directly applicable to unstructured meshes without extra work. Hence this paper is considered a significant step towards automating and extending this method to more general applications.
Numerical results show that accurate and consistent SIFs can be obtained on some pure mode and inclined crack benchmark problems, on structured as well as unstructured meshes. Examples of a crack approaching a hole and two cracks approaching each other are also investigated. The latter illustrate the advantage of this method over a J-integral class of methods, as SIFs can still be calculated when cracks are in close proximity and no remeshing is required. Hence potentially this method can address crack coalescence and branching more rigorously.
Wen Jiang
Chandrasekhar Annavarapu
John Dolbow
Isaac Harari
Department of Mecahnical Engineering and Material Science, Duke University
Durham, NC 27708, USA
Atmospheric, Earth, and Energy Division
Lawrence Livermore National Laboratory
Livermore, CA 94550, USA
Department of Civil and Environmental Engineering, Duke University
Durham, North Carolina 27708, USA
Department of Mechanical Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
The extended finite element method (X-FEM) has proven to be an accurate, robust method for solving embedded interface problems. With a few exceptions, the X-FEM has mostly been used in conjunction with piecewise-linear shape functions and an associated piecewise-linear geometrical representation of interfaces. In the current work, the use of spline-based finite elements is examined along with a Nitsche technique for enforcing constraints on an embedded interface. To obtain optimal rates of convergence, we employ a hierarchical local refinement approach to improve the geometrical representation of curved interfaces. We further propose a novel weighting for the interfacial consistency terms arising in the Nitsche variational form with B-splines. A qualitative dependence between the weights and the stabilization parameters is established with additional element level eigenvalue calculations. An important consequence of this weighting is that the bulk as well as the interfacial fields remain well behaved in the presence of large heterogeneities as well as elements with arbitrarily small volume fractions. We demonstrate the accuracy and robustness of the proposed method through several numerical examples.
Igor Sokolov, Slava Krylov, and Isaac Harari
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Geometrically exact beam theory is extended to allow distortion of the cross section. We present an appropriate set of cross-section basis functions and provide physical insight to the crosssectional distortion from linear elastostatics. The beam formulation in terms of material (backrotated) beam internal force resultants and work-conjugate kinematic quantities emerges naturally from the material description of virtual work of constrained finite elasticity. The inclusion of crosssectional deformation allows straightforward application of three-dimensional constitutive laws in the beam formulation. Beam counterparts of applied loads are expressed in terms of the original threedimensional data. Special attention is paid to the treatment of the applied stress, keeping in mind applications such as hydrogel actuators under environmental stimuli or devices made of electroactive polymers. Numerical comparisons show the ability of the beam model to reproduce finite elasticity results with good efficiency.
Leszek Demkowicz
Isaac Harari
Institute for Computational Engineering and Sciences,
University of Texas at Austin
Austin, TX 78712, USA
School of Mechanical Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Using the reaction-dominated diffusion problem, we investigate various DPG formulations in an attempt to construct a scheme that delivers robust discretization in the norm of our choice. In particular, we generalize the strategy of Heuer and Demkowicz [9] to the primal DPG formulation.
Yixiao Zhang
Paul E. Barbone
Isaac Harari
Assad A. Oberai
Department of Mechanical, Aerospace, and Nuclear Engineering,
Rensselaer Polytechnic Institute
Troy, NY 12180, USA
Department of Mechanical Engineering, Boston University
Boston, MA 02215, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mechanical, Aerospace, and Nuclear Engineering,
Rensselaer Polytechnic Institute
Troy, NY 12180, USA
Elasticity imaging has emerged as a promising medical imaging technique with applications in the detection, diagnosis and treatment monitoring of several types of disease. In elasticity imaging measured displacement fields are used to generate images of elastic parameters of tissue by solving an inverse problem. When the tissue excitation, and the resulting tissue motion is time-harmonic, elasticity imaging can be extended to image the viscoelastic properties of the tissue. This leads to an inverse problem for the complex-valued shear modulus at a given frequency. In this manuscript we have considered the uniqueness of this inverse problem for an incompressible, isotropic linear viscoelastic solid in a state of plane strain. For a single measured displacement field we conclude that the solution is infinite dimensional, and the data required to render it unique is not easily quantified. In contrast, for two independent displacement fields such that the principal directions of the resulting strain fields are different, the space of possible solutions is eight dimensional, and given additional data, like the value of the shear modulus at four locations, or over a calibration region, we may determine the shear modulus everywhere. We have also considered simple analytical examples that verify these results and offer additional insights. The results derived in this paper may be used as guidelines by the practitioners of elasticity imaging in designing more robust and accurate imaging protocols.
Dominik Schillinger
Isaac Harari
Ming-Chen Hsu
David Kamensky
Klaas F.S. Stoter
Yue Yu
Ying Zhao
Department of Civil, Environmental, and Geo-Engineering,
University of Minnesota
Minneapolis, MN 55455, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mechanical Engineering, Iowa State University
Ames, IA 50011, USA
Institute for Computational Engineering and Sciences,
University of Texas at Austin
Austin, TX 78712, USA
Department of Civil, Environmental, and Geo-Engineering,
University of Minnesota
Minneapolis, MN 55455, USA
Department of Mathematics,
Lehigh University
Bethlehem, PA 18015, USA
Graduate School of Computational Engineering,
Technical University of Darmstadt
64287 Darmstadt, Germany
We explore the use of the non-symmetric Nitsche method for the weak imposition of boundary and coupling conditions along interfaces that intersect through a finite element mesh. In contrast to symmetric Nitsche methods, it does not require stabilization and therefore does not depend on the appropriate estimation of stabilization parameters. We first review the available mathematical background, recollecting relevant aspects of the method from a numerical analysis viewpoint. We then compare accuracy and convergence of symmetric and non-symmetric Nitsche methods for a Laplace problem, a Kirchhoff plate, and in 3D elasticity. Our numerical experiments confirm that the non-symmetric method leads to reduced accuracy in the L2 error, but exhibits superior accuracy and robustness for derivative quantities such as diffusive flux, bending moments or stress. Based on our numerical evidence, the non-symmetric Nitsche method is a viable alternative for problems with diffusion-type operators, in particular when the accuracy of derivative quantities is of primary interest.
Yongxiang Wang and Haim Waisman
Isaac Harari
Department of Civil Engineering & Engineering Mechanics,
Columbia University
New York, NY 10027, USA
Department of Mechanical Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
This paper presents a comprehensive study on the use of Irwin's crack closure integral, within a high-order extended finite element method (XFEM), for direct evaluation of mixed-mode stress intensity factors (SIFs) in curved crack problems. The approach employs high-order enrichment functions derived from Williams asymptotic solution and SIFs are computed in closed-form without any special postprocessing requirements. Linear triangular elements are used to discretize the domain, and the crack curvature within an element is represented explicitly. An improved quadrature scheme using high-order isoparametric mapping together with a generalized Duffy transformation is proposed to integrate singular fields in tip elements with curved cracks. Furthermore, since Williams asymptotic solution is derived for straight cracks, an appropriate definition of the angle in the enrichment functions is presented and discussed. This contribution is an important extension of our previous work on straight cracks and illustrates the applicability of the SIF extraction method to curved cracks.
The performance of the method is studied on several circular and parabolic arc crack benchmark examples. With two layers of elements enriched in the vicinity of the crack tip, striking accuracy, even on relatively coarse meshes, is obtained and the method converges to reference SIFs at a rate of . Furthermore, while the popular interaction integral (a variant of the J-integral method) requires special auxiliary fields for curved cracks and also requires cracks to be sufficiently apart from each other in multicrack systems, the proposed approach shows none of those limitations.
Zilong Zou and Wilkins Aquino
Isaac Harari
Department of Civil and Environmental Engineering, Duke University
Durham, NC 27708, USA
Department of Mechanical Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
In this work, we use Nitsche's formulation to weakly enforce kinematic constraints at an embedded interface in Helmholtz problems. Allowing embedded interfaces in a mesh provides significant ease for discretization, especially when material interfaces have complex geometries. We provide analytical results that establish the well-posedness of Helmholtz variational problems and convergence of the corresponding finite element discretizations when Nitsche's method is used to enforce kinematic constraints. As in the analysis of conventional Helmholtz problems, we show that the inf-sup constant remains positive provided that the Nitsche's stabilization parameter is judiciously chosen. We then apply our formulation to several 2D plane-wave examples that confirm our analytical findings. Doing so, we demonstrate the asymptotic convergence of the proposed method and show that numerical results are in accordance with the theoretical analysis.
Isaac Harari
Uri Albocher
Department of Mechanical Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mechanical Engineering, Tel-Aviv University, 69978 Ramat
Aviv, Israel
Afeka, Tel Aviv Academic College of Engineering, Israel
Incompatible discretization methods provide added flexibility in computation by allowing meshes to be unaligned with geometric features and easily accommodating non-interpolatory approximations. Such formulations that are based on Nitsche's approach to enforce surface constraints weakly, which shares features with stabilized methods, combine conceptual simplicity and computational efficiency with robust performance. The basic workings of the method are well understood, in terms of a bound on the parameter. However, its spectral behavior has not been explored in depth. Such investigations can shed light on properties of the operator that effect the solution of boundary-value problems. Furthermore, incompatible discretizations are rarely used for eigenvalue problems. The spectral investigations lead to practical procedures for solving eigenvalue problems that are formulated by Nitsche's approach, with bearing on explicit dynamics.
Bingbing San
Haim Waisman
Isaac Harari
College of Civil and Transportation Engineering,
Hohai University
China 210098
Department of Civil Engineering & Engineering Mechanics,
Columbia University
New York, NY 10027, USA
Department of Mechanical Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
We extend a classical solution concerning the shape optimization of a hanging rope. The well known solution was derived with the objective to find the optimal cross section of a homogeneous rope that minimizes elongation under its own weight and a given applied force, subject to a total volume constraint. Herein, the analytical solution is generalized to materials with variable density and elastic modulus along the rope, subject to a total mass constraint. In the second part, a gradient-based numerical optimization algorithm is developed and used to solve the inverse problem in order to validate the analytical results. The approach is then extended to two dimensional structures by parameterization of the external boundary using nonuniform rational B-Splines (NURBS) functions and solving repeated forward problems with updated meshes. We study three different cases: (i) homogeneous elastic, (ii) homogeneous hyperelastic and (iii) inhomogeneous elastic materials. The results show the differences between optimal shape of one- and two-dimensional models, and the effect of material models on the optimal solutions.
Paul E. Barbone and Navid Nazari
Isaac Harari
Department of Mechanical Engineering, Boston University
Boston, MA 02215, USA
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
The propagation of waves in elastic solids at or near the incompressible limit is of interest in many current and emerging applications. Standard low-order Galerkin finite element discretization struggles with both incompressibility and wave dispersion. Galerkin least-squares stabilization is known to improve computational performance of each of these ingredients separately. A novel approach of combined pressure-curl stabilization is presented, facilitating the use of continuous, equal-order interpolation of displacements and pressure. The pressure stabilization parameter is determined by stability considerations, while the curl stabilization parameter is determined by dispersion considerations. The proposed pressure-curl stabilized scheme provides stable and accurate results on a variety of numerical tests for incompressible and nearly incompressible elastic waves computed with linear elements.
Isaac Harari
Uri Albocher
Faculty of Engineering, Tel-Aviv University
69978 Ramat Aviv, Israel
Department of Mechanical Engineering, Tel-Aviv University, 69978 Ramat
Aviv, Israel
Afeka, Tel Aviv Academic College of Engineering, Israel
Embedded methods that are based on Nitsche's approach can facilitate the task of mesh generation in many configurations. The basic workings of the method are well understood, in terms of a bound on the stabilization parameter. However, its spectral behavior has not been explored in depth. In addition to the eigenpairs which approximate the exact ones, as in the standard formulation, Nitsche's method gives rise to mesh-dependent complementary pairs. The dependence of the eigenvalues on the Nitsche parameter is related to a boundary quotient of the eigenfunctions, explaining the manner in which stabilization engenders coercivity without degrading the accuracy of the discrete eigenpairs. The boundary quotient proves to be useful for separating the two types of solutions. The quotient space is handy for determining the number of eigenpairs and complementary pairs. The complementary solutions approximate functions in the orthogonal complement of the kinematically admissible subspace. A global result for errors in the Galerkin approximation of the eigenvalue problem that pertains to all modes of the discretization, is extended to the Nitsche formulation. Numerical studies on non-conforming aligned meshes confirm the dependence of the eigenvalues on the parameter, in line with the corresponding boundary quotients. The spectrum of a reduced system obtained by algebraic elimination is free of complementary solutions, warranting its use in the solution of boundary-value problems. The reduced system offers an incompatible discretization of eigenvalue problems that is suitable for engineering applications. Using Irons-Guyan reduction yields a spectrum that is virtually insensitive to stabilization, with high accuracy in both eigenvalues and eigenfunctions.
Please send comments or questions about this page to harari AT eng.tau.ac.il