Physical phenomena are observed in many fields (science and engineering) and are often studied by time-consuming computer codes. These codes are analyzed with statistical models, often called emulators. In many situations, the physical system (computer model output) may be known to satisfy inequality constraints with respect to some or all input variables. The aim is to build a model capable of incorporating both data interpolation and inequality constraints into a Gaussian process emulator. By using a functional decomposition, a finite-dimensional approximation of Gaussian processes such that all conditional simulations satisfy the inequality constraints in the entire domain is proposed. To show the performance of the proposed model, some conditional simulations with inequality constraints such as boundedness, monotonicity or convexity conditions in one and two dimensions are given. A simulation study to investigate the efficiency of the method in terms of prediction is included.