In order to obtain a better understanding of the origin of the hydration force, three molecular dynamic simulations of phospholipid/water multilamellar systems were performed. The simulated systems only differed in the amount of interbilayer water, ranging from the minimum to the maximum amount of swelling in the liquid-crystal phase. The analysis of orientational polarization, hydrogen bonding and diffusion rates of the water molecules between the membranes reveals a strong perturbing effect, which decays smoothly and approximately exponentially (with a decay length of about 0.25 nm) toward the middle of the water layer. The electrostatic potential profiles show that the decay of water ordering is directly correlated with the decay of the interfacial dipolar charges. Therefore, the propagation of the ordering of water molecules is not intrinsic to water but is merely determined by the local field from the interfacial charge distribution. The decay of the electrostatic potentials appears to be as a stretched exponential in all simulated systems. This type of decay can be interpreted as normal exponential but with a varying local decay length. We speculate that it is the fractal nature of the membrane surface which results in the stretched exponential behavior. The hydration force, resulting from the ordering of the water molecules between the membranes, will also exhibit this type of decay and is thus dependent on solvent as well as membrane structure.