空間図形「【昭和大（医2022】」三角錘の体積をFreeCADでやってみた。ソリッド

三角錘の体積の計算です。

オリジナル
youtubeです。

https://www.youtube.com/watch?v=DEHmjKGLRjM

大学入試数学問題集成＞
未登録


# sympyだけで(外積にはいろんな情報が詰まっている)
(勉強中)


# sympyだけで(垂線計算はPlane().distance()を使っています。)
SymPy Liveでできました。
以下サイトに、ソースコードを貼り付けて、改行です。

https://live.sympy.org/

以下
(1)「外積」
(2)「ヘロンの公式」
(3)「ベクトルを用いた三角形の面積の公式」
(4)「行列式」


```python
# 計算
from sympy import *
def Area3(A,B,C):
    return Rational(1,2)*sqrt((A.distance(B)**2*A.distance(C)**2-(B-A).dot(C-A)**2))
def Herons(a,b,c):
    s=(a+b+c)/2
    return sqrt(s*(s-a)*(s-b)*(s-c)).simplify()
def Herons3(p,q,r):
    return Herons(p.distance(q), q.distance(r), r.distance(p))
def SimentaiTaiseki(PTO,PTA,PTB,PTC):
    return Matrix([[PTA.x-PTO.x, PTA.y-PTO.y, PTA.z-PTO.z],
                   [PTB.x-PTO.x, PTB.y-PTO.y, PTB.z-PTO.z],
                   [PTC.x-PTO.x, PTC.y-PTO.y, PTC.z-PTO.z]]
                  ).det()/6
def myCross(A,B):
    return Point(A.y*B.z-A.z*B.y,A.z*B.x-A.x*B.z,A.x*B.y-A.y*B.x)
O=Point3D(0, 0, 0)
A=Point3D(1,-1,-1)
B=Point3D(2, 0, 2)
C=Point3D(2, 2, 0)
print("#(1)",Rational(1,3)*Rational(1,2)*myCross(B-A,C-A).distance(Point3D(0,0,0))*Plane(A,B,C).distance(O))
print("#(2)",Rational(1,3)*Herons3(A,B,C)                                         *Plane(A,B,C).distance(O))
print("#(3)",Rational(1,3)*Area3(A,B,C)*Plane(A,B,C).distance(O))
print("#(4)",abs(SimentaiTaiseki(O,A,B,C)))
#(1) 2
#(2) 2
#(3) 2
#(4) 2
```



# FreeCADのマクロで(体積計算にsympy使用)ワイヤーフレーム
上のソースコードを、FreeCADマクロに追加しました。
図は、そのままではありません。手を加えています。
　
図1.アイソメトリック方向、Origin表示

![png1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/109698/ebe4444a-c6fa-a82d-289a-6b421126e95b.png)

図2.ビューを回転、Origin非表示
Hが△ABC上に見えますか？

![png2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/109698/dfff332a-9983-a809-c8d9-5b93a7c9c66b.png)


