【机器学习】PangRank的c++和python实现

ystraw_ah

已于 2022-04-01 10:52:31 修改

阅读量826

点赞数 1

分类专栏：机器学习图计算文章标签： python

于 2020-07-20 16:35:11 首次发布

本文链接：https://blog.csdn.net/qq_39451578/article/details/107466453

版权

机器学习同时被 2 个专栏收录

4 篇文章 1 订阅

订阅专栏

图计算

4 篇文章 0 订阅

订阅专栏

详细讲解视频推荐

文章目录

原理公式：

PageRank的一般公式：
$\frac {1-d} {n} 1}$

其中d为随机概率([0,1]),一般去0.85，M为基本转移矩阵， ${\frac {1-d} {n} 1}$ 表示完全随机转移矩阵，表示从任意一个节点到任意一个节点的转移概率都是1/n.

0.network x

利用python的第三方库networkx实现的pagerank算法, 这就非常简单,可以用于测试自己写的pagerank算法是否正确.参考博客

import networkx as nx
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

plt.figure(figsize=(10,5))
G = nx.DiGraph()
G.add_nodes_from(['0', '1', '2', '3'])    #加点集合
G.add_edge('0' ,'2' )   #一次添加一条边
G.add_edge('0' ,'1' )   #一次添加一条边
G.add_edge('0' ,'3' )   #一次添加一条边
G.add_edge('1' ,'0' )   #一次添加一条边
G.add_edge('1' ,'3' )   #一次添加一条边
G.add_edge('2' ,'0' )   #一次添加一条边
G.add_edge('3' ,'1' )   #一次添加一条边
G.add_edge('3' ,'2' )   #一次添加一条边
nx.draw(G, with_labels=True , node_size = 2000)
plt.show()

pr = nx.pagerank(G, alpha=0.85, personalization=None,
             max_iter=100, tol=1.0e-6, nstart=None, weight='weight',
             dangling=None)
print(pr)

输出结果:
在这里插入图片描述

{'0': 0.3245609358176831, '1': 0.22514635472743894, '2': 0.22514635472743894, '3': 0.22514635472743894}

1.1 c++实现代码(邻接矩阵)：

/*
Function: c++实现PageRank算法
Author: ystraw
Time: 2020-7-20 13:58
*/
#include<iostream>
#include<fstream>
#include<cmath>
#include <cstring>
using namespace std;
#define MAX 10000
int A[MAX][MAX];   // 关系矩阵 可以改为bool类型
float M[MAX][MAX]; // 权重矩阵
float newPR[MAX];  // 存储pr 
float oldPR[MAX];
float d = 0.80;      // 初始化参数

void getM(int n){
	// 计算初始权重矩阵 
	for(int i = 0; i < n; i++){
		int s = 0; // 统计出度 
		for(int j = 0; j < n; j++){
			s += A[j][i]; // 每列求和 
		}
		if(s == 0){  // 处理Dead Ends: 没有出度 
			for(int j = 0; j < n; j++){
				M[j][i] = 1.0 / n;
			}	
		}
		else{
			float fm = s;
			for(int j = 0; j < n; j++){
				M[j][i] = A[j][i] / fm;
			}	
		}
	}
//	for(int i = 0 ; i < n; i++){
//		for(int j = 0; j < n; j++){
//			cout << M[i][j] << " ";
//		}
//		cout << endl;
//	}
}

void getBaseLev(int n){
	//计算完全随机转移矩阵: bl = (1-d)/n*1
	for(int i = 0; i < n; i++){
		newPR[i] = 1.0/n;
	}
//	for(int i = 0; i < n; i++){
//		cout << newPR[i] << endl;
//	} 
}

void getPR(float *a, float *b, int n){
	// d * M @ oldPR + BL
	for(int i = 0; i < n; i++){
		a[i] = 0;
		for(int j = 0; j < n; j++){
			a[i] += M[i][j] * b[j];
		}
		a[i] = a[i] * d + 1.0/n * (1-d);
	}
//	for(int i = 0; i < n; i++){
//		cout << a[i] << " ";
//	} 
//	cout << endl;
} 

