mersenneforum.org  

Go Back   mersenneforum.org > Great Internet Mersenne Prime Search > Math

Reply
 
Thread Tools
Old 2009-11-10, 02:18   #1
Dougy
 
Dougy's Avatar
 
Aug 2004
Melbourne, Australia

23×19 Posts
Default Counting Latin rectangles

Hi folks,

As some may know already, I'm about to complete my thesis on the enumeration of Latin rectangles. Here's a link discussing some of my work.

Currently I'm thinking about future projects after I submit (in a few weeks). One such idea is extending the tables of R(k,n) (counting k x n reduced Latin rectangles) for small k. Specifically, I have in mind conjecture about divisors of these numbers and more data for R(5,n) would support (or, if I'm wrong, not support) this conjecture.

I've written c code and uploaded it here. It implements Doyle's formula in the cases R(4,n), R(5,n) and R(6,n), which I have used to find previously unknown values of these numbers. It's not very optimised and I'm wondering if anyone here would be able and willing to help out.

It uses GMP, which I'm told is pretty slow. It doesn't use anything too fancy, eg. multiplication, division, factorials. Since I'm only interested in prime-power divisors of these numbers, it might be possible to only run this code modulo 2^N and then again for 3^N, hopefully eliminating the need for GMP (although I'd prefer to know the numbers exactly if possible).

Anyway, I'm just looking for expressions of interest and an idea of what's possible at the moment... I know there's a lot of people here with expertise in coding.
Dougy is offline   Reply With Quote
Old 2009-11-10, 06:22   #2
Primeinator
 
Primeinator's Avatar
 
"Kyle"
Feb 2005
Somewhere near M50..sshh!

37E16 Posts
Default

Well, I don't know anything about coding, nor about your topic of expertise. However, I would like to congratulate you in advance for completing your dissertation and doctorate!
Primeinator is offline   Reply With Quote
Old 2009-11-12, 07:29   #3
Dougy
 
Dougy's Avatar
 
Aug 2004
Melbourne, Australia

23·19 Posts
Default

