Skip to content

Commit

Permalink
Merge branch 'master' of https://github.com/Allinsights/Gravity
Browse files Browse the repository at this point in the history
  • Loading branch information
hhijazi committed Jan 10, 2018
2 parents b8abac6 + d66f557 commit 943dc40
Show file tree
Hide file tree
Showing 4 changed files with 91 additions and 37 deletions.
40 changes: 25 additions & 15 deletions examples/MINLP/Power/SDPOPF/SDPOPF_main.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -322,14 +322,15 @@ int main (int argc, char * argv[]) {
}
what = b->nfp();
node_pairs b_pairs;
// param<> R_Wij_star("R_Wij_star"), R_Wij_hat("I_Wij_hat"),
param<double> R_star("R_star"), I_star("I_star"), W_star("W_star");
param<double> R_diff("R_Diff"), I_diff("I_diff"), W_diff("W_diff");
param<double> R_hat("R_hat"), I_hat("I_hat"), W_hat("W_hat");
// param<> Wii_star("Wii_star"), Wii_hat("Wii_hat"), W_diff("W_Diff");
// the cuts for different dimensions don't have the same form...
Constraint sdpcut("sdpcut_" + to_string(numcuts));
Node *ni;
Arc *aij;
double sdp_cst = 0;
for (int i = 0; i < b->_nodes.size(); i++) {
for (int j = i; j < b->_nodes.size(); j++) {
if (i == j) {
Expand All @@ -339,6 +340,10 @@ int main (int argc, char * argv[]) {

W_diff.set_val(ni->_name,((Bus *) ni)->w - what(namew).eval());
W_hat.set_val(ni->_name,what(namew).eval());
W_star.set_val(ni->_name,((Bus *) ni)->w);

// sdp_cst += what(namew).eval()*what(namew).eval() - ((Bus *) ni)->w*what(namew).eval();
// sdp_cst += ((Bus *) ni)->w - what(namew).eval();
} else {
aij = grid->get_arc(b->_nodes[i]->_name, b->_nodes[j]->_name);
if(aij->_imaginary && !aij->_active) {
Expand All @@ -358,35 +363,40 @@ int main (int argc, char * argv[]) {
I_diff.set_val(aij->_src->_name + "," + aij->_dest->_name,((Line *) aij)->wi - what(namewi).eval());
R_hat.set_val(aij->_src->_name + "," + aij->_dest->_name,what(namewr).eval());
I_hat.set_val(aij->_src->_name + "," + aij->_dest->_name,what(namewi).eval());
R_star.set_val(aij->_src->_name + "," + aij->_dest->_name,((Line *) aij)->wr);
I_star.set_val(aij->_src->_name + "," + aij->_dest->_name,((Line *) aij)->wi);

sdp_cst += what(namewr).eval()*what(namewr).eval() - ((Line *) aij)->wr*what(namewr).eval();
// sdp_cst += what(namewi).eval()*what(namewi).eval() - ((Line *) aij)->wi*what(namewi).eval();
// sdp_cst += ((Line *) aij)->wr - what(namewr).eval();
// sdp_cst += ((Line *) aij)->wi - what(namewi).eval();
}
}
}
sdpcut.print();
SDP.add_constraint(sdpcut <= 0);
sdpcut.print_expanded();
// SDP.add_constraint(sdpcut <= 0);


Constraint lin("lin"+to_string(numcuts));
cout << "\nbpairs size = " << b_pairs._keys.size() << endl;
// lin = (R_diff.in(b_pairs._keys) + R_Wij.in(b_pairs._keys) - R_hat.in(b_pairs._keys));
lin = product(R_diff.in(b_pairs._keys),(R_Wij.in(b_pairs._keys) - R_hat.in(b_pairs._keys)));
// lin += product(I_diff.in(b_pairs._keys),(Im_Wij.in(b_pairs._keys) - I_hat.in(b_pairs._keys)));
// lin += product(W_diff.in(b->_nodes),(Wii.in(b->_nodes) - W_hat.in(b->_nodes)));
SDP.add_constraint(lin <= 0);
lin += product(I_diff.in(b_pairs._keys),(Im_Wij.in(b_pairs._keys) - I_hat.in(b_pairs._keys)));
lin += product(W_diff.in(b->_nodes),(Wii.in(b->_nodes) - W_hat.in(b->_nodes)));
lin.print();
lin.print_expanded();
SDP.add_constraint(lin <= 0);

numcuts++;
}

// if(!bus_pairs_sdp._keys.empty()) {
// Constraint SOC_im("SOC_im");
// SOC_im = power(R_Wij, 2) + power(Im_Wij, 2) - Wii.from() * Wii.to();
// SDP.add_constraint(SOC_im.in(bus_pairs_sdp._keys) <= 0);
// bus_pairs_sdp._keys.clear();
// }

if(!bus_pairs_sdp._keys.empty()) {
Constraint SOC_im("SOC_im");
SOC_im = power(R_Wij, 2) + power(Im_Wij, 2) - Wii.from() * Wii.to();
SDP.add_constraint(SOC_im.in(bus_pairs_sdp._keys) <= 0);
bus_pairs_sdp._keys.clear();
}

// param<> R_Wij_star("R_Wij_star"), R_Wij_hat("I_Wij_hat"), R_diff("R_Diff");
// param<> Wii_star("Wii_star"), Wii_hat("Wii_hat"), W_diff("W_Diff");

for(auto& a: grid->arcs){
if(a->_imaginary && !a->_active) a->_free = true;
Expand Down
12 changes: 12 additions & 0 deletions include/gravity/func.h
Original file line number Diff line number Diff line change
Expand Up @@ -299,6 +299,18 @@ namespace gravity {
_coef = coef;
_p = new pair<param_*, param_*>(make_pair(p1,p2));
_sign = sign;
// if (coef->_is_transposed){
// p1->_is_vector=true;
// p2->_is_vector=true;
// }
// if(p1->_is_vector){
// coef->_is_transposed=true;
// p2->_is_vector=true;
// }
// if(p2->_is_vector){
// coef->_is_transposed=true;
// p1->_is_vector=true;
// }
if (coef->_is_transposed && p1->_is_transposed) {
throw invalid_argument("Check the transpose operator, there seems to be a dimension issue\n");
}
Expand Down
2 changes: 1 addition & 1 deletion src/constraint.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -132,9 +132,9 @@ size_t Constraint::get_id_inst(size_t ind) const{
/* Output */

void Constraint::print_expanded(){
eval();
auto nb_inst = get_nb_instances();
for (unsigned inst = 0; inst<nb_inst; inst++) {
eval(inst);
print(inst);
}
}
Expand Down
74 changes: 53 additions & 21 deletions src/func.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -856,6 +856,12 @@ namespace gravity{
res += poly_eval(_coef,i,j) * poly_eval(_p->first, i,j)* poly_eval(_p->second, i,j);
}
}
else if(_p->first->_is_transposed){
auto dim = _p->first->get_nb_instances(i);
for (int j = 0; j<dim; j++) {
res += poly_eval(_coef,i,j) * poly_eval(_p->first, i,j)* poly_eval(_p->second, i,j);
}
}
else {
res = poly_eval(_coef,i) * poly_eval(_p->first, i) * poly_eval(_p->second, i);
}
Expand Down Expand Up @@ -917,9 +923,15 @@ namespace gravity{
_coef = coef;
_p = p;
_sign = sign;
// if (coef->_is_transposed && p->_is_transposed) {
// throw invalid_argument("Check the transpose operator, there seems to be a dimension issue\n");
// if (coef->_is_transposed){
// _p->_is_vector = true;
// }
// if (_p->_is_vector) {
// _coef->_is_transposed = true;
// }
if (coef->_is_transposed && p->_is_transposed) {
throw invalid_argument("Check the transpose operator, there seems to be a dimension issue\n");
}
if (coef->is_function()) {
assert(((func_*)coef)->is_constant());
}
Expand Down Expand Up @@ -2406,37 +2418,53 @@ namespace gravity{
if (is_constant() && (c.is_sdpvar() || c.is_var() || (c.is_function() && !((func_*)&c)->is_constant()))) {
func_ f(c);
if (!f._cst->is_zero()) {
f._cst = multiply(f._cst, *this);
if (f._cst->is_function()) {
f.embed(*(func_*)f._cst);
auto fc = (func_*)f._cst;
*fc = (*this)* (*fc);
f.embed(*fc);
}
else {
f._cst = multiply(f._cst, *this);
}
}
for (auto &pair:*f._lterms) {
pair.second._coef = multiply(pair.second._coef, *this);
if (pair.second._coef->is_function()) {
auto fc = (func_*)pair.second._coef;
*fc = (*this)* (*fc);
f.embed(*fc);
}
else {
pair.second._coef = multiply(pair.second._coef, *this);
}
if (pair.second._coef->_is_transposed) {
pair.second._p->_is_vector = true;
}
if (pair.second._coef->is_function()) {
f.embed(*(func_*)pair.second._coef);
}
f.update_nb_instances(pair.second);

}
for (auto &pair:*f._qterms) {
pair.second._coef = multiply(pair.second._coef, *this);
if (pair.second._coef->is_function()) {
auto fc = (func_*)pair.second._coef;
*fc = (*this)* (*fc);
f.embed(*fc);
}
else {
pair.second._coef = multiply(pair.second._coef, *this);
}
if (pair.second._coef->_is_transposed) {
pair.second._p->first->_is_vector = true;
pair.second._p->second->_is_vector = true;
}
if (pair.second._coef->is_function()) {
f.embed(*(func_*)pair.second._coef);
}
update_nb_instances(pair.second);
}
for (auto &pair:*f._pterms) {
pair.second._coef = multiply(pair.second._coef, *this);
if (pair.second._coef->is_function()) {
f.embed(*(func_*)pair.second._coef);
auto fc = (func_*)pair.second._coef;
*fc = (*this)* (*fc);
f.embed(*fc);
}
else {
pair.second._coef = multiply(pair.second._coef, *this);
}
// update_nb_instances(pair.second); // TODO
}
Expand Down Expand Up @@ -4866,7 +4894,7 @@ namespace gravity{
return res;
}
else {
return func_(c2) *= c1;
return func_(c1) *= c2;
}
}
else {
Expand All @@ -4878,12 +4906,12 @@ namespace gravity{
//// delete new_c2;
// return func_(c1) *= c2;
// }
if (c2.is_var()) {
return func_(c2) *= c1;
}
else {
// if (c2.is_var()) {
// return func_(c2) *= c1;
// }
// else {
return func_(c1) *= c2;
}
// }
}
}

Expand Down Expand Up @@ -7297,6 +7325,8 @@ namespace gravity{
else {
res = f1.tr()*f2.vec();
}
res._nb_instances = 1;
res._dim[0] = 1;
return res;
}

Expand All @@ -7318,7 +7348,9 @@ namespace gravity{
}
else {
res = (*(var<type>*)&p1).vec()*f.vec();
}
}
res._nb_instances = 1;
res._dim[0] = 1;
return res;
}

Expand Down

0 comments on commit 943dc40

Please sign in to comment.