Finding maximum of function I'm struggling to find the maximum of this function $f:\mathbb{R}^n\times\mathbb{R}^n \rightarrow \mathbb{R}$
$$ f(x,y) = \frac{n+1}{2} \sum_{i=1}^n x_i\,y_i - \sum_{i=1}^n x_i \sum_{i=1}^n y_i,$$
where $x_i,y_i\in[0,1]$ for $i=1,...,n$. It reminded me to Chebyshev's sum inequality, but it didn't help much since I can't sort the variables $x$ and $y$. Any help is welcomed.
 A: Here is a stab.  Let $x_i=y_i=1$ for $i\le k$ and $x_i=y_i=0$ for $i>k$.  Then $f=(n+1)k/2-k^2 = ((n+1)/2-k)k$ which is largest when $k\approx n/4$, giving $f\approx n^2/16$. 
This is only a lower bound on the desired maximum. It does not purport to actually answer the original question.   
ADDED several hours later, and again the next day:  Towards a proof.
Note first that the domain of $f$ is the compact set $[0,1]^{2n}$, so the continuous function $f$ attains its maximal value. The "stab" examples above show this maximum value is positive.
Following zwim, look at the restricted problem where $X=\sum x_i$ and $Y=\sum y_i$ are fixed.  Note that at a maximum both $X<n$ and $Y<n$ because otherwise $f\le0$. Then $f(x,y)=(n+1)\sum x_iy_i/2 - XY$ is maximized when first, the $x_i$ and the $y_i$ are ordered the same (lets say non-increasing, so $1\ge x_1\ge x_2\ge  \cdots \ge x_n\ge0$, and similarly for the $y_i$) and second, as many of the $x_i$ and $y_i$ values are equal to $1$ as can be, consistent with the $X$ and $Y$ sum constraints.  Given that, we have $x_1=x_2=\cdots x_k = 1$ and $ x_{k+1}=u$ with $0\le u<1$ and the rest of the $x_i=0$; similarly $y_1=y_2=\cdots = y_l=1$, $y_{l+1}=v$ with $0\le v <1$, and the rest of the $y_i=0$. (The cases $k=n$ or $l=n$ cannot occur if $X<n$ and $Y<n$.) Clearly $k = \lfloor X\rfloor,$ $l = \lfloor Y\rfloor$, $u = X-k$, $ v=Y-l,$ and $f(x,y) = (n+1)(k+u y_{k+1})/2 - XY$, if we assume $X\le Y$ so $k\le l$.  Denote this last expression by $g(X,Y)$, as it depends only on $X$ and $Y$.  For given $X$ and $Y$, if $l>k$ we can replace the $y_i$ values for $i>k+1$ with $0$, which leaves the positive term $\sum x_i y_i$ unaffected but possibly reduces $Y$ and its  negative effect in the formula for $f$. Thus we may as well assume that $l=k$, and then $g(X,Y') = (n+1) (k+uv)/2 - (k+u)(k+v)$, where $Y'$ is the reduced $Y$.  Now, for given $k$, this expression is bi-linear in $u$ and $v$, which is maximized at one of the 4 corners of the $(u,v)$ square $[0,1)^2$.  Relaxing the $(u,v)$ maximization to the domain $[0,1]^2$, we see the maximal value of $f$ is attained as one of
$$f=(n+1)(k+1)/2-(k+1)^2\tag 1$$ $$f=(n+1)k/2 - (k+1)k\tag 2$$ or $$f=(n+1)k/2 - k^2.\tag 3$$
Clearly (2) is not as big as (3), so it is out of the running.  And (1) is the same as (3) with a different value of $k$.  The conclusion:  it looks like  the bound I showed is in fact the correct answer to the original question.  This is borne out by computer experiments, which suggest that $k=\lfloor (n+2)/4\rfloor$ and $(u,v)=(0,0)$ are always the optimizers.
ADDED  18 Sept, revised 19 Sept: So the optimum is achieved by any integer $k$ maximizing (3) above,  subject to $0\le k \le n$.  If we drop the integrality constraint the optimum is achieved by $k=(n+1)/4$, and with the integrality constraint by the (or more precisely, any) integer closest to $(n+1)/4$.  One formula for such a $k$ is  $k= \lfloor(n+2)/4\rfloor$ and another is $k=\lfloor(n+3)/4\rfloor.$
(Other formulas are possible using other symbols.)
When $n\equiv1 \pmod 4$ these formulas deliver distinct  integer values equally close to $(n+1)/4$.  If $n\equiv 0,2,\text{ or } 3\pmod 4$ the two formulas agree.  Consider the example the cases $n=100,$ $101,$ $102,$ and $103$.  In these cases we need $k$ to be a closest integer to $25.25,$ $25.5,$$ 25.75,$ and $26,$ respectively. Clearly the corresponding values of $k$ are: $25$,  $25$ or $26$, $26,$ and $26$.  The $\lfloor(n+2)/4\rfloor$ formula delivers $25,$ $25,$ $26,$ $ 26$, and the $\lfloor(n+3)/4\rfloor$ formula delivers $25,$ $26,$ $26,$ $26$, respectively.
Thank you, user480959 for a most interesting puzzle, and zwim for suggesting a line of attack.
A: Let's have $\bar x=\frac 1n\sum\limits_{i=1}^n x_i$ and $\bar y=\frac 1n\sum\limits_{i=1}^n y_i$.
$\sum\limits_{i=1}^n x_iy_i=\sum\limits_{i=1}^n (x_i-\bar x)(y_i-\bar y)+\overbrace{\sum\limits_{i=1}^n x_i\bar y}^{n\bar x\bar y}+\overbrace{\sum\limits_{i=1}^n \bar xy_i}^{n\bar x\bar y}-\overbrace{\sum\limits_{i=1}^n \bar x\bar y}^{n\bar x\bar y}=\sum\limits_{i=1}^n (x_i-\bar x)(y_i-\bar y)+n\bar x\bar y$
$\sum\limits_{i=1}^n x_i\sum\limits_{i=1}^n y_i=n^2\bar x\bar y$
So $f(x,y)=\frac{n+1}2\sum\limits_{i=1}^n (x_i-\bar x)(y_i-\bar y)+\bar x\bar y\overbrace{(\frac{n(n+1)}2-n^2)}^{\frac{-n(n-1)}2}$
When fixing the averages, this is maximized when $\sum\limits_{i=1}^n (x_i-\bar x)(y_i-\bar y)$ is maximized. 
By Cauchy-Schwarz it has a higher bound when $(y_i-\bar y)=\alpha(x_i-\bar x),\alpha>0$ which is $\alpha\sum\limits_{i=1}^n (x_i-\bar x)^2$.
Note: that we can now apply the Chebyshev inequality since we don't lose in generality by reordering the $x_i$, and the condition above reorders the $y_i$ in the same way.
By convexity of $x^2$ this is maximized when $x_i=0$ or $x_i=1$ (both possible since some $x_i\le \bar x$ and some $x_i\ge \bar x$ by definition).
So intuitively we fall back to kimchi lover solution, where $\alpha=1$ and the $x_i,y_i$ take only $0,1$ values, but I'm not sure if we didn't miss something by fixing $\bar x,\bar y$ like this, if $n$ is large then constraining $\bar x$ to take rationnal values $\frac kn$ is quite lax, but for small $n$ this troubles me...
