In Anabatic, correctly set the wire width according to the routing gauge.

* New: Anabatic::AutoContact::setLayerAndWidth() to set both layer and
    VIA width/side according to the RoutingGauge. If the delta in zero,
    use the metal gauge, and the VIA gauge otherwise.
* Bug: In Anabatic::AutoContactTerminal, Anabatic::AutoContactTurn,
    Anabatic::AutoContactHTee & Anabatic::AutoContactVTee, in the
    updateTopology() method, set both the layer and the VIA with when
    there is a change of layer. Note that this default size may be
    overriden later by updateSeize() in the case of non-default width
    segments.
* New: In Anabatic::AutoSegment, new overload setLayer(size_t) to set
    both layer and segment with according to the routing gauge.
* Bug: In Anabatic::AutoHorizontal and Anabatic::AutoVertical, in method
    _makeDogleg(), make use of the new setLayer() to correctly set up
    the wire width.
      Idem for Anabatic::AutoSegment::changeDepth() and ::makeDogleg(),
    and in Anabatic::LayerAssign.
This commit is contained in:
Jean-Paul Chaput 2019-02-12 12:43:09 +01:00
parent 2f9a7076cb
commit ea4fdfd230
32 changed files with 76 additions and 935015 deletions

View File

@ -20,6 +20,7 @@
find_package(PythonLibs 2 REQUIRED)
find_package(PythonSitePackages REQUIRED)
find_package(FLUTE REQUIRED)
find_package(VLSISAPD REQUIRED)
find_package(HURRICANE REQUIRED)
find_package(CORIOLIS REQUIRED)

View File

