WOLFRAM

Wolfram言語 ＆ システムドキュメントセンター

外科用固定プレート

形状メッシュ

このノートブックでは，OpenCascadeLinkを使って，パラメータ設計に基づいた外科用固定プレートを作成する方法について説明する．

外科用固定プレートとは，折れた骨が治癒する過程において固定しサポートするために使われる医療器具である．チタンの固定プレートでも特に前腕用に設計されたものは，その強度，生体適合性，耐腐食性のために好まれる．これらのプレートは，骨折骨にネジで外科的に取り付けられ，骨が正しく一直線になるようにし，そして効果的な回復を促進する．パラメータ設計を使うことによって，患者のニーズに基づいて最適な治療を行い，器具を容易にカスタマイズすることができる．

固定プレートの設計は，固定プレートの長さと幅についてのパラメトリック形状関数につながる分析的記述に基づく．OpenCascadeを使った作成手順は，平面の記述から始め，それを後から押し出して3Dオブジェクトを形成する．

4.gif

前腕のとう骨上に置かれた外科用固定プレート（茶色で表示）．

パッケージをロードする：
Wolfram Language code: Needs["OpenCascadeLink`"] Needs["NDSolve`FEM`"]

形状

固定プレートは，中央部分と両側の2つの固定部分に分けられる．固定部分と中央部分はそれぞれ3つの穴がある．

サイズを希望する長さと幅の関数として表現した外科用固定プレートの平面の二次元工業的ダイアグラム．固定プレートは，それぞれの長さが長さ/3である3つの部分に分けられる．

3D形状を作成するために，2D設計図が作られ，それを押し出して3Dにする．この2D設計図を作るために，辺はすべてパラメータ化され，つなげられて，平面の2D面を形成する．その面から穴の部分を削除し，最後に穴を開けられた面を押し出して3Dにする．

最初のステップは2Dプレートの上の辺を作成することである．この長い辺は，それぞれの領域に希望数の穴をフィットするようにスケールされた正弦関数で表される3つの領域に分割することができる．形状特性はセンチメートルで表される．

形状パラメータを指定する：
Wolfram Language code: length = 16; thickness = 0.3; width = 1; waveAmplitude = 0.15;
穴の特性を指定する：
Wolfram Language code: numberFixationHoles = 3; radiusFixationHoles = 0.3;

正弦関数を使って，周期のそれぞれの正の波動サイクル内に固定プレートの穴を配置する．周期は固定の穴の数の2倍で使用可能な長さを割って求める．

固定の辺の周期について倍率 ϕ を計算する：
Wolfram Language code: ϕ = (length/3)(1/numberFixationHoles * 2)
辺の固定部分の穴について補助パラメトリック関数を作成する：
Wolfram Language code: fixationEdge[x_] := Sin[(π x/ϕ)]
辺の中央部分の穴についてスケールされた補助関数を作成する：
Wolfram Language code: centralEdge[x_] := fixationEdge[(x/6)]
上の辺の半分を表す補助関数を組み立てる：
Wolfram Language code: edge[x_] := (width/2) + waveAmplitude * Which[ x <= (length/3), fixationEdge[x], (length/3) < x <= (length/2), -centralEdge[x]]
上の辺を可視化する：
Wolfram Language code: edgePlot = Plot[edge[x + length / 2], {x, -length / 2, 0}, AspectRatio -> Automatic, ImageSize -> Large, Ticks -> {True, False}]

次のステップは，この辺のプロットをスプラインに変換することである．そのために辺が離散化，ダウンサンプリングされ，BSplineCurveに変換される．

辺を離散化する：
Wolfram Language code: discretedEdge = DiscretizeGraphics[edgePlot]
15番目の座標ごとに座標を選んで，座標数を減らす：
Wolfram Language code: topEdgeCoords = MeshCoordinates[discretedEdge][[1 ;; -1 ;; 15]];
一番右の点が において対称軸上にあるようにする：
Wolfram Language code: topEdgeCoords[[-1, 1]] = 0.;Show[discretedEdge, Graphics[{Red, Point[topEdgeCoords]}]]
上の辺の半分を表すスプラインを作成して可視化する：
Wolfram Language code: topEdgeSpline = BSplineCurve[topEdgeCoords]; Graphics[topEdgeSpline]

