function [R]=SymICrateOEN(g,snr,K)
%This function calculates the achievable symmetric rate from the paper
%"The Approximate Sum Capacity of the Symmetric Gaussian K-User Interference Channel"
%by Or Ordentlich, Uri Erez and Bobak Nazer

%The assumed channel is a (real) Gaussian symmetric K-user interference channel
%y_k=x_k+g\sum_{j\neq k}x_j+z_k
%where z_k is AWGN with aero mean and variance 1, and all K users are
%subject to a power constraint of snr.

%Functions inputs:
%g - interference gain
%snr - signal-to-noise ratio in LINEAR scale (not dB!)
%K - number of users

%The functions output is a lower bound on the symmetric capacity.

h=[1 g]';

%First coding strategy - single stream per user (no rate-splitting)

%finding the two optimal equations and computation rates
B=diag([1 K-1]);
[R,A]=CoFtransform(h,snr,B);
R_comp1=R(1);
R_comp2=R(2);
%finding the achievable symmetric rate
if A(1,2)==0%best equation to decode is a_1=[1 0], in this case decoding this equation is sufficient in order to get the desired codeword
    R_strat1=R_comp1;
else%must decode two equations in order to get the desired codeword
    R_strat1=R_comp2;
end

%Second coding strategy - two stream per user (lattice Han-Kobayashi)
gamma=1/sqrt(g^2*snr);
if gamma<1%Otherwise the power allocation we use is not valid
    %calculating effective channel coefficients when lattice Han-Kobayashi is applied
    layVec=[sqrt(1-gamma^2) gamma];
    kappa=1/sqrt(1+snr*gamma^2*g^2*(K-1));
    h=[1*layVec g*layVec(1)]'*kappa;
    %finding the three optimal equations and computation rates
    B=diag([1 1 K-1]);
    [R,A]=CoFtransform(h,snr,B);
    R_comp1=R(1);
    R_comp2=R(2);
    R_comp3=R(3);
    
    R_strat2=R_comp2+R_comp3;
else
    R_strat2=0;
end

if ((g^2>snr^(-1/2))&(g^2<1))
    R=R_strat2;
else
    R=R_strat1;
end


end

function [R,A]=CoFtransform(h,snr,B)

%This function calculates the compute-and-forward transform for the K-user MAC
%channel y=h^T x+z, where the (column) vector h holds the channel
%coefficients, the vector x=[x_1 ... x_K]^T holds the K channel inputs, and
%z is AWGN with zero mean and unit variance. All K user are subject to the
%same power constraint of snr.

%The function can also calaculate the compute-and-forward transform of an
%effective MAC (i.e., a MAC where several users have integer valued
%coefficients and are therefore aligned to one effective user). In this
%case the matrix B is the effective weight matrix.

%Functions inputs are:
%h - a column vector holding the K channel coefficients
%snr - signal-to-noise ratio in LINEAR scale (not dB!)
%B - the effective weight matrix. B should be specified only if the MAC is
%an effective one, otherwise there is no need to give the function this
%parameter.

%outputs are:
%R - a vector that holds the K optimal computation rates R=[R_comp1 ... R_comp_K]
%A - a matrix that holds the optimal coefficient vectors A=[a_1^T ... a_K^T]^T

if ~exist('B')
    B=eye(length(h));
end
    Hmat=inv(snr^(-1)*inv(B)+h*h');
    [U S V]=svd(Hmat);
    H_lat=sqrt(S)*V';
    [Tmp T]=fullLLL(H_lat,1);
    prices=diag(T'*Hmat*T);
    R=0.5*log2(snr./prices);
    idx=find(R<0);
    R(idx)=0;
    A=T';
end


function [Hreduced T]=fullLLL(H,delta)

[r c]=size(H);
[Qs Rs]=qr(H);

Q=Qs;
R=Rs;
T=eye(r);

ctr=1;
k=2;
while k<=c
    for l=[k-1:-1:1]
        
        u=round(R(l,k)/R(l,l));
       
        if u~=0
            R(1:l,k)=R(1:l,k)-u*R(1:l,l);
            T(:,k)=T(:,k)-u*T(:,l);
        end
    end
    
    if (delta*(R(k-1,k-1)^2)>((R(k,k)^2)+(R(k-1,k)^2)))
        tmpR=R(:,k-1);        
        R(:,k-1)=R(:,k);
        R(:,k)=tmpR;
        
        tmpT=T(:,k-1);
        T(:,k-1)=T(:,k);
        T(:,k)=tmpT;
        
        alpha=R(k-1,k-1)/sqrt(R(k-1:k,k-1)'*R(k-1:k,k-1));
        beta=R(k,k-1)/sqrt(R(k-1:k,k-1)'*R(k-1:k,k-1));
        theta=[alpha beta;-beta alpha];
        
        R(k-1:k,k-1:c)=theta*R(k-1:k,k-1:c);
        Q(:,k-1:k)=Q(:,k-1:k)*theta';
        
        k=max(k-1,2);
    else
        k=k+1;
    end
    ctr=ctr+1;
    if ctr>1000*c
        break
    end
end
Hreduced=H*T(1:c,1:c);
end
