The symmetry of the waveform can be exploited for additional precision for a similar memory and computational footprint. Make the 256-entry LUT of half of a sine lobe rather than a full period then do some bookkeeping to track where in the period you are for reflections of LUT index and y value.

Example (1-indexed):

   f  = 1000;
   tPer = 1 / f;
   ts = 1 / 48000;
   LUT = int16( 2^15 * sinpi( (1:256) / 256 / 2 ) );
   
   t = 0;
   while true
       t = mod( t + ts, tPer);
       i = uint16( floor( t * f * 1024 ) );
       
       if i <= 256
           y = LUT(i);
       elseif i <= 512
           y = LUT(257-mod(i,257));
       elseif i <= 768
           y = -LUT(mod(i,257));
       else
           y = -LUT(257-mod(i,257));
       end
   end

I left out interpolation to keep the concept highlighted.

Unsurprisingly: it's faster and more accurate to skip the LUT when using MATLAB. It's faster still to use vectors and batched function calls.