RMQ问题的Fischer-Heun算法(直接方法)

在普通RMQ问题的<O(n),O(1)>算法中,由于需要构造Cartesian Tree和得到Euler tour,两个2*n-1大小的数组E和L使得空间消耗增加O(4*n)。本文介绍的Fischer-Heun算法绕过构建Cartesian Tree的步骤,也不需要将普通RMQ转化为RMQ+1/-1问题再求解。该算法的基本思想也是采用RMQ+1/-1的分组处理, 分别生成组内的Lookup Table和组外的稀疏表。但是在每个组的大小缩小了一倍,成为s=(logn)/4。由于每个组对应一颗笛卡尔树,若两个不同的组对应的两棵笛卡尔树相同的话,那么这两个组的lookup table完全一样;另外由于大小为n的笛卡尔树的个数为Catalan数Cn,因此构建的基本组大小是Cn。在组大小为s=(logn)/4时,根据Catalan数Cn的近似公式得到基本组个数为O(4^s/s^(3/2))数量级,虽然Catalan数Cn增长的很快,但在实践中Cs是一个比较小的数:比如,当s=8时,C8=1430,可以处理的A数组大小为n=2^32,n已经很大了。类似于RMQ+1/-1算法的分析,组内预处理的复杂度为线性:O(4^s/s^(3/2) * s^2) = O(n)!组外预处理的复杂度为O(s*logs)=O(n)。查询方式也是采用RMQ+1/-1算法中的查询做法,复杂度也是O(1)。

 Fischer-Heun算法的关键是如何确定每一个组的Cartesian Tree编号,即如何确定组的类型数组T。在J. Fischer和V. Heun 2006年发表的论文中,比较巧妙地使用ballot数来处理。ballot数对应的三角形中,设当前位置在(p,q),那么ballot数C[p,q]相当于从(p,q)到(0,0)走法总数,如果p=q,C[p,q]就是Catalan数Cp,即对角线上的点对应Catalan数。确定组的类型数组的算法利用构造Cartesian Tree的过程,但是不真正构造树,而且堆栈中存放的是A[i]值,而不是下标值。首先在stack底部放置一个sentinel元素Integer.MIN_VALUE。设置Cartesian Tree编号N=0。在构造笛卡尔树的算法中,从ballot数三角形中的点(s,s)出发,每执行一次循环,对应往左走一步;如果某个元素A[i]小于栈顶元素,那么弹出栈顶元素,这时对应往上走一步,若当前位置在(p,q),N的值需要增加ballot数C[p-1,q];最后将A[i]入栈。由于构造Cartesian Tree算法的特点,按照这样的走法,可以保证在ballot数对应的三角形中不会超出边界范围。最后s个元素处理完时,N的值就是block的类型值,范围为0..Cs-1。


实现:

/**
 * 
 * Fischer-Heun Algorithm 
 * A space-economic general RMQ algorithm without transforming into RMQ+1/-1 
 * time complexity: <O(n),O(1)>
 * 
 * see also: (2006) J. Fischer, V. Heun, Theoretical and Practical Improvements on the RMQ-Problem, with Applications to LCA and LCE
 * 
 * Copyright (c) 2011 ljs (http://blog.csdn.net/ljsspace/)
 * Licensed under GPL (http://www.opensource.org/licenses/gpl-license.php) 
 * 
 * 
 * @author ljs
 * 2011-08-07
 *
 */
public class RMQ_FH {
	 
	private int blocksize;
	private int blockcount;
	
	
	//block type
	private int[] T;
	//block LookUp table
	private int[][] P;
	
	//out-of-block sparse table 
	//the first-dimension indices are block IDs (from 0...blockcount-1)
	private int[][] M;
	
 
	 
