Based on tensor neural networks, we propose a novel machine learning method for solving elliptic partial differential equations with high-dimensional random or parametric coefficients in a bounded physical domain. Utilizing the tensor product structure, we transform the high-dimensional integrals of the tensor neural network functions into one-dimensional integrals, which can be efficiently computed using classical quadrature schemes with high accuracy. This transformation reduces the complexity of high-dimensional integrations from an exponential to a polynomial scale. Numerical examples are provided to demonstrate the accuracy and efficiency of the proposed algorithms.