Parallelization of the distinct lattice spring model

The distinct lattice spring model (DLSM) is a newly developed numerical tool for modeling rock dynamics problems, i.e. dynamic failure and wave propagation. In this paper, parallelization of DLSM is presented. With the development of parallel computing technologies in both hardware and software, parallelization of a code is becoming easier than before. There are many available choices now. In this paper, Open Multi‐Processing (OpenMP) with multicore personal computer (PC) and message passing interface (MPI) with cluster are selected as the environments to parallelize DLSM. Performances of these parallel DLSM codes are tested on different computers. It is found that the parallel DLSM code with OpenMP can reach a maximum speed‐up of 4.68× on a quad‐core PC. The parallel DLSM code with MPI can achieve a speed‐up of 40.886× when 256 CPUs are used on a cluster. At the end of this paper, a high‐resolution model with four million particles, which is too big to handle by the serial code, is simulated by using the parallel DLSM code on a cluster. It is concluded that the parallelization of DLSM is successful. Copyright © 2011 John Wiley & Sons, Ltd.     

Validation study of the distinct lattice spring model (DLSM) on P-wave propagation across multiple parallel joints

A validation study of the distinct lattice spring model (DLSM) for wave propagation problems is performed. DLSM is a microstructure-based numerical model, which is meshless and has advantages in modelling dynamic problems where stress wave propagation is important. To verify the applicability of DLSM to modelling wave propagation through a discontinuous medium, the virtual wave source (VWS) method is used to obtain analytical solutions for wave propagation across a jointed rock mass. Numerical modelling results of the commercial code UDEC are selected as the reference. The effects of particle size and lattice rotation angle on wave propagation are first studied. Then, the results of wave transmission across a single joint with a different joint stiffness and across multiple parallel joints with different joint spacings are derived with DLSM, UDEC and VWS. These results are in good agreement with each other. Therefore, the capability of DLSM to model P-wave propagation across jointed rock mass is verified, which provides confidence for the further application of DLSM to modelling more complex problems.     

Dynamic Fracturing Simulation of Brittle Material using the Distinct Lattice Spring Method with a Full Rate-Dependent Cohesive Law

A full rate-dependent cohesive law is implemented in the distinct lattice spring method (DLSM) to investigate the dynamic fracturing behavior of brittle materials. Both the spring ultimate deformation and spring strength are dependent on the spring deformation rate. From the simulation results, it is found that the dynamic crack propagation velocity can be well predicted by the DLSM through the implemented full rate-dependent cohesive law. Furthermore, a numerical investigation on dynamic branching is also conducted by using the DLSM code.     

Development of the distinct lattice spring model for large deformation analyses

This study develops the distinct lattice spring model (DLSM) for geometrically nonlinear large deformation problems. The formulation of a spring bond deformation under a large deformation is derived under the Lagrange framework using polar decomposition. The results reveal that the DLSM's stiffness matrix under small deformations is the tangent stiffness matrix of the DLSM under large deformations. The formulation of the spring bond internal force under a given configuration is also presented and can be used to calculate the unbalanced force. Using these formulations, three nonlinear solving methods (the Euler method, modified Euler method, and Newton method) are developed for the DLSM with which to tackle large deformation problems. To investigate the performance of the developed model, three numerical examples involving large deformations are presented, the results of which are also in good agreement with the analytical and finite element method solutions. Copyright © 2013 John Wiley & Sons, Ltd.     

A 3D distinct lattice spring model for elasticity and dynamic failure

A 3D distinct lattice spring model (DLSM) is proposed where matter is discretized into individual particles linked by springs. The presented model is different from the conventional lattice spring models where a shear spring is introduced to model the multibody force by evaluating the spring deformation from the local strain rather than the particle displacement. By doing this, the proposed model can represent the diversity of Poisson's ratio without violating the rotational invariance. The local strain of the spring is calculated through a least square method which makes the model possessing meshless properties. Because of this and explicitly representing the microstructure, DLSM is able to model dynamic fracturing problems and can be used to study the microstructure influences. The material parameters inputted in the model is the conventional material parameters, e.g. the elastic modules and the Poisson's ratio. Relationships between microscopic spring parameters and macroscopic material constants are derived based on the Cauchy–Born rules and the hyperelastic theory. Numerical examples are presented to show the abilities and properties of DLSM in modeling elastic and dynamic failure problems. Copyright © 2010 John Wiley & Sons, Ltd.     

Implementation of a coupled plastic damage distinct lattice spring model for dynamic crack propagation in geomaterials

