5
$\begingroup$

I have a deformed 2D region of triangular Elements specified by their Nodes in $x$ and $y$ coordinate system. In addition, I have a list of Displacements in $y$ direction for each node. I want to ContourPlot the displacements over my deformed region of elements. So far my reasoning was to make a list Nodes2 where each sublist has in the first two places the node position and the third is displacement. Then I use ListContourPlot.

Nodes = {{0., 0.}, {10., 0.}, {20., 0.}, {30., 0.}, {-0.101573, 9.81671}, {9.92361, 10.0459}, {19.8914, 10.2668}, {29.7869, 10.5013}, {-0.453558, 19.6276}, {9.57795, 20.086}, {19.5471, 20.5443}, {29.4276, 21.0344}, {-1.06222, 29.4477}, {8.97881, 30.1177}, {18.9679, 30.8116}, {28.8439, 31.5901}, {-1.90436,  39.3141}, {8.12513, 40.1282}, {18.1438, 41.0447}, {28.0142, 42.1991}, {-2.88222, 49.2855}, {7.0891, 50.1568}, {16.9537, 51.1588}, {26.904, 52.9408}};

Elements = {{1, 2, 6}, {6, 5, 1}, {2, 3, 7}, {7, 6, 2}, {3, 4, 8}, {8, 7, 3}, {5, 6, 10}, {10, 9, 5}, {6, 7, 11}, {11, 10, 6}, {7, 8, 12}, {12, 11, 7}, {9, 10, 14}, {14, 13, 9}, {10, 11, 15}, {15, 14, 10}, {11, 12, 16}, {16, 15, 11}, {13, 14, 18}, {18, 17, 13}, {14, 15, 19}, {19, 18, 14}, {15, 16, 20}, {20, 19, 15}, {17, 18, 22}, {22, 21, 17}, {18, 19, 23}, {23, 22, 18}, {19, 20, 24}, {24, 23, 19}};

displacements = {0, 0, 0, 0, -0.0091645, 0.00229398, 0.013342, 0.0250656, -0.0186183, 0.0043011, 0.0272138, 0.0517176, -0.0276138, 0.00588411, 0.0405803, 0.0795059, -0.0342935, 0.00641231, 0.0522328, 0.109954, -0.0357271, 0.00784115, 0.0579424, 0.14704};

stn = Nodes//Length;
Nodes2 = Nodes;
For[i = 1, i <= stn, i++,
 AppendTo[Nodes2[[i]], displacements[[i]]]
 ]
ListContourPlot[Nodes2, Contours -> 20, ColorFunction -> "Rainbow", 
 AspectRatio -> 5/3, PlotRange -> All, PlotRange -> Full]

The result is as shown in the picture. The left one represents the deformed shape made of triangles and the right one my ListContourPlot with highlighted spots where there shouldn't be any color as the region does not exist (refer to the left picture). I have found a similar problem here, but I don't see how it can be generalized for non DelaunayMesh problems (at least with my knowledge). My question is: How to plot DISCRETE values (in this case displacements) over ONLY a SPECIFIED DISCRETE region?

enter image description here

$\endgroup$
8
  • $\begingroup$ I get a completely different image from your code and many error messages about invalid access to list elements. Can you revise your example? $\endgroup$
    – halirutan
    Commented Jun 10, 2018 at 1:34
  • $\begingroup$ My bad, I had a vector of displacements both in $x$ and $y$ but for the sake of simplicity only included $y$ here. Should be fine now. $\endgroup$
    – KeVal
    Commented Jun 10, 2018 at 1:42
  • $\begingroup$ Nice! Now let the answers roll and tell you that you should never ever use AppendTo in a For loop :) $\endgroup$
    – halirutan
    Commented Jun 10, 2018 at 1:44
  • $\begingroup$ :) why is that? $\endgroup$
    – KeVal
    Commented Jun 10, 2018 at 1:46
  • $\begingroup$ Very ineffective on large datasets since it copies and reallocates a lot of memory. There are several easier ways to manipulate lists in a more Mathematica-like way: MapThread[Append, {Nodes, displacements}] or Flatten[{Nodes, List /@ displacements}, {{2}, {1, 3}}]. $\endgroup$
    – halirutan
    Commented Jun 10, 2018 at 1:54

4 Answers 4

4
$\begingroup$

Based on your cells, you can easily create MeshRegion. From this region, you can create a function that tests if a point is inside the region.

What I'm doing is creating your ListDensityPlot and in an after step remove all polygons and lines that are not inside your region. I'm not sure how robust this is, but it seems to work

ldp = ListContourPlot[MapThread[Append, {Nodes, displacements}], 
   Contours -> 20, ColorFunction -> "Rainbow", AspectRatio -> 5/3, 
   PlotRange -> All, PlotRange -> Full];