```python
import FreeCAD
import Part
import Draft
import Mesh
#########################################################################################################
# 計算
from sympy import *
def Area3(A,B,C):
    return Rational(1,2)*sqrt((A.distance(B)**2*A.distance(C)**2-(B-A).dot(C-A)**2))
def Herons(a,b,c):
    s=(a+b+c)/2
    return sqrt(s*(s-a)*(s-b)*(s-c)).simplify()
def Herons3(p,q,r):
    return Herons(p.distance(q), q.distance(r), r.distance(p))
def SimentaiTaiseki(PTO,PTA,PTB,PTC):
    return Matrix([[PTA.x-PTO.x, PTA.y-PTO.y, PTA.z-PTO.z],
                   [PTB.x-PTO.x, PTB.y-PTO.y, PTB.z-PTO.z],
                   [PTC.x-PTO.x, PTC.y-PTO.y, PTC.z-PTO.z]]
                  ).det()/6
O=Point3D(0, 0, 0)
A=Point3D(1,-1,-1)
B=Point3D(2, 0, 2)
C=Point3D(2, 2, 0)
print("#(1)",Rational(1,3)*Herons3(A,B,C)*Plane(A,B,C).distance(O))
print("#(2)",abs(SimentaiTaiseki(O,A,B,C)))
print("#(3)",Rational(1,3)*Area3(A,B,C)*Plane(A,B,C).distance(O))
#(1) 2
#(2) 2
#(3) 2
H=Plane(A,B,C).projection(O)
# #########################################################################################################
# 3D作図
def myXYZ2Txt(A):
    return '(' + str(A.x) + ',' + str(A.y) + ',' + str(A.z) + ')'
def myTxtXYZ(A,myWedgei):
    P5x=float(A.x)
    P5y=float(A.y)
    P5z=float(A.z)
    p5 = FreeCAD.Vector(P5x, P5y, P5z)
    myText = Draft.makeText(myWedgei, p5)
    myText.Label = myWedgei
    FreeCADGui.ActiveDocument.ActiveObject.FontSize = '0.5 mm'
def myTxtXYZ_S(*xy_tx):
    for i in range(1,int(len(xy_tx)/2)+1):
        myTxtXYZ(xy_tx[2*i-2],xy_tx[2*i-1]+myXYZ2Txt(xy_tx[2*i-2]) )
def myLine(A,B):
    Ax,Ay,Az=float(A.x),float(A.y),float(A.z)
    Bx,By,Bz=float(B.x),float(B.y),float(B.z)
    pl = FreeCAD.Placement()
    pl.Rotation.Q = (0.4247081540122249, 0.17592004639554645, 0.33985110062924484, 0.8204732460821097)
    pl.Base = FreeCAD.Vector(-3.9166066876399563, -2.1670824762243774, 1.7495260956243028)
    points = [FreeCAD.Vector(Ax,Ay,Az), FreeCAD.Vector(Bx,By,Bz)]
    line = Draft.make_wire(points, placement=pl, closed=False, face=True, support=None)
    Draft.autogroup(line)
def myLine_S(*args):
    for i in range(1,len(args)):
        myLine(args[i-1],args[i])
    return 0
def myLine_H(*args):
    for i in range(1,len(args)):
        myLine(args[0],args[i])
    return 0
myTxtXYZ_S(O,"O",A,"A",B,"B",C,"C",H,"H")
myLine_S(A,B,C,A)
myLine(O,A)
myLine(O,B)
myLine(O,C)
myLine(O,H)
myLine_H(O, A,B,C,H)
doc = App.activeDocument()
App.ActiveDocument.addObject("App::Origin", "Origin")
App.ActiveDocument.getObject('Origin').Visibility = True
App.ActiveDocument.recompute()
Gui.activeDocument().activeView().viewAxonometric()
Gui.SendMsgToActiveView("ViewFit")
```

# FreeCADのマクロで(体積計算にsympy不使用)FreeCADの三角錐ソリッドのFreeCADの体積を表示
ソリッドで作図。transparency30。
図は、そのままではありません。手を加えています。

![png.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/109698/4d3ca031-3100-b0d3-1d1b-db5fc049e6c2.png)

