「オイラープロジェクト321～330」の編集履歴（バックアップ）一覧はこちら

オイラープロジェクト321～330 - (2012/09/29 (土) 17:51:17) の１つ前との変更点

追加された行は緑色になります。

削除された行は赤色になります。

*問い324
3x3xn の塔を 2x1x1 のブロックで埋める場合の数を f(n) とする。
ブロックは好きな方向に回転することができる；しかし、塔自身を回転・反転等させたものは区別して数える。
たとえば（q=100000007 として）：
f(2) = 229,
f(4) = 117805,
f(10) mod q = 96149360,
f(10^3) mod q = 24806056,
f(10^6) mod q = 30808124.
f(10^10000) mod 100000007 を求めよ。



何か漸化式があるのだろうか？
数式を思いつかないのでとりあえず動的計画法で解こうとしてみる。

２段を作り、２段を２つ結合で4段、４段を２つ結合して８段、８段を２つ結合して１６段、、として10＾10000段まで到達し後は10＾10000を２進数に分解２進数で1のついた段だけを掛け算していけば理論上計算できるはず。
多分相当遅いが自分はこれしか思いつかない。
１０億段くらいの小さい数だったら簡単なのにな。

とりあえず小さな数で正しい答えを出せるところまで持ってきた。
最初とても遅い動的計画法から始め何とか高速化して10＾３や10^6段で正しい答えを出せるところまで前進した。
だがまだ遅いとても遅い。
10＾6段を計算するのに役１０分。
10＾50段程度なら８０分、10＾10000段、２進数で3万数千ケタを取り扱うと16000分も計算時間がいることになる。
時間にすると約２７０時間だ。
これでは話にならない。

色々考えて数値をピーピングしてみた結果16段以上からは組み合わせの集合が5種類に分割されているような感じがする。

5つのまとまりに分かれていて、まとまり内部で組み合わせの和をとる感じである。
一度５つにわかれるとまとまりはまとまり内部でのみ組み合わせが出来て外部とは関係を持たなくなる？
のかな？
何かの構造を持ってるようなのだがそれが何か分からない。