	//Note: this method is not quite different from the one in class MinusOrPlusOne_RMQ:
	//Except that blockTypesCount and computeBlockType are different, this method is 
	//almost the same as the one in class MinusOrPlusOne_RMQ.
	public void preprocess(int[] A) throws Exception {
		int n = A.length;
		
		int paddingsize = 0;	
				
		//block size: (logn)/4
		int s = (int) (Math.log(n) / Math.log(2)) >> 2;
		if(s==0){
			s = (n>4)?4:n;	//small problem
		}	
		int count = (int) Math.ceil(n / (double)s);
		this.blocksize = s;
		this.blockcount = count;
		
		//compute the ballot numbers for determining each block's cartesian trees or block type
		int[] ballotnums = RMQ_FH.makeBallotNumbersTable(s);
		
		// padding the last block
		int endblocksize = n - s * (count - 1);
		if (endblocksize > 0 && endblocksize < s) {
			paddingsize = s - endblocksize;
		}
		
		//step 1: in-block preprocess
		//the size of LU table (one-dimensional)
		int size = s*(s+1)/2;
		
		int start = 0; // j is the index of A
		int end = -1;
		T = new int[count];
		int blockTypesCount = ballotnums[(s+2)*(s+1)/2 - 1]; //catalan number
		P = new int[blockTypesCount][size];
		int[] B = new int[count]; //used in ST algorithm: the min-value array for each block
		boolean[] pDone =new boolean[blockTypesCount];
		
		int fullblockscnt = count;
		if(paddingsize>0){
			fullblockscnt = count - 1;
		}	
		
		for (int i = 0; i < fullblockscnt; i++) {
			start = end+1;
			end = (i+1)*s-1;
			//compute the type of the block 
			int type = computeBlockType(ballotnums,A,start,end);	
			T[i] = type;
			if(!pDone[type]){
				//if LU table is not done yet for this type of block
				P[type] = makeLUTable(s,A,start,size);		
				pDone[type]=true;
			}			
			B[i] = queryLUTable(s,start,start,end,type);			
		}
		
		//the end block
		if(paddingsize>0){
			start = end+1;
			end = count*s-1;
			
			//extend the end block
			int actualsize = n-start;
			int[] D = new int[s];
			System.arraycopy(A, start, D, 0, actualsize); 
			for (int k = actualsize; k < s; k++) {
				D[k] = D[k - 1] + 1;
			}			
			
			int type = computeBlockType(ballotnums,D,0,s-1);	
			T[count-1] = type;
			if(!pDone[type]){
				P[type] = makeLUTable(s,D,0,size);		
				pDone[type]=true;
			}	
			//the min-index from start...n-1, not start...end
			B[count-1] = queryLUTable(s,start,start,n-1,type);		
		}
		
		//step 2: Sparse table algorithm applied to out-of-blocks		
		this.M = outBlockPreprocess(B,A,count);
	}
	
	//Note: this method is the same as the one in class MinusOrPlusOne_RMQ. It is copied from there.
	//return the index
	public int query(int[] A,int p,int q){
		if(q<p){
			//swap  
			int tmp=p;p=q;q=tmp;
		}
		
		int start = 0; 
		int end = -1;
		int s = 0,t=0;  
		int startMin=-1,endMin=-1; //the start block and end block's min index
		for (int i = 0; i < this.blockcount; i++) {
			start = end+1;
			end = (i+1)*this.blocksize-1;
			if(p>=start && q<=end){
				//within a block
				return queryLUTable(blocksize,start,p,q,T[i]);
			}else if(p>=start && p<=end){
				startMin = queryLUTable(blocksize,start,p,end,T[i]);
				s=i+1;
			}else if(q<=end && q>=start){
				endMin = queryLUTable(blocksize,start,start,q,T[i]);
				t=i-1;
				break;
			}
		}
		int minIndex = startMin;
		
		if(s<=t){
			int outBlocksMin = outBlockQuery(A,s,t);
			if(A[startMin]>A[outBlocksMin]){
				minIndex = outBlocksMin;
			}
		}
		if(A[minIndex]>A[endMin]){
			minIndex = endMin;
		}
		return minIndex;
	}
	
	//Note: this method is the same as the one in class MinusOrPlusOne_RMQ. It is copied from there.
	//ST: O(1) for querying
	//precondition: s<=t
	private int outBlockQuery(int[] A,int s,int t){
		int k = (int)(Math.log(t-s+1)/Math.log(2));
		//the first interval
		int mina = M[s][k];
		int minb = M[t-(1<<k)+1][k];
		if(A[mina]<=A[minb])
			return mina;
		else
			return minb;
	}
	
	//Note: this method is the same as the one in class MinusOrPlusOne_RMQ. It is copied from there.
	private int[][] outBlockPreprocess(int[] B,int[] A,int count){
		//floor value
		int maxJ=(int)(Math.log(count)/Math.log(2));
		 
		int[][] M = new int[count][maxJ+1];
		
		//initial condition for dynamic programming: the RMQ for interval length=1
		for (int i = 0; i < count; i++)
	          M[i][0] = B[i];
		
		//dynamic programming: compute values from smaller(j=1) to bigger intervals
	    for (int j = 1; j<=maxJ; j++){
	        for (int i = 0; i + (1 << j) - 1 < count; i++){
	        	int nexti = i + (1 << (j - 1));
	            if (A[M[i][j - 1]] <= A[M[nexti][j - 1]])
	                M[i][j] = M[i][j - 1];
	            else
	                M[i][j] = M[nexti][j - 1];
	        }
	    }
	    return M;
	}
	
