The creation of the finite-element mesh is a very cumbersome process, especially for geometrically complex and large-scale models. On the other hand, the meshfree methods do not require such mesh data in an explicit fashion. Among the meshfree methods, the element-free Galerkin method (EFGM) has attracted interest since the computational implementation of the method is quite similar to that of the finite element procedures. Furthermore, EFG has been reported to be as competent as the finite elements in various types of problems. However, it has been indicated that when compared to FE analyses. EFG computations generally require more time due to the implementation of the moving least-squares approximation scheme. In this study, a hybrid modeling methodology implementing both the FE and EFG procedures for a given domain is reported. The penalty function method is implemented to approximately satisfy the continuity constraints for the displacements at FE-EFG interfaces. A data structure is also proposed for the parallel computation using this hybrid method. Linear equation solver module from GeoFEM software is adopted as the parallel middleware, allowing for efficient and effective parallelization. Effectiveness of the present scheme is shown through the two-dimensional elastostatic problems ranging up to I million DOFs, by comparing the CPU time and parallel efficiency among FE (GeoFEM), EFG and the hybrid methodologies.