My objective function, $f(r,k;s,d)$, takes three different expressions, denoted by $f^1$, $f^2$, and $f^3$, respecitvely, depending on the conditions for the arguments and parameters as follows:
$f^1(r,k;s,d) = \frac{3 s \left(d^2-k^2 r\right)+k^2 r (3 d+k r)-6 d s^2+3 s^3}{6 s^2}$ when $0\leq d\leq 1,s\geq 2 d,0\leq k<d,0\leq r\leq 1$,
$f^2(r,k;s,d) = \frac{d^3+3 d^2 (k (r-1)-s)+3 d \left(k^2 ((r-3) r+1)+2 k r s-2 (r-1) s^2\right)+k^3 (r ((r-1) r+3)-1)-3 k^2 r^2 s+3 (r-1) s^3}{6 (r-1) s^2}$ when $0\leq d\leq 1,s\geq 2 d,d\leq k\leq s,0\leq r<\frac{d}{k}$,
$f^3(r,k;s,d) = \frac{d^3+3 d^2 k r+3 d r \left(k^2 (r-1)-2 s^2\right)+k^3 (r-1)^2 r+3 r s^3}{6 r s^2}$ when $0\leq d\leq 1,s\geq 2 d,d\leq k\leq s,\frac{d}{k}\leq r\leq 1$.
My ultimate goal is to find the optimal value of $r$ and $k$ at which the objective function $f$ gets maximized over all the three cases. It has turned out that working with analytical solutions is quite daunting. So I'm working with numerical simulations limited to the parameter values of $d \in [0,1]$ and $s \in [0,2]$.
My strategy is as follows:
First, construct a code that solves the maximization problem for each of the above three functions and computes the maximized value of the function and the corresponding optimal $r$ and $k$ (Obviously, these will be a function of $s$ and $d$.) For instance, for the first case of $f^1$, let me denote them by $f^{1*}$, $r^{1*}$, and $k^{1*}$; and the same for the cases of $f^2$ and $f^3$.
Second, using Plot3D, plot $f^{1*}$, $f^{2*}$, and $f^{3*}$ in one same diagram against $d \in [0,1]$ and $s \in [0,2]$; similarly, plot $r^{1*}$, $r^{2*}$, and $r^{3*}$ in one same diagram against $d \in [0,1]$ and $s \in [0,2]$, and plot $k^{1*}$, $k^{2*}$, and $k^{3*}$ in one same diagram against $d \in [0,1]$ and $s \in [0,2]$.
Lastly, since the first diagram allows comparing $f^{1*}$, $f^{2*}$, and $f^{3*}$, I can get the maximum value of the objective function $f$, and by looking at the second and third diagrams, I can get the corresponding optimal $r$ and $k$, which are what I'm ultimately looking for.
My Mathematica code for doing this is as follows (which is an extension of Alex Trounev's answer in 3D-Plot optimization results for varying parameter values):
Block[{t = 0}, f1 = (k^2 r (3 d + k r) + 3 (d^2 - k^2 r) s - 6 d s^2 + 3 s^3)/(6 s^2); max1 = Flatten[Table[{d, s, MaxValue[{f1, 0 <= d <= 1, s >= 2 d, 0 <= k < d, 0 <= r <= 1}, {k, r}]}, {d, 0, 1, .1}, {s, 0, 2, .1}], 1]; maxk1 = Flatten[Table[{d, s, k /. Last@Maximize[{f1, 0 <= d <= 1, s >= 2 d, 0 <= k < d, 0 <= r <= 1}, {k, r}]}, {d, 0, 1, .1}, {s, 0, 2, .1}], 1]; maxr1 = Flatten[Table[{d, s, r /. Last@Maximize[{f1, 0 <= d <= 1, s >= 2 d, 0 <= k < d, 0 <= r <= 1}, {k, r}]}, {d, 0, 1, .1}, {s, 0, 2, .1}], 1];] Block[{t = 0}, f2 = 1/(6 (-1 + r) s^2) (d^3 + k^3 (-1 + r (3 + (-1 + r) r)) + 3 d^2 (k (-1 + r) - s) - 3 k^2 r^2 s + 3 (-1 + r) s^3 + 3 d (k^2 (1 + (-3 + r) r) + 2 k r s - 2 (-1 + r) s^2)); max2 = Flatten[Table[{d, s, MaxValue[{f2, 0 <= d <= 1, s >= 2 d, d <= k <= s, 0 <= r < d/k}, {k, r}]}, {d, 0, 1, .1}, {s, 0, 2, .1}], 1]; maxk2 = Flatten[Table[{d, s, k /. Last@Maximize[{f2, 0 <= d <= 1, s >= 2 d, d <= k <= s, 0 <= r < d/k}, {k, r}]}, {d, 0, 1, .1}, {s, 0, 2, .1}], 1]; maxr2 = Flatten[Table[{d, s, r /. Last@Maximize[{f2, 0 <= d <= 1, s >= 2 d, d <= k <= s, 0 <= r < d/k}, {k, r}]}, {d, 0, 1, .1}, {s, 0, 2, .1}], 1];] Block[{t = 0}, f3 = (d^3 + 3 d^2 k r + k^3 (-1 + r)^2 r + 3 r s^3 + 3 d r (k^2 (-1 + r) - 2 s^2))/(6 r s^2); max3 = Flatten[Table[{d, s, MaxValue[{f3, 0 <= d <= 1, s >= 2 d, d <= k <= s, d/k <= r <= 1}, {k, r}]}, {d, 0, 1, .1}, {s, 0, 2, .1}], 1]; maxk3 = Flatten[Table[{d, s, k /. Last@Maximize[{f3, 0 <= d <= 1, s >= 2 d, d <= k <= s, d/k <= r <= 1}, {k, r}]}, {d, 0, 1, .1}, {s, 0, 2, .1}], 1]; maxr3 = Flatten[Table[{d, s, r /. Last@Maximize[{f3, 0 <= d <= 1, s >= 2 d, d <= k <= s, d/k <= r <= 1}, {k, r}]}, {d, 0, 1, .1}, {s, 0, 2, .1}], 1];] {ListPlot3D[max1, max2, max3, AxesLabel -> {"d", "s", "V"}], ListPlot3D[maxk1, maxk2, maxk3, PlotRange -> {0, 2}, AxesLabel -> {"d", "s", "k"}], ListPlot3D[maxr1, maxr2, maxr3, PlotRange -> {0, 1}, AxesLabel -> {"d", "s", "r"}]}
(In the code, you can simply by-pass the meaning of $t=0$.)
When this code is run, I get a list of error messages. Can anyone help? Thanks so much, in advance!
ListPlot3Dtakes all the data in a single argument, not a sequence, i.e.ListPlot3D[{a,b,c}]notListPlot3D[a, b, c]. $\endgroup$