A level set based topology optimization for finite unidirectional acoustic phononic structures using boundary element method

doi:10.1016/j.cma.2021.113776

Computer Methods in Applied Mechanics and Engineering

Volume 381, 1 August 2021, 113776

https://doi.org/10.1016/j.cma.2021.113776 Get rights and content

Highlights

•
Average topological derivative is specified for each design domain.
•
Original and adjoint fields are both computed through a size reduced system matrix.
•
Resluts from different numbers of design domains are obtained.
•
The band gap properties are enhanced with more layers.

Abstract

The paper presents a level set based topology optimization method for unidirectional phononic structures with finite layers of lattice cells. Boundary element method(BEM) is employed as the numerical approach to solve the acoustic problems governed by Helmholtz equation. A sized reduced coefficient matrix is derived due to the iteration forms for input and output quantities on the periodic boundary of unit cells. Topological derivatives are formulated by boundary integral equation combined with adjoint variable method and computed for each layer. An average topological sensitivity of a single design domain is proposed for the updating of the level set function(LSF) governed by an evolution equation. Numerical models with different number of layers are considered and several optimized structures of unit cells are obtained in concerned frequencies. A further investigation into the transmission of acoustic waves is carried out by employing more layers of the periodic structures between the input and output domains. The results demonstrate the effectiveness of the proposed optimization method for the finite unidirectional phononic structures.

Introduction

There are continuing interests in design of phononic crystals and metamaterials which are artificially manufactured periodic structures composed of multi materials [1], [2] . The novel properties of phononic structures such as waveguides [3], [4], negative refraction [5], [6], [7] and negative elastic modulus [8], [9], [10], make them have a promising prospect of applications in manipulation of acoustic or elastic waves. Since the experimental investigation of the propagation characteristics of acoustic waves in a sculpture located in Madrid [11], a great amount of efforts has been devoted to developing numerical methods for simulations of phononic structures [12], [13], [14], [15]. The design of phononic structures is usually carried out by adjusting the configuration of the inclusions and host materials to extend or obtain band gaps in desired ranges [16], [17], [18]. The influence of geometric and material parameters, filling fraction and shapes of inclusions on the width of the lowest band gap is investigated through finite element method (FEM) [16]. A tunable passband is designed by applying external electric field onto the piezoelectric layers in phononic crystals [17]. Unidirectional phononic crystals with functionally graded interlayers are investigated through transfer matrix method [18].

With increased interest in the development of phononic devices, it is necessary to construct effective methodologies for the design of phononic structure/metamaterials. In preliminary studies, an investigation of optimization design for maximizing the band gaps of phononic structures is carried out by using FEM [19], and after that, various methodologies [20], [21], [22] such as genetic algorithm, evolutionary algorithm, and level set method, are adopted for optimization design of 2D phononic crystals. The combination of optimization algorithms and various numerical methods [23], [24], [25] provides more efficient approaches for the creation of unit cells through computational techniques. A topological design on the unit cell optimization is carried out through FEM and bidirectional evolutionary structural method to obtain ultrawide band gap for in-plane and out-of-plane wave modes [23]. Adaptive genetic algorithm is also employed to calculate the mutation and crossover rate for the optimization of topological configurations [24]. In the literature [25], bidirectional evolutionary structure optimization method is applied to the design of phononic crystals for maximizing spatial decay of evanescent waves. The density based-topology optimization methods [26], [27] usually result in a continuous distribution of density which shows grayscale areas. Therefore, it is difficult to specify the boundary conditions in the fixed design domain. The adoption of level set-based structural optimization methods [28], [29] avoids the intermediate state of the boundary definition. However, the boundary conditions on the fixed design domain are approximated through ersatz material approach [30] which employs a domain integral to replace the boundary integral. To avoid the weak material or intermediate material and generate actual boundaries defined by the iso-surface of LSF, we employ BEM as a numerical approach.

