TweetFollow Us on Twitter

Apr 99 Challenge

Volume Number: 15 (1999)
Issue Number: 4
Column Tag: Programmer's Challenge

Apr 99 Challenge

by Bob Boonstra, Westford, MA

Shortest Network

This month's problem was suggested by Michael Kennedy, who wins two Challenge points for making the suggestion. The problem is to find the shortest network of line segments interconnecting a specified set of points. Shortest network algorithms have obvious practical application in constructing transportation and communications networks. In a January 1989, Scientific American article, Marshall Bern and Ronald Graham discussed the shortest network "Steiner" problem as one of a class of NP-hard problems. While no polynomial-time algorithm is known, the article (which, unfortunately, I have not been able to find online) discusses practical algorithms that produce networks slightly longer than the optimal one. Your Challenge for this month is to produce a near-optimal network in minimum time. Fortunately, we have been granted unlimited power of eminent domain, so there are no restrictions on where intermediate nodes may be placed or where connections may be routed.

The prototype for the code you should write is:

#if defined(__cplusplus)
extern "C" {
#endif

typedef struct Node {   /* node coordinates */
   double x;
   double y;
} Node;

typedef struct Connection {
         /* connection between Node[index1] and Node[index2] */
   long index1;
   long index2;
} Connection;

long /* numConnections */ ShortestNetwork(
   long numInitialNodes,         /* number of nodes to connect */
   long *numIntermediateNodes,   /* number of nodes added by ShortestNetwork */
   Node nodes[], 
      /* Nodes 0..numInitialNodes-1 are initialized on entry. */
      /* Nodes numInitialNodes..numInitialNodes+*numIntermediateNodes 
                  are added by ShortestNetwork */
   Connection connections[],   /* connections between nodes */
   long maxNodes,        /* number of entries allocated for nodes */
   long maxConnections   /* number of entries allocated for connections */
);

#if defined(__cplusplus)
}
#endif

Your ShortestNetwork routine will be given a list of numInitialNodes nodes to connect. You may add intermediate nodes to help you form a shorter network, and must produce as output a list of connections between pairs of nodes. The connections must provide a path between any pair of the initial nodes.

Your solution must return the number of intermediate nodes added to the network in *numIntermediateNodes, while storing the location of those nodes in nodes[numInitialNodes+k], k=0..*numIntermediateNodes-1. A connection is specified by storing the indices of the two nodes being connected into the connection array. Your ShortestNetwork routine should return the number of connections created.

The maxNodes and maxConnections parameters indicate how much storage has been allocated for nodes and connections. It is my intention to allocate enough storage for all the nodes and connections your solution might create, but if it turns out that there is not enough storage, your solution should return a value of -1 to indicate that storage was exhausted.

The winner will be the solution that generates the shortest network in the minimum amount of time. Specifically, your solution will be assigned a cost equal to the sum of the distances between nodes in your list of connections, plus a penalty of 10% for each second of execution time. Solutions that do not connect all of the initial nodes will be penalized with a very large cost. The solution with the lowest total cost over a series of networking problems will be the winner.

This will be a native PowerPC Challenge, using the latest CodeWarrior environment. Solutions may be coded in C, C++, or Pascal. Thanks to Michael for suggesting this Challenge.

Three Months Ago Winner

Congratulations to Tom Saxton for submitting the winning solution to the January Sphere Packing Challenge. You may recall that this Challenge was to pack a set of spheres of varying size into a box with minimum volume, and to do so in the shortest amount of time possible. Tom submitted one of only two solutions received for this Challenge, and his was the only one that performed correctly.

Tom's approach is to decide on a footprint for the box to contain the spheres, "drop" the spheres individually into the box until they hit another sphere or the bottom of the box, while attempting to move the dropped sphere around the obstacle without going outside the box footprint. The solution then iterates with random movements to try to converge to a better solution. Tom observed in his submission that the time penalty for this problem (1% per millisecond of execution time) was very severe, making it unproductive to let his algorithm iterate very long. Every tenth of a second of execution time requires a factor of 2 reduction in volume to be productive, a rate of improvement smaller than what Tom was able to achieve.

I evaluated the solutions using six test cases with between 200 and 2000 spheres per test case. As one might expect, execution time grew exponentially with the number of spheres. A test case with 1000 spheres took about 20 times as long to solve as a 200-sphere case, and a 2000-sphere case took about 4 times longer than the 1000-sphere case. Tom's solution generated solutions that, in aggregate, occupied between 1.3 and 3.9 times the volume of individual cubes containing the individual spheres, which suggests that better solutions could be achieved with a more relaxed time penalty.

The table below lists, for each of the solutions submitted, the total volume of the boxes containing the spheres, the total execution time, and the total score including the time penalty, as well as the code and data sizes for each entry. As usual, the number in parentheses after the entrant's name is the total number of Challenge points earned in all Challenges prior to this one.

