| polygon1.mos |
(!*********************************************************************
Mosel NL examples
=================
file polygon1.mos
`````````````````
Maximize the area of polygon of N vertices and diameter of 1.
The position of vertices is indicated as (rho,theta) coordinates
where rho denotes the distance to the base point (vertex with number N)
and theta the angle from the x-axis.
-- Formulation using direct algebraic expressions --
(c) 2008 Fair Issac Corporation
Creation: 2002, rev. Jun. 2023
*********************************************************************!)
model "Polygon 1"
uses "mmxnlp"
parameters
N=5 ! Number of vertices
SOLVER=0 ! 0: SLP, 1: Knitro
end-parameters
declarations
RN = 1..N
Area: nlctr
rho : array(RN) of mpvar ! Distance of vertex from the base point
theta : array(RN) of mpvar ! Angle from x-axis
D: array(RN,RN) of nlctr ! Limit on side length
end-declarations
! Objective: sum of areas
Area:= (sum (i in 2..N-1) (rho(i)*rho(i-1)*sin(theta(i)-theta(i-1)))) * 0.5
! Bounds and start values for decision variables
forall(i in 1..N-1) do
rho(i) >= 0.1
rho(i) <= 1
setinitval(rho(i), 4*i*(N + 1 - i)/((N+1)^2))
setinitval(theta(i), M_PI*i/N)
end-do
! Third side of all triangles <= 1
forall(i in 1..N-2, j in i+1..N-1)
D(i,j):= rho(i)^2 + rho(j)^2 - rho(i)*rho(j)*2*cos(theta(j)-theta(i)) <= 1
! Vertices in increasing order
forall(i in 2..N-1) theta(i) >= theta(i-1) +.01
! Boundary conditions
theta(N-1) <= M_PI ! Last vertex above x-axis
! Optional parameter settings
setparam("xnlp_verbose", true) ! Enable XNLP output log
setparam("xnlp_solver", SOLVER) ! Select the solver
! In this example we will use a local solver
setparam("xprs_nlpsolver", 1)
! Solve the problem
maximise(Area)
! Solution output
writeln("Area = ", getobjval)
forall (i in 1..N-1)
writeln("V", i, ": r=", getsol(rho(i)), " theta=", getsol(theta(i)))
end-model
|
|
| polygon1_graph.mos |
(!*********************************************************************
Mosel NL examples
=================
file polygon1_graph.mos
```````````````````````
Maximize the area of polygon of N vertices and diameter of 1.
The position of vertices is indicated as (rho,theta) coordinates
where rho denotes the distance to the base point (vertex with number N)
and theta the angle from the x-axis.
-- Formulation using direct algebraic expressions, solution graph --
(c) 2013 Fair Issac Corporation
Creation: Feb. 2013, rev. Jun. 2023
*********************************************************************!)
model "Polygon 1 (graph)"
uses "mmxnlp", "mmsvg"
parameters
N=5 ! Number of vertices
SOLVER=0 ! 0: SLP, 1: Knitro
end-parameters
declarations
RN = 1..N
Area: nlctr
rho : array(RN) of mpvar ! Distance of vertex from the base point
theta : array(RN) of mpvar ! Angle from x-axis
D: array(RN,RN) of nlctr ! Limit on side length
end-declarations
! Objective: sum of areas
Area:= (sum (i in 2..N-1) (rho(i)*rho(i-1)*sin(theta(i)-theta(i-1)))) * 0.5
! Bounds and start values for decision variables
forall(i in 1..N-1) do
rho(i) >= 0.1
rho(i) <= 1
setinitval(rho(i), 4*i*(N + 1 - i)/((N+1)^2))
setinitval(theta(i), M_PI*i/N)
end-do
! Third side of all triangles <= 1
forall(i in 1..N-2, j in i+1..N-1)
D(i,j):= rho(i)^2 + rho(j)^2 - rho(i)*rho(j)*2*cos(theta(j)-theta(i)) <= 1
! Vertices in increasing order
forall(i in 2..N-1) theta(i) >= theta(i-1) +.01
! Boundary conditions
theta(N-1) <= M_PI ! Last vertex above x-axis
! Optional parameter settings
setparam("xnlp_verbose", true) ! Enable XNLP output log
setparam("xnlp_solver", SOLVER) ! Select the solver
! In this example we will use a local solver
setparam("xprs_nlpsolver", 1)
! Solve the problem
maximise(Area)
! Solution output
writeln("Area = ", getobjval)
forall (i in 1..N-1)
writeln("V", i, ": r=", getsol(rho(i)), " theta=", getsol(theta(i)))
!**************** Graphical representation of results ****************
declarations
X,Y: array(0..N) of real
end-declarations
X(N):=0; Y(N):=0 ! Position for base vertex N
X(0):=X(N); Y(0):=Y(N)
forall(i in 1..N-1) do ! Calculate vertex positions
X(i):=cos(theta(i).sol)*rho(i).sol+X(N)
Y(i):=sin(theta(i).sol)*rho(i).sol+Y(N)
end-do
svgaddgroup("P", "Polygon")
forall(i in 1..N) do ! Draw the resulting polygon
svgaddpoint(X(i),Y(i))
svgaddtext(X(i)+0.03,Y(i)+0.02, string(i))
end-do
svgaddpolygon(sum(i in 1..N) [X(i),Y(i)])
! Scale the size of the displayed graph
svgsetgraphscale(200)
svgsetgraphpointsize(2)
svgsave("polygon.svg")
svgrefresh
svgwaitclose("Close browser window to terminate model execution.", 1)
end-model
|
|
| polygon2.mos |
(!*********************************************************************
Mosel NL examples
=================
file polygon2.mos
`````````````````
Maximize the area of polygon of N vertices and diameter of 1.
The position of vertices is indicated as (rho,theta) coordinates
where rho denotes the distance to the base point (vertex with number N)
and theta the angle from the x-axis.
-- Formulation using a simple Mosel user function --
(c) 2008 Fair Issac Corporation
Creation: 2002, rev. Feb. 2013
*********************************************************************!)
model "Polygon 2"
uses "mmxnlp"
parameters
N=5 ! Number of vertices
SOLVER=0 ! 0: SLP, 1: Knitro
end-parameters
declarations
RN = 1..N
Area: nlctr
rho : array(RN) of mpvar ! Distance of vertex from the base point
theta : array(RN) of mpvar ! Angle from x-axis
D: array(RN,RN) of nlctr ! Limit on side length
FunctionArg: array(RN,{"rho","theta"}) of nlctr ! User function arguments
AreaFunction : userfunc ! User function definition
end-declarations
! Objective: sum of areas. Definition of a user function
AreaFunction := userfuncMosel("MoselArea")
! Create function arguments
forall (i in 1..N) do
FunctionArg(i, "rho") := rho(i)
FunctionArg(i, "theta") := theta(i)
end-do
! Use the Mosel user function in a formula for the objective
Area := F(AreaFunction,FunctionArg)
! Bounds and start values for decision variables
forall (i in 1..N-1) do
rho(i) >= 0.1
rho(i) <= 1
setinitval(rho(i),4*i*(N + 1 - i)/((N+1)^2))
setinitval(theta(i),M_PI*i/N)
end-do
! Third side of all triangles <= 1
forall(i in 1..N-2, j in i+1..N-1)
D(i,j) := rho(i)^2 + rho(j)^2 - rho(i)*rho(j)*2*cos(theta(j)-theta(i)) <= 1
! Vertices in increasing order
forall (i in 2..N-1) theta(i) >= theta(i-1) +.01
! Boundary conditions (last vertex above x-axis)
theta(N-1) <= M_PI
! Uncomment to display user function info
! userfuncinfo(AreaFunction)
! Optional parameter settings
setparam("xnlp_verbose", true) ! Enable XNLP output log
setparam("xnlp_solver", SOLVER) ! Select the solver
! Solve the problem
maximise(Area)
! Solution output
writeln("Area = ", getobjval)
forall (i in 1..N-1)
writeln("V",i,": r=",getsol(rho(i))," theta=",getsol(theta(i)))
! **** Definition of the Mosel user function ****
public function MoselArea(I: array(Indices: range, Types: set of string) of real): real
returned := (sum (i in 2..N-1) (I(i,"rho")*I(i-1,"rho")*sin(I(i,"theta")-I(i-1,"theta")))) * 0.5
end-function
end-model
|
|
| polygon3.mos |
(!*********************************************************************
Mosel NL examples
=================
file polygon3.mos
`````````````````
Maximize the area of polygon of N vertices and diameter of 1.
The position of vertices is indicated as (rho,theta) coordinates
where rho denotes the distance to the base point (vertex with number N)
and theta the angle from the x-axis.
-- Formulation using a multi-valued Mosel user function --
(c) 2008 Fair Issac Corporation
Creation: 2002, rev. Feb. 2013
*********************************************************************!)
model "Polygon 3"
uses "mmxnlp"
parameters
N=5 ! Number of vertices
SOLVER=0 ! 0: SLP, 1: Knitro
end-parameters
declarations
RN = 1..N
Area: nlctr
rho : array(RN) of mpvar ! Distance of vertex from the base point
theta : array(RN) of mpvar ! Angle from x-axis
D: array(RN,RN) of nlctr ! Limit on side length
FunctionArg: array(RN,{"rho","theta"}) of nlctr ! User function arguments
AreaFunction : userfunc ! User function definition
k: integer
end-declarations
! Objective: sum of areas. Definition of a user function
AreaFunction := userfuncMosel("MoselArea")
! Create function arguments
forall (i in 1..N) do
FunctionArg(i, "rho") := rho(i)
FunctionArg(i, "theta") := theta(i)
end-do
! Use the Mosel user function in a formula for the objective
Area := F(AreaFunction, FunctionArg, 1)
! Bounds and start values for decision variables
forall (i in 1..N-1) do
rho(i) >= 0.1
rho(i) <= 1
setinitval(rho(i),4*i*(N + 1 - i)/((N+1)^2))
setinitval(theta(i),M_PI*i/N)
end-do
! Third side of all triangles <= 1
k := 2;
forall(i in 1..N-2, j in i+1..N-1, k as counter)
D(i,j) := F(AreaFunction, FunctionArg, k) <= 1
! Vertices in increasing order
forall (i in 2..N-1) theta(i) >= theta(i-1) +.01
! Boundary conditions (last vertex above x-axis)
theta(N-1) <= M_PI
! Uncomment to display user function info
! userfuncinfo(AreaFunction)
! Optional parameter settings
setparam("xnlp_verbose", true) ! Enable XNLP output log
setparam("xnlp_solver", SOLVER) ! Select the solver
! Solve the problem
maximise(Area)
! Solution output
writeln("Area = ", getobjval)
forall (i in 1..N-1)
writeln("V",i,": r=",getsol(rho(i))," theta=",getsol(theta(i)))
! **** Definition of the Mosel user function ****
public procedure MoselArea(I: array(Indices: range, Types: set of string) of real,
CalcValues: array(RetIndices: set of integer) of real)
! Area
CalcValues(1) := (sum (i in 2..N-1) (I(i,"rho")*I(i-1,"rho")*sin(I(i,"theta")-I(i-1,"theta")))) * 0.5
! Distances
k := 2
forall (i in 1..N-2,j in i+1..N-1, k as counter)
CalcValues(k) := I(i,"rho")^2 + I(j,"rho")^2 - I(i,"rho")*I(j,"rho")*2*cos(I(j,"theta")-I(i,"theta"))
end-procedure
end-model
|
|
| polygon8.mos |
(!*********************************************************************
Mosel NL examples
=================
file polygon8.mos
`````````````````
Maximize the area of polygon of N vertices and diameter of 1.
The position of vertices is indicated as (rho,theta) coordinates
where rho denotes the distance to the base point (vertex with number N)
and theta the angle from the x-axis.
-- Formulation using a simple Mosel user function redirected into Java --
!!! Before running this model, compile Polygon.java into Polygon.class
(c) 2018 Fair Isaac Corporation
author: Z.Csizmadia, Sep. 2018
*********************************************************************!)
model "Polygon 8 JavaSimpleFunction"
uses "mmxnlp"
uses 'mosjvm'
parameters
N=5 ! Number of vertices
end-parameters
declarations
Vertices = 1..N
Polar = {"rho","theta"}
Area: nlctr
rho: array(Vertices) of mpvar ! Distance of vertex from the base point
theta: array(Vertices) of mpvar ! Angle from x-axis
D: array(Vertices,Vertices) of nlctr ! Limit on side length
FunctionArg: array(Vertices,Polar) of nlctr ! User function arguments
AreaFunction: userfunc ! User function definition
end-declarations
! Objective - sum of areas. Definition of a user function
AreaFunction := userfuncMosel("AreaInJava")
! Create function arguments
forall(i in Vertices) do
FunctionArg(i, "rho") := rho(i)
FunctionArg(i, "theta") := theta(i)
end-do
! Use the Mosel user function in a formula for the objective
Area := F(AreaFunction,FunctionArg)
! Bounds and start values for decision variables
forall(i in 1..N-1) do
rho(i) >= 0.1
rho(i) <= 1
setinitval(rho(i),4*i*(N + 1 - i)/((N+1)^2))
setinitval(theta(i),M_PI*i/N)
end-do
! Third side of all triangles <= 1
forall (i in 1..N-2, j in i+1..N-1)
D(i,j) := rho(i)^2 + rho(j)^2 - rho(i)*rho(j)*2*cos(theta(j)-theta(i)) <= 1
! Vertices in increasing order
forall(i in 2..N-1) theta(i) >= theta(i-1) +.01
! Boundary conditions (last vertex above x-axis)
theta(N-1) <= M_PI
! Abort model if we encounter a Java exception
setparam('jvmabortonexception', true)
! Tell Java to look for classes in working directory
setparam('jvmclasspath', getparam('workdir'))
! Optional parameter settings
setparam("xnlp_verbose", true) ! Enable XNLP output log
setparam("xnlp_solver", 0) ! Select SLP as the solver
! Uncomment to display user function info
userfuncinfo(AreaFunction)
! Solve the problem
maximise(Area)
! Solution output
writeln("Area = ", getobjval)
forall(i in 1..N-1)
writeln("V", i, ": r=", getsol(rho(i)), " theta=", getsol(theta(i)))
! **** Definition of the Mosel user function ****
public function AreaInJava(I: array(Indices: set of integer, Types: set of string) of real): real
declarations
! Structures for communicating the input variable values to Java
IndicesToJava=0..N-1
rhoToJava: array(IndicesToJava) of real
thetaToJava: array(IndicesToJava) of real
end-declarations
! Project "rho" and "theta" values to structures formatted for Java
forall(i in Vertices) rhoToJava(i-1) := I(i,"rho")
forall(i in Vertices) thetaToJava(i-1) := I(i,"theta")
returned := jvmcallreal( "Polygon.CalculateArea", rhoToJava, thetaToJava )
end-function
end-model
|
|
| polygon8_delta.mos |
(!*********************************************************************
Mosel NL examples
=================
file polygon8_delta.mos
```````````````````````
Maximize the area of polygon of N vertices and diameter of 1.
The position of vertices is indicated as (rho,theta) coordinates
where rho denotes the distance to the base point (vertex with number N)
and theta the angle from the x-axis.
-- Formulation using a Mosel user function returning its own partial
derivatives that is redirected to Java --
!!! Before running this model, compile Polygon.java into Polygon.class
(c) 2018 Fair Isaac Corporation
author: Z.Csizmadia, Sep. 2018
*********************************************************************!)
model "Polygon 8 JavaDelta"
uses "mmxnlp"
uses 'mosjvm'
parameters
N=5 ! Number of vertices
end-parameters
declarations
Vertices = 1..N
Polar = {"rho","theta"}
Area: nlctr
rho: array(Vertices) of mpvar ! Distance of vertex from the base point
theta: array(Vertices) of mpvar ! Angle from x-axis
D: array(Vertices,Vertices) of nlctr ! Limit on side length
FunctionArg: array(Vertices,Polar) of nlctr ! User function arguments
AreaFunction: userfunc ! User function definition
end-declarations
! Objective - sum of areas. Definition of a user function
AreaFunction := userfuncMosel("AreaInJavaDelta",XNLP_DELTAS)
! Create function arguments
forall(i in Vertices) do
FunctionArg(i, "rho") := rho(i)
FunctionArg(i, "theta") := theta(i)
end-do
! Use the Mosel user function in a formula for the objective
Area := F(AreaFunction,FunctionArg)
! Bounds and start values for decision variables
forall(i in 1..N-1) do
rho(i) >= 0.1
rho(i) <= 1
setinitval(rho(i),4*i*(N + 1 - i)/((N+1)^2))
setinitval(theta(i),M_PI*i/N)
end-do
! Third side of all triangles <= 1
forall (i in 1..N-2, j in i+1..N-1)
D(i,j) := rho(i)^2 + rho(j)^2 - rho(i)*rho(j)*2*cos(theta(j)-theta(i)) <= 1
! Vertices in increasing order
forall(i in 2..N-1) theta(i) >= theta(i-1) +.01
! Boundary conditions (last vertex above x-axis)
theta(N-1) <= M_PI
! Abort model if we encounter a Java exception
setparam('jvmabortonexception', true)
! Tell Java to look for classes in working directory
setparam('jvmclasspath',getparam('workdir'))
! Solver parameter settings
setparam("xslp_cascade", 0)
setparam("xnlp_solver", 0) ! Select SLP as the solver
! Optional parameter settings
setparam("xnlp_verbose", true) ! Enable XNLP output log
! Uncomment to display user function info
userfuncinfo(AreaFunction)
! Solve the problem
maximise(Area)
! Solution output
writeln("Area = ", getobjval)
forall(i in 1..N-1)
writeln("V", i, ": r=", getsol(rho(i)), " theta=", getsol(theta(i)))
! **** Definition of the Mosel user function ****
public function AreaInJavaDelta(I: array(Indices: set of integer, Types: set of string) of real, ! Input array
D: array(DIndices: set of integer, DTypes: set of string) of real, ! Return partial derivatives array
Delta: array(YIndices: set of integer, YTypes: set of string) of real): real ! Suggested perturbations
declarations
ctk: integer
! Structures for communicating the input values to Java
IndicesToJava=0..N-1
rhoToJava: array(IndicesToJava) of real
thetaToJava: array(IndicesToJava) of real
! Structures for communicating the input deltas to Java
! Important: make all arrays you pass down to the Java invocation local;
! missing global values would just be substituted by zeroes
DeltasJava: range
deltaIndexVToJava: array(DeltasJava) of integer
deltaIndexPToJava: array(DeltasJava) of integer
deltaToJava: array(DeltasJava) of real
! Return array from Java
javaReturn: jvmobject
end-declarations
! Project "rho" and "theta" values to structures formatted for Java
forall(i in Vertices) rhoToJava(i-1) := I(i,"rho")
forall(i in Vertices) thetaToJava(i-1) := I(i,"theta")
! Project delta values to Java data structures:
! For the partial derivative information, we need to pass to Java the index
! pairs for which the solver is requesting partial derivate information,
! i.e. where Delta(k,l) is nonzero. Since we can only pass flat arrays to Java,
! we store these as 3 separate arrays, the first two storing the index pair of
! the delta we are interested in and the third for the value of the delta.
! The value in "Delta" is the perturbation size suggested by the solver.
ctk := -1 ! Start with index value 0
forall(ctk as counter, k in DIndices, l in DTypes | Delta(k,l) <> 0.0) do
deltaIndexVToJava(ctk) := k-1
deltaIndexPToJava(ctk) := if(l="rho", 1, 2)
deltaToJava(ctk) := Delta(k,l)
end-do
! Call a complex Java function implementing both the function and its
! partial derivatives
javaReturn:= jvmcallobj("Polygon.AreaWithDeltas", rhoToJava, thetaToJava,
deltaIndexVToJava, deltaIndexPToJava, deltaToJava)
! Pick up the objective value, placed on the first index
returned:= jvmcallreal("java.lang.reflect.Array.getDouble", javaReturn, 0)
! Retrieve the delta values stored in the array returned from Java
ctk := 0 ! Start with index value 1 (index zero is the area itself)
forall(ctk as counter, k in DIndices, l in DTypes | Delta(k,l) <> 0.0)
D(k,l):= jvmcallreal("java.lang.reflect.Array.getDouble", javaReturn, ctk)
end-function
end-model
|
|
