Nonlinear Programming/ Steiner Weber Modell with barriers

Hello everyone,
I would like to solve a mathematical model with GAMS and unfortunately I can’t find the solution. I would like to solve the following problem with the GAMS programme:

I would like to carry out location planning taking a barrier into account. I will take the federal state of Brandenburg as an example. My demand/requirement locations are all cities in Brandenburg. I get the coordinates of the cities from Excel using geocoordinates and can thus determine the respective distances between the cities. To do this, I form a distance matrix between the cities and insert the values in a table in GAMS. I have one demand point per 200 inhabitants. For example, if I have 1000 inhabitants in a small town, then I have 5 demand points for this town. My barrier should be Berlin. The barrier area is represented by a circle or, even more simply, by the corner points of the rectangle as an area. The corner points of the rectangle represent the area for the barrier zone. The new locations should be shown with X and Y coordinates. If a new location falls on the barrier area, the location is moved to the edge of the barrier. There should also be a capacity limit for the new location. For example, a location can only cover 20 demand points.

My objective function is to minimise the distance between demand points and new locations.There should also be the possibility to limit the number of new locations.

How can I create a formula, that moves locations out of the barrier zone to the edge of the barrier zones?

Many thanks in advance :slight_smile:

This is a straight forward text book optimization problem. Get the math straight and GAMS will follow. If you have detailed questions on how to to particular formulations or tasks in GAMS you will find answers here. Getting the math right is the whole fun and you should not shortcut this experience. Good luck and have fun.

-Michael

Hi Michael,

I know try to program the capacitive multi weber problem with a Cirucums barriers.
The function doesn’t work. I want a forbidden barrier where travelling through it is also forbidden.
Can you please help ?

this is the program:
Sets

j Nachfragestandort /130/
i Neue Standorte /1
15/
;

Parameter

a(j) X-Koordinate Nachfrage Ort also die Breitengrade in Radian /
1 0.914799298
2 0.90339648
3 0.913541692
4 0.914552528
5 0.922110773
6 0.912661755
7 0.901082827
8 0.918043187
9 0.917708665
10 0.920738751
11 0.899135453
12 0.913970752
13 0.923667025
14 0.914169525
15 0.925119257
16 0.900095384
17 0.913679863
18 0.926125031
19 0.916297857
20 0.916297857
21 0.916879634
22 0.916879634
23 0.916588746
24 0.915973032
25 0.917461416
26 0.917461992
27 0.917170522
28 0.915134304
29 0.915716081
30 0.915425193

/

b(j) Y-Koordinate Nachfrage Ort also Längengrade als radiant /
1 0.914799298
2 0.90339648
3 0.913541692
4 0.914552528
5 0.922110773
6 0.912661755
7 0.901082827
8 0.918043187
9 0.917708665
10 0.920738751
11 0.899135453
12 0.913970752
13 0.923667025
14 0.914169525
15 0.925119257
16 0.900095384
17 0.913679863
18 0.926125031
19 0.916297857
20 0.916297857
21 0.916879634
22 0.916879634
23 0.916588746
24 0.915973032
25 0.917461416
26 0.917461992
27 0.917170522
28 0.915134304
29 0.915716081
30 0.915425193


/



d(j) Nachfragewerte pro 100 Einwohner /
1 736.09
2 995.15
3 582.3
4 1857.5
5 1827.6
6 1603.14
7 1026.38
8 1552.26
9 1907.14
10 2112.49
11 1152.12
12 1786.58
13 991.25
14 2043.88
15 775.73
16 1144.29
17 1604.48
18 1208.78
19 3435.92
20 2903.86
21 2942.01
22 2699.67
23 3857.48
24 1672.48
25 571.13
26 806.87
27 2451.97
28 3100.71
29 3509.84
30 2736.89
/

s(i) Kapastandort /1*15 5000/


;


Scalar
sx X-Koordinate vom Kreis Mittelpunkt /52.52/
sy X-Koordinaten vom Kreis Mittelpunkt /13.40495/
R Radius von der Barriere /15/
;


free Variable f;

Positive Variable w(i,j);

