この記事を読むとできるようになること
matplotlibで三角グラフに3変数関数の解が描ける
matplotlibで三角グラフを作成する の続編です．
仕事柄，三角グラフに点をプロットするよりも3変数関数の縮退している解（三角グラフ上で直線になる）を描くことが多いので，その描き方のまとめです．
- 環境
- macOS Mojave 10.14.4
- Python 3.7.2
例として，三角グラフ上に次の3変数関数の解を描くとします．
ax+by+cz = E \quad\cdots (1)\\
x+y+z = 1 \qquad\cdots (2)\\
0 \leq x,y,z \leq 1
3変数で式が2本なので，解は縮退します．
点$(x,y,z)$を，2次元直交座標$(X,Y)$に変換して三角グラフ上に描画することにします．
変換の式は以下です．
X = z +\frac{x}{2}\\
Y = \frac{\sqrt{3}}{2}x
まず，$x$をx = np.arange(0,1.0,0.01)などで用意しておきます．
(1)(2)より$y$を消去すると，
z = ((E-b) -(a-b)x)/(c-b)
となります．
これを使って$(X,Y)$座標に直します．
triangle2.py
#!/usr/bin/env python
# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as plt
import math
x = np.arange(0,1.0,0.01)
fig = plt.figure()
ax1 = fig.add_subplot(111)
ax1.set_aspect('equal', 'datalim')
plt.tick_params(labelbottom=False, labelleft=False, labelright=False, labeltop=False)
plt.tick_params(bottom=False, left=False, right=False, top=False)
plt.gca().spines['bottom'].set_visible(False)
plt.gca().spines['left'].set_visible(False)
plt.gca().spines['right'].set_visible(False)
plt.gca().spines['top'].set_visible(False)
h = np.sqrt(3.0)*0.5
#ax+by+cz = E の係数
a = 2
b = 3
c = 6
E = 4
#####formule for triangle plot######
z = ((E-b) -(a-b)*x)/(c-b)
X = z +x*0.5
Y = h*x
ax1.plot(X,Y, color='red')
################# triangle plot settings #######################
#内側目盛
for i in range(1,10):
ax1.plot([i/20.0, 1.0-i/20.0],[h*i/10.0, h*i/10.0], color='gray', lw=0.5)
ax1.plot([i/20.0, i/10.0],[h*i/10.0, 0.0], color='gray', lw=0.5)
ax1.plot([0.5+i/20.0, i/10.0],[h*(1.0-i/10.0), 0.0], color='gray', lw=0.5)
#外周
ax1.plot([0.0, 1.0],[0.0, 0.0], 'k-', lw=2)
ax1.plot([0.0, 0.5],[0.0, h], 'k-', lw=2)
ax1.plot([1.0, 0.5],[0.0, h], 'k-', lw=2)
#頂点のラベル
ax1.text(0.52, h, 'x', fontsize=16)
ax1.text(-0.05, 0.03, 'y', fontsize=16, rotation=300)
ax1.text(0.95, -0.05, 'z', fontsize=16, rotation=60)
#軸ラベル
for i in range(1,10):
ax1.text(0.5+(10-i)/20.0, h*(1.0-(10-i)/10.0), '%d0' % i, fontsize=10)
ax1.text((10-i)/20.0-0.04, h*(10-i)/10.0+0.02, '%d0' % i, fontsize=10, rotation=300)
ax1.text(i/10.0-0.03, -0.025, '%d0' % i, fontsize=10, rotation=60)
#################################################################
plt.tight_layout()
plt.show()
こんな結果になります．
線がはみ出してしまいますね．．．
これは$y$が負になるときの解も描いてしまっているのが原因です．
なんとかしたいので，無理やりはみ出た分を消してみました．
triangle2_2.py
#!/usr/bin/env python
# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as plt
import math
from pylab import *
#xの範囲を適当に狭めてみる
x = np.arange(0,0.7,0.01)
fig = plt.figure()
ax1 = fig.add_subplot(111)
ax1.set_aspect('equal', 'datalim')
plt.tick_params(labelbottom=False, labelleft=False, labelright=False, labeltop=False)
plt.tick_params(bottom=False, left=False, right=False, top=False)
plt.gca().spines['bottom'].set_visible(False)
plt.gca().spines['left'].set_visible(False)
plt.gca().spines['right'].set_visible(False)
plt.gca().spines['top'].set_visible(False)
h = np.sqrt(3.0)*0.5
#ax+by+cz = E の係数
a = 2
b = 3
c = 6
E = 4
#####formule for triangle plot######
z = ((E-b) -(a-b)*x)/(c-b)
X = z +x*0.5
Y = h*x
#はみ出た部分を外枠で隠すために，むりやり折り返してみる
y1 = np.sqrt(3.0)*X
y2 = np.sqrt(3.0)*(1.0-X)
y = minimum(y1, y2)
Y0 = minimum(Y, y)
ax1.plot(X,Y0, color='red')
################# triangle plot settings #######################
#内側目盛
for i in range(1,10):
ax1.plot([i/20.0, 1.0-i/20.0],[h*i/10.0, h*i/10.0], color='gray', lw=0.5)
ax1.plot([i/20.0, i/10.0],[h*i/10.0, 0.0], color='gray', lw=0.5)
ax1.plot([0.5+i/20.0, i/10.0],[h*(1.0-i/10.0), 0.0], color='gray', lw=0.5)
#外周
ax1.plot([0.0, 1.0],[0.0, 0.0], 'k-', lw=2)
ax1.plot([0.0, 0.5],[0.0, h], 'k-', lw=2)
ax1.plot([1.0, 0.5],[0.0, h], 'k-', lw=2)
#頂点のラベル
ax1.text(0.52, h, 'x', fontsize=16)
ax1.text(-0.05, 0.03, 'y', fontsize=16, rotation=300)
ax1.text(0.95, -0.05, 'z', fontsize=16, rotation=60)
#軸ラベル
for i in range(1,10):
ax1.text(0.5+(10-i)/20.0, h*(1.0-(10-i)/10.0), '%d0' % i, fontsize=10)
ax1.text((10-i)/20.0-0.04, h*(10-i)/10.0+0.02, '%d0' % i, fontsize=10, rotation=300)
ax1.text(i/10.0-0.03, -0.025, '%d0' % i, fontsize=10, rotation=60)
#################################################################
plt.tight_layout()
plt.show()
安直にはみ出た部分を隠してみました．
あまりスマートなやり方ではないと思うので，もっと賢い方法をご存知でしたらご教示ください．