次のステップは，形状の左の辺を作ることであり，これは円状のセクターで表すことができる．難しいのはセクターの厳密なサイズを求めることである．この過程を容易にするため，半円と上の辺，そして半円と交差する上の辺の延長線を可視化すると便利である．

左の辺の円の中心を定義する：
Wolfram Language code: leftCircleCenter = {-length / 2 + width / 2, 0}
π/2と3π/2の間の円を定義する：
Wolfram Language code: circleSegement = Circle[leftCircleCenter, width, {π / 2, π}]
半円を通る，上の辺からの延長線を定義する：
Wolfram Language code: extensionLine = Line[{topEdgeCoords[[1]], {leftCircleCenter[[1]] - width, topEdgeCoords[[1, 2]]}}]
上の辺，延長線，半円を可視化する：
Wolfram Language code: Graphics[{topEdgeSpline, {Red, extensionLine}, circleSegement}]

延長線と半円の交差点が計算されると，左の辺のセクターを決定することができる．

上の辺の延長線と左の辺の交差点を計算する：
Wolfram Language code: intersectionPoint = RegionIntersection[extensionLine, circleSegement]
交差点の座標を抽出する：
Wolfram Language code: leftIntersectionCoord = intersectionPoint[[1, 1]]
左の円状セクターの角度を計算する：
Wolfram Language code: leftIntersectAngle = NSolveValues[leftIntersectionCoord[[1]] == leftCircleCenter[[1]] + width * Cos[θ], θ, MaxRoots -> 1][[1]]
上の辺，延長線の辺，左の辺を可視化する：
Wolfram Language code: leftLine = Line[{leftIntersectionCoord, topEdgeCoords[[1]]}]; leftEdge = Circle[leftCircleCenter, width, {2π - leftIntersectAngle, π}]; Graphics[{topEdgeSpline, {Red, leftLine}, leftEdge}]
y 軸に対して鏡映の変換を作成する：
Wolfram Language code: rt1 = ReflectionTransform[{1, 0}, {0, 0}];
辺を変換し可視化する：
Wolfram Language code: tr1 = TransformedRegion[topEdgeSpline, rt1]; tr2 = TransformedRegion[leftLine, rt1]; tr3 = TransformedRegion[leftEdge, rt1]; tr3 = tr3 /. Circle[c_List, r_, s_] :> Circle[c, r, {0, s[[1]]}]; Graphics[{topEdgeSpline, leftLine, leftEdge, {Blue, tr1}, {Red, tr2}, {Magenta, tr3}}]
辺を OpenCascade に変換する：
Wolfram Language code: e1 = OpenCascadeShape[leftEdge]; e2 = OpenCascadeShape[leftLine]; e3 = OpenCascadeShape[topEdgeSpline]; e4 = OpenCascadeShape[tr1]; e5 = OpenCascadeShape[tr2]; e6 = OpenCascadeShape[tr3];
3つの辺を単一のワイヤにする：
Wolfram Language code: w1 = OpenCascadeShapeWire[{e1, e2, e3, e4, e5, e6}]

次に，上半分全体を下に鏡映させる．

上の半分を鏡映変換して下の辺を作成する：
Wolfram Language code: rt2 = ReflectionTransform[{0, -1}, {0, 0}];
上の半分を下の半分に変換したものを可視化する：
Wolfram Language code: gt1 = GeometricTransformation[topEdgeSpline, rt2]; gt2 = GeometricTransformation[leftLine, rt2]; gt3 = GeometricTransformation[leftEdge, rt2]; Graphics[{topEdgeSpline, leftLine, leftEdge, {Blue, gt1}, {Red, gt2}, {Magenta, gt3}}]

実際の変換を上の辺全体に行う．