	//Note: this method is the same as the one in class MinusOrPlusOne_RMQ. It is copied from there.
	//return the index  
	//precondition: i<=j
	private int queryLUTable(int blocksize,int offset,int i,int j,int type){
		i -= offset; j-=offset;
		int[] L = P[type];
		int index = blocksize*i - (i-1)*i/2 + (j-i);
		return L[index]+offset;
	}
	 
	//Note: this method is the same as the one in class MinusOrPlusOne_RMQ. It is copied from there.
	//use naive method to compute the lookup table for a block
	//the return index is relative to the block itself
	private int[] makeLUTable(int blocksize, int[] D,int offset,int size) {
 		int[][] Q = new int[blocksize][blocksize];
		for (int i = 0; i < blocksize; i++)
			Q[i][i] = i;
		for (int i = 0; i < blocksize; i++)
			for (int j = i + 1; j < blocksize; j++)
				if (D[Q[i][j - 1]+offset] <= D[j+offset])
					Q[i][j] = Q[i][j - 1];
				else
					Q[i][j] = j;
		//convert to one-dimension array		
		int[] L = new int[size];
		int k=0;
		for(int i=0;i<blocksize;i++){
			for(int j=i;j<blocksize;j++,k++){
				L[k] = Q[i][j];
			}
		}
		return L;
	}
	/*
	 Ballot numbers: 
	 see: Knuth's Art of Computer Programming Vol 4A section 7.2.1.6 "Generating all trees"
	 The diagonal numbers are Catalan numbers.
	 
	 Catalan numbers:(0..16) 
	 see: http://en.wikipedia.org/wiki/Catalan_number
	 
	 1, 1, 2, 5, 14, 42, 132, 429, 1430, 
	 4862, 16796, 58786, 208012, 742900, 
	 2674440, 9694845, 35357670 ....
	 */	
	public static int[] makeBallotNumbersTable(int s){
		int k = (s+2)*(s+1)/2; 
		int[] ballotnums = new int[k];
		ballotnums[0] = 1;
		int i = 1;
		for(int q=1;q<=s;q++){
			//p=0, only top element is used, left element is 0
			ballotnums[i] = ballotnums[i-q]; i++;		
			for(int p=1;p<q;p++){				
				ballotnums[i] = ballotnums[i-q] + ballotnums[i-1];
				i++;
			}
			//p=q, only left element is used, top element is 0
			ballotnums[i] = ballotnums[i-1]; i++;			
		}
		return ballotnums;
	}

	private static void printBallotNumbersTriangle(int[] ballotnums,int s){
		//print ballot numbers triangle
		System.out.println("Ballot numbers: ");
		int start = 0, end = -1;
		for(int m=0;m<=s;m++){
			start = end + 1;
			end = (m+2)*(m+1)/2 - 1; 
			for(int n=start;n<=end;n++)
				System.out.format(" %d", ballotnums[n]);
			System.out.println();
		}
	}
	
	private static void printCatalanNumbers(int[] ballotnums,int s){
		//print catalan numbers: 
		System.out.println("Catalan numbers: ");
		for(int m=0;m<=s;m++){
			int j = (m+2)*(m+1)/2 - 1; 
			if(m % 5 == 0){				
				System.out.format("(%d):",m);
			}
			System.out.format(" %d", ballotnums[j]);
			if((m+1) % 5 == 0)
				System.out.println();
		}
	} 
	public static void testBallotNumberTable(){
		int s = 16;		
		int[] ballotnums = RMQ_FH.makeBallotNumbersTable(s);
		RMQ_FH.printCatalanNumbers(ballotnums, s);
		System.out.println();
		RMQ_FH.printBallotNumbersTriangle(ballotnums, s);
	}
	
	//compute the block type, ie. the sequence number of its cartesian tree
	private int computeBlockType(int[] ballotnums, int[] D,int start,int end){
		int s = end - start + 1;
		//the stack is the rightmost path
		int[] stack = new int[s+1];
		stack[0] = Integer.MIN_VALUE; //the first element is a sentinel
		int top = 0; //stack top
		int N=0; //the sequence number of this block's cartesian tree
		int q = s;
		for(int i=0;i<s;i++){
			int p = s-i-1;
			while(stack[top]>D[i+start]){			
				//numbering: go up one level in the ballot number triangle				
				N += ballotnums[q*(q+1)/2+p];
				top--; //remove the element
				q--; //up one level
			}
			//when stack[top]<= D[i+start], we append the current number to the rightmost path
			stack[++top] = D[i+start];
		}
		return N;
	}
	public void testComputeBlockType(int[] D){		
		int s = D.length; //block size
		int[] ballotnums = RMQ_FH.makeBallotNumbersTable(s);
		
		int nr = this.computeBlockType(ballotnums,D, 0, D.length-1);
		System.out.format("block type(=cartesian tree sequence number): %d%n",nr);
	}
	
