Tuesday, July 7, 2009

Farey - Huffman Compress

The algorithm starts with computing a vector with the number of appearance of each character. At each step are combined two nodes to form an binary tree structure: the one with lowest appearance values. At the first step are combined the two characters with lowest appearance. These two node forms a new node with an appearance equal with the sum of the children. Than the vector of appearance is sorted and the process of combining two nodes starts again (it is possible to combine a simple character with a “sum” node or two “sum” nodes).


http://en.wikipedia.org/wiki/Huffman_coding
This want to be a method to optimize the Huffman tree compression and all the other algorithms based on this. For this is used the Farey addition.

The main idea use a correspondence between the Huffman tree and the binary search; for this the structure for the Huffman tree is viewed like a same tree but who stores in the roots the median values of an array in which is supposed to make an virtual search. All this values are after this, researched with a new method; instead calculating the half value of an interval, is computed the median value of that interval. This is known also like Farey addition; let’s suppose that the left and right limits of the interval "a", and "c" are viewed like a/1 and c/1. At each step starting from the fractions a/b and c/d is calculated (a+c)/(b+d). Similar to the Huffman codes - if the search is made in the left subinterval we memo an "0" and for the right an "1"- here if we want to merge the two methods and compare with the median value of binary search and with the Farey median value, we have six cases. But we can memo also with 0 and 1, using for this also codes with 0 and 1 and an array in which we store something similar to an impress; in this vector we memo 0 if the code starting for a position is not optimized and a 1 in other case.