Name Volume (x1.0E12) Time (secs) Score (x1.0e12) Code Size Data Size
Tom Saxton (79)65.3142.310107.25796372
A. D.***820104

Top Contestants

Listed here are the Top Contestants for the Programmer's Challenge, including everyone who has accumulated 20 or more points during the past two years. The numbers below include points awarded over the 24 most recent contests, including points earned by this month's entrants.

  1. Munter, Ernst 200
  2. Saxton, Tom 99
  3. Boring, Randy 56
  4. Mallett, Jeff 50
  5. Rieken, Willeke 47
  6. Maurer, Sebastian 40
  7. Heithcock, JG 37
  8. Cooper, Greg 34
  9. Murphy, ACC 34
  10. Lewis, Peter 31
  11. Nicolle, Ludovic 27
  12. Brown, Pat 20
  13. Day, Mark 20
  14. Higgins, Charles 20
  15. Hostetter, Mat 20

There are three ways to earn points: (1) scoring in the top 5 of any Challenge, (2) being the first person to find a bug in a published winning solution or, (3) being the first person to suggest a Challenge that I use. The points you can win are:

1st place20 points
2nd place10 points
3rd place7 points
4th place4 points
5th place2 points
finding bug2 points
suggesting Challenge2 points

Here is Tom's winning Sphere Packing solution:

Spheres.cpp
Copyright © 1999 Tom Saxton

#include "Spheres.h"
#include "VecUtil.h"

#include <math.h>
#include <stdlib.h>

enum {
   fFalse = 0,
   fTrue = 1
};

typedef unsigned long ulong;

// disable asserts
#define Assert(f)

// hard iteration limit
#define cIterLim   10000

// scoring an accepting solutions
#define _FAccept(volNew, volBest) ((volNew) < (volBest))
#define _Score(vol, dtick)      ((vol) * (1.0 + (dtick)*10.0/60.0))

// define this to ignore the time penalty
// #define KEEP_GOING

// time checking parameters
#define dtickSec         60
#define dtickCheckScore      (dtickSec/30)
#define dtickFirstCheck      (dtickSec/30)
#define dtickLastCheck      (10*dtickSec)

static const Position s_normalX = { 1.0, 0.0, 0.0 };
static const Position s_normalY = { 0.0, 1.0, 0.0 };
static const Position s_normalZ = { 0.0, 0.0, 1.0 };
static const Position s_normalXNeg = { -1.0, 0.0, 0.0 };
static const Position s_normalYNeg = { 0.0, -1.0, 0.0 };
static const Position s_normalZNeg = { 0.0, 0.0, -1.0 };

static void _InitStartingPos(
   const long csphere,
   long aisphere[],
   const double aradius[],
   double baseMin,
   double baseBest,
   double baseMax,
   double *pbase,
   Position aposStart[]);
static void _TweakStartingPos(
   const long csphere,
   long aisphere[],
   const double aradius[],
   double baseMin,
   double baseBest,
   double baseMax,
   double *pbase,
   Position aposStart[]);
static void _DropSpheres(
   long csphere,
   const long *paisphere,
   const double aradius[],
   const Position *paposStart,
   Position apos[],
   double base,
   double *pvolume);
static void _DropOneSphere(
   const Position &posStart,
   double radius,
   int csphere,
   const long *paisphere,
   const double aradius[],
   const Position apos[],
   Position * pposResult,
   long * pisphereHit);
static int _FFindObstruction(
   const Position normalMove,
   int fNear,
   const Position &posStart,
   double radius,
   int csphere,
   const long *paisphere,
   const double aradius[],
   const Position apos[],
   Position * pposResult,
   long * pisphereHit);

