Matching Conditions

The smoothness of the common element spaces in Table 2 on p66 depends on the smoothness of the matching between adjacent triangles:  after all the polynomials themselves are infinitely differentiable. We need to work out how to see how smooth the patched together surface is for the various element spaces.

Courant Triangles: (Subscript[ℳ, 0]^1)

The simplest space is the Courant triangles (Subscript[ℳ, 0]^1)of piecewise linear functions evaluated at the corners of the triangles.  It is clear that since the two linear functions agree at two points on the shared edge (the two corners) that the linear functions agree everywhere along the edge.  This means that the space is at least continuous.  It is not differentiable since the gradients are piecewise constant and can only match (at the border) if they are identical.

The three corners provide the three values neccessary to fix the three coefficients (Subscript[a, 0])+(Subscript[a, 1,1])x+(Subscript[a, 1,2])y.

Nodal Functions on reference triangle 0<X<1-Y and 0≤Y≤1: 
1- X - Y: matched to each vertex

Triangles={ {1,2,3},{2,3,4}};
points={ {0,0},{2,0},{1,1},{3,2}};
Show[
Graphics[{
PointSize[0.02],
Table[ {Hue[i/Length[Triangles]],Polygon[points[[Triangles[[i]]]]]},{i, Length[Triangles]}],
Table[Point[points[[i]]],{i,Length[points]}]
}
],
Frame->True
]

f1[x_,y_]=x
f2[x_,y_]=y
f3[x_,y_]=1-(x+y)
{Plot3D[
f1[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1],
PlotRange->{0,1}
],
Plot3D[
f2[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1],
PlotRange->{0,1}
],
Plot3D[
f3[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1],
PlotRange->{0,1}
]
}

(Subscript[ℳ, 0]^2)

Piecewise quadratic functions evaluated at the corners and midpoints of the triangles.

The three corners plus three midpoint values provide the six values neccessary to determine the six coefficients
 (Subscript[a, 0])+(Subscript[a, 1,1])x+(Subscript[a, 1,2]) y+(Subscript[a, 2,11])(x^2)+(Subscript[a, 2,12])x y+(Subscript[a, 2,22])(y^2)

It is clear that since the two quadratic functions agree at three points on the shared edge (the two corners and the midpoint) that the quadratic functions agree everywhere along the edge. This means that the space is at least continuous.

The space is not a subspace of (𝒞^1).  The gradients are piecewise linear 
({{(Subscript[a, 1,1])+2(Subscript[a, 2,11])x+(Subscript[a, 2,12])y}, {(Subscript[a, 1,2])+2(Subscript[a, 2,22])y+(Subscript[a, 2,12])x}})
and contain a total of 5 constants in each element.  The tangential derivative matches across the boundary:  each is the derivative of an identical quadratic along the boundary. However the normal derivative in each element depends on the non-shared points and do not match.

Nodal Functions on reference triangle 0<X<1-Y and 0≤Y≤1:

f[x_,y_]= (Subscript[a, 0])+(Subscript[a, 1,1])x+(Subscript[a, 1,2]) y+(Subscript[a, 2,11])(x^2)+(Subscript[a, 2,12])x y+(Subscript[a, 2,22])(y^2);
f1[x_,y_]=Factor[
f[x,y]/.Solve[ {f[0,0]==1, f[0,1]==0,f[1,0]==0,f[1/2,0]==0,f[0,1/2]==0,f[1/2,1/2]==0}]
]
f2[x_,y_]=Factor[
f[x,y]/.Solve[ {f[0,0]==0, f[0,1]==0,f[1,0]==0,f[1/2,0]==1,f[0,1/2]==0,f[1/2,1/2]==0}]
]
{Plot3D[
f1[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1],
PlotRange->{-0.5,1}
],
Plot3D[
f2[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1],
PlotRange->{-0.5,1}
]
}

