Building models of a geological system using potential field is almost a routine exercise, where both gravity and magnetic data are used to interpret the architectural framework of a geological setting. Although nonuniqueness in potential field interpretation is inevitably a major issue (Skeels, 1947, Roy, 1962), it is argued by many authors (Pedersen, 1979; Bosch et al., 2006) that the use of data redundancy, for example, joint use of gravity and magnetic data as well as appropriate constraints in the interpretation/inverse modelling overcomes the ambiguity limitation of potential field interpretation. However, joint use of gravity and magnetic data especially for geometric inversion is a non-trivial task. By geometric inversion of potential field data we mean determining the geometric configuration of the modelled body whose physical properties, such as density/susceptibility are known. Such potential field inversions are inherently nonlinear. The presence of nonlinearity between the data and model space often results in an apparent paradoxical situation, offering a less stable and thus a less reliable inverse model while exercising the data redundancy by joint use of gravity and magnetic data. We address the issue of producing a stable solution for a 2D joint geometric inversion of gravity and magnetic data via the ?Pareto optimal? criterion. To this end, we state that the joint inversion of gravity and magnetic data can be realised via a multi-objective optimisation problem (MOP), where the objective functions related to the best fitting gravity and magnetic data are satisfied simultaneously. Very recently, MOP has been introduced by Kozlovskaya et al. (2007) in seismic anisotropy studies of lithosphere by jointly inverting shear wave splitting and P-wave residual data. With regard to the MOP our goal is not to determine a unique, possibly globally optimal, solution, but to determine a class of optimal solutions defining a compact set. Such a class of solutions is known as a Pareto set or Pareto front, where each solution is optimal in a Pareto sense; meaning there exists no other solution which will improve one objective function without degrading the other objective function (Coello Coello et al., 2002). To determine a Pareto front we use the epsilon-constraint method, where one objective is minimised while the other one is constrained with a bound by a pre-selected value of epsilon. To minimise the objective function we use a Particle Swarm Optimization (PSO) scheme, which is a global stochastic optimisation scheme. The rationale in selecting PSO is that it is easy to implement and is computationally inexpensive since its memory and its CPU speed requirements are low (Jones, 2006). PSO has been used for both continuous and discrete optimisation problems (Clerc and Kennedy, 2002). The application of PSO in geophysical inversion is demonstrated by Shaw and Shrivastava, (2007). In this study we discuss the application of MOP and PSO in the joint inversion of gravity/magnetic data as well as testing the method through numerical experiment using synthetically generated gravity and magnetic data for a simple but realistic geological model.