PackSpheres
void PackSpheres(
  long csphere,        /* input: number of spheres to pack */
  double aradius[],    /* input: radius of each of numSpheres spheres */
  Position aposBest[]  /* output: location of center of each sphere */
)
{
   int isphere;
   double volGuess, vol, volBest;
   double base, baseMin, baseMax, baseBest;
   double radiusLarge, radiusSum;
   ulong tickStart, tickCur;
   
   tickStart = LMGetTicks();
   radiusLarge = radiusSum = 0.0;
   for (isphere = 0, volGuess = 0.0; isphere < csphere; ++isphere)
   {
      double radius = aradius[isphere];
      volGuess += 8.0 * radius * radius * radius;
      
      if (radius > radiusLarge)
         radiusLarge = radius;
      radiusSum += radius;
   }
   
   baseMin = 2.0 * radiusLarge;
   baseMax = 2.0 * radiusSum;
   Assert(baseMin <= baseMax);
   
   baseBest = baseMin;
   _DropSpheres(csphere, NULL, aradius, NULL, aposBest, 
            baseBest, &volBest);
   
   base = baseMax;
   _DropSpheres(csphere, NULL, aradius, NULL, aposBest, 
            base, &vol);
   if (vol < volBest)
   {
      volBest = vol;
      baseBest = base;
   }
   
   base = sqrt(baseMin * baseMax);
   _DropSpheres(csphere, NULL, aradius, NULL, aposBest, 
            base, &vol);
   if (vol < volBest)
   {
      volBest = vol;
      baseBest = base;
   }
   
   char * pbBlock = NewPtr(csphere * sizeof(Position) + 
                  csphere * sizeof(Position) + csphere * sizeof(long));
   
   if (pbBlock != NULL)
   {
      long iIter;

      Position * aposStart = (Position *)pbBlock;
      Position * aposEnd = &aposStart[csphere];
      long * aisphere = (long *)&aposEnd[csphere];
      long tickNext = tickStart + dtickCheckScore;
   double scorePrev = _Score(volBest, LMGetTicks() - tickStart);
#ifdef KEEP_GOING
      double scoreBest = scorePrev;
      int iIterBest = 0;
#endif   
      
      for (iIter = 0; iIter < cIterLim; ++iIter)
      {
         tickCur = LMGetTicks();
         if (tickCur >= tickNext)
         {
            ulong dtickCur = tickCur - tickStart;
            if (dtickCur >= dtickFirstCheck)
            {
               if (dtickCur >= dtickLastCheck)
                  break;
                  
               double score = _Score(volBest, dtickCur);
#ifdef KEEP_GOING
               if (score < scoreBest)
               {
                  scoreBest = score;
                  iIterBest = iIter;
               }
#else
               if (score > scorePrev)
                  break;
#endif
               scorePrev = score;
            }
            while (tickNext < tickCur)
               tickNext += dtickCheckScore;
         }
         
         // pick a new scenario
         if (iIter == 0)
            _InitStartingPos(csphere, aisphere, aradius, 
                  baseMin, baseBest, baseMax, &base, aposStart);
         else
            _TweakStartingPos(csphere, aisphere, aradius, 
               baseMin, baseBest, baseMax, &base, aposStart);
         
         // try the new scenario
         _DropSpheres(csphere, aisphere, aradius, aposStart, 
               aposEnd, base, &vol);
         if (_FAccept(vol, volBest))
         {
            volBest = vol;
            baseBest = base;
         BlockMove(aposEnd, aposBest, csphere * sizeof(Position));
         }
         
         // if the largest sphere determined the height, then reduce baseMax
      if (vol <= 2.0 * (radiusLarge + epsilon) * base * base)
         {
            Assert(base <= baseMax);
            baseMax = base;
         }
      }
   }
   
   if (pbBlock != NULL)
      DisposePtr((Ptr) pbBlock);
}

_InitStartingPos
static void _InitStartingPos(
   const long csphere,
   long aisphere[],
   const double aradius[],
   double baseMin,
   double baseBest,
   double baseMax,
   double *pbase,
   Position aposStart[])
{
   long isphereCur;
   
   *pbase = baseBest;
   for (isphereCur = 0; isphereCur < csphere; ++isphereCur)
   {
      Position *ppos = &aposStart[isphereCur];
      double radiusCur = aradius[isphereCur];
      
      aisphere[isphereCur] = isphereCur;
      ppos->coordinate[0] = 
            GRandInRange(radiusCur, *pbase - radiusCur);
      ppos->coordinate[1] = 
            GRandInRange(radiusCur, *pbase - radiusCur);
      ppos->coordinate[2] = csphere * *pbase;
   }
}

_TweakStartingPos
static void _TweakStartingPos(
   const long csphere,
   long aisphere[],
   const double aradius[],
   double baseMin,
   double baseBest,
   double baseMax,
   double *pbase,
   Position aposStart[])
{
   long isphereCur;
   double dbase;
   
   // change the base size?
   if (GRandInRange(0.0, 1.0) < 0.1)
   {
      dbase = GRandInRange(-1.0, 1.0);
      dbase *= fabs(dbase);
      dbase *= 0.25 * (baseMax - baseMin);
      *pbase = baseBest + dbase;
      *pbase = fmax(baseMin, *pbase);
      *pbase = fmin(baseMax, *pbase);
   }
   
   // rearrange the drop order?
   if (GRandInRange(0.0, 1.0) < 4.0)
   {
      for (long index = csphere; - index > 0; )
      {
         long indexSwap;
         long isphereSav;
         
         indexSwap = ((unsigned long)LRand()) % index;
         Assert(0 <= indexSwap && indexSwap < index);
         isphereSav = aisphere[index];
         aisphere[index] = aisphere[indexSwap];
         aisphere[indexSwap] = isphereSav;
      }
   }
   
   // change the starting positions
   for (isphereCur = 0; isphereCur < csphere; ++isphereCur)
   {
      Position *ppos = &aposStart[isphereCur];
      double radiusCur = aradius[isphereCur];
      
      ppos->coordinate[0] = 
            GRandInRange(radiusCur, *pbase - radiusCur);
      ppos->coordinate[1] = 
            GRandInRange(radiusCur, *pbase - radiusCur);
      ppos->coordinate[2] = csphere * *pbase;
   }
}