	private void reportLUTable(int[] A){
		for(int x=0;x<A.length;x++){
			System.out.format("%d..[%d-%d]",x,x,A.length-1);
			for(int y=x;y<A.length;y++){
				int p = query(A,x,y);				
				System.out.format(" %d/%d",A[p],p);
			}
			System.out.println();
		}	
		
	}
	
	public static void main(String[] args) throws Exception {
		//test ballot numbers
		RMQ_FH.testBallotNumberTable();
		
		//test block types
		System.out.println("*************");
		int[] A = new int[]{1,2,3,4,5};
		RMQ_FH rmq = new RMQ_FH();
		rmq.testComputeBlockType(A);
		
		System.out.println("*************");
		A = new int[]{5,4,3,2,1};
		rmq = new RMQ_FH();
		rmq.testComputeBlockType(A);
		
		//RMQ operations			
		System.out.format("***********************%n");
		A=new int[]{1,7,3};
		rmq = new RMQ_FH();
		rmq.preprocess(A);		
		rmq.reportLUTable(A);
		
		System.out.format("***********************%n");
		A=new int[]{2,4,3,1,6,7,8,9,1,7};
		rmq = new RMQ_FH();
		rmq.preprocess(A);		
		rmq.reportLUTable(A);
		
		System.out.format("***********************%n");
		A=new int[]{10,15,34,20,7,5,18,68,29,40, //0..9
				24,3,45,26,7,23,43,12,68,34,  //10..19
				26,34,33,12,80,57,24,42,77,27, //20..29
				56,33,23,32,54,13,79,65,19,33,  //30..39
				15,24,43,73,55,13,63,8,23,17};  //40..49
		rmq = new RMQ_FH();
		rmq.preprocess(A);
		
		
		int i=0;
		int j=49;		
		int min = rmq.query(A,i,j);
		System.out.format("RMQ for A[%d..%d]: A[%d]=%d%n", i,j,min,A[min]);
		
	}

}



测试输出:

Catalan numbers:
(0): 1 1 2 5 14
(5): 42 132 429 1430 4862
(10): 16796 58786 208012 742900 2674440
(15): 9694845 35357670
Ballot numbers:
 1
 1 1
 1 2 2
 1 3 5 5
 1 4 9 14 14
 1 5 14 28 42 42
 1 6 20 48 90 132 132
 1 7 27 75 165 297 429 429
 1 8 35 110 275 572 1001 1430 1430
 1 9 44 154 429 1001 2002 3432 4862 4862
 1 10 54 208 637 1638 3640 7072 11934 16796 16796
 1 11 65 273 910 2548 6188 13260 25194 41990 58786 58786
 1 12 77 350 1260 3808 9996 23256 48450 90440 149226 208012 208012
 1 13 90 440 1700 5508 15504 38760 87210 177650 326876 534888 742900 742900
 1 14 104 544 2244 7752 23256 62016 149226 326876 653752 1188640 1931540 2674440 2674440
 1 15 119 663 2907 10659 33915 95931 245157 572033 1225785 2414425 4345965 7020405 9694845 9694845
 1 16 135 798 3705 14364 48279 144210 389367 961400 2187185 4601610 8947575 15967980 25662825 35357670 35357670
*************
block type(=cartesian tree sequence number): 0
*************
block type(=cartesian tree sequence number): 41
***********************
0..[0-2] 1/0 1/0 1/0
1..[1-2] 7/1 3/2
2..[2-2] 3/2
***********************
0..[0-9] 2/0 2/0 2/0 1/3 1/3 1/3 1/3 1/3 1/3 1/3
1..[1-9] 4/1 3/2 1/3 1/3 1/3 1/3 1/3 1/3 1/3
2..[2-9] 3/2 1/3 1/3 1/3 1/3 1/3 1/3 1/3
3..[3-9] 1/3 1/3 1/3 1/3 1/3 1/3 1/3
4..[4-9] 6/4 6/4 6/4 6/4 1/8 1/8
5..[5-9] 7/5 7/5 7/5 1/8 1/8
6..[6-9] 8/6 8/6 1/8 1/8
7..[7-9] 9/7 1/8 1/8
8..[8-9] 1/8 1/8
9..[9-9] 7/9
***********************
RMQ for A[0..49]: A[11]=3



参考资料:

J. Fischer和V. Heun于2006年发表的论文“Theoretical and Practical Improvements on the RMQ-Problem, with Applications to LCA and LCE“



  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 1
    评论
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值