ABSTRACT Large‐scale three‐dimensional (3D) finite‐discrete element method (FDEM) simulations impose high demands on contact models. Existing contact formulations typically evaluate these contacts by explicitly reconstructing overlap volumes or contact surface patches, which requires expensive geometric clipping and can become the dominant computational bottleneck. In this work, we develop a symmetric boundary‐representation contact model for 3D tetrahedral FDEM (SBC‐FDEM). Based on Stokes' theorem, this approach recasts the contact vector area as a boundary integral evaluated exclusively from oriented boundary segments. For each tetrahedron pair, these segments are identified via a fixed 4 4 sweep of lightweight triangle‐triangle intersection tests. By bypassing computationally expensive overlap volume or surface patch reconstructions, this boundary‐only formulation is considerably simpler and maps naturally onto modern GPU architectures. The same boundary representation also yields a boundary depth factor that removes the flat‐facet force jump, so the contact force grows continuously from zero. Moreover, strict swap invariance is achieved by evaluating both the boundary segments and the depth factor symmetrically. The method is verified using pair‐level swap tests, elastic rebound tests, closed‐loop energy tests, and laboratory‐scale UCS and BTS simulations. SBC‐FDEM is swap‐invariant to machine precision and keeps the shallow‐contact spurious energy small. Large‐scale impact fragmentation and slope collapse simulations further demonstrate its robustness in handling rapidly evolving contact interactions. On NVIDIA V100 and A100 GPUs, SBC‐FDEM accelerates the contact force computation by 1.9–2.9× compared with the state‐of‐the‐art contact algorithm (ENIC‐FDEM) and by 7.6–17.0× compared with the widely used classical potential‐based contact model (CPC‐FDEM). The proposed method can be generalized to FDEM formulations based on polygonal and polyhedral discretizations.
ABSTRACT In geometrically nonlinear topology optimization, elements with low pseudo‐densities, subjected to compression from neighboring stiffer regions, can lead to complications, such as “inside‐out” elements in 2D structures or spurious local buckling in shells and plates. The former typically causes computational divergence, while the latter, though not always divergent, could significantly increase the analysis iterations required for convergence, resulting in significantly higher computational costs. To enhance efficiency, reduced‐order modeling (ROM) could be used, as it employs a small number of generalized variables to approximate the behavior governed by the numerous degrees of freedom (DOFs) in the full‐order model (FOM). However, the presence of spurious buckling can severely compromise ROM efficiency, even requiring more analysis iterations compared to methods relying on FOMs. To address the analysis convergence difficulties caused by spurious buckling in ROMs, two types of approaches are investigated. The first and simplest approach involves removing spurious modes from the ROM basis. The second approach focuses on eliminating spurious modes directly from the underlying physics, which can be extended to standard FOM‐based analyses. The above strategies are primarily evaluated using shell models, providing a detailed analysis of their effectiveness. The objective of this paper is not to propose a complete or definitive solution to spurious modes. Rather, the goal is to highlight and compare the strengths and limitations of different methodologies, and to provide insight into their behavior when applied to the problems considered.
ABSTRACT This paper presents a fractional thermo‐viscoelastic formulation for studying the nonlinear vibrational behavior caused by thermal loads in a functionally graded (FGM) Timoshenko‐type beam exposed to a sudden thermal shock. It is presumed that the beam consists of a ceramic‐metal mixture whose effective properties vary continuously across the thickness following a power‐law gradation. The transient heat transfer is described using the Lord–Shulman generalized thermoelasticity theory extended to incorporate a fractional‐order time derivative of order α . Subsequently, the mechanical behavior follows the fractional‐order Kelvin‐Voigt viscoelastic constitutive law, where the stress–strain relation contains a fractional derivative of order γ . Hamilton's variational principle is employed to derive the weak form of the governing equations, with geometric nonlinearity incorporated through the von Kármán strain relations. Spatial discretization is performed by a combination of the variational differential quadrature (VDQ) and generalized differential quadrature (GDQ) methods. The fractional time derivatives are treated by the Grünwald‐Letnikov (GL) scheme, and the resulting nonlinear system is integrated in time using the Newmark‐β algorithm together with Newton–Raphson iterations. A detailed parametric study investigates the effects of the fractional‐order parameters and , the viscoelastic coefficient, the relaxation time, the power‐law index, and different thermal loading conditions on the transient displacement and temperature fields. The numerical results demonstrate that lower values of increase damping and reduce oscillation amplitudes. Regarding the parameter , increasing it leads to a higher peak temperature and larger oscillation amplitudes. The presented computational approach offers a reliable and accurate tool for designing FGM structures under extreme thermal environments.
ABSTRACT A new version of the multiscale finite element method (MsFEM) for elastic solid bodies, based on the mixed ultra‐weak formulation is presented. The mixed formulation ensures fast stress convergence, while its ultra‐weak variant facilitates the construction of conforming shape functions. Furthermore, the stability guaranteed by the discontinuous Petrov–Galerkin (DPG) methodology, combined with the high approximability inherent to MsFEM, ensures the convergence of the numerical solution. The key innovation of the proposed MsFEM modification lies in the computation of the prolongation operator between fine and coarse finite element meshes. The proposed approach, which accounts for high‐order approximation, has been successfully verified on several representative numerical examples for 3D elasticity problems and can also be applied to direct FEM computations (without upscaling) for heterogeneous materials, when the multigrid method is used to solve the systems of algebraic equations for the ultra‐weak mixed formulation.
ABSTRACT The evolution of strain localization and the corresponding macroscopic response obtained in multiscale modeling of softening heterogeneous materials strongly depend on the boundary conditions imposed on the fine‐scale model. For rectilinear models (e.g., squares or cubes), employing standard Periodic BCs leads to overestimation of energy dissipation and artificially ductile behavior, whenever the localization band orientation deviates from the periodic directions. Recently proposed Tessellation and Percolation‐path‐aligned BCs promise to address this by aligning the periodicity frame with evolving localization bands. Alternatively, spherical/circular models provide an orientation independent response by design. Unfortunately, the standard Periodic BCs do not allow development of proper localization band crossing the spherical model's boundaries. A recently proposed modification addresses this by adding a displacement jump to the spherical periodic BCs. This study evaluates the applicability of these novel BCs to a mesoscale discrete particle model of concrete. Two‐dimensional square and circular models under uniaxial tension with different loading directions are analyzed, with selected approaches extended to three‐dimensional cube models. Percolation‐path‐aligned BCs were found to suffer from several drawbacks. Uneven straining of the weakly and strongly constrained sections can trigger localization into multiple bands, and the weakly constrained section is susceptible to spurious strain localization. Tessellation BCs, on the other hand, consistently allow the development of a single localization band of predictable length, determined only by the model geometry, enabling a straightforward correction of its influence during post‐processing. As a result, variations in dissipated energy with loading direction can be attributed solely to model geometry rather than to artifacts introduced by the boundary conditions. Periodic boundary conditions augmented with a displacement jump applied to a circular model yield inconsistent results. Some cases exhibit correct localization band development, but some produce crack patterns similar to those under the standard Periodic BCs. Which of these occurs is dictated by the internal heterogeneity of the model, which effectively determines whether strain localization is initiated in the boundary segment enhanced by the displacement jump.
