Modeling a Partially Rigid Molecule

I am not aware of any restriction in fix rigid which prevents atoms
from two different rigid bodies from interacting with each other via
bonded interactions (bonds, angle, dihedral, or improper
interactions).
As far as LAMMPS is concerned, I would guess that it just treats these
as other forces acting on these particles (just as it would treat
pairwise repulsive forces when they bump into each other). It seems
you are aware of this.

In your case, the first rigid body could be the
benzene-ring+oxygen+carbon, and the second rigid body could be the 3
remaining hydrogens. (Or you could make the second body the entire
CH3 group.)

However if, in addition to this, you want the length of the bond(s)
connecting them to be fixed (and/ the bond-angle to be frozen), then
that could be a problem. I don't know if you can combine fix shake
with fix rigid to constrain the distance of a bond connected two rigid
bodies (for example). It might not work, but have you tried it?

But the easiest thing to do is us fix rigid, and use a fairly stiff
bond (and bond-angle) connecting the two rigid bodies. You thought
this was inelegant, but how stiff do you really need it be? You can
probably get away with a bond which is 4 times stiffer (K is 4 times
larger) than you normally would use. This does not cost you anything.
This is because fix rigid forces you to use smaller timesteps (0.5fs
instead of 1fs). If you need something stiffer, lower the timestep.

This isn't the answer you wanted, but I doubt LAMMPS can do it another
way (unless you were somehow able to combine fix rigid with shake).

If I'm wrong, hopefully somebody will correct me
Good luck!

andrew