int main(){
	// 迭代终止阈值
    float threshold = 1e-7;
	// 读入/输出数据 
	ifstream fin(".\\edges.txt");
	ofstream fout(".\\output.txt"); 
	int i, a, b;
	int n = 0; // 顶点数据量 
	memset(A,0,sizeof(A));
	// 输出有向边 
	for(i=0; fin>>a>>b; ++i){ 
		A[b][a] = 1; //a->b
		a = a < b ? b : a;
		n = n < a ? a : n; //最大值为顶点数量 
//		cout << a << " " << b << " " << n << endl;
	}
	n += 1; // 数量为最大值+1 
	cout << n << endl;
	// 通过邻接矩阵计算状态转移矩阵M:
    getM(n);
	// 得到完全随机转移矩阵
    getBaseLev(n);
    // 记录迭代次数
	int cnt = 0;
	while(1){
		cnt++;
		if(cnt % 2 == 1){
			getPR(oldPR, newPR, n); // newPR迭代存入oldPR 
		}
		else{
			getPR(newPR, oldPR, n); 
		}
		// 检测是否结束：
		float wc = 0;
		for(int i = 0; i < n; i++){
			wc += abs(oldPR[i] - newPR[i]);
		} 
//		cout << cnt << ": " << wc << " " << threshold << endl; 
		if(wc < threshold){
			break;
		}
	} 
	// 将结果写入 
	if(cnt % 2 == 1){
		for(int i=0; i<n; i++)
			fout<< i << ' ' << oldPR[i] << endl;
	}
	else{
		for(int i=0; i<n; i++)
			fout << i << ' ' << newPR[i] << endl;
	}
	return 0;
}

测试样例：

edges.txt:
0 2
0 1
0 3
1 0
1 3
2 0
3 1
3 2
output.txt:
0 0.33333
1 0.22222
2 0.22222
3 0.22222

edges.txt:
0 2
0 1
0 3
1 0
1 3
2 2
3 1
3 2
output.txt:
0 0.101351
1 0.128378
2 0.641892
3 0.128378

1.2 c++实现代码(邻接链表_直接申请内存)：

对于数据量太大的数据集，申请太大的数组会报错，所以将上面的代码进行了修改：

/*
Function: c++实现PageRank算法, 邻接链表 
Author: ystraw
Time: 2020-7-22 13:58
*/

#include <iostream>
#include <vector>
#include <cstring>
#include<fstream>
#include<cmath>
using namespace std;
#define MAX 1000000

vector <int> in_ver[MAX]; // 存储各点入度的点 
int out_num[MAX]; // 存储各点出度数量 
float newPR[MAX]; // 存储pr 
float oldPR[MAX];
float d = 0.80;         // 初始化参数
float threshold = 1e-7; // 迭代终止阈值

void getPR(float *a, float *b, int n){
	for(int i = 0; i < n; i++){
		a[i] = 0;
		int len = in_ver[i].size();
		for(int j = 0; j < len; j++){
			a[i] += b[in_ver[i][j]] / out_num[in_ver[i][j]]; // 累加：该点入度点的pr值除以该点的出度 
		}
		a[i] = a[i] * d + (1 - d) / n; 
	}
}

float* computePR(int n, char* outpath){ 
	// 处理没有出度的点
	for(int i = 0; i < n; i++){
		if(out_num[i] == 0){ // 如果没有出度，则将其指向另外的所有点 
			out_num[i] = n;
			for(int j = 0; j < n; j++){
				in_ver[j].push_back(i);
			}
		}
		newPR[i] = 1.0/n;
	}  
    // 记录迭代次数
	int cnt = 0;
	while(1){
		cnt++;
		if(cnt % 2 == 1){
			getPR(oldPR, newPR, n); // newPR迭代存入oldPR 
		}
		else{
			getPR(newPR, oldPR, n); 
		}
		// 检测是否结束：
		float wc = 0;
		for(int i = 0; i < n; i++){
			wc += abs(oldPR[i] - newPR[i]);
		} 
		cout << cnt << ": " << wc << " " << threshold << endl; 
		if(wc < threshold){
			break;
		}
	} 
	// 将结果写入 
	if(*outpath != '\0'){
//		ofstream fout(".\\output.txt");
		cout << "outpath:" << outpath << endl;
		ofstream fout(outpath);
		if(cnt % 2 == 1){
			for(int i=0; i<n; i++)
				fout<< i << ' ' << oldPR[i] << endl;
		}
		else{
			for(int i=0; i<n; i++)
				fout << i << ' ' << newPR[i] << endl;
		}
	}
	// 将结果返回
	if(cnt % 2 == 1){
		return oldPR;
	}
	else{
		return newPR;
	}
}