@ -552,6 +552,20 @@ namespace Anabatic {
}
void AutoContact::setLayerAndWidth ( size_t delta, size_t depth )
{
if (delta == 0) {
setLayer( Session::getRoutingLayer(depth) );
setSizes( Session::getWireWidth (depth)
, Session::getWireWidth (depth) );
} else {
setLayer( Session::getContactLayer(depth) );
setSizes( Session::getViaWidth (depth)
, Session::getViaWidth (depth) );
}
}
AutoContact* AutoContact::createFrom ( Contact* hurricaneContact )
{
AutoContact* autoContact = NULL;

View File

@ -309,10 +309,10 @@ namespace Anabatic {
if (depthH1 == depthH2) {
// Dogleg on the vertical.
switch ( delta ) {
case 0: setLayer( rg->getRoutingLayer(minDepth) ); break;
case 1: setLayer( rg->getContactLayer(minDepth) ); break;
case 0:
case 1: setLayerAndWidth( delta, minDepth ); break;
default:
setLayer( rg->getContactLayer( depthH1 + ((depthH1==minDepth)?0:-1) ) );
setLayerAndWidth( delta, depthH1 + ((depthH1==minDepth)?0:-1) );
_vertical1 = static_cast<AutoVertical*>( _vertical1->makeDogleg(this) );
break;
}
@ -322,12 +322,12 @@ namespace Anabatic {
int deltaH2 = (int)depthH2 - (int)depthV1;
if (std::abs(deltaH1) > std::abs(deltaH2)) {
setLayer( rg->getContactLayer( depthH2 + ((depthH2<depthV1)?0:-1) ) );
setLayerAndWidth( 2, depthH2 + ((depthH2<depthV1)?0:-1) );
//_horizontal1 = static_cast<AutoHorizontal*>( _horizontal1->makeDogleg(this) );
_horizontal1->makeDogleg(this);
cdebug_log(145,0) << "New h1:" << _horizontal1 << endl;
} else {
setLayer( rg->getContactLayer( depthH1 + ((depthH1<depthV1)?0:-1) ) );
setLayerAndWidth( 2, depthH1 + ((depthH1<depthV1)?0:-1) );
//_horizontal2 = static_cast<AutoHorizontal*>( _horizontal2->makeDogleg(this) );
_horizontal2->makeDogleg(this);
cdebug_log(145,0) << "New h2:" << _horizontal2 << endl;

View File

@ -550,8 +550,8 @@ namespace Anabatic {
cdebug_log(145,0) << "Update seg: " << _segment << endl;
delta = abssub( anchorDepth, rg->getLayerDepth( _segment->getLayer() ) );
}
else if (delta == 0) setLayer( rg->getRoutingLayer(anchorDepth) );
else if (delta == 1) setLayer( rg->getContactLayer(std::min(anchorDepth,segmentDepth)) );
else if (delta == 0) setLayerAndWidth( delta, anchorDepth );
else if (delta == 1) setLayerAndWidth( delta, std::min(anchorDepth,segmentDepth) );
}
_segment->invalidate( this );

View File

@ -263,7 +263,7 @@ namespace Anabatic {
delta = abssub ( depthH1, depthV1 );
}
setLayer ( (delta == 0) ? rg->getRoutingLayer(depthContact) : rg->getContactLayer(depthContact) );
setLayerAndWidth( delta, depthContact );
}
_horizontal1->invalidate( this );

View File

@ -278,11 +278,11 @@ namespace Anabatic {
cdebug_log(145,0) << "depthV1 == depthV2 (" << depthV1 << ")" << endl;
// Dogleg on the horizontal.
switch ( delta ) {
case 0: setLayer( rg->getRoutingLayer(minDepth) ); break;
case 1: setLayer( rg->getContactLayer(minDepth) ); break;
case 0:
case 1: setLayerAndWidth( delta, minDepth ); break;
default:
cdebug_log(145,0) << "Restore connectivity: dogleg on h1." << endl;
setLayer( rg->getContactLayer( depthV1 + ((depthV1==minDepth)?0:-1) ) );
setLayerAndWidth( delta, depthV1 + ((depthV1==minDepth)?0:-1) );
_horizontal1 = static_cast<AutoHorizontal*>( _horizontal1->makeDogleg(this) );
break;
}
@ -292,11 +292,11 @@ namespace Anabatic {
int deltaV2 = (int)depthV2 - (int)depthH1;
if (std::abs(deltaV1) > std::abs(deltaV2)) {
setLayer( rg->getContactLayer( depthV2 + ((depthV2<depthH1)?0:-1) ) );
setLayerAndWidth( 2, depthV2 + ((depthV2<depthH1)?0:-1) );
//_vertical1 = static_cast<AutoVertical*>( _vertical1->makeDogleg(this) );
_vertical1->makeDogleg(this);
} else {
setLayer( rg->getContactLayer( depthV1 + ((depthV1<depthH1)?0:-1) ) );
setLayerAndWidth( 2, depthV1 + ((depthV1<depthH1)?0:-1) );
//_vertical2 = static_cast<AutoVertical*>( _vertical2->makeDogleg(this) );
_vertical2->makeDogleg(this);
}

View File

@ -804,8 +804,9 @@ namespace Anabatic {
upLayer = (depth+1 <= Session::getConfiguration()->getAllowedDepth());
}
size_t doglegDepth = depth + ((upLayer)?1:-1);
Layer* contactLayer = Session::getRoutingGauge()->getContactLayer( depth + ((upLayer)?0:-1) );
const Layer* doglegLayer = Session::getRoutingGauge()->getRoutingLayer( depth + ((upLayer)?1:-1) );
const Layer* doglegLayer = Session::getRoutingGauge()->getRoutingLayer( doglegDepth );
Session::dogleg( this );
targetDetach();
@ -814,7 +815,7 @@ namespace Anabatic {
AutoContact* dlContact1 = AutoContactTurn::create( doglegGCell, _horizontal->getNet(), contactLayer );
AutoContact* dlContact2 = AutoContactTurn::create( doglegGCell, _horizontal->getNet(), contactLayer );
AutoSegment* segment1 = AutoSegment::create( dlContact1 , dlContact2, Flags::Vertical );
segment1->setLayer( doglegLayer );
segment1->setLayer( doglegDepth );
segment1->_setAxis( doglegAxis );
segment1->setFlags( SegDogleg|SegSlackened|SegCanonical|SegNotAligned );
@ -825,7 +826,7 @@ namespace Anabatic {
targetAttach( dlContact1 );
AutoSegment* segment2 = AutoSegment::create( dlContact2 , autoTarget, Flags::Horizontal );
autoTarget->cacheAttach( segment2 );
segment2->setLayer( getLayer() );
segment2->setLayer( depth );
segment2->_setAxis( getY() );
segment2->setFlags( (isSlackened()?SegSlackened:0) );
Session::dogleg( segment2 );

View File

@ -1632,7 +1632,7 @@ namespace Anabatic {
const Layer* newLayer = Session::getRoutingGauge()->getRoutingLayer(depth);
if (getLayer() != newLayer) {
setLayer( newLayer );
setLayer( depth );
getAutoSource()->invalidate( Flags::Topology|Flags::NoCheckLayer );
getAutoTarget()->invalidate( Flags::Topology|Flags::NoCheckLayer );
}
@ -2110,18 +2110,18 @@ namespace Anabatic {
if (getLayer()->above(from->getLayer())) {
cdebug_log(149,0) << "Go Down from depth " << segmentDepth << endl;
doglegs[ index + 1 ]->setLayer( rg->getRoutingLayer(segmentDepth-1) );
doglegs[ index + 1 ]->setLayer( segmentDepth-1 );
cdebug_log(149,0) << "doglegs[i+1]: " << doglegs[index+1] << endl;
if (isSource) {
doglegs[ index + 0 ]->setLayer( rg->getRoutingLayer(segmentDepth-2) );
doglegs[ index + 0 ]->setLayer( segmentDepth-2 );
doglegs[ index + 1 ]->getAutoSource()->setLayer( rg->getContactLayer(segmentDepth-2) );
doglegs[ index + 1 ]->getAutoTarget()->setLayer( rg->getContactLayer(segmentDepth-1) );
cdebug_log(149,0) << "doglegs[i+0]: " << doglegs[index+0] << endl;
cdebug_log(149,0) << "doglegs[i+1]: " << doglegs[index+1]->getAutoSource() << endl;
cdebug_log(149,0) << "doglegs[i+1]: " << doglegs[index+1]->getAutoTarget() << endl;
} else {
doglegs[ index + 2 ]->setLayer( rg->getRoutingLayer(segmentDepth-2) );
doglegs[ index + 2 ]->setLayer( segmentDepth-2 );
doglegs[ index + 1 ]->getAutoTarget()->setLayer( rg->getContactLayer(segmentDepth-2) );
doglegs[ index + 1 ]->getAutoSource()->setLayer( rg->getContactLayer(segmentDepth-1) );
cdebug_log(149,0) << "doglegs[i+2]: " << doglegs[index+2] << endl;
@ -2131,18 +2131,18 @@ namespace Anabatic {
} else {
cdebug_log(149,0) << "Go Up from depth " << segmentDepth << endl;
doglegs[ index + 1 ]->setLayer( rg->getRoutingLayer(segmentDepth+1) );
doglegs[ index + 1 ]->setLayer( segmentDepth+1 );
cdebug_log(149,0) << "doglegs[i+1]: " << doglegs[index+1] << endl;
if (isSource) {
doglegs[ index + 0 ]->setLayer( rg->getRoutingLayer(segmentDepth+2) );
doglegs[ index + 0 ]->setLayer( segmentDepth+2 );
doglegs[ index + 1 ]->getAutoSource()->setLayer( rg->getContactLayer(segmentDepth+1) );
doglegs[ index + 1 ]->getAutoTarget()->setLayer( rg->getContactLayer(segmentDepth ) );
cdebug_log(149,0) << "doglegs[i+0]: " << doglegs[index+0] << endl;
cdebug_log(149,0) << "doglegs[i+1]: " << doglegs[index+1]->getAutoSource() << endl;
cdebug_log(149,0) << "doglegs[i+1]: " << doglegs[index+1]->getAutoTarget() << endl;
} else {
doglegs[ index + 2 ]->setLayer( rg->getRoutingLayer(segmentDepth+2) );
doglegs[ index + 2 ]->setLayer( segmentDepth+2 );
doglegs[ index + 1 ]->getAutoTarget()->setLayer( rg->getContactLayer(segmentDepth+1) );
doglegs[ index + 1 ]->getAutoSource()->setLayer( rg->getContactLayer(segmentDepth ) );
cdebug_log(149,0) << "doglegs[i+2]: " << doglegs[index+2] << endl;

View File

@ -712,8 +712,9 @@ namespace Anabatic {
upLayer = (depth+1 <= Session::getConfiguration()->getAllowedDepth());
}
size_t doglegDepth = depth + ((upLayer)?1:-1);
Layer* contactLayer = Session::getRoutingGauge()->getContactLayer ( depth + ((upLayer)?0:-1) );
const Layer* doglegLayer = Session::getRoutingGauge()->getRoutingLayer ( depth + ((upLayer)?1:-1) );
const Layer* doglegLayer = Session::getRoutingGauge()->getRoutingLayer ( doglegDepth );
Session::dogleg( this );
targetDetach();
@ -725,7 +726,7 @@ namespace Anabatic {
cdebug_log(149,0) << dlContact2 << endl;
AutoSegment* segment1 = AutoSegment::create( dlContact1 , dlContact2, Flags::Horizontal );
cdebug_log(149,0) << segment1 << endl;
segment1->setLayer( doglegLayer );
segment1->setLayer( doglegDepth );
segment1->_setAxis( doglegAxis );
segment1->setFlags( SegDogleg|SegSlackened|SegCanonical|SegNotAligned );
cdebug_log(149,0) << "New " << dlContact1->base() << "." << endl;
@ -735,7 +736,7 @@ namespace Anabatic {
targetAttach( dlContact1 );
AutoSegment* segment2 = AutoVertical::create ( dlContact2, autoTarget, Flags::Vertical );
autoTarget->cacheAttach( segment2 );
segment2->setLayer( getLayer() );
segment2->setLayer( depth );
segment2->_setAxis( getX() );
segment2->setFlags( (isSlackened()?SegSlackened:0) );
Session::dogleg( segment2 );

View File

@ -8,6 +8,7 @@ endif ( CHECK_DETERMINISM )
${CORIOLIS_INCLUDE_DIR}
${HURRICANE_INCLUDE_DIR}
${CONFIGURATION_INCLUDE_DIR}
${FLUTE_INCLUDE_DIR}
${Boost_INCLUDE_DIRS}
${QtX_INCLUDE_DIR}
${PYTHON_INCLUDE_PATH}
@ -76,6 +77,7 @@ endif ( CHECK_DETERMINISM )
${CONFIGURATION_LIBRARY}
${CIF_LIBRARY}
${AGDS_LIBRARY}
${FLUTE_LIBRARIES}
${LEFDEF_LIBRARIES}
${OA_LIBRARIES}
${QtX_LIBRARIES}

View File

@ -117,8 +117,14 @@ namespace Anabatic {
}
global++;
if ((*isegment)->getLayer() == Session::getRoutingLayer(1)) (*isegment)->setLayer( Session::getRoutingLayer(3) );
if ((*isegment)->getLayer() == Session::getRoutingLayer(2)) (*isegment)->setLayer( Session::getRoutingLayer(4) );
if ((*isegment)->getLayer() == Session::getRoutingLayer(1)) {
(*isegment)->setLayer( Session::getRoutingLayer(3) );
(*isegment)->setWidth( Session::getWireWidth (3) );
}
if ((*isegment)->getLayer() == Session::getRoutingLayer(2)) {
(*isegment)->setLayer( Session::getRoutingLayer(4) );
(*isegment)->setWidth( Session::getWireWidth (4) );
}
}
}
@ -172,8 +178,14 @@ namespace Anabatic {
netGlobal++;
cdebug_log(145,0) << "Migrate to M4/M5: " << autoSegment << endl;
if (autoSegment->isHorizontal()) autoSegment->setLayer( Session::getRoutingLayer(3) );
if (autoSegment->isVertical ()) autoSegment->setLayer( Session::getRoutingLayer(4) );
if (autoSegment->isHorizontal()) {
autoSegment->setLayer( Session::getRoutingLayer(3) );
autoSegment->setWidth( Session::getWireWidth (3) );
}
if (autoSegment->isVertical()) {
autoSegment->setLayer( Session::getRoutingLayer(4) );
autoSegment->setWidth( Session::getWireWidth (4) );
}
}
}
}

View File

@ -178,6 +178,7 @@ namespace Anabatic {
, Flags flags=Flags::WarnOnError );
void restoreNativeConstraintBox ();
void migrateConstraintBox ( AutoContact* other );
void setLayerAndWidth ( size_t delta, size_t depth );
void destroy ();
// Inspector Management.
Record* _getRecord () const;

View File

@ -178,6 +178,8 @@ namespace Anabatic {
inline DbU::Unit getTargetY () const;
inline void invert ();
inline void setLayer ( const Layer* );
inline void setLayer ( size_t depth );
inline void setWidth ( DbU::Unit );
// Predicates.
inline bool isHorizontal () const;
inline bool isVertical () const;
@ -538,6 +540,7 @@ namespace Anabatic {
inline void AutoSegment::incReduceds () { if (_reduceds<3) ++_reduceds; }
inline void AutoSegment::decReduceds () { if (_reduceds>0) --_reduceds; }
inline void AutoSegment::setLayer ( const Layer* layer ) { base()->setLayer(layer); _depth=Session::getLayerDepth(layer); _flags|=SegInvalidatedLayer; }
inline void AutoSegment::setWidth ( DbU::Unit width ) { base()->setWidth(width); }
inline void AutoSegment::setOptimalMin ( DbU::Unit min ) { _optimalMin = (unsigned int)DbU::getLambda(min-getOrigin()); }
inline void AutoSegment::setOptimalMax ( DbU::Unit max ) { _optimalMax = (unsigned int)DbU::getLambda(max-getOrigin()); }
inline void AutoSegment::mergeNativeMin ( DbU::Unit min ) { _nativeConstraints.getVMin() = std::max( min, _nativeConstraints.getVMin() ); }
@ -547,6 +550,17 @@ namespace Anabatic {
inline void AutoSegment::resetUserConstraints () { _userConstraints = Interval(false); }
inline void AutoSegment::setLayer ( size_t depth )
{
RoutingLayerGauge* layerGauge = Session::getLayerGauge( depth );
base()->setLayer( layerGauge->getLayer () );
base()->setWidth( layerGauge->getWireWidth() );
_depth = depth;
_flags|=SegInvalidatedLayer;
}
inline DbU::Unit AutoSegment::getContactWidth () const
{ return getWidth() + Session::getViaWidth(getLayer()) - Session::getWireWidth(getLayer()); }

View File

@ -1,50 +0,0 @@
READ THIS LICENSE AGREEMENT CAREFULLY BEFORE USING THIS PRODUCT. BY USING
THIS PRODUCT YOU INDICATE YOUR ACCEPTANCE OF THE TERMS OF THE FOLLOWING
AGREEMENT. THESE TERMS APPLY TO YOU AND ANY SUBSEQUENT LICENSEE OF THIS
PRODUCT.
License Agreement for FLUTE
Copyright (c) 2004 by Dr. Chris C. N. Chu
All rights reserved
ATTRIBUTION ASSURANCE LICENSE (adapted from the original BSD license)
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the conditions below are
met. These conditions require a modest attribution to Dr. Chris C. N. Chu
(the "Author").
1. Redistributions of the source code, with or without modification (the
"Code"), must be accompanied by any documentation and, each time
the resulting executable program or a program dependent thereon is
launched, a prominent display (e.g., splash screen or banner text) of
the Author's attribution information, which includes:
(a) Dr. Chris C. N. Chu ("AUTHOR"),
(b) Iowa State University ("PROFESSIONAL IDENTIFICATION"), and
(c) http://home.eng.iastate.edu/~cnchu/ ("URL").
2. Users who intend to use the Code for commercial purposes will notify
Author prior to such commercial use.
3. Neither the name nor any trademark of the Author may be used to
endorse or promote products derived from this software without
specific prior written permission.
4. Users are entirely responsible, to the exclusion of the Author and any
other persons, for compliance with (1) regulations set by owners or
administrators of employed equipment, (2) licensing terms of any other
software, and (3) local, national, and international regulations
regarding use, including those regarding import, export, and use of
encryption software.
THIS FREE SOFTWARE IS PROVIDED BY THE AUTHOR "AS IS" AND ANY EXPRESS OR
IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
IN NO EVENT SHALL THE AUTHOR OR ANY CONTRIBUTOR BE LIABLE FOR ANY DIRECT,
INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, EFFECTS OF UNAUTHORIZED OR MALICIOUS
NETWORK ACCESS; PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

File diff suppressed because it is too large Load Diff

File diff suppressed because it is too large Load Diff

View File

@ -1,44 +0,0 @@
#include "knik/global.h"
/*********************************************************************/
/*
Return the Manhattan distance between two points
*/
long dist(
Point p,
Point q
)
{
long dx, dy;
dx = (p.x) - (q.x);
if( dx < 0 ) dx = -dx;
dy = (p.y) - (q.y);
if( dy < 0 ) dy = -dy;
return dx + dy;
}
/*********************************************************************/
/*
Return the Manhattan distance between two points
*/
long dist2(
Point* p,
Point* q
)
{
long dx, dy;
dx = (p->x) - (q->x);
if( dx < 0 ) dx = -dx;
dy = (p->y) - (q->y);
if( dy < 0 ) dy = -dy;
return dx + dy;
}
/*********************************************************************/
/*********************************************************************/

View File

@ -1,161 +0,0 @@
#include "knik/dl.h"
#include <assert.h>
#include <stdio.h>
dl_t dl_alloc()
{
dl_t dl = (dl_t)malloc(sizeof(dl_s));
if (!dl) {
printf("Out of memory!!\n");
} else {
dl->first = dl->last = 0; dl->count = 0;
}
return dl;
}
void dl_delete(dl_t dl, dl_el *el)
{
if (dl->first == el) {
dl->first = el->next;
}
if (dl->last == el) {
dl->last = el->prev;
}
if (el->next) {
el->next->prev = el->prev;
}
if (el->prev) {
el->prev->next = el->next;
}
free(el); dl->count--;
}
void dl_clear(dl_t dl)
{
dl_el *el, *next;
if (dl->count > 0) {
for (el=dl->first; el; el=next) {
next = el->next;
free(el);
}
}
dl->first = dl->last = 0;
dl->count = 0;
}
void dl_concat(dl_t first_list, dl_t second_list)
{
if (first_list->count <= 0) {
*first_list = *second_list;
} else if (second_list->count > 0) {
first_list->last->next = second_list->first;
second_list->first->prev = first_list->last;
first_list->last = second_list->last;
first_list->count += second_list->count;
}
free(second_list);
}
static void dl_insertion_sort(dl_t dl, size_t el_size,
int(*compar)(void *, void *))
{
char *buf;
void *curr_d, *srch_d;
dl_el *curr, *srch;
if (dl_length(dl) <= 1) {
return;
}
buf = (char*)malloc(el_size);
for (curr=dl->first; curr!=dl->last; curr=curr->next) {
curr_d = (void*)(((dl_el*)curr)+1);
for (srch=dl->last; srch!=curr; srch=srch->prev) {
srch_d = (void*)(((dl_el*)srch)+1);
if (compar(curr_d, srch_d) > 0) {
memcpy((void*)buf, curr_d, el_size);
memcpy(curr_d, srch_d, el_size);
memcpy(srch_d, (void*)buf, el_size);
}
}
}
free(buf);
}
void dl_sort(dl_t dl, size_t el_size, int(*compar)(void *, void *))
{
dl_el *el, *first_head, *second_head;
dl_s first_list, second_list;
void *first_item, *second_item;
int i, len;
if (dl_length(dl) <= 25) {
dl_insertion_sort(dl, el_size, compar);
return;
}
len = dl_length(dl)/2;
for (i=0, el=dl->first; i<len; i++) {
el = el->next;
}
first_list.first = dl->first;
first_list.last = el->prev;
first_list.count = len;
first_list.last->next = 0;
second_list.first = el;
second_list.last = dl->last;
second_list.count = dl_length(dl)-len;
second_list.first->prev = 0;
dl_sort(&first_list, el_size, compar);
dl_sort(&second_list, el_size, compar);
/* in-place merging */
first_head = first_list.first;
second_head = second_list.first;
first_item = (void*)(((dl_el*)first_head)+1);
second_item = (void*)(((dl_el*)second_head)+1);
if (compar(first_item, second_item) <= 0) {
dl->first = el = first_head;
first_head = first_head->next;
} else {
dl->first = el = second_head;
second_head = second_head->next;
}
while (1) {
first_item = (void*)(((dl_el*)first_head)+1);
second_item = (void*)(((dl_el*)second_head)+1);
if (compar(first_item, second_item) <= 0) {
el->next = first_head;
first_head->prev = el;
el = first_head;
first_head = first_head->next;
if (!first_head) {
el->next = second_head;
second_head->prev = el;
dl->last = second_list.last;
break;
}
} else {
el->next = second_head;
second_head->prev = el;
el = second_head;
second_head = second_head->next;
if (!second_head) {
el->next = first_head;
first_head->prev = el;
dl->last = first_list.last;
break;
}
}
}
}

View File

@ -1,28 +0,0 @@
#include <stdio.h>
#include <stdlib.h>
/**************************************************************************/
/*
print error message and continue
*/
void err_msg(
char* msg
)
{
fprintf(stderr, "%s\n", msg);
}
/**************************************************************************/
/*
print error message and exit
*/
void err_exit(
char* msg
)
{
fprintf(stderr, "%s\n", msg);
exit(1);
}

File diff suppressed because it is too large Load Diff

File diff suppressed because it is too large Load Diff

View File

@ -1,177 +0,0 @@
/****************************************************************************/
/*
Binary heap routines for use in Prim's algorithm,
with points are numbered from 0 to n-1
*/
#include <stdlib.h>
#include "knik/heap.h"
#include "knik/err.h"
Heap* _heap = (Heap*)NULL;
long _max_heap_size = 0;
long _heap_size = 0;
/****************************************************************************/
/*
*/
void allocate_heap( long n )
{
if( _max_heap_size < n )
{
_heap = (Heap*)realloc( (void*)_heap, (size_t)(n+1)*sizeof(Heap) );
if( ! _heap )
{
err_exit( (char*)"Cannot reallocate memory in allocate_heap!" );
}
_max_heap_size = n;
}
}
/****************************************************************************/
/*
*/
void deallocate_heap()
{
_max_heap_size = 0;
if( _heap )
{
free( (void*)_heap );
_heap = (Heap*)NULL;
}
}
/****************************************************************************/
void heap_init( long n )
{
register long p;
allocate_heap( n );
_heap_size = 0;
for( p = 0; p < n; p++ )
{
heap_idx( p ) = 0;
}
} /* END heap_init() */
/****************************************************************************/
void heap_insert(
long p,
long key
)
{
register long k; /* hole in the heap */
register long j; /* parent of the hole */
register long q; /* heap_elt(j) */
heap_key( p ) = key;
if( _heap_size == 0 )
{
_heap_size = 1;
heap_elt( 1 ) = p;
heap_idx( p ) = 1;
return;
}
k = ++ _heap_size;
j = k >> 1; /* k/2 */
while( (j > 0) && (heap_key(q=heap_elt(j)) > key) ) {
heap_elt( k ) = q;
heap_idx( q ) = k;
k = j;
j = k>>1; /* k/2 */
}
/* store p in the position of the hole */
heap_elt( k ) = p;
heap_idx( p ) = k;
} /* END heap_insert() */
/****************************************************************************/
void heap_decrease_key
(
long p,
long new_key
)
{
register long k; /* hole in the heap */
register long j; /* parent of the hole */
register long q; /* heap_elt(j) */
heap_key( p ) = new_key;
k = heap_idx( p );
j = k >> 1; /* k/2 */
if( (j > 0) && (heap_key(q=heap_elt(j)) > new_key) ) { /* change is needed */
do {
heap_elt( k ) = q;
heap_idx( q ) = k;
k = j;
j = k>>1; /* k/2 */
} while( (j > 0) && (heap_key(q=heap_elt(j)) > new_key) );
/* store p in the position of the hole */
heap_elt( k ) = p;
heap_idx( p ) = k;
}
} /* END heap_decrease_key() */
/****************************************************************************/
long heap_delete_min()
{
long min, last;
register long k; /* hole in the heap */
register long j; /* child of the hole */
register long l_key; /* key of last point */
if( _heap_size == 0 ) /* heap is empty */
return( -1 );
min = heap_elt( 1 );
last = heap_elt( _heap_size -- );
l_key = heap_key( last );
k = 1; j = 2;
while( j <= _heap_size ) {
if( heap_key(heap_elt(j)) > heap_key(heap_elt(j+1)) )
j++;
if( heap_key(heap_elt(j)) >= l_key)
break; /* found a position to insert 'last' */
/* else, sift hole down */
heap_elt(k) = heap_elt(j); /* Note that j <= _heap_size */
heap_idx( heap_elt(k) ) = k;
k = j;
j = k << 1;
}
heap_elt( k ) = last;
heap_idx( last ) = k;
heap_idx( min ) = -1; /* mark the point visited */
return( min );
} /* END heap_delete_min() */
/****************************************************************************/

View File

@ -1,16 +0,0 @@
#ifndef _DIST_H_
#define _DIST_H_
#include "global.h"
long dist(
Point p,
Point q
);
long dist2(
Point* p,
Point* q
);
#endif

View File

@ -1,180 +0,0 @@
#ifndef DL_H
#define DL_H
#include <string.h>
#include <stdlib.h>
typedef struct dl_el_s {
struct dl_el_s *prev, *next;
} dl_el;
typedef struct {
dl_el *first, *last;
unsigned int count;
} dl_s, *dl_t;
dl_t dl_alloc(void);
void dl_delete(dl_t dl, dl_el *el);
void dl_clear(dl_t dl);
void dl_concat(dl_t list1, dl_t list2);
void dl_sort(dl_t dl, size_t el_size, int(*compar)(void *, void *));
#define dl_length(dl) (dl)->count
#define dl_empty(dl) ((dl)->count <= 0)
#define dl_data(type, el) \
*(type*)(((dl_el*)(el))+1)
#define dl_data_p(type, el) \
((type*)(((dl_el*)(el))+1))
#define dl_forall(type, dl, data) \
{ \
dl_el *_el, *_next; \
dl_t _curr_dl = (dl); \
for (_el=_curr_dl->first; _el; _el=_next) { \
_next = _el->next; \
(data) = dl_data(type, _el);
#define dl_forall_p(type, dl, data_p) \
{ \
dl_el *_el, *_next; \
dl_t _curr_dl = (dl); \
for (_el=_curr_dl->first; _el; _el=_next) { \
_next = _el->next; \
(data_p) = dl_data_p(type, _el);
#define dl_current() _el
#define dl_delete_current() dl_delete(_curr_dl, _el)
#define dl_endfor \
} \
}
#define dl_forall_reverse(type, dl, data) \
{ \
dl_el *_el, *_next; \
dl_t _curr_dl = (dl); \
for (_el=_curr_dl->last; _el; _el=_next) { \
_next = _el->prev; \
(data) = dl_data(type, _el);
#define dl_forall_reverse_p(type, dl, data_p) \
{ \
dl_el *_el, *_next; \
dl_t _curr_dl = (dl); \
for (_el=_curr_dl->last; _el; _el=_next) { \
_next = _el->prev; \
(data_p) = dl_data_p(type, _el);
#define dl_first(type, dl) \
dl_data(type, (dl)->first)
#define dl_first_element(dl) (dl)->first
#define dl_last(type, dl) \
dl_data(type, (dl)->last)
#define dl_pop_first(type, dl, data) \
{ \
(data) = dl_first(type, dl); \
dl_delete((dl), (dl)->first); \
}
#define dl_pop_last(type, dl, data) \
{ (data) = dl_last(type, dl); dl_delete((dl), (dl)->last); }
#define dl_insert_before(type, dl, element, data) \
{ \
if ((element) == (dl)->first) { \
dl_prepend(type, dl, data); \
} else { \
dl_el *_el = (dl_el*) malloc(sizeof(dl_el)+sizeof(type)); \
if (!_el) { \
printf("Out of memory!!\n"); \
} else { \
memcpy(_el+1, &(data), sizeof(type)); \
_el->prev = (element)->prev; _el->next = (element); \
(element)->prev->next = _el; (element)->prev = _el; \
(dl)->count++; \
} \
} \
}
#define dl_insert_after(type, dl, element, data) \
{ \
if ((element) == (dl)->last) { \
dl_append(type, dl, data); \
} else { \
dl_el *_el = (dl_el*) malloc(sizeof(dl_el)+sizeof(type)); \
if (!_el) { \
printf("Out of memory!!\n"); \
} else { \
memcpy(_el+1, &(data), sizeof(type)); \
_el->next = (element)->next; _el->prev = (element); \
(element)->next->prev = _el; (element)->next = _el; \
(dl)->count++; \
} \
} \
}
#define dl_append(type, dl, data) \
{ \
dl_el *_el = (dl_el*) malloc(sizeof(dl_el)+sizeof(type)); \
if (!_el) { \
printf("Out of memory!!\n"); \
} else { \
memcpy(_el+1, &(data), sizeof(type)); \
_el->next = 0; \
if ((dl)->count <= 0) { \
_el->prev = 0; \
(dl)->first = (dl)->last = _el; \
(dl)->count = 1; \
} else { \
_el->prev = (dl)->last; \
(dl)->last->next = _el; \
(dl)->last = _el; \
(dl)->count++; \
} \
} \
}
#define dl_prepend(type, dl, data) \
{ \
dl_el *_el = (dl_el*) malloc(sizeof(dl_el)+sizeof(type)); \
if (!_el) { \
printf("Out of memory!!\n"); \
} else { \
memcpy(_el+1, &(data), sizeof(type)); \
_el->prev = 0; \
if ((dl)->count <= 0) { \
_el->next = 0; \
(dl)->first = (dl)->last = _el; \
(dl)->count = 1; \
} else { \
_el->next = (dl)->first; \
(dl)->first->prev = _el; \
(dl)->first = _el; \
(dl)->count++; \
} \
} \
}
#define dl_free(dl) \
{ \
dl_clear(dl); free(dl); dl = 0; \
}
#define dl_duplicate(dest, src, type) \
{ \
dest = dl_alloc(); \
type _data_el; \
dl_forall(type, src, _data_el) { \
dl_append(type, dest, _data_el); \
} dl_endfor; \
}
#endif

View File

@ -1,12 +0,0 @@
#ifndef _ERR_H_
#define _ERR_H_
void err_msg(
char* msg
);
void err_exit(
char* msg
);
#endif

View File

@ -1,99 +0,0 @@
#ifndef _KNIK_FLUTE_H
#define _KNIK_FLUTE_H
/*****************************/
/* User-Defined Parameters */
/*****************************/
#define MAXD 150 // max. degree that can be handled
#define ACCURACY 3 // Default accuracy
#define ROUTING 1 // 1 to construct routing, 0 to estimate WL only
#define LOCAL_REFINEMENT 0 // Suggestion: Set to 1 if ACCURACY >= 5
#define REMOVE_DUPLICATE_PIN 0 // Remove dup. pin for flute_wl() & flute()
#ifndef DTYPE // Data type for distance
#define DTYPE int
#endif
/*****************************/
/* User-Callable Functions */
/*****************************/
// void readLUT();
// DTYPE flute_wl(int d, DTYPE x[], DTYPE y[], int acc);
// DTYPE flutes_wl(int d, DTYPE xs[], DTYPE ys[], int s[], int acc);
// FTree flute(int d, DTYPE x[], DTYPE y[], int acc);
// FTree flutes(int d, DTYPE xs[], DTYPE ys[], int s[], int acc);
// DTYPE wirelength(FTree t);
// void printtree(FTree t);
// void plottree(FTree t);
/*************************************/
/* Internal Parameters and Functions */
/*************************************/
#define POWVFILE "POWV9.dat" // LUT for POWV (Wirelength Vector)
#define POSTFILE "POST9.dat" // LUT for POST (Steiner FTree)
#define D 9 // LUT is used for d <= D, D <= 9
#define TAU(A) (8+1.3*(A))
#define D1(A) (25+120/((A)*(A))) // flute_mr is used for D1 < d <= D2
#define D2(A) ((A)<=6 ? 500 : 75+5*(A))
typedef struct
{
DTYPE x, y; // starting point of the branch
int n; // index of neighbor
} Branch;
struct FTree
{
int deg; // degree
DTYPE length; // total wirelength
Branch *branch; // array of tree branches
};
// User-Callable Functions
extern void readLUT();
extern DTYPE flute_wl(int d, DTYPE x[], DTYPE y[], int acc);
//Macro: DTYPE flutes_wl(int d, DTYPE xs[], DTYPE ys[], int s[], int acc);
extern FTree flute(int d, DTYPE x[], DTYPE y[], int acc);
//Macro: FTree flutes(int d, DTYPE xs[], DTYPE ys[], int s[], int acc);
extern DTYPE wirelength(FTree t);
extern void printtree(FTree t);
extern void plottree(FTree t);
// Other useful functions
extern void init_param();
extern DTYPE flutes_wl_LD(int d, DTYPE xs[], DTYPE ys[], int s[]);
extern DTYPE flutes_wl_MD(int d, DTYPE xs[], DTYPE ys[], int s[], int acc);
extern DTYPE flutes_wl_RDP(int d, DTYPE xs[], DTYPE ys[], int s[], int acc);
extern FTree flutes_LD(int d, DTYPE xs[], DTYPE ys[], int s[]);
extern FTree flutes_MD(int d, DTYPE xs[], DTYPE ys[], int s[], int acc);
extern FTree flutes_HD(int d, DTYPE xs[], DTYPE ys[], int s[], int acc);
extern FTree flutes_RDP(int d, DTYPE xs[], DTYPE ys[], int s[], int acc);
#if REMOVE_DUPLICATE_PIN==1
#define flutes_wl(d, xs, ys, s, acc) flutes_wl_RDP(d, xs, ys, s, acc)
#define flutes(d, xs, ys, s, acc) flutes_RDP(d, xs, ys, s, acc)
#else
#define flutes_wl(d, xs, ys, s, acc) flutes_wl_ALLD(d, xs, ys, s, acc)
#define flutes(d, xs, ys, s, acc) flutes_ALLD(d, xs, ys, s, acc)
#endif
#define flutes_wl_ALLD(d, xs, ys, s, acc) flutes_wl_LMD(d, xs, ys, s, acc)
#define flutes_ALLD(d, xs, ys, s, acc) \
(d<=D ? flutes_LD(d, xs, ys, s) \
: (d<=D1(acc) ? flutes_MD(d, xs, ys, s, acc) \
: flutes_HD(d, xs, ys, s, acc)))
#define flutes_wl_LMD(d, xs, ys, s, acc) \
(d<=D ? flutes_wl_LD(d, xs, ys, s) : flutes_wl_MD(d, xs, ys, s, acc))
#define flutes_LMD(d, xs, ys, s, acc) \
(d<=D ? flutes_LD(d, xs, ys, s) : flutes_MD(d, xs, ys, s, acc))
#define max(x,y) ((x)>(y)?(x):(y))
#define min(x,y) ((x)<(y)?(x):(y))
#define abs(x) ((x)<0?(-x):(x))
#define ADIFF(x,y) ((x)>(y)?(x-y):(y-x)) // Absolute difference
#endif

View File

@ -1,19 +0,0 @@
#ifndef _GLOBAL_H_
#define _GLOBAL_H_
#include <stdio.h>
#define TRUE 1
#define FALSE 0
#define MAXLONG 0x7fffffffL
struct point
{
long x, y;
};
typedef struct point Point;
typedef long nn_array[8];
#endif /* _GLOBAL_H_ */

View File

@ -1,31 +0,0 @@
#ifndef _HEAP_H_
#define _HEAP_H_
#include "global.h"
struct heap_info
{
long key;
long idx;
long elt;
};
typedef struct heap_info Heap;
extern Heap* _heap;
#define heap_key( p ) ( _heap[p].key )
#define heap_idx( p ) ( _heap[p].idx )
#define heap_elt( k ) ( _heap[k].elt )
#define in_heap( p ) ( heap_idx(p) > 0 )
#define never_seen( p ) ( heap_idx(p) == 0 )
void allocate_heap( long n );
void deallocate_heap();
void heap_init( long n );
void heap_insert( long p, long key );
void heap_decrease_key( long p, long new_key );
long heap_delete_min();
#endif /* _HEAP_H_ */

View File

@ -1,11 +0,0 @@
#ifndef _MST2_H_
#define _MST2_H_
#include "global.h"
void mst2_package_init( long n );
void mst2_package_done();
void mst2( long n, Point* pt, long* parent );
#endif

View File

@ -1,19 +0,0 @@
#include "global.h"
void allocate_nn_arrays( long n );
void deallocate_nn_arrays();
void brute_force_nearest_neighbors
(
long n,
Point* pt,
nn_array* nn
);
void dq_nearest_neighbors
(
long n,
Point* pt,
nn_array* nn
);

View File

@ -1,92 +0,0 @@
#include <stdlib.h>
#include <stdio.h>
#include <assert.h>
#include "knik/global.h"
#include "knik/neighbors.h"
#include "knik/dist.h"
#include "knik/heap.h"
#include "knik/err.h"
void mst2_package_init( long n )
{
allocate_heap( n );
allocate_nn_arrays( n );
}
/****************************************************************************/
/*
*/
void mst2_package_done()
{
deallocate_heap();
deallocate_nn_arrays();
}
/****************************************************************************/
/*
*/
void mst2
(
long n,
Point* pt,
long* parent
)
{
long i, k, nn1;
long d;
long oct;
long root = 0;
extern nn_array* nn;
// brute_force_nearest_neighbors( n, pt, nn );
dq_nearest_neighbors( n, pt, nn );
/*
Binary heap implementation of Prim's algorithm.
Runs in O(n*log(n)) time since at most 8n edges are considered
*/
heap_init( n );
heap_insert( root, 0 );
parent[root] = root;
for( k = 0; k < n; k++ ) /* n points to be extracted from heap */
{
i = heap_delete_min();
if (i<0) break;
#ifdef DEBUG
assert( i >= 0 );
#endif
/*
pt[i] entered the tree, update heap keys for its neighbors
*/
for( oct = 0; oct < 8; oct++ )
{
nn1 = nn[i][oct];
if( nn1 >= 0 )
{
d = dist( pt[i], pt[nn1] );
if( in_heap(nn1) && (d < heap_key(nn1)) )
{
heap_decrease_key( nn1, d );
parent[nn1] = i;
}
else if( never_seen(nn1) )
{
heap_insert( nn1, d );
parent[nn1] = i;
}
}
}
}
}
/****************************************************************************/
/****************************************************************************/

View File

@ -1,527 +0,0 @@
#include <assert.h>
#include <string.h>
#include <stdlib.h>
#include "knik/global.h"
#include "knik/err.h"
#include "knik/dist.h"
long octant
(
Point from,
Point to
);
static Point* _pt;
/***************************************************************************/
/*
For efficiency purposes auxiliary arrays are allocated as globals
*/
long max_arrays_size = 0;
nn_array* nn = (nn_array*)NULL;
Point* sheared = (Point*)NULL;
long* sorted = (long*)NULL;
long* aux = (long*)NULL;
/***************************************************************************/
/*
resize the auxiliary arrays to fit the specified number of points
*/
void allocate_nn_arrays( long n )
{
if( max_arrays_size < n )
{
nn = (nn_array*)realloc( (void*)nn, (size_t)n*sizeof(nn_array) );
sheared = (Point*)realloc( (void*)sheared, (size_t)n*sizeof(Point) );
sorted = (long*)realloc( (void*)sorted, (size_t)n*sizeof(long) );
aux = (long*)realloc( (void*)aux, (size_t)n*sizeof(long) );
if( !nn || !sheared || !sorted || !aux )
{
err_exit( (char*)"Cannot allocate memory in allocate_nn_arrays!" );
}
max_arrays_size = n;
}
}
/***************************************************************************/
/*
free memory used by auxiliary arrays
*/
void deallocate_nn_arrays()
{
max_arrays_size = 0;
if( nn )
{
free( (void*)nn );
nn = (nn_array*)NULL;
}
if( sheared )
{
free( (void*)sheared );
sheared = (Point*)NULL;
}
if( sorted )
{
free( (void*)sorted );
sorted = (long*)NULL;
}
if( aux )
{
free( (void*)aux );
aux = (long*)NULL;
}
}
/***************************************************************************/
/*
comparison function for use in quicksort
*/
static int compare_x
(
const void* i,
const void* j
)
{
/*
points with the same x must appear in increasing order of y
*/
if( sheared[*((long*)i)].x == sheared[*((long*)j)].x)
{
return sheared[*((long*)i)].y - sheared[*((long*)j)].y;
}
else
{
return sheared[*((long*)i)].x - sheared[*((long*)j)].x;
}
}
/***************************************************************************/
/*
Combine step of the Guibas-Stolfi divide-and-conquer NE nearest neighbor
algorithm. For efficiency purposes SW nearest neighbors are computed
at the same time.
*/
void ne_sw_combine
(
long left,
long mid,
long right,
Point* pt,
long* sorted,
long* aux,
long oct,
nn_array* nn
)
{
long i, j, k, y2;
long i1;
long i2;
long best_i2; /* index of current best nearest-neighbor */
long best_dist; /* distance to best nearest-neighbor */
long d;
#ifdef DEBUG
assert( right > mid );
assert( mid > left );
#endif
/*
update north-east nearest neighbors accross the mid-line
*/
i1 = left;
i2 = mid; y2 = pt[ sorted[i2] ].y;
while( (i1 < mid) && (pt[ sorted[i1] ].y >= y2) )
{
i1++;
}
if( i1 < mid )
{
best_i2 = i2;
best_dist = dist2( pt + sorted[i1], pt + sorted[best_i2] );
i2++;
while( (i1 < mid) && (i2 < right) )
{
if( pt[ sorted[i1] ].y < pt[ sorted[i2] ].y )
{
d = dist2( pt + sorted[i1], pt + sorted[i2] );
if( d < best_dist )
{
best_i2 = i2;
best_dist = d;
}
i2++;
}
else
{
if( (nn[ sorted[i1] ][oct] == -1) ||
( best_dist < dist2( pt + sorted[i1], pt + nn[ sorted[i1] ][oct]) )
)
{
nn[ sorted[i1] ][oct] = sorted[best_i2];
}
i1++;
if( i1 < mid )
{
best_dist = dist2( pt + sorted[i1], pt + sorted[best_i2] );
}
}
}
while( i1 < mid )
{
if( (nn[ sorted[i1] ][oct] == -1) ||
( dist2( pt + sorted[i1], pt + sorted[best_i2] ) <
dist2( pt + sorted[i1], pt + nn[ sorted[i1] ][oct]) )
)
{
nn[ sorted[i1] ][oct] = sorted[best_i2];
}
i1++;
}
}
/*
repeat for south-west nearest neighbors
*/
oct = (oct + 4) % 8;
i1 = right - 1;
i2 = mid - 1; y2 = pt[ sorted[i2] ].y;
while( (i1 >= mid) && (pt[ sorted[i1] ].y <= y2) )
{
i1--;
}
if( i1 >= mid )
{
best_i2 = i2;
best_dist = dist2( pt + sorted[i1], pt + sorted[best_i2] );
i2--;
while( (i1 >= mid) && (i2 >= left) )
{
if( pt[ sorted[i1] ].y > pt[ sorted[i2] ].y )
{
d = dist2( pt + sorted[i1], pt + sorted[i2] );
if( d < best_dist )
{
best_i2 = i2;
best_dist = d;
}
i2--;
}
else
{
if( (nn[ sorted[i1] ][oct] == -1) ||
( best_dist < dist2( pt + sorted[i1], pt + nn[ sorted[i1] ][oct]) )
)
{
nn[ sorted[i1] ][oct] = sorted[best_i2];
}
i1--;
if( i1 >= mid )
{
best_dist = dist2( pt + sorted[i1], pt + sorted[best_i2] );
}
}
}
while( i1 >= mid )
{
if( (nn[ sorted[i1] ][oct] == -1) ||
( dist2( pt + sorted[i1], pt + sorted[best_i2] ) <
dist2( pt + sorted[i1], pt + nn[ sorted[i1] ][oct]) )
)
{
nn[ sorted[i1] ][oct] = sorted[best_i2];
}
i1--;
}
}
/*
merge sorted[left..mid-1] with sorted[mid..right-1] by y-coordinate
*/
i = left; /* first unprocessed element in left list */
j = mid; /* first unprocessed element in right list */
k = left; /* first free available slot in output list */
while( (i < mid) && (j < right) )
{
if( pt[ sorted[i] ].y >= pt[ sorted[j] ].y )
{
aux[k++] = sorted[i++];
}
else
{
aux[k++] = sorted[j++];
}
}
/*
copy leftovers
*/
while( i < mid ) { aux[k++] = sorted[i++]; }
while( j < right ) { aux[k++] = sorted[j++]; }
/*
now copy sorted points from 'aux' to 'sorted'
*/
for( i = left; i < right; i++ ) { sorted[i] = aux[i]; }
#if 0
memcpy( (void*)(sorted+left), /* destination */
(void*)(aux+left), /* source */
(size_t)(right-left)*sizeof(long) /* number of bytes */
);
#endif
}
/***************************************************************************/
/*
compute north-east and south-west nearest neighbors for points indexed
by {sorted[left],...,sorted[right-1]}
*/
void ne_sw_nearest_neighbors
(
long left,
long right,
Point* pt,
long* sorted,
long* aux,
long oct,
nn_array* nn
)
{
long mid;
#ifdef DEBUG
assert( right > left );
#endif
if( right == left + 1 )
{
nn[ sorted[left] ][oct] = nn[ sorted[left]][(oct+4) % 8] = -1;
}
else
{
mid = (left + right) / 2;
ne_sw_nearest_neighbors( left, mid, pt, sorted, aux, oct, nn );
ne_sw_nearest_neighbors( mid, right, pt, sorted, aux, oct, nn );
ne_sw_combine( left, mid, right, pt, sorted, aux, oct, nn );
}
}
/***************************************************************************/
/*
Guibas-Stolfi algorithm for computing nearest NE neighbors
*/
void dq_nearest_neighbors
(
long n,
Point* pt,
nn_array* nn
)
{
long i, oct;
void check_nn( long, Point*, nn_array* );
long shear[4][4] = {
{1, -1, 0, 2},
{2, 0, -1, 1},
{1, 1, -2, 0},
{0, 2, -1, -1}
};
_pt = pt;
for( oct = 0; oct < 4; oct++ )
{
for( i = 0; i < n; i++ )
{
sheared[i].x = shear[oct][0]*pt[i].x + shear[oct][1]*pt[i].y;
sheared[i].y = shear[oct][2]*pt[i].x + shear[oct][3]*pt[i].y;
sorted[i] = i;
}
qsort( sorted, n, sizeof(long), compare_x );
ne_sw_nearest_neighbors( 0, n, sheared, sorted, aux, oct, nn );
}
#ifdef DEBUG
check_nn( n, pt, nn );
#endif
}
/***************************************************************************/
/***************************************************************************/
/*
Brute-force nearest-neighbor computation for debugging purposes
*/
/***************************************************************************/
/*
Half-open octants are numbered from 0 to 7 in anti-clockwise order
starting from ( dx >= dy > 0 ).
*/
#define sgn(x) ( x>0 ? 1 : (x < 0 ? -1 : 0) )
long octant
(
Point from,
Point to
)
{
long dx = to.x - from.x;
long dy = to.y - from.y;
long sgn1 = sgn(dx)*sgn(dy);
long sgn2 = sgn(dx+dy)*sgn(dx-dy);
long oct = 0x0;
if( (dy < 0) || ((dy==0) && (dx>0)) ) oct += 4;
if( (sgn1 < 0) || (dy==0) ) oct += 2;
if( (sgn1*sgn2 < 0) || (dy==0) || (dx==0) ) oct += 1;
return oct;
}
/***************************************************************************/
/*
O(n^2) algorithm for computing all nearest neighbors
*/
void brute_force_nearest_neighbors
(
long n,
Point* pt,
nn_array* nn
)
{
long i, j, oct;
long d;
/*
compute nearest neighbors by inspecting all pairs of points
*/
for( i = 0; i < n; i++ )
{
for( oct = 0; oct < 8; oct++ )
{
nn[i][oct] = -1;
}
}
for( i = 0; i < n; i++ )
{
for( j = i+1; j < n; j++ )
{
d = dist(pt[i], pt[j]);
oct = octant( pt[i], pt[j] );
if( ( nn[i][oct] == -1 ) ||
( d < dist(pt[i], pt[ nn[i][oct] ]) )
)
{
nn[i][oct] = j;
}
oct = (oct + 4) % 8;
if( ( nn[j][oct] == -1 ) ||
( d < dist(pt[j], pt[ nn[j][oct] ]) )
)
{
nn[j][oct] = i;
}
}
}
}
/***************************************************************************/
/*
compare nearest neighbors against those computed by brute force
*/
void check_nn
(
long n,
Point* pt,
nn_array* nn
)
{
long i, j, oct;
nn_array* nn1;
nn1 = (nn_array*)calloc( (size_t)n, (size_t)sizeof(nn_array) );
brute_force_nearest_neighbors( n, pt, nn1 );
for( i = 0; i < n; i++ )
{
for( oct = 0; oct < 8; oct++ )
{
if( nn[i][oct] == -1 )
{
assert( nn1[i][oct] == -1 );
}
else
{
assert( nn1[i][oct] != -1 );
if( octant(pt[i], pt[ nn[i][oct] ]) != oct )
{
printf( "WRONG OCTANT!\noct=%ld\n", oct );
printf( "i=%ld, x=%ld, y=%ld\n", i, pt[i].x, pt[i].y );
j = nn[i][oct];
printf( "nn=%ld, x=%ld, y=%ld, dist = %ld\n", j, pt[j].x, pt[j].y,
dist(pt[i], pt[j ]) );
}
// assert( octant(pt[i], pt[ nn[i][oct] ]) == oct );
assert( octant(pt[i], pt[ nn1[i][oct] ]) == oct );
if( dist(pt[i], pt[ nn[i][oct] ]) !=
dist(pt[i], pt[ nn1[i][oct] ]) )
{
printf( "NNs DON'T MATCH!\noct=%ld\n", oct );
printf( "i=%ld, x=%ld, y=%ld\n", i, pt[i].x, pt[i].y );
j = nn[i][oct];
printf( "nn=%ld, x=%ld, y=%ld, dist = %ld\n", j, pt[j].x, pt[j].y,
dist(pt[i], pt[j ]) );
j = nn1[i][oct];
printf( "nn1=%ld, x=%ld, y=%ld, dist = %ld\n", j, pt[j].x, pt[j].y,
dist(pt[i], pt[ j ]) );
}
// assert( dist(pt[i], pt[ nn[i][oct] ]) ==
// dist(pt[i], pt[ nn1[i][oct] ]) );
}
}
}
free( nn1 );
}
/***************************************************************************/
/***************************************************************************/