The occurrence and properties of natural faults and fractures in geothermal reservoirs are key in determining reservoir flow properties, and thereby the performance of geothermal doublets placed in fractured reservoirs or in the vicinity of fault zones. In this paper, an analytical model is presented that describes the 3D non-isotropic permeability of a geothermal reservoir around a major fault zone, taking into account typical fault architectures consisting of a fault core, a damage zone and surrounding intact reservoir matrix. The sensitivity of model predictions to orientations of sedimentary layers, damage zone fractures and fault core is analysed for typical permeability contrasts and dimensions of intact reservoir, damage zone fractures and fault core. The model can be used to determine optimum orientation of geothermal doublets around fault zones, taking into account the distribution and characteristics of faults, fractures and sedimentary layering. Implications for optimizing the design of geothermal doublets placed in the vicinity of large fault zones are given.