int main(){
	// 读入/输出数据 
	ifstream fin(".\\edges.txt");
	int n = 0; // 顶点数据量 
	int a, b; 
	// 输入有向边, 获取邻接矩阵 
	for(int i=0; fin>>a>>b; ++i){ 
		in_ver[b].push_back(a); //a->b, 记录每个点的入度点 
		out_num[a]++; // 记录每个点的出度数量 
		a = a < b ? b : a;
		n = n < a ? a : n; //最大值为顶点数量 
//		cout << a << " " << b << " " << n << endl;
	}
	n++; 
	// 进行计算 
	float *pr = computePR(n, ".\\output.txt");
	for(int i = 0; i <= n; i++){
		cout << pr[i] << " ";
	}
 
	return 0;
}

1.2 c++实现代码(邻接链表_动态申请内存)：

#include <iostream>
#include <vector>
#include <cstring>
#include<fstream>
#include<cmath>
#include <cstdlib>
using namespace std;
#define MAX 100000
 
vector <vector<int> > in_ver;// 存储各点入度的点 
int *out_num; // 存储各点出度数量 
float *newPR; // 存储pr 
float *oldPR;
float d = 0.80;         // 初始化参数
float threshold = 1e-4; // 迭代终止阈值

void getPR(float *a, float *b, int n){
	for(int i = 0; i < n; i++){
		a[i] = 0;
		int len = in_ver[i].size();
		for(int j = 0; j < len; j++){
			a[i] += b[in_ver[i][j]] / out_num[in_ver[i][j]]; // 累加：该点入度点的pr值除以该点的出度 
		}
		//a[i] = a[i] * d + (1 - d) / n;
		a[i] = a[i] * d + 1 - d;  // (1 - d) 不除 n 
	}
}

float* computePR(int n, char* outpath){ 
	// 处理没有出度的点
	for(int i = 0; i < n; i++){
		if(out_num[i] == 0){ // 如果没有出度，则将其指向另外的所有点 
			out_num[i] = n;
			for(int j = 0; j < n; j++){
				in_ver[j].push_back(i);
			}
		}
		//newPR[i] = 1.0/n; // 初始化为 1.0/n 
		newPR[i] = 1.0; // 初始化1
	}  
    // 记录迭代次数
	int cnt = 0;
	while(cnt < 500){  // 改为固定迭代500次 
		cnt++;
		if(cnt % 2 == 1){
			getPR(oldPR, newPR, n); // newPR迭代存入oldPR 
		}
		else{
			getPR(newPR, oldPR, n); 
		}
		// 检测是否结束：
		float wc = 0;
		for(int i = 0; i < n; i++){
			wc += abs(oldPR[i] - newPR[i]);
		} 
//		cout << cnt << ": " << wc << " " << threshold << endl; 
		//if(wc < threshold){ // 收敛判断
		//	break;
		//}
	} 
	// 将结果写入 
	if(*outpath != '\0'){
//		ofstream fout(".\\output.txt");
		cout << "outpath:" << outpath << endl;
		ofstream fout(outpath);
		if(cnt % 2 == 1){
			for(int i=0; i<n; i++)
				fout<< i << ' ' << oldPR[i] << endl;
		}
		else{
			for(int i=0; i<n; i++)
				fout << i << ' ' << newPR[i] << endl;
		}
	}
	// 将结果返回
	if(cnt % 2 == 1){
		return oldPR;
	}
	else{
		return newPR;
	}
}

