A three-dimensional constraint-driven dynamic rigid-link numerical model of a flapping wing structure with compliant joints (CJs) called the dynamic spar numerical model is introduced and implemented. CJs are modeled as spherical joints with distributed mass and spring-dampers with coupled nonlinear spring and damping coefficients, which models compliant mechanisms spatially distributed in the structure while greatly reducing computation time compared to a finite element model. The constraints are established, followed by the formulation of a state model used in conjunction with a forward time integrator, an experiment to verify a rigid-link assumption and determine a flapping angle function, and finally several example runs. Modeling the CJs as coupled bi-linear springs shows the wing is able to flex more during upstroke than downstroke. Coupling the spring stiffnesses allows an angular deformation about one axis to induce an angular deformation about another axis, where the magnitude is proportional to the coupling term. Modeling both the leading edge and diagonal spars shows that the diagonal spar changes the kinematics of the leading edge spar verses only considering the leading edge spar, causing much larger axial rotations in the leading edge spar. The kinematics are very sensitive to CJ location, where moving the CJ toward the wing root causes a stronger response, and adding multiple CJs on the leading edge spar with a CJ on the diagonal spar allows the wing to deform with larger magnitude in all directions. This model lays a framework for a tool which can be used to understand flapping wing flight.