How does Overlap-Add work for IIR filter?
According to wiki, it seems that overlap add technique usually works for FIR filters. I hasn't read the link you posted and the matrix form for IIR before. But when it comes to buffer for IIR, I have an alternative (hope simpler) method which might works in your case.


  If I want to process a buffer of audio with the IIR filter coming from the LPC analysis, the IR should be infinite, so the only way I can think of using a block convolution is by convolving 


The concept is filtering in sections and utilizing the memory in filters. Let's first look at the testing code and usage with following steps,


Dividing the x in to x1 and x2, 
Filtering x1 and x2 seprately (with memories in filterszf)
Combining the ouput together. 


The results confirm that filtering in section can be the same as x directly.

%% Filter in section
x = randn(10000,1); % psuido signal

x1 = x(1:5000);
x2 = x(5001:end);

b = [2,3]; 
a = [1,0.2]; % IIR filter coeffs

[y1,zf] = myFilter(b,a,x1); % The zf is the memory in IIR filter!
y2 = myFilter(b,a,x2,zf); % assign the zf when filter the next section.

y = myFilter(b,a,x); % filtering x directly.

isequal(y,[y1;y2]) % Verification 


Here is the implementation of myFilter Function. In fact, myFilter imitates(implements) the built-in function fitler in matlab. The implementation architecture is Direct form II, which has advantages over least usage of delay operator and disadvantages over not free for numerical overflow for fixed points.

function [Y, z] = myFilter(b, a, X, z)
    n = max(length(a),length(b));
    a = [a zeros(1,n-length(a))];
    b = [b zeros(1,n-length(b))];
    z(n) = 0;     
    b = b / a(1); 
    a = a / a(1);
    Y    = zeros(size(X));
    for m = 1:length(Y)
       Y(m) = b(1) * X(m) + z(1);
       for i = 2:n
          z(i - 1) = b(i) * X(m) + z(i) - a(i) * Y(m);
       end
    end
    z = z(1:n - 1);
end