heuristica.c

001   #include "heuristica.h"
002   #include <limits.h>
003   #include <string.h>
004   #include "priorqueue.h"
005   
006   
007   #include <stdio.h>
008   #include <stdlib.h>
009   
010   //------------------------------------------------------------------------------
011   
012   typedef struct {
013       int nVert;
014       int nPai;
015       int nDist;
016       int nPqIndex;
017   } TVerticeData;
018   
019   
020   typedef struct {
021       TVerticeData *data[MAX_VERTICES];
022   } TVerticeDataPriorQueue;
023   
024   
025   // Todas as rotinas seguintes implementam a abstração da manipulação
026   // de vetores de vértices, para o uso em filas de prioridades.
027   
028   int VdpqCompare(void *pq, int nIndex1, int nIndex2)
029   {
030       int k1 = ((TVerticeDataPriorQueue *) pq)->data[nIndex1]->nDist;
031       int k2 = ((TVerticeDataPriorQueue *) pq)->data[nIndex2]->nDist;
032   
033       return (k1 < k2 ? -1 : k1 > k2 ? 1 : 0);
034   }
035   
036   
037   int VdpqCompareKey(void *pq, int nIndex, int nKey)
038   {
039       int k = ((TVerticeDataPriorQueue *) pq)->data[nIndex]->nDist;
040   
041       return (k < nKey ? -1 : k > nKey ? 1 : 0);
042   }
043   
044   
045   int VdpqGetId(void *pq, int nIndex)
046   {
047       return ((TVerticeDataPriorQueue *) pq)->data[nIndex]->nVert;
048   }
049   
050   
051   void VdpqQueueIndexChanged(void *pq, int nIndex)
052   {
053       ((TVerticeDataPriorQueue *) pq)->data[nIndex]->nPqIndex = nIndex;
054   }
055   
056   
057   void VdpqSetKey(void *pq, int nIndex, int nKey)
058   {
059       ((TVerticeDataPriorQueue *) pq)->data[nIndex]->nDist = nKey;
060   }
061   
062   
063   void VdpqSwap(void *pq, int nIndex1, int nIndex2)
064   {
065       TVerticeData *aux = ((TVerticeDataPriorQueue *) pq)->data[nIndex1];
066       ((TVerticeDataPriorQueue *) pq)->data[nIndex1] =
067           ((TVerticeDataPriorQueue *) pq)->data[nIndex2];
068       ((TVerticeDataPriorQueue *) pq)->data[nIndex2] = aux;
069   }
070   
071   
072   // Estrutura que representa o controle de vetores de vértices
073   TPriorityQueueControl pqcVdpq = {
074       &VdpqCompare,
075       &VdpqCompareKey,
076       &VdpqGetId,
077       &VdpqQueueIndexChanged,
078       &VdpqSetKey,
079       &VdpqSwap,
080   };
081   
082   //------------------------------------------------------------------------------
083   
084   typedef struct {
085       int u, v;
086       int nDist;
087       char cCor;
088   } TArestaData;
089   
090   
091   typedef struct {
092       TArestaData *data[MAX_VERTICES * MAX_VERTICES / 2];
093   } TArestaDataPriorQueue;
094   
095   
096   // Todas as rotinas seguintes implementam a abstração da manipulação
097   // de vetores de arestas, para o uso em filas de prioridades.
098   
099   int AdpqCompare(void *pq, int nIndex1, int nIndex2)
100   {
101       int k1 = ((TArestaDataPriorQueue *) pq)->data[nIndex1]->nDist;
102       int k2 = ((TArestaDataPriorQueue *) pq)->data[nIndex2]->nDist;
103   
104       return (k1 < k2 ? -1 : k1 > k2 ? 1 : 0);
105   }
106   
107   
108   int AdpqCompareKey(void *pq, int nIndex, int nKey)
109   {
110       int k = ((TArestaDataPriorQueue *) pq)->data[nIndex]->nDist;
111   
112       return (k < nKey ? -1 : k > nKey ? 1 : 0);
113   }
114   
115   
116   int AdpqGetId(void *pq, int nIndex)
117   {
118       return ((TArestaDataPriorQueue *) pq)->data[nIndex]->u * MAX_VERTICES +
119        ((TArestaDataPriorQueue *) pq)->data[nIndex]->v;
120   }
121   
122   
123   void AdpqQueueIndexChanged(void *pq, int nIndex)
124   {
125   }
126   
127   
128   void AdpqSetKey(void *pq, int nIndex, int nKey)
129   {
130       ((TArestaDataPriorQueue *) pq)->data[nIndex]->nDist = nKey;
131   }
132   
133   
134   void AdpqSwap(void *pq, int nIndex1, int nIndex2)
135   {
136       TArestaData *aux = ((TArestaDataPriorQueue *) pq)->data[nIndex1];
137       ((TArestaDataPriorQueue *) pq)->data[nIndex1] =
138           ((TArestaDataPriorQueue *) pq)->data[nIndex2];
139       ((TArestaDataPriorQueue *) pq)->data[nIndex2] = aux;
140   }
141   
142   
143   // Estrutura que representa o controle de vetores de arestas
144   TPriorityQueueControl pqcAdpq = {
145       &AdpqCompare,
146       &AdpqCompareKey,
147       &AdpqGetId,
148       &AdpqQueueIndexChanged,
149       &AdpqSetKey,
150       &AdpqSwap,
151   };
152   
153   //------------------------------------------------------------------------------
154   
155   // Dado o conjunto de arestas 'ad', de tamanho 'nSizead', esta rotina
156   // procura por um ciclo que parte e chega ao vértice 'u'. O ciclo encontrado
157   // poderá incluir sub-ciclos.
158   void heuristProcuraCiclo(TPcvData *PcvData, TArestaData *ad, int nSizead, int u,
159       int *anCicloEncontrado, int *nSizeCiclo)
160   {
161       int a, v, aAdot, vAdot, nPrev = u;
162       int nMinDist;
163   
164       *nSizeCiclo = 1;
165       anCicloEncontrado[0] = u;
166   
167       for (;;) {
168           vAdot = -1;
169           nMinDist = INT_MAX;
170   
171           for (a = 0; a < nSizead; a++) {
172               if (ad[a].cCor != COR_BRANCO) continue;
173               if (ad[a].u == u) {
174                   v = ad[a].v;
175               } else if (ad[a].v == u) {
176                   v = ad[a].u;
177               } else {
178                   continue;
179               }
180   
181               if (PcvData->anDist[nPrev][v] < nMinDist) {
182                   nMinDist = PcvData->anDist[nPrev][v];
183                   vAdot = v;
184                   aAdot = a;
185               }
186           }
187   
188           if (vAdot == -1) break;
189           v = vAdot;
190           ad[aAdot].u = u;
191           ad[aAdot].v = v;
192           ad[aAdot].cCor = COR_CINZA;
193   
194           if (PcvData->acCorVert[v] != COR_PRETO) {
195               PcvData->acCorVert[v] = COR_PRETO;
196               anCicloEncontrado[*nSizeCiclo] = v;
197               (*nSizeCiclo)++;
198               nPrev = v;
199           }
200           u = v;
201       }
202   }
203   
204   
205   // Dado o conjunto de arestas 'ad', de tamanho 'nSizead', esta rotina
206   // procura por um ciclo que parte e chega ao vértice 'u' e que inclui
207   // todas as arestas de 'ad'. Portanto, se houver mais arestas do que
208   // vértices, o ciclo encontrado necessariamente terá sub-ciclos.
209   void heuristProcuraCicloTotal(TPcvData *PcvData, TArestaData *ad, int nSizead,
210       int u)
211   {
212       int a, v;
213       int anSubCiclo[MAX_VERTICES], anSizeSubCiclo;
214   
215       PcvData->acCorVert[u] = COR_PRETO;
216       heuristProcuraCiclo(PcvData, ad, nSizead, u, PcvData->anCaminho,
217           &PcvData->nVisitados);
218   
219       while (PcvData->nVisitados < PcvData->nNumVertices) {
220           u = -1;
221           for (a = 0; a < nSizead; a++) {
222               if (ad[a].cCor != COR_BRANCO) continue;
223               if (PcvData->acCorVert[ad[a].u] == COR_PRETO) {
224                   u = ad[a].u;
225                   break;
226               } else if (PcvData->acCorVert[ad[a].v] == COR_PRETO) {
227                   u = ad[a].v;
228                   break;
229               }
230           }
231   
232           if (u == -1) break;
233   
234           heuristProcuraCiclo(PcvData, ad, nSizead, u, anSubCiclo, &anSizeSubCiclo);
235   
236           for (v = 0; v < PcvData->nVisitados; v++) {
237               if (PcvData->anCaminho[v] != u) continue;
238   
239               anSizeSubCiclo--;
240               memmove(&PcvData->anCaminho[v + anSizeSubCiclo],
241                   &PcvData->anCaminho[v], (PcvData->nVisitados - v) * sizeof(int));
242               memmove(&PcvData->anCaminho[v + 1], &anSubCiclo[1],
243                   anSizeSubCiclo * sizeof(int));
244               PcvData->nVisitados += anSizeSubCiclo;
245               break;
246           }
247       }
248   }
249   
250   
251   // Calcula a distância do ciclo armazenado no vetor PcvData->anDist,
252   // possivelmente armazenando-o como sendo o melhor caminho conseguido
253   // se for este o caso
254   int calcDistCaminho(TPcvData *PcvData)
255   {
256       int i;
257   
258       // Calcula o custo do caminho descoberto
259       PcvData->nDistTotal = PcvData->anDist
260           [PcvData->anCaminho[PcvData->nNumVertices - 1]]
261           [PcvData->anCaminho[0]];
262       for (i = 0; i < PcvData->nNumVertices - 1; i++) {
263           PcvData->nDistTotal += PcvData->anDist[PcvData->anCaminho[i]]
264               [PcvData->anCaminho[i + 1]];
265       }
266       PcvData->nDistTotal -= PcvData->nNumVertices * PcvData->nDistDelta;
267   
268       // Se for o menor até agora, então armazena-o
269       if (PcvData->nDistTotal < PcvData->nMinDistTotal) {
270           PcvData->nMinDistTotal = PcvData->nDistTotal;
271           memcpy(PcvData->anMinCaminho, PcvData->anCaminho,
272               sizeof(PcvData->anCaminho));
273       }
274   
275       return PcvData->nDistTotal;
276   }
277   
278   
279   // Corpo principal da heurística proposta
280   void procCaminhoHeuristica(TPcvData *PcvData)
281   {
282       int a, i, j, k, u, v, nIndexMin, nVertStart = 0;
283       int nChave_i, nChave_j, nChave_k;
284       TVerticeData vd[MAX_VERTICES], *vdMin;
285       int anGrauVert[MAX_VERTICES];
286       TVerticeDataPriorQueue vdpq;
287       TArestaData *ad, *adMatching;
288       TArestaDataPriorQueue adpq;
289       int nNumVertImpares, nIdAresta;
290   
291       // Inicializa o menor caminho descoberto até agora como "infinito"
292       PcvData->nMinDistTotal = INT_MAX;
293   
294       // Laço principal: Todo o algoritmo é percorrido iniciando em
295       // cada um dos vértices
296       for (nVertStart = 0; nVertStart < PcvData->nNumVertices; nVertStart++) {
297           // Inicializa memória
298           memset(&PcvData->acCorVert, COR_BRANCO, PcvData->nNumVertices
299               * sizeof(char));
300           memset(&anGrauVert, 0, sizeof(anGrauVert));
301   
302           pqcVdpq.nElems = PcvData->nNumVertices;
303           for (i = 0; i < PcvData->nNumVertices; i++) {
304               vd[i].nVert = i;
305               vd[i].nPai = -1;
306               vd[i].nDist = INT_MAX;
307               vd[i].nPqIndex = i;
308               vdpq.data[i] = &vd[i];
309           }
310   
311           // Assegura que o primeiro vértice a ser selecionado é o
312           // de índice nVertStart
313           PriorQueueDecreaseKey(&vdpq, &pqcVdpq, nVertStart, 0);
314           nNumVertImpares = 0;
315   
316           for (;;) {
317               // Obtém o vértice mais próximo ao conjunto de vértices já descobertos.
318               // Caso esta seja a primeira iteração deste laço, o vértice obtido
319               // é o de índice nVertStart
320               nIndexMin = PriorQueueExtractMin(&vdpq, &pqcVdpq);
321               if (nIndexMin < 0) break;
322               vdMin = &vd[nIndexMin];
323               u = vdMin->nVert;
324   
325               // Pinta o vértice de cinza (torna-o já descoberto)
326               PcvData->acCorVert[u] = COR_CINZA;
327   
328               // Iteração para cada um dos vértices adjacentes a 'u',
329               // para descobrir uma Árvore Geradora Mínima (algoritmo de Prim)
330               for (v = 0; v < PcvData->nNumVertices; v++) {
331                   // Apenas prossegue se 'v' não tiver sido descoberto
332                   if (PcvData->acCorVert[v] != COR_BRANCO) continue;
333   
334                   // Se a distância de 'u' a 'v' for menor do que a distância
335                   // aos vértices anteriormente descobertos:
336                   if (PcvData->anDist[u][v] < vd[v].nDist) {
337                       // Atualiza a informação de paridade de grau do antigo vértice pai
338                       if (vd[v].nPai == -1) {
339                           anGrauVert[v] ^= 1;
340                           nNumVertImpares += (anGrauVert[v] ? 1 : -1);
341                       } else {
342                           anGrauVert[vd[v].nPai] ^= 1;
343                           nNumVertImpares += (anGrauVert[vd[v].nPai] ? 1 : -1);
344                       }
345                       // Marca o novo pai
346                       vd[v].nPai = u;
347                       // Atualiza a informação de paridade de grau do novo vértice pai
348                       anGrauVert[u] ^= 1;
349                       nNumVertImpares += (anGrauVert[u] ? 1 : -1);
350   
351                       // Decrementa a chave do vértice 'v' para a nova distância
352                       PriorQueueDecreaseKey(&vdpq, &pqcVdpq, vd[v].nPqIndex,
353                           PcvData->anDist[u][v]);
354                   }
355               }
356           }
357   
358           // Aloca e inicializa memória para descobrir um "matching"
359           adMatching = (TArestaData *) malloc(nNumVertImpares
360               * (nNumVertImpares - 1) / 2 * sizeof(TArestaData));
361           memset(&PcvData->acCorVert, COR_BRANCO, PcvData->nNumVertices
362               * sizeof(char));
363           i = 0;
364   
365           // Para cada par de vétrices de grau ímpar, insere a
366           // aresta correspondente no vetor adMatching
367           for (u = 0; u < PcvData->nNumVertices; u++) {
368               if (!anGrauVert[u]) continue;
369   
370               for (v = u + 1; v < PcvData->nNumVertices; v++) {
371                   if (!anGrauVert[v]) continue;
372                   adMatching[i].u = u;
373                   adMatching[i].v = v;
374                   adMatching[i].nDist = PcvData->anDist[u][v];
375                   adpq.data[i] = &adMatching[i];
376                   i++;
377               }
378           }
379   
380           // Constrói um heap mínimo (adpq) baseado no vetor adMatching,
381           // para ser usado como fila de prioridades
382           pqcAdpq.nElems = i;
383           PriorQueueBuildMinHeap(&adpq, &pqcAdpq);
384           a = 0;
385   
386           // Aloca memória para o conjunto de arestas consideradas
387           ad = (TArestaData *) malloc((PcvData->nNumVertices - 1 + nNumVertImpares / 2)
388               * sizeof(TArestaData));
389   
390           // Inclui na lista principal de arestas todas as arestas da AGM
391           for (u = 0; u < PcvData->nNumVertices; u++) {
392               if (vd[u].nPai < 0) continue;
393   
394               ad[a].u = vd[u].nPai;
395               ad[a].v = u;
396               ad[a].nDist = PcvData->anDist[vd[u].nPai][u];
397               ad[a].cCor = COR_BRANCO;
398               a++;
399           }
400   
401           // Realiza o "matching" por escolha gulosa, usando o heap adpq
402           // como base para a fila de prioridades. Cada aresta descoberta
403           // é adicionada na lista principal de arestas
404           for (i = 0; i < nNumVertImpares; ) {
405               nIdAresta = PriorQueueExtractMin(&adpq, &pqcAdpq);
406               u = nIdAresta / MAX_VERTICES;
407               v = nIdAresta % MAX_VERTICES;
408               if (PcvData->acCorVert[u] != COR_BRANCO) continue;
409               if (PcvData->acCorVert[v] != COR_BRANCO) continue;
410   
411               PcvData->acCorVert[u] = COR_CINZA;
412               PcvData->acCorVert[v] = COR_CINZA;
413               ad[a].u = u;
414               ad[a].v = v;
415               ad[a].nDist = PcvData->anDist[u][v];
416               ad[a].cCor = COR_BRANCO;
417               a++;
418               i += 2;
419           }
420   
421           // Libera memória alocada
422           free(adMatching);
423   
424           // Realiza o caminhamento com curto-circuitos
425           heuristProcuraCicloTotal(PcvData, ad, a, nVertStart);
426   
427           // Calcula o custo do caminho descoberto
428           calcDistCaminho(PcvData);
429   
430           // Algoritmo de melhoramento: Dados os índices de três vértices no caminho,
431           // 'i', 'j' e 'k', o algoritmo realiza todas as permutações diferentes dos
432           // três vértices para procurar por uma solução melhor
433           for (i = 0; i < PcvData->nNumVertices - 2; i++) {
434               nChave_i = PcvData->anCaminho[i];
435               for (j = i + 1; j < PcvData->nNumVertices - 1; j++) {
436                   nChave_j = PcvData->anCaminho[j];
437                   for (k = j + 1; k < PcvData->nNumVertices; k++) {
438                       nChave_k = PcvData->anCaminho[k];
439   
440                       // A ordem original no vetor é: ... i ... j ... k ...
441   
442                       // Altera para ... j ... i ... k ... e testa
443                       PcvData->anCaminho[i] = nChave_j;
444                       PcvData->anCaminho[j] = nChave_i;
445                       calcDistCaminho(PcvData);
446   
447                       // Altera para ... j ... k ... i ... e testa
448                       PcvData->anCaminho[j] = nChave_k;
449                       PcvData->anCaminho[k] = nChave_i;
450                       calcDistCaminho(PcvData);
451   
452                       // Altera para ... i ... k ... j ... e testa
453                       PcvData->anCaminho[i] = nChave_i;
454                       PcvData->anCaminho[k] = nChave_j;
455                       calcDistCaminho(PcvData);
456   
457                       // Altera para ... k ... i ... j ... e testa
458                       PcvData->anCaminho[i] = nChave_k;
459                       PcvData->anCaminho[j] = nChave_i;
460                       calcDistCaminho(PcvData);
461   
462                       // Altera para ... k ... j ... i ... e testa
463                       PcvData->anCaminho[j] = nChave_j;
464                       PcvData->anCaminho[k] = nChave_i;
465                       calcDistCaminho(PcvData);
466   
467                       // Restaura a ordem original: ... i ... j ... k ...
468                       PcvData->anCaminho[i] = nChave_i;
469                       PcvData->anCaminho[k] = nChave_k;
470                   }
471               }
472           }
473   
474           // Libera memória alocada
475           free(ad);
476       }
477   }