0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 1 year has passed since last update.

四面体「2023名古屋工業大学前期【3】」をsympyとFreeCADでやってみた。

Last updated at Posted at 2023-04-03

FreeCADのマクロで作図しました。
(ワイヤーフレームです。ソリッドを勉強中。)

オリジナル 公式ホームページ

解答例 公式ホームページ

上と同じです。大学入試数学問題集成>【3】テキスト

sympyで(みなさんの模範解答を待ってます。)!

sympyで(sympy的?安易なやり方)

座標計算です。
・点Bの座標は、円と円の交点計算です。
・点Cの座標は、3球面の交点計算です。
・(4)について、角度が難しかったので、垂線の足を2倍に伸ばしました。
・教えて下さい。sympyのcoordinatesのかっこの外し方がわかりませんでした。

from sympy import *
var('x y z')
var('u v w')
def myABC(OQ, OA,OB,OC):
    return solve([
                  Eq(OQ.x,u*OA.x+v*OB.x+w*OC.x),
                  Eq(OQ.y,u*OA.y+v*OB.y+w*OC.y),
                  Eq(OQ.z,u*OA.z+v*OB.z+w*OC.z)
                 ],
                 [u,v,w])
def mySphereFormula(myP,myR):
    return ((x-myP.x)**2+(y-myP.y)**2+(z-myP.z)**2-myR**2).expand()
def myMenseki3D(O,A,B):
    P=A-O
    Q=B-O
    return Rational(1, 2) * sqrt(P.distance(Point(0, 0, 0)) ** 2 * Q.distance(Point(0, 0, 0)) ** 2 \
                                 - P.dot(Q) ** 2).simplify()
OA=BC=5
OB=AC=7
OC=AB=8
O=Point3D(0    ,0    ,0)
A=Point3D(OA   ,0    ,0)
Bxy=Circle(O,OB).intersection(Circle(A,AB))[1]
#B =Point(Bxy.coordinates[0],Bxy.coordinates[1],0)
B  =Point(Bxy.x,Bxy.y,0)
dum,Cxyz=solve([mySphereFormula(O,OC),
                mySphereFormula(A,AC),
                mySphereFormula(B,BC)],
               [x,y,z]
              )
C=Point(Cxyz)
print("#(1)",(A-O).dot(B-O),
             (B-O).dot(C-O),
             (C-O).dot(A-O)
     )
H=Plane(O,A,B).projection(C)
print("#(2)",myABC( (H-O), (A-O),(B-O),(C-O)) )
K=Line(O,A).projection(B)
print("#(3)",O.distance(K)/O.distance(A)) 
B_OA=Line(O,A).projection(B)
D=Line(A,B+(B_OA-B)*2).intersection(Line(O,H))[0]
print("#(4)",O.distance(D)/O.distance(H)) 
#(1) 5 44 20
#(2) {u: 19/30, v: 5/6, w: 0}
#(3) 1/5
#(4) 15/2

①私の環境は,pycharmです。
②よく聞くのは、Jupyterです。
③web上で、上記のソースを「SymPy Live shell」に、コピー貼り付けでもできました。
黒背景の右上に、マウスを移動すると、コピーマークが発生します。
計算時間は約1分です。資源の無駄使い。ごめんなさい。

FreeCADのマクロで作図

import FreeCAD
import Part
import Draft
import Mesh
#########################################################################################################
# 計算
from sympy import *
var('x y z')
var('u v w')
def myABC(OQ, OA,OB,OC):
    return solve([
                  Eq(OQ.x,u*OA.x+v*OB.x+w*OC.x),
                  Eq(OQ.y,u*OA.y+v*OB.y+w*OC.y),
                  Eq(OQ.z,u*OA.z+v*OB.z+w*OC.z)
                 ],
                 [u,v,w])
def mySphereFormula(myP,myR):
    return ((x-myP.x)**2+(y-myP.y)**2+(z-myP.z)**2-myR**2).expand()
def myMenseki3D(O,A,B):
    P=A-O
    Q=B-O
    return Rational(1, 2) * sqrt(P.distance(Point(0, 0, 0)) ** 2 * Q.distance(Point(0, 0, 0)) ** 2 \
                                 - P.dot(Q) ** 2).simplify()
OA=BC=5
OB=AC=7
OC=AB=8
O=Point3D(0    ,0    ,0)
A=Point3D(OA   ,0    ,0)
Bxy=Circle(O,OB).intersection(Circle(A,AB))[1]
#B =Point(Bxy.coordinates[0],Bxy.coordinates[1],0)
B  =Point(Bxy.x,Bxy.y,0)
dum,Cxyz=solve([mySphereFormula(O,OC),
                mySphereFormula(A,AC),
                mySphereFormula(B,BC)],
               [x,y,z]
              )
C=Point(Cxyz)
print("#(1)",(A-O).dot(B-O),
             (B-O).dot(C-O),
             (C-O).dot(A-O)
     )
H=Plane(O,A,B).projection(C)
print("#(2)",myABC( (H-O), (A-O),(B-O),(C-O)) )
K=Line(O,A).projection(B)
print("#(3)",O.distance(K)/O.distance(A)) 
B_OA=Line(O,A).projection(B)
D=Line(A,B+(B_OA-B)*2).intersection(Line(O,H))[0]
print("#(4)",O.distance(D)/O.distance(H)) 
##########################################################################################################
### 作図用
##########################################################################################################
### 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.4 mm'
    return
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]) )
    return
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)
    return
def myLine_S(*args):
    for i in range(1,len(args)):
        myLine(args[i-1],args[i])
    return
def myLine_C(*args):
    for i in range(1,len(args)):
        myLine(args[i-1],args[i])
    myLine(args[i],args[0])
    return
def myLine_H(*args):
    for i in range(1,len(args)):
        myLine(args[0],args[i])
    return
#
myLine_C  (O,A,B)
myLine_H  (C, O,A,B)
myLine    (B,K)
myLine    (C,H)
myLine    (O,D)
myLine    (O,H)
myLine    (O,K)
myTxtXYZ_S(O,"O",A,"A",B,"B",C,"C",D,"D",H,"H",K,"K")
#
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")

isometric方向です。

0png.png

上図の拡大図です。
Originは、非表示です。マウス操作で、線分と文字を赤色にしました。

1png.png

参考

この名古屋工業大学の問題と岐阜大学の問題は、似ていると思いました。

The intersection of this circle with another geometrical entity.

Project the given point onto the plane along the plane normal.

Returns the two coordinates of the Point.

0
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?