Being capable of predicting seakeeping capabilities in the time domain is of great interest for the marine and offshore industries. However, most computer programs used work in the frequency domain, with the subsequent limitation in the accuracy of their model predictions. The main objective of this work is the development of a time domain solver based on the finite element method capable of solving multi-body seakeeping problems on unstructured meshes. In order to achieve this objective, several techniques are combined: the use of an efficient algorithm for the free surface boundary conditions, the use of deflated conjugate gradients, and the use of a graphic processing unit for speeding up the linear solver. The results obtained by the developed model showed good agreement with analytical solutions, experimental data for an actual offshore structure model, as well as numerical solutions obtained by other numerical method. Also, a simulation with sixteen floating bodies was carried out in an affordable CPU time, showing the potential of this approach for multi-body simulation.