insideQ = RegionMember[MeshRegion[Nodes, Triangle /@ Elements]];

DeleteCases[Normal[ldp], (Line | Polygon)[pts_, ___] /; 
  Not[And @@ (insideQ /@ pts)], Infinity]

Mathematica graphics

For reference, please see my answer here. As a side note: the visible white borders have nothing to do with my approach. It is an artifact that is introduced by using Normal (which turns a GraphicsComplex into a list of graphics primitives). We have several questions and solutions about this here.

You can get rid of them by explicitly including a color for the edges:

DeleteCases[
  Normal[ldp], (Line | Polygon)[pts_, ___] /; 
   Not[And @@ (insideQ /@ pts)], 
  Infinity] /. {EdgeForm[], 
   r_?(MemberQ[{RGBColor, Hue, CMYKColor, GrayLevel}, Head[#]] &), 
   i___} :> {EdgeForm[r], r, i}

Mathematica graphics

See this post for reference.

$\endgroup$
3
$\begingroup$

You can create an ElementMesh and then an InterpolatingFunction which you can use for a plot:

Needs["NDSolve`FEM`"]
m = ToElementMesh["Coordinates" -> Nodes, 
   "MeshElements" -> {TriangleElement[Elements]}];
eif = ElementMeshInterpolation[{m}, displacements];


ContourPlot[eif[x, y], Element[{x, y}, m], Contours -> 20, 
 ColorFunction -> "Rainbow", AspectRatio -> 5/3, PlotRange -> All, 
 PlotRange -> Full]

enter image description here

You could also use the ElementMeshDeformation function to create the deformation:

ElementMeshDeformation[m, {eif, eif}, 
  "ScalingFactor" -> 1]["Wireframe"]

enter image description here

In this case I used the same deformation in x and y direction.

From version 13.0 onward there is also VectorDisplacementPLot:

VectorDensityPlot[{eif[x, y], eif[x, y]}, {x, y} \[Element] m]

enter image description here

$\endgroup$
2
$\begingroup$

Here is a completely different approach that circumvents ListDensityPlot. I assume you have properly defined your Nodes, Elements, and displacements.

Since you already have nicely constructed cells defined by your "point references" in Elements, there is absolutely no reason not to use them. What if we draw a polygon for each cell and color the vertices according to your displacement? The good thing is that inside each polygon, colors are linearly interpolated. That doesn't look as good as in ListContourPlot, but the code is very easy to understand so let's give this a try.

A GraphicsComplex[pts, graphx] is a way to define points pts only once and refer to them by using their position. How handy it is, that you already have exactly this structure. What we need is to turn each displacement into a color. Here is a small function that takes an element index and returns a rainbow color:

With[{dis = Rescale[displacements]},
 toColor[n_Integer] := List @@ ColorData["Rainbow", dis[[n]]];
 SetAttributes[toColor, {Listable}]
]

We need Listable because we want to be able to call toColor[{{1,2,3},...] and get {{toColor[1], toColor[2], toColor[3]}, ...}. That's all. Now look at this

Graphics@GraphicsComplex[Nodes, 
  Polygon[Elements, VertexColors -> toColor[Elements]]
]

Mathematica graphics

Adding lines to the cells is easy too

Graphics@GraphicsComplex[Nodes,
  {
   Polygon[Elements, VertexColors -> toColor[Elements]],
   GrayLevel[.2], Line[Elements]
   }
  ]

Mathematica graphics

And maybe a bit of antialiasing if this looks as pixelized on your machine

Graphics@GraphicsComplex[Nodes,
  {
   Polygon[Elements, VertexColors -> toColor[Elements]],
   GrayLevel[.2], Opacity[.2], Thick, Line[Elements],
   Thin, Opacity[1], Line[Elements]
   }
  ]

Mathematica graphics

$\endgroup$
0
$\begingroup$

This is a bit of a hack.

I generate a mesh from the points you have, then subtract it from a rectangle sized as the plot range of ListContourPlot to leave a transparent hole and a white mask, then superimposing this mask on the contour plot using Epilog:

points = MapThread[Join[#1, {#2}] &, {Nodes, displacements}];

mask =
  MeshPrimitives[#, 2] &@
   DiscretizeRegion@
    RegionDifference[
     Rectangle[{-3, -0.5}, {30, 53}],
     MeshRegion[Nodes, Triangle /@ Elements]
    ];

ListContourPlot[
  points, Contours -> 20, ColorFunction -> "Rainbow",
  AspectRatio -> Automatic,
  Epilog -> {FaceForm[White], EdgeForm[None], mask}
]

Mathematica graphics

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.