Three-dimensional elastoplastic finite element (FE) analysis models considering a damage zone and cleavage were developed, including plastic behavior of embedment tests parallel and perpendicular to the grain at ambient and high temperatures, in order to study the embedment behavior of timber. In a numerical analysis of compression tests to define material properties of timber, the post-yield behavior of the X-axis in the axial direction was larger than the test values when the post-yield slope of the test result were set to the plastic modulus. The reason was that the stress on the Y-axis and Z-axis increased after the yield point; thus, the equivalent yield stress decreased. In the numerical analysis of the embedment tests, the effect of the damage zone must be considered in the analytical models because the initial stiffness of the numerical results in compression parallel to the grain was considerably higher than that of the test results when the damage zone was not considered. In the numerical analysis of the embedment tests in compression perpendicular to the grain, the two models reproduced the stiffness and behavior after the yield point of the test results. The model that considered cleaving showed the stress distribution after cleavage.