[Vari] Contest 8: Compattamento Barionico

[gugoXX]

1) Sia data una lista di punti contenuti in uno spazio cubico, con coordinate tridimensionali in virgola mobile.
Il formato della lista e'

NP   Size
X1 Y1 Z1
X2 Y2 Z2
...
Xn Yn Zn

Dove NP e' il numero di punti (quindi le linee che seguono)
e Size e' un numero floating point pari alla dimensione del cubo di contenimento (quindi nessuna coordinata superera' tale valore)

Trovare ed enunciare la coppia di punti a distanza minima e la coppia di punti a distanza massima (e le distanze relative)

PS: Per quanto riguarda la formula della distanza tridimensionale in uno spazio euclideo, dati i punti
Pa di coordinate (Xa, Ya, Za) e Pb di coordinate (Xb, Yb, Zb)


(Da Wiki, sotto distanza Euclidea: )
E non ditemi che vi arenate sulla distanza 3D per favore. E' pieno di gente qui che vorrebbe fare giochi 3D, e chi si arena su una semplice distanza...

Lista dei punti: http://www.usaupload.net/d/cjepq9exi1j

2) Teorema: Sia dato un cubo di lato d ed un numero di particelle N
Esiste almeno una configurazione spaziale casuale delle N particelle tale per cui non esista nemmeno una coppia di particelle
con distanza c inferiore a
d / ³√N




La semplice dimostrazione è lasciata come esercizio al lettore (Ho sempre sognato di scriverlo)

In pratica questo teorema ci assicura che tirando a caso molte volte la posizione di ciascuna delle N particelle, c'e' almeno una configurazione tale per cui le particelle non sono troppo vicine e distano tra loro almeno c.

Problema:
Si trovi almeno una configurazione spaziale (casuale) delle particelle nei seguenti casi:
d = 10 N = 5
d = 10 N = 10
d = 10 N = 50
d = 10 N = 100
d = 10 N = 200

PS: La soluzione di questo problema è in uso al CERN, per le condizioni iniziali di alcune simulazioni.

Domani posterò anche l'algoritmo di confronto a forza bruta.
Il primo esercizio serve come test di verifica per il secondo, non mi aspetto un'ottimizzazione particolare, anche se magari potrebbe essere interessante.

[banryu79]

Originariamente inviato da gugoXX

Il primo esercizio serve come test di verifica per il secondo, non mi aspetto un'ottimizzazione particolare, anche se magari potrebbe essere interessante.