_DropSpheres
static void _DropSpheres(
   const long csphere,
   const long *paisphere,
   const double aradius[],
   const Position *paposStart,
   Position aposEnd[],
   double base,
   double *pvol)
{
   long csphereDone;
   
   for (csphereDone = 0; csphereDone < csphere; ++csphereDone)
   {
      Position posStart, posLand;
      double radiusCur;
      long isphereHit;
      long isphereCur;

      isphereCur = paisphere == NULL ? csphereDone : 
            paisphere[csphereDone];

      radiusCur = aradius[isphereCur];
      
      // pick a starting point for the current sphere
      Assert(base >= radiusCur*2.0);
      if (paposStart == NULL)
      {
         posStart.coordinate[0] = 
            GRandInRange(radiusCur, base - radiusCur);
         posStart.coordinate[1] = 
            GRandInRange(radiusCur, base - radiusCur);
      }
      else
      {
         posStart.coordinate[0] = 
               paposStart[isphereCur].coordinate[0];
         posStart.coordinate[1] = 
               paposStart[isphereCur].coordinate[1];
      }
      
      // drop it
   _DropOneSphere(posStart, radiusCur, csphereDone, paisphere, 
         aradius, aposEnd, &aposEnd[isphereCur], &isphereHit);
      
      // try to move it around the sphere it hit
   for (int cIter = 0; isphereHit != -1 && cIter < isphereCur; 
            ++cIter)
      {
         Position vecMove, vecMoveXY, normalMove;
         Position posHit;
         double distH, distMove;
         int icoord;
         
         posHit = aposEnd[isphereHit];
         SubVec(aposEnd[isphereCur], posHit, &vecMove);
         vecMoveXY = vecMove;
         vecMoveXY.coordinate[2] = 0;
         distH = LengthVec(vecMoveXY);
         
         if (distH < epsilon)
            break;
            
         ScaleVec(1.0/distH, vecMoveXY, &normalMove);
         distMove = radiusCur + aradius[isphereHit];
         Assert(distMove > distH - epsilon);
         
         // don't move out of the box
         for (icoord = 0; icoord <= 1; ++icoord)
         {
            if (normalMove.coordinate[icoord] < -epsilon)
            {
               if (posHit.coordinate[icoord] + 
         distMove * normalMove.coordinate[icoord] < radiusCur)
            distMove = (radiusCur - posHit.coordinate[icoord]) / 
                     normalMove.coordinate[icoord];
            }
            else if (normalMove.coordinate[icoord] > epsilon)
            {
               if (posHit.coordinate[icoord] + distMove * 
               normalMove.coordinate[icoord] > base - radiusCur)
                  distMove = (base - radiusCur - 
                                             posHit.coordinate[icoord]) / 
                                          normalMove.coordinate[icoord];
            }
         }
         
         Assert(distMove >= distH - epsilon);
         if (distMove < distH + epsilon)
            break;
            
         AddScaleVec(posHit, distMove, normalMove, &posStart);
         
   _DropOneSphere(posStart, radiusCur, csphereDone, paisphere, 
                  aradius, aposEnd, &posLand, &isphereHit);
         
         if (posLand.coordinate[2] > 
                  aposEnd[isphereCur].coordinate[2] - epsilon)
            break;
         
         aposEnd[isphereCur] = posLand;
      
      }

      // try move it toward the edges
      int fImproved, cIter;
      for (fImproved = fTrue, cIter = 1; fImproved; ++cIter)
      {
         Assert(cIter < 15);
         fImproved = fFalse;
         for (int dir = 0; dir < 4; ++dir)
         {
            Position normalMove;
            int fHit;
            double sEdge;
            Position aposStart[2];
            int cposStart;
            
            switch (dir)
            {
            case 0:
               normalMove = s_normalX;
               sEdge = base - radiusCur;
               break;
            case 1:
               normalMove = s_normalY;
               sEdge = base - radiusCur;
               break;
            case 2:
               normalMove = s_normalXNeg;
               sEdge = -radiusCur;
               break;
            case 3:
               normalMove = s_normalYNeg;
               sEdge = -radiusCur;
               break;
            }
            
            fHit = _FFindObstruction(
                     normalMove,
                     fTrue/*fNear*/,
                     aposEnd[isphereCur],
                     radiusCur,
                     csphereDone,
                     paisphere,
                     aradius,
                     aposEnd,
                     &posLand,
                     &isphereHit);
            
            cposStart = 0;
            if (!fHit || DotVec(posLand, normalMove) > sEdge)
            {
               posLand = aposEnd[isphereCur];
               AddScaleVec(posLand, sEdge - 
                        DotVec(posLand, normalMove), normalMove, 
                                    &aposStart[cposStart++]);
               cposStart = 1;
            }
            else
            {
         LinearComboVec(0.5, posLand, 0.5, aposEnd[isphereCur], 
                                 &aposStart[cposStart++]);
               aposStart[cposStart++] = posLand;
            }
            
for (int iposStart = 0; iposStart < cposStart; ++iposStart)
            {
               _DropOneSphere(aposStart[iposStart], radiusCur, 
            csphereDone, paisphere, aradius, aposEnd, &posLand, 
                     &isphereHit);
               
               if (posLand.coordinate[2] < 
                     aposEnd[isphereCur].coordinate[2] + epsilon)
               {
                  if (aposEnd[isphereCur].coordinate[2] - 
                           posLand.coordinate[2] > radiusCur * 0.05)
                     fImproved = fTrue;
                  aposEnd[isphereCur] = posLand;
            
               }
            }
         }
      }
   }
   
   ComputeVol(csphere, NULL, aradius, aposEnd, base, pvol);
}