Triangles={ {1,2,3},{2,3,4}};
points={ {0,1},{2,0},{1,1},{4,2}};
Show[
Graphics[{
PointSize[0.02],
Table[
Vertices=points[[Triangles[[i]]]];
 {
{Hue[i/Length[Triangles]],Polygon[Vertices]},
Point[(Vertices[[1]]+Vertices[[2]])/2],
Point[(Vertices[[2]]+Vertices[[3]])/2],
Point[(Vertices[[3]]+Vertices[[1]])/2]
},
{i, Length[Triangles]}],
Table[Point[points[[i]]],{i,Length[points]}],
MidPt=(points[[2]]+points[[3]])/2;
TanVec=(points[[2]]-points[[3]])/2;
NormVec={
{0,1},
{-1,0}
}.TanVec;
Thickness[0.01],
Arrow[{MidPt,MidPt+0.8TanVec}],
Arrow[{MidPt,MidPt+0.8NormVec}],
Arrow[{MidPt,MidPt-0.8NormVec}]
}
],
Frame->True
]

(Subscript[ℳ, 0]^3)

Piecewise cubic functions evaluated at the corners, thirdpoints of the triangles, and center of the triangles.

The three corners plus six midpoint values plus the center provide the ten values neccessary to determine the ten coefficients (Subscript[a, 0])+(Subscript[a, 1,1])x+(Subscript[a, 1,2]) y+(Subscript[a, 2,11])(x^2)+(Subscript[a, 2,12])x y+(Subscript[a, 2,22])(y^2)+(Subscript[a, 3,111])(x^3)+(Subscript[a, 3,112])((x )^2)y+(Subscript[a, 3,122])x (y^2)+(Subscript[a, 3,222])(x^3)

It is clear that since the two cubice functions agree at four points on the shared edge (the two corners and the two third point) that the cubic functions agree everywhere along the edge. This means that the space is at least continuous.

The space is not a subspace of (𝒞^1).  The gradients are piecewise quadratic.  The tangential derivative matches across the boundary:  each is the derivative of an identical cubic along the boundary. However the normal derivatives depend on the non-shared points and do not match.

Triangles={ {1,2,3},{2,3,4}};
points={ {0,1},{2,0},{1,1},{4,2}};
Show[
Graphics[{
PointSize[0.02],
Table[
Vertices=points[[Triangles[[i]]]];
 {
{Hue[i/Length[Triangles]],Polygon[Vertices]},
Point[(Vertices[[1]]+2Vertices[[2]])/3],
Point[(Vertices[[2]]+2Vertices[[3]])/3],
Point[(Vertices[[3]]+2Vertices[[1]])/3],
Point[(2Vertices[[1]]+Vertices[[2]])/3],
Point[(2Vertices[[2]]+Vertices[[3]])/3],
Point[(2Vertices[[3]]+Vertices[[1]])/3],
Point[(Vertices[[1]]+Vertices[[2]]+Vertices[[3]])/3]
},
{i, Length[Triangles]}],
Table[Point[points[[i]]],{i,Length[points]}],
MidPt=(points[[2]]+points[[3]])/2;
TanVec=(points[[2]]-points[[3]])/2;
NormVec={
{0,1},
{-1,0}
}.TanVec;
Thickness[0.01],
Arrow[{MidPt,MidPt+0.8TanVec}],
Arrow[{MidPt,MidPt+0.8NormVec}],
Arrow[{MidPt,MidPt-0.8NormVec}]
}
],
Frame->True
]

Clear[f,x,y]
f[x_,y_]= (∑_(i=0)^3)((∑_(j=0)^i)(Subscript[a, i,j,i-j])(x^(i-j))(y^j))
f1[x_,y_]=Factor[
f[x,y]/.Solve[ {
f[0,0]==1, 
f[1/3,0]==0,
f[2/3,0]==0,
f[1,0]==0,
f[2/3,1/3]==0,
f[1/3,2/3]==0,
f[0,1]==0,
f[0,2/3]==0,
f[0,1/3]==0,
f[1/3,1/3]==0
}]
]

