In this chapter we present Kriging— also known as a Gaussian process (GP) model— which is a mathematical interpolation method. To select the input combinations to be simulated, we use Latin hypercube sampling (LHS); we allow uniform and non-uniform distributions of the simulation inputs. Besides deterministic simulation we discuss random simulation, which requires adjusting the design and analysis. We discuss sensitivity analysis of simulation models, using "functional analysis of variance" (FANOVA)— also known as Sobol sensitivity indexes. Finally, we discussoptimization of the simulated system, including "robust" optimization.