_DropOneSphere
static void _DropOneSphere(
   const Position &posStart,
   double radius,
   int csphere,
   const long *paisphere,
   const double aradius[],
   const Position apos[],
   Position * pposResult,
   long * pisphereHit)
{
   Position posBase;
   int fHit;
   
   posBase = posStart;
   posBase.coordinate[2] = 0.0;
   
   *pposResult = posBase;

   fHit = _FFindObstruction(
            s_normalZ,
            fFalse, /* fNear */
            posBase,
            radius,
            csphere,
            paisphere,
            aradius,
            apos,
            pposResult,
            pisphereHit
            );
   
   if (!fHit || pposResult->coordinate[2] < radius)
   {
      *pisphereHit = -1;
      pposResult->coordinate[2] = radius;
   }

   // add some fudge
   pposResult->coordinate[2] += epsilon;
   
#ifdef DEBUG
   for (long csphereChecked = 0; csphereChecked < csphere; 
               ++csphereChecked)
   {
      Position vecT;
      double dist, distGap;
      int isphere;
      
      isphere = paisphere == NULL ? csphereChecked : 
                                                paisphere[csphereChecked];
      
      SubVec(apos[isphere], *pposResult, &vecT);
      dist = LengthVec(vecT);
      distGap = dist - (radius + aradius[isphere]);
      Assert(distGap >= 0.0);
   }
#endif
}

_FFindObstruction
// moving a sphere with specifed radius from posStart in the direction normalMove,
// find the nearest or farthest obstruction
// If there is an obstruction, return the index to the obstructing sphere
// and the position to which the object can move.
static int _FFindObstruction(
   const Position normalMove,
   int fNear,
   const Position &posStart,
   double radius,
   int csphere,
   const long *paisphere,
   const double aradius[],
   const Position apos[],
   Position * pposResult,
   long * pisphereHit)
{
   double zBest;
   
   *pisphereHit = -1;

   for (int csphereChecked = 0; csphereChecked < csphere; 
               ++csphereChecked)
   {
      Position vecToOther, vecPerp, vecParallel;
      double distPerpSq, distSep, distSepSq;
      double z, dz;
      int isphere;
      
      isphere = paisphere == NULL ? csphereChecked : 
                                                paisphere[csphereChecked];
      SubVec(apos[isphere], posStart, &vecToOther);
      ProjectVec(vecToOther, normalMove, &vecParallel);
      SubVec(vecToOther, vecParallel, &vecPerp);
      
      distPerpSq = DotVec(vecPerp, vecPerp);
      distSep = radius + aradius[isphere];
      distSepSq = distSep * distSep;
      
      if (distPerpSq > distSepSq)
         continue;
      
      dz = sqrt(distSepSq - distPerpSq);
      if (fNear)
         dz = -dz;
      z = DotVec(vecParallel, normalMove) + dz;
      
      if (z >= 0.0 && (*pisphereHit == -1 || 
                                    (fNear ? z < zBest : z > zBest)))
      {
         zBest = z;
         *pisphereHit = isphere;
      }
   }
   
   if (*pisphereHit == -1)
      return fFalse;
      
   *pposResult = posStart;
   AddScaleVec(posStart, zBest, normalMove, pposResult);
   
   return fTrue;
}

VecUtil.cpp
#include "Spheres.h"
#include "VecUtil.h"

