Go to online
Centers-of-Gravity Calculator

How to calculate Centers‑of‑Gravity

It is common practice to run a Center‑of‑Gravity analysis in supply chain network design.
This page explains the calculations behind Centers‑of‑Gravity Calculator, a web app that determines Centers‑of‑Gravity - indicative warehouse locations that minimize transport costs.
It explains how to determine a single Center‑of‑Gravity (two alternatives). This involves more than applying a simple formula, as there exists no formula that brings the optimal position.
It explains how to determine multiple Centers‑of‑Gravity. This takes various iterative calculations that are embedded in two separate algorithms.


 

Centers‑of‑Gravity are indicative warehouse locations that minimize transport costs

In fact, Centers‑of‑Gravity are those locations that minimize the sum of weighted distances. Weighted distance is the distance from warehouse to customer multiplied by customer demand. If a customer has a demand of 10 and is located 25 kilometers from its warehouse, then its weighted distance is 250 kilometers. The sum of weighted distances of all customers acts as an indicator for transport costs.

Towards the algorithms - Introductory thoughts

Upfront: there exists no mathematical formula that can (simply) calculate the optimal warehouse position given customer locations and demand. So...

Imagine a case with two customers only. Customer A has a demand of 10 and customer B of 1. Where is the Center‑of‑Gravity? Somewhere on line A-B, closer to A?


Well, you could say that customer A pulls 10 times harder than customer B. If the 'Center‑of‑Gravity' * moves a distance d towards A the sum of weighted distances decreases (10 × d) – (1 × d) = 9 × d.
* If put between quotes then it refers to a location under investigation which is not yet the true Center‑of‑Gravity.


So, the Center‑of‑Gravity is right on top of customer A, not somewhere in between A and B!



The resultant pull force - sum of all customer pull forces - is leading. It always points in the direction that the 'Center‑of‑Gravity' should move into to decrease the sum of weighted distances. This resultant pull force is constructed by linking all customer pull forces together, head to tail.

However, the sum of weighted distances will only decrease if the 'Center‑of‑Gravity' moves the right distance! The resultant pull force does not tell how far to move. Move this far?


Or this far?



The smaller the move distance, the less chance of 'overshooting' the Center‑of‑Gravity, but the many more moves to be made. Therefore, it is better start with big moves. In case of 'extreme overshooting' the sum of weighted distances will have increased instead of decreased, and the resultant force arrow changed - or even completely reversed - its direction. By then, the move size needs to be reduced in order to get closer to the Center‑of‑Gravity. Else, the next move would end up at the previous position, followed by getting stuck in a loop of going back-and-forth without ever getting any closer. By reducing the move size each time the sum of weighted distances increased the 'Center‑of‑Gravity' ends up on top of the Center‑of‑Gravity.



Of course, the above was an extreme example. Usually, there is no single dominant customer, and the Center‑of‑Gravity will be somewhere in the center of all customers. Also then, the Center‑of‑Gravity is found by moving in the direction of the resultant pull force, but not too far. The resultant pull force becomes smaller the closer the location gets to the Center‑of‑Gravity (which implicates that the length of the resultant pull force somehow does relate to the distance to be moved). It becomes zero* on top of the Center‑of‑Gravity, with all underlying pull forces (F) in a state of equilibrium in both x and y direction: ΣFx=0, ΣFy=0. All forces linked together start and finish at exactly the same location = zero resultant force.



* Except if one dominant customer outweighs all other customers. Then the resultant force never becomes zero. Another exceptional case, rather theoretical, is one in which there are only two customers, A and B, with equal demand. Then the resultant force is zero for all points between A and B on line A-B. Each of these points is Center-of-Gravity.

Single-Center‑of‑Gravity algorithm

The Single-Center‑of‑Gravity algorithm works like this on an x,y plane:
  1. Initialization: locate the 'Center‑of‑Gravity' at an initial (random) position, and give the move distance a large value.
  2. Calculate the resultant pull force. Move the 'Center‑of‑Gravity' in that direction, at move distance.
  3. Calculate the sum of weighted distances. If it increased, then divide the move distance by two.
  4. Repeat steps 2 and 3 until the move distance becomes very small.

