Finally I found the problem(s), there were many of them!
1) The biggest one were the initial conditions (thanks E-design!). After selecting "Skip the initial transient bias point calculation" (SKIPBP) in PSpice simulation's options everything got better.
2) Another big problem was in the math: I didn't size reactances well. I did it well, but just barely. In fact assuming A*B(jw0)=1 you barely get a persistent oscillation and a small rounging is enough to let it fade. So I assumed A*B(jw0)=100 when sizing reactances.
3) I had to put another small resistor in the project (R2=1f) because otherwise the simulation didn't converge, but it seems it was too small and I still had convergence problems in some circumstances. So I changed it to 1 mOhm instead.
4) The output capacitor took a long time to load (which wasn't a problem), but also lead to convergence problems in some circumstances. So I changed it to 1nF.
5) I also increased the current generator I1 to 1mA to get the oscillations sooner.
6) Finally I sized the reactances once again, assuming L=C1*10^3 to get more realistic values.
7) Now THD is ~1% for the first 50 harmonics!
Thanks to everybody, I would have never made it if it wasn't for your help!