#include <math.h>
#include <stdlib.h>

enum {
   fFalse = 0,
   fTrue = 1
};

// disable asserts
#define Assert(f)

// math utilities

double GRandInRange(double gLow, double gHigh)
{
   double g;
   
   g = gLow + rand() * (gHigh - gLow) / RAND_MAX;
   Assert(gLow <= g && g <= gHigh);
   return g;
}

// return a long's worth of randomness
long LRand()
{
   long lw;
   
   Assert(RAND_MAX > 256);
   
   lw = 0;
   for (int ib = 0; ib < sizeof(long); ++ib)
      lw = (lw << 8) + (rand() & 0xFF);
   return lw;
}

// vector utilities

void SubVec(const Position &pos1, const Position &pos2, 
      Position *pposResult)
{
   for (int i = 0; i < 3; ++i)
      pposResult->coordinate[i] = pos1.coordinate[i] - 
                                                         pos2.coordinate[i];
}

double DotVec(const Position &pos1, const Position &pos2)
{
   double g = 0.0;
   for (int i = 0; i < 3; ++i)
      g += pos1.coordinate[i] * pos2.coordinate[i];
   return g;
}

double LengthVec(const Position &pos)
{
   return sqrt(DotVec(pos, pos));
}

void ScaleVec(double g, const Position &pos, 
   Position *pposResult)
{
   for (int i = 0; i < 3; ++i)
      pposResult->coordinate[i] = g * pos.coordinate[i];
}

void AddScaleVec(const Position &posBase, double g, 
   const Position &posAdd, Position *pposResult)
{
   for (int i = 0; i < 3; ++i)
      pposResult->coordinate[i] = posBase.coordinate[i] + 
                                                g * posAdd.coordinate[i];
}

void LinearComboVec(double g1, const Position &pos1, double g2, 
      const Position &pos2, Position *pposResult)
{
   for (int i = 0; i < 3; ++i)
      pposResult->coordinate[i] = g1 * pos1.coordinate[i] + 
                                                g2 * pos2.coordinate[i];
}

// project "vec" onto a "normal" vector
void ProjectVec(const Position &vec, const Position &normal, 
      Position *pvecResult)
{
   ScaleVec(DotVec(vec, normal), normal, pvecResult);
}

// sphere stuff

void ComputeVol(
   const long csphere,
   const long *paisphere,
   const double aradius[],
   const Position apos[],
   double base,
   double *pvol)
{
   Position posMin, posMax;
   long index;
   int icoord;
   double radius;
   const Position * ppos;

   posMin = posMax = apos[0];

   for (index = 0; index < csphere; ++index)
   {
      long isphere;
      
      isphere = paisphere == NULL ? index : paisphere[index];
      ppos = &apos[isphere];
      radius = aradius[isphere];
      for (icoord = 0; icoord < 3; ++icoord)
      {
         if (ppos->coordinate[icoord] - radius < 
                  posMin.coordinate[icoord])
            posMin.coordinate[icoord] = ppos->coordinate[icoord] - 
                                                            radius;

         if (ppos->coordinate[icoord] + radius > 
                  posMax.coordinate[icoord])
            posMax.coordinate[icoord] = ppos->coordinate[icoord] + 
                                                            radius;
      }
   }

   *pvol = 1.0;
   
   for (icoord = 0; icoord < 3; ++icoord)
   {
      Assert(posMin.coordinate[icoord] >= -epsilon);
   Assert(base == 0 || posMax.coordinate[icoord] <= base+epsilon 
                                       || icoord == 2);
      *pvol *= posMax.coordinate[icoord] - 
                        posMin.coordinate[icoord];
   }
}

Spheres.h

#if defined(__cplusplus)
extern "C" {
#endif

typedef struct Position {
  double coordinate[3];  /* coordinate[0]==X position, [1]==Y, [2]==Z */
} Position;

void PackSpheres(
  long numSpheres,       /* input: number of spheres to pack */
  double radius[],       /* input: radius of each of numSpheres spheres */
  Position location[]    /* output: location of center of each sphere */
);

#if defined (__cplusplus)
}
#endif

VecUtil.h
// error tolerance

const double epsilon (1.0e-10);

// math utilities

double GRandInRange(double gLow, double gHigh);
long LRand();

// vector utilities

void SubVec(const Position &pos1, const Position &pos2, 
   Position *pposResult);
double DotVec(const Position &pos1, const Position &pos2);
double LengthVec(const Position &pos);
void ScaleVec(double g, const Position &pos, 
   Position *pposResult);
void AddScaleVec(const Position &posBase, double g, 
   const Position &posAdd, Position *pposResult);
void ProjectVec(const Position &vec, const Position &normal, 
   Position *pvecResult);
void LinearComboVec(double g1, const Position &pos1, 
   double g2, const Position &pos2, Position *pposResult);

