The following code draws the intersection of a surface N: eqNcar=0 with three other surfaces defined by eq1=0, eq2=0, eq3=0.

I would like to shade with dots a specific region on $$N$$ namely the region defined by $$\{eq1>0\}, \{eq2>0\} \cap \{eq3<0\}$$. I would like to use the option MeshShading but I don't understand what tensor I should put. The documentation says an array of depth $$k$$ (here $$k=3$$ because I have $$3$$ mesh functions) but such specification raises an error (like "bad tensor dimension"). Anyway I am fighting with what to put in the structure of MeshShading.

Clear["Global*"];

nn          = 3;
eps         = Table[Subscript[ep,i],{i,nn}];
xs          = Table[Subscript[x,i],{i,nn}];
alp         = 1/2;
bet         = 2;
Subscript[ep,2] = -alp Subscript[ep,1]/bet;
Subscript[ep,1] = 1;
Subscript[ep,3] = 0;

eqNcar      = -(1 + alp + bet)^(-1) - ((1 + alp + bet)^(-1) - Subscript[x, 1])^2 - ((1 + alp + bet)^(-1) - Subscript[x, 2])^2 + Subscript[x, 3];
eqNpar      = {w,s,w^2+s^2} + 1/(1+ alp+bet) {1,1,1};

eq1         = (-((1 + (1 + alp + bet)*s^2 + (1 + alp + bet)*w^2)*\
(bet*s + s^2 + w*(alp + w))) + 2*s*(1 + (1 + alp + bet)*s)*\
(s + alp*s^2 + w*(bet + alp*w)) + 2*w*(1 + (1 + alp + bet)*w)*\
(alp*s + w + bet*(s^2 + w^2)))/(1 + alp + bet);

eq2         = (alp*(-s + w)*((1 + alp + bet)^(-1) + s^2 + w^2) -
(2*w*((1 + alp + bet)^(-1) + w)*(-(alp^2*((1 + alp + bet)^(-1) + s)) +
bet*((1 + alp + bet)^(-1) + w)))/bet -
(2*s*((1 + alp + bet)^(-1) + s)*(-(alp*((1 + alp + bet)^(-1) + s)) +
bet^2*((1 + alp + bet)^(-1) + w)))/bet)*Subscript[ep, 1];

eq3 = (-3/(1 + alp + bet) + s + 2*w)*Subscript[ep, 1] +
((-2 + alp + bet)*(3 + (1 + alp + bet)^2*s^2 + (1 + alp + bet)*s*
(-3 + (1 + alp + bet)*w) + (1 + alp + bet)*w*
(-3 + (1 + alp + bet)*w)) + (1 + alp + bet)*(-3 + 2*(1 + alp + bet)*s +
(1 + alp + bet)*w)*Subscript[ep, 2])/(1 + alp + bet)^2;

pN = ParametricPlot3D[eqNpar,{w,-1,1},{s,-1,1},\
MeshStyle -> {{Blue},{Red},{Black}},\
{{Automatic,Automatic,Automatic}},\
{{Automatic,Automatic,Automatic}},\
},
Mesh -> {{0}},PlotStyle->Orange,PlotPoints->200,
BoundaryStyle->{1->None},\
AxesLabel->xs
];

cm              = 72/2.54;
image       = Rasterize[Show[pN,ImageSize->10 cm],"Image",ImageResolution->800];
Export["ex.pdf", Show[image, ImageSize -> 10 cm]]


In my example, MeshShading is an array $$M$$ of size $$(1,3,3)$$ and of depth $$4$$ (note the documentation would expect a array of depth $$3$$) and I was expecting the coefficient $$m_{1ij}$$ to be the color of the region between the ith meshcurves and jth mesh curve but it is not the case.

ps: I would like the method to be easily adaptable for a number of mesh functions up to $$10$$.

• What is xf in the code? Jul 11, 2023 at 14:24
• I corrected, it's eqNpar Jul 11, 2023 at 14:26
• Since your three MeshFunctions are difference from the three functions eq1,eq2,eq3, we can only use RegionFunction -> Function[{x, y, z, w, s}, eq1 > 0 && eq2 > 0 && eq3 < 0]. to get the region eq1 > 0 && eq2 > 0 && eq3 < 0. and the rest region by RegionFunction -> Function[{x, y, z, w, s}, ! (eq1 > 0 && eq2 > 0 && eq3 < 0)]. Jul 11, 2023 at 14:54
• If ms is the MeshShading tensor, then its depth should equal the number of mesh functions. At level j, ms[[..., j0, ...]] is a directive that applies when the $j$-th mesh function lies between j0-1 and j0 values of Mesh division points (with the indices ... at the other levels depending on the values of the other mesh functions). Notes: 1. If j0 == 1, then the division point corresponding to j0-1 is -Infinity. 2. If the dimension of the mesh shading tensor at level j is less than 1 plus the number of division points, then the mesh shading tensor is cyclically extended. Jul 11, 2023 at 22:22

I think you are misunderstanding how MeshShading is supposed to work. The way it works is as follows: You give a set of MeshFunctions $$\{f_1, f_2, \dots\}$$, and a set of mesh divisions $$\{\{d_{1,1},d_{1,2},\dots\},\{d_{2,1},d_{2,2},\dots\},\dots\}$$ (effectively the setting of Mesh, i.e. the values of $$f_1,\dots$$ that where mesh contours should be drawn). The mesh divisions effectively define intervals as $$\{\{I_{1,1},I_{1,2},\dots\},\{I_{2,1},I_{2,2},\dots\},\dots\}=\{\{(-\infty,c_{1,1}],(c_{1,1},c_{1,2}],\dots\},\{(-\infty,c_{2,1}],(c_{2,1},c_{2,2}],\dots\},\dots\}.$$

This means that every point can be assigned to a set of intervals $$\{I_{1,i},I_{2,j},\dots\}$$, whose style is given by the entry $$m_{i,j,\dots}$$ from the MeshShading option.

To illustrate this, consider the following:

With[{plot =
ParametricPlot3D[eqNpar, {w, -1, 1}, {s, -1, 1},
MeshFunctions -> {Function[{x, y, z, w,
s}, (eqNcar /. Thread[xs -> eqNpar]) - eq1],
Function[{x, y, z, w, s}, (eqNcar /. Thread[xs -> eqNpar]) -
eq2], Function[{x, y, z, w,
s}, (eqNcar /. Thread[xs -> eqNpar]) - eq3]},
MeshStyle -> {{Blue}, {Red}, {Black}},
SparseArray[{{i_, j_, k_} :> style[i, j, k]}, {2, 2, 2},
Directive[]], Mesh -> {{0}}, PlotStyle -> Orange,
PlotPoints -> 200, BoundaryStyle -> {1 -> None}, AxesLabel -> xs]},
Manipulate[
plot /.
style[i_, j_, k_] :>
If[{i, j, k} == {i0, j0, k0}, Green, Directive[]],
{i0, 1, 2, 1},
{j0, 1, 2, 1},
{k0, 1, 2, 1}
]
]
`

As you can see, the first slider selects which "side" of the black mesh line you are on, the second which side of the red mesh line, and the third slider which side of the blue mesh line you are on.