A new isogeometric density based approach for the topology optimization of multi-material structures is presented. In this method, the density fields of multiple material phases are represented using the isogeometric NURBS-based parameterization leading to exact modeling of the geometry, removing numerical artifacts and full analytical computation of sensitivities in a cost effective manner. An extension of the perimeter control technique is introduced where restrictions are imposed on the perimeters of density fields of all phases. Consequently, not only can one control the complexity of the optimal design but also the minimal lengths scales of all material phases. This leads to optimal designs with significantly enhanced manufacturability and comparable performance. Unlike the common element-wise or nodal based density representations, owing to higher order continuity of density fields in this method, their gradients required for perimeter control restrictions are calculated exactly without additional computational cost. The problem is formulated with constraints on either: (a) volume fractions of different material phases, or (b) the total mass of the structure. The proposed method is applied for the minimal compliance design of two-dimensional structures consisting of multiple distinct materials as well as functionally graded ones.

Digital imaging technologies such as X-ray scans and ultrasound provide a convenient and non-invasive way to capture high resolution images. The colour intensity of digital images provides information on the geometrical features and material distribution which can be utilised for stress analysis. The proposed approach employs an automatic and robust algorithm to generate quadtree (2D) or octree (3D) meshes from digital images. The use of polygonal elements (2D) or polyhedral elements (3D) constructed by the scaled boundary finite element method (SBFEM) avoids the issue of hanging nodes (mesh incompatibility) commonly encountered by finite elements on quadtree or octree meshes. The computational effort is reduced by considering the small number of cell patterns occurring in a quadtree or an octree mesh. Examples with analytical solutions in 2D and 3D are provided to show the validity of the approach. Other examples including the analysis of 2D and 3D microstructures of concrete specimens as well as of a domain containing multiple spherical holes are presented to demonstrate the versatility and the simplicity of the proposed technique. This article is protected by copyright. All rights reserved.

This paper proposes an efficient metamodeling approach for UQ of complex system based on Gaussian process model (GPM). The proposed GPM-based method is able to efficiently and accurately calculate the mean and variance of model outputs with uncertain parameters specified by arbitrary probability distributions. Due to the use of GPM, the closed form expressions of mean and variance can be derived by decomposing high-dimensional integrals into one-dimensional integrals. This paper details on how to efficiently compute the one-dimensional integrals. When the parameters are either uniformly- or normally distributed, the one-dimensional integrals can be analytically evaluated, while when parameters do not follow normal or uniform distributions, this paper adopts the effective Gaussian quadrature technique for the fast computation of the one-dimensional integrals. As a result, the developed GPM method is able to calculate mean and variance of model outputs in an efficient manner independent of parameter distributions. The proposed GPM method is applied to a collection of examples. And its accuracy and efficiency is compared with Monte Carlo simulation (MCS), which is used as benchmark solution. Results show that the proposed GPM method is feasible and reliable for efficient UQ of complex systems in terms of the computational accuracy and efficiency. Copyright © 2016 John Wiley & Sons, Ltd.

We propose a numerical method for a fluid-structure interaction problem. The material of the structure is homogeneous, isotropic and it can be described by the compressible Neo-Hookean constitutive equation, while the fluid is governed by the Navier-Stokes equations. Our study does not use turbulence model. Updated Lagrangian method is used for the structure and fluid equations are written in Arbitrary Lagrangian Eulerian coordinates. One global moving mesh is employed for the fluid-structure domain, where the fluid-structure interface is an “interior boundary” of the global mesh. At each time step, we solve a monolithic system of unknowns velocity and pressure defined on the global mesh. The continuity of velocity at the interface is automatically satisfied, while the continuity of stress does not appear explicitly in the monolithic fluid-structure system. This method is very fast because at each time step we solve only one linear system. This linear system was obtained by the linearization of the structure around the previous position in the updated Lagrangian formulation and by the employment of a linear convection term for the fluid. Numerical results are presented. This article is protected by copyright. All rights reserved.

The variational approach to fracture is effective for simulating the nucleation and propagation of complex crack patterns, but is computationally demanding. The model is a strongly nonlinear non-convex variational inequality that demands the resolution of small length scales. The current standard algorithm for its solution, alternate minimization, is robust but converges slowly and demands the solution of large, ill-conditioned linear subproblems. In this paper, we propose several advances in the numerical solution of this model that improve its computational efficiency. We reformulate alternate minimization as a nonlinear Gauss-Seidel iteration and employ over-relaxation to accelerate its convergence; we compose this accelerated alternate minimization with Newton's method, to further reduce the time to solution; and we formulate efficient preconditioners for the solution of the linear subproblems arising in both alternate minimization and in Newton's method. We investigate the improvements in efficiency on several examples from the literature; the new solver is 5–6× faster on a majority of the test cases considered. This article is protected by copyright. All rights reserved.

Accurate sizing functions are crucial for efficient generation of high-quality meshes, but to define the sizing function is often the bottleneck in complicated mesh generation tasks due to the tedious user interaction involved. We present a novel algorithm to automatically create high-quality sizing functions for surface mesh generation. First, the tessellation of a CAD model is taken as the background mesh, in which an initial sizing function is defined by considering geometrical factors and user-specified parameters. Then, a convex nonlinear programming problem is formulated and solved efficiently to obtain a smoothed sizing function that corresponds to a mesh satisfying necessary gradient constraint conditions and containing a significantly reduced element number. Finally, this sizing function is applied in an advancing front mesher. With the aid of a walk-through algorithm, an efficient sizing-value query scheme is developed. Meshing experiments of some very complicated geometry models are presented to demonstrate that the proposed sizing-function approach enables accurate and fully automatic surface mesh generation.

A new higher-order accurate method is proposed that combines the advantages of the classical *p*-version of the FEM on body-fitted meshes with embedded domain methods. A background mesh composed by higher-order Lagrange elements is used. Boundaries and interfaces are described implicitly by the level set method and are within elements. In the elements cut by the boundaries or interfaces, an automatic decomposition into higher-order accurate sub-elements is realized. Therefore, the zero level sets are detected and meshed in a first step which is called reconstruction. Then, based on the topological situation in the cut element, higher-order sub-elements are mapped to the two sides of the boundary or interface. The quality of the reconstruction and the mapping largely determines the properties of the resulting, automatically generated conforming mesh. It is found that optimal convergence rates are possible although the resulting sub-elements are not always well-shaped. Copyright © 2016 John Wiley & Sons, Ltd.

This article revisits a topology optimization design approach for micro-manufacturing and extends it to optical microlithography with partially coherent illumination. The solution is based on a combination of two technologies, the topology optimization and the proximity error correction in micro-/nano- lithography. The key steps include (i) modeling the physical inputs of the fabrication process, including the UV-light illumination source and the mask, as the design variables in optimization and (ii) applying physical filtering and heaviside projection for topology optimization, which correspond to the aerial image formulation and the pattern development processes, respectively. The proposed approach results in an effective source and a binary design mask, which can be sent directly to fabrication without additional post-processing steps for proximity error correction. Meanwhile, the performance of the device is optimized and robust w.r.t. process variations, such as dose / photo-resist variations and lens defocus. A compliant micro-gripper design example is considered to demonstrate the applicability of this approach. Copyright © 2016 John Wiley & Sons, Ltd.

This paper develops and demonstrates a strategy for identifying the lack of dynamic interaction, or coupling, between potential design modifications to a vibrating base structure. Such decoupling may be exploited to efficiently compare the performance of competing engineering designs. In particular, it is shown how different design modifications may be represented as the addition or removal of substructures. When these substructures are uncoupled according to the metric developed here, the computational cost of determining the optimal system design can be greatly reduced. For example, if a designer considers seven possible modifications and wishes to examine all possible combinations of the modifications, 128 possible structures must be analyzed. However, if all modifications are dynamically uncoupled, significant computational effort need only be spent on eight of the possible structures in order to generate the responses for all remaining designs. Example problems demonstrate this cost reduction and illustrate cases where dynamic uncoupling occurs. This article is protected by copyright. All rights reserved.

A non-gradient based approach for topology optimization using a genetic algorithm is proposed in this paper. The genetic algorithm used in this paper is assisted by the Kriging surrogate model to reduce computational cost required for function evaluation. To validate the non-gradient based topology optimization method in flow problems, this research focuses on two single-objective optimization problems, where the objective functions are to minimize pressure loss and to maximize heat transfer of flow channels, and one multi-objective optimization problem, which combines above two single-objective optimization problems. The shape of flow channels is represented by the level set function. The pressure loss and the heat transfer performance of the channels are evaluated by the Building-Cube Method (BCM) code, which is a Cartesian-mesh CFD solver. The proposed method resulted in an agreement with previous study in the single-objective problems in its topology, and achieved global exploration of non-dominated solutions in the multi-objective problems.

