# solve L_inf--L_inf norm iteratively by solving
# simplified problem first


option solver loqo;
param segments := 8;
param stopSteps := 128;
param passSteps := 64;
param stopBands := 3;

var square;
var solution {i in 1..segments};

param prev_a {i in 1..segments};
param stopPosition {i in 1..stopSteps, j in 1..segments, k in 
1..stopBands};
param stopValue {i in 1..stopSteps, j in 1..segments, k in 1..stopBands};
let {i in 1..stopSteps, j in 1..segments, k in 1..stopBands} 
stopPosition[i,j,k] := 3.1415926536*j*(i+64+128+256*(k-1))/1024;
let {i in 1..stopSteps, j in 1..segments, k in 1..stopBands} 
stopValue[i,j,k] := 
sin(stopPosition[i,j,k])*sin(stopPosition[i,j,k])
/(stopPosition[i,j,k]*stopPosition[i,j,k]);
param passPosition {i in 1..passSteps, j in 1..segments};
param passValue {i in 1..passSteps, j in 1..segments};
let {i in 1..passSteps, j in 1..segments} passPosition[i,j] := 
3.1415926536*j*(i)/1024;
let {i in 1..passSteps, j in 1..segments} passValue[i,j] := 
sin(passPosition[i,j])*sin(passPosition[i,j])
/(passPosition[i,j]*passPosition[i,j]);

var a {i in 1..segments};
var alpha;

maximize maxMin: alpha;

subject to stop1 {i in 1..(stopSteps/2 -1), k in 1..stopBands} :
 -1  <= 
 sum {m in 1..segments} (a[m]*stopValue[i,m,k])
/ sum {j in 1..segments} (prev_a[j]*passValue[64-i,j] ) 
 <= 1;

subject to stop2 {i in (stopSteps/2 +1)..stopSteps, k in 1..stopBands} :
 -1  <= 
  sum {m in 1..segments} (a[m]*stopValue[i,m,k]) 
/  sum {j in 1..segments} (prev_a[j]*passValue[i-64,j] ) 

<= 1;


subject to pass {i in 1..passSteps} :
 alpha <= sum {j in 1..segments} (a[j]*passValue[i,j] ) ;

subject to positive: a[segments] >= 0;



data;
param prev_a :=

1  -133.576
2  -255.048
3   123.333
4   721.332
5  2051.37
6  2322.58
7  1732.38
8  1500.89
;

 
 
solve;
 

display alpha / (sum {i in 1..segments} abs(a[i]));

display {j in 1..segments} a[j] / (sum {i in 1..segments} abs(a[i]));
display a;
display alpha;
 
display {j in 1..segments} a[j] / (a[1]);