More detailed

  1. Initialization: locate the 'Center‑of‑Gravity' at an initial (random) position, and give the move distance a large value.
    • A good initial position is the weighted average x and y coordinates of its customers. It is often near-optimal, but not the Center-of-Gravity, though you may have read elsewhere it is.
    • Often, the sum of weighted distances of the weighted average x,y position versus that of the true Center‑of‑Gravity will be 0-5% worse, but more than 5% worse is also possible, as sometimes happens in the simulation further below.

      Theoretically even 100% worse is possible! Imagine a case with only two customers, customer A with demand 1000 at x,y position (0, 0) and B with demand 1 at position (100, 0). The weighted average x,y position is (0.0999, 0), with 0.0999 calculated as (1000×0 + 1×100)/(1000+1). If the 'Center‑of‑Gravity' moves a distance d towards customer A the goal value improves 1000×d (closer to A) − 1×d (further from B) = 999×d. So the Center‑of‑Gravity is on top of customer A, not at the initial weighted average x,y position! The Center‑of‑Gravity has a sum of weighted distances of 1000×0 (customer A) + 1×100 (customer B) = 100, whereas the initial position has a value of 1000×0.0999 (customer A) + 1×99.9 (customer B) = 199.8, which is 99.8% worse!

    • A random location within the customer cloud can be used just as well, as the remainder of the algorithm will always find the Center‑of‑Gravity, in almost equal time.

    • It is better to start with a bit too large move distance (100 kilometers) than with a bit too small, because if too small then it will take many iterations before the Center‑of‑Gravity is reached.

  2. Calculate the resultant pull force. Move the 'Center‑of‑Gravity' in that direction, at move distance.
    • Sum up all pull forces (vectors) to get the resultant pull force (vector) with a certain size (less relevant) and direction (most relevant) by summing up x and y components separately. For each customer pull force F: Fx = length demand vector × (customer_x − CoG_x)/(distance from customer to CoG), or shorter: Fx = demand × cosinus(angle F). Fy = demand × sinus(angle F).
    • Normalize this resultant pull force (a vector with origin at 0,0 and end point at x,y) by dividing x by its length √(x2 + y2) and y by its length. The vector keeps pointing in the same direction, but its length becomes 1, which makes it a unit vector.
    • Multiply x and y values of this unit vector with the step size to get the move vector representing move distance in x and y direction.
    • If the move vector's tail starts at the current 'Center‑of‑Gravity' position, then its head ends at the next 'Center‑of‑Gravity' position. So, add the move distance in x direction to the current x position of the 'Center‑of‑Gravity', and add distance in y direction to its current y position. This brings the next 'Center‑of‑Gravity' position.

  3. Calculate the sum of weighted distances. If it increased, then divide the move distance by two.
    • This approach originates from the binary search algorithm - a.k.a. half interval search - to find the position of a target value within a sorted array. See half-interval search (wikipedia).
    • Optional: move the 'Center‑of‑Gravity' back to its previous position if the sum increased, so with each 'accepted move' the solution always improves (or remains equal, but never worsens).

  4. Repeat steps 2 and 3 until the move distance becomes very small.
    • Then the Center‑of‑Gravity hardly changes anymore, and the sum of weighted distances converges.
    • The Center‑of‑Gravity may be in the middle of a lake or on top of a mountain.
This algorithm will always find the global optimum with a single run. With this type of problem, there is no chance of getting stuck in local optimum. The algorithm relates to the gradient descent method (wikipedia), also used by Excel Solver's method GRG Nonlineair (GRG stands for Generalized Reduced Gradient).

Determining the Center‑of‑Gravity is like finding the x and y coordinate of the lowest point (x,y,z) in a bowl-shaped 3D‑chart - visualizing the sum of weighted distances z = f(x,y) for each x and y - by going into the direction of steepest descent as pointed to by the resultant pull force. During this descent the direction gets adjusted after each move made/step taken.