This paper builds on recent work developed by the authors for the numerical analysis of large strain solid dynamics, by introducing an upwind cell centred hexahedral Finite Volume framework implemented within the open source code OpenFOAM [http://www.openfoam.com/www.openfoam.com]. In Lee, Gil and Bonet [1], a first order hyperbolic system of conservation laws was introduced in terms of the linear momentum and the deformation gradient tensor of the system, leading to excellent behaviour in two dimensional bending dominated nearly incompressible scenarios. The main aim of this paper is the extension of this algorithm into three dimensions, its tailor-made implementation into OpenFOAM and the enhancement of the formulation with three key novelties. First, the introduction of two different strategies in order to ensure the satisfaction of the underlying involutions of the system, that is, that the deformation gradient tensor must be curl-free throughout the deformation process. Second, the use of a discrete angular momentum projection algorithm and a monolithic Total Variation Diminishing Runge-Kutta time integrator combined in order to guarantee the conservation of angular momentum. Third, and for comparison purposes, an adapted Total Lagrangian version of the Hyperelastic-GLACE nodal scheme of Kluth and Després [2] is presented. A series of challenging numerical examples are examined in order to assess the robustness and accuracy of the proposed algorithm, benchmarking it against an ample spectrum of alternative numerical strategies developed by the authors in recent publications. Copyright © 2016 John Wiley & Sons, Ltd.

By separation of scales and the homogenization of a flow through porous media, a two-scale problem arises where a Darcy-type flow is present on the macroscale and a Stokes flow on the subscale. In this paper, the problem is given as the minimization of a potential. Additional constraints imposing periodicity in a weak sense are added using Lagrange multipliers. In particular, the upper and lower energy bounds for the corresponding strongly periodic problem are produced, quantifying the accuracy of the weakly periodic boundary conditions. A numerical example demonstrates the evaluation of energy bounds and the performance of weakly periodic boundary conditions on an Representative Volume Element. Copyright © 2016 John Wiley & Sons, Ltd.

The contribution of this work is the implementation of a new elastic solution method for thick laminated composites and sandwich structures based on a Generalized Unified Formulation (GUF) using Finite Elements. A quadrilateral 4-node element was developed and evaluated using in-house Finite Element program. The C-1 continuity requirements are fulfilled for the transversal displacement field variable. This method is tagged as Caliri's Generalized Formulation (CGF). The results employing the proposed solution method yielded coherent results with deviations as low as 0.05% for a static simply supported symmetric laminate and 0.5% for the modal analyses of a soft core sandwich structure.

A new nonconforming brick element is introduced, which only has 13 degrees of freedom locally and takes as its shape functions space. The vector-valued version generates, together with a discontinuous approximation, an inf-sup stable finite element pair of order 2 for the Stokes problem in the energy norm. This article is protected by copyright. All rights reserved.

This work discusses a discontinuous Galerkin (DG) discretization for two-phase flows. The fluid interface is represented by a level-set and the DG approximation space is adapted such that jumps and kinks in pressure and velocity fields can be approximated sharply. This adaption of the original DG space, which can be performed “on-the-fly” for arbitrary interface shapes, is referred to as extended discontinuous Galerkin (XDG). By combining this ansatz with a special quadrature technique, one can regain spectral convergence properties for low-regularity solutions, which is demonstrated by numerical examples. This work focuses on the aspects of spatial discretization, and special emphasis is devoted on how to overcome problems related to quadrature, small cut-cells and condition number of linear systems. Temporal discretization will be discussed in future works. This article is protected by copyright. All rights reserved.

A robust and efficient strategy is proposed to simulate mechanical problems involving cohesive fractures. This class of problems is characterized by a global structural behavior that is strongly affected by localized non-linearities at relatively small-sized critical regions. The proposed approach is based on the division of a simulation into a suitable number of sub-simulations where adaptive mesh refinement is performed only once based on refinement window(s) around crack front process zone(s). The initialization of Newton-Raphson non-linear iterations at the start of each sub-simulation is accomplished by solving a linear problem based on a secant stiffness, rather than a volume mapping of non-linear solutions between meshes. The secant stiffness is evaluated using material state information stored/read on crack surface facets which are employed to explicitly represent the geometry of the discontinuity surface independently of the volume mesh within the generalized finite element method framework. Moreover, a simplified version of the algorithm is proposed for its straightforward implementation into existing commercial software. Data transfer between sub-simulations is not required in the simplified strategy. The computational efficiency, accuracy, and robustness of the proposed strategies are demonstrated by an application to cohesive fracture simulations in 3-D. Copyright © 2016 John Wiley & Sons, Ltd.

Three new contributions to the field of multisurface plasticity are presented for general situations with an arbitrary number of nonlinear yield surfaces with hardening or softening. A method for handling linearly-dependent flow directions is described. A residual that can be used in a line search is defined. An algorithm that has been implemented and comprehensively tested is discussed in detail.

Examples are presented to illustrate the computational cost of various components of the algorithm. The overall result is that a single Newton-Raphson iteration of the algorithm costs between 1.5 and 2 times that of an elastic calculation.

Examples also illustrate the successful convergence of the algorithm in complicated situations. For example, without using the new contributions presented here, the algorithm fails to converge for approximately 50% of the trial stresses for a common geomechanical model of sedementary rocks, while the current algorithm results in complete success.

Since it involves no approximations, the algorithm is used to quantify the accuracy of an efficient, pragmatic, but approximate, algorithm used for sedimentary-rock plasticity in a commercial software package.

The main weakness of the algorithm is identified as the difficulty of correctly choosing the set of initially active constraints in the general setting. This article is protected by copyright. All rights reserved.

Reduced order models are useful for accelerating simulations in many-query contexts, such as optimization, uncertainty quantification, and sensitivity analysis. However, offline training of reduced order models can have prohibitively expensive memory and floating-point operation costs in high-performance computing applications, where memory per core is limited. To overcome this limitation for proper orthogonal decomposition, we propose a novel adaptive selection method for snapshots in time that limits offline training costs by selecting snapshots according an error control mechanism similar to that found in adaptive time-stepping ordinary differential equation solvers. The error estimator used in this work is related to theory bounding the approximation error in time of proper orthogonal decomposition-based reduced order models, and memory usage is minimized by computing the singular value decomposition using a single-pass incremental algorithm. Results for a viscous Burgers' test problem demonstrate convergence in the limit as the algorithm error tolerances go to zero; in this limit, the full order model is recovered to within discretization error. The resulting method can be used on supercomputers to generate proper orthogonal decomposition-based reduced order models, or as a subroutine within hyperreduction algorithms that require taking snapshots in time, or within greedy algorithms for sampling parameter space. Copyright © 2016 John Wiley & Sons, Ltd.

This paper describes a novel methodology that combines smoothed discrete particle hydrodynamics (SDPH) and finite volume method (FVM) to enhance the effective performance in solving the problems of gas-particle multiphase flow. To describe the collision and fluctuation of particles, this method also increases a new parameter, namely, granular temperature, according to the kinetic theory of granular flow. The coupled framework of SDPH–FVM has been established, in which the drag force and pressure gradient act on the SDPH particles and the momentum sources of drag force are added back onto the FVM mesh. The proposed technique is a coupled discrete–continuum method based on the two-fluid model. To compute for the discrete phase, its SDPH is developed from smoothed particle hydrodynamics (SPH), in which the properties of SPH are redefined with some new physical quantities added into the traditional SPH parameters, so that it is more beneficial for SDPH in representing the particle characteristics. For the continuum phase, FVM is employed to discretize the continuum flow field on a stationary grid by capturing fluid characteristics. The coupled method exhibits strong efficiency and accuracy in several 2-D numerical simulations. This article is protected by copyright. All rights reserved.

In this paper, we present a generalized prismatic hybrid meshing method for viscous flow simulations. One major difficulty in implementing a robust prismatic hybrid meshing tool is to handle boundary layer mesh collisions and normally an extra data structure (e.g. quadtree in 2D and octree in 3D) is required. The proposed method overcomes this difficulty via an heuristic approach and it only relies on Constrained Delaunay Triangulation/Tetrahedralization(CDT). No extra data structures are required. Geometrical reasoning is used to approximate the maximum marching distance of each point by walking through the CDT. This is combined with post-processing of marching vectors and distance and prohibition of multilevel differences to form an automatic and robust mechanism to remove boundary layer mesh collisions. Benefiting from the matureness of CDT techniques, the proposed method is robust, efficient and simple to implement. Its capability is demonstrated by generating quality prismatic hybrid meshes for industrial models with complex geometries. The proposed method is believed to be able considerably reduce the effort to implement a robust hybrid prismatic mesh generator for viscous flow simulations. This article is protected by copyright. All rights reserved.

A computational framework for assisting in the development of novel textiles is presented. Electronic textiles are key in the rapidly growing field of wearable electronics for both consumer and military uses. There are two main challenges to the modeling of electronic textiles: the discretization of the textile microstructure and the interaction between electromagnetic and mechanical fields. A director-based beam formulation with an assumed electrical current is used to discretize the fabric at the level of individual fibrils. The open source package FEniCS was used to implement the finite element model. Contact integrals were added into the FEniCS framework so that multiphysics contact laws can be incorporated in the same framework, leveraging the code generation and automated differentiation capabilities of FEniCS to produce the tangents needed by the implicit solution method. The computational model is used to construct and determine the mechanical, thermal, and electrical properties of a representative volume elements of a plain woven textile. Dynamic relaxation to solve the mechanical fields, and the electrical and thermal fields are solved statically for a given mechanical state. The simulated electrical responses are fit to a simplified Kirchhoff network model to determine effective resistances of the textile. Copyright © 2016 John Wiley & Sons, Ltd.

It is common in solving topology optimization problems to replace an integer-valued characteristic function design field with the material volume fraction field, a real-valued approximation of the design field that permits “fictitious” mixtures of materials during intermediate iterations in the optimization process. This is reasonable so long as one can interpolate properties for such materials and so long as the final design is integer valued. For this purpose, we present a method for smoothly thresholding the volume fractions of an arbitrary number of material phases which specify the design. This method is trivial for two-material design problems, for example the canonical topology design problem of specifying the presence or absence of a single material within a domain, but it becomes more complex when three or more materials are used, as often occurs in material design problems. We take advantage of the similarity in properties between the volume fractions and the barycentric coordinates on a simplex to derive a thresholding method which is applicable to an arbitrary number of materials. As we show in a sensitivity analysis, this method has smooth derivatives, allowing it to be used in gradient-based optimization algorithms. We present results which show synergistic effects when used with SIMP and RAMP material interpolation functions, popular methods of ensuring integerness of solutions. This article is protected by copyright. All rights reserved.

In this paper we develop a dual-horizon peridynamics (DH-PD) formulation that naturally includes varying horizon sizes and completely solves the “ghost force” issue. Therefore, the concept of dual-horizon is introduced to consider the unbalanced interactions between the particles with different horizon sizes. The present formulation fulfills both the balances of linear momentum and angular momentum exactly. Neither the “partial stress tensor” nor the “slice” technique are needed to ameliorate the ghost force issue in [1]. We will show that the traditional peridynamics can be derived as a special case of the present DH-PD. All three peridynamic formulations, namely bond based, ordinary state based and non-ordinary state based peridynamics can be implemented within the DH-PD framework. Our DH-PD formulation allows for *h*-adaptivity and can be implemented in any existing peridynamics code with minimal changes. A simple adaptive refinement procedure is proposed reducing the computational cost. Both two- and three- dimensional examples including the Kalthoff-Winkler experiment and plate with branching cracks are tested to demonstrate the capability of the method. This article is protected by copyright. All rights reserved.

In this paper, a 4-node quadrilateral flat shell element is proposed for geometrically nonlinear analysis based on updated Lagrangian formulation with the co-rotational kinematics concept. The flat shell element combines the membrane element with drilling degrees of freedom and the plate element with shear deformation. By means of these linearized elements, a simplified nonlinear analysis procedure allowing for warping of the flat shell element and large rotation is proposed. The tangent stiffness matrix and the internal force recovery are formulated in this paper. Several classic benchmark examples are presented to validate the accuracy and efficiency of the proposed new and more proficient element for practical engineering analysis of shell structures. This article is protected by copyright. All rights reserved.

In this work a mixed variational formulation to simulate quasi-incompressible electro- or magneto-active polymers immersed in the surrounding free space is presented. A novel domain decomposition is used to disconnect the primary coupled problem and the arbitrary free space mesh update problem. Exploiting this decomposition we describe a block iterative approach to solving the linearised multiphysics problem, and a physically and geometrically based, three-parameter method to update the free space mesh. Several application-driven example problems are implemented to demonstrate the robustness of the mixed formulation for both electro-elastic and magneto-elastic problems involving both finite deformations and quasi-incompressible media. This article is protected by copyright. All rights reserved.

Uniform grid solvers of the periodic Lippmann–Schwinger equation have been introduced by Moulinec and Suquet for the numerical homogenization of heterogeneous materials. Based on the fast Fourier transform, these methods use the strain as main unknown and usually do not produce displacement fields. While this is generally not perceived as a restriction for homogenization purposes, some tasks might require kinematically admissible displacement fields. In this paper, we show how the numerical solution to the periodic Lippmann–Schwinger equation can be post-processed to reconstruct a displacement field. Our procedure applies to any variant of the Moulinec–Suquet solver. The reconstruction is formulated as an auxiliary elastic equilibrium problem of a homogeneous material, which is solved with displacement-based finite elements. Taking advantage of periodicity, uniformity of the grid and homogeneity of the material, the resulting linear system is formulated and solved efficiently in Fourier space. The cost of our procedure is lower than that of one iteration of the Lippmann–Schwinger solver. Two applications are proposed, in two and three dimensions. In the first application, the reconstructed displacement field is used to compute a rigorous upper bound on the effective shear modulus. In the second application, the quality of the reconstruction is assessed quantitatively. Copyright © 2016 John Wiley & Sons, Ltd.

Mechanical computations in multiphase domains raise numerous difficulties from the generation of the initial mesh to its adaptation throughout the simulation. All alternatives to mesh adaptation, such as level-set methods, have the well-known drawback of inducing volume conservation issues. In this paper, a moving mesh method is coupled to a topological mesh adaptation technique in order to track moving and deforming interfaces in multiphase simulations, with a robust control of mesh quality. Level-set functions are used as intermediaries to enhance the mesh adaptation technique with a volume conservation constraint, which is compatible both with implicit and with body-fitted interfaces. Results show that this method has the same advantage of permitting important displacements, deformations, and topological changes (coalescence of interfaces, for example) as a standard level-set method, while volume diffusion is drastically reduced. Copyright © 2016 John Wiley & Sons, Ltd.

This study concerns the development of a numerical methodology for initializing immersed interface-based CFD solvers for using complex computer-aided design (CAD) geometry. CFD solvers with higher-order discretization stencils require larger stencil widths, which become problematic in regions of space where insufficient mesh resolution is available. This problem becomes especially challenging when convoluted triangulated surface meshes generated from complex solid models are used to initialize the cut-cells. A pragmatic balance between desired local geometry resolution and numerical accuracy is often required to find a practical solution. Here, a robust iterative fill algorithm is presented that balances geometry resolution with numerical accuracy (*via* stencil size). Several examples are presented to illustrate the use of this initialization procedure that employs both the original CAD generated triangulated surface mesh, along with a level set representation of the surface to initialize cut-cells and boundary proximity measures for creation of CFD stencils. Convergence error analysis of surface area and enclosed volumes is first presented to show the effects of fill on the geometry as a function of desired stencil size and grid resolution. The algorithm is then applied to geometrically complex problems using large eddy simulation. Two problems are considered. The first is flow around the Eiffel Tower. The second is a combustion swirler in the context of a design problem. Copyright © 2016 John Wiley & Sons, Ltd.

A comprehensive study of the two sub-steps composite implicit time integration scheme for the structural dynamics is presented in this paper. A framework is proposed for the convergence accuracy analysis of the generalized composite scheme. The local truncation errors of the acceleration, velocity, and displacement are evaluated in a rigorous procedure. The presented and proved accuracy condition enables the displacement, velocity, and acceleration achieving second-order accuracy simultaneously, which avoids the drawback that the acceleration accuracy may not reach second order. The different influences of numerical frequencies and time step on the accuracy of displacement, velocity, and acceleration are clarified. The numerical dissipation and dispersion and the initial magnitude errors are investigated physically, which measure the errors from the algorithmic amplification matrix's eigenvalues and eigenvectors, respectively. The load and physically undamped/damped cases are naturally accounted. An optimal algorithm-Bathe composite method (Bathe and Baig, 2005; Bathe, 2007; Bathe and Noh, 2012) is revealed with unconditional stability, no overshooting in displacement, velocity, and acceleration, and excellent performance compared with many other algorithms. The proposed framework also can be used for accuracy analysis and design of other multi-sub-steps composite schemes and single-step methods under physical damping and/or loading. Copyright © 2016 John Wiley & Sons, Ltd.

This contribution presents a numerical strategy to evaluate the effective properties of image-based microstructures in the case of random material properties. The method relies on three points: (1) a high-order fictitious domain method; (2) an accurate spectral stochastic model; and (3) an efficient model-reduction method based on the proper generalized decomposition in order to decrease the computational cost introduced by the stochastic model. A feedback procedure is proposed for an automatic estimation of the random effective properties with a given confidence. Numerical verifications highlight the convergence properties of the method for both deterministic and stochastic models. The method is finally applied to a real 3D bone microstructure where the empirical probability density function of the effective behaviour could be obtained. Copyright © 2016 John Wiley & Sons, Ltd.

In this contribution, a novel method for a fail-safe optimal design of structures is proposed, which is a coupled approach of optimization employing a genetic algorithm, the structural analysis conducted in the framework of fracture mechanics and uncertainty analysis. The idea of fail-safe structures is to keep their functionality and integrity even under damage conditions, for example, a local failure of substructures. In the present work, a design concept of a substructure exhibiting a damage accumulating function due to the application of crack arresters is introduced. If such a substructure is integrated within a system of coupled substructures, it will accumulate the damage arising from the boundary conditions change induced by the failure of certain neighbouring structural elements and hinder further damage escalation. The investigation of failure of the damage accumulating substructure is introduced within a finite element framework by a combination of discrete fracturing and configurational mechanics based criteria. In order to design a structure, which will fail safely according to a predefined scenario, uncertainties are taken into account. The developed approach optimizes the configuration of crack arresters within the damage accumulating substructure so that the uncertain crack propagation is hindered and only a local failure of this element occurs. Copyright © 2016 John Wiley & Sons, Ltd.

In this paper, a multibody meshfree framework is proposed for the simulation of granular materials undergoing deformations at the grain scale. This framework is based on an implicit solving of the mechanical problem based on a weak form written on the domain defined by all the grains composing the granular sample. Several technical choices, related to the displacement field interpolation, to the contact modelling, and to the integration scheme used to solve the dynamic equations, are explained in details. A first implementation is proposed, under the acronym Multibody ELement-free Open code for DYnamic simulation (MELODY), and is made available for free download. Two numerical examples are provided to show the convergence and the capability of the method. Copyright © 2016 John Wiley & Sons, Ltd.

Parametric and implicit methods are traditionally thought to be two irrelevant approaches in structural shape optimization. Parametric method works as a Lagrangian approach and often uses the parametric boundary representation (B-rep) of curves/surfaces, for example, Bezier and B-splines in combination with the conformal mesh of a finite element model, while implicit method relies upon level-set functions, that is, implicit functions for B-rep, and works as an Eulerian approach in combination with the fixed mesh within the scope of extended finite element method or finite cell method. The original contribution of this work is the unification of both methods. First, a new shape optimization method is proposed by combining the features of the parametric and implicit B-reps. Shape changes of the structural boundary are governed by parametric B-rep on the fixed mesh to maintain the merit in computer-aided design modeling and avoid laborious remeshing. Second, analytical shape design sensitivity is formulated for the parametric B-rep in the framework of fixed mesh of finite cell method by means of the Hamilton–Jacobi equation. Numerical examples are solved to illustrate the unified methodology. Copyright © 2016 John Wiley & Sons, Ltd.

This paper studies the behaviour of the error incurred when numerically integrating the elasto-plastic mechanical relationships of a constitutive model for soils using an explicit substepping formulation with automatic error control. The correct update of all the variables involved in the numerical integration of the incremental stress–strain relationships is central to the computational performance of the integration scheme, and, although often missed in the literature, all variables (including specific volume) should be explicitly considered in the algorithmic formulation. This is demonstrated in the paper by studying, in the context of the Cam clay formulations for saturated soils, the influence that the updating of the specific volume has on the accuracy of the numerical solution. The fact that the variation of the local error with the size of the integrated strain depends on the order of local accuracy of the numerical method is also used in the paper to propose a simple and powerful strategy to verify the correctness of the implemented mathematical formulation. Copyright © 2016 John Wiley & Sons, Ltd.

An efficient method is proposed to estimate the first-order global sensitivity indices based on failure probability and variance by using maximum entropy theory and Nataf transformation. The computational cost of this proposed method is quite small, and the proposed method can efficiently overcome the ‘dimensional curse’ due to dimensional reduction technique. Ideas for the estimation of higher-order sensitivity indices are discussed. Copyright © 2016 John Wiley & Sons, Ltd.

This paper describes a novel method for mapping between basis representation of a field variable over a domain in the context of numerical modelling and inverse problems. In the numerical solution of inverse problems, a continuous scalar or vector field over a domain may be represented in different finite-dimensional basis approximations, such as an unstructured mesh basis for the numerical solution of the forward problem, and a regular grid basis for the representation of the solution of the inverse problem. Mapping between the basis representations is generally lossy, and the objective of the mapping procedure is to minimise the errors incurred. We present in this paper a novel mapping mechanism that is based on a minimisation of the *L*^{2} or *H*^{1} norm of the difference between the two basis representations. We provide examples of mapping in 2D and 3D problems, between an unstructured mesh basis representative of an FEM approximation, and different types of structured basis including piecewise constant and linear pixel basis, and blob basis as a representation of the inverse basis. A comparison with results from a simple sampling-based mapping algorithm shows the superior performance of the method proposed here. © 2016 The Authors. International Journal for Numerical Methods in Engineering Published by John Wiley & Sons Ltd.

In finite element analysis of volume coupled multiphysics, different meshes for the involved physical fields are often highly desirable in terms of solution accuracy and computational costs. We present a general methodology for volumetric coupling of different meshes within a monolithic solution scheme. A straightforward collocation approach is compared to a mortar-based method for nodal information transfer. For the latter, dual shape functions based on the biorthogonality concept are used to build the projection matrices, thus further reducing the evaluation costs. We give a detailed explanation of the integration scheme and the construction of dual shape functions for general first-order and second-order Langrangian finite elements within the mortar method, as well as an analysis of the conservation properties of the projection operators. Moreover, possible incompatibilities due to different geometric approximations of curved boundaries are discussed. Numerical examples demonstrate the flexibility of the presented mortar approach for arbitrary finite element combinations in two and three dimensions and its applicability to different multiphysics coupling scenarios. Copyright © 2016 John Wiley & Sons, Ltd.

This paper presents an eight-node nonlinear solid-shell element for static problems. The main goal of this work is to develop a solid-shell formulation with improved membrane response compared with the previous solid-shell element (MOS2013), presented in . Assumed natural strain concept is implemented to account for the transverse shear and thickness strains to circumvent the curvature thickness and transverse shear locking problems. The enhanced assumed strain approach based on the Hu–Washizu variational principle with six enhanced assumed strain degrees of freedom is applied. Five extra degrees of freedom are applied on the in-plane strains to improve the membrane response and one on the thickness strain to alleviate the volumetric and Poisson's thickness locking problems. The ensuing element performs well in both in-plane and out-of-plane responses, besides the simplicity of implementation. The element formulation yields exact solutions for both the membrane and bending patch tests. The formulation is extended to the geometrically nonlinear regime using the corotational approach, explained in . Numerical results from benchmarks show the robustness of the formulation in geometrically linear and nonlinear problems. Copyright © 2016 John Wiley & Sons, Ltd.

Reliability analysis with both aleatory and epistemic uncertainties is investigated in this paper. The aleatory uncertainties are described with random variables, and epistemic uncertainties are tackled with evidence theory. To estimate the bounds of failure probability, several methods have been proposed. However, the existing methods suffer the dimensionality challenge of epistemic variables. To get rid of this challenge, a so-called random-set based Monte Carlo simulation (RS-MCS) method derived from the theory of random sets is offered. Nevertheless, RS-MCS is also computational expensive. So an active learning Kriging (ALK) model that only rightly predicts the sign of performance function is introduced and closely integrated with RS-MCS. The proposed method is termed as ALK-RS-MCS. ALK-RS-MCS accurately predicts the bounds of failure probability using as few function calls as possible. Moreover, in ALK-RS-MCS, an optimization method based on Karush–Kuhn–Tucker conditions is proposed to make the estimation of failure probability interval more efficient based on the Kriging model. The efficiency and accuracy of the proposed approach are demonstrated with four examples. Copyright © 2016 John Wiley & Sons, Ltd.

This paper deals with the numerical analysis of instabilities for elastic-plastic materials undergoing large deformations in non-isothermal conditions. The considered isotropic model is fully thermomechanically coupled and includes temperature-induced softening, which is another source of strain localization next to geometrical effects. Due to complexity of the model, a symbolic-numerical tool *Ace* is used for the preparation of user-supplied subroutines for the finite element method. The computational verification is performed using two benchmark tests: necking of circular bar in tension and shear banding of elongated rectangular plate in plain strain conditions. The attention is focused on mesh dependence of the numerical results and the regularizing effect of heat conduction. The research reveals that the conductivity influences the shear band width and ductility of the material response; however, for the adiabatic case, the results are discretization sensitive, and another regularization is needed. A new gradient-enhanced thermomechanical model is developed that introduces an internal length parameter governing the size of the shear band caused by thermal softening. The numerical verification of the non-local model is performed for the adiabatic case. Subsequently, the simultaneous application of the gradient enhancement and heat conduction in the model is analyzed, which reproduces an evolving shear band. Copyright © 2016 John Wiley & Sons, Ltd.

Computational aspects of a recently developed gradient elasticity model are discussed in this paper. This model includes the (Aifantis) strain gradient term along with two higher-order acceleration terms (micro-inertia contributions). It has been demonstrated that the presence of these three gradient terms enables one to capture the dispersive wave propagation with great accuracy. In this paper, the discretisation details of this model are thoroughly investigated, including both discretisation in time and in space. Firstly, the critical time step is derived that is relevant for conditionally stable time integrators. Secondly, recommendations on how to choose the numerical parameters, primarily the element size and time step, are given by comparing the dispersion behaviour of the original higher-order continuum with that of the discretised medium. In so doing, the accuracy of the discretised model can be assessed a priori depending on the selected discretisation parameters for given length-scales. A set of guidelines can therefore be established to select optimal discretisation parameters that balance computational efficiency and numerical accuracy. These guidelines are then verified numerically by examining the wave propagation in a one-dimensional bar as well as in a two-dimensional example. Copyright © 2016 John Wiley & Sons, Ltd.

In this article, a black-box higher-order fast multipole method for solving boundary integral equations on parametric surfaces in three spatial dimensions is proposed. Such piecewise smooth surfaces are the topic of recent studies in isogeometric analysis. Due to the exact surface representation, the rate of convergence of higher-order methods is not limited by approximation errors of the surface. An element-wise clustering strategy yields a balanced cluster tree and an efficient numerical integration scheme for the underlying Galerkin method. By performing the interpolation for the fast multipole method directly on the reference domain, the cost complexity in the polynomial degree is reduced by one order. This gain is independent of the application of either - or -matrices. In fact, several simplifications in the construction of -matrices are pointed out, which are a by-product of the surface representation. Extensive numerical examples are provided in order to quantify and qualify the proposed method. Copyright © 2016 John Wiley & Sons, Ltd.

A differential quadrature hierarchical finite element method (DQHFEM) is proposed by expressing the hierarchical finite element method matrices in similar form as in the differential quadrature finite element method and introducing interpolation basis on the boundary of hierarchical finite element method elements. The DQHFEM is similar as the fixed interface mode synthesis method but the DQHFEM does not need modal analysis. The DQHFEM with non-uniform rational B-splines elements were shown to accomplish similar destination as the isogeometric analysis. Three key points that determine the accuracy, efficiency and convergence of DQHFEM were addressed, namely, (1) the Gauss–Lobatto–Legendre points should be used as nodes, (2) the recursion formula should be used to compute high-order orthogonal polynomials, and (3) the separation variable feature of the basis should be used to save computational cost. Numerical comparison and convergence studies of the DQHFEM were carried out by comparing the DQHFEM results for vibration and bending of Mindlin plates with available exact or highly accurate approximate results in literatures. The DQHFEM can present highly accurate results using only a few sampling points. Meanwhile, the order of the DQHFEM can be as high as needed for high-frequency vibration analysis. Copyright © 2016 John Wiley & Sons, Ltd.

This manuscript presents the development of novel high-order complete shape functions over star-convex polygons based on the scaled boundary finite element method. The boundary of a polygon is discretised using one-dimensional high order shape functions. Within the domain, the shape functions are analytically formulated from the equilibrium conditions of a polygon. These standard scaled boundary shape functions are augmented by introducing additional bubble functions, which renders them high-order complete up to the order of the line elements on the polygon boundary. The bubble functions are also semi-analytical and preserve the displacement compatibility between adjacent polygons. They are derived from the scaled boundary formulation by incorporating body force modes. Higher-order interpolations can be conveniently formulated by simultaneously increasing the order of the shape functions on the polygon boundary and the order of the body force mode. The resulting stiffness-matrices and mass-matrices are integrated numerically along the boundary using standard integration rules and analytically along the radial coordinate within the domain. The bubble functions improve the convergence rate of the scaled boundary finite element method in modal analyses and for problems with non-zero body forces. Numerical examples demonstrate the accuracy and convergence of the developed approach. Copyright © 2016 John Wiley & Sons, Ltd.

In this contribution, we propose a dynamic gradient damage model as a phase-field approach for studying brutal fracture phenomena in quasi-brittle materials under impact-type loading conditions. Several existing approaches to account for the tension–compression asymmetry of fracture behavior of materials are reviewed. A better understanding of these models is provided through a uniaxial traction experiment. We then give an efficient numerical implementation of the model in an explicit dynamics context. Simulations results obtained with parallel computing are discussed both from a computational and physical point of view. Different damage constitutive laws and tension–compression asymmetry formulations are compared with respect to their aptitude to approximate brittle fracture. Copyright © 2016 John Wiley & Sons, Ltd.

Many challenging engineering and scientific problems involve the response of nonlinear solid materials to high-rate dynamic loading. Accompanying hydrodynamic effects are crucial, where shock-driven pressures may dominate material response. In this work, a hydrodynamic meshfree formulation is developed under the Lagrangian reproducing kernel particle method framework. The volumetric stress divergence is enhanced using a Rankine–Hugoniot-enriched Riemann solution that introduces the essential physics; oscillation control is introduced through appropriate state and field variable approximations that define the Riemann problem initial conditions. Consequently, non-physical numerical parameters and length scales required in the traditional artificial viscosity technique for shock modeling are avoided here. Several numerical examples are provided to verify the formulation accuracy across a range of shock loading conditions. Copyright © 2016 John Wiley & Sons, Ltd.

The production of new composite laminates with variable stiffness within the surface of plies was enabled by tow-placement machines. Because of the variation of stiffness, these materials are called variable stiffness composite laminates (VSCL). Recently, many attempts were made to investigate their structural behaviour. In this contribution, a first-order shear deformation theory is selected to model the multilayered composite laminates. The adopted theory is enhanced by the extended finite element method (XFEM) to describe discontinuities at element level of any interface of interest. To predict the location of the delamination onset, a traction–separation law is developed that is consistent with the XFEM topology. An exponential softening behaviour is implemented within the interface to model the delamination growth in a mixed-mode direction. In order to solve the non-linear equations of the delamination propagation, an arc-length method is applied. The effect of the curvilinear fibre orientation on the location of the delamination onset is investigated. Subsequently, the structural response of the laminates is computed. According to the simplicity of the new approach using the XFEM; and based on the computational cost for calculating the stiffness of VSCL, the method is able to determine structural response of VSCL with less computational effort. Copyright © 2016 John Wiley & Sons, Ltd.

Modal derivative is an approach to compute a reduced basis for model order reduction of large-scale nonlinear systems that typically stem from the discretization of partial differential equations. In this way, a complex nonlinear simulation model can be integrated into an optimization problem or the design of a controller, based on the resulting small-scale state-space model. We investigate the approximation properties of modal derivatives analytically and thus lay a theoretical foundation of their use in model order reduction, which has been missing so far. Concentrating on the application field of structural mechanics and structural dynamics, we show that the concept of modal derivatives can also be applied as nonlinear extension of the Craig–Bampton family of methods for substructuring. We furthermore generalize the approach from a pure projection scheme to a novel reduced-order modeling method that replaces all nonlinear terms by quadratic expressions in the reduced state variables. This complexity reduction leads to a frequency-preserving nonlinear quadratic state-space model. Numerical examples with carefully chosen nonlinear model problems and three-dimensional nonlinear elasticity confirm the analytical properties of the modal derivative reduction and show the potential of the proposed novel complexity reduction methods, along with the current limitations. Copyright © 2016 John Wiley & Sons, Ltd.

Computational modelling of fracture has been attempted in the past with a range of numerical approaches including finite element, extended finite element and meshless methods. The cracking particle method (CPM) of Rabczuk is a pragmatic alternative to explicit modelling of crack surfaces in which a crack is represented by a set of cracking particles that can be easily updated when the crack propagates. The change of cracking angle is recorded in discrete segments of broken lines, which makes this methodology suitable to model discontinuous cracks. In this paper, a new CPM is presented that improves on two counts: firstly, crack path curvature modelling is improved by the use of bilinear segments centred at each particle and secondly, efficiency for larger problems is improved via an adaptive process of both refinement and recovery. The system stiffness is calculated and stored in local matrices, so only a small influenced domain should be recalculated for each step while the remainder can be read directly from storage, which greatly reduces the computational expense. The methodology is applied to several 2D crack problems, and good agreement to analytical solutions and previous work is obtained. Copyright © 2016 John Wiley & Sons, Ltd.

This paper describes a selective mass scaling method which is designed for the analysis of wave propagation problems in nearly incompressible materials. The incompressibility of materials leads to a high value of the compressional wave speed, which makes the time step extremely small in explicit time integration method. The proposed selective mass scaling method selects the eigenfrequencies related to volumetric deformation modes to decrease them, while it keeps the shear eigenmodes unchanged. This makes the time step no longer limited by the compressional wave speed but by the shear wave speed. A significant reduction of CPU time is obtained with a good accuracy for transient problems in small strains on free or largely prestressed media. Copyright © 2016 John Wiley & Sons, Ltd.

This paper presents a 3D formulation for quasi-kinematic limit analysis, which is based on a radial point interpolation meshless method and numerical optimization. The velocity field is interpolated using radial point interpolation shape functions, and the resulting optimization problem is cast as a standard second-order cone programming problem. Because the essential boundary conditions can be only guaranteed at the position of the nodes when using radial point interpolation, the results obtained with the proposed approach are not rigorous upper bound solutions. This paper aims to improve the computing efficiency of 3D upper bound limit analysis and large problems, with tens of thousands of nodes, can be solved efficiently. Five numerical examples are given to confirm the effectiveness of the proposed approach with the von Mises yield criterion: an internally pressurized cylinder; a cantilever beam; a double-notched tensile specimen; and strip, square and rectangular footings. Copyright © 2016 John Wiley & Sons, Ltd.

The numerical errors was used to verify the correctness of key results. The truncation errors, which are larger than the round-off errors by orders of magnitude, have a superlinear relationship with both the simulation time-step and the interparticle collision speed. This remains the case regardless of the simulation details including the chosen contact model, particle size distribution, particle density or stiffness. Hence, the total errors can usually be reduced by choosing a smaller time-step. Increasing the polydispersity in a simulation by including smaller particles necessitates choosing a smaller time-step to maintain simulation stability and reduces the truncation errors in most cases. The truncation errors are increased by the dissipation of energy by frictional sliding or by the inclusion of damping in the system. The number of contacts affects the accuracy, and one can deduce that because 2D simulations contain fewer interparticle contacts than the equivalent 3D simulations, they therefore have lower accrued simulation errors. Copyright © 2016 John Wiley & Sons, Ltd.

A computational framework for scale-bridging in multi-scale simulations is presented. The framework enables seamless combination of at-scale models into highly dynamic hierarchies to build a multi-scale model. Its centerpiece is formulated as a standalone module capable of fully asynchronous operation. We assess its feasibility and performance for a two-scale model applied to two challenging test problems from impact physics. We find that the computational cost associated with using the framework may, as expected, become substantial. However, the framework has the ability of effortlessly combining at-scale models to render complex multi-scale models. The main source of the computational inefficiency of the framework is related to poor load balancing of the lower-scale model evaluation We demonstrate that the load balancing can be efficiently addressed by recourse to conventional load-balancing strategies. Copyright © 2016 John Wiley & Sons, Ltd.

An optimization framework is developed for surface texture design in hydrodynamic lubrication. The microscopic model of the lubrication interface is based on the Reynolds equation, and the macroscopic response is characterized through homogenization. The microscale setting assumes a unilateral periodic texture but implicitly accounts for the bilateral motion of the surfaces. The surface texture in a unit cell is described indirectly through the film thickness, which is allowed to vary between prescribed minimum and maximum values according to a morphology variable distribution that is obtained through the filtering of a design variable. The design and morphology variables are discretized using either element-wise constant values or through first-order elements. In addition to sharp textures, which are characterized by pillars and holes that induce sudden transitions between extreme film thickness values, the framework can also attain a variety of non-standard smoothly varying surface textures with a macroscopically isotropic or anisotropic response. Copyright © 2016 John Wiley & Sons, Ltd.

The macroscopic behavior of materials is affected by their inner micro-structure. Elementary considerations based on the arrangement, and the physical and mechanical features of the micro-structure may lead to the formulation of elastoplastic constitutive laws, involving hardening/softening mechanisms and non-associative properties. In order to model the non-linear behavior of micro-structured materials, the classical theory of time-independent multisurface plasticity is herein extended to Cosserat continua. The account for plastic relative strains and curvatures is made by means of a robust quadratic-convergent projection algorithm, specifically formulated for non-associative and hardening/softening plasticity. Some important limitations of the classical implementation of the algorithm for multisurface plasticity prevent its application for any plastic surfaces and loading conditions. These limitations are addressed in this paper, and a robust solution strategy based on the singular value decomposition technique is proposed. The projection algorithm is then implemented into a finite element formulation for Cosserat continua. A specific finite element is considered, developed for micropolar plates. The element is validated through illustrative examples and applications, showing able performance. Copyright © 2016 John Wiley & Sons, Ltd.

The main feature of the geometrically exact theory of contact is that all objects, which are necessary for computation, weak form and residual, linearized weak form, and tangent matrices are given in a covariant closed form in the local coordinate system corresponding to the geometry of contact pairs. This allows easily to construct computational algorithms for the normal and tangential follower forces as an inverse contact algorithm. In this case, following the definition of the follower forces as given and not changing in the local coordinate system, we have to modify all objects for the contact taking into account the definition of follower forces instead of constitutive relationships for the contact interfaces. The main feature is that the tangent matrices for both normal and tangential part being split into the rotational and the curvature parts are symmetric for any order of approximation. The following numerical examples are selected in the current article to illustrate the effectiveness of implementation: (1) the modeling of a pure bending with a moment applied either as a pair of single forces or a distributed follower forces (pressure) in both 2D and 3D cases; (2) modeling of inflation of a plate as application of the distributed follower normal forces (pressure); and (3) modeling of twisting of a beam with the rectangular cross-section as application of the tangential follower forces algorithm. Copyright © 2016 John Wiley & Sons, Ltd.

This article deals with the computation of guaranteed lower bounds of the error in the framework of finite element and domain decomposition methods. In addition to a fully parallel computation, the proposed lower bounds separate the algebraic error (due to the use of a domain decomposition iterative solver) from the discretization error (due to the finite element), which enables the steering of the iterative solver by the discretization error. These lower bounds are also used to improve the goal-oriented error estimation in a substructured context. Assessments on 2D static linear mechanic problems illustrate the relevance of the separation of sources of error and the lower bounds' independence from the substructuring. We also steer the iterative solver by an objective of precision on a quantity of interest. This strategy consists in a sequence of solvings and takes advantage of adaptive remeshing and recycling of search directions. Copyright © 2016 John Wiley & Sons, Ltd.

A new boundary element formulation is developed to analyze two-dimensional size-dependent piezoelectric response in isotropic dielectric materials. The model is based on the recently developed consistent couple stress theory, in which the couple-stress tensor is skew-symmetric. For isotropic materials, there is no classical piezoelectricity, and the size-dependent piezoelectricity or flexoelectricity effect is solely the result of coupling of polarization to the skew-symmetric mean curvature tensor. As a result, the size-dependent effect is specified by one characteristic length scale parameter *l*, and the electromechanical effect is specified by one flexoelectric coefficient *f*. Interestingly, in this size-dependent multi-physics model, the governing equations are decoupled. However, the problem is coupled, because of the existence of a flexoelectric effect in the boundary couple-traction and normal electric displacement. We discuss the boundary integral formulation and numerical implementation of this size-dependent piezoelectric boundary element method, which provides a boundary-only formulation involving displacements, rotations, force-tractions, couple-tractions, electric potential, and normal electric displacement as primary variables. Afterwards, we apply the resulting BEM formulation to several computational problems to confirm the validity of the numerical implementation and to explore the physics of the flexoelectric coupling. Copyright © 2016 John Wiley & Sons, Ltd.

When using particle methods to simulate water–air flows with compressible air pockets, a major challenge is to deal with the large differences in physical properties (e.g., density and viscosity) between water and air. In addition, the accurate modeling of air compressibility is essential. To this end, a new two-phase strategy is proposed to simulate incompressible and compressible fluids simultaneously without iterations between the solvers for incompressible and compressible flows. Water is modeled by the recently developed 2-phase Consistent Particle Method for incompressible flows. For air modeling, a new compressible solver is proposed based on the ideal gas law and thermodynamics. The formulation avoids the problem of determining the actual sound speed that is dependent on the temperature and is therefore not necessarily constant. In addition, the compressible air solver is seamlessly integrated with the incompressible solver 2-phase Consistent Particle Method because they both use the same predictor–corrector scheme to solve the governing equations. The performance of the proposed method is demonstrated by three benchmark problems as well as an experimental study of sloshing impact with entrapped air pockets in an oscillating tank. Copyright © 2016 John Wiley & Sons, Ltd.

Stress-related problems have not been given the same attention as the minimum compliance topological optimization problem in the literature. Continuum structural topological optimization with stress constraints is of wide engineering application prospect, in which there still are many problems to solve, such as the stress concentration, an equivalent approximate optimization model and etc. A new and effective topological optimization method of continuum structures with the stress constraints and the objective function being the structural volume has been presented in this paper. To solve the stress concentration issue, an approximate stress gradient evaluation for any element is introduced, and a total aggregation normalized stress gradient constraint is constructed for the optimized structure under the *r*−th load case. To obtain stable convergent series solutions and enhance the control on the stress level, two p-norm global stress constraint functions with different indexes are adopted, and some weighting p-norm global stress constraint functions are introduced for any load case. And an equivalent topological optimization model with reduced stress constraints is constructed,being incorporated with the rational approximation for material properties, an active constraint technique, a trust region scheme, and an effective local stress approach like the qp approach to resolve the stress singularity phenomenon. Hence, a set of stress quadratic explicit approximations are constructed, based on stress sensitivities and the method of moving asymptotes. A set of algorithm for the one level optimization problem with artificial variables and many possible non-active design variables is proposed by adopting an inequality constrained nonlinear programming method with simple trust regions, based on the primal-dual theory, in which the non-smooth expressions of the design variable solutions are reformulated as smoothing functions of the Lagrange multipliers by using a novel smoothing function. Finally, a two-level optimization design scheme with active constraint technique, i.e. varied constraint limits, is proposed to deal with the aggregation constraints that always are of loose constraint (non active constraint) features in the conventional structural optimization method. A novel structural topological optimization method with stress constraints and its algorithm are formed, and examples are provided to demonstrate that the proposed method is feasible and very effective. © 2016 The Authors. International Journal for Numerical Methods in Engineering published by John Wiley & Sons Ltd.

Fibrous soft biological tissues such as skin, ligaments, tendons, and arteries are non-homogeneous composite materials composed of fibers embedded in a ground substance. Cyclic tensile tests on these type of materials usually show a hysteretic stress–strain behavior in which strain rate dependence (viscoelasticity) and softening (Mullins' effect) play a coupled role. The main contribution of the present paper is to present unified variational approach to model both coupled phenomena: nonlinear viscoelasticity and Mullins-like softening behavior. The approach is labeled as variational because viscous-strain and damage internal variables are updated based on the minimization of a hyperelastic-like potential that takes a renewed value at each time step. Numerical examples explores (a) the versatility of the proposed model to account for the two described phenomena according to the chosen functions for the free-energy and dissipative potentials, (b) the ability of the time-integration scheme embedded in the incremental potential definition to allow for large time increments, and (c) the capability of the model to mimic experimentally obtained stress–strain cyclic curves of soft tissues. The model implementation on standard finite elements is also tested in which symmetric analytic tangent matrices are used as a natural consequence of the variational nature of the proposed approach. Copyright © 2016 John Wiley & Sons, Ltd.

This paper presents the study on non-deterministic problems of structures with a mixture of random field and interval material properties under uncertain-but-bounded forces. Probabilistic framework is extended to handle the mixed uncertainties from structural parameters and loads by incorporating interval algorithms into spectral stochastic finite element method. Random interval formulations are developed based on K–L expansion and polynomial chaos accommodating the random field Young's modulus, interval Poisson's ratios and bounded applied forces. Numerical characteristics including mean value and standard deviation of the interval random structural responses are consequently obtained as intervals rather than deterministic values. The randomised low-discrepancy sequences initialized particles and high-order nonlinear inertia weight with multi-dimensional parameters are employed to determine the change ranges of statistical moments of the random interval structural responses. The bounded probability density and cumulative distribution of the interval random response are then visualised. The feasibility, efficiency and usefulness of the proposed interval spectral stochastic finite element method are illustrated by three numerical examples. Copyright © 2016 John Wiley & Sons, Ltd.

A stochastic thermo-mechanical model for strength prediction of concrete is developed, based on the two-scale asymptotic expressions, which involves both the macroscale and the mesoscale of concrete materials. The mesoscale of concrete is characterized by a periodic layout of unit cells of matrix-aggregate composite materials, consisting of randomly distributed aggregate grains and cement matrix. The stochastic second-order and two-scale computational formulae are proposed in detail, and the maximum normal stress is assumed as the strength criterion for the aggregates, and the cement paste, in view of their brittle characteristics. Numerical results for the strength of concrete obtained from the proposed model are compared with those from known experiments. The comparison shows that the proposed method is validated for strength prediction of concrete. The proposed thermo-mechanical model is also employed to investigate the influence of different volume fraction of the aggregates on the strength of concrete. Copyright © 2016 John Wiley & Sons, Ltd.

In this work, the authors formulate a 2-D linearized ordinary state-based peridynamic model of elastic deformations and compute the stiffness matrix for 2-D plane stress/strain conditions. This model is then verified by testing the recovery of elastic properties for given Poisson's ratios in the range 0.1–0.45. The convergence behavior of peridynamic solutions in terms of the size of the nonlocal region by comparison with the classical (local) mechanics model is also discussed. The degree to which the peridynamic surface effect influences the recovery of elastic properties is examined, and stress/strain recovery values are found to have a definite influence on the results. The technique used here can provide the basis for applying 2-D peridynamic models to the study of fatigue failure and quasi-static fracture problems. Copyright © 2016 John Wiley & Sons, Ltd.

A new method for the real-time simulation of surgical cutting in haptic environments is presented. It is based on the intensive use of computational vademecums, that is, a sort of computational parametric meta-model, which is computed offline and only evaluated online. Therefore, the necessary time savings are obtained, allowing for feedback responses on the order of kilohertz. Such a high-dimensional, parametric solution of the problem is computed by employing proper generalized decomposition for the offline phase of the method, along with X-FEM techniques for the incorporation of the discontinuities in the displacement field after cutting, in the online phase. A thorough description of the proposed method, along with examples of its performance in the simulation of corneal surgery, are provided. Copyright © 2016 John Wiley & Sons, Ltd.

Thin-walled structures, when compressed, are prone to buckling. To fully utilize the capabilities of such structures, the post-buckling response should be considered and optimized in the design process. This work presents a novel method for gradient-based design optimization of the post-buckling performance of structures. The post-buckling analysis is based on Koiter's asymptotic method. To perform gradient-based optimization, the design sensitivities of the Koiter factors are derived, and new design optimization formulations based on the Koiter factors are presented. The proposed optimization formulations are demonstrated on a composite square plate and a curved panel where the post-buckling stability is optimized. Copyright © 2016 John Wiley & Sons, Ltd.

Predicting localized, nonlinear, thermoplastic behavior and residual stresses and deformations in structures subjected to intense heating is a prevalent challenge in a range of modern engineering applications. The authors present a generalized finite element method targeted at this class of problems, involving the solution of intrinsically parallelizable local boundary value problems to capture localized, time-dependent thermo-elasto-plastic behavior, which is embedded in the coarse, structural-scale approximation via enrichment functions. The method accommodates approximation spaces that evolve in between time or load steps while maintaining a fixed global mesh, which avoids the need to map solutions and state variables on changing meshes typical of traditional adaptive approaches. Representative three-dimensional examples exhibiting localized, transient, nonlinear thermal and thermomechanical effects are presented to demonstrate the advantages of the method with respect to available approaches, especially in terms of its flexibility and potential for realistic future applications in this area. Parallelism of the approach is also discussed. Copyright © 2016 John Wiley & Sons, Ltd.

A numerical method based on the movable cellular automata (MCA) is presented. The proposed method aims to investigate the tribological characteristics of polytetrafluoroethylene, polyimide, and polyetheretherketone sealing polymer materials in the microscale. The microscale frictional behavior and wear processes of four sealing materials are vividly shown through this grid method. The formation of a mechanically mixed layer is visualized, and the rotation angles of the automata and number of worn cellular automata are obtained. Friction coefficients of the four sealing materials are calculated and compared on a microscale mesh. Tests are conducted on a special test rig of the seal. Scanning electron microscopy investigations reveal the wear patterns of the four sealing composites. Results indicate that the MCA simulations are consistent with the scanning electron microscopy investigations in terms of the wear mechanism. The proposed MCA method can be used to investigate and compare the tribological characteristics of different sealing materials, which will significantly improve the efficiency of sealing material selection. Copyright © 2016 John Wiley & Sons, Ltd.

An isogeometric model is developed for the analysis of fluid transport in pre-existing faults or cracks that are embedded in a fluid-saturated deformable porous medium. Flow of the interstitial fluid in the porous medium and fluid transport in the discontinuities are accounted for and are coupled. The modelling of a fluid-saturated porous medium in general requires the interpolation of the displacements of the solid to be one order higher than that of the pressure of the interstitial fluid. Using order elevation and Bézier projection, a consistent procedure has been developed to accomplish this in an isogeometric framework. Particular attention has also been given to the spatial integration along the isogeometric interface element in order to suppress traction oscillations that can arise for certain integration rules when a relatively high dummy stiffness is used in a poromechanical model. © 2016 The Authors. International Journal for Numerical Methods in Engineering Published by John Wiley & Sons Ltd.

This paper presents the first method that enables the fully automatic generation of triangular meshes suitable for the so-called non-uniform rational B-spline (NURBS)-enhanced finite element method (NEFEM). The meshes generated with the proposed approach account for the computer-aided design boundary representation of the domain given by NURBS curves. The characteristic element size is completely independent of the geometric complexity and of the presence of very small geometric features. The proposed strategy allows to circumvent the time-consuming process of de-featuring complex geometric models before a finite element mesh suitable for the analysis can be produced. A generalisation of the original definition of a NEFEM element is also proposed, enabling to treat more complicated elements with an edge defined by several NURBS curves or more than one edge defined by different NURBS. Three examples of increasing difficulty demonstrate the applicability of the proposed approach and illustrate the advantages compared with those of traditional finite element mesh generators. Finally, a simulation of an electromagnetic scattering problem is considered to show the applicability of the generated meshes for finite element analysis. ©2016 The Authors. International Journal for NumericalMethods in Engineering published by John Wiley & Sons Ltd.

The finite element methods (FEMs) are important techniques in engineering for solving partial differential equations, but they depend heavily on element shape quality for stability and good performance. In this paper, we introduce the *Adaptive Extended Stencil Finite Element Method* (AES-FEM) as a means for overcoming this dependence on element shape quality. Our method replaces the traditional basis functions with a set of *generalized Lagrange polynomial basis functions*, which we construct using local weighted least-squares approximations. The method preserves the theoretical framework of FEM and allows imposing essential boundary conditions and integrating the stiffness matrix in the same way as the classical FEM. In addition, AES-FEM can use higher-degree polynomial basis functions than the classical FEM, while virtually preserving the sparsity pattern of the stiffness matrix. We describe the formulation and implementation of AES-FEM and analyze its consistency and stability. We present numerical experiments in both 2D and 3D for the Poisson equation and a time-independent convection–diffusion equation. The numerical results demonstrate that AES-FEM is more accurate than linear FEM, is also more efficient than linear FEM in terms of error versus runtime, and enables much better stability and faster convergence of iterative solvers than linear FEM over poor-quality meshes. Copyright © 2016 John Wiley & Sons, Ltd.

The mechanics of the interaction between a fluid and a soft interface undergoing large deformations appear in many places, such as in biological systems or industrial processes. We present an Eulerian approach that describes the mechanics of an interface and its interactions with a surrounding fluid via the so-called Navier boundary condition. The interface is modeled as a curvilinear surface with arbitrary mechanical properties across which discontinuities in pressure and tangential fluid velocity can be accounted for using a modified version of the extended finite element method. The coupling between the interface and the fluid is enforced through the use of Lagrange multipliers. The tracking and evolution of the interface are then handled in a Lagrangian step with the grid-based particle method. We show that this method is ideal to describe large membrane deformations and Navier boundary conditions on the interface with velocity/pressure discontinuities. The validity of the model is assessed by evaluating the numerical convergence for a axisymmetrical flow past a spherical capsule with various surface properties. We show the effect of slip length on the shear flow past a two-dimensional capsule and simulate the compression of an elastic membrane lying on a viscous fluid substrate. Copyright © 2015 John Wiley & Sons, Ltd.

This paper is concerned with numerical solution of the transient acoustic–structure interaction problems in three dimensions. An efficient and higher-order method is proposed with a combination of the exponential window technique and a fast and accurate boundary integral equation solver in the frequency-domain. The exponential window applied to the acoustic–structure system yields an artificial damping to the system, which eliminates the wrap-around errors brought by the discrete Fourier transform. The frequency-domain boundary integral equation approach relies on accurate evaluations of relevant singular integrals and fast computation of nonsingular integrals via the method of equivalent source representations and the fast Fourier transform. Numerical studies are presented to demonstrate the accuracy and efficiency of the method. Copyright © 2016 John Wiley & Sons, Ltd.

The basic principles of the discrete duality and nonlinear monotone finite volume methods are combined in order to obtain a new monotone nonlinear finite volume method for the approximation of diffusion operators on general meshes. Numerical results highlight both the second-order accuracy of this method on general meshes and its capability to deal with challenging anisotropic diffusion problems on various computational domains. Copyright © 2015 John Wiley & Sons, Ltd.

In this paper, based on the general stress–strain relationship, displacement and stress boundary-domain integral equations are established for single medium with varying material properties. From the established integral equations, single interface integral equations are derived for solving general multi-medium mechanics problems by making use of the variation feature of the material properties. The displacement and stress interface integral equations derived in this paper can be applied to solve non-homogeneous, anisotropic, and non-linear multi-medium problems in a unified way. By imposing some assumptions on the derived integral equations, detailed expressions for some specific mechanics problems are deduced, and a few numerical examples are given to demonstrate the correctness and robustness of the derived displacement and stress interface integral equations. Copyright © 2015 John Wiley & Sons, Ltd.

A fully discrete second-order decoupled implicit/explicit method is proposed for solving 3D primitive equations of ocean in the case of Dirichlet boundary conditions on the side, where a second-order decoupled implicit/explicit scheme is used for time discretization, and a finite element method based on the *P*_{1}(*P*_{1}) − *P*_{1}−*P*_{1}(*P*_{1}) elements for velocity, pressure and density is used for spatial discretization of these primitive equations. Optimal *H*^{1}−*L*^{2}−*H*^{1} error estimates for numerical solution
and an optimal *L*^{2} error estimate for
are established under the convergence condition of 0 < *h*≤*β*_{1},0 < *τ*≤*β*_{2}, and *τ*≤*β*_{3}*h* for some positive constants *β*_{1},*β*_{2}, and *β*_{3}. Furthermore, numerical computations show that the *H*^{1}−*L*^{2}−*H*^{1} convergence rate for numerical solution
is of *O*(*h* + *τ*^{2}) and an *L*^{2} convergence rate for
is *O*(*h*^{2}+*τ*^{2}) with the assumed convergence condition, where *h* is a mesh size and *τ* is a time step size. More practical calculations are performed as a further validation of the numerical method. Copyright © 2016 John Wiley & Sons, Ltd.

A number of methods have been developed for solving the dynamics of saturated porous media. However, most solutions are based on the finite element method, and only a few employ finite differences (FDM). One problem with the FDM is the difficulty in fulfilling the inf-sup (Ladyženskaja-Babuška-Brezzi) condition. This paper explores solutions with the FDM, including the development of new schemes aiming at stabilised formulations. The efficiency, accuracy and stability of several FDM and finite element method algorithms are thoroughly investigated as well. A combination of primary variables from the theory of porous media is considered, including the so-called *up* and *uvp* formulations. Six numerical schemes are produced and quantitatively studied. Simulations of 1D and 2D wave propagation problems are performed in order to reveal the advantages and drawbacks of all schemes. Copyright © 2016 John Wiley & Sons, Ltd.

This paper proposes a first step towards a framework to develop shell elements applicable to any deformation regime. Here, we apply it to the large and moderate deformations of, respectively, plates and shells, showing with standard benchmarks that the resulting low-order discretization is competitive against the best elements for either membrane-dominated or bending-dominated scenarios. Additionally, we propose a new test for measuring membrane locking, which highlights the mesh-independence properties of our element. Our strategy is based on building a discrete model that mimics the smooth behavior by construction, rather than discretizing a smooth energy. The proposed framework consists of two steps: (i) defining a discrete kinematics by means of constraints and (ii) formulating an energy that vanishes on such a constraint manifold. We achieve (i) by considering each triangle as a tensegrity structure, constructed to be unstretchable but bendable isometrically (in a discrete sense). We then present a choice for (ii) based on assuming a linear strain field on each triangle, using tools from differential geometry for coupling the discrete membrane energy with our locking-free kinematic description. We argue that such a locking-free element is only a member of a new family that can be created using our framework (i) and (ii). Copyright © 2016 John Wiley & Sons, Ltd.

It is known that the explicit time integration is conditionally stable. The very small time step leads to increase of computational time dramatically. In this paper, a mass-redistributed method is formulated in different numerical schemes to simulate transient quasi-harmonic problems. The essential idea of the mass-redistributed method is to shift the integration points away from the Gauss locations in the computation of mass matrix for achieving a much larger stable time increment in the explicit method. For the first time, it is found that the stability of explicit method in transient quasi-harmonic problems is proportional to the softened effect of discretized model with mass-redistributed method. With adjustment of integration points in the mass matrix, the stability of transient models is improved significantly. Numerical experiments including 1D, 2D and 3D problems with regular and irregular mesh have demonstrated the superior performance of the proposed mass-redistributed method with the combination of smoothed finite element method in terms of accuracy as well as stability. Copyright © 2016 John Wiley & Sons, Ltd.

The equivalent static load (ESL) method is a powerful approach to solve dynamic response structural optimization problems. The method transforms the dynamic response optimization into a static response optimization under multiple load cases. The ESL cases are defined based on the transient analysis response whereupon all the standard techniques of static response optimization can be used. In the last decade, the ESL method has been applied to perform the structural optimization of flexible components of mechanical systems modeled as multibody systems (MBS). The ESL evaluation strongly depends on the adopted formulation to describe the MBS and has been initially derived based on a floating frame of reference formulation. In this paper, we propose a method to derive the ESL adapted to a nonlinear finite element approach based on a Lie group formalism for two main reasons. Firstly, the finite element approach is completely general to analyze complex MBS and is suitable to perform more advanced optimization problems like topology optimization. Secondly, the selected Lie group formalism leads to a formulation of the equations of motion in the local frame, which turns out to be a strong practical advantage for the ESL evaluation. Examples are provided to validate the proposed method. Copyright © 2016 John Wiley & Sons, Ltd.

This paper proposes a generalized finite element method based on the use of parametric solutions as enrichment functions. These parametric solutions are precomputed off-line and stored in memory in the form of a computational vademecum so that they can be used on-line with negligible cost. This renders a more efficient computational method than traditional finite element methods at performing simulations of processes. One key issue of the proposed method is the efficient computation of the parametric enrichments. These are computed and efficiently stored in memory by employing proper generalized decompositions. Although the presented method can be broadly applied, it is particularly well suited in manufacturing processes involving localized physics that depend on many parameters, such as welding. After introducing the vademecum-generalized finite element method formulation, we present some numerical examples related to the simulation of thermal models encountered in welding processes. Copyright © 2016 John Wiley & Sons, Ltd.

The variational theory of complex rays (VTCR) is an indirect Trefftz method designed to study systems governed by Helmholtz-like equations. It uses wave functions to represent the solution inside elements, which reduces the dispersion error compared with classical polynomial approaches, but the resulting system is prone to be ill-conditioned. This paper gives a simple and original presentation of the VTCR using the discontinuous Galerkin framework, and it traces back the ill-conditioning to the accumulation of eigenvalues near zero for the formulation written in terms of wave amplitude. The core of this paper presents an efficient solving strategy that overcomes this issue. The key element is the construction of a search subspace where the condition number is controlled at the cost of a limited decrease of attainable precision. An augmented LSQR solver is then proposed to solve efficiently and accurately the complete system. The approach is successfully applied to different examples. Copyright © 2015 John Wiley & Sons, Ltd.

Incompressible free-surface flow is a common assumption for the modelling of water waves. Connected with the aim to develop very large floating platforms, air chamber supported floating structures have attracted considerable research interest in the past. Such structures are carried by air entrapped in chambers formed by stiff, vertical walls. In order to model these types of structures, the interactions between surface gravity waves and compressible air must be taken into account. If the payload requirements for air chamber supported structures are low enough, the air chambers may be formed by flexible membrane cylinders. In such systems, pressure variations can lead to considerable changes in chamber volume. Therefore, the flexibility of the bounding structures must be taken into account. We present a modelling strategy to tackle the fully coupled problem of compressible gas in a flexible chamber and incompressible free-surface flow in an unbounded domain. The governing equations and boundary conditions are described and solved by the finite element method. A perfectly matched layer is used to obtain a solution for an unbounded domain. Finally, the numerical implementation is validated by various test cases. Copyright © 2015 John Wiley & Sons, Ltd.

The problem of robust optimal Robin boundary control for a parabolic partial differential equation with uncertain input data is considered. As a measure of robustness, the variance of the random system response is included in two different cost functionals. Uncertainties in both the underlying state equation and the control variable are quantified through random fields. The paper is mainly concerned with the numerical resolution of the problem. To this end, a gradient-based method is proposed considering different functional costs to achieve the robustness of the system. An adaptive anisotropic sparse grid stochastic collocation method is used for the numerical resolution of the associated state and adjoint state equations. The different functional costs are analysed in terms of computational efficiency and its capability to provide robust solutions. Two numerical experiments illustrate the performance of the algorithm. Copyright © 2016 John Wiley & Sons, Ltd.

We propose a method for efficient evaluation of surface integrals arising in boundary element methods for three-dimensional Helmholtz problems (with real positive wavenumber *k*), modelling wave scattering and/or radiation in homogeneous media. To reduce the number of degrees of freedom required when *k* is large, a common approach is to include in the approximation space oscillatory basis functions, with support extending across many wavelengths. A difficulty with this approach is that it leads to highly oscillatory surface integrals whose evaluation by standard quadrature would require at least O(*k*^{2}) quadrature points. Here, we use equivalent contour integrals developed for aperture scattering in optics to reduce this requirement to O(*k*), and possible extensions to reduce it further to O(1) are identified. The contour integral is derived for arbitrary shaped elements, but its application is limited to planar elements in many cases. In addition, the transform regularises the singularity in the surface integrand caused by the Green's function, including for the hyper-singular case under appropriate conditions. An open-source Matlab™code library is available to demonstrate our routines. © 2016 The Authors International Journal for Numerical Methods in Engineering Published by John Wiley & Sons Ltd.

No abstract is available for this article.

]]>No abstract is available for this article.

]]>In transient finite element analysis, reducing the time-step size improves the accuracy of the solution. However, a lower bound to the time-step size exists, below which the solution may exhibit spatial oscillations at the initial stages of the analysis. This numerical ‘shock’ problem may lead to accumulated errors in coupled analyses. To satisfy the non-oscillatory criterion, a novel analytical approach is presented in this paper to obtain the time-step constraints using the *θ*-method for the transient coupled analysis, including both heat conduction–convection and coupled consolidation analyses. The expressions of the minimum time-step size for heat conduction–convection problems with both linear and quadratic elements reduce to those applicable to heat conduction problems if the effect of heat convection is not taken into account. For coupled consolidation analysis, time-step constraints are obtained for three different types of elements, and the one for composite elements matches that in the literature. Finally, recommendations on how to handle the numerical ‘shock’ issues are suggested. Copyright © 2015 John Wiley & Sons, Ltd.

