/*************************************************************************************************************************************/
/* */
/* M A X I M I S A T I O N D E L A P L U S P E T I T E D I S T A N C E */
/* A L ' I N T E R I E U R D ' U N N U A G E D E P O I N T S T R I D I M E N S I O N N E L */
/* E N C O O R D O N N E E S C A R T E S I E N N E S */
/* P A R U N E M E T H O D E D U T Y P E " R E C U I T S I M U L E " */
/* A V E C T E M P E R A T U R E : */
/* */
/* */
/* Author of '$xtc/recuit_si.43$c' : */
/* */
/* Jean-Francois COLONNA (LACTAMME, AAAAMMJJhhmmss). */
/* */
/*************************************************************************************************************************************/
#include "INCLUDES.01.I"
/* Introduit le 20051116102340... */
#define nPOINTS \
27 \
/* Nombre de points du nuage... */
#define nITERATIONS \
90000 \
/* Nombre d'iterations maximum... */
#define DELTA_INITIAL \
80.0
#define DELTA_FINAL \
1.0
/* Extrema des variations possibles des coordonnees {x,y}. */
#define PROBABILITE_DE_PERTURBER_UN_POINT \
0.1 \
/* Probabilite de perturber le point courant... */
#define PERIODISER_L_UNIVERS \
0 \
/* L'univers n'est pas periodique et donc limite (0) ou est periodique et donc infini (1). */
#define PROBABILITE_D_ACCEPTER_UNE_MAUVAISE_SOLUTION \
0.02 \
/* Probabilite d'accepter une mauvaise solution... */
#define TEMPS \
(((double)(iterations-1))/((double)(nITERATIONS-1))) \
/* Simulation d'un temps dans [0,1]. */
extern double drand48();
#define PERTURBATION_ALEATOIRE \
((drand48() - 0.5)*((delta_initial*(1-temps))+(delta_final*temps))) \
/* Plus grande variation possible des coordonnees {x,y}, en rappelant que : */ \
/* */ \
/* drand48() E [0,1] */ \
/* */ \
/* d'ou la translation de 0.5 qui permet d'avoir un deplacement algebrique dans n'importe */ \
/* quelle direction et ce de facon isotrope... */
typedef struct point {
double x;
double y;
double z;
} point;
/* Structure de definition d'un point... */
#define X(P) \
(P.x)
#define Y(P) \
(P.y)
#define Z(P) \
(P.z)
/* Definition de l'acces aux coordonnees d'un point. */
extern double sqrt();
#define Carre(x) \
((x)*(x)) \
/* Definition du carre de 'x'. */
#define DistanceEuclidienne(A,B) \
sqrt(Carre(X(B)-X(A)) + Carre(Y(B)-Y(A)) + Carre(Z(B)-Z(A))) \
/* Definition de la distance entre deux points 'A' et 'B'. */
#define PERTURBATION(coordonnee,minimum,maximum) \
{ \
int perturber=1; \
double coordonnee_perturbee; \
double temps=TEMPS; \
\
while (perturber == 1) \
{ \
coordonnee_perturbee = coordonnee + PERTURBATION_ALEATOIRE; \
/* Perturbation de la coordonnee courante... */ \
\
if (PERIODISER_L_UNIVERS == 0) \
{ \
} \
else \
{ \
if (coordonnee_perturbee < minimum) \
{ \
coordonnee_perturbee = coordonnee_perturbee - minimum + maximum; \
} \
else \
{ \
} \
\
if (coordonnee_perturbee > maximum) \
{ \
coordonnee_perturbee = coordonnee_perturbee - maximum + minimum; \
} \
else \
{ \
} \
} \
\
if ((coordonnee_perturbee < minimum) || (coordonnee_perturbee > maximum)) \
{ \
/* Les coordonnees hors-ecran sont rejetees... */ \
} \
else \
{ \
coordonnee = coordonnee_perturbee; \
perturber=0; \
/* On s'arrete sur la premiere coordonnee perturbee situee dans l'ecran... */ \
} \
} \
}
#define DISTANCE_MINIMALE(distance_minimale) \
{ \
distance_minimale=+1e+300; \
\
for (i=0 ; i<nPOINTS ; i++) \
{ \
for (j=0 ; j<nPOINTS ; j++) \
{ \
if (i > j) \
{ \
double distance_courante=DistanceEuclidienne(Lpoints[i],Lpoints[j]); \
distance_minimale=MIN2(distance_minimale,distance_courante); \
/* Recherche de la distance minimale courante du nuage de points... */ \
} \
else \
{ \
} \
} \
} \
}
extern void *malloc();
static int dimX=0;
#define Xmin 0
#define Xmax (Xmin + (dimX-1))
/* Definition des abscisses. */
static int dimY=0;
#define Ymin 0
#define Ymax (Ymin + (dimY-1))
/* Definition des ordonnees. */
static int dimZ=0;
#define Zmin 0
#define Zmax (Zmin + (dimZ-1))
/* Definition des profondeurs. */
#define IMAGE(x,y) \
(*(image + (((((int)y)-Ymin)*dimX) + (((int)x)-Xmin)))) \
/* Acces a un point de l'image. */
#define PROJECT 0.0
#define store(n,x,y,z) \
{ \
int xp=x+((int)(PROJECT*((double)(z)))); \
int yp=y+((int)(PROJECT*((double)(z)))); \
/* "Pseudo-projection" du point {x,y,z}. */ \
\
if (((xp>=Xmin) && (xp<=Xmax)) && ((yp>=Ymin) && (yp<=Ymax))) \
{ \
IMAGE(xp,yp) = n; \
} \
else \
{ \
} \
} \
/* Rangement d'un point valide d'une image. */
#define initialisation(x,y) \
{ \
if (((x>=Xmin) && (x<=Xmax)) && ((y>=Ymin) && (y<=Ymax))) \
{ \
IMAGE(x,y) = NOIR; \
} \
else \
{ \
} \
} \
/* Initialisation de l'image. */
main()
{
int generer_les_positions_finales=0;
/* Pour controler la generation des positions finales uniquement (=1) ou des positions */
/* finales et intermediaires (=0). */
int sortir_l_image=1;
/* Pour controler la sortie (=1) ou pas (=0) de l'image... */
int x,y;
/* Definition des coordonnees. */
unsigned char *image;
/* Definition de l'image a generer... */
int i,j,n;
int iterations;
/* Index divers... */
point Lpoints[nPOINTS];
point SLpoints[nPOINTS];
/* Liste des points et sa Sauvegarde. */
double delta_initial,delta_final=DELTA_FINAL;
/* Bornes perturbation des coordonnees. */
double distance_minimale_avant;
double distance_minimale_apres;
int distance_minimale_apres_valide=0;
Get(dimX,"dimX");
Get(dimY,"dimY");
Get(dimZ,"dimZ");
dimZ = 4*dimZ;
/* Recuperation des dimensions en 'X', en 'Y' et en 'Z' de l'image a generer. On notera la */
/* multiplication par 4 de 'dimZ' car il est en general 4 fois plus petit que les autres... */
image=malloc(dimY*dimX);
/* Definition de l'image a generer... */
for (y=Ymin ; y<=Ymax ; y++)
{
for (x=Xmin ; x<=Xmax ; x++)
{
initialisation(x,y);
/* Initialisation de l'image a generer... */
}
}
delta_initial = DELTA_INITIAL*(((double)MIN2(dimX,dimY))/((double)512));
for (n=0 ; n<nPOINTS ; n++)
{
double aveugle1=drand48();
double aveugle2=drand48();
double aveugle3=drand48();
double aveugle4=drand48();
double aveugle5=drand48();
double aveugle6=drand48();
/* Necessaire (c'est l'experience qui le montre...). */
X(Lpoints[n]) = Xmin + (0.5*dimX);
Y(Lpoints[n]) = Ymin + (0.5*dimY);
Z(Lpoints[n]) = Zmin + (0.5*dimZ);
/* Initialisation arbitraire des points. L'experience montre qu'une configuration */
/* initiale tres fortement symetrique est preferable a une situation aleatoire... */
}
for (iterations=1 ; iterations<=nITERATIONS ; iterations++)
{
for (n=0 ; n<nPOINTS ; n++)
{
X(SLpoints[n]) = X(Lpoints[n]);
Y(SLpoints[n]) = Y(Lpoints[n]);
Z(SLpoints[n]) = Z(Lpoints[n]);
/* Sauvegarde de l'etat courant des points... */
}
if (distance_minimale_apres_valide == 0)
{
if (iterations == 1)
{
DISTANCE_MINIMALE(distance_minimale_avant);
/* Recherche de la distance minimale avant la perturbation. */
}
else
{
/* Dans ce cas 'distance_minimale_avant' est encore valide (calculee lors d'une iteration */
/* anterieure...). */
}
}
else
{
distance_minimale_avant = distance_minimale_apres;
}
if (iterations == 1)
{
fprintf(stderr,"distance minimale initiale=%f\n",distance_minimale_avant);
}
else
{
}
for (n=0 ; n<nPOINTS ; n++)
{
if (drand48() < PROBABILITE_DE_PERTURBER_UN_POINT)
{
double X_perturbe=X(Lpoints[n]);
double Y_perturbe=Y(Lpoints[n]);
double Z_perturbe=Z(Lpoints[n]);
PERTURBATION(X_perturbe,Xmin,Xmax);
PERTURBATION(Y_perturbe,Ymin,Ymax);
PERTURBATION(Z_perturbe,Zmin,Zmax);
X(Lpoints[n]) = X_perturbe;
Y(Lpoints[n]) = Y_perturbe;
Z(Lpoints[n]) = Z_perturbe;
/* Perturbation aleatoire des points. */
}
else
{
}
}
DISTANCE_MINIMALE(distance_minimale_apres);
/* Recherche de la distance minimale apres la perturbation. */
if ( (distance_minimale_apres < distance_minimale_avant)
|| (drand48() < PROBABILITE_D_ACCEPTER_UNE_MAUVAISE_SOLUTION)
)
/* On cherche a maximiser la distance minimale. On notera le test '<' (alors que la */
/* logique voudrait '<=') a cause du cas ou les points occuperaient initialement tous */
/* la meme position (d'ou : distance_minimale_avant=0) : dans ce cas, comme un certain */
/* nombre de points ne sont pas perturbes, on a donc distance_minimale_apres=0, d'ou avec */
/* avec le test '<=' le rejet de cette solution ; ainsi, on risque d'etre bloque */
/* indefiniment puisque l'on restaure alors l'etat precedent... */
{
distance_minimale_apres_valide = 0;
/* L'etat anterieur est inchange... */
for (n=0 ; n<nPOINTS ; n++)
{
X(Lpoints[n]) = X(SLpoints[n]);
Y(Lpoints[n]) = Y(SLpoints[n]);
Z(Lpoints[n]) = Z(SLpoints[n]);
/* Restauration de l'etat anterieur lorsqu'il y a degradation des performances... */
}
}
else
{
distance_minimale_apres_valide = 1;
/* On conserve un etat perturbe qui maximise la distance minimale... */
}
if (generer_les_positions_finales == 1)
{
}
else
{
for (n=0 ; n<nPOINTS ; n++)
{
store(MAX2(NOIR+1,(BLANC*iterations)/nITERATIONS),X(Lpoints[n]),Y(Lpoints[n]),Z(Lpoints[n]));
/* Marquage des points... */
}
}
if (iterations == nITERATIONS)
{
fprintf(stderr,"distance minimale finale..=%f\n",distance_minimale_avant);
/* ATTENTION : c'est la distance minimale a l'avant-derniere iteration qui est editee car, */
/* en effet, il n'est pas sur que 'distance_minimale_apres' etait maximale... */
}
else
{
}
}
if (generer_les_positions_finales == 1)
{
for (n=0 ; n<nPOINTS ; n++)
{
store(BLANC,X(Lpoints[n]),Y(Lpoints[n]),Z(Lpoints[n]));
/* Marquage des points... */
}
}
else
{
}
if (sortir_l_image == 1)
{
write(1,image,dimX*dimY);
/* Sortie de l'image... */
}
else
{
}
}