根据http://theory.stanford.edu/~virgi/cs267/lecture1.pdf论文第二部分的伪代码实现出来的C++实现代码,有改动,原文第5步没有还原。
#include<iostream>
#include<vector>
#include<time.h>
using namespace std;
#define MATRIX vector<vector<int>>
#define DEBUG
vector<int> getDegree(MATRIX m){
int p = m.size();
vector<int> degree(p,0);
for (int i = 0; i < p; i++){
for (int j = 0; j < p; j++){
degree[i] += m[i][j];
}
}
return degree;
}
MATRIX multi(MATRIX& M1, MATRIX& M2){
int rowN = M1.size();
int columnN = M2[0].size();
MATRIX result;
for (int i = 0; i < rowN; i++){
result.push_back(vector<int>(columnN , 0));
for (int j = 0; j < columnN; j++){
for (int k = 0; k < M2.size(); k++){
result[i][j] += M1[i][k] * M2[k][j];
}
}
}
return result;
}
MATRIX transposition(MATRIX& M){
int rowN = M.size();
int columnN = M[0].size();
MATRIX result;
for (int j = 0; j < columnN; j++){
result.push_back(vector<int>(rowN));
for (int i = 0; i < rowN; i++){
result[j][i] = M[i][j];
}
}
return result;
}
void printMatrix(MATRIX& M){
for (MATRIX::iterator it = M.begin(); it != M.end(); it++){
for (vector<int>::iterator itN = (*it).begin(); itN != (*it).end(); itN++){
cout << *itN << " ";
}
cout << endl;
}
}
bool checkM(MATRIX& M,MATRIX& MA,MATRIX& MB){
MATRIX C = multi(M, transposition(multi(M, MB)));
int rowN = MA.size();
int columnN = MA[0].size();
bool flag = true;
for (int i = 0; i < rowN; i++){
for (int j = 0; j < columnN; j++){
if (MA[i][j] == 1 && C[i][j] != 1){
flag = false;
break;
}
}
if (!flag)
break;
}
return flag;
}
void ullmann(MATRIX &ma, MATRIX &mb){
int pa = ma.size(), qa = ma[0].size();
int pb = mb.size(), qb = mb[0].size();
vector<int> degreeA = getDegree(ma);
vector<int> degreeB = getDegree(mb);
vector<vector<int>> M0;
for (int i = 0; i < pa; i++){
M0.push_back(vector<int>(pb));
for (int j = 0; j < pb; j++){
M0[i][j] = degreeB[j] >= degreeA[i] ? 1 : 0;
}
}
vector<int> F(pb, 0);
vector<int> H(pa, -1);
vector<MATRIX> M(pa);
//step 1
MATRIX CM = M0;
int d = 0;
H[0] = -1;
int k = -1;
step2:
bool existJ = false;
for (int j = 0; j < pb; j++){
if (CM[d][j] == 1 && F[j] == 0){
existJ = true;
break;
}
}
if (existJ){
M[d] = CM;
k = d == 0 ? H[0] : -1;
}
else{
goto step7;
}
step3:
k = k + 1;
if (CM[d][k] == 0 || F[k] == 1)
goto step3;
H[d] = k;//
for (int j = 0; j < pb; j++){
if (j != k){
CM[d][j] = 0;
}
}
step4:
if (d < pa - 1 ){
goto step6;
}else{
if (checkM(CM, ma, mb)){
#ifdef DEBUG
//printMatrix(CM);
static int sum = 0;//the number of answers
cout << "---#" << ++sum << "---" << endl;
for (int n = 0; n < H.size(); n++){
cout << H[n]+1 << " ";
}
cout << endl;
#endif
}
}
step5:
F[k] = 0;
CM = M[d];//复原
existJ = false;
for (int j = k + 1; j < pb; j++){
if (CM[d][j] == 1 && F[j] == 0){
existJ = true;
break;
}
}
if (existJ){
goto step3;
}
else{
goto step7;
}
step6:
H[d] = k;
F[k] = 1;
//set F by H is right,or reset F[k] to 0 when search a new k in step 7.
/*for (int i = 0; i < F.size(); i++)
F[i] = 0;
for (int i = 0; i <= d; i++){
F[H[i]] = 1;
}*/
d = d + 1;
goto step2;
step7:
if (d == 0){
goto END;
}
else{
H[d] = -1;
F[k] = 0;
d = d - 1;
// CM = M[d];
k = H[d];
goto step5;
}
END:
#ifdef DEBUG
cout <<"---end---"<< endl;
#endif
return;
}
MATRIX getMatrix(int** a, int m, int n){
MATRIX result;
for (int i = 0; i < m; i++){
result.push_back(vector<int>((int*)a + n*i, (int*)a + n*i + n));
}
return result;
}
int main(){
int a[3][3] = { { 0, 1, 1 }, { 1, 0, 0 }, { 1, 0, 0 } };
int b[4][4] = { { 0, 1, 1, 0 }, { 1, 0, 1, 1 }, { 1, 1, 0, 1 }, { 0, 1, 1, 0 } };
MATRIX ma = getMatrix((int **)a, 3, 3);
MATRIX mb = getMatrix((int **)b, 4, 4);
clock_t t1 = clock();
ullmann(ma, mb);
clock_t t2 = clock();
cout << "ull time is : " << t2 - t1 << endl;
return 0;
}