Thanks (:. They say the thesis will be available online, so I can post a link.

Well not much of a race to respond... I'm not sure if it's because it's quite off-topic for this forum. The code, however, should be able to be improved without knowledge of Latin rectangles. Perhaps I'll try to explain how it works a bit better. I'm mostly interested in R(5,n) since for R(4,n) I have lots of data already and R(6,n) seems too hard.

Overall, it's a sum over "partitions" of n=A[0]+A[1]+...+A[15] (where the A[i] can take the values of 0...n). For each partition of n the variable "internal" is added to "count."

Code:
#include <stdio.h>
#include <gmp.h>
#include <time.h>

int n,x,y;
time_t start_time;
mpz_t count,internal,temp;
int A[16],B[16];

int DoyleF5_1(int A[16]) { return A[1]+A[2]+A[3]+A[4]+A[5]+A[6]+A[7]+A[8]; }
int DoyleF5_12(int A[16]) { return A[1]+A[2]+A[3]+A[4]; }
int DoyleF5_123(int A[16]) { return A[1]+A[2]; }
int DoyleF5_1234(int A[16]) { return A[1]; }
int DoyleF5_124(int A[16]) { return A[1]+A[3]; }
int DoyleF5_13(int A[16]) { return A[1]+A[2]+A[5]+A[6]; }
int DoyleF5_134(int A[16]) { return A[1]+A[5]; }
int DoyleF5_14(int A[16]) { return A[1]+A[3]+A[5]+A[7]; }
int DoyleF5_2(int A[16]) { return A[1]+A[2]+A[3]+A[4]+A[9]+A[10]+A[11]+A[12]; }
int DoyleF5_23(int A[16]) { return A[1]+A[2]+A[9]+A[10]; }
int DoyleF5_234(int A[16]) { return A[1]+A[9]; }
int DoyleF5_24(int A[16]) { return A[1]+A[3]+A[9]+A[11]; }
int DoyleF5_3(int A[16]) { return A[1]+A[2]+A[5]+A[6]+A[9]+A[10]+A[13]+A[14]; }
int DoyleF5_34(int A[16]) { return A[1]+A[5]+A[9]+A[13]; }
int DoyleF5_4(int A[16]) { return A[1]+A[3]+A[5]+A[7]+A[9]+A[11]+A[13]+A[15]; }
After the front matter, the above is required for Doyle's formula. In each case, it's just addition of some collection of the A[i]'s.

Code:
int DoyleG5(int A[16]) { return DoyleF5_1(A)*DoyleF5_2(A)*DoyleF5_3(A)*DoyleF5_4(A)-DoyleF5_1(A)*DoyleF5_2(A)*DoyleF5_34(A)-DoyleF5_1(A)*DoyleF5_23(A)*DoyleF5_4(A)+2*DoyleF5_1(A)*DoyleF5_234(A)-DoyleF5_1(A)*DoyleF5_24(A)*DoyleF5_3(A)-DoyleF5_12(A)*DoyleF5_3(A)*DoyleF5_4(A)+DoyleF5_12(A)*DoyleF5_34(A)+2*DoyleF5_123(A)*DoyleF5_4(A)-6*DoyleF5_1234(A)+2*DoyleF5_124(A)*DoyleF5_3(A)-DoyleF5_13(A)*DoyleF5_2(A)*DoyleF5_4(A)+DoyleF5_13(A)*DoyleF5_24(A)+2*DoyleF5_134(A)*DoyleF5_2(A)-DoyleF5_14(A)*DoyleF5_2(A)*DoyleF5_3(A)+DoyleF5_14(A)*DoyleF5_23(A); }
This is the guts. But it's just multiplication, addition and subtraction of the earlier functions.

Code:
int MinusOneOrNot(int A[16]) { if((A[2]+A[3]+A[5]+A[8]+A[9]+A[12]+A[14]+A[15])%2==0) { return 1; } else { return -1; } }
This will later decide the sign of "internal."

Code:
  start_time=time(NULL);
  mpz_init(count);
  mpz_init(temp);
  mpz_init(internal);
  for(n=20;n<=100;n++) { mpz_set_si(count,0);
Initialisation... I've set it to run for n=20 to 100.

Code:
  for(A[1]=0;A[1]<=n;A[1]++) {
  for(A[2]=0;A[2]<=n-A[1];A[2]++) {
  for(A[3]=0;A[3]<=n-A[1]-A[2];A[3]++) {
  for(A[4]=0;A[4]<=n-A[1]-A[2]-A[3];A[4]++) {
  for(A[5]=0;A[5]<=n-A[1]-A[2]-A[3]-A[4];A[5]++) {
  for(A[6]=0;A[6]<=n-A[1]-A[2]-A[3]-A[4]-A[5];A[6]++) {
  for(A[7]=0;A[7]<=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6];A[7]++) {
  for(A[8]=0;A[8]<=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6]-A[7];A[8]++) {
  for(A[9]=0;A[9]<=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6]-A[7]-A[8];A[9]++) {
  for(A[10]=0;A[10]<=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6]-A[7]-A[8]-A[9];A[10]++) {
  for(A[11]=0;A[11]<=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6]-A[7]-A[8]-A[9]-A[10];A[11]++) {
  for(A[12]=0;A[12]<=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6]-A[7]-A[8]-A[9]-A[10]-A[11];A[12]++) {
  for(A[13]=0;A[13]<=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6]-A[7]-A[8]-A[9]-A[10]-A[11]-A[12];A[13]++) {
  for(A[14]=0;A[14]<=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6]-A[7]-A[8]-A[9]-A[10]-A[11]-A[12]-A[13];A[14]++) {
  for(A[15]=0;A[15]<=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6]-A[7]-A[8]-A[9]-A[10]-A[11]-A[12]-A[13]-A[14];A[15]++) {
    A[0]=n-A[1]-A[2]-A[3]-A[4]-A[5]-A[6]-A[7]-A[8]-A[9]-A[10]-A[11]-A[12]-A[13]-A[14]-A[15];
Chooses the partition of n=A[0]+A[1]+...+A[15]. Is there a better way?

Code:
    mpz_fac_ui(internal,n);
    for(x=0;x<=15;x++) { mpz_fac_ui(temp,A[x]); mpz_div(internal,internal,temp); }
    mpz_mul_si(internal,internal,MinusOneOrNot(A));
Initialising "internal" to be the multinomial coefficient {n \choose A[0],A[1],\dots,A[15]}=n!/(A[0]! \cdot A[1]! \cdots A[15]!). Then MinusOneOrNot(A) determines the appropriate sign.

Code:
    for(x=0;x<=15;x++) {
      for(y=0;y<=15;y++) { B[y]=A[y]; }
      B[x]--;
      mpz_set_si(temp,0);
      mpz_set_si(temp,DoyleG5(B));
      mpz_pow_ui(temp,temp,A[x]);
      mpz_mul(internal,internal,temp);
    }
This is the product required for Doyle's formula. for each 0<=x<=15 the B[i]=A[i] except B[x] is reduced by one. "internal" is then created via a "temp" variable which is \prod_{0 \leq x \leq 15} DoyleG5(B)^{A[x]} as a result of this loop.

Code:
    }
    mpz_add(count,count,internal);
  } } } } } } } } } } } } } } }
  if(n>=5) { mpz_div_ui(count,count,(n-1)*(n-2)*(n-3)*(n-4)); }
  gmp_printf("R(5,%d)=%Zd   %ld sec.\n",n,count,time(NULL)-start_time);
  }
}
Adds "internal" to the variable "count" for each partition. At the end it divides by (n-1)(n-2)(n-3)(n-4) when it can and prints out the result.

Hope this helps. But if there's anything I could explain more clearly, I'm happy to.
Dougy is offline   Reply With Quote
Old 2010-02-16, 10:20   #4
Dougy
 
Dougy's Avatar
 
Aug 2004
Melbourne, Australia

23·19 Posts
Default

I'm surprised how little interest there is in this... ah well. Anyway, my thesis has now been accepted so I will have my Ph.D. soon (after some admin hurdles). I uploaded the thesis to this website.
Dougy is offline   Reply With Quote
Reply

Thread Tools


Similar Threads
Thread Thread Starter Forum Replies Last Post
Counting ECM work twice alpertron PrimeNet 1 2009-12-01 23:33
Counting duplicates fivemack Factoring 1 2009-08-24 02:26
Does anyone know latin well? crash893 Lounge 3 2008-04-23 23:17
Prime squares/rectangles roger Puzzles 10 2007-05-04 16:07
counting bricks michael Puzzles 8 2004-01-14 16:27

All times are UTC. The time now is 22:51.

Sat Oct 31 22:51:48 UTC 2020 up 51 days, 20:02, 2 users, load averages: 2.30, 2.41, 2.36

Powered by vBulletin® Version 3.8.11
Copyright ©2000 - 2020, Jelsoft Enterprises Ltd.

This forum has received and complied with 0 (zero) government requests for information.

Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.2 or any later version published by the Free Software Foundation.
A copy of the license is included in the FAQ.