Dato che quello che interessa non è la reale distanza 3D ma solo il confronto di valori riguardanti le varie distanze dei vari punti per trovare la coppia con la minima distanza e la coppia con la massima distanza, nella formula per calcolare la distanza di un punto non dovrebbe essere neccessario fare la radice quadrata, si può così evitarne la sua chiamata.
(Sì lo so che parlavi di ottimizzazioni a livello algoritmico, però mi è saltata all'occhio questa cosa).

[gugoXX]

Dai cos'e'? Passiamo dalle cose quadrate a quelle sferiche e ci "Impalliamo" subito?

[banryu79]

Mah, non saprei, non è che sia perchè il Forum è un po' smorto in quest'ultimo periodo?
Magari gioverebbe anche capire meglio che cosa si sia alla fine deciso/accolto nella discussione circa i nuovi contest... link

[Tommo]

Non capisco... ma per "casuale" si intende che va bene una configurazione qualsiasi, oppure che i valori iniziali devono essere tirati a caso?
C'è una bella differenza...

Nel primo caso la soluzione più veloce è:
nextPos = pos + c;
Sarebbe ben facile

Nel secondo invece c'è un problema... la configurazione casuale migliore potrebbe uscire alla seconda iterazione come dopo 1000... rendendo difatti insensato il contest...
credo sia il caso di fare un file che contenga le varie configurazioni iniziali.

[Vincenzo1968]

Ohé,

per me che, in matematica, sono deboluccio: è corretto il seguente output?:



Questo è il codice con l'array di prova:

#include <stdio.h>
#include <math.h>

typedef struct tagPunti
{
	double A;
	double B;
	double C;
} Punti;

typedef struct tagRes
{
	Punti P1;
	Punti P2;
	double Distanza;
	int index1;
	int index2;
} Res;

void Calcola(Punti P[], int dim, Res *pMin, Res *pMax)
{
	int x, y;
	Punti P1, P2;
	double distanza, distanzaMin, distanzaMax;

	distanzaMin = distanzaMax = -1;

	x = y = 0;

	while( x < dim - 1 )
	{
		P1 = P[x];
		y = x + 1;
		while( y < dim )
		{
			P2 = P[y];
			distanza = sqrt( pow(P1.A - P2.A, 2.0) + pow(P1.B - P2.B, 2.0) + pow(P1.C - P2.C, 2.0) );
			//printf("Distanza -> %.15lf\n", distanza);
			if ( distanza < distanzaMin || distanzaMin == -1 )
			{
				distanzaMin = distanza;
				pMin->P1 = P[x];
				pMin->P2 = P[y];
				pMin->Distanza = distanza;
				pMin->index1 = x;
				pMin->index2 = y;
			}
			if ( distanza > distanzaMax )
			{
				distanzaMax = distanza;
				pMax->P1 = P[x];
				pMax->P2 = P[y];
				pMax->Distanza = distanza;
				pMax->index1 = x;
				pMax->index2 = y;
			}
			++y;
		}
		++x;
	}
}

int main()
{
	Punti P[5];
	Res r1, r2;

	P[0].A = 1.5;
	P[0].B = 2.3;
	P[0].C = 5.1;

	P[1].A = 2.8;
	P[1].B = 4.3;
	P[1].C = 7.5;

	P[2].A = 8.5;
	P[2].B = 2.3;
	P[2].C = 4.1;

	P[3].A = 5.5;
	P[3].B = 8.1;
	P[3].C = 6.1;

	P[4].A = 7.5;
	P[4].B = 9.3;
	P[4].C = 2.6;

	Calcola(P, 5, &r1, &r2);

	printf("\nDistanza Min: P[%d] - P[%d] -> %.15lf\n", r1.index1, r1.index2, r1.Distanza);
	printf("\nDistanza Max: P[%d] - P[%d] -> %.15lf\n", r2.index1, r2.index2, r2.Distanza);

	return 0;
}

[gugoXX]

Originariamente inviato da Tommo

Non capisco... ma per "casuale" si intende che va bene una configurazione qualsiasi, oppure che i valori iniziali devono essere tirati a caso?
C'è una bella differenza...

Nel primo caso la soluzione più veloce è:
nextPos = pos + c;
Sarebbe ben facile

Nel secondo invece c'è un problema... la configurazione casuale migliore potrebbe uscire alla seconda iterazione come dopo 1000... rendendo difatti insensato il contest...
credo sia il caso di fare un file che contenga le varie configurazioni iniziali.

Casuale perche' o riesci a trovare una soluzione a priori "cristallina" a questo problema, e non e' semplicissimo quando non si hanno numeri di particelle "belli", oppure mi sa tanto che quella da te proposta non e' sufficientemente esaustiva.
Mi spiego, a forza di aggiungere +c alla coordinata da te scelta, prima o poi cozzerai contro la scatola, e allora cosa farai?

L'idea sarebbe quella di trovare un metodo per la ricerca delle soluzioni casuali, non un'unica soluzione deterministica e sempre uguale a fronte dei valori scelti, ma anche cosi' va bene, c'e' da divertirsi.

[Tommo]

Originariamente inviato da gugoXX

Casuale perche' o riesci a trovare una soluzione a priori "cristallina" a questo problema, e non e' semplicissimo quando non si hanno numeri di particelle "belli", oppure mi sa tanto che quella da te proposta non e' sufficientemente esaustiva.
Mi spiego, a forza di aggiungere +c alla coordinata da te scelta, prima o poi cozzerai contro la scatola, e allora cosa farai?

L'idea sarebbe quella di trovare un metodo per la ricerca delle soluzioni casuali, non un'unica soluzione deterministica e sempre uguale a fronte dei valori scelti, ma anche cosi' va bene, c'e' da divertirsi.

vabè quella non era una soluzione seria

Era solo per dire che il problema non è affatto difficile se posso mettermi le particelle come mi pare... basta pensarci un pò di più per capire che per occupare il volume minore devi affiancare tutti tetraedroni, come nei cristalli... in quella maniera ogni particella è al centro di un poliedro di raggio c... in 2d sarebbero tutti esagoni.

E' un pò difficile a farsi ma credo non si possano compattare di più...

[marco.r]

Originariamente inviato da gugoXX

Casuale perche' o riesci a trovare una soluzione a priori "cristallina" a questo problema, e non e' semplicissimo quando non si hanno numeri di particelle "belli", oppure mi sa tanto che quella da te proposta non e' sufficientemente esaustiva.
Mi spiego, a forza di aggiungere +c alla coordinata da te scelta, prima o poi cozzerai contro la scatola, e allora cosa farai?

Continuo con la fila successiva . E da quattro conti che ho fatto sembra funzionare per le varie scelte di N.

In un cubo mi ci stanno ( floor(d/c)+1 )^3 sfere (nell'ipotesi che una molecola possa stare "appicicata al bordo").
Quindi diciamo almeno (d/c)^3
Sostituendo la definizione di c ottieni (d/d * pow(N,1/3))^3 molecole, che guardacaso coincide proprio con N.
Dimostra contemporaneamente il teorema sopra e mostra un modo per impacchettare le molecole (anzi, scommetto che e' il modo in cui e' stato ricavato quel valore).
Quindi la soluzione "cristallina" e' effettivamente banale.

[gugoXX]

Esatto! Questa e' la dimostrazione che pero' funziona per "cubi", nel senso quando le particelle sono 1,8,27, etc.

Quando sono... 17, impaccarle in quel modo non e' semplice (e mi sa che non funziona). Pero' il teorema e' sempre valido.

[Vincenzo1968]

Mbe'? Nessuno sa dirmi se ho eseguito correttamente i calcoli?

Data la seguente lista di punti:

P[0] -> 1.5 2.3 5.1
P[1] -> 2.8 4.3 7.5
P[2] -> 8.5 2.3 4.1
P[3] -> 5.5 8.1 6.1
P[4] -> 7.5 9.3 2.6

è corretto dire che la coppia di punti a distanza minima è (P[0],P[1]) e quella a distanza massima è, invece, (P[0],P[4])?

[marco.r]

Originariamente inviato da Tommo

vabè quella non era una soluzione seria

Era solo per dire che il problema non è affatto difficile se posso mettermi le particelle come mi pare... basta pensarci un pò di più per capire che per occupare il volume minore devi affiancare tutti tetraedroni, come nei cristalli... in quella maniera ogni particella è al centro di un poliedro di raggio c... in 2d sarebbero tutti esagoni.

E' un pò difficile a farsi ma credo non si possano compattare di più...

Non penso sia quella ottimale, resta molto spazio tra le sfere.
Si potrebbe provare a vedere quante sfere si rescono effettivamente a compattare al massimo

[Tommo]

Originariamente inviato da marco.r

Non penso sia quella ottimale, resta molto spazio tra le sfere.
Si potrebbe provare a vedere quante sfere si rescono effettivamente a compattare al massimo

E' ottimale

Semplicemente, TUTTE le coppie di sfere sono distanti c, che è il minimo.
E se tutte quante le sfere sono fra loro distanti c, non si possono compattare di più...

prendi il Geomag per esempio: la struttura più densa possibile è quella composta da piramidi a base triangolare... se uno prende le palline come particelle, e le stecche come la distanza minima c.
Non c'è maniera di compattarlo di più, e credo valga lo stesso per i barioni, qualunque cosa siano

[Torav]

Originariamente inviato da Tommo

E' ottimale

Semplicemente, TUTTE le coppie di sfere sono distanti c, che è il minimo.
E se tutte quante le sfere sono fra loro distanti c, non si possono compattare di più...

prendi il Geomag per esempio: la struttura più densa possibile è quella composta da piramidi a base triangolare... se uno prende le palline come particelle, e le stecche come la distanza minima c.
Non c'è maniera di compattarlo di più, e credo valga lo stesso per i barioni, qualunque cosa siano

In verità i reticoli con il packing più alto (74.0% del volume occupato dalle pallozze ) sono diversi: ad esempio quello a base esagonale e il reticolo cubico a facce centrate. Per convincersi della cosa basta pensare di "stendere" delle palline su un piano in modo che siano il più vicino possibile. Ora mettiamo il secondo strato, le palline finiranno proprio negli avvallamenti del primo strato. Se ora mettiamo un terzo strato nella stessa posizione del primo otteniamo un reticolo esagonale, se invece lo mettiamo ancora girato di 60 gradi otteniamo l'FCC. Se provate con le biglie si vede chiaramente che più strette di così non ci possono stare
se domani ho tempo voglio vedere cosa riesco a tirare fuori

[Tommo]

Originariamente inviato da Torav

In verità i reticoli con il packing più alto (74.0% del volume occupato dalle pallozze ) sono diversi: ad esempio quello a base esagonale e il reticolo cubico a facce centrate. Per convincersi della cosa basta pensare di "stendere" delle palline su un piano in modo che siano il più vicino possibile. Ora mettiamo il secondo strato, le palline finiranno proprio negli avvallamenti del primo strato. Se ora mettiamo un terzo strato nella stessa posizione del primo otteniamo un reticolo esagonale, se invece lo mettiamo ancora girato di 60 gradi otteniamo l'FCC. Se provate con le biglie si vede chiaramente che più strette di così non ci possono stare
se domani ho tempo voglio vedere cosa riesco a tirare fuori

Molto interessante... ma mi sfugge una cosa: se ruoti un esagono di 60 gradi, non cambia niente...
Per cui non capisco come fa la rotazione a compattarlo di più

[gugoXX]

Per stimolare i giovani virgulti ho aggiornato il punto1, migliorandone la definizione e estendendolo con un file di prova.
Ora cerco di mettere il primo risultato forza bruta di base.

[gugoXX]

Ecco

Max Dist: 1.70445415447397 Indexes: 48483-72832
Min Dist: 0.000366913641159708 Indexes: 22503-72620

234827 (Millisecondi)

Qui il codice.
Vi chiederei quanto piu' possibile di tenere separate la fase di input da quella di elaborazione, qualunque essa sia.
Il risultato della fase di input vorrebbe essere una semplice lista o array in memoria, da elaborare successivamente.

PS: C'e' gia' una prima banalissima ottimizzazione, come gia' proposto, ovvero di non confrontare le vere e proprie distanze, ma il loro quadrato
in modo tale da evitare tante, tantissime radici quadrate.
Nel file di prova si risparmiano la bellezza di 4999950000 radici quadrate.

class Program
{    
    static void Main(string[] args)
    {
        DistCalc dc = new DistCalc();
        dc.Load(@"C:\temp\Coords.dat");
        Stopwatch sw = new Stopwatch();
        sw.Start();
        var res = dc.Research();
        sw.Stop();
        Console.WriteLine(res);
        Console.WriteLine(sw.ElapsedMilliseconds);
        Console.ReadKey();
    }
}

public class DistCalc
{
    List<Point3d> allpoints;
    public DistCalc()
    {
        allpoints = new List<Point3d>();
    }

    public void Load(string fname)
    {
        StreamReader sr = new StreamReader(fname);
        string stt = sr.ReadLine();
        string[] first = stt.Split(new char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries);
        int howmany = int.Parse(first[0]);
        InputSize = double.Parse(first[1]);
        string line;
        while ((line = sr.ReadLine()) != null)
        {
            string[] coords = line.Split(new char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries);
            allpoints.Add(new Point3d(double.Parse(coords[0]),double.Parse(coords[1]),double.Parse(coords[2])));
        }
        sr.Close();
    }

    public void Save(string fname)
    {
        StreamWriter sw = new StreamWriter(fname);
        sw.WriteLine("{0} {1}", allpoints.Count, InputSize);
        foreach (Point3d p3d in allpoints)
        {
            sw.WriteLine(p3d.ToString());
        }
        sw.Close();
    }

    public double InputSize;
    public void Generate(int npoints,double dim)
    {
        InputSize = dim;
        allpoints.Clear();
        for (int t = 0; t < npoints; t++)
        {
            allpoints.Add(Point3d.Rand(dim));
        }
    }

    public class Result
    {
        public Result()
        {
            curdimmin = double.MaxValue;
            curdimmax = 0d;
        }

        public double curdimmin;
        public int i1dimmin, i2dimmin;

        public double curdimmax;
        public int i1dimmax, i2dimmax;

        public override string ToString()
        {
            StringWriter sw = new StringWriter();
            sw.WriteLine("Max Dist: {0} Indexes: {1}-{2}", Math.Sqrt(curdimmax), i1dimmax, i2dimmax);
            sw.WriteLine("Min Dist: {0} Indexes: {1}-{2}", Math.Sqrt(curdimmin), i1dimmin, i2dimmin);
            return sw.ToString();
        }
    }

    public Result Research()
    {
        Result res = new Result();
        int ln = allpoints.Count;
        for (int t = 0; t < ln - 1; t++)
        {
            for (int u = t+1; u < ln; u++)
            {
                double cd = allpoints[t] - allpoints[u];
                if (cd < res.curdimmin)
                {
                    res.curdimmin = cd;
                    res.i1dimmin = t;
                    res.i2dimmin = u;
                }
                else if (cd > res.curdimmax)
                {
                    res.curdimmax = cd;
                    res.i1dimmax = t;
                    res.i2dimmax = u;
                }
            }
        }
        return res;
    }
}

public struct Point3d
{
    public double x;
    public double y;
    public double z;
    public Point3d(double ix,double iy,double iz)
    {
        x = ix;
        y = iy;
        z = iz;
    }

    public static Random rnd = new Random(18);
    public static Point3d Rand(double dim)
    {
        return new Point3d(rnd.NextDouble() * dim, rnd.NextDouble() * dim, rnd.NextDouble() * dim);
    }

    public static double operator -(Point3d d1, Point3d d2)
    {
        double dx=d2.x-d1.x;
        double dy=d2.y-d1.y;
        double dz=d2.z-d1.z;
        return (dx * dx + dy * dy + dz * dz);
    }

    public override string ToString()
    {
        return string.Format("{0} {1} {2}", x, y, z);
    }
}

[Vincenzo1968]

Originariamente inviato da gugoXX

Per stimolare i giovani virgulti ho aggiornato il punto1, migliorandone la definizione e estendendolo con un file di prova.
Ora cerco di mettere il primo risultato forza bruta di base.

Cosa rappresenta quell'uno dopo il diecimila nella prima riga del file?

[Torav]

Originariamente inviato da Tommo

Molto interessante... ma mi sfugge una cosa: se ruoti un esagono di 60 gradi, non cambia niente...
Per cui non capisco come fa la rotazione a compattarlo di più

ups! Era tardi..quello che intendo dire è che il terzo strato si può mettere in due modi diversi riempiendo esattamente la stessa % di spazio! Ora sto uscendo, se trovo qualche disegno più chiaro lo posto perchè altrimenti ci vuole una bella dose di immaginazione!

[Torav]

ecco un link che fa vedere decentemente la cosa, e con questo chiudo il leggero OT, scusate