This paper reviews the application of the algorithm particle swarm optimization (PSO) to perform stochastic inverse modeling of geophysical data. The main features of PSO are summarized, and the most important contributions in several geophysical fields are analyzed. The aim is to indicate the fundamental steps of the evolution of PSO methodologies that have been adopted to model the Earth’s subsurface and then to undertake a critical evaluation of their benefits and limitations. Original works have been selected from the existing geophysical literature to illustrate successful PSO applied to the interpretation of electromagnetic (magnetotelluric and time-domain) data, gravimetric and magnetic data, self-potential, direct current and seismic data. These case studies are critically described and compared. In addition, joint optimization of multiple geophysical data sets by means of multi-objective PSO is presented to highlight the advantage of using a single solver that deploys Pareto optimality to handle different data sets without conflicting solutions. Finally, we propose best practices for the implementation of a customized algorithm from scratch to perform stochastic inverse modeling of any kind of geophysical data sets for the benefit of PSO practitioners or inexperienced researchers.

We implement the particle swarm optimization (PSO) algorithm for the two-dimensional (2D) magnetotelluric (MT) inverse problem. We first validate PSO on two synthetic models of different complexity and then apply it to an MT benchmark for real-field data, the COPROD2 data set (Canada). We pay particular attention to the selection of the PSO input parameters to properly address the complexity of the 2D MT inverse problem. We enhance the stability and convergence of the solution of the geophysical problem by applying the hierarchical PSO with time-varying acceleration coefficients (HPSO-TVAC). Moreover, we parallelize the code to reduce the computation time because PSO is a computationally demanding global search algorithm. The inverse problem was solved for the synthetic data both by giving a priori information at the beginning and by using a random initialization. The a priori information was given to a small number of particles as the initial position within the search space of solutions, so that the swarming behavior was only slightly influenced. We have demonstrated that there is no need for the a priori initialization to obtain robust 2D models because the results are largely comparable with the results from randomly initialized PSO. The optimization of the COPROD2 data set provides a resistivity model of the earth in line with results from previous interpretations. Our results suggest that the 2D MT inverse problem can be successfully addressed by means of computational swarm intelligence.

We have developed an innovative, simultaneous 1D optimization of electromagnetic (EM) data. Our scheme is suitable for the simultaneous analysis of magnetotelluric (MT) and time-domain EM (TDEM) data based on the probabilistic and evolutionary particle swarm optimization (PSO) algorithm. The simultaneous optimization also identifies and removes the static shift from the MT data. In our scheme, the static shift of the MT apparent resistivity curve is considered as an additional parameter [ ] to be optimized. We tested the suggested method on the synthetic data and then applied it to the data from an EM geophysical study carried out in the geothermal area of Larderello-Travale (Tuscany, Italy). Apart from the novelty of using the PSO algorithm to estimate the model parameters by joint analysis, the simultaneous optimization of the static shift parameter addresses a major problem in MT, i.e., how to define and remove the galvanic effects on MT curves according to independent information, such as that provided by TDEM data. The procedure is expected to strongly influence the application of MT, particularly in geothermal exploration, which commonly relies extensively on EM methods.

We tested synthetic seismic reflection modelling along a seismic line (CROP-18A) in the geothermal field at Larderello (Italy). This seismic line is characterized by a discontinuous but locally very bright seismic marker, named K-horizon, which has been associated with various geological processes, including the presence of fluids at supercritical conditions. Geological and geophysical data were integrated in order to develop a 3D subsurface model of a portion of the Larderello field, where extremely high heat flow values have been recorded. In the study area, the K-horizon is particularly shallow and supercritical deep conditions were accessed at depth. The 2D model of the main geological units up to the K-horizon was extracted from the 3D model along the CROP-18A and used to generate the synthetic TWT stacked seismic sections which were then compared with the observed stacked CROP-18A seismic section. To build the synthetic sections, generated through the exploding reflector approach, a 2D velocity model was created assigning to each pixel of the model a constant P-wave velocity corresponding to the related geological unit. The geophysical parameters and the geological model reconstructions used in the modelling process derive from a multidisciplinary integration process including geological outcrop analogues, core samples, and geophysical and laboratory information. Two geophysical models were used to test the seismic response of the K-horizon, which is associated with (1) a lithological discontinuity or (2) a physically perturbed layer, represented by a randomized velocity distribution in a thin layer. For the latter geophysical model (i.e., the physically perturbed layer), we have tested three different scenarios changing the shape and the thickness of the modelled layer. Despite the reliable calibration implied by the use of homogeneous units, the seismic modelling clearly shows that the physically perturbed layer provides a better explanation of the reflectivity features associated with the K-horizon.