Variable
*c(i,j)
x(i) X-Koordinate der Standorte
y(i) Y-Koordinate der Standorte
;



equations
Gesamtdistanz Die Gesamtdistanz zwischen Nachfrager und neuem Standort
Kapa(i) Die maximale Kapazität der neuen Standorte muss kleiner gleich der Nachgefragtenmengen sein
Nachfragebed(j) Die Nachfrage menge entspricht der Transpoortierten Menge
Barriere(i) Die Restriktion für die Bairre
haversine_formula
;

Function haversine_distance /lat1, lon1, lat2, lon2/;
haversine_distance(lat1, lon1, lat2, lon2) =
2 * 6371 * asin(sqrt(sqr(sin((lat2-lat1)/2)) + cos(lat1) * cos(lat2) * sqr(sin((lon2-lon1)/2))));




Gesamtdistanz…
f =e= sum((i,j), haversine_distance(x(i), y(i), a(j), b(j)) * w(i,j));

Nachfragebed(j)…
sum(i,w(i,j))=E=d(j);

Kapa(i)…
sum(j,w(i,j))=L=s(i);

Barriere(i)…
haversine_distance(x(i), y(i), sx, sy) =g= R
;


Model MultiWeber /all/;

Solve MultiWeber using NLP minimizing f;

GAMS has no simple user defined functions (only extrinsic function implemented in a shared library where you need to provide derivatives etc). Anyhow, you can define your function as a macro. Moreover, arcus sinus is called arcsin in GAMS not asin. I have made the changes for you in the attached model. Conopt is the only solver that works okay. The numerics of the model are not very good and hence Conopt terminates with a feasible but not (locally) optimal solution.

-Michael
new8.gms (2.82 KB)

Hello Michael,

First of all, thank you very much for your help. Unfortunately, I had to change the programme after consulting my professor.
My professor now wants to model the Multi Steiner Weber model with me, but it should be a linear problem. We have added some conditions and linearised the program. He also wants me to use the Manhattan distance instead of the Euclidean distance. However, the attached programme has an extremely long calculation time. After more than 35 hours, the programme did not produce any results. So I turned to my professor again, who suggested the following:
‘About your GAMS programme: have you tried debugging it? It seems to me to be just a numerical instability. Try reducing M significantly (e.g. to 100 for the current numbers).
In general, the input data is rather unsuitable for optimisation. I would recommend selecting a desired precision (how many decimal places are sensible/necessary) and then only using whole numbers to represent the coordinates and the demand. M must then be scaled to the smallest possible value (maximum difference that can occur).
In addition, symmetry should be removed from the problem (currently any location can be assigned to any combination of demand locations. It is better to assign location 1 to demand location 1 (z_11 = 1) and then require that location i (>1) can only be assigned to demand location j if its predecessor (j-1) has already been assigned to one of i’s predecessors:

z_{ij} < \sum_{k=1}^{i-1} z_{k,j-1}

Unfortunately, I can’t change the input data as it is the latitude and longitude of the earth, so the number of decimal places must be correct.
Unfortunately, I can’t solve the problem with the symmetry.
Can you perhaps insert my professor’s line into my Gams programme?
I don’t understand how I can say in the sum that k=1 to i-1 and it sums over z(k,j-1). I always get error messages in Gams.
Thank you very much in advance
Multi Weber Final2.gms (2.87 KB)

You need to comment the definition of k as a new set and just make it an alias. Otherwise, the solution is relatively simple (if you know your basic GAMS):

alias(i,k);
sym(i,j)$(ord(i)>1 and ord(j)>1)..
z(i,j)=L=sum(k$(ord(k)<ord(i)),z(k,j-1));

I have attached the entire model. The symmetry breaking constraint help, but that won’t be enough to solve to optimality. Running MIPs for many hours only helps if you see progress in the run. You need to watch the increase of the lower bound and the decrease of the primal bound. Only if the solver makes progress then you can hope for a reasonable finite running time. Don’t forget MIP problems are in general hard and there is no guarantee that you get an optimal solution for a problem of this size in reasonable time.

-Michael
Multi Weber Final2.gms (2.91 KB)