// sphere stuff

void ComputeVol(
   const long csphere,
   const long *paisphere,
   const double aradius[],
   const Position apos[],
   double base,
   double *pvol);
 

Community Search:
MacTech Search:

Software Updates via MacUpdate

WhiteCap 6.7 - Visual plug-in for iTunes...
WhiteCap is a sleek and sophisticated music visualizer and screensaver that features futuristic, wireframe mesh visuals with dynamic backgrounds and colors. WhiteCap contains thousands of visual... Read more
Dropbox 24.4.16 - Cloud backup and synch...
Dropbox is an application that creates a special Finder folder that automatically syncs online and between your computers. It allows you to both backup files and keep them up-to-date between systems... Read more
Amazon Chime 4.2.5645 - Amazon-based com...
Amazon Chime is a communications service that transforms online meetings with a secure, easy-to-use application that you can trust. Amazon Chime works seamlessly across your devices so that you can... Read more
Notion 0.1.8 - A unified workspace for m...
Notion is the unified workspace for modern teams. Features: Integration with Slack Documents Wikis Tasks Note: This application contains in-app and/or external module purchases. Version 0.1.8:... Read more
Google Chrome 58.0.3029.81 - Modern and...
Google Chrome is a Web browser by Google, created to be a modern platform for Web pages and applications. It utilizes very fast loading of Web pages and has a V8 engine, which is a custom built... Read more
Notion 0.1.8 - A unified workspace for m...
Notion is the unified workspace for modern teams. Features: Integration with Slack Documents Wikis Tasks Note: This application contains in-app and/or external module purchases. Version 0.1.8:... Read more
WhiteCap 6.7 - Visual plug-in for iTunes...
WhiteCap is a sleek and sophisticated music visualizer and screensaver that features futuristic, wireframe mesh visuals with dynamic backgrounds and colors. WhiteCap contains thousands of visual... Read more
Dropbox 24.4.16 - Cloud backup and synch...
Dropbox is an application that creates a special Finder folder that automatically syncs online and between your computers. It allows you to both backup files and keep them up-to-date between systems... Read more
Google Chrome 58.0.3029.81 - Modern and...
Google Chrome is a Web browser by Google, created to be a modern platform for Web pages and applications. It utilizes very fast loading of Web pages and has a V8 engine, which is a custom built... Read more
Amazon Chime 4.2.5645 - Amazon-based com...
Amazon Chime is a communications service that transforms online meetings with a secure, easy-to-use application that you can trust. Amazon Chime works seamlessly across your devices so that you can... Read more

Latest Forum Discussions

See All

ChordFlow (Music)
ChordFlow 1.0.0 Device: iOS Universal Category: Music Price: $6.99, Version: 1.0.0 (iTunes) Description: ChordFlow is a chord sequencer with a unique 4-track polyphonic arpeggiator, extensive chord library, MIDI out and Ableton Link... | Read more »
The Walking Dead: A New Frontier is out...
The newest season of Telltale Games'The Walking Dead is well underway. After the release of the third episode, "Above the Law" about a month ago, episode four, "Thicker Than Water" is hot and ready for more zombies and gut-wrenching emotional... | Read more »
Best games we played this week
Another week, another new wave of mobile games do dive into. We've dug through the list of apps that came out this week to tell you which apps are worth your sweet time. And while there weren't too many games this week, there were some big ones.... | Read more »
Vignettes (Games)
Vignettes 1.0.1 Device: iOS Universal Category: Games Price: $2.99, Version: 1.0.1 (iTunes) Description: Vignettes is a casual but unique exploration game without text or characters, where objects shapeshift as you spin them around... | Read more »
Get Me Outta Here is an 80s retro shoote...
Are you ready to fight some aliens? Because Crescent Moon Games has released the retro shooter Get Me Outta Here on iOS devices today. [Read more] | Read more »
Get a bunch of Apple productivity apps f...
If you're an Apple Mac owner, you're probably aware of the host of Apple productivity apps the company includes in all new Mac purchases. Apps like iMovie, Keynote, and of course, GarageBand. While you used to be able to also buy these apps... | Read more »
Terra Mystica (Games)
Terra Mystica 1.03 Device: iOS Universal Category: Games Price: $9.99, Version: 1.03 (iTunes) Description: Short Summary:≈≈≈≈≈≈≈≈≈≈≈≈≈ | Read more »
Ms. Spell (Games)
Ms. Spell 1.0 Device: iOS Universal Category: Games Price: $.99, Version: 1.0 (iTunes) Description: Cast spells and battle monsters in this turn based game, that has you delving into ever the changing Dreadwood to retrieve the lost... | Read more »
Invert - A Minimal Puzzle Game (Games)
Invert - A Minimal Puzzle Game 1.0 Device: iOS Universal Category: Games Price: $2.99, Version: 1.0 (iTunes) Description: Invert is a minimalist puzzle game for fans of brain teasers, logic puzzles, and Rubik's Cube fiddlers. | Read more »
Evergrow: Paper Forest (Games)
Evergrow: Paper Forest 1.0 Device: iOS Universal Category: Games Price: $2.99, Version: 1.0 (iTunes) Description: Solve puzzles in the handcrafted forests of Evergrow through the eyes of an imaginative child and his parents. Discover... | Read more »

