In the notation of my "answer" to my still-unanswered question Almost simple Hermite interpolation, we can compute a quintic polynomial $l_a(x)$ such that $l_a(0) = l_a(1) = 0,$ $l_a(a) = 1,$ and $l_a'(0) = l_a'(1) = l_a'(a) = 0,$ where $a = \frac1{c + 1},$ so that we can consider, as a possible solution to the problem, at least for some values of $c$:
$$
f(x) = x + \left(\tfrac12 - a\right)l_a(x) \quad (0 \leqslant x \leqslant 1).
$$
After much simplification, we arrive at the formula
\begin{equation}
\label{3339606:eq:1}\tag{$1$}
\boxed{f(x) = x + \frac{(c^2 - 1)(c + 1)^2x^2(1 - x)^2[(3c - 2) - 2(c^2 - 1)x]}{2c^3}.}
\end{equation}
The appendix gives a range of values of $c$ for which this polynomial function satisfies the conditions of the question. For the moment, I'll just give two examples (excluding the trivial case $c = 1$):
When $c = 2,$
$$
f(x) = x + \frac{27x^2(1 - x)^2(2 - 3x)}8.
$$
From Wolfram Alpha:

When $c = \tfrac12,$
$$
f(x) = x + \frac{27x^2(1 - x)^2(1 - 3x)}8.
$$
From Wolfram Alpha:

The reason why these graphs are $180^\circ$ rotated images of one another is as follows:
By the uniqueness of the Hermite interpolating polynomial,
$$
l_{1 - a}(x) = l_a(1 - x)
\quad (0 < a < 1, \ 0 \leqslant x \leqslant 1).
$$
Write $c = (1 - a)/a,$ i.e. $a = 1/(c + 1).$ Then $a$ is related to
$c$ as $1 - a$ is to $1/c,$ and
$$
f_{1/c}(x) = x + \left(\tfrac12 - (1 - a)\right)l_{1 - a}(x)
= x - \left(\tfrac12 - a\right)l_a(1 - x)
= 1 - f_c(1 - x),
$$
where, for all $c > 0,$ $f_c(x)$ denotes the function in
\eqref{3339606:eq:1} with parameter $c.$ $\ \square$
Dropping the constraint $f'\left(\frac1{c + 1}\right) = 1$ simplifies the formula somewhat:
$$
f(x) = x + \frac{(c^2 - 1)(c + 1)^2x^2(1 - x)^2}{2c^2},
$$
but this doesn't greatly increase the range of usable values of $c.$ Also, the graphs take on a squashed appearance when $x$ approaches $1,$ as this example for $c = 2$ illustrates:

The graph for $c = 5/2,$ although still monotonic - unlike \eqref{3339606:eq:1}, in this case - is even worse:

