构造三元网格，在Matlab中评估网格和轮廓图上的函数 [英] Construct ternary grid, evaluate a function on the grid and contour plot in Matlab

问题描述

我需要评估一个函数(例如) Fxy = 2 * x.^ 2 +3 * y.^ 2; 在三元网格x范围(0-1)，y范围(0-1)和1-x-y(0-1)上. 我无法构造需要评估上述功能的三元网格.同样，一旦评估，我就需要在三元轮廓图中绘制函数.理想情况下，我需要使轴在(x-> y->(1-x-y))的意义上逆时针旋转.

I need to evaluate a function (say) Fxy = 2*x.^2 +3 *y.^2; on a ternary grid x-range (0 - 1), y-range (0-1) and 1-x-y (0 - 1). I am unable to construct the ternary grid on which I need to evaluate the above function. Also, once evaluated I need to plot the function in a ternary contour plot. Ideally, I need the axes to go counter clockwise in the sense (x -> y--> (1-x-y)).

我已经尝试过该功能

function tg = triangle_grid ( n, t )

  ng = ( ( n + 1 ) * ( n + 2 ) ) / 2;
  tg = zeros ( 2, ng );

  p = 0;

  for i = 0 : n
    for j = 0 : n - i
      k = n - i - j;
      p = p + 1;
      tg(1:2,p) = ( i * t(1:2,1) + j * t(1:2,2) + k * t(1:2,3) ) / n;
    end
  end

  return
end

获取三角形边缘坐标之间子间隔的数量

for the number of sub intervals between the triangle edge coordinates

n = 10 (say)

以及等边三角形的边缘坐标

and for the edge coordinates of an equilateral triangle

t = tcoord = [0.0, 0.5,           1.0;
              0.0, 1.0*sqrt(3)/2, 0.0];

这生成了一个x轴为0-1的三角形网格，但其他两个都不是0-1.

This generated a triangular grid with the x-axis from 0-1 but the other two are not from 0-1.

我需要这样的东西:

...，轴范围为0-1( 0-100也可以).

... with the axes range 0-1 (0-100 would also do).

此外，我需要知道三角形网格内所有交叉点的坐标点.一旦有了这个，我就可以继续评估该网格中的功能了.

In addition, I need to know the coordinate points for all intersections within the triangular grid. Once I have this I can proceed to evaluate the function in this grid.

我的最终目标是得到这样的东西.与我现在删除的上一个情节相比，这是我需要实现的更好的表示.

My final aim is to get something like this. This is a better representation of what I need to achieve (as compared to the previous plot which I have now removed)

请注意，两个三元图的等值线轮廓大小不同.在我的情况下，差异是一个数量级，两个非常不同的Fxy.

Note that the two ternary plots have iso-value contours which are different in in magnitude. In my case the difference is an order of magnitude, two very different Fxy's.

如果我可以将两个三元图相互绘制在一起，然后在三元平面上两个等值轮廓的交点处评估构图.成分应从三元图读取，而不是从其上定义了三角形的矩形网格中读取. 当前存在问题(如评论部分中突出显示，一旦问题更接近解决方案，它将进行更新).

If I can plot the two ternary plots on top of each other then and evaluate the compositions at the intersection of two iso-value contours on the ternary plane. The compositions should be as read from the ternary plot and not the rectangular grid on which triangle is defined. Currently there are issues (as highlighted in the comments section, will update this once the problem is closer to solution).

推荐答案

我是 ternplot .正确推测后，ternpcolor不会做您想要的事情，因为它是自动构建的，可对数据进行网格化.回想起来，这不是一个特别明智的决定，我已经作了记录更改设计.同时，此代码应执行您想要的操作:

I am the author of ternplot. As you have correctly surmised, ternpcolor does not do what you want, as it is built to grid data automatically. In retrospect, this was not a particularly wise decision, I've made a note to change the design. In the mean time this code should do what you want:

我已经更改了代码以查找两条曲线的交点，而不仅仅是一条曲线的交点.

I've changed the code to find the intersection of two curves rather than just one.

N = 10;
x = linspace(0, 1, N);
y = x;

% The grid intersections on your diagram are actually rectangularly arranged,
% so meshgrid will build the intersections for us
[xx, yy] = meshgrid(x, y);
zz = 1 - (xx + yy);

% now that we've got the intersections, we can evaluate the function
f1 = @(x, y) 2*x.^2 + 3*y.^2 + 0.1;
Fxy1 = f1(xx, yy);
Fxy1(xx + yy > 1) = nan;

f2 = @(x, y) 3*x.^2 + 2*y.^2;
Fxy2 = f2(xx, yy);
Fxy2(xx + yy > 1) = nan;

f3 = @(x, y) (3*x.^2 + 2*y.^2) * 1000; % different order of magnitude
Fxy3 = f3(xx, yy);
Fxy3(xx + yy > 1) = nan;

subplot(1, 2, 1)
% This constructs the ternary axes
ternaxes(5);

% These are the coordinates of the compositions mapped to plot coordinates
[xg, yg] = terncoords(xx, yy);
% simpletri constructs the correct triangles
tri = simpletri(N);

hold on
% and now we can plot
trisurf(tri, xg, yg, Fxy1);
trisurf(tri, xg, yg, Fxy2);
hold off
view([137.5, 30]);

subplot(1, 2, 2);
ternaxes(5)
% Here we plot the line of intersection of the two functions
contour(xg, yg, Fxy1 - Fxy2, [0 0], 'r')
axis equal

如果要查找两个轮廓之间的交点，则可以有效地求解两个联立方程.这段额外的代码将为您解决这一问题(请注意，我现在也在上面的代码中使用了一些匿名函数):

EDIT 2: If you want to find the point of intersection between two contours, you are effectively solving two simultaneous equations. This bit of extra code will solve that for you (notice I've used some anonymous functions in the code above now, as well):

f1level = 1;
f3level = 1000;
intersection = fsolve(@(v) [f1(v(1), v(2)) - f1level; f3(v(1), v(2)) - f3level], [0.5, 0.4]);
% if you don't have the optimization toolbox, this command works almost as well
intersection = fminsearch(@(v) sum([f1(v(1), v(2)) - f1level; f3(v(1), v(2)) - f3level].^2), [0.5, 0.4]);

ternaxes(5)
hold on
contour(xg, yg, Fxy1, [f1level f1level]);
contour(xg, yg, Fxy3, [f3level f3level]);
ternplot(intersection(1), intersection(2), 1 - sum(intersection), 'r.');
hold off

这篇关于构造三元网格，在Matlab中评估网格和轮廓图上的函数的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！