The fiber orientation in composite materials is highly dependent on the manufacturing process and plays a crucial role in determining the mechanical, thermal, and electrical properties of the fabricated parts. Several simulation frameworks have been developed to predict the fiber orientation in the printed strands for Material EXtrusion Additive Manufacturing (MEX-AM) process, but studies have been limited to planar and 2D model reductions. Consequently, it is currently impossible to predict the full spatial variation of fiber orientations in a three-dimensional printed strand, and thus the effects of printing conditions remain not well understood. This work seeks to address this issue by introducing the first three-dimensional model capable of simulating the MEX-AM process with fiber-reinforcements. A fully coupled model is developed in this work, which is based on the finite-volume method and solved in the open-source multiphase solver OpenFOAM. Using this framework, we explored the effects of different printing conditions on fiber orientation and strand shape morphology, and we compared our results to experimental observations when possible. We found that the extrusion and nozzle velocities significantly affected the fiber orientation while altering the gap distance between the nozzle exit and substrate had a limited impact. Moreover, increasing anisotropy using longer fibers and higher volume fractions had little influence on the fiber orientation, but their impact on the rheology altered the strand shape considerably. Prediction of fiber orientation for different printing conditions within strands will open the possibility of manufacturing products that has locally engineered properties and tailored anisotropic behavior.