```python
import FreeCAD
import Part
import DraftTools
import Draft
import Mesh
from sympy import *
##########################################################################################################
# 計算
from sympy import *
def Area3(A,B,C):
    return Rational(1,2)*sqrt((A.distance(B)**2*A.distance(C)**2-(B-A).dot(C-A)**2))
def Herons(a,b,c):
    s=(a+b+c)/2
    return sqrt(s*(s-a)*(s-b)*(s-c)).simplify()
def Herons3(p,q,r):
    return Herons(p.distance(q), q.distance(r), r.distance(p))
def SimentaiTaiseki(PTO,PTA,PTB,PTC):
    return Matrix([[PTA.x-PTO.x, PTA.y-PTO.y, PTA.z-PTO.z],
                   [PTB.x-PTO.x, PTB.y-PTO.y, PTB.z-PTO.z],
                   [PTC.x-PTO.x, PTC.y-PTO.y, PTC.z-PTO.z]]
                  ).det()/6
def myCross(A,B):
    return Point(A.y*B.z-A.z*B.y,A.z*B.x-A.x*B.z,A.x*B.y-A.y*B.x)
O=Point3D(0, 0, 0)
A=Point3D(1,-1,-1)
B=Point3D(2, 0, 2)
C=Point3D(2, 2, 0)
print("#(1)",Rational(1,3)*Rational(1,2)*myCross(B-A,C-A).distance(Point3D(0,0,0))*Plane(A,B,C).distance(O))
print("#(2)",Rational(1,3)*Herons3(A,B,C)                                         *Plane(A,B,C).distance(O))
print("#(3)",abs(SimentaiTaiseki(O,A,B,C)))
print("#(4)",Rational(1,3)*Area3(A,B,C)*Plane(A,B,C).distance(O))
#(1) 2
#(2) 2
#(3) 2
#(4) 2
H=Plane(A,B,C).projection(O)
##########################################################################################################
def Tetrahedron(O,A,B,C):
    P1x,P1y,P1z=float(O.x),float(O.y),float(O.z)
    P2x,P2y,P2z=float(A.x),float(A.y),float(A.z)
    P3x,P3y,P3z=float(B.x),float(B.y),float(B.z)
    P4x,P4y,P4z=float(C.x),float(C.y),float(C.z)
    xDummy=0
    yDummy=0
    zDummy=0
    a =[[xDummy,yDummy,zDummy],  #dummy
        [P1x,P1y,P1z],
        [P2x,P2y,P2z],
        [P3x,P3y,P3z],
        [P4x,P4y,P4z]]
    m1 = Mesh.Mesh([a[1],a[2],a[3]])     
    m2 = Mesh.Mesh([a[1],a[2],a[4]])
    m3 = Mesh.Mesh([a[2],a[3],a[4]])
    m4 = Mesh.Mesh([a[3],a[1],a[4]])
    mesh=Mesh.Mesh()
    mesh.addMesh(m1)
    mesh.addMesh(m2)
    mesh.addMesh(m3)
    mesh.addMesh(m4)
    shape = Part.Shape()        
    shape.makeShapeFromMesh(mesh.Topology,0.05)
    solid = Part.makeSolid(shape).removeSplitter()
    myWedgei='myTetrahedron_'+'x' + str(P1x) + 'y' + str(P1y)+ 'z' + str(P1z)
    f = doc.addObject("Part::Feature",myWedgei)
    f.Shape = solid
    # gui_doc.getObject('myTetrahedron_x0_0y0_0z5_0').ShapeColor = (1.00,1.00,0.50)
    target_object=doc.getObjectsByLabel("myTetrahedron_x0_0y0_0z0_0")[0]
    myVolume=target_object.Shape.Volume
    print("myVolume",myVolume)
    return myVolume
def myXYZ2Txt(A):
    return '(' + str(A.x) + ',' + str(A.y) + ',' + str(A.z) + ')'
def myTxtXYZ(A,myWedgei):
    P5x=float(A.x)
    P5y=float(A.y)
    P5z=float(A.z)
    p5 = FreeCAD.Vector(P5x, P5y, P5z)
    myText = Draft.makeText(myWedgei, p5)
    myText.Label = myWedgei
    FreeCADGui.ActiveDocument.ActiveObject.FontSize = '0.2 mm'
def myTxtXYZ_S(*xy_tx):
    for i in range(1,int(len(xy_tx)/2)+1):
        myTxtXYZ(xy_tx[2*i-2],xy_tx[2*i-1]+myXYZ2Txt(xy_tx[2*i-2]) )
def myLine(A,B):
    Ax,Ay,Az=float(A.x),float(A.y),float(A.z)
    Bx,By,Bz=float(B.x),float(B.y),float(B.z)
    pl = FreeCAD.Placement()
    pl.Rotation.Q = (0.4247081540122249, 0.17592004639554645, 0.33985110062924484, 0.8204732460821097)
    pl.Base = FreeCAD.Vector(-3.9166066876399563, -2.1670824762243774, 1.7495260956243028)
    points = [FreeCAD.Vector(Ax,Ay,Az), FreeCAD.Vector(Bx,By,Bz)]
    line = Draft.make_wire(points, placement=pl, closed=False, face=True, support=None)
    Draft.autogroup(line)
O=Point3D(0, 0, 0)
A=Point3D(1,-1,-1)
B=Point3D(2, 0, 2)
C=Point3D(2, 2, 0)
doc = App.activeDocument()
myTxtXYZ_S(A,"A",B,"B",C,"C",H,"H")
myLine(O,H)
myV=Tetrahedron(O,A,B,C)
myTxtXYZ(O,myXYZ2Txt(O)+"V={:.3f}".format(myV))
doc = App.activeDocument()
App.ActiveDocument.addObject("App::Origin", "Origin")
App.ActiveDocument.getObject('Origin').Visibility = True
App.ActiveDocument.recompute()
Gui.activeDocument().activeView().viewAxonometric()
Gui.SendMsgToActiveView("ViewFit")
```
# 参考

https://qiita.com/mrrclb48z/items/15beeb08dcd74e8bc183

https://forum.freecadweb.org/viewtopic.php?style=1&f=16&t=59757