Graph-Maker-Other
view release on metacpan or search on metacpan
devel/wiener.pl view on Meta::CPAN
#!/usr/bin/perl -w
# Copyright 2015, 2016, 2017, 2021 Kevin Ryde
#
# This file is part of Graph-Maker-Other.
#
# Graph-Maker-Other is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by the
# Free Software Foundation; either version 3, or (at your option) any later
# version.
#
# Graph-Maker-Other is distributed in the hope that it will be useful, but
# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
# or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
# for more details.
#
# You should have received a copy of the GNU General Public License along
# with Graph-Maker-Other. If not, see <http://www.gnu.org/licenses/>.
use strict;
use 5.010;
use File::Slurp;
use List::Util 'min','max';
use POSIX 'ceil';
use FindBin;
use lib "$FindBin::Bin/../devel/lib";
use MyGraphs;
# GP-DEFINE nearly_equal(x,y,delta=1-e10) = abs(x-y) < delta;
$|=1;
# uncomment this to run the ### lines
# use Smart::Comments;
{
# average_path_length() / diameter of various graph types
#
#
#--------------------
# linear
# *--*--*--*
# n=4; sum(i=1,n-1, i) == n*(n-1)/2
# 2*sum(n=2,N, n*(n-1)/2) = (N-1)*(N)*(N+1)/3 \\ 2*tetrahedral
# (N-1)*(N)*(N+1)/3 / (N-1) / N^2 -> 1/3
#
# GP-DEFINE linear_diameter(k) = k-1;
# GP-DEFINE linear_W_by_sum(k) = sum(i=1,k, sum(j=i+1,k, j-i));
# vector(20,k,k--; linear_W_by_sum(k))
# GP-DEFINE linear_W(k) = 1/6*k*(k^2-1);
# \\ A000292 tetrahedral n*(n+1)*(n+2)/6 = (n+1)*((n+1)^2-1)/6
# GP-Test vector(100,k,k--; linear_W(k)) == vector(100,k,k--; linear_W_by_sum(k))
# GP-DEFINE linear_binomial(k) = k*(k-1)/2;
# vector(100,k,k--; linear_binomial(k)) == vector(100,k,k--; binomial(k,2))
# GP-DEFINE linear_mean_by_formula(k) = linear_W(k) / (linear_diameter(k) * binomial(k,2));
# linear_mean(k) = 1/6*k*(k^2-1) / (k*(k-1)/2) / (k-1);
# linear_mean(k) = 1/3*(k+1)/(k-1);
# linear_mean(k) = my(n=k-1); (n+2)/(3*n);
# linear_mean(k) = my(n=k+1); n/(3*(n-2));
# GP-DEFINE linear_mean(k) = 1/3 + 2/3/(k-1);
# GP-Test vector(100,k,k++; linear_mean(k)) == vector(100,k,k++; linear_mean_by_formula(k))
# GP-Test my(k=1000000); nearly_equal(linear_mean(k), 1/3, 1e-3) \\ -> 1/3
# GP-DEFINE A060789(n) = n/gcd(n,2)/gcd(n,3); \\ numerator
# GP-Test vector(100,k,k++; numerator(linear_mean(k))) == vector(100,k,k++; A060789(k+1))
# vector(30,k,k++; denominator(linear_mean(k))) \\ *1.0
# can divide out possible 2 if n even, possible 3 too,
# otherwise n,n-1 no common factor
# denominator not
#
#--------------------
# star of k vertices
# *
# |
# *--*--*
# GP-DEFINE star_W_by_terms(k) = if(k==0,0, 2*binomial(k-1,2) + k-1);
# star_W(0) == 0
# star_W(4) == 1+2+2 + 1+1 + 2
# GP-DEFINE star_W(k) = if(k==0,0, (k-1)^2);
# GP-Test vector(100,k,k--; star_W(k)) == vector(100,k,k--; star_W_by_terms(k))
# GP-DEFINE star_diameter(k) = if(k<=1,0, k==2,1, 2);
# GP-DEFINE star_mean_by_formula(k) = star_W(k) / star_diameter(k) / binomial(k,2);
# GP-DEFINE star_mean(k) = if(k==2,1, 1 - 1/k); \\ -> 1
# GP-Test vector(100,k,k++; star_mean(k)) == vector(100,k,k++; star_mean_by_formula(k))
# vector(20,k,k++; star_mean(k))
#
#----------
# Graph::Maker::BalancedTree, rows 0 to n inclusive
# A158681 Wiener of complete binary tree = 4*(n-2)*4^n + 2*(n+4)*2^n
# = 4, 48, 368, 2304 starting n=1
# * (1+2 + 2 + 1+2)/2 = 4
# / \ has 2^(n+1)-1 vertices
# * *
# diameter 2*n
# GP-DEFINE complete_binary_W(n) = 4*(n-2)*4^n + 2*(n+4)*2^n;
# GP-Test complete_binary_W(1) == 4
# GP-Test complete_binary_W(2) == 48
# GP-DEFINE complete_binary_N(n) = 2^(n+1)-1;
# GP-Test complete_binary_N(1) == 3
# GP-DEFINE complete_binary_diameter(n) = 2*n;
# GP-Test complete_binary_diameter(1) == 2
# GP-DEFINE complete_binary_mean(n) = 2*complete_binary_W(n) \
# GP-DEFINE / complete_binary_N(n)^2 \
# GP-DEFINE / complete_binary_diameter(n);
# GP-DEFINE complete_binary_mean_formula(n) = \
# GP-DEFINE 1 - 2/n + 3*1/(2*2^n - 1) + 2*(1 + 1/n)*1/(2*2^n - 1)^2;
# GP-Test vector(100,n, complete_binary_mean_formula(n)) == vector(100,n, complete_binary_mean(n))
# complete_binary_mean -> 1 slightly slowly
#
#----------
# BinomialTree
# A192021 Wiener of binomial tree = 2*n*4^n + 2^n
# W Binomial tree = (k-1)*2^(2k-1) + 2^(k-1)
# (in POD)
#
#----------
# Cycle
# GP-DEFINE cycle_W_by_sum(k) = sum(i=1,k, sum(j=i+1,k, min(j-i, (i-j)%k)));
# cycle_W_by_sum(0) == 0
# cycle_W_by_sum(1) == 0
# cycle_W_by_sum(2) == 1
# cycle_W_by_sum(3) == 3
# cycle_W_by_sum(4) == 1+1+2 + 1+2 + 1
# vector(20,k,k--; cycle_W_by_sum(k))
# GP-DEFINE cycle_W(k) = 1/2 * k * floor((k/2)^2);
# GP-DEFINE cycle_W(k) = 1/2*k*((k/2)^2 - if(k%2==1,1/4));
# GP-Test vector(100,k,k--; cycle_W(k)) == vector(100,k,k--; cycle_W_by_sum(k))
# GP-DEFINE cycle_diameter(k) = floor(k/2);
# GP-DEFINE cycle_mean_by_formula(k) = cycle_W(k) / cycle_diameter(k) / binomial(k,2);
# GP-DEFINE cycle_mean(k) = 1/2 + (k^2/4 - 1/2*(k-1)*(k/2-if(k%2==1,1/2)) - if(k%2==1,1/4)) / (k-1) / floor(k/2);
# GP-DEFINE cycle_mean(k) = 1/2 + (k^2/4 - (k-1)*(k/4-if(k%2==1,1/4)) - if(k%2==1,1/4)) / (k-1) / floor(k/2);
# GP-DEFINE cycle_mean(k) = 1/2 + 1/4*(k*if(k%2==0,1,2) - if(k%2==0,0,1) - if(k%2==0,0,1)) / (k-1) / floor(k/2);
# GP-DEFINE cycle_mean(k) = 1/2 + 1/4*(if(k%2==0,k,2*k-2)) / (k-1) / if(k%2==0,k/2,(k-1)/2);
# GP-DEFINE cycle_mean(k) = 1/2 + 1/2*if(k%2==0,k / (k-1) / if(k%2==0,k,k-1),(2*k-2) / (k-1) / if(k%2==0,k,k-1));
# GP-DEFINE cycle_mean(k) = 1/2 + 1/2*if(k%2==0, 1/k + 1/(k*(k-1)), 2/(k-1));
# GP-Test vector(100,k,k++; cycle_mean(k)) == vector(100,k,k++; cycle_mean_by_formula(k))
# GP-Test my(k=1000000); nearly_equal(cycle_mean(k), 1/2, 1e-3) \\ -> 1/2
# vector(30,k,k++; numerator(cycle_mean(k)))
# vector(30,k,k++; denominator(cycle_mean(k)))
#
#----------
# Hypercube
# GP-DEFINE hypercube_W_by_sum(k) = sum(i=0,2^k-1, sum(j=i+1,2^k-1, hammingweight(bitxor(i,j))));
# hypercube_W_by_sum(0) == 0
# hypercube_W_by_sum(1) == 1
# hypercube_W_by_sum(2) == 8
# vector(8,k,k--; hypercube_W_by_sum(k))
#
# GP-DEFINE hypercube_W(k) = k*4^(k-1);
# GP-Test vector(10,k,k--; hypercube_W(k)) == \
# GP-Test vector(10,k,k--; hypercube_W_by_sum(k))
#
# GP-DEFINE hypercube_diameter(k) = k;
# GP-DEFINE hypercube_mean_by_formula(k) = \
# GP-DEFINE hypercube_W(k) / hypercube_diameter(k) / binomial(2^k-1,2);
# GP-DEFINE hypercube_mean(k) = 1/2 * 4^k / (2^k-1)/(2^k-2);
# GP-Test vector(20,k,k++; hypercube_mean_by_formula(k)) == \
# GP-Test vector(20,k,k++; hypercube_mean(k))
# GP-Test vector(20,k,k++; hypercube_mean_by_formula(k)) == \
# GP-Test vector(20,k,k++; 1/2 + 1/2*(3*2^k - 2) / (2^k-1)/(2^k-2))
# GP-Test vector(20,k,k++; hypercube_mean_by_formula(k)) == \
# GP-Test vector(20,k,k++; 1/2 + 1/2*(3*2^k-3 +3 - 2) / (2^k-1)/(2^k-2))
# GP-Test vector(20,k,k++; hypercube_mean_by_formula(k)) == \
# GP-Test vector(20,k,k++; 1/2 + 3/2/(2^k-2) + 1/2/(2^k-1)/(2^k-2))
# GP-Test my(k=100); nearly_equal(hypercube_mean(k), 1/2, 1e-3) \\ -> 1/2
#
# vector(30,k,k++; numerator(hypercube_mean(k))) \\ 4^k
# vector(30,k,k++; denominator(hypercube_mean(k))) \\ binomial(2^k-1,2)
# vector(10,k,k++; hypercube_mean(k)*1.0)
require Graph::Maker::Star;
require Graph::Maker::Linear;
my @values;
foreach my $k (1..4) {
# my $graph = Graph::Maker->new('star', N=>$k, undirected=>1);
# my $graph = Graph::Maker->new('linear', N=>$k, undirected=>1);
# require Graph::Maker::BalancedTree;
# my $graph = Graph::Maker->new('balanced_tree',
# fan_out => 2, height => $k,
# );
# require Graph::Maker::BinomialTree;
# my $graph = Graph::Maker->new('binomial_tree',
# order => $k,
# undirected => 1,
# );
# require Graph::Maker::Cycle;
# my $graph = Graph::Maker->new('cycle', N => $k);
# require Graph::Maker::Hypercube;
# my $graph = Graph::Maker->new('hypercube', N => $k);
require Graph::Maker::Keller;
my $graph = Graph::Maker->new('Keller', N => $k);
# require Graph::Maker::FibonacciTree;
# my $graph = Graph::Maker->new('fibonacci_tree',
# height => $k,
# # leaf_reduced => 1,
# # series_reduced => 1,
# );
# my $total = Graph_total_path_length($graph);
my $W = MyGraphs::Graph_Wiener_index($graph);
my $vertices = $graph->vertices;
my $diameter = $graph->diameter || 0;
my $div = $vertices*($vertices-1)/2;
my $average = ($div == 0 || $diameter == 0 ? 'undef'
: $W / $div / $diameter);
# my $average = ($diameter == 0 ? 'undef'
# : $total / $diameter / $vertices**2);
if ($k==1) { print $graph->get_graph_attribute ('name'),"\n"; }
print "$k W=$W diam=$diameter $average\n";
push @values, $W;
}
require Math::OEIS::Grep;
Math::OEIS::Grep->search(array => \@values, verbose=>1);
exit 0;
sub Graph_total_path_length {
my ($graph) = @_;
( run in 3.679 seconds using v1.01-cache-2.11-cpan-5837b0d9d2c )