Price Scanner via MacPrices.net

15-inch 2.7GHz Touch Bar MacBook Pros on sale...
Amazon has 2016 15″ 2.7GHz Apple Touch Bar MacBook Pros in stock today and on sale for $150-$200 off MSRP. Shipping is free: - 15″ 2.7GHz Touch Bar MacBook Pro Space Gray (sku MLH42LL/A): $2599 $200... Read more
Apple now offering Certified Refurbished 13-i...
Apple is now offering Certified Refurbished 2016 13″ Touch Bar MacBook Pros for $270-$300 off original MSRP. An Apple one-year warranty is included with each model, and shipping is free: - 13″ 2.9GHz... Read more
MyGiHealth Digestive Symptom Tracker Version...
My Total Health, Inc. has announced the release of MyGiHealth 2.1, an important update to their digestive symptom tracker developed exclusively for iPhone, iPad and iPod touch devices. MyGiHealth is... Read more
Galaxy S8 Materials Costs Highest by Far Comp...
The new Samsung Galaxy S8 equipped with 64 gigabytes (GB) of NAND flash memory carries a bill of materials (BOM) cost that comes out to US$301.60, much higher than for previous versions of the... Read more
iCarMode 4.0 Car Dashboard App For iOS Integr...
Indie developer Diego Resnik has announced the release of iCarMode 4.0, an update to his productivity app developed for iOS devices. iCarMode has positioned itself as a true car dashboard app,... Read more
How to save $150+ on Apple’s 13-inch 2.0GHz n...
Apple Authorized Reseller B&H Photo has non-Touch Bar 13″ 2.0GHz MacBook Pros on sale for $150 off MSRP for a limited time. Shipping is free, and B&H charges NY sales tax only: - 13″ 2.0GHz... Read more
15-inch 2.2GHz Retina MacBook Pro, Apple refu...
Apple has Certified Refurbished 2015 15″ 2.2GHz Retina MacBook Pros available for $1699. That’s $300 off MSRP, and it’s the lowest price available for a 15″ MacBook Pro. An Apple one-year warranty is... Read more
Apple Certified Refurbished iMacs available f...
Apple has Certified Refurbished 2015 21″ & 27″ iMacs available for up to $350 off MSRP. Apple’s one-year warranty is standard, and shipping is free. The following models are available: - 21″ 3.... Read more
Save up to $160 with Apple refurbished 9-inch...
Apple has Certified Refurbished 9″ and 12″ Apple iPad Pros available for up to $160 off the cost of new iPads. An Apple one-year warranty is included with each model, and shipping is free: - 32GB 9″... Read more
27-inch Apple iMacs on sale for $200 off MSRP
Amazon has 27″ iMacs on sale for $200-$201 off MSRP, each including free shipping: - 27″ 3.3GHz iMac 5K: $2099.99 $200 off MSRP - 27″ 3.2GHz/1TB Fusion iMac 5K: $1798 $201 off MSRP - 27″ 3.2GHz/1TB... Read more

Jobs Board

*Apple* Media Products - Commerce Engineerin...
Apple Media Products - Commerce Engineering Manager Job Number: 57037480 Santa Clara Valley, California, United States Posted: Apr. 18, 2017 Weekly Hours: 40.00 Job Read more
*Apple* Mac Computer Technician - GeekHampto...
…complex computer issues over the phone and in person? GeekHampton, Long Island's Apple Premium Service Provider, is looking for you! Come work with our crew Read more
*Apple* Mobile Master - Best Buy (United Sta...
**493714BR** **Job Title:** Apple Mobile Master **Location Number:** 001024-Weatherford-Store **Job Description:** **What does a Best Buy Apple Mobile Master Read more
Best Buy *Apple* Computing Master - Best Bu...
**496963BR** **Job Title:** Best Buy Apple Computing Master **Location Number:** 001061-Marina-Store **Job Description:** **What does a Best Buy Apple Computing Read more
*Apple* Mobile Master - Best Buy (United Sta...
What does a Best Buy Apple Mobile Master do? At Best Buy, our mission is to leverage the unique talents and passions of our employees to inspire, delight, and enrich Read more
All contents are Copyright 1984-2011 by Xplain Corporation. All rights reserved. Theme designed by Icreon.