So I won't consider this simplification any further.
If $f$ need not be analytic, and if continuous differentiability is enough, and if the values of $f'(0)$ and $f'(1)$ do not matter so long as they are equal, then it is easy to solve the problem using cubic splines. For example:
$$
f(x) =
\begin{cases}
\tfrac12g\left[(c + 1)x\right]
& \text{if } 0 \leqslant x \leqslant \frac1{c + 1}, \\
\tfrac12\left\{1 + g\left[\frac{(c + 1)x - 1}{c}\right]\right\}
& \text{if } \frac1{c + 1} \leqslant x \leqslant 1,
\end{cases}
$$
where
\begin{gather*}
g(t) = 3t^2 - 2t^3, \ g'(t) = 6t(1 - t)
\ \, (0 \leqslant t \leqslant 1), \\
g(0) = 0, \ g(1) = 1, \ g'(0) = g'(1) = 0, \\
g'(t) > 0 \ \, (0 < t < 1).
\end{gather*}
This has $f'(0) = f'(1) = f'\left(\frac1{c + 1}\right) = 0.$ It is valid for all $c > 0.$
But we can do a lot better than that. I delayed looking at this
possibility, wrongly imagining that it would work only for a
restricted range of values of $c,$ like the Hermite interpolation
solution. In fact, it works for all values of $c$ (the value of $c$
must be strictly positive, of course), without exception.
We continue to use the same "cardinal" cubic spline function $g,$
but now we define
$$
f(x) = x + \left(\tfrac12 - a\right)s_a(x)
\quad (0 \leqslant x \leqslant 1),
$$
where
$$
s_a(x) =
\begin{cases}
g\left(\frac{x}{a}\right)
& \text{ if } 0 \leqslant x \leqslant a,\\
g\left(\frac{1 - x}{1 - a}\right)
& \text{ if } a \leqslant x \leqslant 1.
\end{cases}
$$
Differentiating,
$$
s_a'(x) =
\begin{cases}
\frac1{a}g'\left(\frac{x}{a}\right)
& \text{ if } 0 \leqslant x \leqslant a,\\
-\frac1{1 - a}g'\left(\frac{1 - x}{1 - a}\right)
& \text{ if } a \leqslant x \leqslant 1.
\end{cases}
$$
Because $g'(t) \geqslant 0$ for all $t \in [0, 1],$ and
$$
\max_{0 \leqslant t \leqslant 1}g'(t) = \tfrac32,
$$
we have
\begin{gather*}
\min_{0 \leqslant x \leqslant 1}s_a'(x) = -\frac3{2(1 - a)}, \\
\max_{0 \leqslant x \leqslant 1}s_a'(x) = \frac3{2a}.
\end{gather*}
Therefore, if $a \leqslant \tfrac12,$ i.e. $c \geqslant 1,$
$$
\min_{0 \leqslant x \leqslant 1}f'(x) =
1 + \left(\frac12 - a\right)\left(-\frac3{2(1 - a)}\right) =
1 - \frac{3(1 - 2a)}{4(1 - a)} =
\frac{1 + 2a}{4(1 - a)} > 0.
$$
On the other hand, if $a \geqslant \tfrac12,$ i.e. $c \leqslant 1,$
$$
\min_{0 \leqslant x \leqslant 1}f'(x) =
1 + \left(\frac12 - a\right)\left(\frac3{2a}\right) =
1 - \frac{3(2a - 1)}{4a} =
\frac{3 - 2a}{4a} > 0.
$$
In all cases, therefore, $f$ is strictly increasing on $[0, 1].$
In terms of the constant $c = (1 - a)/a,$ the definition of
the function $f$ is:
\begin{equation}
\label{3339606:eq:2}\tag{$2$}
\boxed{f(x) =
\begin{cases}
x + \frac{c - 1}{2(c + 1)}g[(c + 1)x]
& \text{ if } 0 \leqslant x \leqslant \frac1{c + 1}, \\
x + \frac{c - 1}{2(c + 1)}g\left[\frac{(c + 1)(1 - x)}{c}\right]
& \text{ if } \frac1{c + 1} \leqslant x \leqslant 1.
\end{cases}
}
\end{equation}
For example, when $c = 4,$ \eqref{3339606:eq:2} becomes:
$$
f(x) =
\begin{cases}
x + \frac{15}2x^2(3 - 10x)
& \text{ if } 0 \leqslant x \leqslant \frac15, \\
x + \frac{15}{64}(1 - x)^2(1 + 5x)
& \text{ if } \frac15 \leqslant x \leqslant 1.
\end{cases}
$$
From Wolfram Alpha:

Here is a closer look at the knot of that spline function:

