imos法.html

<title>いもす法</title>

<!--一新、改変OK,ここから-->
原本
<a href="https://imoz.jp/algorithms/imos_method.html">いもす法</a>

<pre>
想定
初期値0の縦R行、横C行の二次元配列Aに
以下のクエリをQ回行う
クエリ：矩形(R1,C1)-(R2,C2)に定数dを足す。
この計算をO(Q + R*C * R*C)で計算する。

元々二次元配列Aに、0以外の値が設定されている場合、
別のアルゴリズムを考えた方が良い。
(Aから勾配配置にO(R*C*R*C)で戻すなど工夫が必要)

--------------------

2次元配列
矩形領域に定数dを足す。
矩形領域は(R1,C1)-(R2,C2)で表すことにします。
R1<=R2 && C1<=C2
ナイーブな実装
<code>
for(int r=R1;r<=R2;r=r+1){
	for(int c=C1;c<=C2;c=c+1){
		A[r][c] = A[r][c] + d;
	}
}
</code>

いもす法は大きく分けて
二つのステップから実行される
１．勾配の配置
２．累積和

例1
縦R行、横C列の２次元配列Aにおいて
矩形(1,1)-(3,3)に1を加算する
<code>
pos1[0] = 1;//row
pos1[1] = 1;//col
pos2[0] = 3;//row
pos2[1] = 3;//col

int add = 1;//加算する値
for(int s=0;s<(1 << 2);++s){
	int point[2];//(r,c):勾配を設置する位置

	int count = 0;
	for(int i=0;i < 2;++i){
		if( ( s&(1 << i) ) == 0){
			point[i] = pos1[i];
		}else{
			++count;
			point[i] = pos2[i]+1;//要:配列レンジチェック
		}
	}

	if( count%2 == 0 ){//パリティチェック
		A[ point[0] ][ point[1] ] += add;
	}else{
		A[ point[0] ][ point[1] ] -= add;
	}
}
</code>
/*
やっていること
A[1][1]に1加算
A[3+1][1]に-1加算
A[1][3+1]に-1加算
A[3+1][3+1]に1加算
*/
初期値
 0  0  0  0  0
 0  0  0  0  0
 0  0  0  0  0
 0  0  0  0  0
 0  0  0  0  0

勾配の配置
 0  0  0  0  0
 0  1  0  0 -1
 0  0  0  0  0
 0  0  0  0  0
 0 -1  0  0  1
 
 //c+方向への累積和:A[r][c] += A[r][c-1]
<code>
for(int r=0;r < R;r=r+1){
	for(int c=1;c < C;c=c+1){
		A[r][c] += A[r][c-1];
	}
}
</code>
 0  0  0  0  0
 0  1  1  1  0
 0  0  0  0  0
 0  0  0  0  0
 0 -1 -1 -1  0
 
 //r+方向への累積和:A[r][c] += A[r-1][c]
 <code>
for(int c=0;c < C;c=c+1){
	for(int r=1;r < R;r=r+1){
		A[r][c] += A[r-1][c];
	}
}
</code>
 0  0  0  0  0
 0  1  1  1  0
 0  1  1  1  0
 0  1  1  1  0
 0  0  0  0  0

矩形(1,1)-(3,3)に1を加算出来た。

-------
例2　複数の矩形をまとめて処理できる。
(1,1)-(3,3)に1加算、(2,2)-(3,4)に2加算
勾配の配置,
 0  0  0  0  0
 0  1  0  0 -1
 0  0  2  0  0 (-2:-2は配置されていない)
 0  0  0  0  0
 0 -1 -2  0  1 (+2:+2は配置されていない)

c+方向の累積和
 0  0  0  0  0
 0  1  1  1  0
 0  0  2  2  2
 0  0  0  0  0
 0 -1 -3 -3 -2
 
r+方向の累積和
 0  0  0  0  0
 0  1  1  1  0
 0  1  3  3  2
 0  1  3  3  2
 0  0  0  0  0
 
計算できた

---------------------------
 
勾配配置を求める
例1でのいもす法の手順は
1.勾配の配置
2.c+方向の累積和
3.r+方向の累積和
で行った。
累積和を逆に行うことで、勾配の配置を求めることができる。

