# Find maximum with a binary and linear restriction

I'm trying to find the maximum of a function $\texttt{vm1}$, depending of a variable $\texttt{delta} \in [0,1)$.

The function $\texttt{f}$ was created with the restriction that the files of the resulting matrix sum 1. The function $\texttt{g}$ was created with the restriction that only binary variable could be in the resulting matrix. (thanks to @m_goldberg for helped me to write $\texttt{f}$ and $\texttt{g}$)

My code is:

um1 = {0, 0, 0, 1, 1, 2, 2};
um2 = {0, 2, 1, 0, 1, 0, 2};
uw1 = {0, 0, 2, 0, 2, 1, 1};
uw2 = {0, 1, 0, 2, 2, 0, 1};

f[p : {_?NumericQ ..}, rows_Integer: 7, cols_Integer: 3] :=
Module[{m}, m = ArrayReshape[p, {rows, cols}];
If[Total[m, {2}] == ConstantArray[1, rows], m,
Defer[f[p, rows, cols]]]]

g[p : {_Integer ..}, rows_Integer: 7, cols_Integer: 3] :=
Module[{m}, m = Pick[p, p, 0 | 1];
If[m != p, Defer[g[p, rows, cols]], ArrayReshape[p, {rows, cols}]]]

qm1[f_] := {{f[[1, 3]], 0, 0, f[[1, 2]], 0, f[[1, 1]], 0},
{0, f[[2, 3]], 0, f[[2, 2]], 0, 0, f[[2, 1]]},
{0, 0, f[[3, 3]], 0, f[[3, 2]], f[[3, 1]], 0},
{f[[4, 3]], 0, 0, f[[4, 2]], 0, f[[4, 1]], 0},
{0, 0, f[[5, 3]], 0, f[[5, 2]], f[[5, 1]], 0},
{f[[6, 3]], 0, 0, f[[6, 2]], 0, f[[6, 1]], 0},
{0, f[[7, 3]], 0, f[[7, 2]], 0, 0, f[[7, 1]]}}

qm2[f_] := {{f[[1, 3]], f[[1, 2]], f[[1, 1]], 0, 0, 0, 0},
{f[[2, 3]], f[[2, 2]], f[[2, 1]], 0, 0, 0, 0},
{f[[3, 3]], f[[3, 2]], f[[3, 1]], 0, 0, 0, 0},
{0, f[[4, 2]], 0, f[[4, 3]], f[[4, 1]], 0, 0},
{0, f[[5, 2]], 0, f[[5, 3]], f[[5, 1]], 0, 0},
{0, 0, f[[6, 1]], 0, 0, f[[6, 3]], f[[6, 2]]},
{0, 0, f[[7, 1]], 0, 0, f[[7, 3]], f[[7, 2]]}}

qw1[f_] := {{f[[1, 3]], 0, f[[1, 2]], 0, 0, f[[1, 1]], 0},
{0, f[[2, 3]], f[[2, 2]], 0, 0, 0, f[[2, 1]]},
{f[[3, 3]], 0, f[[3, 2]], 0, 0, f[[3, 1]], 0},
{0, 0, 0, f[[4, 3]], f[[4, 2]], f[[4, 1]], 0},
{0, 0, 0, f[[5, 3]], f[[5, 2]], f[[5, 1]], 0},
{f[[6, 3]], 0, f[[6, 2]], 0, 0, f[[6, 1]], 0},
{0, f[[7, 3]], f[[7, 2]], 0, 0, 0, f[[7, 1]]}}

qw2[f_] := {{f[[1, 3]], f[[1, 2]], 0, f[[1, 1]], 0, 0, 0},
{f[[2, 3]], f[[2, 2]], 0, f[[2, 1]], 0, 0, 0},
{0, f[[3, 2]], f[[3, 3]], 0, f[[3, 1]], 0, 0},
{f[[4, 3]], f[[4, 2]], 0, f[[4, 1]], 0, 0, 0},
{0, f[[5, 2]], f[[5, 3]], 0, f[[5, 1]], 0, 0},
{0, 0, 0, f[[6, 1]], 0, f[[6, 3]], f[[6, 2]]},
{0, 0, 0, f[[7, 1]], 0, f[[7, 3]], f[[7, 2]]}}

am1[gw1_, gw2_] := {{1, 0, 0, gw2[[1, 1]], 0, gw1[[1, 1]], 0},
{0, 1, 0, gw2[[2, 1]], 0, 0, gw1[[2, 1]]},
{0, 0, 1, 0, gw2[[3, 1]], gw1[[3, 1]], 0},
{1, 0, 0, gw2[[4, 1]], 0, gw2[[4, 1]], 0},
{0, 0, 1, 0, gw2[[5, 1]], gw1[[5, 1]], 0},
{1, 0, 0, gw2[[6, 1]], 0, gw1[[6, 1]], 0},
{0, 1, 0, gw2[[7, 1]], 0, 0, gw1[[7, 1]]}}

