(!*********************************************************************
Mosel NL examples
=================
file sphere.mos
```````````````
Find the eqilibrium state distribution of electrons
positioned on a conducting sphere.
Nonconvex NLP problem.
- Alternative objective functions selected via parameter ALG -
Based on AMPL model fekete.mod
Source: http://www.orfe.princeton.edu/~rvdb/ampl/nlmodels/fekete/
(c) 2013 Fair Issac Corporation
author: S. Heipcke, Mar. 2013
*********************************************************************!)
model "Points on sphere"
uses "mmxnlp"
parameters
NP = 50
MDIM = 3
ALG = 2 ! 0: minimize z, 1: maximize Separation, 2: minimize Potential
end-parameters
declarations
POINTS = 1..NP ! Set of points to place
DIMS = 1..MDIM ! Dimensions
x: array(POINTS, DIMS) of mpvar ! Coordinate tuples
z: mpvar ! Objective variable
end-declarations
! Set initial values (= randomly distributed points on the sphere)
setrandseed(3)
if MDIM=2 then
forall(i in POINTS) do
theta(i):= 2*M_PI*random
setinitval(x(i,1), cos(theta(i)) )
setinitval(x(i,2), sin(theta(i)) )
end-do
elif MDIM=3 then
forall(i in POINTS) do
theta(i):= 2*M_PI*random
phi(i):= M_PI*random
setinitval(x(i,1), cos(theta(i))*sin(phi(i)) )
setinitval(x(i,2), sin(theta(i))*sin(phi(i)) )
setinitval(x(i,3), cos(phi(i)) )
end-do
end-if
forall (i in 1..NP, j in DIMS) x(i,j) is_free
! Alternative objective function definitions
forall (i,j in POINTS | iXNLP_STATUS_LOCALLY_OPTIMAL and nlstat<>XNLP_STATUS_OPTIMAL then
writeln("No solution found.")
exit(0)
end-if
! Display the results
writeln("Value of z = ", z.sol)
writeln("Potential = ", Potential.sol)
writeln("Separation = ", Separation.sol)
forall(i in POINTS) do
write(" x(", i, ") = [")
forall(k in DIMS) write(x(i,k).sol, if(k |
(!*********************************************************************
Mosel NL examples
=================
file sphere_graph.mos
`````````````````````
Find the eqilibrium state distribution of electrons
positioned on a conducting sphere.
Nonconvex NLP problem.
- Alternative objective functions selected via parameter ALG -
Based on AMPL model fekete.mod
Source: http://www.orfe.princeton.edu/~rvdb/ampl/nlmodels/fekete/
- Graphical representation of results -
(c) 2013 Fair Issac Corporation
author: S. Heipcke, Mar. 2013, Sep. 2017
*********************************************************************!)
model "Points on sphere"
uses "mmxnlp", "mmsvg"
parameters
NP = 50
MDIM = 3
ALG = 2 ! 0: minimize z, 1: maximize Separation, 2: minimize Potential
end-parameters
declarations
POINTS = 1..NP ! Set of points to place
DIMS = 1..MDIM ! Dimensions
x: array(POINTS, DIMS) of mpvar ! Coordinate tuples
z: mpvar ! Objective variable
end-declarations
! Set initial values (= randomly distributed points on the sphere)
setrandseed(3)
if MDIM=2 then
forall(i in POINTS) do
theta(i):= 2*M_PI*random
setinitval(x(i,1), cos(theta(i)) )
setinitval(x(i,2), sin(theta(i)) )
end-do
elif MDIM=3 then
forall(i in POINTS) do
theta(i):= 2*M_PI*random
phi(i):= M_PI*random
setinitval(x(i,1), cos(theta(i))*sin(phi(i)) )
setinitval(x(i,2), sin(theta(i))*sin(phi(i)) )
setinitval(x(i,3), cos(phi(i)) )
end-do
end-if
forall (i in 1..NP, j in DIMS) x(i,j) is_free
! Alternative objective function definitions
forall (i,j in POINTS | iXNLP_STATUS_LOCALLY_OPTIMAL and nlstat<>XNLP_STATUS_OPTIMAL then
writeln("No solution found.")
exit(0)
end-if
! Display the results
writeln("Value of z = ", z.sol)
writeln("Potential = ", Potential.sol)
writeln("Separation = ", Separation.sol)
forall(i in POINTS) do
write(" x(", i, ") = [")
forall(k in DIMS) write(x(i,k).sol, if(k=3 then ! 3D projection
svgaddgroup("S", "Sphere")
forall(sz in union(i in 0..50) {i*1/50}) do
offset:=0.2*sz
! Overlay circles of different size+colour to obtain 3D effect
svgaddellipse(0,0, 1-sz+offset, 1-sz+offset/2)
col:=svgcolor(50,50+round(100*sz),75+round(125*sz))
svgsetstyle(svggetlastobj, SVG_FILL,col)
svgsetstyle(svggetlastobj, SVG_STROKE,col)
end-do
svgaddgroup("Sol","Solution", SVG_GOLD)
forall(i in POINTS | x(i,3).sol<=0 or
(x(i,1).sol+x(i,1).sol/2>0 and x(i,3).sol<=0.1)) do
offset:=0.2*x(i,3).sol
svgaddpoint(getsol(x(i,1))+offset, getsol(x(i,2))+offset/2 )
end-do
end-if
svgsetgraphscale(100)
svgsetgraphlabels("x1","x2")
! Alternatively change individual sizes:
! svgsetgraphstyle(SVG_STROKEWIDTH,0.01)
! svgsetgraphpointsize(0.02)
svgsave("sphere.svg")
svgrefresh
svgwaitclose("Close browser window to terminate model execution.", 1)
end-model
|