今のところ次の２倍段を求めるのに必要な計算は全て行ってるし必要ない計算は全て省いている。
動的計画法ではここから先に進めないような感じもする？
それとも何かまだ纏めて計算できる何かがあるのか？
残るは纏めて計算できる何かを探すくらいだろうか？
一番いいのは漸化式だろうけど問題が複雑すぎてそのアプローチは中々手が出ない。



 #include<stdio.h>
 #include<map>
 #include<set>
 #include<string.h>
 #include<iostream>
 #include<vector>
 const __int64 base=100000007;
 struct S{
	int blocks[2][3][3];//上面、下面のみをもってブロックを判別する
	//立てた棒は11、横か縦の棒は1
	//上面は隙間を許さない
	bool operator<(const S& s)const{
		for(int h=0;h<2;h++){
			for(int r=0;r<3;r++){
				for(int c=0;c<3;c++){
					if(blocks[h][r][c]!=s.blocks[h][r][c])return blocks[h][r][c]<s.blocks[h][r][c];
				}
			}
		}
		return false;
	}
 };
 std::map<S,__int64> memo,nextMemo,ansMemo,ansNext;
 std::map<S,std::vector<S> > okPerms,nextOKPerms;
 std::vector<S> vec;
 bool isConnectOK(S& s1,S& s2){
	//ｓ1が下のブロック、ｓ２が上のブロック
	int upSide,downSide;
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			upSide=s1.blocks[1][r][c];//接合部下面
			downSide=s2.blocks[0][r][c];//接合部上面
			if(downSide==0&&upSide!=11)return false;//結合部上面に隙間があった
			if(downSide!=0&&upSide==11)return false;//結合部上面と下面に出っ張りがあった
		}
	}
	//ぴったり結合出来た
	return true;
 }
 S connectBlock(S& s1,S& s2){
	S s3;
	//結合可能なので結合し上面と下面をもって種類を判別する
	memcpy(s3.blocks[1],s2.blocks[1],sizeof(s2.blocks[1]));
	memcpy(s3.blocks[0],s1.blocks[0],sizeof(s1.blocks[0]));
	return s3;
 }
 void first(int deep,int perm[2][3][3],int sukima,int hBar){
	//ブロックの組み合わせ２段を作る
	//上面に隙間は許されない
	if(deep>=18){
		S s;
		memcpy(s.blocks,perm,sizeof(s.blocks));
		if(memo.find(s)==memo.end()){
			memo[s]=1;
		}else{
			memo[s]++;
		}
		return ;
	}else{
		int h=deep/9;
		int r=(deep%9)/3;
		int c=deep%3;
		//この部分にはもうブロックがあった
		if(perm[h][r][c]!=0){
			first(deep+1,perm,sukima,hBar);
			return ;
		}
		//ここを隙間にする、隙間が許されるのは下面のみ
		if(h==0){
			first(deep+1,perm,sukima+1,hBar);
		}
		
		//横
		perm[h][r][c]=1;
		if(c+1<3&&perm[h][r][c+1]==0){
			perm[h][r][c+1]=1;
			first(deep+1,perm,sukima,hBar);
			perm[h][r][c+1]=0;
		}
		//縦
		
		if(r+1<3&&perm[h][r+1][c]==0){
			perm[h][r+1][c]=1;
			first(deep+1,perm,sukima,hBar);
			perm[h][r+1][c]=0;
		}
		
		//高さ
		perm[h][r][c]=11;
		if(h==0){
			perm[h][r][c]=1;
			perm[h+1][r][c]=1;
		}
		first(deep+1,perm,sukima,hBar+1);
		
		if(h==0)perm[h+1][r][c]=0;
		perm[h][r][c]=0;
	}
	
 }
 bool downCheck(S& s1,S& s2){
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			if(s1.blocks[0][r][c]!=s2.blocks[0][r][c])return false;
		}
	}
	return true;
 }
 S createNext(S& s1){
	S s2;
	memset(s2.blocks,0,sizeof(s2.blocks));
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			if(s1.blocks[1][r][c]==11)s2.blocks[0][r][c]=0;
			else s2.blocks[0][r][c]=1;
		}
	}
	return s2;
 }
 void ansCalc(char c){
	if(c=='0')return;
	if(ansMemo.empty()){
		ansMemo.insert(memo.begin(),memo.end());
	}else{
		ansNext.clear();
		std::map<S,__int64>::iterator it1,it2;
		S s1,s2,s3,s4;
		for(it1=ansMemo.begin();it1!=ansMemo.end();it1++){
			s1=(*it1).first;
			s2=createNext(s1);
			for(it2=memo.lower_bound(s2);it2!=memo.end()&&downCheck((*it2).first,s2);it2++){
				s4=connectBlock(s1,(*it2).first);
				if(ansNext.find(s4)==ansNext.end()){
					ansNext[s4]=((*it1).second*(*it2).second)%base;
				}else{
					ansNext[s4]=(ansNext[s4]+((*it1).second*(*it2).second)%base)%base;;
				}
			}
		}
		ansMemo.clear();
		ansMemo.insert(ansNext.begin(),ansNext.end());
	}
 }
 bool lastCheck(S& s1){
	int t;
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			if(s1.blocks[0][r][c]==0)return false;
			if(s1.blocks[1][r][c]!=1)return false;
		}
	}
	return true;
 }
 int main(){
	int perm[2][3][3];
	memset(perm,0,sizeof(perm));
	first(0,perm,0,0);
	std::map<S,__int64>::iterator it1,it2,itD,itU;
	S s1,s2,s3,s4;
	//char text[30]="01111100010000101111";//10^6の逆転
	char text[30]="0001011111";//10＾３の逆転
	int len=strlen(text);
	std::cout<<len<<"\n";
	int keta;
	for(keta=1;keta<len;keta++){
		ansCalc(text[keta]);
		nextMemo.clear();
		for(it1=memo.begin();it1!=memo.end();it1++){
			s1=(*it1).first;
			s2=createNext(s1);
			for(it2=memo.lower_bound(s2);it2!=memo.end()&&downCheck((*it2).first,s2);it2++){
				s4=connectBlock(s1,(*it2).first);
				if(nextMemo.find(s4)==nextMemo.end()){
					nextMemo[s4]=((*it1).second*(*it2).second)%base;
				}else{
					nextMemo[s4]=(nextMemo[s4]+((*it1).second*(*it2).second)%base)%base;;
				}
			}
		}
		memo.clear();
		memo.insert(nextMemo.begin(),nextMemo.end());
		std::cout<<keta<<" "<<ansMemo.size()<<" "<<memo.size()<<"\n";
	}
	
	//蓋はいらなくなった
	
	__int64 ans=0,t;
	for(it1=ansMemo.begin();it1!=ansMemo.end();it1++){
		s1=(*it1).first;
		if(lastCheck(s1)==true)ans=(ans+(*it1).second)%base;
	}
	std::cout<<"\nans="<<ans;
 }

