Make the stationary center of your coordinate system the stationary center of mass of the two bodies (called the "barycenter"). Then by conservation of momentum, the total momentum of the system will always be zero. If the net momentum of the two bodies is not zero, then first subtract the velocity of their center of mass from both to make it stationary.
What you normally see as the orbit solution is in reality the vector difference between the positions of the two bodies. If one of the bodies is effectively not moving because it is so much larger than the other, then you simply have that solution as the motion of the smaller body.
For the general case of the two-body problem, you solve for the vector difference of the two positions using the sum of the two masses as the central mass, getting the usual Kepler solution for that vector difference. Then to get the motions of the two actual bodies, you solve for two equations in two unknowns, where the difference of the two positions is your solution, and the sum of the momentums is zero.
Try this simple example in two dimensions:
$$\vec{r_1}=\left(2,0\right),\,\,\vec{r_2}=\left(-1,0\right)$$
$$\vec{v_1}=\left(0,2\right),\,\,\vec{v_2}=\left(0,-1\right)$$
$$\mu_1=6,\,\,\mu_2=12$$
where $\mu$ is the notation for the mass times Newton's gravitational constant $G$. Verify that the CG position and velocity are at the origin. Then determine the orbit for:
$$\vec{r}=\vec{r_1}-\vec{r_2}$$
$$\vec{v}=\vec{v_1}-\vec{v_2}$$
$$\mu=\mu_1+\mu_2$$
where $\vec{r}$ and $\vec{v}$ are the initial position and velocity of a massless object orbiting a central body with mass times $G$ of $\mu$.
You get $a=6$, $e=1/2$, $\omega=0$, and $\nu=0$. We then have:
$$\vec{r}=\left(6\left(-{1\over 2}+\cos\tau\right),6{\sqrt{3}\over 2}\sin\tau\right)$$
where $\tau$ is the eccentric anomaly. The time equation is the same for both the zero mass and general two-body cases:
$$t=2\sqrt{3}\left(\tau-{1\over 2}\sin\tau\right)$$
Solving the simultaneous equations, we get in general:
$$\vec{r_1}={\mu_2\over\mu}\vec{r},\,\,\vec{r_2}=-{\mu_1\over\mu}\vec{r}$$
That gives:
$$\vec{r_1}=\left(-2+4\cos\tau,\,\,2\sqrt{3}\sin\tau\right)$$
$$\vec{r_2}=\left(1-2\cos\tau,\,\,-\sqrt{3}\sin\tau\right)$$
Now we can plot it:

Update for updated question:
Here is the same Newtonian solution for PSR J0737-3039, where this presentation provides the masses, period, and eccentricity:

The scale on the plot is km, and the period is 2.45 hours (the animated gif runs a little faster than that).
The eccentricity is small, 0.0878, and the masses are close to equal, where A is 1.338 and B is 1.249 solar masses. As a result, they are in close to a common circular orbit about their center of mass. (A is blue and B is orange.) The positions in km and time in seconds are:
$$\vec{r_A}=\left(424000\,(-0.0878+\cos\tau),\,\,422000\sin\tau\right)$$
$$\vec{r_B}=\left(-454000\,(-0.0878+\cos\tau),\,\,-452000\sin\tau\right)$$
$$t=1404\,(\tau-0.0878\sin\tau)$$
where $\tau$ goes from $0$ to $2\pi$.
Their speed in their orbits is an impressive 0.1% the speed of light.