Breeze and GRB data types link by chi sq distribution
AnsweredDear all,
I am trying to translate the following problem (written in Python)
def chi_sq(x, mod):
UVM = np.zeros(np.shape(mod),dtype=np.complex128)
for m in range(len(mod)):
UVM[m] = np.sum( x*np.exp(2.*np.pi*1.j )
return np.sum(np.abs(UVM-mod)**2)
So, a kind of a chi^2 distribution (np.sum(np.abs(UVM-mod)**2) in a Fourier Domain, subject to a non linear constraint which I know how to deal with.
Having thousands of variables, my naive problem is to write this problem
min_x chi^2(x,some_data); s.t. f(x)=0
in a way that gurobi, using Java, can read it.
My code is
val env: GRBEnv = new GRBEnv("test.log")
val model: GRBModel = new GRBModel(env)
val x: List[GRBVar] = List.tabulate(nPix)(i=>model.addVar(1e-3/nPix, GRB.INFINITY, 0.0, GRB.CONTINUOUS, "x"))
val y: List[GRBVar] = List.tabulate(nPix)(i=>model.addVar(-GRB.INFINITY, GRB.INFINITY, 1.0, GRB.CONTINUOUS, "y"))
val obj: GRBLinExpr = new GRBLinExpr()
obj.addTerms(Array.fill(nPix)(1.0), y.toArray)
model.setObjective(obj, GRB.MINIMIZE)
def createVar(stPoint: Double): LazyList[Double] = LazyList.cons(stPoint, createVar(stPoint + 0.01))
val xs: LazyList[Double] = createVar(0.01).take(1000)
val ys: List[Double] = xs.view.map(el => f(el))
model.addConstr(eq, GRB.EQUAL, 1.75, "eq")
(x zip y).map(pair => model.addGenConstrPWL(pair._1, pair._2, xs.toArray, ys.toArray, "pwl"))
model.optimize()
Where the chisq I've just defined to the breeze Dense vectors. What I am solving here is
min y
s.t. y(x) := f(x) + log(chiSq(x) - 0.05)
x,y being vectors.
Since chisq is defined as DenseVectors, I don't how to link it with GRB data types. I will love to have it in the Objective functions.
My main concern is how to write (data-x)^2, where data is a DenseVector and x an array of GRBVal. In other words, in how to rewrite chisq as a function of DenseVectors and GRBVal.
I would be very happy if you could help me,
Best regards,
Alejandro
-
Official comment
-
Gurobi Staff
This post is more than three years old. Some information may not be up to date. For current information, please check the Gurobi Documentation or Knowledge Base. If you need more help, please create a new post in the community forum. Or why not try our AI Gurobot?. -
-
-
Hi Alejandro,
Could you elaborate more on what exactly you are trying to achieve? You say you have a \(\texttt{data}\) vector of double values and an array of GRBVar objects (you have written GRBVal, but I assume you actually meant GRBVar, right?)
Are you trying to model something like
\[ \sum_i (d_i - x_i)^2 \]
or rather
\[ \left ( \sum_i (d_i - x_i) \right )^2 \]
or something else (\(d_i\) are the data entries)?
Best regards,
Jaromił0 -
-
Hi Jaromil,
Thank you very much for your answer.
What I am trying to get is precissely
\[ \sum_i (d_i - x_i)^2 \]
where x_i are (indeed, sorry for the confusion) GRBVar objects, with coefficients in \mathbb{R}, and d_i are complex values. Thus, actually is not only (d_i - x_i)^2 but norm(d_i - x_i)^2.
The python code --using numpy and not Gurobi objects-- of the function is:
def chi_sq(x, mod):
UVM = np.zeros(np.shape(mod),dtype=np.complex128)
for m in range(len(mod)):
UVM[m] = np.sum( x*np.exp(2.*np.pi*1.j )
return np.sum(np.abs(UVM-mod)**2)
(NOTE: mod variable is our d) and my optimization problem is
min_x x*log(x) + log(\sum_i (d_i - x_i)^2 - constant).
Since we need to rewrite the log in a form of GenConstraints, I am solving
min y
s.t. y(x) := x*log(x) + log(\sum_i (d_i - x_i)^2 - constant).
So my problems are:
- How would be written \sum_i (d_i - x_i)^2 , in terms of GRBVar objects. Indeed, how would be the function chi_sq implemented if the inputs are x (list of GRBVar) and d (list of complex values).
- How could be written x*log(x) + log(\sum_i (d_i - x_i)^2 - constant) in terms of GenConstraint.
For the log I know, but since the chi_sq that I have done in Python has as input a vector, how would it fit? Or can chi_sq be included directly on the cost function and added the log in a GenConstraint form?
Thank you,
Best,
Alejandro Mus
0 -
-
-
Hi Alejandro,
How would be written \sum_i (d_i - x_i)^2 , in terms of GRBVar objects. Indeed, how would be the function chi_sq implemented if the inputs are x (list of GRBVar) and d (list of complex values).
Gurobi does not offer the functionality of complex numbers. So as long as your data values have an imaginary part, i.e., \(d_i= a+ i\cdot b\), there is currently no way to implement it with the \(\texttt{gurobipy}\) module. However, if you can find a suitable way to use real number values for your data, you can write \(\sum_i (d_i - x_i)^2\) as
import gurobipy as gp
from gurobipy import GRB
m = gp.Model("test")
d = [1,2,3]
c = [4,5,6]
cst = 42
x = m.addVars(3,name="x")
aux1 = m.addVar(lb=-GRB.INFINITY,vtype=GRB.CONTINUOUS,name="auxVar1")
qexpr = gp.QuadExpr(0)
for i in range(3):
qexpr.add(d[i]*d[i] - 2*d[i]*c[i]*x[i] + c[i]*c[i]*x[i]*x[i])
qexpr.add(-cst)
m.addConstr(qexpr - aux1 ==0, name="qConstr")
m.write("myLP.lp")The above code constructs the constraint
\[\sum_i (d_i^2 - 2d_ic_ix_i + c_i^2x_i^2) - \text{cst} = \text{auxVar1}\]
How could be written x*log(x) + log(\sum_i (d_i - x_i)^2 - constant) in terms of GenConstraint.
You can then use addGenConstrLog to construct your objective
y = m.addVar(name="y")
log1 = m.addVar(lb=-GRB.INFINITY,vtype=GRB.CONTINUOUS,name="log1_var")
log2 = m.addVar(lb=-GRB.INFINITY,vtype=GRB.CONTINUOUS,name="log2_var")
m.addGenConstrLog(y,log1,name="log(y)")
m.addGenConstrLog(aux1,log2,name="log(sum)")
m.setObjective(y*log1 + log2, GRB.MINIMIZE)
m.write("myLP.lp")Note that I used the write function to analyze the LP file and check whether the final problem looks what I expect it to be.
Best regards,
Jaromił0 -
-
Thanks a lot for your answer!!
I leave the code with your toy example translated to Scala and working with DenseVectors, in case someone have the same problem in other language (because, for example y*log1 + log2 I was not able to write it directly in Java due to typation limitations) and can be inspired with it.
What I've done to overcome y*log1 + log2 is to solve
Minimize
y + log1_var
Subject To
qConstr >= cst
General Constraints
log(chi_sq): log1_var = LOG ( auxiliaryVariable1 )
y(x) := x*log(x) = PWL(x)In Scala code:
val aux1: GRBVar = model.addVar(1.0, GRB.INFINITY, 1.0, GRB.CONTINUOUS, "auxiliaryVariable1")
val quadExpr: GRBQuadExpr = new GRBQuadExpr()
V.mapValues(Complex.complexNorm(_)).mapValues(el=>quadExpr.addConstant(el*el))
//SAME: (V*:*V).mapValues(Complex.complexNorm(_)).mapValues(quadExpr.addConstant)
quadExpr.addTerms(V.mapValues(_.real).toArray,x.toArray)
quadExpr.addTerms(Array.fill(n)(1.0),x.toArray,x.toArray)
quadExpr.addConstant(cst)
model.addQConstr(quadExpr,GRB.GREATER_EQUAL,aux1,"qConstr")
val log1: GRBVar = model.addVar(1e-3/n, GRB.INFINITY, 1.0, GRB.CONTINUOUS, "log1_var")
model.addGenConstrLog(aux1,log1,"log(chi_sq)","")
def createVar(stPoint: Double): LazyList[Double] = LazyList.cons(stPoint, createVar(stPoint + 0.01))
val xs: LazyList[Double] = createVar(0.01).take(100)
val ys: Array[Double] = xs.view.map(el => el*log(el)).toArray
val y: GRBVar = model.addVar(-GRB.INFINITY, GRB.INFINITY, 0.0, GRB.CONTINUOUS, "y")
x.map(el=>model.addGenConstrPWL(el,y,xs.toArray,ys,"pwl"))
val obj: GRBLinExpr = new GRBLinExpr()
obj.addTerm(1.0, y)
obj.addTerm(1.0,log1)
model.setObjective(obj, GRB.MINIMIZE)
model.write("myLP.lp")
model.optimize()Indeed, it is a pity that can't deal with complex values, but I will figure out on how to rewrite the problem.
Thanks a lot again for your time and your help,
Bests,
Alejandro Mus
0 -
Post is closed for comments.
Comments
5 comments