Percolation properties and permeability of heterogeneous and anisotropic three-dimensional fracture networks are studied numerically. Both models composed of mono-disperse plane hexagons. In heterogeneous fracture networks, the fracture centres distributed non-uniformly according to a Gaussian correlated field in space, and the fracture orientations are isotropic. In anisotropic fracture networks, the fractures are directed in an anisotropic way according to Fisher distribution with uniform distribution of their fracture centres. Finite-size-scaling method is used to extrapolate the percolation thresholds of infinite networks. The influences of degree of heterogeneity and angular dispersion parameter on percolation thresholds of fracture networks are analyzed. In heterogeneous fracture networks, increased degree of heterogeneity and correlation length increases the percolation threshold of the networks. For anisotropic networks percolation thresholds of three models of anisotropic three-dimensional (3D) fracture networks have been calculated numerically. The first model is when the fracture networks are comprised of one family of fractures that are distributed in an anisotropic manner about Z-direction. We call this model mono-polar anisotropic fracture network (MFN). The second model is when the fracture networks are comprised of two groups of fractures that are distributed in an anisotropic manner about two orthogonal mean directions, i.e., Z- and X-directions. We call this model bipolar anisotropic fracture network (BFN). The second model is when three groups of fractures are distributed about three orthogonal mean directions, that is Z-, X-, and Y-directions. We call this model tripolar anisotropic fracture network (TFN). The effect of anisotropicity on percolation thresholds in X-, Y-, and Z-directions is investigated. We have revealed that increased anisotropicity in MFN leads into decreased percolation thresholds in both X- and Y-directions, and in these two directions percolation thresholds in anisotropic networks demonstrate a declining trend as anisotropicity goes up. However, in the Z-direction the trend is the opposite. In BFN and TFN as the anisotropicity of networks increases, the percolation thresholds in X-, Y-, and Z-directions span the range of 2.3 to 2.0, where 2.3 and 2.0 are extremums of percolation thresholds for isotropic and non-isotropic orthogonal fracture networks, respectively. The fracture networks are triangulated via an advancing front technique and the macroscopic permeability of both networks are determined by solving the two-dimensional Darcy equation in each fracture. The macroscopic pressure gradient of the system is related to flow problem by an effective permeability; this permeability is averaged over a large number of independent statistical realizations of the fracture .....