A challenge common to weather, climate and seasonal numerical prediction is the need to simulate accurately long range transport and reversible isentropic processes in combination with appropriate determination of sources/sinks of energy and entropy. With respect to the entropy of matter, this includes the distribution, transport and transformation of internal, gravitational and kinetic energies, the energies of water substances in all forms, and the related thermodynamic processes of phase changes involved with clouds including condensation, evaporation, precipitation, and cloud radiation interaction. A means to study a model's accuracy in simulating internal hydrologic processes is to determine its capability to simulate the appropriate conservation of potential and equivalent potential temperature as surrogates of dry and moist entropy under reversible moist adiabatic processes in which clouds form, evaporate and precipitate. An analysis of variance of model error from numerics is developed in which the sums of squares of the simulated differences between equivalent potential temperature and its proxy is partitioned in three components: the square of the deviations of differences from area mean difference, the square of the area mean deviation differences from the global mean difference and the square of the global difference. The three components provide information on the conservation of moist entropy and energy in relation to growth of bias and random components of model error horizontally, vertically and globally. The analysis of variance of the error sums of squares may also be applied to any regional climate domain. Results from a series of six global ten day simulations employing various versions of CCM2 and 3 --- all Eulerian spectral numerics, all semi-Lagrangian numerics, mixed Eulerian spectral and semi-Lagrangian numerics, --- and two UW hybrid isentropic models over a ten day period will be presented to detail differences in the simulation of reversible processes that develop from nonlinear numerical inaccuracies involving the exchange of energy and water substances. The partitioning of the global sum of squares isolates a structure of bias differences associated with the spurious vertical exchange of entropy and energy in all of the models except the UW isentropic-eta models. Furthermore, the empirical relative frequency distributions of the differences from the UW models were symmetric, triangular in form and equilibrated as expected from statistical theory in which the random variate is given by the difference of two variates, each of which is drawn from a rectangular distribution of rounding errors. With judicious modification to the Lagrangian source for the trace continuity equation, the strategy may be extended to study the accuracy of a model in simulating reversibility as different dry and moist convective parameterization algorithms are introduced and tested in model development. The strategy developed herein extends to potential vorticity and total energy.