A subsurface hydrocarbon reservoir is simulated by sequential solution of reservoir and well equations to simulate fluid flow inside the reservoir and well production rates. Sequential solution of reservoir and well equations treats wells as specified bottom hole pressure wells. This avoids solving large matrices resulting from the simultaneous solution of the reservoir and well equations which can be computationally very expensive for large number of unknowns and require special sparse matrix solvers. Such sequential solution involves regular reservoir system solvers complemented by small matrices for the numerical solution of the well bottom hole pressures. The solution is performed on tridiagonal matrices for the adjacent reservoir cells to the well cells at the perforated well intervals; and a vector of the unknown reservoir potentials for the adjacent reservoir cells.