= 垂心： r
6.4 Plain_Basic を使ったプログラミング例
========================================================
10 rem 三辺の長さから求める三角形の諸量、Triang1 20 data 350,400,450
30 defdbl L 40 read La,Lb,Lc
50 rem --- 面積 F 60 S=(La+Lb+Lc)/2
70 F=sqr(S*(S-La)*(S-Lb)*(S-Lc)) 80 print "面積 F= ";F
90 rem --- 頂角 A,B,C 100 SS=(La*La+Lb*Lb+Lc*Lc)/2
110 sinA=2*F/(Lb*Lc):cosA=(SS-La*La)/(Lb*Lc) 120 sinB=2*F/(Lc*La):cosB=(SS-Lb*Lb)/(Lc*La) 130 sinC=2*F/(La*Lb):cosC=(SS-Lc*Lc)/(La*Lb) 140 angA=atn2(sinA,cosA):print "angA= ";angA 150 angB=atn2(sinB,cosB):print "angB= ";angB 160 angC=atn2(sinc,cosc):print "angC= ";angC 170 print "angA+angB+angC = ";(angA+angB+angC) 180 rem --- 内接円半径 Ri 190 Ri=F/S: print "内接円半径 Ri= ";Ri
200 rem --- 外接円半径 Ro 210 Ro=La*Lb*Lc/(4*F): print "外接円半径 Ro= ";Ro 220 rem --- 作図
230 dperas
240 Xb= -La/2:Yb=-Lc/2:dpmove Xb,Yb 250 Xc=Xb+La:Yc=Yb:dpdraw Xc,Yc 260 Xa=Xb+Lc*cosB:Ya=Yb+Lc*sinB 270 dpdraw Xa,Ya: dpdraw Xb,Yb 280 dptext Xa,Ya, "A"
290 dptext Xb,Yb, "B"
300 dptext Xc,Yc, "C"
run
面積 F= 67082.0393249937 angA= 48.1896851042214 angB= 58.4118644947988 angC= 73.3984504009798 angA+angB+angC = 180
内接円半径 Ri= 111.803398874989 外接円半径 Ro= 234.787137637478
========================================================
10 rem 三辺の長さから求める三角形の諸量、Triang2 20 rem 内接円と外接円
30 dperas
40 data 60, 160, -175, -225, 175, -225 50 read Xa,Ya,Xb,Yb,Xc,Yc
60 dpmark Xa,Ya,3: dptext Xa,Ya, "A"
70 dpmark Xb,Yb,3: dptext Xb,Yb, "B"
80 dpmark Xc,Yc,3: dptext Xc,Yc, "C"
90 rem --- 以上は入力データ
100 print "--- 三辺の長さ La, Lb, Lc ---"
110 defdbl L
120 La=sqr((Xb-Xc)^2+(Yb-Yc)^2):print "La=";La 130 Lb=sqr((Xc-Xa)^2+(Yc-Ya)^2):print "Lb=";Lb 140 Lc=sqr((Xa-Xb)^2+(Ya-Yb)^2):print "Lc=";Lc 150 dpmove Xa,Ya
160 dpdraw Xb,Yb 170 dpdraw Xc,Yc 180 dpdraw Xa,Ya
190 rem --- 200 S=(La+Lb+Lc)/2
210 F=sqr(S*(S-La)*(S-Lb)*(S-Lc)) 220 SS=(La*La+Lb*Lb+Lc*Lc)/2 230 defdbl K
240 rem --- 内接円 250 K1=La/(2*s):K2=Lb/(2*s):k3=Lc/(2*s) 260 Xi=K1*Xa+K2*Xb+K3*Xc
270 Yi=K1*Ya+K2*Yb+K3*Yc
280 Ri=F/S: print "内接円半径 Ri= ";Ri 290 dpmark Xi,Yi,6:dptext Xi,Yi,"Ri"
300 dpcirc Xi,Yi,Ri
310 rem --- 外接円 320 K1=La*La*(SS-La*La)/(8*F*F)
330 K2=Lb*Lb*(SS-Lb*Lb)/(8*F*F) 340 K3=Lc*Lc*(SS-Lc*Lc)/(8*F*F) 350 Xo=K1*Xa+K2*Xb+K3*Xc 360 Yo=K1*Ya+K2*Yb+K3*Yc
370 Ro=La*Lb*Lc/(4*F): print "外接円半径 Ro= ";Ro 380 dpmark Xo,Yo,6:dptext Xo,Yo,"Ro"
390 dpcirc Xo,Yo,Ro run
--- 三辺の長さ La, Lb, Lc --- La= 350
Lb= 401.808412057289 Lc= 451.054320453756
内接円半径 Ri= 112.024420042262 外接円半径 Ro= 235.373273055978
========================================================
10 rem 三辺の長さから求める三角形の諸量、Triang3 20 rem 重心と垂心
30 dperas
40 data 60, 160, -175, -225, 175, -225 50 read Xa,Ya,Xb,Yb,Xc,Yc
60 dpmark Xa,Ya,3: dptext Xa,Ya, "A"
70 dpmark Xb,Yb,3: dptext Xb,Yb, "B"
80 dpmark Xc,Yc,3: dptext Xc,Yc, "C"
90 rem --- 以上は入力データ
100 print "--- 三辺の長さ La, Lb, Lc ---"
110 defdbl L
120 La=sqr((Xb-Xc)^2+(Yb-Yc)^2):print "La=";La 130 Lb=sqr((Xc-Xa)^2+(Yc-Ya)^2):print "Lb=";Lb 140 Lc=sqr((Xa-Xb)^2+(Ya-Yb)^2):print "Lc=";Lc 150 dpmove Xa,Ya
160 dpdraw Xb,Yb 170 dpdraw Xc,Yc 180 dpdraw Xa,Ya
190 rem --- 200 S=(La+Lb+Lc)/2
210 F=sqr(S*(S-La)*(S-Lb)*(S-Lc)) 220 SS=(La*La+Lb*Lb+Lc*Lc)/2 230 defdbl K
240 rem --- 重心 250 Xg=(Xa+Xb+Xc)/3
260 Yg=(Ya+Yb+Yc)/3 270 print "重心= ";Xg,Yg
280 dpmark Xg,Yg,6:dptext Xg,Yg,"G"
290 rem --- 垂心 300 K1=(SS-Lb*Lb)*(SS-Lc*Lc)/(4*F*F) 310 K2=(SS-Lc*Lc)*(SS-La*La)/(4*F*F) 320 K3=(SS-La*La)*(SS-Lb*Lb)/(4*F*F) 330 Xh=K1*Xa+K2*Xb+K3*Xc
340 Yh=K1*Ya+K2*Yb+K3*Yc 350 print "垂心= ";Xh,Yh
360 dpmark Xh,Yh,6:dptext Xh,Yh,"H"
run
--- 三辺の長さ La, Lb, Lc --- La= 350
Lb= 401.808412057289 Lc= 451.054320453756 重心= 20 -96.6666666666667 垂心= 60 -154.805194805195
========================================================
10 rem 三頂点座標から求める三角形の諸量、Triang4 20 rem 慣性モーメントと主軸の計算
30 dperas
40 data 125, 220, -260, -170, 40, -240 50 read Xa,Ya,Xb,Yb,Xc,Yc
60 rem --- 以上は入力データ 70 dpmark Xa,Ya,3: dptext Xa,Ya, "A"
80 dpmark Xb,Yb,3: dptext Xb,Yb, "B"
90 dpmark Xc,Yc,3: dptext Xc,Yc, "C"
100 rem --- 面積 110 F=((Xc-Xb)*(Ya-Yb)-(Yc-Yb)*(Xa-Xb))/2 120 print "面積= ";F
130 dpmove Xa,Ya 140 dpdraw Xb,Yb 150 dpdraw Xc,Yc
160 dpdraw Xa,Ya
170 rem --- 重心 180 Xg=(Xa+Xb+Xc)/3
190 Yg=(Ya+Yb+Yc)/3 200 print "重心= ";Xg,Yg
210 dpmark Xg,Yg,6:dptext Xg,Yg,"G"
220 rem --- 慣性モーメント 230 defdbl J:dim J[2,2]
240 J[1,1]=(Xa*Xa+Xb*Xb+Xc*Xc)*F/12 + 3*F*Xg*Xg/4 250 J[1,2]=(Xa*Ya+Xb*Yb+Xc*Yc)*F/12 + 3*F*Xg*Yg/4 260 J[2,2]=(Ya*Ya+Yb*Yb+Yc*Yc)*F/12 + 3*F*Yg*Yg/4 270 print "Jxx,Jxy,Jyy=";J[1,1],J[1,2],j[2,2]
280 rem --- 主軸の計算 290 J0=(J[1,1]+J[2,2])/2
300 Jx=(J[1,1]-J[2,2])/2 310 Jy=J[1,2]
320 Sign=1:if Jx<0 then Sign=-1 330 Jr=sqr(Jx*Jx + Jy*Jy) 340 Ang=atn2(Sign*Jy,Sign*Jx)/2 350 print "主軸の回転角度=";Ang 360 J11=J0+Sign*Jr
370 J22=J0-Sign*Jr
380 print "主軸回りの慣性モーメント=";J11,J22 390 rem --- 主軸の作図 400 ccc=cos(Ang):sss=sin(Ang)
410 X1=Xg+300*ccc:Y1=Yg+300*sss:dpmove X1,Y1 420 X1=Xg-300*ccc:Y1=Yg-300*sss:dpdraw X1,Y1 430 X1=Xg+300*sss:Y1=Yg-300*ccc:dpmove X1,Y1 440 X1=Xg-300*sss:Y1=Yg+300*ccc:dpdraw X1,Y1 run
面積= 71975
重心= -31.6666666666667 -63.3333333333333 Jxx,Jxy,Jyy= 562904479.166667
480733020.833333 1025643750 主軸の回転角度= -32.1495908225377
主軸回りの慣性モーメント= 260761091.769097 1327787137.39757
========================================================
10 rem 多角形の幾何学的諸量 "PcSection1"
20 rem --- 頂点座標を左回りに準備する
30 dim P[2,20]: REM reserve an array for points 40 data -0.665, 0.00, -0.665, -0.17, -0.385, -0.17 50 data -0.085, -0.27, -0.085, -0.63, -0.205, -0.75 60 data -0.205, -1.00, 0.205, -1.00, 0.205, -0.75 70 data 0.085, -0.63, 0.085, -0.27, 0.385, -0.17 80 data 0.665, -0.17, 0.665, 0.00, 9000 / 90 read P
100 rem: --- 頂点数を勘定する 110 Np=1
120 IF P[1,NP+1]=9000 THEN GOTO 140 130 NP=NP+1: goto 120
140 P[1,NP+1]=P[1,1]: P[2,NP+1]=P[2,1]
150 rem --- データ確認の作図をする 160 dperas: dpwind 0,-0.5,2 170 dpmove P[1,1],P[2,1]
180 for K=2 to NP+1 190 dpdraw P[1,K],P[2,K]
200 next
210 rem --- 三角形要素の集合を計算する 220 defdbl I,J
230 for K=1 to NP
240 dA=(P[1,K]*P[2,K+1]-P[2,K]*P[1,K+1])/2 250 Area=Area+dA
260 dGx=(P[1,K]+P[1,K+1])/3:Gx=Gx+dA*dGx 270 dGy=(P[2,K]+P[2,K+1])/3:Gy=Gy+dA*dGy
280 dIxx=(P[1,K]*P[1,K]+P[1,K+1]*P[1,K+1])*dA/12 290 dIxx=dIxx+3*dA*dGx*dGx/4: Ixx=Ixx+dIxx 300 dIyy=(P[2,K]*P[2,K]+P[2,K+1]*P[2,K+1])*dA/12 310 dIyy=dIyy+3*dA*dGy*dGy/4: Iyy=Iyy+dIyy 320 dIxy=(P[1,K]*P[2,K]+P[1,K+1]*P[2,K+1])*dA/12 330 dIxy=dIxy+3*dA*dGx*dGy/4: Ixy=Ixy+dIxy 340 next K
350 rem --- 重心位置
360 Xg=Gx/Area:Yg=Gy/Area:rem center of gravity 370 print "Xg, Yg = ";Xg,Yg
380 dpmark Xg,Yg,3:dptext Xg,Yg, "G"
390 rem --- 断面積、慣性モーメント 400 print "Area = ";Area
410 Ixx=Ixx-Area*Xg*Xg 420 Iyy=Iyy-Area*Yg*Yg 430 Ixy=Ixy-Area*Xg*Yg
440 print "Ixx,Ixy,Iyy=";Ixx,Ixy,Iyy 450 rem --- 主軸の計算 460 J0=(Ixx+Iyy)/2
470 Jx=(Ixx-Iyy)/2 480 Jy=Ixy
490 Sign=1:if Jx<0 then Sign=-1 500 Jr=sqr(Jx*Jx + Jy*Jy)
510 Ang=atn2(Sign*Jy,Sign*Jx)/2 520 print "主軸の回転角度=";Ang 530 J11=J0+Sign*Jr
540 J22=J0-Sign*Jr
550 print "主軸回りの慣性モーメント=";J11,J22 560 rem --- 主軸の作図 570 ccc=cos(-Ang):sss=sin(-Ang)
580 X1=Xg+ccc:Y1=Yg+sss:dpmove X1,Y1 590 X1=Xg-ccc:Y1=Yg-sss:dpdraw X1,Y1 600 X1=Xg+sss:Y1=Yg-ccc:dpmove X1,Y1 610 X1=Xg-sss:Y1=Yg+ccc:dpdraw X1,Y1 run
rest variables ignored
Xg, Yg = 0 -0.361717557251908 Area = 0.4716
Ixx,Ixy,Iyy= 0.03641561 2.16840434497101E-19 0.0516399687786259
主軸の回転角度= -8.16063382709388E-16 主軸回りの慣性モーメント= 0.03641561 0.0516399687786259
========================================================
_______________________________________________________________________________