int main(){
	// 读入/输出数据 
	ifstream fin(".\\edges.txt");
//	ifstream fin(".\\google_90w_edge.txt");
	int n = 0; // 顶点数据量 
	int a, b; 
	// 初始化in_ver:
	for(int i = 0; i < MAX; i++){
		in_ver.push_back(vector<int>()); //添加一个空行 
	} 
	// 动态申请各个数组存储空间，长度为MAX：
	out_num = (int*)malloc(sizeof(int) * MAX);  
	memset(out_num, 0, MAX * sizeof(int));
	newPR = (float*)malloc(sizeof(float) * MAX);  
	oldPR = (float*)malloc(sizeof(float) * MAX);  
	// 输入有向边, 获取邻接矩阵 
	for(int i=0; fin>>a>>b; ++i){ 
		in_ver[b].push_back(a); //a->b, 记录每个点的入度点 
		out_num[a]++; // 记录每个点的出度数量 
		a = a < b ? b : a;
		n = n < a ? a : n; //最大值为顶点数量 
//		cout << a << " " << b << " " << n << endl;
	}
	n++; 
	// 进行计算 
	float *pr = computePR(n, ".\\output.txt");
	
	// 输出结果 
	for(int i = 0; i < n; i++){
		cout << pr[i] << " ";
	}
 
 	// 使用完需要回收:由编程人员手动申请，手动释放，若不手动释放，程序结束后由系统回收
	free(out_num);	//释放
	free(newPR);	//释放
	free(oldPR);	//释放
 	 
 	
	return 0;
}

2.python代码：

import numpy as np 
import networkx as nx
import matplotlib.pyplot as plt
'''
    实现PageRank算法
'''

def getM(A):
    '''
    由图的邻接矩阵来计算基本状态转移矩阵
    '''
    # 转置
    A = A.copy().T
    # 行求和
    sum1 = np.sum(A, axis=1)
    # 状态矩阵：
    M = []
    # 节点个数：
    n = len(A)
    for i in range(len(A)):
        # 处理Dead Ends
        if sum1[i] != 0: # 如果有出度
            M.append(A[i] / sum1[i])
        else:  # 没有出度,就假设都可以到任何点
            M.append(np.ones(n)/n)
    return np.array(M).T

def getBaseLev(n, d):
    '''
    计算完全随机转移矩阵: bl = (1-d)/n*1
    '''
    bl = np.ones(n) * (1-d)/n
#     bl = np.ones(n)
    return bl
    
if __name__ == '__main__':
    # 初始化参数
    d = 0.80 # 一个加权系数，通常取 0.85 左右,按照超链接进行浏览的概率，有效解决：Spider Traps(只链向自己)问题
    # 给定一个网页链接图的邻接矩阵，每一列表示一个网页的出度
    # eq1:
#     A =  np.array([[0,1,1,0,1,1,0],
#                    [1,0,1,1,0,0,0],
#                    [1,0,0,1,1,0,0],
#                    [1,0,0,0,1,0,0],
#                    [1,0,0,1,0,1,1],
#                    [0,0,0,0,1,0,0],
#                    [1,0,0,0,0,0,0]])
    # eq2:
    A = np.array([[0, 1, 0, 0],
                [1, 0, 0, 1],
                [1, 0, 1, 1],
                [1, 1, 0, 0]])
    # eq3:
#     A = np.array([[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, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 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, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 1, 0, 1, 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, 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, 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, 0, 0, 0, 0, 0, 0, 0, 0, 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, 1, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 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, 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, 0, 0, 0, 0, 0, 1, 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, 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, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 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], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]])
    # 通过邻接矩阵计算状态转移矩阵M:
    M = getM(A)
#     print(M)
    # 得到完全随机转移矩阵
    BL = getBaseLev(len(A), d)
    newPR = BL
    print('初始概率：\n', newPR)
    # 迭代终止阈值
    threshold = 1e-7
    # 进行迭代计算PR:
    cnt = 0
    while 1:
        cnt += 1
        oldPR = newPR
        # 根据公式计算：
        newPR = d * M @ oldPR + BL
        # 输出迭代结果：
        print('第 %d 次迭代：' % cnt)
        print(newPR)
        # 如果趋于稳定，则退出
        if np.linalg.norm(oldPR - newPR, ord=1) < threshold:
            break

3.c++拆分写：

文件：pagerank_c++.h

/*
Function: c++实现PageRank算法
Author: ystraw
Time: 2020-7-20 13:58
*/
#include<iostream>
#include<fstream>
#include<cmath>
#include <cstring>
using namespace std;
#define MAX 10000
int A[MAX][MAX];   // 关系矩阵 可以改为bool类型
float M[MAX][MAX]; // 权重矩阵
float newPR[MAX];  // 存储pr 
float oldPR[MAX];
float d = 0.85;      // 初始化参数
float threshold = 1e-4; // 迭代终止阈值