上の形状全体を下の部分に鏡映させる：
Wolfram Language code: w2 = OpenCascadeShapeTransformation[w1, rt2]
上と下の形状から面を作成する：
Wolfram Language code: face = OpenCascadeShapeFace[{w1, w2}]
境界メッシュを可視化する：
Wolfram Language code: bmesh = OpenCascadeShapeSurfaceMeshToBoundaryMesh[face]; bmesh["Wireframe"]

次のステップは，固定プレートに穴を開けることである．3つのグループの穴がある．2つのグループの固定部分の穴と中央の穴である．固定部分の穴は，正弦の空間的周期 , と一直線に並ぶ円で表され，中央の穴は，下に説明される超楕円境界で表される．

正弦周期に基づいて穴の距離を計算する：
Wolfram Language code: holesCentersAlongX = Table[-length / 2 + (2i - 3 / 2)ϕ, {i, 1, 3, 1}];
対象条件を使って他の3つの穴を追加する：
Wolfram Language code: holesCentersAlongX = Join[holesCentersAlongX, -#& /@ holesCentersAlongX];
固定部分の穴の形状を定義する：
Wolfram Language code: sideHoles = OpenCascadeShape[Disk[{#, 0}, radiusFixationHoles]]& /@ holesCentersAlongX;

超楕円は方程式 で表される．ここで と は半軸長を表す．さらに穴は45°回転させる．

超楕円パラメータを定義する：
Wolfram Language code: a = 0.2; b = 0.5; n = 3;
超楕円関数を定義する：
Wolfram Language code: hyperellipseFun[a_, b_, n_] := Abs[(x/a)]^n + Abs[(y/b)]^n == 1
超楕円を離散化する：
Wolfram Language code: discretizedHyperellipse = DiscretizeRegion[Rotate[Scale[Region@ImplicitRegion[hyperellipseFun[a, b, n], {x, y}], (width/1.25){1, 1}], -π / 4]];
超楕円の座標を抽出する：
Wolfram Language code: hyperellipsePts = MeshCoordinates[discretizedHyperellipse];
超楕円の点を正しく中心に置くための補助関数を定義する：
Wolfram Language code: centralHoleFun[l_] := {#[[1]] + l, #[[2]]}& /@ hyperellipsePts;
中央の穴の形状を定義する：
Wolfram Language code: centralHoles = OpenCascadeShape[Polygon[centralHoleFun[#]]]& /@ Table[((length/10))i, {i, -1, 1, 1}];
横と中央の穴を繋げて単一の形状にする：
Wolfram Language code: holes = OpenCascadeShapeUnion[Join[sideHoles, centralHoles]];

これで平面から穴の形を切り取ることが可能になった．

平面から穴を切り取る：
Wolfram Language code: base = OpenCascadeShapeDifference[face, holes]
穴を含む平面形状を可視化する：
Wolfram Language code: OpenCascadeShapeSurfaceMeshToBoundaryMesh[base]["Wireframe"["MeshElementStyle" -> Directive[FaceForm[LightBlue], EdgeForm[]]]]

メッシュ

最後に線形スイープを使って，平面を3Dオブジェクトに変換することができる．その後OpenCascade形状は，適切な大きさの希望の形状を表す境界メッシュと中空ではないメッシュに変換することができる．

線形スイープを使って3Dオブジェクトを押し出す：
Wolfram Language code: shape3D = OpenCascadeShapeLinearSweep[base, {{0, 0, 0}, {0, 0, thickness}}];
境界メッシュを作成する：
Wolfram Language code: (bmesh = OpenCascadeShapeSurfaceMeshToBoundaryMesh[shape3D])["Wireframe"]
中空ではない要素メッシュを作成する：
Wolfram Language code: (mesh = ToElementMesh[bmesh])["Wireframe"]
単位をセンチメートルにする：
Wolfram Language code: SetLengthUnit[mesh, "Centimeters"]
標準測定単位でプレートをスケールし直す：
Wolfram Language code: plateMeshed = ElementMeshCoordinateRescale[mesh, "Meters"];
最終的なメッシュ付きのプレートを可視化する：
Wolfram Language code: plateMeshed["Wireframe"["MeshElementStyle" -> Directive[EdgeForm[], FaceForm[Gray]]]]