This search produce smallest codes in approx. 30% of total cases in comparison with binary search.
Why we can use this method instead of using an other division of interval? The answers can be:- if we used an logx(n) the complexity of divide et impera algorithm grows; the steps are smaller but the length of the codes increase if we want to represent using 0 and 1 or if want to use an impress array the impress will grow because that can be cases in which an same binary code to represent multiple searched values.- the Farey bisection of interval produce an proportion -related to the length of intervals- equal to gold section number phi if we calculate at limit. - the number of "questions" to find a fraction represent log2(1+sum(fi(n)) approx. log2(3*n^2/pi^2) < log2(n).- on subintervals the envelope of the distribution looks like the Gauss bell. The generation of the Farey numbers can use some matrixes which are sub matrix of the Hadamard matrix/transform which is used for the generating of the Gauss distribution due to the central limit theorem.

I observed that if the "impress" is compressed again with UHBC algorithm (I found it on the internet) than the result is smaller, lets note this by "UhbcA".
The size of Farey codes + "UhbcA" << The size of Huffman codes.

PS. The impress can be memo using something like sparse matrix memorisation.

Saturday, July 4, 2009

Farey Fourier Window

The classic method use different method of windowing in wich is better to compute more armonics of the signal for obtaining an accurate aproximation:
http://en.wikipedia.org/wiki/Window_function
This technique can be used in the algorithms for: detection, selection, rejection of targets/clutter.
We can see from the programs from below that the leakage of the spectrum density is better for the Farey method vs. rectangular method. The leakage of the spectrum density shows how much each armonics from the Fourier series contribute to the energy of spectrum. The windowing methods want to eliminate the big errors of calculus that can appear at the limits of the sample analyzed. We can see also in the graphics of the spectum density of energy that the values near these boundary are smaller – their armonics do not contribute so much in the final result of the Fourier approximation; only the values near the middle of the graphics on the OX axis.

A program for generating in a file the Farey numbers.

//
//Optimizat cu Programare Dinamica
//
// N maxim = 328 => 32751 termeni
//

# include ;
# include ;
# include ;
# include ;

int N;
int v_N[3280], v_D[3280];
int v_Y[3280], v_X[3280];

FILE * pf;
int x,y,k,NT,i,j;

//
int Y(int n);
int X(int n)
{
if(n==0) {v_X[0]=0;return 0;}
else if(n==1) {v_X[1]=1;return 1;}
else
{
if(v_X[n]!=0) return v_X[n];
else
{
v_X[n]=floor( (v_Y[n-2]+N)/v_Y[n-1] )*v_X[n-1] - v_X[n-2];
return ( floor( (Y(n-2)+N)/Y(n-1) )*X(n-1) - X(n-2) );
}
}
}
int Y(int n)
{
if(n==0) {v_Y[0]=1;return 1;}
else if(n==1) {v_Y[1]=N;return N;}
else
{
if(v_Y[n]!=0) return v_Y[n];
else
{
v_Y[n]=floor( (v_Y[n-2]+N)/v_Y[n-1] )*v_Y[n-1] - v_Y[n-2];
return ( floor( (Y(n-2)+N)/Y(n-1) )*Y(n-1) - Y(n-2) );
}
}
}
//

void main()
{

//clrscr();
cout<<"Intr n= ";
cin>>N;
pf=fopen("FAREY.txt","w+");

k=0;
v_N[k]=x=X(k);
v_D[k]=y=Y(k);

while( (x!=y) && ((float)x/y<=0.5) )
{
cout< fprintf(pf,"%d/%d, ",x,y);
k++;
v_N[k]=x=X(k);
v_D[k]=y=Y(k);
}
NT=k;

//Calculez ceilalti termeni

if(N!=1)
k--;

i=1;
for(j=k-1;j>0;j--)
{
v_N[k+i]=v_D[j]-v_N[j];
v_D[k+i]=v_D[j];

cout< fprintf(pf,"1-%d/%d, ",v_N[k+i],v_D[k+i]);

i++;
NT++;
}

//

cout<<"1"<<"/"<<"1";
fprintf(pf,"%d/%d ",1,1);

v_N[NT]=1;
v_D[NT]=1;
NT++;

cout<<"\nNr de termeni = "<
fclose(pf);

getche();
}


The numbers that are bigger than the ½ are substracted from 1 and the final result is multiplied with 2 to respect the propriety of a window function: it has to have an amplitude equal with 1 and at the limits of interval it must be near 0.


A program in MATLAB with smaller modification from the programs of wikipedia.
The vector “w” can be generated fastest inside the program; it contains the values obtained from Farey series to generate an window function.

N=33;
k=0:N-1;

dr = 60;
w = [0/1, 1/10, 1/9, 1/8, 1/7, 1/6, 1/5, 2/9, 1/4, 2/7, 3/10, 1/3, 3/8, 2/5, 3/7, 4/9, 1/2, 1-5/9, 1-4/7, 1-3/5, 1-5/8, 1-2/3, 1-7/10, 1-5/7, 1-3/4, 1-7/9, 1-4/5, 1-5/6, 1-6/7, 1-7/8, 1-8/9, 1-9/10, 0 ] * 2;
%w = ones(1,N)

B = N*sum(w.^2)/sum(w)^2; % noise bandwidth (bins)

H = abs(fft([w zeros(1,7*N)]));
H = fftshift(H);
H = H/max(H);
H = 20*log10(H);
H = max(0,dr+H);

figure
area(k,w,'FaceColor', [0 .4 .6])
xlim([0 N-1])
ylim([0 2])
set(gca,'XTick', [0 : 1/8 : 1]*(N-1))
set(gca,'XTickLabel','0| | | | | | | |N-1')
grid on
ylabel('amplitude')
xlabel('samples')
title('Window function (Farey)')

figure
stem(([1:(8*N)]-1-4*N)/8,H,'-');
set(findobj('Type','line'),'Marker','none','Color',[.871 .49 0])
xlim([-4*N 4*N]/8)
ylim([0 dr])
set(gca,'YTickLabel','-60|-50|-40|-30|-20|-10|0')
grid on
ylabel('decibels')
xlabel('DFT bins')
title('Spectral "leakage" from a Farey windowed sinusoid')




The advantages of using the Farey windowing are that only the first armonics must be computed and this is obtained for rectangular window or others only for big values of samples (N big).
This window function not use sin of other function that use more steps in the computer to be aproximated.
The Farey nunbers can be generated faster or can be use an Farey extension/generalisation fore each interval and computed than.
Only for big values of samples (N == 1500 here) the first armonics aproximate in more steps for Fourier calculation the signal.
In this way the clutter or the detection of the targets can be done faster.

Wednesday, July 1, 2009

Farey Search vs Fibonacci Search

The Farey search algorithm is a method of searching a sorted array using the Farey addition. Compared to Fibonacci search, this algorithm is better regarding the number of steps.
For an example, for the array [1-23] there are 6 cases when Farey number of steps are lesser than Binary search; in Fibonacci search there are 4 cases when this search produce number of steps lesser than Binary search.

Here is a code that was used for Farey search vs. Binary search:

# include ;
# include ;
# include ;
# include ;
# include ;
# include ;
//# include ;

# define MAX 100000

long int a,a_aux,b_aux,b,mij,mij_v,v_N[3],v_D[3];
//[0] pt a, [1] pt mij, [2] pt b
long int a2,b2,mij2;
int opt;
long int r,scb,scbf;
unsigned long int c=0;
long int i,val;
long int huge v[MAX];
time_t t1,t2;
struct time t3,t4;

long int ncb,ncbf,o;



void main()
{
clrscr();
opt=0;

while(opt!=3)
{
clrscr();
cout<<"Please choose an option :\n";
cout<<"\t1. Binary search\n";
cout<<"\t2. Binary search using the FAREY addition\n";
cout<<"\t3. Enter the interval\n";
cout<<"\t4. END of the PROGRAM\n";
cin>>opt;
if(opt!=3) cout<<"\a";
}

while(opt<=3)
{
c=1;
if(opt==3)
{
cout<<"Enter a,b :\n";
cin>>a>>b;
a_aux=a;
b_aux=b;
for(i=a;i<=b+1;i++)
v[i-a]=i;
cout<<"The searched value will be each value of the interval... ";
//cin>>val;
getche();
}

if(opt<=2)
{
r=0;
scb=0;
scbf=0;
for(i=a_aux+1;i<=b_aux;i++)
{
ncb=1;
ncbf=1;

val=i;
//if(val==176) exit(1);

for(o=1;o<=2;o++)
{
opt=o;


a=a_aux;
b=b_aux;
mij2=(a+b)/2;


//initializari
if(opt==1) mij=(a+b)/2;
else if(opt==2)
{
v_N[0]=a;
v_D[0]=1;
v_N[2]=b;
v_D[2]=1;
v_N[1]=v_N[0]+v_N[2];
v_D[1]=v_D[0]+v_D[2];
mij=v_N[1]/v_D[1];
}

c=0;

a2=a;
b2=b;
mij2=(a2+b2)/2;

t1=time(NULL);
gettime(&t3);


while(
(val!=v[mij2])
&&(val!=v[mij])
&&(val!=v[a])
&&(val!=v[b])
)
{
if(opt==1)
{
if(val else if(val>v[mij2]) a=mij2;
mij2=(a+b)/2;
}
else if(opt==2)
{

if(val {
if(val>v[mij])
{
a=mij;
b=mij2;
v_N[0]=v_N[1];
v_D[0]=v_D[1];
v_N[2]=(a+b);
v_D[2]=2;
}
else //val {
if(mij {
b=mij;
v_N[2]=v_N[1];
v_D[2]=v_D[1];
}
else
{
b=mij2;
v_N[2]=a2+b2;
v_D[2]=2;
}
}
b2=mij2;
}
else if(val>v[mij2])
{
if(val {
a=mij2;
b=mij;
v_N[0]=a+b;
v_D[0]=2;
v_N[2]=v_N[1];
v_D[2]=v_D[1];
}
else//val>v[mij] and val>v[mij2], mij ? mij2
{
if(mij>mij2)
{
a=mij;
v_N[0]=v_N[1];
v_D[0]=v_D[1];
}
else
{
a=mij2;
v_N[0]=a2+b2;
v_D[0]=2;
}
}
a2=mij2;
}
v_N[1]=v_N[0]+v_N[2];
v_D[1]=v_D[0]+v_D[2];
mij=v_N[1]/v_D[1];
mij2=(a2+b2)/2;
}
c++;

}



if(o==1) {ncb=c;}
else {ncbf=c;}

t2=time(NULL);
gettime(&t4);

}
if(ncbelse cout<

if(ncb>ncbf)
{
r++;
scb+=ncb;
scbf+=ncbf;
}

}
float pc=0.0,pro=0.0;
pro=((float)scb)/scbf;
pc=(r*100.0)/(b_aux-a_aux+1);
cout<<"\nCBF < CB for "<printf("\nIn %g%% from cases, in total CBF < CB for %g times",pc,pro);
getche();
}

clrscr();
cout<<"Please choose an option :\n";
cout<<"\t1. Binary search\n";
cout<<"\t2. Binary search using the FAREY addition\n";
cout<<"\t3. Enter the interval\n";
cout<<"\t4. END of the PROGRAM\n";
cin>>opt;

}

}

Fibonacci search program

///////////////////////////////////////////////////////////////////////////////////////
//////
////// Fibonaccian search in an ordered array
////// This program can be downloaded from http://www.ics.forth.gr/~lourakis/fibsrch/
////// Copyright (C) 2005 Manolis Lourakis (lourakis AT .DeleteThis.ics.forth.gr)
////// Institute of Computer Science, Foundation for Research & Technology - Hellas
////// Heraklion, Crete, Greece.
//////
////// This program is free software; you can redistribute it and/or modify
////// it under the terms of the GNU General Public License as published by
////// the Free Software Foundation; either version 2 of the License, or
////// (at your option) any later version.
//////
////// This program is distributed in the hope that it will be useful,
////// but WITHOUT ANY WARRANTY; without even the implied warranty of
////// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
////// GNU General Public License for more details.
//////
///////////////////////////////////////////////////////////////////////////////////////

#include
#include
#include
#include

//#include "fibsrch.h"

int s1 = 0;
int s2 = 0;

int data[]={1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23};
int i, x, n, val;


/*
* If val is found in arr, return the index of its location in arr.
* Otherwise, return the index of the smallest element greater than val
*/
static int binsrch_geq(const int *arr, int n, int val)
{
register int low, high, mid;
int geq;

low=0; high=n-1; geq=-1;

/* binary search for finding the element with value val */
while(low<=high){

s2++;

mid=(low+high)>>1; //(low+high)/2;
if(val high=mid-1;
geq=mid;
}
else if(val>arr[mid])
low=mid+1;
else
return mid; /* found */
}

return geq; /* not found */
}

/*
Fibonaccian search for locating the index of "val" in an array "arr" of size "n"
that is sorted in ascending order. See http://doi.acm.org/10.1145/367487.367496

Algorithm description
-----------------------------------------------------------------------------
Let Fk represent the k-th Fibonacci number where Fk+2=Fk+1 + Fk for k>=0 and
F0 = 0, F1 = 1. To test whether an item is in a list of n = Fm ordered numbers,
proceed as follows:
a) Set k = m.
b) If k = 0, finish - no match.
c) Test item against entry in position Fk-1.
d) If match, finish.
e) If item is less than entry Fk-1, discard entries from positions Fk-1 + 1 to n.
Set k = k - 1 and go to b).
f) If item is greater than entry Fk-1, discard entries from positions 1 to Fk-1.
Renumber remaining entries from 1 to Fk-2, set k = k - 2 and go to b)

If Fm>n then the original array is augmented with Fm-n numbers larger
than val and the above algorithm is applied.
*/

int fibsrch(const int *arr, int n, int val)
{
register int k, idx, offs;
static int prevn=-1, prevk=-1;

/* Precomputed Fibonacci numbers F0 up to F46. This implementation assumes that the size n
* of the input array fits in 4 bytes. Note that F46=1836311903 is the largest Fibonacci
* number that is less or equal to the 4-byte INT_MAX (=2147483647). The next Fibonacci
* number, i.e. F47, is 2971215073 and is larger than INT_MAX, implying that it does not
* fit in a 4 byte integer. Note also that the last array element is INT_MAX rather than
* F47. This ensures correct operation for n>F46.
*/
const static int Fib[47+1]={0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181, 6765,
10946, 17711, 28657, 46368, 75025, 121393, 196418, 317811, 514229, 832040, 1346269, 2178309, 3524578,
5702887, 9227465, 14930352, 24157817, 39088169, 63245986, 102334155, 165580141, 267914296,
433494437, 701408733, 1134903170, 1836311903, INT_MAX};

/* find the smallest fibonacci number that is greater or equal to n. Store this
* number to avoid recomputing it in the case of repetitive searches with identical n.
*/
if(n!=prevn){
#if 1
k=(n>1)? binsrch_geq(Fib, sizeof(Fib)/sizeof(int), n) : 1;
#else /* the above binary search call can be substituted by the following code fragment: */
{
register int f0, f1, t;

for(f0=0, f1=1, k=1; f1 ;
}
#endif
prevk=k;
prevn=n;
}
else k=prevk;

/* If the sought value is larger than the largest Fibonacci number less than n,
* care must be taken top ensure that we do not attempt to read beyond the end
* of the array. If we do need to do this, we pretend that the array is padded
* with elements larger than the sought value.
*/
for(offs=0; k>0; ){

//s1++;

idx=offs+Fib[--k];

/* note that at this point k has already been decremented once */
if(idx>=n || val continue;
else if (val>arr[idx]){ // val in 2nd part
s1++;
offs=idx;
--k;
}
else // val==arr[idx], found
return idx;
}

return -1; // not found
}

#if 1
/* Sample driver program for fibsrch() */

void main(void)
{

for(val=1;val<=n;val++)
{
x=val; n=sizeof(data)/sizeof(int);

i=fibsrch(data, n, x);
if(i>=0)
printf("%d found at index %d\n", x, i);
else
printf("%d was not found\n", x);
printf("\nFIB search steps = %d",s1);

/****************************************/

i=binsrch_geq(data, n, x);
if(i>=0)
printf("%d found at index %d\n", x, i);
else
printf("%d was not found\n", x);
printf("\nBIN search steps = %d",s2);

printf("\n--------------------------------\n");
s1 = 0;
s2 = 0;
}

getche();
}
#endif



References

http://en.wikipedia.org/wiki/Fibonacci_search_technique
http://www.ics.forth.gr/~lourakis/fibsrch/#Algo

FAREY Search vs Binary Search

Some values can be found faster usinf the Farey addition than using the classic algorithm.

# include ;
# include ;
# include ;
# include ;
# include ;
# include ;

# define MAX 10000

int a,a_aux,b_aux,b,mij,mij_v,v_N[3],v_D[3];
//[0] pt a, [1] pt mij, [2] pt b
int a2,b2,mij2;
int opt;
unsigned long int c=0;
int i,val,v[MAX];
struct time t3,t4;


void main()
{
clrscr();
opt=0;

while(opt!=3)
{
clrscr();
cout<<"Please choose an option :\n";
cout<<"\t1. Binary search\n";
cout<<"\t2. Binary search using the FAREY addition\n";
cout<<"\t3. Enter the interval\n";
cout<<"\t4. END THE PROGRAM\n";
cin>>opt;
if(opt!=3) cout<<"\a";
}

while(opt<=3)
{
c=1;
if(opt==3)
{
cout<<"Please enter a,b :\n";
cin>>a>>b;
a_aux=a;
b_aux=b;
for(i=a;i<=b+1;i++)
v[i-a]=i;
cout<<"Please enter the shearched value: ";
cin>>val;
}

if(opt<=2)
{
a=a_aux;
b=b_aux;
mij=(a+b)/2;

//initializari
if(opt==1) mij=(a+b)/2;
else if(opt==2)
{
v_N[0]=a;
v_D[0]=1;
v_N[2]=b;
v_D[2]=1;
v_N[1]=v_N[0]+v_N[2];
v_D[1]=v_D[0]+v_D[2];
mij=v_N[1]/v_D[1];
}

c=0;
mij_v=0;

a2=a;
b2=b;
mij2=(a2+b2)/2;

time_t t1,t2;
t1=time(NULL);
gettime(&t3);

int q,r;

while(
v[mij]!=val
)
{

r=0;

if(v[a]==val) mij=a;
else if(v[b]==val) mij=b;
else if(v[mij2]==val) mij=mij2;
else
{

q=0;
if(opt==2)
{
if((mij {
a2=a=mij2;
v_N[0]=v_N[1];
v_D[0]=v_D[1];
b2=b=mij;
v_N[2]=mij;
v_D[2]=1;
q=1;
}
else
if((mij2 {
b2=b=mij2;
v_N[2]=v_N[1];
v_D[2]=v_D[1];
a2=a=mij;
v_N[0]=mij;
v_D[0]=1;
q=1;
}
}
if(q==0)
{
/*
if(opt==2)
{
a=a2;b=b2;
if(mij2 }
*/

if(v[(a2+b2)/2]==val)
{
mij=(a2+b2)/2;
v_N[0]=a2;
v_D[0]=1;
v_N[2]=b2;
v_D[2]=1;
v_N[1]=a2+b2;
v_D[1]=2;
}
else
{
if(opt==2)
if(v[mij]>val)
if(mij2 if((v[a]<=val)&&(val<=v[mij]))
{
if(opt==1)
b=mij;
else
{
if(mij2 b=b2=mij;
else
b=b2=mij2;
if(v[mij2]==val)
{
mij=mij2;
}
else
{
//Test if mij_FAREY < mij_Bisectie
if(v[mij]<=v[(a+b)/2])
{
//F(a) * F(mij) < 0
v_N[2]=v_N[1];
v_D[2]=v_D[1];
b=v_N[2]/v_D[2];
}
// mij_FAREY >= mij_Bisectie
else
{
if((v[mij]<=val)&&(val<=v[(a+b)/2]))
{
v_N[0]=a+b;
v_D[0]=2;
v_N[2]=v_N[1];
v_D[2]=v_D[1];
a2=a=v_N[0]/v_D[0];
b2=b=v_N[2]/v_D[2];
}
else if((v[a]<=val)&&(val<=v[(a+b)/2]))
{
v_N[2]=a+b;
v_D[2]=2;
b2=b=v_N[2]/v_D[2];
}
}
}
}
}

else if((v[mij]<=val)&&(val<=v[b]))
{
if(opt==1)
a=mij;
else
{
if(mij2 a=a2=mij2;
else
a=a2=mij;
if(v[mij2]==val)
mij=mij2;
else
{
//Testez daca mij_FAREY < mij_Bisectie
if(mij<=(a+b)/2)
{
if((v[mij]<=val)&&(val<=v[(a+b)/2]))
{
// F(mij) * F((a+b)/2) < 0
v_N[0]=v_N[1];
v_D[0]=v_D[1];
v_N[2]=a+b;
v_D[2]=2;
a2=a=v_N[0]/v_D[0];
b2=b=v_N[2]/v_D[2];
}
else // F(mij) * F(b) < 0
{
v_N[0]=a+b;
v_D[0]=2;
a2=a=v_N[0]/v_D[0];
}
}
// mij_FAREY >= mij_Bisectie
else
{
// F(mij) * F(b) < 0
v_N[0]=v_N[1];
v_D[0]=v_D[1];
a2=a=v_N[0]/v_D[0];
}
}
}
}
}
}

if(opt==1) {mij_v=mij; mij=(a+b)/2;}
else
{
mij2=(a2+b2)/2;
if(mij2 else
if(v[(a2+b2)/2]==val)
{
v_N[0]=a2;
v_D[0]=1;
v_N[2]=b2;
v_D[2]=1;
v_N[1]=a2+b2;
v_D[1]=2;
mij=(a2+b2)/2;
}
else
{
if(r==0)
{
mij_v=mij;
v_N[1]=v_N[0]+v_N[2];
v_D[1]=v_D[0]+v_D[2];
mij=(v_N[0]+v_N[2])/(v_D[0]+v_D[2]);
}
}
}
c++;

}
}
cout<<"\nNo. of iteration of the program: "< cout<
t2=time(NULL);
gettime(&t4);

cout<<"\n\n";

printf("Time of execution in second : %f\n",difftime(t2,t1));

printf("System clock before execution is: %2d:%02d:%02d.%02d\n",
t3.ti_hour, t3.ti_min, t3.ti_sec, t3.ti_hund);
printf("System clock after execution is: %2d:%02d:%02d.%02d\n",
t4.ti_hour, t4.ti_min, t4.ti_sec, t4.ti_hund);


getche();
}

clrscr();
cout<<"Please choose an option :\n";
cout<<"\t1. Binary search\n";
cout<<"\t2. Binary search using the FAREY addition\n";
cout<<"\t3. Enter the interval\n";
cout<<"\t4. END THE PROGRAM\n";
cin>>opt;

}

}


Sparse matrix using the Stern-Brocot Tree

Sparse/Rare Matrix representation using the Stern-Brocot Tree

1/2
1/3 2/3
1/4 3/4 2/5 3/5
1/5 4/5 3/7 4/7 2/7 5/7 3/8 5/8
1/6 5/6 4/9 5/9 3/10 7/10 4/11 7/11 2/9 7/9 5/12 7/12 3/11 8/11 5/13 8/13
1/7 6/7 5/11 6/11 4/13 9/13 5/14 9/14 3/13 10/13 7/17 10/17 4/15 11/15 7/18 11/18 2/11 9/11 7/16 9/16 5/17 12/17 7/19 ...
...
All the fraction in each row of the Stern-Brocot tree can be viewed like coordinates in a matrix(row/column).
For the generation of the Stern-Brocot tree can be used a method that produce the Farey fraction but not in ascending order(see[1]). If we visit the Stern-Brocot tree by rows we obtain at the end of each row a complete set of the Farey sequence. In this way, at the end of the n'th row we obtained {Fn}. Some of the fractions with the denomitor equal to n are stored above in the Stern-Brocot tree.
This is the ideea of the representation of the sparse/rare matrix using the Stern-Brocot tree.
Starting from the fraction 1/1 we can associate to each fraction a binary code; left = = 0, right = = 1.
If a sparse/rare matrtix contain a non-zero value on the position (x,y), we can associate for this position two decimal numbers obtained by the binary code described above - spliting the code by half. Let's note the code by z; than in some cases, z will be smaller than the binary representation of x and of y .
A better representation of a sparse/rare matrix than the classic method - storing a compact memory space - use dynamic memory - circular dynamic lists. In these lists are memo some values related to the next row and column in which is stored a non-zero value.

Let's take an example, the matix:

0 1 0 6 0
0 0 0 1 0
4 0 9 0 0
0 0 0 0 0
0 13 0 2 0
0 0 0 0 1
The positions (1, 4) - "0 0"and (5, 2) "11 01" (starting from 1/1) can be optimized like described above to (0, 0), (3, 1).
When allocating the nodes in the circular dynamic lists, these two nodes can use an other struct definition for the node: for the type of the row and column to use less bits - memo only one byte: the decimal representation of "1"+z .
(*) The codes obtained from the Stern-Brocot tree, are smaller than x and y only for the values of the interior of the tree - in general the values that do not appear in {Fn}.

References:
Pages maintained by Mr. Alexander Bogomolny :
1) http://www.cut-the-knot.org/blue/b-tree.shtml
2) http://www.cut-the-knot.org/blue/Stern.shtml
2) http://www.cut-the-knot.org/blue/encoding.shtml
3) http://www.cut-the-knot.org/blue/chaos_game.shtml




/************************************************************************************
Program that generates all the fractions from a level "h", of the Stern-Brocot Tree @ For a detailed description of the algorithm, please visit the site: @ http://www.cut-the-knot.org/blue/b-tree.shtml ************************************************************************************/
/************************************************************************************INCLUDES************************************************************************************/# include "stdio.h"# include "malloc.h"# include "math.h"
/************************************************************************************TYPES & VARIABLES************************************************************************************/
typedef unsigned long int UL_int;typedef unsigned char ui8_t;
UL_int h;
struct nod{ UL_int num,den; struct nod * s, * d;};
typedef struct nod NOD;
NOD * rad;
/* Output file wich stores the optimized positions */FILE * pf;
/* Variables that store the number of nodes */UL_int N, c_init;/* Variable that stores the number of optimized positions */UL_int c_opt;
/************************************************************************************FUNCTIONS************************************************************************************/
/*----------------------------------------------------*//* Functions for comparing the binary representations *//*----------------------------------------------------*/
ui8_t BIN_comp(UL_int a, UL_int b, UL_int h){ UL_int nr_cif_a, nr_cif_b; nr_cif_a = (UL_int) (log(a)/log(2.0)) + 1; nr_cif_b = (UL_int) (log(b)/log(2.0)) + 1; if(ceil((nr_cif_a + nr_cif_b)/8.0)+ 1 > ceil((h-2)/8.0) + 1) return 1; else return 0;}
/*----------------------------------------------------------------------*//* Function for the generation of the Stern-Brocot tree with "h" levels *//*----------------------------------------------------------------------*/
NOD * SB_tree(NOD * r, UL_int x, UL_int y, UL_int c){ NOD * p, * q; if(c==h-1) { /* Print all the fractions from a read level */ printf("%lu/%lu ",r->num,r->den); /* Check for optimized position */ if(BIN_comp(r->num,r->den,h) == 1) { c_opt++; fprintf(pf,"%lu/%lu\n",r->num,r->den); } return r; } else { p =(NOD *)malloc(sizeof(NOD)); p->s = NULL; p->d = NULL; p->num = x; p->den = x+y; r->s = SB_tree(p,x,x+y,c+1); q =(NOD *)malloc(sizeof(NOD)); q->s = NULL; q->d = NULL; q->num = y; q->den = x+y; r->d = SB_tree(q,y,x+y,c+1); }}
/************************************************************************************MAIN ()************************************************************************************/
void main(void){ pf=fopen("out.txt","w+"); /* Read the level */ printf("Enter the level in the tree: "); scanf("%lu",&h); /* The root is treated separatly */ rad =(NOD *)malloc(sizeof(NOD)); rad->s = NULL; rad->d = NULL; rad->num = 1; rad->den = 2; /* Generation of the STERN-BROCOT tree fractions */ /* http://www.cut-the-knot.org/blue/b-tree.shtml */ c_init = 1; SB_tree(rad,1,2,c_init); /* Print the number of optimised positions */ N = (UL_int) pow(2,h); printf("\nThe number of optimised positions is: %lu / %lu", c_opt, N); printf("\nPercent = %.2f %% optimised",(c_opt*100/(float)N));
fclose(pf);}/************************************************************************************END************************************************************************************/


/********************************************************************************************************************************************************************************/
/************************************************************************************INCLUDES************************************************************************************/
# include "stdio.h"
# include "malloc.h"
# include "math.h"
# include "conio.h"
/************************************************************************************TYPES & VARIABLES************************************************************************************/
typedef unsigned __int64 UL_int;
typedef unsigned char ui8_t;
UL_int h;
struct nod{ UL_int num,den; struct nod * s, * d;};
typedef struct nod NOD;
NOD * rad;

/* Output file wich stores the optimized positions */
FILE * pf;
/* Variables that store the number of nodes */
UL_int N, c_init;
/* Variable that stores the number of optimized positions */
UL_int c_opt;
/************************************************************************************FUNCTIONS************************************************************************************/
/*----------------------------------------------------*//* Functions for comparing the binary representations *//*----------------------------------------------------*/
ui8_t BIN_comp(UL_int a, UL_int b, UL_int h)
{
UL_int nr_cif_a, nr_cif_b;

nr_cif_a = (UL_int) (log(a)/log(2.0)) + 1;
nr_cif_b = (UL_int) (log(b)/log(2.0)) + 1;

if(ceil((nr_cif_a + nr_cif_b)/8.0)+ 1 > ceil((h-2)/8.0) + 1)
    return 1;
else
    return 0;
}
/*----------------------------------------------------------------------*//* Function for the generation of the Stern-Brocot tree with "h" levels *//*----------------------------------------------------------------------*/
NOD * SB_tree(NOD * r, UL_int x, UL_int y, UL_int c)
{
    NOD * p, * q;
    if(c==h-1)
    { /* Print all the fractions from a read level */
        //printf("%lu/%lu ",r->num,r->den);
        /* Check for optimized position */
        if(BIN_comp(r->num,r->den,h) == 1)
        {
            c_opt++;
            //fprintf(pf,"%lu/%lu\n",r->num,r->den);
        }
            return r;
    }
    else
    {
        p =(NOD *)malloc(sizeof(NOD));
        p->s = NULL;
        p->d = NULL;
        p->num = x;
        p->den = x+y;
        r->s = SB_tree(p,x,x+y,c+1);
        q =(NOD *)malloc(sizeof(NOD));
        q->s = NULL;
        q->d = NULL;
        q->num = y;
        q->den = x+y;
        r->d = SB_tree(q,y,x+y,c+1);
    }
}
/************************************************************************************MAIN ()************************************************************************************/
int main(void)
{
    pf=fopen("out.txt","w+");
    /* Read the level */
    printf("Enter the level in the tree: ");
    scanf("%llu",&h);
    /* The root is treated separatly */
    rad =(NOD *)malloc(sizeof(NOD));
    rad->s = NULL;
    rad->d = NULL;
    rad->num = 1;
    rad->den = 2;
    /* Generation of the STERN-BROCOT tree fractions */ /* http://www.cut-the-knot.org/blue/b-tree.shtml */
    c_init = 1;

    SB_tree(rad,1,2,c_init);

    /* Print the number of optimised positions */
    N = (UL_int) pow(2,h);
     printf("\nThe number of optimised positions is: %llu / %llu", c_opt, N);
     printf("\nPercent = %.2f %% optimised",(c_opt*100/(float)N));


fclose(pf);

getche();

return 1;
}
/************************************************************************************END************************************************************************************/