Segédlet BMEGEMMAGM5 “Végeselem módszer alapjai” tárgy laborgyakorlatához.
Készítette: Dr. Kossa Attila, BME Műszaki Mechanikai Tanszék, 2013. szeptember 18.

Határozzuk meg az alábbi szelvény súlypontját, a főmásodrendű nyomatékait és főirányait!

felkor=Graphics[{EdgeForm[Thick],Yellow,Disk[{0,0},2,{Pi/2,Pi}]}];
haromszog=Graphics[{EdgeForm[Thick],Yellow,Polygon[{{0,2},{0,0},{4,0}}]}];
szelveny=Show[{felkor,haromszog},AxesOrigin->{0,0},Axes->True,PlotRangePadding->1,AxesStyle->Arrowheads[0.05],AxesLabel->{x,y}]

A szelvény két részből áll: egy negyedkörből (r=2) és egy derékszögű háromszögből (2r alaphosszúságú, r magasságú). Ezen részkeresztmetszetek keresztmetszeti jellemzői segítségével számítjuk a keresett mennyiségeket.
Az integrálási sorok egyszerűsítése végett itt most megadhatjuk, hogy az r változó 0-nál nagyobb értékű. Ezzel gyorsíthatjuk a számításokat azon részeknél, ahol a Mathematica általánosabban próbálná megoldani a feladatot.

$Assumptions=r>0;

A továbbiakban az 1-es rész alatt a negyedkört, 2-es alatt pedig a háromszöget értjük.
Az 1-es rész területe:

Területek számítása

A negyedkör terültét megkapjuk ha kiszámítjuk a negyed körív alatti területet, ami az alábbi integrállal adható meg:

A1=Integrate[Integrate[1,{y,0,(Sqrt[(r^2)-(x^2)])}],{x,-r,0}]

A háromszög (r magasságú és 2r alaphosszúságó) területe a felső kontúr általi függvény alatti területet jelenti::

A2=Integrate[Integrate[1,{y,0,r-x/2}],{x,0,2r}]

A teljes szelvény területe:

A=FullSimplify[A1+A2]

Numerikusan:

N[A]

Súlypontok számítása:

Az 1-es keresztmetszet súlypontjának x koordinátája:

xs1=(Integrate[Integrate[x,{y,0,(Sqrt[(r^2)-(x^2)])}],{x,-r,0}])/A1

A 2-es keresztmetszet súlypontjának x koordinátája:

xs2=(Integrate[Integrate[x,{y,0,r-x/2}],{x,0,2r}])/A2

Az 1-es keresztmetszet súlypontjának y koordinátája:

ys1=(Integrate[Integrate[y,{y,0,(Sqrt[(r^2)-(x^2)])}],{x,-r,0}])/A1

A 2-es keresztmetszet súlypontjának y koordinátája:

ys2=(Integrate[Integrate[y,{y,0,r-x/2}],{x,0,2r}])/A2

A teljes szelvény súlypontjának koordinátái:

xs=(Integrate[Integrate[x,{y,0,(Sqrt[(r^2)-(x^2)])}],{x,-r,0}]+Integrate[Integrate[x,{y,0,r-x/2}],{x,0,2r}])/A

ys=(Integrate[Integrate[y,{y,0,(Sqrt[(r^2)-(x^2)])}],{x,-r,0}]+Integrate[Integrate[y,{y,0,r-x/2}],{x,0,2r}])/A

Ábrázolás: teljes szelvény, a részkeresztmetszetek súlypontjai (kék), az eredő súlypont (piros), a részkeresztmetszetek súlypontjait összekötő egyenes:

S1=Graphics[{PointSize[Medium],Blue,Point[{xs1,ys1}/.{r->2}]}];
S2=Graphics[{PointSize[Medium],Blue,Point[{xs2,ys2}/.{r->2}]}];
sline=Graphics[{Dashed,Line[{({xs1,ys1}/.{r->2}),({xs2,ys2}/.{r->2})}]}];
S=Graphics[{PointSize[Large],Red,Point[{xs,ys}/.{r->2}]}];
Show[{szelveny,sline,S1,S2,S}]

Másodrendű nyomatékok számítása:

A részkeresztmetszetek másodrendű nyomatékai az x tengelyre:

Ix1=Integrate[Integrate[(y^2),{y,0,(Sqrt[(r^2)-(x^2)])}],{x,-r,0}]
Ix2=Integrate[Integrate[(y^2),{y,0,r-x/2}],{x,0,2r}]