f2[x_,y_]=Factor[
f[x,y]/.Solve[ {
f[0,0]==0, 
f[1/3,0]==1,
f[2/3,0]==0,
f[1,0]==0,
f[2/3,1/3]==0,
f[1/3,2/3]==0,
f[0,1]==0,
f[0,2/3]==0,
f[0,1/3]==0,
f[1/3,1/3]==0
}]
]
f3[x_,y_]=Factor[
f[x,y]/.Solve[ {
f[0,0]==0, 
f[1/3,0]==0,
f[2/3,0]==0,
f[1,0]==0,
f[2/3,1/3]==0,
f[1/3,2/3]==0,
f[0,1]==0,
f[0,2/3]==0,
f[0,1/3]==0,
f[1/3,1/3]==1
}]
]
{Plot3D[
f1[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1],
PlotRange->{-0.5,1}
],
Plot3D[
f2[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1],
PlotRange->{-0.5,1}
],
Plot3D[
f3[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1],
PlotRange->{-0.5,1}
]
}

Argyris Triangles

The first three spaces should convince you that the problem is matching the normal derivatives.

The three corners provide the three values neccessary to fix the three coefficients (Subscript[a, 0])+(Subscript[a, 1,1])x+(Subscript[a, 1,2])y.

Clear[f,x,y]
f[x_,y_]= (∑_(i=0)^5)((∑_(j=0)^i)(Subscript[a, i,j,i-j])(x^(i-j))(y^j))
f1[x_,y_]=Factor[
f[x,y]/.Solve[ {
(* (0,0) equations *)
f[0,0]==1, 
Derivative[1,0][f][0,0]==0,
Derivative[0,1][f][0,0]==0,
Derivative[2,0][f][0,0]==0,
Derivative[1,1][f][0,0]==0,
Derivative[0,2][f][0,0]==0,
(* (0,1) equations *)
f[0,1]==0, 
Derivative[1,0][f][0,1]==0,
Derivative[0,1][f][0,1]==0,
Derivative[2,0][f][0,1]==0,
Derivative[1,1][f][0,1]==0,
Derivative[0,2][f][0,1]==0,
(* (1,0) equations *)
f[1,0]==0, 
Derivative[1,0][f][1,0]==0,
Derivative[0,1][f][1,0]==0,
Derivative[2,0][f][1,0]==0,
Derivative[1,1][f][1,0]==0,
Derivative[0,2][f][1,0]==0,
(* MidPoint Normal Equations *)
1/(Sqrt[2]) (Derivative[1,0][f][1/2,1/2]+ Derivative[0,1][f][1/2,1/2]) ==0,
-Derivative[1,0][f][0,1/2]==0,
-Derivative[0,1][f][1/2,0]==0
}]
]
f2[x_,y_]=Factor[
f[x,y]/.Solve[ {
(* (0,0) equations *)
f[0,0]==0, 
Derivative[1,0][f][0,0]==0,
Derivative[0,1][f][0,0]==0,
Derivative[2,0][f][0,0]==0,
Derivative[1,1][f][0,0]==1,
Derivative[0,2][f][0,0]==0,
(* (0,1) equations *)
f[0,1]==0, 
Derivative[1,0][f][0,1]==0,
Derivative[0,1][f][0,1]==0,
Derivative[2,0][f][0,1]==0,
Derivative[1,1][f][0,1]==0,
Derivative[0,2][f][0,1]==0,
(* (1,0) equations *)
f[1,0]==0, 
Derivative[1,0][f][1,0]==0,
Derivative[0,1][f][1,0]==0,
Derivative[2,0][f][1,0]==0,
Derivative[1,1][f][1,0]==0,
Derivative[0,2][f][1,0]==0,
(* MidPoint Normal Equations *)
1/(Sqrt[2]) (Derivative[1,0][f][1/2,1/2]+ Derivative[0,1][f][1/2,1/2]) ==0,
-Derivative[1,0][f][0,1/2]==0,
-Derivative[0,1][f][1/2,0]==0
}]
]
f3[x_,y_]=Factor[
f[x,y]/.Solve[ {
(* (0,0) equations *)
f[0,0]==0, 
Derivative[1,0][f][0,0]==0,
Derivative[0,1][f][0,0]==0,
Derivative[2,0][f][0,0]==0,
Derivative[1,1][f][0,0]==0,
Derivative[0,2][f][0,0]==0,
(* (0,1) equations *)
f[0,1]==0, 
Derivative[1,0][f][0,1]==0,
Derivative[0,1][f][0,1]==0,
Derivative[2,0][f][0,1]==0,
Derivative[1,1][f][0,1]==0,
Derivative[0,2][f][0,1]==0,
(* (1,0) equations *)
f[1,0]==0, 
Derivative[1,0][f][1,0]==0,
Derivative[0,1][f][1,0]==0,
Derivative[2,0][f][1,0]==0,
Derivative[1,1][f][1,0]==0,
Derivative[0,2][f][1,0]==0,
(* MidPoint Normal Equations *)
1/(Sqrt[2]) (Derivative[1,0][f][1/2,1/2]+ Derivative[0,1][f][1/2,1/2]) ==1,
-Derivative[1,0][f][0,1/2]==0,
-Derivative[0,1][f][1/2,0]==0
}]
]