Wikipedia presents complex methods to determine a proper move distance (a.k.a. step size) per iteration. Luckily, there is no urgent need for those, because the method of step 3 works quite good!

Visual simulation

Press button Next process step to start the simulation.

Earth is a globe - Required adjustments

X,y on a flat plane (Chartesian coordinate system) becomes latitude,longitude on a globe (Spherical coordinate system). The algorithm involves summing up vectors with a length and an angle. The distance between two points on a globe is calculated with the Haversine formula, and the angle with a Bearing formula. Demand vectors are mapped to a flat tangent space before summing up, and the result is translated back to the globe - see also the picture of section 'Weiszfeld's algorithm as alternative'.
"No, it's not
a flat plane."
 
Latitudes and longitudes in both formulas need to be expressed in radians = degrees × 𝜋/180.

Haversine formula - in 2 steps
  1. a = sinus2((latitude_customer − latitude_cog)/2) + cos(latitude_cog) × cos(latitude_customer) × sinus2((longitude_customer − longitude_cog)/2)
  2. As-the-crow-flies-distance = 2 × atan2(√a , √(1-a)) × 6371. 6371 kilometers is Earth's mean radius.
See Excel - As-the-crow-flies distance function for an Excel implementation. Excel's atan2 function takes its arguments in reversed order!