void getM(int n){
	// 计算初始权重矩阵 
	for(int i = 0; i < n; i++){
		int s = 0; // 统计出度 
		for(int j = 0; j < n; j++){
			s += A[j][i]; // 每列求和 
		}
		if(s == 0){  // 处理Dead Ends: 没有出度 
			for(int j = 0; j < n; j++){
				M[j][i] = 1.0 / n;
			}	
		}
		else{
			float fm = s;
			for(int j = 0; j < n; j++){
				M[j][i] = A[j][i] / fm;
			}	
		}
	}
//	for(int i = 0 ; i < n; i++){
//		for(int j = 0; j < n; j++){
//			cout << M[i][j] << " ";
//		}
//		cout << endl;
//	}
}

void getBaseLev(int n){
	//计算完全随机转移矩阵: bl = (1-d)/n*1
	for(int i = 0; i < n; i++){
		newPR[i] = 1.0/n;
	}
//	for(int i = 0; i < n; i++){
//		cout << newPR[i] << endl;
//	} 
}

void getPR(float *a, float *b, int n){
	// d * M @ oldPR + BL
	for(int i = 0; i < n; i++){
		a[i] = 0;
		for(int j = 0; j < n; j++){
			a[i] += M[i][j] * b[j];
		}
		a[i] = a[i] * d + 1.0/n * (1-d);
	}
//	for(int i = 0; i < n; i++){
//		cout << a[i] << " ";
//	} 
//	cout << endl;
} 

float* computePR(int n, char* outpath){ 
	
	//	cout << n << endl;
	// 通过邻接矩阵计算状态转移矩阵M:
    getM(n);
	// 得到完全随机转移矩阵
    getBaseLev(n);
    // 记录迭代次数
	int cnt = 0;
	while(1){
		cnt++;
		if(cnt % 2 == 1){
			getPR(oldPR, newPR, n); // newPR迭代存入oldPR 
		}
		else{
			getPR(newPR, oldPR, n); 
		}
		// 检测是否结束：
		float wc = 0;
		for(int i = 0; i < n; i++){
			wc += abs(oldPR[i] - newPR[i]);
		} 
//		cout << cnt << ": " << wc << " " << threshold << endl; 
		if(wc < threshold){
			break;
		}
	} 
	// 将结果写入 
	if(*outpath != '\0'){
//		ofstream fout(".\\output.txt");
		cout << "outpath:" << outpath << endl;
		ofstream fout(outpath);
		if(cnt % 2 == 1){
			for(int i=0; i<n; i++)
				fout<< i << ' ' << oldPR[i] << endl;
		}
		else{
			for(int i=0; i<n; i++)
				fout << i << ' ' << newPR[i] << endl;
		}
	}
	// 将结果返回
	if(cnt % 2 == 1){
		return oldPR;
	}
	else{
		return newPR;
	}
}

pagerankmain.cpp

/*
Function: c++实现PageRank算法
Author: ystraw
Time: 2020-7-20 13:58
*/
#include <iostream>
#include <fstream>
#include <cmath>
#include <cstring>
#include "pagerank_c++.h"
using namespace std;

int main(){
	// 读入/输出数据 
	ifstream fin(".\\edges.txt");
	int n = 0; // 顶点数据量 
	memset(A,0,sizeof(A));
	int a, b; 
	// 输入有向边, 获取邻接矩阵 
	for(int i=0; fin>>a>>b; ++i){ 
		A[b][a] = 1; //a->b
		a = a < b ? b : a;
		n = n < a ? a : n; //最大值为顶点数量 
//		cout << a << " " << b << " " << n << endl;
	}
	n++; 
	// 进行计算 
	float *pr = computePR(n, ".\\output.txt");
	for(int i = 0; i <= n; i++){
		cout << pr[i] << " ";
	}
	return 0;
}

测试样例

edges.txt:
0 2
0 1
0 3
1 0
1 3
2 0
3 1
3 2
output.txt(除n):
0 0.33333
1 0.22222
2 0.22222
3 0.22222

edges.txt:
0 2
0 1
0 3
1 0
1 3
2 2
3 1
3 2
output.txt:
0 0.101351
1 0.128378
2 0.641892
3 0.128378

Pagerank收敛性验证

