3D Piecewise Plot with varying bounds in Tikz/PGF

Question

The function to be plotted reads as

In it, these constrains hold true:

If plotted it looks like

Please excuse the different symbols, I will explain: p/p_0 is equivalent to \Pi and \gamma is equivalent to \kappa.

The plot was generated using data from Matlab, fed into LaTeX/PGF via a CSV-file. This works fine, but of course the data cannot be modified anymore, so I wanted to have it all in Latex naturally. Plus, doing stuff in Latex is fun.

I defined a function using the math library of Tikz/PGF such that

\usepackage{tikz}
\usepackage{pgfplots}
\usetikzlibrary{math}
\tikzmath{%
function prcrit(\x) {%
    return (2/(\x+1))^(\x/(\x-1));
};
function psimax(\x) {%
    return (sqrt(\x/(\x+1))) * (2/(\x+1))^(1/(\x-1));   
};
function psi(\x,\y) {%
    \psimax = psimax(\x);%maximum psi only dependent on kappa / k
    \prcrit = prcrit(\x);%critical pressure ratio only dependent on kappa / k
        if \y<\prcrit then {\efflux = \psimax;} else {\efflux = sqrt(\x/(\x-1))*(\y)^(1/\x)*sqrt(1-\y^((\x-1)/\x));};
    return \efflux;
};
}%

in which prcrit (pressure ratio Pi) stands for $\Pi_{mathrm{crit}}$ from the above image, and psimax for the maximum value for \Psi (see first of the two cases of the initial equation). Thus, the function \Psi is set to its respective maximum \Psi_{\mathrm{max}} for any pressure ratio below critical.

The issue with this is that the boundaries that define the piecewise 3D function are themselves functions of one of the two input variables, namely \kappa, see the second equation.

Therefore, this thread and this one did not help me. This answer seems quite close to what my attempt is, but it does not seem to work with two variables. It also does not make use of the axis environment and merely draws a regular tikz picture. This is also the approach I found in the package's documentation.

Ideally, the end result could look similar to this:

\begin{tikzpicture}
    \begin{axis}
         \addplot3[domain=1.3:1.7,y domain=0:1,surf] {psi(\x,\y)};
    \end{axis}
\end{tikzpicture}

but that I could not get to work. The issue is essentially marrying the thus far functional function declaration in the preamble and a regular surf plot using the axis environment. I could not find such a solution yet.

The function also needs to be contained in quite narrow boundaries for \Pi and \kappa, namely 0 to 1 and roughly 1.3 to 1.7 respectively. This is what I would try to achieve by assigning corresponding domain boundaries. Towards 0 and 1 for \Pi, the function also has an infinite gradient (can be seen quite clearly in the function plot above), so that pgfplots also had some hiccups there in the past; but that is not so bad.

Answer 1

With declare function it is arguably easier this way (since you cannot use tikzmath easily with fpu). The value y=0 is a bit tricky, but since the true dependence on y is not too strong, I'd like to argue that it is good enough to start plotting at y=0.01 or so.

\documentclass[border=3.14mm,tikz]{standalone}
\usepackage{pgfplots}
\pgfplotsset{compat=1.16}
\begin{document}
\begin{tikzpicture}[declare function={%
prcrit(\x)=(2/(\x+1))^(\x/(\x-1));
psimax(\x)=(sqrt(\x/(\x+1))) * (2/(\x+1))^(1/(\x-1));
psi(\x,\y)=ifthenelse(\y<prcrit(\x),psimax(\x),
sqrt(\x/(\x-1))*((\y)^(1/\x))*sqrt(1-\y^((\x-1)/\x)));}]
\begin{axis}[view={45}{45},xlabel={$\gamma$},ylabel={$p/p_0$},zlabel={$\Psi(\gamma,p/p_0)$},
point meta min=-1]
  \addplot3[opacity=0.6,colormap/viridis,domain=1.3:1.7,y domain=0.01:1,surf] {psi(\x,\y)};
\end{axis}
\end{tikzpicture}
\end{document}