Appendix
For any $a$ such that $0 < a < 1,$ we define the quintic polynomial
function
$$
l_a(x) = \frac{x^2(1 - x)^2[a(3 - 5a) - 2(1 - 2a)x]}{a^3(1 - a)^3}.
$$
Its derivative is given by
$$
l_a'(x) =
\frac{2x(1 - x)(x - a)[5(1 - 2a)x - (3 - 5a)]}{a^3(1 - a)^3}.
$$
It satisfies (and indeed it is uniquely determined by) six
constraints
\begin{gather*}
l_a(0) = l_a(1) = 0, \ l_a(a) = 1, \\
l_a'(0) = l_a'(1) = l_a'(a) = 0.
\end{gather*}
Writing $c = (1 - a)/a,$ or equivalently $a = 1/(c + 1),$ where $c$
is any strictly positive number, we define
$$
f_c(x) = x + \left(\tfrac12 - a\right)l_a(x).
$$
Then
\begin{gather*}
f_c(0) = 0, \ f_c(1) = 1, \
f_c\left(\frac1{c + 1}\right) = \frac12, \\
f_c'(0) = f_c'(1) = f_c'\left(\frac1{c + 1}\right) = 1.
\end{gather*}
I shall determine a set of values of $c$ such that $f_c'(x) > 0$ for
all $x \in [0, 1].$ (I shall not try to determine all such values of
$c.$) It was shown above that for all $c > 0,$ if either of
$f_c',$ $f_{1/c}'$ is strictly positive on $[0, 1],$ then so is the
other. Because $f_1(x) = x,$ it suffices to consider only the case
$c > 1,$ i.e., $a < \tfrac12.$
Differentiating:
\begin{equation}
\label{3339606:eq:3}\tag{$3$}
1 - f_c'(x) = -\left(\tfrac12 - a\right)l_a'(x) =
\frac{20\left(\tfrac12 - a\right)^2x(1 - x)(x - a)(b - x)}
{a^3(1 - a)^3},
\end{equation}
where
$$
b = \frac{3 - 5a}{5(1 - 2a)} =
\tfrac12\cdot\frac{\tfrac35 - a}{\tfrac12 - a} =
\tfrac12\left(1 + \frac{\tfrac1{10}}{\tfrac12 - a}\right),
$$
i.e.,
$$
\left(\frac12 - a\right)\left(b - \frac12\right) = \frac1{20},
$$
so we can rewrite \eqref{3339606:eq:3} as
\begin{equation}
\label{3339606:eq:4}\tag{$4$}
1 - f_c'(x) =
\frac{\left(\tfrac12 - a\right)x(1 - x)(x - a)(b - x)}
{a^3(1 - a)^3\left(b - \tfrac12\right)}.
\end{equation}
We are interested in determining $c > 1$ such that $1 - f_c'(x) < 1$
for all $x \in [0, 1].$ By \eqref{3339606:eq:4}, we only need to
consider $x$ such that $a < x < \min\{1, b\}.$
Case 1:
$$
\boxed{c \leqslant \tfrac32 \iff a \geqslant \tfrac25 \iff
\tfrac12 - a \leqslant \tfrac1{10} \iff b \geqslant 1.}
$$
Looking at the factors in \eqref{3339606:eq:4}, we have:
\begin{gather*}
\frac{b - x}{b - \tfrac12} = 1 + \frac{\tfrac12 - x}{b - \tfrac12}
\leqslant 1 + \frac{\tfrac12 - x}{1 - \tfrac12} = 2(1 - x)
\leqslant \frac65, \\
\frac12 - a \leqslant \frac1{10}, \\
x \leqslant 1, \\
(1 - x)(x - a) \leqslant \left(\frac{1 - a}2\right)^2
\leqslant \left(\frac3{10}\right)^2, \\
a(1 - a) = \frac14 - \left(\frac12 - a\right)^2
\geqslant \frac6{25},
\end{gather*}
therefore
$$
1 - f_c'(x) \leqslant
\frac{\tfrac65\cdot\tfrac1{10}\cdot\left(\tfrac3{10}\right)^2}
{\left(\tfrac6{25}\right)^3} = \frac{5^2}{2^5} = \frac{25}{32} < 1.
$$
This completes the proof that $f_c'(x) > 0$ for $x \in [0, 1]$
and $c \in \left[\tfrac23, \tfrac32\right].$ $\ \square$
Case 2:
$$
\boxed{c \geqslant \tfrac32 \iff a \leqslant \tfrac25 \iff
\tfrac12 - a \geqslant \tfrac1{10} \iff b \leqslant 1.}
$$
From \eqref{3339606:eq:3}, using the inequalities
$x(1 - x) \leqslant \tfrac14$ and
$(x - a)(b - x) \leqslant ((b - a)/2)^2,$
$$
1 - f_c'(x) \leqslant
\frac{5\left(\tfrac12 - a\right)^2(b - a)^2}{4a^3(1 - a)^3}.
$$
Reparameterising in terms of $p,$ where
$$
a = \frac12 - p,
\quad 1 - a = \frac12 + p,
\quad b - \frac12 = \frac1{20p}
\qquad \left(\frac1{10} \leqslant p < \frac12\right),
$$
we have
$$
1 - f_c'(x) \leqslant
\frac{5p^2(p + 1/20p)^2}{4\left(\frac14 - p^2\right)^3} =
\frac{5(p^2 + 1/20)^2}{4\left(\frac14 - p^2\right)^3}
= \frac{5\left(\frac3{10} - q\right)^2}{4q^3}
= \frac{(3 - 10q)^2}{80q^3},
$$
where
$$
q = \frac14 - p^2 = a(1 - a) \in \left(0, \, \frac6{25}\right)\!.
$$
According to Wolfram Alpha, the cubic equation
$80q^3 = (3 - 10q)^2$ has a single real root,
$$
q_0 \bumpeq 0.212428328248244.
$$
We therefore have $f_c'(x) > 0$ for all $x \in [0, 1]$ if any of the
following list of equivalent conditions is satisfied:
\begin{align*}
q > q_0 & \iff \frac{c}{(c + 1)^2} > q_0 \\
& \iff c^2 - 2\left(\frac1{2q_0} - 1\right)c + 1 < 0 \\
& \iff c < c_0 = \left(\frac1{2q_0} - 1\right) +
\sqrt{\left(\frac1{2q_0} - 1\right)^2 - 1}
\bumpeq 2.266203431.
\end{align*}
Finally, then: $f_c'(x) > 0$ for all $x \in [0, 1]$ if
$c \in [0.4413, 2.2662]$. $\ \square$
The graph of $f_c$ for $c = c_0$ looks like this:

Here is a close-up view of the flattest part of that graph:

This estimated value, $c_0,$ is evidently quite close to the least
upper bound of the set (presumably a closed interval) of all values
of $c$ for which $f_c$ is monotonic.
By "non-piecewise" I probably meant to say that I don't want it to be non-smooth at any point, but is that just implied by calling the function smooth?
– Acreol Aug 31 '19 at 00:29