Does the length of the sidereal day vary systematically? I'm confused about some properties of the sidereal day, in particular whether its duration varies systematically over the course of the year.1 It seems to me that that must be the case, but the details are a bit confusing to me.
I understand that length of a sidereal day is the time interval between successive transits of the line of equatorial longitude corresponding to a specified right ascension, $R$ (usually $0^\text{h}$). Each sidereal day, this line shifts against the earth's direction of rotation by some amount, $\Delta\alpha_{SID}$, to a new position, so that the sidereal day is always shorter than the earth's rotational period on its axis (the "stellar day"). 
What's confusing me is that the length of the sidereal day is often described as if it were constant; but it can't be, nor can the corresponding changes in the equatorial longitude of $R$, $\Delta\lambda_{SID}$, even assuming a constant rate of precession:2 If the sidereal day has a constant length, $\Delta\alpha_{SID}$ must be constant, but then $\Delta\lambda_{SID}$ will vary (between $\Delta\alpha_{SID}/\cos\varepsilon$ and $\Delta\alpha_{SID}\cdot\cos\varepsilon$) as a consequence of the varying equatorial latitude at which the the line of equatorial longitude corresponding to $R$ intersects the ecliptic. But this is impossible, since equatorial longitude should change uniformly, so that $\Delta\lambda_{SID}$ should be the same for each (same-length) sidereal day. The reverse can't be the case either: that each sidereal day, the position of the line corresponding to $R$ shifts by the same fixed amount, $\Delta\lambda_{SID}$, along the ecliptic. If that were the case, $\Delta\alpha_{SID}$ and thus the length of the sidereal day would vary (between $\Delta\lambda_{SID}\cdot\cos\varepsilon$ and $\Delta\lambda_{SID}/\cos\varepsilon$) resulting in a sidereal day of varying length; but if the length of the sidereal day varies $\Delta\lambda_{SID}$ can't be constant.
So the only conclusion I can come to is that both $\Delta\lambda_{SID}$ and the duration of the sidereal day, and thus $\Delta\alpha_{SID}$, vary, subject to the constraints imposed by spherical trigonometry: 
$$\Delta\lambda_{SID}\cdot\cos\varepsilon<\Delta\alpha_{SID}<\Delta\lambda_{SID}/\cos\varepsilon$$
and
$$\langle\Delta\lambda_{SID}\rangle=\langle\Delta\alpha_{SID}\rangle$$
Is this correct?3 Do both of these quantities vary systematically in this way, with these properties.4 Is "the" sidereal day really a mean sidereal day?5

(1) Though a distinction between a "mean sidereal day" and an "apparent sidereal day" is sometimes mentioned, the terminology I've come across is confusing, and descriptions that do discuss variations seem to be referring to different, unspecified, phenomena. Here, I am asking about the properties of an idealized "mean sidereal day", specifically whether even this "mean" day has systematic variation in duration. 
(2) My thinking went like this. Take the plane of the ecliptic as reference. Then precession is just the a rotation of the equatorial coordinate system around an axis perpendicular to that plain. The earth's rotational axis is not collinear with that axis, so that the great circle of the equator intersects the ecliptic plane at two points, one of which can be used as a reference to "track the progress" of precession. Clearly, the progress of that point along the ecliptic — i.e., the change in the intersection's ecliptic longitude, $\Delta\lambda$ — has a direct one-to-one relationship with the rate of precession. What's not clear to me is how this relates to the corresponding change in right ascension, $\Delta\alpha$. But what's being measured in that case is the distance between the reference point and some initial point on the ecliptic, within the equatorial system, an arc that passes through different declinations and thus $\Delta\alpha$ values that vary from $\Delta\lambda$ as described above. (Moreover, this seems to me the only way that $\langle\Delta\lambda\rangle$ and $\langle\Delta\lambda\alpha\rangle$ to be the same, which must be the case for the method described in (3) to work.)
(3) I'm reasonably confident that the first is true, but I'm not sure how to arrive at the second; though it must also be true, since otherwise it would not be possible to determine the precessional period (in sidereal days) from $24^\text{h}/\langle\Delta\alpha_{SID}\rangle$, as is commonly done.
(4) The corresponding (by definition, assuming constant rotational rate for the earth) constant quantities for the stellar day must thus be
$$\Delta\lambda_{ST}=\langle\Delta\lambda_{SID}\rangle\cdot\frac{d_{ST}}{d_{SID}}=\Delta\alpha_{ST}=\langle\Delta\alpha_{SID}\rangle\cdot\frac{d_{ST}}{d_{SID}}$$
where ${d_{SID}}$ is the length of a mean sidereal day, and ${d_{ST}}$ the length of a stellar day.
(5) Or perhaps a mean mean sidereal day. (A very mean sidereal day?)

 A: The sidereal day is a mean quantity, derived from many years of observations (it goes back to Halley times). But there are additional variations due to changes in the orbital parameters of the Earth. This changes are partly well known and predictable by celestial mechanics, but a significant part is unpredictable.
