A method for modeling flow properties over a series of time increments of a reservoir in an earth formation comprising a plurality of fractures, the method comprising: obtaining information about the plurality of fractures using a data acquisition tool; using a processor for: building a three-dimensional stress field representing stresses in the reservoir; building a three-dimensional natural fracture model, NFM, comprising fracture flow properties of the plurality of fractures using the obtained information or a change to the three-dimensional stress field; running a flow simulation of the reservoir for a time increment using the NFM to model the flow properties of the reservoir for that time increment; computing a latest change in the three-dimensional stress field from the flow simulation; and incrementing the time increment and iterating the building the NFM using the latest change in the three-dimensional stress field, the running of the flow simulation using the latest NFM to model the flow properties, and the computing a latest change in the three-dimensional stress field from the latest running of the flow simulation, wherein the incrementing continues for the series of time increments; wherein building the three-dimensional NFM further comprises establishing a three-dimensional mesh comprising a plurality of cells to represent the reservoir; and the method is characterised by calculating a fracture intensity at frictional failure for each fracture in each cell using the obtained information and the latest change in the three-dimensional stress field, wherein the fracture flow properties in the NFM comprise the fracture intensity at frictional failure for each fracture, and wherein the fracture intensity at frictional failure of a cell is the frequency of fractures in that cell weighted by each of the fracture's proximity to frictional failure. obtaining information about the plurality of fractures using a data acquisition tool; using a processor for: building a three-dimensional stress field representing stresses in the reservoir; building a three-dimensional natural fracture model, NFM, comprising fracture flow properties of the plurality of fractures using the obtained information or a change to the three-dimensional stress field; running a flow simulation of the reservoir for a time increment using the NFM to model the flow properties of the reservoir for that time increment; computing a latest change in the three-dimensional stress field from the flow simulation; and incrementing the time increment and iterating the building the NFM using the latest change in the three-dimensional stress field, the running of the flow simulation using the latest NFM to model the flow properties, and the computing a latest change in the three-dimensional stress field from the latest running of the flow simulation, wherein the incrementing continues for the series of time increments; wherein building the three-dimensional NFM further comprises establishing a three-dimensional mesh comprising a plurality of cells to represent the reservoir; and the method is characterised by calculating a fracture intensity at frictional failure for each fracture in each cell using the obtained information and the latest change in the three-dimensional stress field, wherein the fracture flow properties in the NFM comprise the fracture intensity at frictional failure for each fracture, and wherein the fracture intensity at frictional failure of a cell is the frequency of fractures in that cell weighted by each of the fracture's proximity to frictional failure.