*問い324
3x3xn の塔を 2x1x1 のブロックで埋める場合の数を f(n) とする。
ブロックは好きな方向に回転することができる；しかし、塔自身を回転・反転等させたものは区別して数える。
たとえば（q=100000007 として）：
f(2) = 229,
f(4) = 117805,
f(10) mod q = 96149360,
f(10^3) mod q = 24806056,
f(10^6) mod q = 30808124.
f(10^10000) mod 100000007 を求めよ。



何か漸化式があるのだろうか？
数式を思いつかないのでとりあえず動的計画法で解こうとしてみる。

２段を作り、２段を２つ結合で4段、４段を２つ結合して８段、８段を２つ結合して１６段、、として10＾10000段まで到達し後は10＾10000を２進数に分解２進数で1のついた段だけを掛け算していけば理論上計算できるはず。
多分相当遅いが自分はこれしか思いつかない。
１０億段くらいの小さい数だったら簡単なのにな。

とりあえず小さな数で正しい答えを出せるところまで持ってきた。
最初とても遅い動的計画法から始め何とか高速化して10＾３や10^6段で正しい答えを出せるところまで前進した。
だがまだ遅いとても遅い。
10＾6段を計算するのに役１０分。
10＾50段程度なら８０分、10＾10000段、２進数で3万数千ケタを取り扱うと16000分も計算時間がいることになる。
時間にすると約２７０時間だ。
これでは話にならない。

色々考えて数値をピーピングしてみた結果16段以上からは組み合わせの集合が5種類に分割されているような感じがする。

5つのまとまりに分かれていて、まとまり内部で組み合わせの和をとる感じである。
一度５つにわかれるとまとまりはまとまり内部でのみ組み合わせが出来て外部とは関係を持たなくなる？
のかな？
何かの構造を持ってるようなのだがそれが何か分からない。


今のところ次の２倍段を求めるのに必要な計算は全て行ってるし必要ない計算は全て省いている。
動的計画法ではここから先に進めないような感じもする？
それとも何かまだ纏めて計算できる何かがあるのか？
残るは纏めて計算できる何かを探すくらいだろうか？
一番いいのは漸化式だろうけど問題が複雑すぎてそのアプローチは中々手が出ない。