This paper studies dynamic crack propagation by employing the distinct lattice spring model (DLSM) and 3‐dimensional (3D) printing technique. A damage‐plasticity model was developed and implemented in a 2D DLSM. Applicability of the damage‐plasticity DLSM was verified against analytical elastic solutions and experimental results for crack propagation. As a physical analogy, dynamic fracturing tests were conducted on 3D printed specimens using the split Hopkinson pressure bar. The dynamic stress intensity factors were recorded, and crack paths were captured by a high‐speed camera. A parametric study was conducted to find the influences of the parameters on cracking behaviors, including initial and peak fracture toughness, crack speed, and crack patterns. Finally, selection of parameters for the damage‐plasticity model was determined through the comparison of numerical predictions and the experimentally observed cracking features.     

应用图形处理器快速计算逆时偏移

逆时偏移是目前精度最高的地震数据叠前深度偏移方法,但高强度的计算需求限制了其在工业生产领域的大规模应用。可编程图形处理器的发展为逆时偏移的快速计算提供了一种新的计算选择。围绕如何在图形处理器上开展逆时偏移计算展开,总结了图形处理器计算的优化关键,并根据逆时偏移的特点着重介绍了两个优化环节：一个是应用随机边界条件,以计算换存储,减少数据在主机和图形处理器间的传输;二是应用共享存储器来存储正演计算的波场,相比全局存储器,提高了数据读取的带宽。应用Marmousi模型数据对经过上述优化后的程序进行了测试,结果表明,图形处理器逆时偏移程序得到了很好的优化,提高了计算效率。     

A SPH framework for dynamic interaction between soil and rigid body system with hybrid contact method

A smoothed particle hydrodynamics (SPH) framework for three-dimensional dynamic soil-multibody interaction modeling is presented, where both soils and rigid bodies are discretized using SPH particles. In the framework, soils are modeled using the Drucker-Prager model, while rigid bodies are considered with a multibody dynamics solver. A hybrid contact method suitable for three-dimensional simulations is developed to model the soil-body and body-body frictionless and frictional contacts, where contact forces are calculated based on ideal plastic collision and the unit normal/tangential vectors of the actual surface. Owing to its simplicity in contact detection and accuracy in contact force calculation, the hybrid contact method can be easily incorporated into SPH. Furthermore, graphics processing unit (GPU) parallelization is utilized to improve efficiency. The presented numerical framework and the hybrid contact method are validated using several examples. Numerical results are compared with analytical solutions and results from the literature. Furthermore, two three-dimensional simulations involving dynamic soil-multibody interaction are included to demonstrate the application.     

Parallel finite element analysis of seismic soil structure interaction using a PC cluster

This paper describes an implementation of a highly scalable parallel computational facility with high speedup efficiency using relatively low-cost hardware, which consists of a cluster of desktop personal computers (PCs) connected via a 10-Gigabit Ethernet. Two-levels of parallelization were implemented. Communication between different PCs was achieved using message passing interface (MPI) protocol. Domain decomposition was automated and based on element numbering. Domain continuity was assured largely by re-numbering the elements using a “front squasher” code prior to decomposition. Within each PC, the shared memory parallelization was implemented using either the open-multiprocessing (OpenMP) or the MPI protocol. Analysis of three different problems with number of degrees-of-freedom ranging from about 129,000 to about 2,260,000 shows a speedup efficiency generally above 70%. Super-linear speedup was achieved in several of the cases examined in this study, with the hybrid MPI-OpenMP approach generally performing better compared to the pure MPI method for parallelization. The results demonstrate the feasibility of acquiring a parallel computing facility with relatively modest outlay that is within the reach of consulting or engineering offices.     

On the linear elastic responses of the 2D bonded discrete element model

The bonded discrete element model (DEM) is a numerical tool that is becoming widely used when studying fracturing, fragmentation, and failure of solids in various disciplines. However, its abilities to solve elastic problems are usually overlooked. In this work, the main features of the 2D bonded DEM which influence Poisson's ratio and Young's modulus, and accuracy when solving elastic boundary value problems, are investigated. Outputs of numerical simulations using the 2D bonded DEM, the finite element method, a hyper elasticity analysis, and the distinct lattice spring model (DLSM) are compared in the investigation. It is shown that a shear interaction (local) factor and a geometric (global) factor are two essential elements for the 2D bonded DEM to reproduce a full range of Poisson's ratios. It is also found that the 2D bonded DEM might be unable to reproduce the correct displacements for elastic boundary value problems when the represented Poisson's ratio is close to 0.5 or the long-range interaction is considered. In addition, an analytical relationship between the shear stiffness ratio and the Poisson's ratio, derived from a hyper elasticity analysis and applicable to discontinuum-based models, provides good agreement with outputs from the 2D bonded DEM and DLSM. Finally, it is shown that the selection of elastic parameters used the 2D bonded DEM has a significant effect on fracturing and fragment patterns of solids.     

Parallel processing for a discrete element program