{Plot3D[
f1[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1],
PlotRange->{-0.5,1}
],
Plot3D[
f2[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1]
],
Plot3D[
f3[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1]
]
}

Bell Triangles

There are more degrees fo freedom (and complicated ones at that) in the Argyris triangle than are required. The Bell triangle requires creating polynomials (fifth order) that have lower order polynomials for their normal derivatives on the edges.

Clear[f,x,y]
f[x_,y_]= (∑_(i=0)^5)((∑_(j=0)^i)(Subscript[a, i,j,i-j])(x^(i-j))(y^j))

On the first edge (y=0 and 0≤x≤1) the normal derivative is

Derivative[0,1][f][x,0]

This is a cubic iff (Subscript[a, 5,1,4])=0.  Along the second edge (0≤x≤1 and y=1-x) the normal derivative is complicated but the coefficient of (x^4) is  computed below.

1/(Sqrt[2])(Derivative[1,0][f][x,1-x]+Derivative[0,1][f][x,1-x])
Coefficient[1/(Sqrt[2])(Derivative[1,0][f][x,1-x]+Derivative[0,1][f][x,1-x]),(x^4)]

The normal derivative is a polynomial on the edge iff this is zero. On the last edge (x=0 and 0≤y≤1) the normal derivative is

-Derivative[1,0][f][0,y]

This is a cubic iff (Subscript[a, 5,4,1])=0.  Combining these we get the conditions that the quintic is a cubic along the three edges iff (Subscript[a, 5,1,4])=(Subscript[a, (Subscript[a, 5,4,1])])=0  and 5 (Subscript[a, 5,0,5])+(Subscript[a, 5,2,3])+(Subscript[a, 5,3,2])+5 (Subscript[a, 5,5,0])=0.We can use this condition to plot the basis functions for the Bell triangular elements.

