In order to avoid selection of the type of functions when using traditional copula methods to fit the joint probability distribution, the copula function and the maximum entropy principle were jointly used(in a method also called the copula entropy method), and the two-dimensional joint distribution functions of flood variates were obtained by solving the copula function with the maximum entropy. With the solved copula function, the two-dimensional joint distribution functions of flood variates(peak discharge, flood volume, and flood duration)were mutually constructed. The Gibbs sampling method and copula functions were used to stochastically simulate trivariate flood events. A case study was conducted using the observed flood data from the Lutaizi Hydrological Station on the Huaihe River. A goodness-of-fit analysis shows that the copula entropy method is effective in fitting multivariate probability distributions and the Gibbs sampling method is effective in simulating trivariate flood events.