Reconfiguration of a discrete element code for parallel operation provided the opportunity to compare processing speeds on various hardware platforms. The program employed in this comparison, NURBM3DP, is a three dimensional, distinct element code employed to calculate dynamic response of a cavern in a jointed rock mass. On a 16 processor IBM SP2, it is capable of calculating dynamic response with 1000's of explicit time steps of jointed rock masses with up to 2,000,000 blocks. Comparison of single instruction multiple data stream (SIMD) and multiple-instruction multiple-data stream (MIMD) operation showed MIMD processing to provide the best overall parallelization. The full report of the comparisons of operation on different hardware with different data streaming configurations can be found at the research section of the Northwestern University Computational Mechanics site: http://www.tam.nwu.edu/compmech.html. In addition, a color movie of dynamic response of a million block model of a cavern responding to dynamic excitation can be seen at: http://geotech.civen.okstate.edu/ejge/ppr9801/index.htm.     

Development of a GPGPU-parallelized hybrid finite-discrete element method for modeling rock fracture

The hybrid finite-discrete element method (FDEM) is widely used for engineering applications, which, however, is computationally expensive and needs further development, especially when rock fracture process is modeled. This study aims to further develop a sequential hybrid FDEM code formerly proposed by the authors and parallelize it using compute unified device architecture (CUDA) C/C++ on the basis of a general-purpose graphics processing unit (GPGPU) for rock engineering applications. Because the contact detection algorithm in the sequential code is not suitable for GPGPU parallelization, a different contact detection algorithm is implemented in the GPGPU-parallelized hybrid FDEM. Moreover, a number of new features are implemented in the hybrid FDEM code, including the local damping technique for efficient geostatic stress analysis, contact damping, contact friction, and the absorbing boundary. Then, a number of simulations with both quasi-static and dynamic loading conditions are conducted using the GPGPU-parallelized hybrid FDEM, and the obtained results are compared both quantitatively and qualitatively with those from either theoretical analysis or the literature to calibrate the implementations. Finally, the speed-up performance of the hybrid FDEM is discussed in terms of its performance on various GPGPU accelerators and a comparison with the sequential code, which reveals that the GPGPU-parallelized hybrid FDEM can run more than 128 times faster than the sequential code if it is run on appropriate GPGPU accelerators, such as the Quadro GP100. It is concluded that the GPGPU-parallelized hybrid FDEM developed in this study is a valuable and powerful numerical tool for rock engineering applications.     

Study on the particle breakage of ballast based on a GPU accelerated discrete element method

Breakage of particles will have greatly influence on mechanical behavior of granular material(GM)under external loads,such as ballast,rockfill and sand.The discrete element method(DEM)is one of the most popular methods for simulating GM as each particle is represented on its own.To study breakage mechanism of particle breakage,a cohesive contact mode is developed based on the GPU accelerated DEM code-Blaze-DEM.A database of the 3D geometry model of rock blocks is established based on the 3D scanning method.And an agglomerate describing the rock block with a series of non-overlapping spherical particles is used to build the DEM numerical model of a railway ballast sample,which is used to the DEM oedometric test to study the particles’breakage characteristics of the sample under external load.Furthermore,to obtain the meso-mechanical parameters used in DEM,a black-analysis method is used based on the laboratory tests of the rock sample.Based on the DEM numerical tests,the particle breakage process and mechanisms of the railway ballast are studied.All results show that the developed code can better used for large scale simulation of the particle breakage analysis of granular material.     

GPU通用计算模式在岩土工程中的应用

由于岩土工程地质条件的复杂性及其规模的不断增大,对大规模数值计算速度的要求越来越高。显卡核心单元（GPU）由于其硬件构造特殊,有着并行计算上的独特优势、高速浮点运算性能和超高的内存带宽,可以很好地解决大规模的科学计算速度问题。文中介绍了GPU与CPU的硬件构架差异,总结了多核CPU、工作站等方式发展的局限性及GPU在并行运算方面的优势,详细阐述了GPU各类计算模式的发展特点及其成果,展示了其在坝区渗透特性中随机微分方程加速求解过程中的优越性,探讨了采用GPU进行大规模岩土工程数值计算的应用前景。     

OpenMP在CO2地质储存数值模拟并行计算中的应用

为减少CO2地质储存数值模拟计算时间和增强计算规模,文章提出基于OpenMP和动态内存分配的方式,重构TOUGH2-ECO2N数值模拟器。通过耗时评估可知,模拟器的主要耗时部分为状态方程计算、组建矩阵方程和方程求解。基于此,在遵循OpenMP并行化原则下,采用动态内存分配、处理跳转语句和算法内部的相关性,以及函数内部参数优化等措施,完成了多相流模拟器的并行化。计算试验结果表明,并行化模拟器算法正确、执行效率高,且具有加速效果良好的特点。针对中小规模的模型对比试验,4核的加速比可以达到2.28倍。     