f1[x_,y_]=Factor[
f[x,y]/.Solve[ {
(* Cubic Conditions *)
(Subscript[a, 5,1,4])==0,
(Subscript[a, 5,4,1])==0,
5 (Subscript[a, 5,0,5])+(Subscript[a, 5,2,3])+(Subscript[a, 5,3,2])+5 (Subscript[a, 5,5,0])==0,
(* (0,0) equations *)
f[0,0]==0, 
Derivative[1,0][f][0,0]==0,
Derivative[0,1][f][0,0]==0,
Derivative[2,0][f][0,0]==1,
Derivative[1,1][f][0,0]==0,
Derivative[0,2][f][0,0]==0,
(* (0,1) equations *)
f[0,1]==0, 
Derivative[1,0][f][0,1]==0,
Derivative[0,1][f][0,1]==0,
Derivative[2,0][f][0,1]==0,
Derivative[1,1][f][0,1]==0,
Derivative[0,2][f][0,1]==0,
(* (1,0) equations *)
f[1,0]==0, 
Derivative[1,0][f][1,0]==0,
Derivative[0,1][f][1,0]==0,
Derivative[2,0][f][1,0]==0,
Derivative[1,1][f][1,0]==0,
Derivative[0,2][f][1,0]==0
}]
]
f2[x_,y_]=Factor[
f[x,y]/.Solve[ {
(* Cubic Conditions *)
(Subscript[a, 5,1,4])==0,
(Subscript[a, 5,4,1])==0,
5 (Subscript[a, 5,0,5])+(Subscript[a, 5,2,3])+(Subscript[a, 5,3,2])+5 (Subscript[a, 5,5,0])==0,
(* (0,0) equations *)
f[0,0]==0, 
Derivative[1,0][f][0,0]==0,
Derivative[0,1][f][0,0]==0,
Derivative[2,0][f][0,0]==0,
Derivative[1,1][f][0,0]==1,
Derivative[0,2][f][0,0]==0,
(* (0,1) equations *)
f[0,1]==0, 
Derivative[1,0][f][0,1]==0,
Derivative[0,1][f][0,1]==0,
Derivative[2,0][f][0,1]==0,
Derivative[1,1][f][0,1]==0,
Derivative[0,2][f][0,1]==0,
(* (1,0) equations *)
f[1,0]==0, 
Derivative[1,0][f][1,0]==0,
Derivative[0,1][f][1,0]==0,
Derivative[2,0][f][1,0]==0,
Derivative[1,1][f][1,0]==0,
Derivative[0,2][f][1,0]==0
}]
]
f3[x_,y_]=Factor[
f[x,y]/.Solve[ {
(* Cubic Conditions *)
(Subscript[a, 5,1,4])==0,
(Subscript[a, 5,4,1])==0,
5 (Subscript[a, 5,0,5])+(Subscript[a, 5,2,3])+(Subscript[a, 5,3,2])+5 (Subscript[a, 5,5,0])==0,
(* (0,0) equations *)
f[0,0]==1, 
Derivative[1,0][f][0,0]==0,
Derivative[0,1][f][0,0]==0,
Derivative[2,0][f][0,0]==0,
Derivative[1,1][f][0,0]==0,
Derivative[0,2][f][0,0]==0,
(* (0,1) equations *)
f[0,1]==0, 
Derivative[1,0][f][0,1]==0,
Derivative[0,1][f][0,1]==0,
Derivative[2,0][f][0,1]==0,
Derivative[1,1][f][0,1]==0,
Derivative[0,2][f][0,1]==0,
(* (1,0) equations *)
f[1,0]==0, 
Derivative[1,0][f][1,0]==0,
Derivative[0,1][f][1,0]==0,
Derivative[2,0][f][1,0]==0,
Derivative[1,1][f][1,0]==0,
Derivative[0,2][f][1,0]==0
}]
]

{Plot3D[
f1[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1]
],
Plot3D[
f2[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1]
],
Plot3D[
f3[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1]
]
}

Macro Elements

The Hsieh-Clough-Tougher macro element will wind up looking a bit different.  Lets pick an interior point and call it {(Subscript[x, 0]),(Subscript[y, 0])} The we are looking for cubics on the three sub triangles that match along the interior sub-domain edges.

{x0,y0}={1/3,1/3};
Triangles={ {1,2,4},{2,3,4},{1,3,4}};
points={ {0,0},{1,0},{0,1},{x0,y0}};
Show[
Graphics[{
PointSize[0.02],
Table[ {Hue[i/Length[Triangles]],Polygon[points[[Triangles[[i]]]]]},{i, Length[Triangles]}],
Table[Point[points[[i]]],{i,Length[points]}]
}
],
Frame->True
]

Clear[f,f1,f2,f3,x,y,a,b,c]
f1[x_,y_]= (∑_(i=0)^3)((∑_(j=0)^i)(Subscript[a, i,j,i-j])(x^(i-j))(y^j));
f2[x_,y_]= (∑_(i=0)^3)((∑_(j=0)^i)(Subscript[b, i,j,i-j])(x^(i-j))(y^j));
f3[x_,y_]= (∑_(i=0)^3)((∑_(j=0)^i)(Subscript[c, i,j,i-j])(x^(i-j))(y^j));

f[x_,y_]:=Piecewise[
{f1[x,y],And[0≤y≤(Subscript[y, 0]), y ((Subscript[x, 0])/(Subscript[y, 0]))≤x≤]},
{f2[x,y], And[0≤x≤(Subscript[x, 0]),((Subscript[y, 0])/(Subscript[x, 0]))≤y≤1-x()]},
{f3[x,y],Cond3},
{0, True}
];