とりあえず目的の形を絶対に生成できない組み合わせを排除するコードを書いてみた。
少し早くなったがコードは膨らむし思ったほど速度が出ないしであまり楽しくはない。



 #include<stdio.h>
 #include<map>
 #include<set>
 #include<string.h>
 #include<iostream>
 #include<vector>
 #include<queue>
 
 const __int64 base=100000007;
 
 struct S{
	int blocks[2][3][3];//上面、下面のみをもってブロックを判別する
	//立てた棒は11、横か縦の棒は1
	//上面は隙間を許さない
	bool operator<(const S& s)const{
		for(int h=0;h<2;h++){
			for(int r=0;r<3;r++){
				for(int c=0;c<3;c++){
					if(blocks[h][r][c]!=s.blocks[h][r][c])return blocks[h][r][c]<s.blocks[h][r][c];
				}
			}
		}
		return false;
	}
 };
 struct S2{
	int blocks[2][3][3];//上面、下面のみをもってブロックを判別する
	//立てた棒は11、横か縦の棒は1
	//上面は隙間を許さない
	bool operator<(const S2& s)const{
		for(int h=1;h>=0;h--){
			for(int r=0;r<3;r++){
				for(int c=0;c<3;c++){
					if(blocks[h][r][c]!=s.blocks[h][r][c])return blocks[h][r][c]<s.blocks[h][r][c];
				}
			}
		}
		return false;
	}
 };
 

 std::map<S,__int64> memo,nextMemo,ansMemo,ansNext;
 bool isConnectOK(S& s1,S& s2){
	//ｓ1が下のブロック、ｓ２が上のブロック
	int upSide,downSide;
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			upSide=s1.blocks[1][r][c];//接合部下面
			downSide=s2.blocks[0][r][c];//接合部上面
			if(downSide==0&&upSide!=11)return false;//結合部上面に隙間があった
			if(downSide!=0&&upSide==11)return false;//結合部上面と下面に出っ張りがあった
		}
	}
	//ぴったり結合出来た
	return true;
 }
 S connectBlock(S& s1,S& s2){
	S s3;
	//結合可能なので結合し上面と下面をもって種類を判別する
	memcpy(s3.blocks[1],s2.blocks[1],sizeof(s2.blocks[1]));
	memcpy(s3.blocks[0],s1.blocks[0],sizeof(s1.blocks[0]));
	return s3;
 }
 void first(int deep,int perm[2][3][3],int sukima,int hBar){
	//ブロックの組み合わせ２段を作る
	//上面に隙間は許されない
	
	
	
	//答えに到達できない組み合わせは無視する
	if(deep==9){
		if(sukima%2==1)return ;
		hBar=0;
	}
	if(deep==18){
		if(hBar%2==1)return ;
	}
 	
 	
	if(deep>=18){
		S s;
		memcpy(s.blocks,perm,sizeof(s.blocks));
		if(memo.find(s)==memo.end()){
			memo[s]=1;
		}else{
			memo[s]++;
		}
		return ;
	}else{
		int h=deep/9;
		int r=(deep%9)/3;
		int c=deep%3;
		//この部分にはもうブロックがあった
		if(perm[h][r][c]!=0){
			first(deep+1,perm,sukima,hBar);
			return ;
		}
		//ここを隙間にする、隙間が許されるのは下面のみ
		if(h==0){
			first(deep+1,perm,sukima+1,hBar);
		}
 		
		//横
		perm[h][r][c]=1;
		if(c+1<3&&perm[h][r][c+1]==0){
			perm[h][r][c+1]=1;
			first(deep+1,perm,sukima,hBar);
			perm[h][r][c+1]=0;
		}
		//縦
 		
		if(r+1<3&&perm[h][r+1][c]==0){
			perm[h][r+1][c]=1;
			first(deep+1,perm,sukima,hBar);
			perm[h][r+1][c]=0;
		}
 		
		//高さ
		perm[h][r][c]=11;
		if(h==0){
			perm[h][r][c]=1;
			perm[h+1][r][c]=1;
		}
		first(deep+1,perm,sukima,hBar+1);
 		
		if(h==0)perm[h+1][r][c]=0;
		perm[h][r][c]=0;
	}
	
 }
 bool downCheck(S& s1,S& s2){
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			if(s1.blocks[0][r][c]!=s2.blocks[0][r][c])return false;
		}
	}
	return true;
 }
 bool upCheck(S2& s1,S2& s2){
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			if(s1.blocks[1][r][c]!=s2.blocks[1][r][c])return false;
		}
	}
	return true;
 }
 S createNext(S& s1){
	S s2;
	memset(s2.blocks,0,sizeof(s2.blocks));
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			if(s1.blocks[1][r][c]==11)s2.blocks[0][r][c]=0;
			else s2.blocks[0][r][c]=1;
		}
	}
	return s2;
 }
 S2 createNext(S2& ss){
	S2 ss2;
	memset(ss2.blocks,0,sizeof(ss2.blocks));
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			if(ss.blocks[0][r][c]==0)ss2.blocks[1][r][c]=11;
			else ss2.blocks[1][r][c]=1;
		}
	}
	return ss2;
 }
 void ansCalc(char c){
	if(c=='0')return;
	if(ansMemo.empty()){
		//ansMemo.insert(memo.begin(),memo.end());
		std::map<S,__int64>::iterator it1;
		S s1;
		for(it1=memo.begin();it1!=memo.end();it1++){
			s1=(*it1).first;
			bool ok=true;
			for(int r=0;r<3;r++){
				for(int c=0;c<3;c++){
					if(s1.blocks[0][r][c]==0)ok=false;
				}
			}
			if(ok==true)ansMemo[s1]=(*it1).second;
		}
	}else if(!memo.empty()){
		ansNext.clear();
		std::map<S,__int64>::iterator it1,it2;
		S s1,s2,s3,s4;
		for(it1=ansMemo.begin();it1!=ansMemo.end();it1++){
			s1=(*it1).first;
			s2=createNext(s1);
			for(it2=memo.lower_bound(s2);it2!=memo.end()&&downCheck((*it2).first,s2);it2++){
				s4=connectBlock(s1,(*it2).first);
				if(ansNext.find(s4)==ansNext.end()){
					ansNext[s4]=((*it1).second*(*it2).second)%base;
				}else{
					ansNext[s4]=(ansNext[s4]+((*it1).second*(*it2).second)%base)%base;;
				}
			}
		}
		ansMemo.clear();
		ansMemo.insert(ansNext.begin(),ansNext.end());
	}
 }
 bool lastCheck(S& s1){
	int t;
	for(int r=0;r<3;r++){
		for(int c=0;c<3;c++){
			if(s1.blocks[0][r][c]==0)return false;
			if(s1.blocks[1][r][c]!=1)return false;
		}
	}
	return true;
 }
 void badPermDell(){
	std::map<S,__int64>::iterator it1;
	std::queue<S2> qu;
	std::set<S2> allPerm,okPerms,dellPerms;
	std::set<S2>::iterator sIt;
	S2 ss,ss2,ss3;
	S s1;
	for(it1=memo.begin();it1!=memo.end();it1++){
		s1=(*it1).first;
		memcpy(ss.blocks,s1.blocks,sizeof(ss.blocks));
		allPerm.insert(ss);
		dellPerms.insert(ss);
		bool ok=true;
		for(int r=0;r<3;r++){
			for(int c=0;c<3;c++){
				int k=ss.blocks[1][r][c];
				if(k==0||k==11)ok=false;
			}
		}
		if(ok==true)qu.push(ss);
	}
	while(!qu.empty()){
		ss=qu.front();
		qu.pop();
		if(okPerms.find(ss)!=okPerms.end())continue;
		okPerms.insert(ss);
		dellPerms.erase(ss);
		ss2=createNext(ss);
		for(sIt=allPerm.lower_bound(ss2);sIt!=allPerm.end()&&upCheck(ss2,(*sIt));sIt++){
			ss3=(*sIt);
			if(okPerms.find(ss3)==okPerms.end()){
				qu.push(ss3);
			}
		}
	}
	for(sIt=dellPerms.begin();sIt!=dellPerms.end();sIt++){
		ss=(*sIt);
		memcpy(s1.blocks,ss.blocks,sizeof(s1.blocks));
		nextMemo.erase(s1);
	}
 }
 int main(){
	int perm[2][3][3];
	memset(perm,0,sizeof(perm));
	first(0,perm,0,0);
	std::map<S,__int64>::iterator it1,it2;
	S s1,s2,s3,s4;
	//char text[30]="01111100010000101111";//10^6の逆転
	char text[30]="0001011111";//10＾３の逆転
	int len=strlen(text);
	std::cout<<len<<"\n";
	int keta;
 	
	for(keta=1;keta<len;keta++){
		ansCalc(text[keta]);
		nextMemo.clear();
		int count=0;
		for(it1=memo.begin();it1!=memo.end();it1++){
			s1=(*it1).first;
			s2=createNext(s1);
			for(it2=memo.lower_bound(s2);it2!=memo.end()&&downCheck((*it2).first,s2);it2++){
				s4=connectBlock(s1,(*it2).first);
				if(nextMemo.find(s4)==nextMemo.end()){
					nextMemo[s4]=((*it1).second*(*it2).second)%base;
				}else{
					nextMemo[s4]=(nextMemo[s4]+((*it1).second*(*it2).second)%base)%base;;
				}
				count++;
			}
		}
		if(keta==3){
			badPermDell();
		}
 		
		memo.clear();
		memo.insert(nextMemo.begin(),nextMemo.end());
		std::cout<<keta<<" "<<count<<" "<<" ans="<<ansMemo.size()<<" memo="<<memo.size()<<"\n";
 	}
 	
 	
	__int64 ans=0,t;
 	for(it1=ansMemo.begin();it1!=ansMemo.end();it1++){
		s1=(*it1).first;
 		if(lastCheck(s1)==true)ans=(ans+(*it1).second)%base;
	}
 	std::cout<<"\nans="<<ans;
 }