In structural analysis and design it is important to consider the effects of uncertainties in loading and material properties in a rational way. Uncertainty in material properties such as heterogeneity in elastic and mass properties can be modeled as a random field. For computational purpose, it is essential to discretize and represent the random field. For a field with known second order statistics, such a representation can be achieved by Karhunen-Lo`eve (KL) expansion. Accordingly, the random field is represented in a truncated series expansion using a few eigenvalues and associated eigenfunctions of the covariance function, and corresponding random coefficients. The eigenvalues and eigenfunctions of the covariance kernel are obtained by solving a second order Fredholm integral equation. A closed-form solution for the integral equation, especially for arbitrary domains, may not always be available. Therefore an approximate solution is sought. While finding an approximate solution, it is important to consider both accuracy of the solution and the cost of computing the solution. This work is focused on exploring a few numerical methods for estimating the solution of this integral equation. Three different methods:(i)using finite element bases(Method1),(ii) mid-point approximation(Method2), and(iii)by the Nystr¨om method(Method3), are implemented and numerically studied. The methods and results are compared in terms of accuracy, computational cost, and difficulty of implementation. In the first method an eigenfunction is first represented in a linear combination of a set of finite element bases. The resulting error in the integral equation is then minimized in the Galerkinsense, which results in a generalized matrix eigenvalue problem. In the second method, the domain is partitioned into a finite number of subdomains. The covariance function is discretized by approximating its value over each subdomain locally, and thereby the integral equation is transformed to a matrix eigenvalue problem. In the third method the Fredholm integral equation is approximated by a quadrature rule, which also results in a matrix eigenvalue problem. The methods and results are compared in terms of accuracy, computational cost, and difficulty of implementation. The first part of the numerical study involves comparing these three methods. This numerical study is first done in one dimensional domain. Then for study in two dimensions a simple rectangular domain(referred toasDomain1)is taken with an uncertain material property modeled as a Gaussian random field. For the chosen covariance model and domain, the analytical solutions are known, which allows verifying the accuracy of the numerical solutions. There by these three numerical methods are studied and are compared for a chosen target accuracy and different correlation lengths of the random field. It was observed that Method 2 and Method 3 are much faster than the Method 1. On the other hand, for Method 2 and 3, additional cost for discretizing the domain into nodes should be considered whereas for a mechanics-related problem, Method 1 can use the available finite element mesh used for solving the mechanics problem. The second part of the work focuses on studying on the effect of the geometry of the model on realizations of the random field. The objective of the study is to see the possibility of generating the random field for a complicated domain from the KL expansion for a simpler domain. For this purpose, two KL decompositions are obtained: one on the Domain1, and another on the same rectangular domain modified with a rectangular hole (referredtoasDomain2) inside it. The random process is generated and realizations are compared. It was observed from the studies that probability density functions at the nodes on both the domains, that is, on Domain 1 and Domain 2, are similar. This observation leads to a possibility that a complicated domain can be replaced by a corresponding simpler domain, thereby reducing the computational cost.