Reactive fluid flow and geochemical transport in unsaturated fractured rock s have received increasing attention for studies of contaminant transport, ground-water quality, waste disposal, acid mine drainage remediation, miner al deposits, sedimentary diagenesis, and fluid-rock interactions in hydroth ermal systems, This paper presents methods for modeling geochemical systems that emphasize: (1) involvement of the gas phase in addition to liquid and solid phases m fluid flow, mass transport, and chemical reactions; (2) tre atment of physically and chemically heterogeneous and fractured rocks, (3) the effect of heat on fluid flow and reaction properties and processes, and (4) the kinetics of fluid-rock interaction. The physical and chemical proc ess model is embodied in a system of partial differential equations for flo w and transport, coupled to algebraic equations and ordinary differential e quations for chemical interactions. For numerical solution, the continuum e quations are discretized in space and time. Space discretization is based o n a flexible integral finite difference approach that can use irregular gri dding to model geologic structure; time is discretized fully implicitly as a first-order finite difference. Heterogeneous and fractured media are trea ted with a general multiple interacting continua method that includes doubl e-porosity, dual-permeability, and multi-region models as special cases, A sequential iteration approach is used to treat the coupling between fluid f low and mass transport on the one hand, chemical reactions on the other. Ap plications of the methods developed here to variably saturated geochemical systems are presented in a companion paper (part 2, this issue).