正多面体の切断

正多面体を平面で切断したときの切断面を描くCindyアプリケーションを作ります。

わかりやすいところで、中学でも学ぶ、立方体の断面を考えます。次のような図です。

まず設計をしましょう。

設計

(1) 正多面体は rpolyhedron() で面データを取得して描く。（色塗りなし）

(2) 辺のリストを作る。

(3) 正多面体の各辺と切断する平面の交点を求める。

(4) (3)で求めた交点を結んだ図をplate3d() で描く。

では，順にやっていきましょう。

(1) 正多面体の面データを取得して描く。

rpolyhedron() は面を塗った多面体を描くので，オプション "nodisp" をつけて，戻り値の面データを取得し，その第１要素の頂点リストを使ってVertexEdgeFace() で描きます。Nohiddenbyfaces() で陰線処理もしましょう。なお，後の都合があるので，１辺を２としたら，原点が中に入らないよう，中心を(2,2,2)になるように平行移動しておきます。

fd=rpolyhedron("1",6,2*sqrt(3),["nodisp"]);

pd=transfer3d(fd_1,["Trans=[2,2,2]"]);

VertexEdgeFace("1",[pd,fd_2]);

Nohiddenbyfaces("1","phe3d1","phf3d1",["dr,2"],["do"]);

(2) 辺のリストを作る。

println(pd);

とすると，頂点リストが確かめられます。

[[0,0,4],[4,0,4],[4,4,4],[0,4,4],[0,0,0],[4,0,0],[4,4,0],[0,4,0]]

とコンソールに表示されます。

fd の第２要素は面を表すデータです。これも println(fd_2) で表示してみましょう。

[[4,3,7,8],[4,8,5,1],[4,1,2,3],[6,2,1,5],[6,5,8,7],[6,7,3,2]]

立方体の６つの面を表しています。それぞれの数は，頂点の番号です。

この面を表すデータを用いて，１２個の辺のリストを作ります。

まず，面ごとに隣り合う２つづつをとって辺データを作ります。

辺のリストを sidelist とします。１つの面について４つの辺データを作り，これを６つの面についておこないます。

sidelist=[];

repeat(length(fd_2),s,

sd=append(fd_2_s,fd_2_s_1);

repeat(length(sd)-1,t,

v1=sd_t;

v2=sd_(t+1);

sidelist=append(sidelist,[pd_v1,pd_v2]);

);

ここで，３行目で sd=append(fd_2_s,fd_2_s_1); としているのは，始点を末尾に追加しているのです。

たとえば，面を表すデータが [1,2,3,4] のとき，辺は [1,2],[2,3],[3,4],[4,1] の４つになります。この４番目のペアを取得するためです。

これで２４個のデータができます。（println(sidelist) で表示できます。）辺の数は１２ですが、各辺とも２回ずつ数えているので２４個になるのです。

これを整理してもいいですが、面倒なので、そのままにしましょう。平面との交点を調べるときに、同じものが２つづつできることになりますが、そのあとで処理します。

(3) 正多面体の各辺と切断する平面の交点を求める。

切断する平面Sは、適当に３点を決めて作ります。ひとまず，[6,0,0],[0,6,0],[0,0,6] としておきます。

この３点で決まる平面と線分の交点は、interps()で求められます。辺リスト sidelist のそれぞれに対して、この関数で交点を求め、それをリストにします。sidelist の要素のそれぞれに処理を適用するので apply()が便利です。交点リストは intpt としましょう。(intersection point)

interps()は交点がない場合は[i,i,i]を返してきますので、[i,i,i]でなければ intpt に追加していきます。

p1=[6,0,0];

p2=[0,6,0];

p3=[0,0,6];

intpt=[];

repeat(length(sidelist),s,

pt=interps(p1,p2,p3,sidelist_s_1,sidelist_s_2);

if(pt!=[i,i,i],intpt=append(intpt,pt));

);

intpt=sort(intpt);

println(intpt); で結果を表示してみましょう。

[[2,4,0],[0,4,2],[0,4,2],[0,2,4],[0,2,4],[2,0,4],[4,0,2],[2,0,4],[2,4,0],[4,2,0],[4,2,0],[4,0,2]]

となっています。同じものが２つずつあります。そこで，sort() で並べ替えて，１つおきにとっていくことにします。たとえば次のようにします。

plist=[];

repeat(length(intpt),s,

if(isodd(s),plist=append(plist,intpt_s));

);

isodd(s) は sが偶数かどうかの判定です。sはrepeatで繰り返されるときの番号ですから、偶数ならplistに追加(append)していくことになります。

さて、これで、断面の頂点リスト

[[0,2,4],[0,4,2],[2,0,4],[2,4,0],[4,0,2],[4,2,0]]

ができました。交点６つの座標のリストです。

(4) (3)で求めた交点を結んだ図をplate3d() で描く。

Cindyscript の convexhull3d() を使うと，この点のリストに対して面データを作ることができます。convexhull3d()は，凸型多面体の頂点リストから面データを作るものです。いまの３個の頂点は同一平面上にありますが、これに１つ点を追加して凸型多面体にすると、それに対して面データを作ってくれます。原点が中に入らないようにしましたので，追加する点を原点にします。

pd0=convexhull3d(prepend([0,0,0],plist));

pd=[];

forall(pd0_2,if(!contains(#,1),pd=append(pd,#)));

先ほど得た plist の前に prepend() で [0,0,0] を追加します。こうすれば、交点がいくつでも，常に原点が１番になります。

convexfull3d() の戻り値をpd0 とします。pd0の第２要素が面データのリストです。

こんどは、この面データの全てを走査するのに forall() を使います。リスト pd0_2 の要素を最初から順に走査します。調べている要素は # で表されます。その中に 1 が含まれていなければ (!contains() !は否定）ppに追加します。実際にできるのは１つだけです。

最後にこの面を表す頂点の順序に従って、多角形の頂点リストを作ります。plist から指定された番号の点を順に拾っていきます。

pp=[];

repeat(length(pd_1),s,

nn=pd_1_s;

pp=append(pp,plist_(nn-1));

);