Google for a curve called Polar Motion. It depicts the movement of the north pole over the surface of the Earth through the years. Look for both the theoretical and the real one (build from measures). You will see the theoretical one is spiral-alike, but the measured "real" curve has a strongly erratic component superposed. It is even affected by earthquakes, and so is the exact position of the vernal equinox...
It has no sense to fine tune the determination of each sidereal day. It would be unrealistic, since the erratic, unpredictable contributions are in most cases bigger than the calculated corrections to the mean sidereal day.
Apart from that, I think there are some points in your reasoning that may be improved:
"Transit" is a well defined thing in astronomy. A line cannot transit. Moreover, the sidereal day is defined as the mean time between two transits of the vernal equinox. That point is, by definition, the zero of right ascension, always.
What you refer as "equatorial longitude" and "equatorial latitude" is very probably right ascension and declination, but I suggest that you carefully check that definitions, especially the sense (clockwise/counterclockwise) of your equatorial longitude vs right ascension, or you might run into a lot of problems when comparing your results with the literature.
In your calculations a factor cosε appears and, although you may be doing something strictly correct from the geometrical point of view, I sincerely don't understand how you came to include the obliquity of the ecliptic in that calculations. You have probably projected some non constant movement from the ecliptic into the equator or something similar. Nobody can take you from doing that, but I think it leads unnecessarily to difficulties.
The ecliptic has nothing to do here. As long as we are not examining the movement of the real Sun, the ecliptic plays no role at all. It is the celestial equator the one that matters, since it is defined as the normal plane to the Earth rotation axis. There, in the celestial equator, all hour angles (and not "equatorial longitudes" or right ascensions, which are fixed) move with constant angular speed (as long as you neglect the 26,000 year periodic equinox precession and the regular + irregular movements of the north pole registered in the polar curve).
I kindly suggest you carefully learn the basic definitions and do lots of drawings. Thereafter, deriving relations will be safer.
A: I've just finished writing scripts for GMST and GAST, and they checked against the GAVO and JPL models. YAY! Perhaps what I learnt may be of interest.
GMST is a function of UT1, and could be described as the smoothed (filtered) measure of a sidereal day wrt UT1 (a mean solar day). It includes precession. There are t^2 t^3 t^4 t^5 terms. So, yes, in a tiny amount, the length of a sidereal day, compared to a mean solar day, is decreasing....
GST or GAST or "Angle of Aries" could be described as the unfiltered measure of an actual sidereal day wrt UT1. It includes precession and nutation. There are the precession t^2 t^3 t^4 t^5 terms, and a Fourier series to account for interactions with the moon and sun. There are also expressions for planets, but they have a tiny effect. My Fourier series used the dominant 10 frequencies.
There are 2 dominant frequencies in nutation.
The biggest is due to the regression of the Moon's nodal line, and has a period of 6798.383 days. The second biggest has a period of 182.621 days.
So I guess the answer to your question is yes again, there are many sinusoidal variations acting together to vary the length of a sidereal day.
This webpage has an excellent explanation of time
https://www.cfa.harvard.edu/~jzhao/times.html
If anyone wants to go down the nutation rabbit hole, I recommend the Meeus Astronomical Algorithms text, which looks at the 1980 precession and nutation models.
Meeus also explores a dramatic simplification to nutation modelling, by picking the dominant four frequencies, and neglecting t^2 t^3 t^4 terms in the nutation fundamental arguments (known as Delaunay variables).
If you're really interested, look at the IERS conventions from 1992, 2003, and 2010. Also called "technical notes" no's 21, 32, and 36. 1992 is good because it gives a commentary on the development of these models. 2003 is good because it has Tables 5.3a and 5.3b set out to point you towards dominant frequencies. The 2010 is terrible to read, and doesn't contain ALL of the information you need to write a nutation script.
You can download the matrix for a nutation Fourier series from the IERS or SOFA websites. There is an IAU 2000B Abridged model, with 77 equations. I found it produced a mean absolute error less of 0.6 milliarcseconds, and a maximum error of 2.5 milliarcseconds.
It surprised me how quickly these models converge, and that with only 5 frequencies you can achieve a mean absolute error of 60 milliarcseconds, and a maximum error of 200 milliarcseconds.
And if you want to be rough and ready, an optimised linear model with one mode of vibration will give a mean absolute error less than 1", and a maximum error of less than 2". Try this formula...
(T is days since J2000.0)
ARIES= 280.460573384848 + 360.985647373871.T - 15.7623".sin(125.344555 - T.(360/6798.383))
A: The sidereal day, $1/k$, is a mean quantity derived from many of observations. Apart from the consideration you are deriving, there are variations due to changes in the orbital parameters of the Earth. These changes are partly well known and predictable by celestial mechanics, but a significant part is unpredictable.
One thing that, critically, does not vary, is the "specified right ascension" you refer to in your question: that is always $0^\text{h}$, the vernal equinox. I believe this is the source of your confusion.
The basics of your reasoning are sound. It is in fact true that at any given moment, the rate at which the meridian sweeps through ecliptic longitude will be different from the rate at which it sweeps through right ascension, depending on the declination of the segment of the ecliptic the meridian is passing through at that time. And you're right that that rate will vary systematically between 
$\dot{\alpha}\cdot\cos\varepsilon$ and $\dot{\alpha}/\cos\varepsilon$. So, as a result, as you state, the omission of different "chunks" of the ecliptic from a complete circuit, will result in different day lengths.
However, this is irrelevant to the length of a sidereal day. The "chunk" of the ecliptic by which a sidereal day differs from a complete ecliptic circuit is always the same: in particular, it is always at the same declinations, those just west of the vernal equinox. As a result, the effect of precession on the length of the day is the same each day, and the systematic variations you describe have no effect.
The variation you describe does indeed occur for the time it takes for the meridian to return to the spot occupied by the sun at the start of each day, because the declination of that spot, and thus of the segment of the ecliptic omitted by precession, varies. But even in that case, those segments follow the sun through precisely the full circuit of the ecliptic over the course of a tropical year, so that the average duration of these "days" is $1/k$. Note that the same process applies to the time between successive transits of the sun. This effect is much more significant (since the sun several orders of magnitude faster than the equinox precesses), and is in fact the component of the equation of time due to the obliquity of the ecliptic.  