giving a total of 3×10 = 30 coefficients.  We specify three conditions (function values and first derivatives) at the outside corners:  this represents 6 conditions (each sub triangle has two exterior corners) on each of the cubics for a total of 18 conditions. We add conditions (f1==f2=f3, ∇f1=∇f2=∇f3; and (∇^2)f1 = (∇^2)f2 = (∇^2)f3)  at the internal nodes  for a total of (2 function conditions; 2×2=4 first derivative conditions; 2×3=6 second derivative conditions) to resolve the remaining 12 degrees of freedom.  There are three duplicate equations in the matching conditions and we can enforce the

Note: All the cubics match their values and normal derivatives along the internal edges.  The same arguments we have been making all along ensure that the big triangles provide a (C^1) space.

Clear[f1s,f2s,f3s,x,y]
{f1s[x_,y_],f2s[x_,y_],f3s[x_,y_]}={f1[x,y],f2[x,y],f3[x,y]}/.Solve[ {
(* (0,0) equations for f1 *)
f1[0,0]==1, 
Derivative[1,0][f1][0,0]==0,
Derivative[0,1][f1][0,0]==0,
(* (0,0) equations for f2 *)
f2[0,0]==f1[0,0], 
Derivative[1,0][f2][0,0]==Derivative[1,0][f1][0,0],
Derivative[0,1][f2][0,0]==Derivative[0,1][f1][0,0],
(* (0,1) equations for f2 *)
f2[0,1]==0, 
Derivative[1,0][f2][0,1]==0,
Derivative[0,1][f2][0,1]==0,
(* (0,1) equations for f3 *)
f3[0,1]==f2[0,1], 
Derivative[1,0][f3][0,1]==Derivative[1,0][f2][0,1],
Derivative[0,1][f3][0,1]==Derivative[0,1][f2][0,1],
(* (1,0) equations for f3 *)
f3[1,0]==0, 
Derivative[1,0][f3][1,0]==0,
Derivative[0,1][f3][1,0]==0,
(* (1,0) equations for f1 *)
f1[1,0]==f3[1,0], 
Derivative[1,0][f1][1,0]==Derivative[1,0][f3][1,0],
Derivative[0,1][f1][1,0]==Derivative[0,1][f3][1,0],
(* Matching Conditions Function *)
f1[x0,y0]==f2[x0,y0],
f2[x0,y0]==f3[x0,y0],
(* Matching Conditions Gradient *)
Derivative[1,0][f1][x0,y0]==Derivative[1,0][f2][x0,y0],
Derivative[0,1][f1][x0,y0]==Derivative[0,1][f2][x0,y0],
Derivative[1,0][f2][x0,y0]==Derivative[1,0][f3][x0,y0],
Derivative[0,1][f2][x0,y0]==Derivative[0,1][f3][x0,y0],
(* Matching Conditions Second Gradient *)
Derivative[2,0][f1][x0,y0]==Derivative[2,0][f2][x0,y0],
Derivative[1,1][f1][x0,y0]==Derivative[1,1][f2][x0,y0],
Derivative[0,2][f1][x0,y0]==Derivative[0,2][f2][x0,y0],
Derivative[2,0][f3][x0,y0]==Derivative[2,0][f2][x0,y0],
Derivative[1,1][f3][x0,y0]==Derivative[1,1][f2][x0,y0],
Derivative[0,2][f3][x0,y0]==Derivative[0,2][f2][x0,y0],
(* Normal Derivative Conditions *)
-Derivative[0,1][f1][1/2,0]==0,
-Derivative[1,0][f2][0,1/2]==0,
1/(Sqrt[2]) (Derivative[1,0][f1][1/2,1/2]+Derivative[0,1][f1][1/2,1/2])==0
}
][[1]]

fs[x_,y_]:=Piecewise[
{
{f1s[x,y],And[0≤y≤y0, y (x0/y0)≤x≤x0+((1-x0)/y0)(y0-y)]},
{f2s[x,y], And[0≤x≤x0,x(y0/x0)≤y≤1+((y0-1)/x0)x]},
{f3s[x,y],And[0≤x≤1,0≤y≤1-x]},
{0, True}
}
];

fs[0.2,0.3]

Plot3D[
fs[x,y],{x,0,1},{y,0,1},
RegionFunction->Function[{x,y},x+y<1]
]