Boundary-only discretization provided by BEM [31] leads to dimension-reduced elements. Therefore, BEM simplifies the modeling process and reduces the DOFs. BEM has been applied to sensitivity analysis for shape optimization problems [32], [33], [34]. With the development of isogeometric analysis (IGA), BEM provides promising applications as a numerical approach which is conveniently combined with the IGA due to the feature of boundary-only discretization. Many research works have been carried out to develop the IGA for various applications [35], [36], [37], [38]. The isogeometric BEM(IGABEM) [39], [40] is developed for acoustic problems and introduced to real-world engineering applications. Moreover, IGABEM is implemented for fracture modeling and crack propagation simulations [41], [42], since knot insertion in B-splines can provide discontinuities in the geometry. For 2D and 3D multi-patch coupling in statics, modal analysis and contact problems, the skew-symmetric Nitsche’s formulation is introduced by the Ref. [43]. Moreover, the IGABEM is also developed for optimization problems [44], [45], [46], [47], [48], [49]. For sound absorbing problems, the IGABEM is adopted for both 2D and 3D cases to optimize the distribution of materials [49], [50]. The research work in the Ref. [51] employs coarse meshes provided by Computer Aided Design (CAD) to reduce the number of design variables and fine meshes for accurate evaluation of sensitivities. Furthermore, Geometry Independent Field approximaTion (GIFT) proposed in the Ref. [52] allows to reduce the number of control points while acquiring a refinement of the variable space for interior regions, which means GIFT can be potentially applied to optimization problems to improve the efficiency.

The combination of BEM and level set method is investigated in various problems [53], [54], [55]. The adoption of level set method leads to clear definition of boundaries which can be directly discretized into boundary elements conveniently. Moreover, due to the features of numerical model from BEM and the properties of fundamental solutions, BEM is also suitable for acoustic or electromagnetic problems defined in infinite domains. The literatures [53], [55] present effective applications of the combination methodology for heat conduction problems and the literature [54] extends it to the optimization design for a coupling problem. The phononic structure composed of periodically distributed scatterers in the infinite domain is investigated through BEM by considering the scattering frequencies [56], [57] and the topology optimization of the scatterers is carried out through scattering matrix method combined with contour integral method [58], [59].

Optimization issues involving eigenvalues arise in many applications related to acoustic or elastodynamic problems [60], [61], [62]. The combination of the contour integral method with BEM requires repeat computations for the numerical evaluation of the contour integrals. The topological optimization algorithm usually leads to high computation cost since the evolution of topology demands the iteration of the solving process for eigenvalues.

However, due to the previous research [63], the finite periodic structures present band gap properties which is determined by the employed unit cell and the reduction of transmission of acoustic waves is enhanced as the number of cells increases. To avoid the high computation cost resulted by the Bloch eigenvalue value analysis which is usually carried out to illustrate the band structure of the unit cell, the topological optimization algorithm for finite unidirectional phononic structures is proposed in this work. Finite layers of design domains are considered and the average topological derivatives which are formulated using adjoint variable method, are adopted for each design domain. The original and adjoint fields are both solved through BEM with a size reduced coefficient matrix by formulating local system matrix of the unit cell. The source code which accompanies this paper is available for downloading from https://sourceforge.net/projects/topology-opt-bem-lsf/ for Linux and Mac OS computers.

The rest parts of the paper are organized as follows: Section 2 explains the derivation of the average topological derivatives using the adjoint variable method. In Section 3, the BEM formulas for the finite unidirectional phononic structures are presented. Several numerical examples with different conditions are investigated to demonstrate the effectiveness of the proposed optimization algorithm. The final parts consist of the conclusions and acknowledgment.

Section snippets

Level set method