step = 100
p = 1
d = 0.85
sum = 0
for i in range(step):
  sum += p * (1 - d)
  p = p * d
  print('step=%2d' % i, 'p=%.10f' % p, 'sum=%.10f' % sum, 'wc=%.10f' % (1 - sum))

output:
运行85轮，和1的误差将小于1e-6.

step= 0 p=0.8500000000 sum=0.1500000000 wc=0.8500000000
step= 1 p=0.7225000000 sum=0.2775000000 wc=0.7225000000
step= 2 p=0.6141250000 sum=0.3858750000 wc=0.6141250000
step= 3 p=0.5220062500 sum=0.4779937500 wc=0.5220062500
step= 4 p=0.4437053125 sum=0.5562946875 wc=0.4437053125
step= 5 p=0.3771495156 sum=0.6228504844 wc=0.3771495156
step= 6 p=0.3205770883 sum=0.6794229117 wc=0.3205770883
step= 7 p=0.2724905250 sum=0.7275094750 wc=0.2724905250
step= 8 p=0.2316169463 sum=0.7683830537 wc=0.2316169463
step= 9 p=0.1968744043 sum=0.8031255957 wc=0.1968744043
step=10 p=0.1673432437 sum=0.8326567563 wc=0.1673432437
step=11 p=0.1422417571 sum=0.8577582429 wc=0.1422417571
step=12 p=0.1209054936 sum=0.8790945064 wc=0.1209054936
step=13 p=0.1027696695 sum=0.8972303305 wc=0.1027696695
step=14 p=0.0873542191 sum=0.9126457809 wc=0.0873542191
step=15 p=0.0742510862 sum=0.9257489138 wc=0.0742510862
step=16 p=0.0631134233 sum=0.9368865767 wc=0.0631134233
step=17 p=0.0536464098 sum=0.9463535902 wc=0.0536464098
step=18 p=0.0455994483 sum=0.9544005517 wc=0.0455994483
step=19 p=0.0387595311 sum=0.9612404689 wc=0.0387595311
step=20 p=0.0329456014 sum=0.9670543986 wc=0.0329456014
step=21 p=0.0280037612 sum=0.9719962388 wc=0.0280037612
step=22 p=0.0238031970 sum=0.9761968030 wc=0.0238031970
step=23 p=0.0202327175 sum=0.9797672825 wc=0.0202327175
step=24 p=0.0171978099 sum=0.9828021901 wc=0.0171978099
step=25 p=0.0146181384 sum=0.9853818616 wc=0.0146181384
step=26 p=0.0124254176 sum=0.9875745824 wc=0.0124254176
step=27 p=0.0105616050 sum=0.9894383950 wc=0.0105616050
step=28 p=0.0089773642 sum=0.9910226358 wc=0.0089773642
step=29 p=0.0076307596 sum=0.9923692404 wc=0.0076307596
step=30 p=0.0064861457 sum=0.9935138543 wc=0.0064861457
step=31 p=0.0055132238 sum=0.9944867762 wc=0.0055132238
step=32 p=0.0046862402 sum=0.9953137598 wc=0.0046862402
step=33 p=0.0039833042 sum=0.9960166958 wc=0.0039833042
step=34 p=0.0033858086 sum=0.9966141914 wc=0.0033858086
step=35 p=0.0028779373 sum=0.9971220627 wc=0.0028779373
step=36 p=0.0024462467 sum=0.9975537533 wc=0.0024462467
step=37 p=0.0020793097 sum=0.9979206903 wc=0.0020793097
step=38 p=0.0017674132 sum=0.9982325868 wc=0.0017674132
step=39 p=0.0015023012 sum=0.9984976988 wc=0.0015023012
step=40 p=0.0012769561 sum=0.9987230439 wc=0.0012769561
step=41 p=0.0010854127 sum=0.9989145873 wc=0.0010854127
step=42 p=0.0009226008 sum=0.9990773992 wc=0.0009226008
step=43 p=0.0007842106 sum=0.9992157894 wc=0.0007842106
step=44 p=0.0006665790 sum=0.9993334210 wc=0.0006665790
step=45 p=0.0005665922 sum=0.9994334078 wc=0.0005665922
step=46 p=0.0004816034 sum=0.9995183966 wc=0.0004816034
step=47 p=0.0004093629 sum=0.9995906371 wc=0.0004093629
step=48 p=0.0003479584 sum=0.9996520416 wc=0.0003479584
step=49 p=0.0002957647 sum=0.9997042353 wc=0.0002957647
step=50 p=0.0002514000 sum=0.9997486000 wc=0.0002514000
step=51 p=0.0002136900 sum=0.9997863100 wc=0.0002136900
step=52 p=0.0001816365 sum=0.9998183635 wc=0.0001816365
step=53 p=0.0001543910 sum=0.9998456090 wc=0.0001543910
step=54 p=0.0001312324 sum=0.9998687676 wc=0.0001312324
step=55 p=0.0001115475 sum=0.9998884525 wc=0.0001115475
step=56 p=0.0000948154 sum=0.9999051846 wc=0.0000948154
step=57 p=0.0000805931 sum=0.9999194069 wc=0.0000805931
step=58 p=0.0000685041 sum=0.9999314959 wc=0.0000685041
step=59 p=0.0000582285 sum=0.9999417715 wc=0.0000582285
step=60 p=0.0000494942 sum=0.9999505058 wc=0.0000494942
step=61 p=0.0000420701 sum=0.9999579299 wc=0.0000420701
step=62 p=0.0000357596 sum=0.9999642404 wc=0.0000357596
step=63 p=0.0000303956 sum=0.9999696044 wc=0.0000303956
step=64 p=0.0000258363 sum=0.9999741637 wc=0.0000258363
step=65 p=0.0000219608 sum=0.9999780392 wc=0.0000219608
step=66 p=0.0000186667 sum=0.9999813333 wc=0.0000186667
step=67 p=0.0000158667 sum=0.9999841333 wc=0.0000158667
step=68 p=0.0000134867 sum=0.9999865133 wc=0.0000134867
step=69 p=0.0000114637 sum=0.9999885363 wc=0.0000114637
step=70 p=0.0000097441 sum=0.9999902559 wc=0.0000097441
step=71 p=0.0000082825 sum=0.9999917175 wc=0.0000082825
step=72 p=0.0000070401 sum=0.9999929599 wc=0.0000070401
step=73 p=0.0000059841 sum=0.9999940159 wc=0.0000059841
step=74 p=0.0000050865 sum=0.9999949135 wc=0.0000050865
step=75 p=0.0000043235 sum=0.9999956765 wc=0.0000043235
step=76 p=0.0000036750 sum=0.9999963250 wc=0.0000036750
step=77 p=0.0000031237 sum=0.9999968763 wc=0.0000031237
step=78 p=0.0000026552 sum=0.9999973448 wc=0.0000026552
step=79 p=0.0000022569 sum=0.9999977431 wc=0.0000022569
step=80 p=0.0000019184 sum=0.9999980816 wc=0.0000019184
step=81 p=0.0000016306 sum=0.9999983694 wc=0.0000016306
step=82 p=0.0000013860 sum=0.9999986140 wc=0.0000013860
step=83 p=0.0000011781 sum=0.9999988219 wc=0.0000011781
step=84 p=0.0000010014 sum=0.9999989986 wc=0.0000010014
step=85 p=0.0000008512 sum=0.9999991488 wc=0.0000008512
step=86 p=0.0000007235 sum=0.9999992765 wc=0.0000007235
step=87 p=0.0000006150 sum=0.9999993850 wc=0.0000006150
step=88 p=0.0000005227 sum=0.9999994773 wc=0.0000005227
step=89 p=0.0000004443 sum=0.9999995557 wc=0.0000004443
step=90 p=0.0000003777 sum=0.9999996223 wc=0.0000003777
step=91 p=0.0000003210 sum=0.9999996790 wc=0.0000003210
step=92 p=0.0000002729 sum=0.9999997271 wc=0.0000002729
step=93 p=0.0000002319 sum=0.9999997681 wc=0.0000002319
step=94 p=0.0000001972 sum=0.9999998028 wc=0.0000001972
step=95 p=0.0000001676 sum=0.9999998324 wc=0.0000001676
step=96 p=0.0000001424 sum=0.9999998576 wc=0.0000001424
step=97 p=0.0000001211 sum=0.9999998789 wc=0.0000001211
step=98 p=0.0000001029 sum=0.9999998971 wc=0.0000001029
step=99 p=0.0000000875 sum=0.9999999125 wc=0.0000000875