GPU提速叠前时间体偏移技术

为了进一步提高叠前时间体偏移的计算效率,实现了在GPU＼CPU协同并行计算模式下Kirchhoff叠前时间体偏移技术,并进行优化。经在Nvida Tesla C1060GPU上的测试表明,GPU（Graphic Processing Unit）的处理速度是CPU（单核）的四十倍左右。同时表明,CUDA（Cornpute Unified Device Architectarc）编程为CPU向GPU的转化提供了一个较为方便的语言环境。     

Reduction of computing time for least-squares migration based on the Helmholtz equation by graphics processing units

In geophysical applications, the interest in least-squares migration (LSM) as an imaging algorithm is increasing due to the demand for more accurate solutions and the development of high-performance computing. The computational engine of LSM in this work is the numerical solution of the 3D Helmholtz equation in the frequency domain. The Helmholtz solver is Bi-CGSTAB preconditioned with the shifted Laplace matrix-dependent multigrid method. In this paper, an efficient LSM algorithm is presented using several enhancements. First of all, a frequency decimation approach is introduced that makes use of redundant information present in the data. It leads to a speedup of LSM, whereas the impact on accuracy is kept minimal. Secondly, a new matrix storage format Very Compressed Row Storage (VCRS) is presented. It not only reduces the size of the stored matrix by a certain factor but also increases the efficiency of the matrix-vector computations. The effects of lossless and lossy compression with a proper choice of the compression parameters are positive. Thirdly, we accelerate the LSM engine by graphics cards (GPUs). A GPU is used as an accelerator, where the data is partially transferred to a GPU to execute a set of operations or as a replacement, where the complete data is stored in the GPU memory. We demonstrate that using the GPU as a replacement leads to higher speedups and allows us to solve larger problem sizes. Summarizing the effects of each improvement, the resulting speedup can be at least an order of magnitude compared to the original LSM method.     

基于GPU的任意三维复杂形体重磁异常快速计算

提出了基于图形处理单元的任意三维复杂形体的重磁异常快速正演计算方法。将地下半空间剖分为大小相等规则排列的一组长方体单元,任意三维复杂形体可以表示成很多不同体积和密度(磁性)的长方体的近似组合。用解析方法计算出所有这些长方体在计算点的重力(磁力)异常,并累加求和,就可以得到整个模型体在计算点引起的重(磁)异常值。为了提高近似程度,需将地下半空间剖分得很细,用传统的CPU串行程序计算相当耗时。GPU在处理能力和存储器带宽上相对CPU有明显优势,采用GPU并行算法,可大大提高计算速度。相关试验结果表明,用GPU实现的正演快速算法计算结果正确,效率明显提高,为重磁异常三维物性反演提供了基础。     

GPU-accelerated iterative solutions for finite element analysis of soil–structure interaction problems

Soil–structure interaction problems are commonly encountered in engineering practice, and the resulting linear systems of equations are difficult to solve due to the significant material stiffness contrast. In this study, a novel partitioned block preconditioner in conjunction with the Krylov subspace iterative method symmetric quasiminimal residual is proposed to solve such linear equations. The performance of these investigated preconditioners is evaluated and compared on both the CPU architecture and the hybrid CPU–graphics processing units (GPU) computing environment. On the hybrid CPU–GPU computing platform, the capability of GPU in parallel implementation and high-intensity floating point operations is exploited to accelerate the iterative solutions, and particular attention is paid to the matrix–vector multiplications involved in the iterative process. Based on a pile-group foundation example and a tunneling example, numerical results show that the partitioned block preconditioners investigated are very efficient for the soil–structure interaction problems. However, their comparative performances may apparently depend on the computer architecture. When the CPU computer architecture is used, the novel partitioned block symmetric successive over-relaxation preconditioner appears to be the most efficient, but when the hybrid CPU–GPU computer architecture is adopted, it is shown that the inexact block diagonal preconditioners embedded with simple diagonal approximation to the soil block outperform the others.     

A GPU parallel computing strategy for the material point method

The material point method (MPM), which is a combination of the finite element and meshfree methods, suffers from significant computational workload due to the fine mesh that is required in spite of its advantages in simulating large deformations. This paper presents a parallel computing strategy for the MPM on the graphics processing unit (GPU) to boost the method’s computational efficiency. The interaction between a structural element and soil is investigated to validate the applicability of the parallelisation strategy. Two techniques are developed to parallelise the interpolation from soil particles to nodes to avoid a data race; the technique that is based on workload parallelisation across threads over the nodes has a higher computational efficiency. Benchmark problems of surface footing penetration and a submarine landslide are analysed to quantify the speedup of GPU parallel computing over sequential simulations on the central processing unit. The maximum speedup with the GPU used is ∼30 for single-precision calculations and decreases to ∼20 for double-precision calculations.