Let us introduce a scalar function $ϕ (x)$ , called level set function defined in a design domain $D$ . With a given range $[- 1, 1]$ , $ϕ (x)$ is defined as follows: $\{\begin{matrix} 0 < ϕ (x) \leq 1, & x \in Ω ∖ Γ, \\ ϕ (x) = 0, & x \in Γ, \\ - 1 \leq ϕ (x) < 0, & x \in D ∖ Ω, \end{matrix}$ where $Ω$ is material region and $Γ$ denotes the boundary of $Ω$ . Obviously it can be seen that the distribution of the material is related to the value of the level set function $ϕ (x)$ . Thus an evolution of $ϕ (x)$ can be carried out to achieve the optimum topology configuration of the material. As shown in Fig. 1

BEM for finite unidirectional phononic structures

To formulate the system equations for the finite periodic structures, Eq. (2) can be transformed into expression of integral equation forms with the Green identity: $p (y) + \int_{Γ} q^{*} (x, y) p (x) d Γ (x) = \int_{Γ} p^{*} (x, y) q (x) d Γ (x) in Ω,$ and its normal directional derivative with respect to point $y$ is written as $q (y) + \int_{Γ} {\hat{q}}^{*} (x, y) p (x) d Γ (x) = \int_{Γ} {\hat{p}}^{*} (x, y) q (x) d Γ (x) in Ω,$ where the $p^{*} (x, y)$ is the fundamental solution for Eq. (2), $y$ , $x$ are collocation point and source point, respectively, and $r = | x - y |$ is the distance between $x$ and $y$ . $p^{*}$

Verification of topological derivative

A model with input and output domains is depicted in Fig. 5. The left boundary of input domain $D_{I}$ is specified with Dirichlet boundary condition $p = {\bar{p}}_{0} = 1$ [Pa] and the upper and lower boundaries are considered as rigid wall with $q = 0$ (Neumann boundary condition). The upper, lower and right boundaries of the output domain $D_{O}$ are also considered as rigid wall with $q = 0$ . Furthermore, the objective function $f = {| p - p_{obj} |}^{2}$ is defined and evaluated on the 2 elements in the center of right boundary of $D_{O}$ . To

Conclusions

This paper presents a level set based topology optimization method for finite unidirectional phononic structures using BEM. The adjoint method is employed to derive the topological derivative associated with the objective functional defined on a portion of boundary in the output domain and the acoustic problems formulated in adjoint and original fields are solved through boundary element method. The boundary elements are generated by discretizing the zero contour line of the level set function.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgment

This work is supported by the Applied Basic Research Programs of Shanxi Province of China [grant number 201801D221022], the National Natural Science Foundation of China [grant numbers 11872168, 52075361], the Key R & D Projects in Shanxi Province of China [grant number 201903D421030], the Major Special Projects of Science and Technology in Shanxi Province of China [grant number 20201102003], and the China Scholarship Council [grant number 201908140142].

References (66)

HuR. et al.
Spectral variational multiscale model for transient dynamics of phononic crystals and acoustic metamaterials
Comput. Methods Appl. Mech. Engrg.
(2020)
PereiraF. et al.
Phononic crystal investigation using a fluid-structure circular cylindrical shell spectral element
Mech. Syst. Signal Process.
(2021)
Guillén-GallegosC. et al.
Phononic band structure of an acoustic waveguide that behaves as a phononic crystal
Results Phys.
(2019)
LiC. et al.
Engineering of band gap and cavity mode in phononic crystal strip waveguides
Phys. Lett. A
(2013)
KeM. et al.
Flat superlens by using negative refraction in two-dimensional phononic crystals
Solid State Commun.
(2007)
LiJ. et al.
Negative refraction imaging of solid acoustic waves by two-dimensional three-component phononic crystal
Phys. Lett. A
(2008)
MokhtariA.A. et al.
On the emergence of negative effective density and modulus in 2-phase phononic crystals
J. Mech. Phys. Solids
(2019)
LiZ. et al.
Modelling of elastic metamaterials with negative mass and modulus based on translational resonance
Int. J. Solids Struct.
(2019)
WangX.
Dynamic behaviour of a metamaterial system with negative mass and modulus
Int. J. Solids Struct.
(2014)
KnippP.A. et al.
Boundary-element calculations of electromagnetic band-structure of photonic crystals
Physica E
(1998)

LiF.L. et al.

Boundary element method for band gap calculations of two-dimensional solid phononic crystals

Eng. Anal. Bound. Elem.

(2013)

KrushynskaA. et al.

Towards optimal design of locally resonant acoustic metamaterials

J. Mech. Phys. Solids

(2014)

FomenkoS. et al.

In-plane elastic wave propagation and band-gaps in layered functionally graded phononic crystals

Int. J. Solids Struct.

(2014)

HyunC.J. et al.

A genetic algorithm for multiple objective sequencing problems in mixed model assembly lines

Comput. Oper. Res.

(1998)

Garcia-LopezN. et al.

An improved robust topology optimization approach using multiobjective evolutionary algorithms

Comput. Struct.

(2013)

JiangY. et al.

Topology optimization under design-dependent loads with the parameterized level-set method based on radial-basis functions

Comput. Methods Appl. Mech. Engrg.

(2020)

WangK. et al.

Ultrawide band gap design of phononic crystals based on topological optimization

Physica B

(2019)

HanX. et al.

Bandgap design of three-phase phononic crystal by topological optimization

Wave Motion

(2020)

ChenY. et al.

Maximizing spatial decay of evanescent waves in phononic crystals by topology optimization

Comput. Struct.

(2017)

AllaireG. et al.

Structural optimization using sensitivity analysis and a level-set method

J. Comput. Phys.

(2004)

BrebbiaC.A.

The birth of the boundary element method from conception to application

Eng. Anal. Bound. Elem.

(2017)

AithalR. et al.

Shape sensitivity analysis in thermal problems using BEM

Eng. Anal. Bound. Elem.

(1995)

SianoD. et al.

Sensitivity analysis and correlation experimental/numerical FEM-bem for noise reduction assessment of an engine beauty cover

Energy Procedia

(2015)

ModaveA. et al.

Corner treatments for high-order local absorbing boundary conditions in high-frequency acoustic scattering

J. Comput. Phys.

(2020)

WuH. et al.

Isogeometric finite element analysis of interior acoustic problems

Appl. Acoust.

(2015)

SimpsonR. et al.

Acoustic isogeometric boundary element analysis

Comput. Methods Appl. Mech. Engrg.

(2014)

CooxL. et al.

A performance study of NURBS-based isogeometric analysis for interior two-dimensional time-harmonic acoustics

Comput. Methods Appl. Mech. Engrg.

(2016)

PengX. et al.

Isogeometric boundary element methods for three dimensional static fracture and fatigue crack growth

Comput. Methods Appl. Mech. Engrg.

(2017)

ChenL.L. et al.

Structural shape optimization of three dimensional acoustic problems with isogeometric boundary element methods

Comput. Methods Appl. Mech. Engrg.

(2019)

LianH. et al.

Shape optimization directly from cad: An isogeometric boundary element approach using t-splines

Comput. Methods Appl. Mech. Engrg.

(2017)

LiuC. et al.

Shape optimization of sound barrier using an isogeometric fast multipole boundary element method in two dimensions - sciencedirect

Eng. Anal. Bound. Elem.

(2017)

ChenL. et al.

Acoustic topology optimization of sound absorbing materials directly from subdivision surfaces with isogeometric boundary element methods

Comput. Methods Appl. Mech. Engrg.

(2020)

ChenL. et al.

An isogeometric approach of two dimensional acoustic design sensitivity analysis and topology optimization analysis for absorbing material distribution

Comput. Methods Appl. Mech. Engrg.

(2018)

Cited by (13)

WBEM-based analysis of band structures of solid-solid and fluid-fluid phononic crystals with frequency-independent fundamental solutions
2023, Engineering Analysis with Boundary Elements
A novel method is presented to determine the band structures of two-dimensional solid-solid and fluid-fluid phononic crystals (PCs) with square and triangular lattices using the wavelet-based boundary element method (WBEM). Both sweeping frequency technique and sweeping wave vector technique are used to determine the band structures of PCs. The fundamental solutions for integral equations established in primitive cell are independent of angular frequency, resulting in an avoidance of solving nonlinear eigenvalue problems and much reduction of computing time. The radial integration method and dual reciprocity method are respectively applied to handle the domain integral terms arising from the use of frequency-independent fundamental solutions, and to make the dimension reduction. The physical fields are approximated by B-spline wavelet on the interval and wavelet coefficients in wavelet space. It is proved that the wavelet coefficients still meet the relationship of their corresponding physical quantities, such as Bloch theorem. By means of the vanishing moment characteristic of wavelets, some small entries are produced in the system matrices which are further compressed into sparse ones, and the influence of sparse matrices on the results is discussed. Finally, the high computing speed and accuracy of the proposed method are shown in examples.
An exact volume constraint method for topology optimization via reaction–diffusion equation
2023, Computers and Structures
A novel approach is proposed to determine the exact Lagrangian multiplier $Λ$ for the volume constraint of topology optimization via the reaction–diffusion equation (RDE). Such $Λ$ enables an exact volume constraint. The mathematical approach involves splitting the density function (or level-set function) and their RDEs into two parts. Both become independent of $Λ$ , and the exact value of $Λ$ can then be determined by superposition to satisfy the constraint at the current step. Because only a few iterations are required on average, our proposed volume constraint method is computationally inexpensive. In the existing volume constraint method (the augmented Lagrangian method for constrained problems), its $Λ$ is determined from the previous optimization step, and hence the volume constraint at the current step is inexact. Such inexactness leads to a large fluctuation in the as-constrained volume fraction, which can jeopardize the appearance of minor geometrical features and hence alter the topology of the optimized structure. The exactness of our proposed volume constraint method can ensure convergence for the minimum compliance problem by only reconstructing and remeshing the material domain, whereas the existing volume constraint method fails to do so. The proposed exact volume constraint during the entire optimization process differs from merely adjusting the volume by choosing a level set other than zero. The latter is merely an engineering trick that neither affects the convergence nor corresponds to the obtained objective functional.
A BEM-based topology optimization for acoustic problems considering tangential derivative of sound pressure
2022, Computer Methods in Applied Mechanics and Engineering
Citation Excerpt :
Jiang et al. [45] propose the topology optimization method based on IGABEM and the piecewise constant level method for the distribution design of the sound absorbing materials. Gao et al. [46]derive the average topological derivative for the topology optimization problem for finite unidirectional phononic structures using the BEM and level set method. Oliveira et al. [47] employ IGABEM and the level set method for topology optimization.
This paper employs a new tangential derivative of boundary integral equation for the optimization problems in the acoustic field with a objective function involving tangential derivatives of the sound pressure on the boundary. The level set method is adopted to generate the topological structure by updating the level set function which defines the boundary of the material domain with its zero contour line. The hyper singular integral is directly derived and singular terms are canceled due to the form of the tangential derivative at the boundary. The topological derivative is derived through the adjoint variable method(AVM) and the most of the unknowns in the variation of objective function can be canceled by evaluating the adjoint field. However, one of the terms which includes the variation of the tangential derivative of the sound pressure is evaluated using integration by parts. The remaining part having unknown variation of sound pressure is neglected by extending the objective function defined boundary by one elements at its start and end points. Numerical implementations demonstrate the effectiveness and correctness of the proposed method for topology optimization problems with the objective function involving tangential derivative quantities.
Enhanced ultrasonic wave generation using energy-localized behaviors of phononic crystals
2022, International Journal of Mechanical Sciences
Citation Excerpt :
One mechanism of phononic band-gap formation is Bragg-scattering; its periodic nature brings about repeated destructive interferences [6]. Several research topics on the phononic band gap include the modeling [7], experimental validation [8], unit-cell designs for lowering [9] or broadening [10] frequency ranges, and design-optimization execution [11]. Promising engineering applications have been presented, including wave mitigation [12], isolation [13], and attenuation [14].
This paper aims to bridge the energy-localized behaviors of a phononic crystal (PnC) and an ultrasonic-wave-generation (UWG) system. A defective PnC manifests an intriguing phenomenon in which input ultrasonic-wave energy is confined near the defect in several defect-band frequencies that appear in the phononic band-gap range. Using the inverse of this energy localization, remarkably amplified ultrasonic waves can be generated compared to a UWG system without the PnC. This work proposes an analytical model that predicts the UWG system's sensitivities and electrical admittance. The mechanical equation of motions and the electrical circuit equation of the piezoelectric-patch-bonded defect are derived using the Rayleigh-Love rod theory and the extended Hamilton's principle. Through the use of Green's functions, these electroelastically coupled equations are solved. With two assumptions, specifically, (i) continuity conditions throughout the PnC and (ii) no incident waves toward the UWG system, the S-parameter-method-based analytical approach presents the system performances explicitly. Furthermore, the modified analytical model for application to engineering situations with finite boundary conditions is also proposed. The predicted results of the proposed analytical model match well with those of the finite element model and show that the UWG performances are amplified to several times their prior value. The main contributions of this work are as follows: (i) this is the first attempt to utilize a defective PnC as UWG systems, (ii) the analytical model is newly proposed and dramatically shortens the computational time compared to the finite-element-method, and (iii) the proposed defective-PnC-incorporated UWG system and its analytical model are effective, regardless of the boundary conditions.
Partitioned path-following strategy for nonlinear structural analyses using the boundary element method
2022, Computer Methods in Applied Mechanics and Engineering
Citation Excerpt :
The Boundary Element Method (BEM) is recognized as a robust and efficient numerical tool, characterized by a one-order reduction in mesh dimensionality and its ability to accurately reproduce high gradient stress fields [9–11]. It finds applications in several areas of engineering, including isogeometric analysis [12,13], biomechanics [14,15], fracture mechanics [16–19], acoustics [20–22], shape [23,24] and topology optimization [25–29]. Some reviews [30–32] provide more details on the development and application of the BEM.
The appearance of snapping phenomena is frequent in structural analyses involving nonlinearities such as the propagation of localized damage zones, or cracks, as encountered in quasi-brittle materials. In the presence of snapping effects, the standard Newton-based algorithms cannot predict the equilibrium curves correctly by usual displacement or force control. This problem can be circumvented by using path-following techniques. This study proposes a numerical approach for nonlinear structural analysis based on the Boundary Element Method (BEM) coupled with a partitioned path-following strategy. In this strategy, the acting loads are divided into two groups: (a) those perfectly known and (b) those in which only the direction is known, while the magnitude (or load factor/multiplier) is determined to respect the equilibrium of the boundary fields (indirectly controlled) and an additional problem equation, the path-following constraint equation. The resulting nonlinear system is solved using an incremental iterative scheme. For each iteration, the corrections to the boundary fields are obtained in a partitioned manner. First, the portion coming from the unbalanced boundary fields is calculated. Then the correction resulting from a unit load in the direction of the load to be indirectly controlled is obtained. Finally, the load factor is calculated independently using the path-following constraint equation, which can be linear, nonlinear, differentiable, or not. In the numerical tests, the Control of the Nodal Displacement Increment (CNDI) is chosen as path-following constraint equation for the sake of illustration. The results show that the proposed approach can efficiently capture the equilibrium curve even in the presence of severe snap-backs. The approach is designed to be implemented in existing BEM codes straightforwardly, provided that the influence matrices can be accessed. Moreover, no additional pre-conditioning technique is required.
Analytically-integrated radial integration PBEM for solving three-dimensional steady heat conduction problems
2021, Engineering Analysis with Boundary Elements
Citation Excerpt :
First, just the 2D boundary of the physical model should be discretized which facilitates the modeling of complex geometries. Second, semi-analytical solutions to problems make the results have high accuracy [13–16]. These merits promote the application of BEM in heat conduction problems, so BEM becomes one of the most popular numerical methods.
Recently, a polygonal boundary element method (PBEM) has been developed for solving heat conduction problems. In the PBEM, the radial integration method (RIM) is employed to shift the domain and surface integrals in the boundary-domain integral equation into equivalent line integrals, and all the radial integrals are computed by Gaussian quadrature, which may obtain an accurate solution. However, the computation costs much. Therefore, analytical expressions of radial integrals are derived in this paper, concerning four kinds of varying thermal conductivities and two kinds of heat generation functions, aiming at improving the efficiency. Then all the radial integrals in the boundary-domain integral equation can be analytically calculated. Three examples are employed to examine the analytically-integrated PBEM for solving steady heat conduction problems. The results show that the proposed method is accurate, and it is more efficient than traditional PBEM.

View all citing articles on Scopus

View full text

A level set based topology optimization for finite unidirectional acoustic phononic structures using boundary element method

Highlights

Abstract

Introduction

Section snippets

Level set method

BEM for finite unidirectional phononic structures

Verification of topological derivative

Conclusions

Declaration of Competing Interest

Acknowledgment

Comput. Methods Appl. Mech. Engrg.

Mech. Syst. Signal Process.

Results Phys.

Phys. Lett. A

Solid State Commun.

Phys. Lett. A

J. Mech. Phys. Solids

Int. J. Solids Struct.

Int. J. Solids Struct.

Physica E

Eng. Anal. Bound. Elem.

J. Mech. Phys. Solids

Int. J. Solids Struct.

Comput. Oper. Res.

Comput. Struct.

Comput. Methods Appl. Mech. Engrg.

Physica B

Wave Motion

Comput. Struct.

J. Comput. Phys.

Eng. Anal. Bound. Elem.

Eng. Anal. Bound. Elem.

Energy Procedia

J. Comput. Phys.

Appl. Acoust.

Comput. Methods Appl. Mech. Engrg.

Comput. Methods Appl. Mech. Engrg.

Comput. Methods Appl. Mech. Engrg.

Comput. Methods Appl. Mech. Engrg.

Comput. Methods Appl. Mech. Engrg.

Eng. Anal. Bound. Elem.

Comput. Methods Appl. Mech. Engrg.

Comput. Methods Appl. Mech. Engrg.