Bearing formula
Bearing = atan2(sinus(longitude_customer − longitude_cog) × cosinus(latitude_customer) , cosinus(latitude_cog) × sinus(latitude_customer) − sinus(latitude_cog) × cosinus(latitude_customer) × cosinus(longitude_customer − longitude_cog)

Bearing is measured against the vertical longitude axis, not the horizontal latitiude axis, so you need to swap the use of cosinus/sinus when calculating forces in x and y direction.

Weiszfeld's algorithm as alternative

Weiszfeld's algorithm - see also Geometric median (wikipedia) - works like this on an x,y plane:
  1. Initialization: locate the 'Center‑of‑Gravity' at an initial (random) position.
  2. Calculate each customer weight as demand divided by distance to 'Center‑of‑Gravity'.
  3. Calculate the weighted average x,y of customers as the next 'Center‑of‑Gravity'.
  4. Repeat steps 2 and 3 until the sum of weighted distances converges.
Applying Weiszfeld's algorithm to Earth is also done via a flat tangent space.

Source: section 2.3 Weiszfeld Algorithm on a Riemannian Manifold of article Generalized Weiszfeld Algorithms for Lq Optimization by K. Aftab, R. Hartley, J. Trumpf.

Differences between both algorithms

Major difference between Single-Center‑of‑Gravity algorithm and Weiszfeld's algorithm is that the former uses an explicitly controlled move distance whereas - you could say - the latter has an implicit variable move distance of 1 ∕ Σ(demandi ∕ distancei), that is applied to a non-normalized resultant pull force vector. The controlled move distance enables a feature of Centers‑of‑Gravity Calculator: allowing predefined warehouses to move within a limited range. They may need to end up somewhere at the border of this range and Weiszfeld's algorithm as such cannot handle this. But this controlled move distance also means that, whereas under Weiszfeld's algorithm the sum of weighted distances decreases with each iteration, under the Single-Center‑of‑Gravity algorithm it decreases only if the move distance did not result in 'extreme overshooting'.

Weiszfeld's algorithm converges approx. 10% faster. After 15 iterations the sum of weighted distances will have become less than 0.01% above absolute minimum. Weiszfeld's algorithm does not use an explicit move distance, so no need to worry about its initialization value. Under the Single‑Center‑of‑Gravity algorithm: if the initial move distance is set too high, then the first few iterations will become useless 'extreme overshoots' (as can be seen in the chart), but if set too low it will take many more iterations to move to the Center‑of‑Gravity, so it's better to initialize it a bit too high.

Centers‑of‑Gravity Calculator uses both algorithms: the faster Weiszfeld's algorithm as its default, and Single‑Center‑of‑Gravity algorithm in case of predefined warehouses with limited move range.

Multiple‑Centers‑of‑Gravity algorithm

Below, 'Center-of-Gravity' and warehouse are used interchangeably.

A single run of the Multiple‑Centers‑of‑Gravity algorithm does the following:
  1. Initialization: locate each warehouse at a random location within the customer cloud.
  2. (Re)assign each customer to its closest warehouse.
  3. Move each warehouse to the Center‑of‑Gravity of its assigned customers.
    = apply Single‑Center‑of‑Gravity algorithm / Weiszfeld's algorithm per warehouse
  4. Repeat steps 2 and 3 until the sum of weighted distances converges.
Multiple runs are done to find the global optimum, as a run may get stuck in a local optimum. The more runs, the more likely the global optimum is found. Usually, this global optimum will then have been found multiple times as best solution. The warehouse locations of the best solution found are considered to be the Centers‑of‑Gravity. But - especially with a higher number of warehouses - it may happen that two different warehouse structures bring an almost similar sum of weighted distances.

The algorithm resembles the k-means clustering algorithm (wikipedia) used for cluster analysis in data mining and machine learning, as depicted on the left.
K-means Clustering creates clusters that are purely distance based, just like each customer gets linked to its closest warehouse.
Expectation–Maximization (EM) algorithm creates clusters differently. With Centers‑of‑Gravity Calculator you can group customers beforehand - like all dots in the yellow oval - that need to be served by a single warehouse.

Animations

The animations show what happens during a single run of this algorithm. Customers are assigned to their closest warehouse, warehouses move to the center‑of‑gravity of their assigned customers, customer are reassigned, warehouses move, et cetera, until the final situation is reached where none of the customer assignments change anymore, and none of the warehouse positions.

Although each run starts with different random locations, the solution found (global optimum with a minimal sum of weighted distances) is the same in the first three run examples.

 Run 1 - global minimum  Run 2 - global minimum  Run 3 - global minimum

The example below shows a run that got stuck in a local optimum. Therefore, multiple runs need to be done to find the global optimum.

 Run 4 - local minimum
 

A completely different approach

A completely different approach is one in which a set of hundreds or thousands of possible warehouse locations are pre-generated, after which the k best ones are selected (see Combination - wikipedia). Advantage is that the warehouse-customer cost matrix can be based on various transport rate structures and can capture warehousing costs. Disadvantages are that it requires many more inputs, and that the number of possible warehouse combinations - thus runtime - may grow exponentially.
Our Warehouse Selection Optimizer web app uses this approach.

Underlying assumptions

A major underlying assumption is that transport cost = rate/kilometer × distance. This assumption is only partly valid. Parcel rates are often distance independent within a region, FTL pallet rate/kilometer is lower than LTL pallet rate/kilometer, macro-economic imbalances cause direction-dependent rates.

Besides, as-the-crow-flies distances × circuitry factor are used as approximation of road distances. What is meant with road distance? The distance of the fastest or the shortest route? The fastest route is usually more relevant, as driver and truck cost are more time-related than distance-related. But the fastest route may depend on the time of the day, and on the day of the week. Reality is complex. As the estimation errors will be more or less the same for each direction, approximations are acceptable, and commonly used in supply chain network design software.

Broader supply chain network design perspective

Of course, transport costs are only a part of supply chain costs. The optimal number of warehouses and their locations are driven by many quantitative and qualitative factors such as (future) transport and warehousing rates, future demand (and supply), lead time requirements, inventory effects, supply chain risk/redundancy, contractual obligations, distance to parcel/logistics hubs, workforce availability/public transport connections, et cetera.

Nevertheless, it is common practice to run a Center‑of‑Gravity-analysis to get a view on which logical warehouse locations to consider in a supply chain network design.

Centers‑of‑Gravity Calculator

Centers‑of‑Gravity Calculator uses a more advanced version of weighted distances, defined as distance × demand (1st weight) × transport rate (2nd weight). The transport rate enables differentiation between a demand volume transported in small shipments versus the same volume transported in large shipments, which is less expensive.

This web app provides many more advanced functionalities. Global logistic service providers, multinationals, and supply chain consultancy firms are using it all over the globe, and have been using the software for years.