The present work addresses shape sensitivity analysis and optimization in two-dimensional elasticity with a regularized isogeometric boundary element method (IGABEM). Non-uniform rational B-splines are used both for the geometry and the basis functions to discretize the regularized boundary integral equations. With the advantage of tight integration of design and analysis, the application of IGABEM in shape optimization reduces the mesh generation/regeneration burden greatly. The work is distinct from the previous literatures in IGABEM shape optimization mainly in two aspects: (1) the structural and sensitivity analysis takes advantage of the regularized form of the boundary integral equations, eliminating completely the need of evaluating strongly singular integrals and jump terms and their shape derivatives, which were the main implementation difficulty in IGABEM, and (2) although based on the same Computer Aided Design (CAD) model, the mesh for structural and shape sensitivity analysis is separated from the geometrical design mesh, thus achieving a balance between less design variables for efficiency and refined mesh for accuracy. This technique was initially used in isogeometric finite element method and was incorporated into the present IGABEM implementation. Copyright © 2016 John Wiley & Sons, Ltd.

In this paper, the extended finite element method (X-FEM) formulation for the modeling of arbitrary crack propagation in coupled shell/solid structures is developed based on the large deformation continuum-based (CB) shell theory. The main features of the new method are as follows: (1) different kinematic equations are derived for different fibers in CB shell elements, including the fibers enriched by shifted jump function or crack tip functions and the fibers cut into two segments by the crack surface or connecting with solid elements. So the crack tip can locate inside the element, and the crack surface is not necessarily perpendicular to the middle surface. (2) The enhanced CB shell element is developed to realize the seamless transition of crack propagation between shell and solid structures. (3) A revised interaction integral is used to calculate the stress intensity factor (SIF) for shells, which avoids that the auxiliary fields for cracks in Mindlin–Reissner plates cannot satisfy exactly the equilibrium equations. Several numerical examples, including the calculation of SIF for the cracked plate under uniform bending and crack propagation between solid and shell structures are presented to demonstrate the performance of the developed method. Copyright © 2015 John Wiley & Sons, Ltd.