Unfortunately none of those are going to do what you need.
what I find works well for Functional analysis is to model the output stage resistors, as resistors, and to control the current thru them from the behavior of the model..
but Generally speaking I use integer code (at least in Verilog-A) to model the function of the gate..
so while one can get decent accuracy of the transistion rise time etc from resistor and any line capacitance,
to talk of "harmonics" with "CML gates" doesn't work well, the gate will have a delay, and so you may get some harmonics thru multiple gates connected to the same power, or sharing a common bias node.. (signal to bias coupling I usually ignore..) But the gate itself cannot generate any, as it cant change state until the inputs cross a threshold..
For high speed functional checks I use Verilog-AMS, and am writing a paper on this very topic for the next BMAS conference.
www.bmas-conf.org..
Jonathan