3D color map of F=F(x,y,z) where data is given as a 2D array of the form [x y z F]

The output of a CFD calculation is usually given in the form of a 2D array [x y z F] where F is a function such as pressure or velocity that is calculated for the given points xyz in the 3D space. The result is then given as a color map, as shown below as an example for a relatively simple channel geomtry. How can this be done in MATLAB giventhe [x y z F] output format?

[Saurabh Chaudhary]

you can try to use contour function.

[Torsten]

Use "slice" for the three visible planes in one plot.

Answer 1

[Karim]

Hello , you can create those plots both for both 2D and 3D elements using the patch command. Essentially you need a connecitivty matrix (i.e. the faces of the patch) which you have since these are the faces of the elements used in the mesh. Hense next to the coordindates of the nodes, you also need the connecivity. If you can share a set of sample data (i.e. the mesh, the nodes and the value for plotting) i can adjust the example to meet you data. You can use the paperclip symbol to attach data (of a small example) to you question.

See below for a demonstration of the concept:

% the data here will differ depending on the used FE package, these are
% outputs from ansys i used as a demo in another matlab answer

% this file contains the xyz coordinates of the grid points
coordinates = readmatrix("coordinates.csv");

% this file contains the deformation values, which we will use to create the 'F' variable
displacements = readmatrix("displacements.csv"); 
% let's use the magnitude of the deformation as color
F = sqrt( sum(displacements(:,2:end).^2,2));

% this file containts the node conectivity (i.e. the nodes that make up the elements)
nodes = readmatrix("nodes.csv"); 


% now create the plot
figure
patch('Faces',nodes(:,2:end),'Vertices',coordinates(:,2:end),'FaceVertexCData',F,'FaceColor','interp',...
    'EdgeAlpha',0.15) % put the edge alpha to 0 if you went to hide the element edges
view([-15 45])
axis equal
colormap jet
colorbar
title("Deformed")

[Karim]

Hi @Saeid, the file you provided is not enough. The mesh itself is also needed. The concectivity between the nodes is required for creating the plot.

[Saeid]

Hi Karim,

I have attached the exported mesh to this message. To make sure that the mesh file corresponds with the solution, I have added the solution file once again, so ignore the one I sent to you in my previous message.

[Torsten]

@Saeid

If you only want a graphical presentation of the results, why don't you use those pictures supplied by the CFD program ?

[Saeid]

Correct, theoretically that is possible, but there are two reasons for the need for a graphical presentation outside of the CFD software:

1. I am writing a post-processing app using MATLAB that reads the data and saves a reduced or edited version into an output file (reasons are given below). It will be much easier to have a preview of the existing solution before getting to edit the results. In addition, this App may be used by other users who do not have a CFD program to view the results and the only possible way for them to have a closer look at the results is to plot the solution in another medium.
2. In order to plot the reduced (.i.e. fewer data points, mentioned above) in a proper graphical program such as Origin, I first need to have a view of the 3D solution in the App to see if everything is OK. You may argue that the CFD program offers that, but MATLAB gives you the opportunity to have a close look and edit parts of the database. The reason why I need a 3rd program such as Origin for presentation is that for scientific papers I usually see people doing what you just said (directly use the graphics from the CFD program) but in many publications these images turn out to have a poor quality and I'd rather have them at the same, unified format in which I paste my other plots into the document. This may sound like a aesthetics issue but as I said, MATLAB works as a very good intermediate here.

[Karim]

@Saeid, the reults like you have them now would require a lot of insight into the problem and a lot of work to map. Notice on the second figure the red points, these corespond with the result grid in your text file. But these points do not match the grid points of your mesh. You even have less results then you have elements, if you would have a result at the center of the element you could use the shapefuncions of the element to determine the results at the nodes. However, the data like you presented it will require some special extrapolation. I typicaly work with Nastran and Samcef, these provide the results either at directly at the nodes or at the element centers. I'm a bit surprised that comsol gives so few result points spread over the mesh.

Below you can get an idea on how to plot the mesh and the points, if you want to plot the mesh with colors for some values, like i showed in my answer, you will need to somehow extrapolate your results to values at the nodes.

Mesh = readmatrix('3D_LPTT_ThinSlitEntrance_CoarseMesh.txt');
idx = find(all(isnan(Mesh),2));
Grid = Mesh((1+idx(end-1)):(idx(end)-1),1:3);
Elem = Mesh((1+idx(end)):end,:);

Results = readmatrix('3D_LPTT_ThinSlit_Coarse_P.txt');
ResultGrid = Results(:,1:3);

% plot the mesh
figure
patch('Faces',Elem,'Vertices',Grid,'FaceColor','r','EdgeColor','k')
grid on
view([-210 40])
axis equal

% plot the grid of the results on top of the edges of the mesh
figure
hold on
patch('Faces',Elem,'Vertices',Grid,'FaceColor','none','EdgeColor','k')
scatter3(ResultGrid(:,1),ResultGrid(:,2),ResultGrid(:,3),'r','filled')
hold off
grid on
view([-210 40])
axis equal

[Saeid]

Hi Karim,

this still looks very promising, thank you for the help!

I may still have sent you two files for mesh and results that correspond to slightly different geometries. Let me check this in my files and, before bothering you again, I will also try to run your routine on a simpler geometry (channel with just a round or slit profile and NOT both of them). Also keep in mind that this is a CFD problem where a viscoelastic (LPTT) fluid enters form the round end of this channel and exits from the narrow end. The BCs are no fluid movement at the outside borders + a certain flow rate (or fluid velocity) at the entrance and the solution is given in the form of fluid velocity at every node, which will then be converted to pressure values at each node using the said viscoelastic constitutive equaion. I don't know whether this then leads to the differences you have observed, but thought it would be helpful to mention.

As to the mesh, I have just used a coarse mesh for speed and to reduce the size of the files I am sending here, and that is probably the reason why the number of nodes is so low, but I can assure you that COMSOL can go up with node numbers.

One final questions about the colormap of the solution: is it not possible because of the discrepancy you mentioned above, or is it generally not posible?