初期値
 0  0  0  0  0
 0  1  1  1  0
 0  1  1  1  0
 0  1  1  1  0
 0  0  0  0  0
 
r+方向の差分
<code>
for(int c=0;c < C;c=c+1){//cの順番は不問
	for(int r=R-1;r > 0;r=r-1){
		A[r][c] -= A[r-1][c];
	}
}
</code>
 0  0  0  0  0
 0  1  1  1  0
 0  0  0  0  0
 0  0  0  0  0
 0 -1 -1 -1  0

c+方向の差分
<code>
for(int r=0;r < R;r=r+1){//rの順番は不問
	for(int c=C-1;c > 0;c=c-1){
		A[r][c] -= A[r][c-1];
	}
}
</code>
 0  0  0  0  0
 0  1  0  0 -1
 0  0  0  0  0
 0  0  0  0  0
 0 -1  0  0  1
 
(1,1)-(3,3)の勾配配置が得られた。
未知の勾配配置に対しては
実験的考察が必要だと思います。
うまーく差分を求めていって、
少ない個数の勾配配置を得られれば成功です。
今のところ思いついている累積方向は
A[r][c] += A[r-1][c];
A[r][c] += A[r][c-1];
A[r][c] += A[r-1][c-1];
A[r][c] += A[r+1][c-1];
周囲8近傍/2の4パターンです。
きっと他にもある。

注意1
加算する矩形が(1,1)-(3,3)なのか
(1,1)-(2,2)なのかはA[r][c]の意味に依存します。
注意2
累積方向によって、配置する勾配は変わります。

--------------------------------
どこまでがいもす法なのか、これがわからない。

一次元配列の累積和を用いて
連続区間の和がO(1)で求められた。
同様に、
二次元配列の累積和を用いて
矩形内の和がO(1)で求められる。
＊（三次元の立方区間の和もO(1)で求められる。
(x1,y1,z1)-(x2,y2,z2):パリティ計算は同じ。)
＊(特殊な座標や特殊な連続区間に対しては、
同様の事が出来るかはわかりません。）

例2で求めた
 0  0  0  0  0
 0  1  1  1  0
 0  1  3  3  2
 0  1  3  3  2
 0  0  0  0  0
この二次元配列に対して
c+方向の累積とr+方向の累積を実行する。
c+方向の累積
 0  0  0  0  0
 0  1  2  3  3
 0  1  4  7  9
 0  1  4  7  9
 0  0  0  0  0

r+方向の累積
 0  0  0  0  0
 0  1  2  3  3
 0  2  6 10 12
 0  3 10 17 21
 0  3 10 17 21

出来た累積和配列を使って矩形(2,2)-(3,3)の総和を求めてみます。

B[3][3] + (-B[2-1][3]) + (-B[3][2-1]) +B[2-1][2-1]
= 17 + (-3) + (-3) + 1 = 12

 0  0  0  0
 0  1  1  1
 0  1  3  3
 0  1  3  3
-
 0  0  0  0
 0  1  1  1
-
 0  0
 0  1
 0  1
 0  1
+
 0  0
 0  1	(1回余分に引いているので足す)
=
 _  _  _  _
 _  _  _  _
 _  _  3  3
 _  _  3  3
矩形(2,2)-(3,3)の総和が累積和から求められました。
矩形(r1,c1)-(r2,c2)の領域の総和を
累積和配列Bを以下のように用いて求められます。
B[r2][c2]-B[r1-1][c2]-B[r2][c1-1] +B[r1-1][c1-1]
ただし、累積和の方向はc+,r+
</pre>
矩形の総和は、こっちのが詳しい
<a href="http://paiza.hatenablog.com/entry/2014/05/28/%E3%82%82%E3%81%97%E5%A5%B3%E5%AD%90%E5%A4%A7%E7%94%9F%E3%83%97%E3%83%AD%E3%82%B0%E3%83%A9%E3%83%9E%E3%81%AB%E3%80%8E%E3%82%A2%E3%83%AB%E3%82%B4%E3%83%AA%E3%82%BA%E3%83%A0%E3%80%8F%E3%82%92%E5%9B%B3">ぱいざ</a>
の「O(H^2 W^2 ) の解法」を参照
<!--ここまで-->