am2[gw1_, gw2_] := {{1, gw2[[1, 2]], gw1[[1, 2]], 0, 0, 0, 0},
{1, gw2[[2, 2]], gw1[[2, 2]], 0, 0, 0, 0},
{1, gw2[[3, 2]], gw1[[3, 2]], 0, 0, 0, 0},
{0, gw2[[4, 2]], 0, 1, gw1[[4, 2]], 0, 0},
{0, gw2[[5, 2]], 0, 1, gw1[[5, 2]], 0, 0},
{0, 0, gw1[[6, 2]], 0, 0, 1, gw2[[6, 2]]},
{0, 0, gw1[[7, 1]], 0, 0, 1, gw2[[7, 2]]}}

aw1[gm1_, gm2_] := {{1, 0, gm2[[1, 1]], 0, 0, gm1[[1, 1]], 0},
{0, 1, gm2[[2, 1]], 0, 0, 0, gm1[[2, 1]]},
{1, 0, gm2[[3, 1]], 0, 0, gm1[[3, 1]], 0},
{0, 0, 0, 1, gm2[[4, 1]], gm1[[4, 1]], 0},
{0, 0, 0, 1, gm2[[5, 1]], gm1[[5, 1]], 0},
{1, 0, gm2[[6, 1]], 0, 0, gm1[[6, 1]], 0},
{0, 1, gm2[[7, 1]], 0, 0, 0, gm1[[7, 1]]}}

aw2[gm1_, gm2_] := {{1, gm2[[1, 2]], 0, gm1[[1, 2]], 0, 0, 0},
{1, gm2[[2, 2]], 0, gm1[[2, 2]], 0, 0, 0},
{0, gm2[[3, 2]], 1, 0, gm1[[3, 2]], 0, 0},
{1, gm2[[4, 2]], 0, gm1[[4, 2]], 0, 0, 0},
{0, gm2[[5, 2]], 1, 0, gm1[[5, 2]], 0, 0},
{0, 0, 0, gm1[[6, 2]], 0, 1, gm2[[6, 2]]},
{0, 0, 0, gm1[[7, 2]], 0, 1, gm2[[7, 2]]}}

qmatrix[fm1_, fm2_, fw1_, fw2_, gm1_, gm2_, gw1_, gw2_] :=
qm1[fm1]*am1[gw1, gw2]*1/4 + qm2[fm2]*am2[gw1, gw2]*1/4 +
qw1[fw1]*aw1[gm1, gm2]*1/4 + qw2[fw2]*aw2[gm1, gm2]*1/4

vm1[am1_, am2_, aw1_, aw2_, bm1_, bm2_, bw1_, bw2_, delta_] :=
Inverse[IdentityMatrix[7] -
delta*
qmatrix[f[am1, 7, 3], f[am2, 7, 3], f[aw1, 7, 3], f[aw2, 7, 3],
g[bm1, 7, 3], g[bm2, 7, 3], g[bw1, 7, 3],
g[bw2, 7, 3]]].qmatrix[f[am1, 7, 3], f[am2, 7, 3], f[aw1, 7, 3],
f[aw2, 7, 3], g[bm1, 7, 3], g[bm2, 7, 3], g[bw1, 7, 3],
g[bw2, 7, 3]].um1


All the previous definitions were tested, but when I tried to find the maximum of $\texttt{v1m}$, fixing other variables, it appears an error.

For example:

ss={0, 1, 1, 0, 1, 1, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0}
dd= ConstantArray[1/3, {1, 21}][[1]]
FindMaximum[
vm1[x, dd, dd, dd, ss, ss, ss, ss, delta], x]


This give me the error: "Part specification f[x,7,3][[1,3]] is longer than depth of object". I don't know why.

I apologize if there's anything grossly written.

I would be very very grateful if someone would help or advise me. Thanks

• Before anybody wastes his or her time on this: Why are you sure that this somewhat unconventional system should have solutions? What are you doing there? What's the context? – Henrik Schumacher Apr 7 '18 at 21:59
• Hi, I will put context on the question, sorry about that. I know that the system has at least one solution by some theorem that I checked before. – hllspwn Apr 7 '18 at 22:02
• Moreover, Solve has no chance to do anything because f evaluates only for numerical vectors as arguments. Solve is a symbolic solver, primarily for polynomial equations. Maybe FindRoot would be a better choice. Moreover, I have my doubts that having FindMaximum in the equations is helpful. Instead, one could use the first order optimality conditions as equations. (It would be good to know if the respective objective is convex with respect to the optimization variables.) – Henrik Schumacher Apr 7 '18 at 22:02
• Is there a paper on this topic? You might include a link to it. – Henrik Schumacher Apr 7 '18 at 22:04
• There is no paper since this is personal research. But I use a book from Putterman ("Markov Desicion Processes") to see Markov systems. I try to explain what I'm doing, sorry if I didn't express myself correctly. – hllspwn Apr 7 '18 at 22:14