A részkeresztmetszetek másodrendű nyomatékai az y tengelyre:

Iy1=Integrate[Integrate[(x^2),{y,0,(Sqrt[(r^2)-(x^2)])}],{x,-r,0}]
Iy2=Integrate[Integrate[(x^2),{y,0,r-x/2}],{x,0,2r}]

A részkeresztmetszetek másodrendű nyomatékai az xy tengelypárra:

Ixy1=Integrate[Integrate[x*y,{y,0,(Sqrt[(r^2)-(x^2)])}],{x,-r,0}]
Ixy2=Integrate[Integrate[x*y,{y,0,r-x/2}],{x,0,2r}]

A Steiner-tétel segítségével átszámíthatjuk a fenti mennyiségeket a saját súlyponton átmenő és az x-szel párhuzamos tengelyekre (AZÉRT SZEREPEL MÍNUSZ JEL AZ ELTOLÁSNÁL, MERT A SÚLYPONT FELÉ KÖZELEDÜNK ÉS NEM TÁVOLODUNK!!!):

Ixs1=FullSimplify[Ix1-(ys1^2)A1]
Ixs2=FullSimplify[Ix2-(ys2^2)A2]

A Steiner-tétel segítségével átszámíthatjuk a fenti mennyiségeket a saját súlyponton átmenő és az y-nal párhuzamos tengelyekre (AZÉRT SZEREPEL MÍNUSZ JEL AZ ELTOLÁSNÁL, MERT A SÚLYPONT FELÉ KÖZELEDÜNK ÉS NEM TÁVOLODUNK!!!):

Iys1=FullSimplify[Iy1-(xs1^2)A1]
Iys2=FullSimplify[Iy2-(xs2^2)A2]

Hasonlóképpen a tengelypárra számított mennyiségek átszámítása:

Ixys1=FullSimplify[Ixy1-ys1*xs1*A1]
Ixys2=FullSimplify[Ixy2-ys2*xs2*A2]

A teljes szelvény súlypontjára átszámított mennyiségek:

IxS=FullSimplify[(Ixs1+(((ys1-ys))^2)A1)+(Ixs2+(((ys2-ys))^2)A2)]
IyS=FullSimplify[(Iys1+(((xs1-xs))^2)A1)+(Iys2+(((xs2-xs))^2)A2)]
IxyS=FullSimplify[(Ixys1+(ys-ys1)A1)+(Ixys2+(ys-ys2)(xs-xs2)A2)]

Másodrendű nyomatékok mátrixa:

A mátrix felírásnál a tengelypárra számított mennyiséget negatív előjellel kell beírni!

Imat={{IxS,-IxyS},{-IxyS,IyS}};
MatrixForm[Imat]

Ennek a mátrixnak a sajátétékei:

Eigenvalues[Imat]

A fenti paraméteres megoldásba behelyettesítve az r=2-t, megkapjuk a feladatra vonatkozó numerikus értékeket:

Imatnum=FullSimplify[Imat/.{r->2}];
MatrixForm[Imatnum]

Amennyiben a 2-t, nem egész számként helyettesítjük be, hanem “numerikus” 2-ként:

Imatnum=FullSimplify[Imat/.{r->2.}];
MatrixForm[Imatnum]

A keresett sajátértékek:

Eigenvalues[Imatnum]

Ezek a keresett fő másodrendű nyomatékok. Ha nem lennének sorbarendezve akkor az alábbi paranccsal megtehetjük:

Sort[%,Greater]

A sajátvektorok irányát megadó vektorok:

{n1,n2}=Eigenvectors[Imatnum]

3D esetben NEM egységvektorként adja ki a Mathematica! Ekkor célszerű normalizálni:

n1=Normalize[n1]
n2=Normalize[n2]

A sajátértékek és a hozzá tartozó sajátvektorok egy paranccsal is számíthatóak:

Eigensystem[Imatnum]

A sajátvektorok irányának szemléltetése:

linen1=Graphics[{Red,Thick,Arrowheads[0.05],Arrow[{({xs,ys}/.{r->2}),(({xs,ys}/.{r->2})+n1*3)}]}];
linen2=Graphics[{Red,Thick,Arrowheads[0.05],Arrow[{({xs,ys}/.{r->2}),(({xs,ys}/.{r->2})+n2*3)}]}];
Show[{szelveny,S,linen1,linen2}]

Az 1-es főirány x-tengellyel bázárt szöge:

φ=ArcTan[n1[[2]]/n1[[1]]]

Szögben kifejezve:

φ/Degree