Physics of the Earth and Planetary Interiors, manuscript No. (will be inserted by the editor) 3D seismic modeling in geothermal reservoirs with a distribution of steam patch sizes, permeabilities and saturations, including ductility of the rock frame Short title: Seismic wave simulation in geothermal reservoirs Jos´e M. Carcione · Flavio Poletto · Biancamaria Farina · Cinzia Bellezza Received: date / Accepted: date Abstract Seismic propagation in the upper part of the crust, where geothermal reser- voirs are located, shows generally strong velocity dispersion and attenuation due to varying permeability and saturation conditions and is affected by the brittleness and/or ductility of the rocks, including zones of partial melting. From the elastic-plastic as- pect, the seismic properties (seismic velocity, quality factor and density) depend on effective pressure and temperature. We describe the related effects with a Burgers me- chanical element for the shear modulus of the dry-rock frame. The Arrhenius equation combined to the octahedral stress criterion define the Burgers viscosity responsible of the brittle-ductile behaviour. The effects of permeability, partial saturation, varying porosity and mineral compo- sition on the seismic properties is described by a generalization of the White mesoscopic- loss model to the case of a distribution of heterogeneities of those properties. White model involves the wave-induced fluid flow attenuation mechanism, by which seismic waves propagating through small-scale heterogeneities, induce pressure gradients be- tween regions of dissimilar properties, where part of the energy of the fast P-wave is converted to slow P (Biot)-wave. We consider a range of variations of the radius and size of the patches and thin layers whose probability density function is defined by different distributions. The White models used here are that of spherical patches (for partial saturation) and thin layers (for permeability heterogeneities). The complex bulk modulus of the composite medium is obtained with the Voigt-Reuss-Hill average. Effective pressure effects are taken into account by using exponential functions. We then solve the 3D equation of motion in the space-time domain, by approximat- ing the White complex bulk modulus with that of a set of Zener elements connected in series. The Burgers and generalized Zener models allows us to solve the equations with a direct grid method by the introduction of memory variables. The algorithm uses the Fourier pseudospectral method to compute the spatial derivatives. It is tested against an analytical solution obtained with the correspondence principle. We consider two main cases, namely the same rock frame (uniform porosity and permeability) satu- Istituto Nazionale di Oceanografia e di Geofisica Sperimentale (OGS), Borgo Grotta Gigante 42c, 34010 Sgonico, Trieste, Italy. E-mail: jcarcione@inogs.it