진동학 및 실습Project과 목 명 : 진동학 및 실습교 수 님 : 이수훈설계 문제에 대해서 적합한 ‘회전하는 편심질량’ 물체는 CD-ROM Drive으로 설정하였습니다. CD-ROM Drive의 경우 회전 시 발생하는 편심에 의해 발생하는 진동을 줄이기 위해 추가적으로 동흡진기를 설치하여 설치 전과 설치 후의 진동을 비교해보고자 합니다.1. 회전하는 편심질량을 가지는 적절한 대상을 선정하여 1자유도 시스템으로모델링하고, 어떻게 모델링 하였는지 설명하시오.CD-ROM Drive는 CD를 넣은 후 회전함에 따라 약간의 무게 차이로 인해 작지만 편심질량을 가지게 되므로 이번 프로젝트의 적합한 주제라고 생각했습니다. 또한 회전하는 편심질량에 의한 힘의 작용 방향은 위-아래 방향만 있다고 가정하였습니다.모델을 1자유도계로 Base에 질량m _{1}으로, 모델을 단순화하여c _{1}의 감쇠상수,k _{1}의 강성이 작용하는 것으로 편심질량을m _{3}으로 동흡진기를 설치하였을 때 질량을m _{2}라 하고 감쇠상수c _{2}, 스프링상수k _{2}로 모델링 할 수 있습니다.2. 모델을 구성하는 각 요소의 적절한 값을 선정하고 그 이유를 설명하라.편심질량을m _{3}으로 0.001 kg이라 가정하였고, 그 질량 중심까지 거리r을 0.01m이라 가정하였습니다.k _{1} 값은 CD-ROM 내의 Rubber mount가 5개가 있고, 일반적인 값k _{1} ``=`2215 TIMES 5=11075N/m으로 가정,c _{1} 값은 0.27kg으로 가정하였습니다. CD-ROM Drive 내의 CD가 회전할 때 각속도w _{max} =1200`Hz로 보편적인 가정 할 수 있습니다.3-1. 모델의 운동에 대하여 기술하시오.m _{1} ddot{x _{1}} `+c _{1} dot{x} +k _{1} x _{1} =m _{3} ew _{s} ^{2} sinwt 이 운동방정식의 해를x _{1} (t)=X _{1} sin(wt+ phi _{`1} )이라고 가정하면` dot{x _{1}} `(t)=wX _{1} cos(wt+ phi _{`1} ),``` ddot{x _{1}} `(t)=-w ^{2} X _{1} sin(wt+ phi _{`1} ) 이를 운동방정식에 대입합니다.-m _{1} `w ^{2} X _{1} sin(wt+ phi _{`1} )+wc _{1} X _{1} cos(wt+ phi _{`1} )c _{1} +k _{1} X _{1} sin(wt+ phi _{`1} )=m _{3} ew _{s} ^{2} sinwtX _{1}을 묶어주면X _{1} [(k _{1} -m _{1} w ^{2} )sin(wt+ phi _{`1} )+cwcos(wt+ phi _{`1} )]=mew _{s} ^{2} sinwtX _{1} [(k _{1} -m _{1} w ^{2} )(sinwtcos phi _{`1} +coswtsin phi _{`1} )+c _{1} w(coswtcos phi _{`1} -sinwtsin phi _{`1} )]=m _{3} ew _{s} ^{2} sinwtsin`wt를 제외한 남은 항들의 값은 0입니다.X _{1} [(k _{1} -m _{1} w ^{2} )sin phi _{`1} +c _{1} wcos phi _{`1} ]=m _{3} ew _{s} ^{2} sinwtX _{1} [(k _{1} -m _{1} w ^{2} )cos phi _{`1} +c _{1} wsin phi _{`1} ]=0위의 두 식을 이용하여X _{1}과phi _{1}를 구하면THEREFORE ``X _{1} = {m _{3} ew _{s} ^{2}} over {sqrt {(k _{1} -m _{1} w ^{2} ) ^{2} +c _{1} w ^{2}}},phi _{1} =-tan ^{-1} ( {c _{1w}} over {k _{1} -m _{1} w ^{2}} )3-2. 모델의 운동을 컴퓨터를 이용하여 풀어보시오.Vibration.m 파일%%%%%%%%%%%%%%%%%%%%%%% 함수설정 %%%%%%%%%%%%%%%%%%%%%%%function u_prime = vib(t,u)%%%%%%%%%%%%%%%%%%%%%%% 초기값 설정 %%%%%%%%%%%%%%%%%%%%%%%m1=0.25; % Sled Base(CD 올려 놓는 곳) 질량c=0.27; % 감쇠상수k=11075; % 스프링상수(강성)w=1200/60/2/pi; % 각속도(Hz단위를 rad/s로 바꿔줍니다.)m2=0.001; % 동흡진기의 질량r=0.01; % 편심질량에서 중심까지의 거리f0=m2*r*w^2; % 편심질량에서 중심방향으로 작용하는 원심력force=f0*sin(w*t); % Sled Base에 가해지는 외력u_prime=zeros(2,1); % u_prime(1)과 u_prime(2)가 들어가는 함수 2행1열 zeros 생성%m ddot{x} +c acute{x} +kx=0을u _{1} (t)=x',u _{2} (t)=x으로 변경한 결과u_prime(1)=-(c/m1)*u(1)-(k/m1)*u(2)+force/m1;u_prime(2)=u(1);%%%%%%%%%%%%%%%%%%%%%%% 함수 ODE45를 이용한 진동해석 %%%%%%%%%%%%%%%%%%%%%%%clc;clear;%%%%%%%%%%%%%%%%%%%%%%% 초기값 설정 %%%%%%%%%%%%%%%%%%%%%%%time=0:0.05:25; % 시간범위 0~25초, 시간 간격을 0.05초로 설정w=1200/60/2/pi; % 각속도(Hz단위를 rad/s로 바꿔줍니다.)m2=0.001; % 동흡진기 질량r=0.01; % 편심질량에서 중심까지의 거리f0=m2*r*w^2; % 편심질량에서 중심방향으로 작용하는 원심력force=f0*sin(w*time); % Sled Base에 작용하는 원심력, 변수 t대신 시간범위 time적용%%%%%%%%%%%%%%%%%%%%%%% 실습시간에 배운 미분방정식 풀이 %%%%%%%%%%%%%%%%%%%%%%%initial_condition=[0 0]; % 위 문제의 진동의 초기조건 설정[t, Out] = ode45('vib', time, initial_condition); % ODE45함수를 이용하여 해석%%%%%%%%%%%%%%%%%%%%%%% 그래프 작성 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 입력 외력 표시 그래프(붉은색) %%%%%%%%%%%%%%%%%%%%%%%figure, plot(t, force, 'b'); % axis([15 time(end) -1.5 1.5]);hold ontitle('Input Signal'); xlabel('Time'); ylabel('Amplitude');% 그래프의 제목은‘Input Siganl', x축은 ‘Time', y축은 '진폭’%%%%%%%%%%%%%%%%%%%%%%% 출력 진동 표시 그래프(푸른색) %%%%%%%%%%%%%%%%%%%%%%%figure, plot(t, Out(:,1), 'r'); % axis([15 time(end) -1.5 1.5]);hold ontitle('Output Signal'); xlabel('Time'); ylabel('Amplitude');% 그래프의 제목은‘Output Siganl', x축은 ‘Time', y축은 '진폭’%%%%%%%%%%%%%%%%%%%%Output Signal Graph와 Input Signal Graph%%%%%%%%%%%%%%%%%%%%figure,subplot(2,1,1); plot(t,force, 'r');%axis([30 40 -1.5 1.5]);hold ontitle('Input Signal'); xlabel('Time'); ylabel('Amplitude');subplot(2,1,2); plot(t,Out(:,1), 'b');%axis([30 40 -0.3 0.3]);hold ontitle('Output Signal'); xlabel('Time'); ylabel('Amplitude');Solution.m 파일Wnwn.m 파일%%%%%%%%%%%%%%%%%%%%%%% Frenquency Response 함수 설정 %%%%%%%%%%%%%%%%%%%%%%%function wnwn%%%%%%%%%%%%%%%%%%%%%%% 초기값 입력 %%%%%%%%%%%%%%%%%%%%%%%m1=0.25;c1=0.27;k1=11075;m3=0.001;r=0.01;Frequency=zeros(4000,1); % 각각의 주파수를 넣기 위한 포맷 형성xp=zeros(4000,1); % 각각의 위상을 넣기 위한 포맷 형성%%%%%%%%%%%%%%%%%%%%%%% Frequency, xp(크기), th(위상)계산 %%%%%%%%%%%%%%%%%%%%%%%for i=1:4000f=m3*r*(i/10)^2;Frequency(i)=(i/10)/(2*pi);xp(i)=f/sqrt((k1-m1*(i/10)^2)^2+(c1*(i/10))^2);if (i/10)
유체역학 및 실습보고서목 차1프로젝트 목적2이론 및 실제 계산2-1 베르누이 방정식질량 보존 법칙Newton의 제 2법칙 - 선형운동량2-2 실제 계산 (수식정리)2-3 Matlab를 통한 그래프3설계 시 고려된 재료의 물성4상세 설계도5기타5-1 MATLAB graph m-file5-2 고찰1. 프로젝트 목적아래의 그림과 같이 중앙으로 공급된 유체(물)가 위쪽과 아래쪽 원판사이를 흐를 때 아래쪽 판이 떨어지지 않고 위쪽 판에 흡착하게 됩니다. 이러한 현상이 일어나는 원리를 유체역학적으로 분석해보고, 아랫판이 윗판에 대한 흡착력이 최대가 되도록 설계를 구상해보는 것이 목적입니다. 조건은 아래의 도표와 같습니다.항 목조 건디스크지름 50cm 이하입구관(연결관)길이 10-30cm외경 3cm작동유체물 (최대유량 = 40 LPM)2. 이론 및 실제 계산2-1 베르누이 방정식베르누이 방정식은 유체에 가해지는 일이 없는 경우에 유체의 속도와 압력, 위치에너지 사이의 관계를 나타낸 식입니다. 일반적으로 나타내면 아래와 같습니다.{P _{}} over {gamma } ``+` {V _{} ^{2}} over {2g} ``+`z````=`H (Constant)[압력 수두, 속도 수두, 위치수두, 전수두]또한 베르누이 방정식이 성립하기 위해서는 아래와 같은 조건이 필요합니다.① 유체가 유선을 따라 움직입니다.② 유체는 점성력이 없습니다. (비점성유동)③ 유동은 정상상태입니다.④ 비압축성 유동을 합니다.질량보존의 법칙시스템의 질량의 시간 변화율은 검사체적의 질량의 시간 변화율과 검사표면을 통한 질량의 유동률로 나타낼 수 있는 레이놀드 수송 정리에서 비롯되었습니다.{Partial } over {Partial t} int _{sys} ^{} {rho ``dV```=`} {Partial } over {Partial t} int _{cv} ^{} {rho ``dV} ````+` int _{cs} ^{} {rho ` vec{v} ` BULLET vec{n} ``dA````}이에 대해서 검사체척을 고정되고 변형이 없는 검사체적에 대해 정리하면 다음과 같은 질량 보존의 법칙 (연속 방정식)을 도출할 수 있습니다.{PARTIAL } over {PARTIAL t} int _{cv} ^{} {rho dv+ int _{cs} ^{} {rho vec{V} BULLET vec{n} dA=0}}이를 질량유동dot{m}에 대해서 나타내면 다음과 같습니다.{Partial } over {Partial t} int _{cv} ^{} {rho dv+` sum _{} ^{} dot{m _{out}} ```-} ` sum _{} ^{} dot{m _{in}} ```=`0`Newton의 제 2법칙 - 선형 운동량시스템의 선형운동량의 시간에 대한 변화율은 시스템에 작용하는 외력의 합으로 나타낼 수 있습니다. 운동량이란 질량과 속도의 곱으로 Newton의 제 2운동법칙에 따라 나타낼 수 있으며, 이를 고정되고 변형이 없는 검사체적에 대하여 나타내면 다음과 같습니다.{PARTIAL } over {PARTIAL t} int _{cv} ^{} {vec{V} `` rho ```dv+ int _{cs} ^{} {vec{V`} `` rho ` vec{V} BULLET vec{n} ```dA= sum _{} ^{} F}}2-2 실제계산 (수식정리)hQP, V Pa, VaDrR간단하게 나타내면 위와 같은 형태라 할 수 있습니다.{P _{i`n}} over {gamma } ``+` {V _{i`n} ^{2}} over {2g} ``+`z _{i`n} ````=` {P _{out}} over {gamma } ``+` {V _{out} ^{2}} over {2g} ``+`z _{o`u`t} (베르누이 방정식에 적용)들어가고, 나가는 유체의 위치에너지는 동일하므로z _{i``n} ``=`z _{out}{P _{i`n}} over {gamma } ``+` {V _{i`n} ^{2}} over {2g} ```=` {P _{out}} over {gamma } ``+` {V _{out} ^{2}} over {2g} `이를 정리하면P _{i``n} ```=` {gamma } over {2g} (V _{out} ^{2} ``-`V _{i``n} ^{2} )여기서 유량Q _{i`n`} ``=`Q _{out}이므로V``=` {Q} over {A} ``=` {Q} over {2 pi rh}정리하게 되면P _{i```n} ``=` {rho } over {2} ( {Q} over {2 pi h} ) ^{2} ( {1} over {R ^{2}} ``-`` {1} over {r ^{2}} )위의 식을 잘 살펴보면 왜 흡착이 일어나는지 알 수 있습니다.R > r 이므로P _{i``n}
기구학 및 실습Project 2과 목 명 : 기구학 및 실습교 수 님 : 송봉섭목 차1.문제 및 설계 조건2.이 론2-1 Normalized Polynomial (General Polynomial Curve)2-2 예제 6.12-3 Swinging Roller Follower3.문제 해결 과정3-1 A번3-2 B번3-3 C번3-4 D번3-5 E번4.참 고 문 헌1. 문제 및 설계 조건설계해야 할 캠의 모델은 ‘요동롤러종동절’와 같은 핀에 연결된 ‘4-Bar-Link’입니다.주어진 조건으로는 종동절의 변위는 교과서 예제 6.1의 운동과 동일하다는 것입니다.2. 이 론2-1 Normalized Polynomial (General Polynomial Curve)Polynomial Curve는 그 차수를 마음대로 조절 할 수 있기 때문에 많은 설계요구조건을 만족 시킬 수 있습니다.S( theta )``=`A _{0} `+A _{1} theta +A _{2} theta ^{2} +A _{3} theta ^{3} +` CDOTS `+A _{n} theta ^{n}(n : Boundary Condition의 개수 - 1,A _{i} : Polynomial 상수,theta : Cam Rotation Angle)사용자가 정의해주는 특정 구간에서의 양끝 또는 어느 부분에서의 Boundary Condition에 따라 식의 차수가 결정될 수 있으며, 아래와 같은 여러 종류의 조건을 주어, 그 식을 결정할 수 있습니다. 즉, Curve 양끝에서 위치, 속도,가속도, 저크 등의 조건을 줄 수도 잇고, 또한 최대 속도의 값, 최대 가속도의 값 등도 조건 값이 될 수 있습니다. 어떠한 조건이 주어지더라도 이를 만족하는 Polynomial Curve 식을 구할 수 있습니다.①theta = theta _{0} ```` RARROW `````S,` {dS} over {d theta } ,` {d ^{2} S} over {d theta ^{2}} ,`` {d ^{3} S} over {d thet치x _{p} ``=`F _{x} cos theta +F _{y} sin theta #y _{p`} ``=`-F _{x} sin theta +F _{y} cos theta 종동절 로울러 중심좌표 (x _{0} ,`y _{0})x _{0} ``=`x _{p} +L _{f} cos(S( theta )+ phi _{0} - theta )#y _{0} ``=`y _{p} +L _{f} sin(S( theta )+ phi _{0} - theta ) 이를 정리하면x``=`F _{x} cos theta ``+F _{y} sin theta +L _{f} cos beta `+ sigma {R _{f}} over {sqrt {1+[ {M} over {N} ] ^{2}}}#y``=`-F _{x} sin theta ``+F _{y} cos theta +L _{f} sin beta `+ sigma {( {M} over {N} )R _{f}} over {sqrt {1+[ {M} over {N} ] ^{2}}}M``=`F _{x} sin theta -F _{y} cos theta +L _{f} sin beta ( {dS} over {d theta } -1),N``=`-F _{x} cos theta -F _{y} cos theta +L _{f} cos beta ( {dS} over {d theta } -1)beta =S( theta )+ phi _{0} - theta ,``=` phi - theta ``` RARROW `` {d beta } over {d theta } = {d phi } over {d theta } -13. 문제 해결 과정3-1. 변위곡선의 각 구간별로, 거리(X)를 각도theta (캠 회전각)의 함수로 표시하라.Matlab으로 구현했을 때, 그래프는 아래와 같습니다.3-2. 각 구간별로, 변위, 속도, 가속도 곡선을 시간 t(sec)의 함수로 표시하고 이를 그래프로 그리고 최대 가속도를 구하라.일반적으로 속도V= {dS} over {dt}이며, 이는 Chain-Rulea1-theta0;dP(x)=(30*Pi(x)^2)-(60*Pi(x)^3)+(30*Pi(x)^4);dS(x)=dP(x)*ds(x)/dtheta(x);d2P(x)=(60*Pi(x))-(180*Pi(x)^2)+(120*Pi(x)^3);d2S(x)=d2P(x)*ds(x)/dtheta(x);endif c==1figure(1);subplot(311)theta=(a:1:b)*pi/180;time=(0.5*theta)/(2*pi); % 0.5theta/2pi의 경우 1회전 각도 2pi를 시간으로 변환과정plot((0.5*theta)/(2*pi),S); grid ontitle('Displacement Graph'); xlabel('Time(sec)'); ylabel('Displacement(mm)');axis([0 0.5 0 8]);subplot(312 )% 120rpm을 rad/s로 변환하면 4pi, 속도의 경우 각속도를 곱해주어야 합니다.plot((0.5*theta)/(2*pi),dS*4*pi); grid ontitle('Velocity Graph'); xlabel('Time(sec)'); ylabel('Velocity(mm/s)'); % 0.5theta/2pi의 경우 1회전 각도 2pi를 시간으로 변환과정axis([0 0.5 -200 200]);subplot(313) % 120rpm을 rad/s로 변환하면 4pi, 가속도의 경우 각속도의 제곱을 곱해주어야 합니다.plot((0.5*theta)/(2*pi),d2S*(4*pi)^2); grid ontitle('Aceelation Graph'); xlabel('Time(sec)'); ylabel('Acceleration(mm/s^2)'); % 0.5theta/2pi의 경우 1회전 각도 2pi를 시간으로 변환과정axis([0 0.5 -10000 10000]);%%% 최대 가속도 계산 %%%%Acceleration_max=max(d2S*(4*pi)^2);disp(Acceleration_max)elseif c==0theta=(a:ver {2a sqrt {1-X ^{2}}},{d phi } over {d theta } = {dX} over {d theta } {d phi } over {dX}pi LEQ theta LEQ {5pi } over {4}{dX} over {d theta } `=`0,{d phi } over {d theta } =0{5pi } over {4} LEQ theta LEQ {3pi } over {2}X( theta )=X _{0} +30( {theta - {5 pi } over {4}} over {{pi } over {4}} ) ^{3} -45( {theta - {5 pi } over {4}} over {{pi } over {4}} ) ^{4} +18( {theta - {5 pi } over {4}} over {{pi } over {4}} ) ^{5}{dX} over {d theta } =90( {theta - {5 pi } over {4}} over {{pi } over {4}} ) ^{2} -180( {theta - {5 pi } over {4}} over {{pi } over {4}} ) ^{3} +90( {theta - {5 pi } over {4}} over {{pi } over {4}} ) ^{4}{d phi } over {dX} ``=`- {1} over {2a sqrt {1-X ^{2}}},{d phi } over {d theta } = {dX} over {d theta } {d phi } over {dX}{3pi } over {2} LEQ theta LEQ {7pi } over {4}{dX} over {d theta } `=`0,{d phi } over {d theta } =0{7pi } over {4} LEQ theta LEQ 2 pi X( theta )=X _{0} -30( {theta - {7 pi } over {4}} over {{pi } over {4}} ) ^{3} +45( {theta - {7 pi } over {4}} over {{pi 225도까지 경우s0=0; s1=0; % dwell이므로 s는 0으로 유지theta0=180*pi/180; theta1=225*pi/180;elseif (225*pi/180=theta) % polynomial, 225도부터 270도까지 경우s0=0; s1=3; % s는 225도에 0, 270도에 3mmtheta0=225*pi/180; theta1=270*pi/180;elseif (270*pi/180=theta) % dwell, 270도부터 315도까지 경우s0=3; s1=3; % s는 270도에 3mm, 315도에 3mmtheta0=270*pi/180; theta1=315*pi/180;elseif (315*pi/180=theta) % polynomial, 315도부터 360도까지 경우s0=3; s1=0; % s는 315도에 3mm, 360에 0theta0=315*pi/180; theta1=360*pi/180;end% Polynomial을 이용해 변위 S를 나타냅니다.Pi(x)=(theta-theta0)/(theta1-theta0);P(x)=(10*Pi(x)^3)-(15*Pi(x)^4)+(6*Pi(x)^5); % 문제 6.1의 Polynomial 함수.S(x)=s0+(P(x)*(s1-s0));% 변위를 미분하면 속도 함수를 얻을 수 있습니다.ds(x)=s1-s0; dtheta(x)=theta1-theta0;dP(x)=(30*Pi(x)^2)-(60*Pi(x)^3)+(30*Pi(x)^4);dS(x)=dP(x)*ds(x)/dtheta(x);endif c==1 % 그래프를 나타내는 경우.figure(1);subplot(211)theta=(a:1:b)*pi/180;plot(theta*180/pi,S); grid ontitle('변위(S) 그래프'); xlabel('theta(degree)'); ylabel('변위(S) (mm)');axis([0 360 0 8]);subplot(212)plot(theta*180/pi,dS); grid ontitle('속도(V) 그래프');-479
기구학 및 실습Project 1과 목 명 : 기구학 및 실습교 수 님 : 송봉섭[Project 1]아래와 같은 Mechanism에서 Link가 일정속도로 회전할 때 이 Mechanism이 움직이는 모양을 Computer Graphics를 이용하여, 화면을 통해 보고자 한다. 별첨 Disk에 있는 라는 컴퓨터 프로그램을 이용하여 아래 문제를 해결하라. (이의 사용법은 부록 B에 있음.) 아래 Mechanism에서 조건은 다음과 같습니다.(Dx,`Dy)``=`(6,`-2),`L _{2} ``=`3`cm,`L _{3} ``=`6cm,``L _{4} ``=`6cm,`BE`=`4cm,`EF``=`2cm,`` alpha =60 DEG1) 이 프로그램이 제공하는 [STEP MOVEMENT] MENU를 이용하여, Link 2의 각도를 변화시키면서, Mechanism 의 운동을 관찰하라.2) Link 2 가 1 회전하는 동안, 커플러 점 F 의 이동 경로를 관찰하고, 이 화면을 Copy 하여 제출하라.3) Link 2 가 1 회전하는 동안 Link 3 및 Link 4 의 각도가 어떻게 변화하는지를 [POSITION ANALYSIS] MENU를 이용하여 관찰하고, 이 그림을 Copy 하여 제출하라.4) Velocity analysis, Animation 등을 수행해보고, 이 프로그램의 사용법을 충분히 익히도록 하라.5) Link 2 가 1 회전하는 동안, 커플러 점 F 의 이동 경로를 확인 할 수 있는MATLAB(또는 유사 프로그램 언어) 코드를 작성하고, 문제 2의 결과와 비교하라.그림#1-1 A 4-Bar Mechanism for Project 11) 이 프로그램이 제공하는 [STEP MOVEMENT] MENU를 이용하여, Link 2의 각도를 변화시키면서, Mechanism 의 운동을 관찰하라.theta `=`0 DEG (링크L _{2}의 각도theta )theta `=`90 DEG (링크L _{2}의 각도theta )theta `=`180 DEG (링크L _{2}의 각도theta )theta `=`270 DEG (링크L _{2}의 각도theta )대표적으로theta 가 0˚, 90˚, 180˚, 270˚ 등으로 바꾸어 보면서 운동을 관찰했습니다.이 때 점 F의 경로를 표시하도록 Trace Point를 이용하여, 흰 점으로 나타내었습니다.2), 3) Link 2 가 1 회전하는 동안 Link 2, Link 3 및 Link 4 의 각도가 어떻게 변화하는지를 [POSITION ANALYSIS] MENU를 이용하여 관찰하고, 이 그림을 Copy 하여 제출하라.Link 2 각도(theta )(˚)Link 3 각도(psi )(˚)Link 4 각도(phi )(˚)Link 2 각도(theta )(˚)Link 3 각도(psi )(˚)Link 4 각도(phi )(˚)1038.82573.7952019031.481129.76121031.14869.6202120036.018131.37832024.96867.8902221040.806132.53843020.18668.1862322045.753133.25754016.56870.0312423050.765133.54265013.88773.0122524055.731133.38776011.96176.8042625060.530132.76987010.65781.1552726065.015131.6399809.88585.8722827069.005129.92010909.58890.8002928072.271127.492111009.73795.8103029074.512124.1911211010.320100.7893130075.342119.7991312011.341105.6363231074.298114.0781413012.813110.2613332070.920106.8631514014.752114.5813433064.98898.2901615017.172118.5283534056.88289.0711716020.076122.0473635047.70980.4851817023.452125.1013736038.82573.7951918027.270127.672점 F의 경로theta (˚)FxFytheta (˚)FxFy058.09544.8411904.05235.6721063.37745.8952002.06533.1512066.20647.0682100.54730.7853066.93348.510220-0.50128.6194065.96650.143230-1.07426.7025063.66551.796240-1.16025.0926060.32553.288250-0.74223.8617056.19554.4592600.21023.0938051.49455.1912701.74422.8869046.41655.4072803.93823.35610041.14155.0722906.92124.61911035.82854.18830010.89526.77512030.61852.79531016.14129.84013025.63050.95932022.95733.63114020.96048.76533031.42537.63415016.67646.31134041.00241.08016012.82543.69235050.37743.4321709.43240.99836058.09544.8411806.50738.305※ 이 값은 10배 크게 표현되었으며, 이유는 마지막장에 언급되었듯이, 초기 데이터를 10배 크게 주었기 때문입니다. 실제 문제에서 얻을 수 있는 값은 위 표의 값에서 0.1을 곱한 값입니다.10˚단위로 각도(theta )의 변화를 주어, 위에 나타난 표로 정리하였습니다.프로그램으로 구현한 화면은 30˚ 단위로 나타내면 다음과 같습니다.theta `=3`0 DEG (링크L _{2}의 각도theta )theta `=6`0 DEG (링크L _{2}의 각도theta )theta `=`90 DEG (링크L _{2}의 각도theta )theta `=12`0 DEG (링크L _{2}의 각도theta )theta `=150 DEG (링크L _{2}의 각도theta )theta `=18`0 DEG (링크L _{2}의 각도theta )theta `=`210 DEG (링크L _{2}의 각도theta )theta `=24`0 DEG (링크L _{2}의 각도theta )theta `=`270 DEG (링크L _{2}의 각도theta )theta `=`300 DEG (링크L _{2}의 각도theta )theta `=330 DEG (링크L _{2}의 각도theta )theta `=360 DEG (링크L _{2}의 각도theta )5) Link 2 가 1 회전하는 동안, 커플러 점 F 의 이동 경로를 확인 할 수 있는MATLAB(또는 유사 프로그램 언어) 코드를 작성하고, 문제 2의 결과와 비교하라.※ 스캔 원본 파일은 같이 첨부하였습니다.MATLAB M-File Codeclear,clc;theta=[0:pi/180:2*pi]'; % theta의 회전각 범위는 0~2pi이며 radian 단위입니다.Dx=6; Dy=-2; L2=3; L3=6; L4=6; BE=4; EF=2;% Dx, Dy, L2, L3, L4, BE, EF 문제에 주어진 값을 입력합니다.P=2*L4*(Dx-L2*cos(theta));Q=2*L4*(Dy-L2*sin(theta));R=(Dx^2)+(Dy^2)+(L2^2)+(L4^2)-(L3^2)-(2*L2*(Dx*cos(theta)+Dy*sin(theta)));cospi=((-P.*R)-(Q.*sqrt((P.^2)+(Q.^2)-(R.^2))))./((P.^2)+(Q.^2));sinpi=((-Q.*R)+(P.*sqrt((P.^2)+(Q.^2)-(R.^2))))./((P.^2)+(Q.^2));% P와 R은 theta에 따라서 달라지므로, 열 벡터(1*N 행렬)을 이루게 됩니다. 따라서 각각의 성분에 대한 값을 고려하기% 위해서, .을 찍어서 나타내었습니다.pi2=atan2(sinpi,cospi); % pi라고 입력시 상수 pi와 중복되므로 pi2으로 지정했습니다.% atan2(높이, 밑변)을 입력시 각도를 출력해주는 명령어입니다. 이를 통해 Link4가 이루는 각, pi2를 구합니다.cosbeta=(Dx+(L4*cos(pi2))-(L2*cos(theta)))/L3;sinbeta=(Dy+(L4*sin(pi2))-(L2*sin(theta)))/L3;% 직각삼각형을 만들어(스캔본 참조), Link3이 이루는 각, beta를 구합니다.beta=atan2(sinbeta,cosbeta);Fx=(L2*cos(theta))+(BE*cos(beta))+(EF*cos(beta+(pi/3)));Fy=(L2*sin(theta))+(BE*sin(beta))+(EF*sin(beta+(pi/3)));plot(Fx,Fy);xlabel('F의 x좌표'), ylabel('F의 y좌표');grid onbeta3=beta.*(180/pi); % 위의 과정은 모두 radian으로 나타내었으므로, 이를 degree로 바꿔주는 과정입니다.pi23=pi2.*(180/pi);F=[Fx,Fy,beta3,pi23];문제2와 문제5의 결과를 비교하면 아래와 같습니다.4BarSim으로 구현한 결과MATLAB으로 구현한 결과문제에서 언급한대로 비교하기 위해 유사한 그래프 개형에서 최대값, 최소값을비교해 볼 수 있습니다. 먼저 MATLAB 구현 결과 최대, 최소는 아래와 같습니다.Fx 최소값, -0.1186Fx 최대값, 6.6906Fy 최소값, 2.2877Fy 최대값, 5.54104BarSim을 1회 회전각도를 1˚로 설정 후 1회전 하는 동안 관찰되는 Fx, Fy의 최대, 최소는 아래와 같습니다.Fx 최소값, -1.186Fx 최대값, 66.943Fy 최소값, 22.877Fy 최대값, 55.410MATLAB을 통해 얻은 결과와 비교시, 4BarSim으로 얻은 결과는 MATLAB 결과의 10배입니다.(물론 Fx 최대의 경우 10.006배이지만 대략 10배입니다.) 이유는 초기에 4BarSim을 실행 시 입력되는 초기값을 문제에서 주어진 값의 10배 곱하여 대입했기 때문입니다.이를 고려해 볼 때, MATLAB을 통해 얻은 결과와 4Barsim을 이용해 얻은 결과는 일치한다고 할 수 있습니다.
고체역학 및 실습Project학부 : 기계공학부학번 : 200921590목차1문제 설명2하중 선도, 전단력 선도, 굽힘 모멘트 선도, 처짐선도(1) 불연속 함수방법(Discontinuity Function)(2) 2차 미분방정식 법(3) MATLAB을 이용하여 선도 표현3최대처짐(Max Deflection) & 최대 유연응력(Max Flexural Stress)4설계의 적합성 판단 & 대안 제시5고 찰6MATLAB m 파일 설명1. 문제 설명(1) 문제에서 수정 사항은 “W 460 × 110"이 아닌 ”W 460 × 113"입니다.(2) 문제에서 고려해야 할 사항은 보의 프레임(보)의 하중 또한 고려해야 한다는 것입니다.(3) 설계의 적합성을 판단하기 위해서S _{design`} ``=` {M _{max}} over {sigma _{allow}}을 사용합니다.-S`` GEQ `S _{design}로 단면을 선택하면, 굽힘응력의 크기는 보의 어느 곳에서도sigma _{allow}를 초과하지 않는 것으로 보장됩니다.2. 하중 선도, 전단력 선도, 굽힘 모멘트 선도, 처짐선도(1) 불연속 함수방법※ 잘 안보이는 부분이 있을 수도 있으므로, 원본 파일을 복사하여 첨부하겠습니다.2. 하중 선도, 전단력 선도, 굽힘 모멘트 선도, 처짐선도(1) 2차 미분방정식 법※ 잘 안보이는 부분이 있을 수도 있으므로, 원본 파일을 복사하여 첨부하겠습니다.2. 하중 선도, 전단력 선도, 굽힘 모멘트 선도, 처짐선도(3) MATLAB을 이용하여 선도 표현기본적으로 Macaulay Functions.을 MATLAB를 이용하여 나타내면 다음과 같습니다. ^{0} ^{1} ^{2} ^{3} ^{3}1) 하중선도① 하중은 아랫방향으로 작용하므로 그래프에서 하중의 방향을 고려해서 음수로 표현됩니다.보의 하중w _{b} , 트레일러 무게에 의한 하중w _{0}입니다.②0 ^{+} `` LEQ `x` LEQ 5.5 ^{-} ``m` 까지의 하중은 보의 하중과 트레일러 무게에 의한 하중의 합으로 크기_{b} )L``-`w _{b} ````(L LEQ x LEQ L+a)※Ay``=`9154.513487`N,``By=9759.166213`N 입니다.3) 굽힘 모멘트 선도< 굽힘 모멘트 선도 >{dM} over {dx} ``=`V(x)을 이용하여 굽힘 모멘트M(x)를 나타낼 수 있습니다.V(x)``=`Ay``-`(w _{0} `+w _{b} )x````(0 LEQ x LEQ L) 범위에서M _{max}를 가지는 것을 그래프를 통해 알 수 있습니다.M _{max}값을 가질 때 위치를 알기 위해M``'(x)```=`V(x)``=`0인 점을 찾으면x _{m} `=`2.742452993`m 임을 알 수 있습니다.따라서M _{max} ``=`M(x _{m} )``=`12552.91335`N` BULLET m 입니다.4) 처짐 선도< 처짐 선도 >M`(x)를 두 번 적분하게 되면EIv(x)를 나타내며, 양변을 EI로 나눠주면 다음과 같은 식을 얻을 수 있습니다.v(x)`=` {1} over {EI} `[ {Ay} over {6} x ^{3} ``-` {w _{0} ``+`w _{b}} over {24} x ^{4} ``+`EIv prime (0)](0` LEQ `x` LEQ `L)v(x)`=` {1} over {EI} `[ {Ay} over {6} x ^{3} ``-` {w _{0} ``+`w _{b}} over {24} [x ^{4} -(x-L) ^{4} ]`+` {By} over {6} (x-L) ^{3}#````````````````````````- {w _{b}} over {24} (x-L) ^{4} ``+`EIv prime (0)x]```(L` LEQ `x` LEQ `L+a) ※EI`v prime (0)``=`-23013.49322 입니다.또한 문제에서 주어진 재질 ASTM A-36의 W 460 × 113의 보의 탄성계수는 관성모멘트는 다음과 같습니다.E``=`200GPa``=`200 TIMES 10 ^{9} ``N/m ^{2} ,I``=556` TIMES 10)이S _{design}보다 크다면 이는 올바른 설계라 할 수 있습니다.S _{min} ``=` {I} over {C(보의`제일`윗부분)} ``=` {0.0000556`m ^{4}} over {0.2315`m} ``=`0.000240172`m ^{3}이므로 이 값은S _{design}보다 크므로 적합한 설계라고 할 수 있습니다. 또한sigma _{allow} 또한 탄성변형 범위를 유지하기 위해sigma _{YS} ``=250MPa를 사용하였습니다. 그리고 문제에서 구한 최대유연응력은 52.26617699 MPa이므로, 250MPa보다 작으므로 이를 통해서도 설계는 적합하다고 판단할 수 있습니다. 따라서 이에 대한은 별도로 없습니다. 추가적으로 경제적인 요소를 고려하여 8개의 보가 아닌 보의 숫자를 1개 정도 줄여도 적합 할 것이라고 판단됩니다.5. 고 찰이번 고체역학 및 실습 프로젝트는 1학기 때 MATLB을 이용하여서, 2차 미분방정식 방법과 불연속함수 방법을 사용하여 문제를 풀고, 이에 대한 설계의 적합성 여부를 판단하는 프로젝트이었습니다. 두 방법 다 많은 시간이 소요 되었지만, MATLAB를 이용함에 있어서 불연속함수 방법을 사용할 경우 쉽게 그래프를 얻을 수 있었으며, 최대 굽힘 모멘트, 최대 처짐을 손쉽게 알 수 있었습니다. 간단한 명령어 UnitStep 등을 잘 이용하면 더 어려운 문제도 간단히 해결할 수 있을 것이라고 생각됩니다. 이 문제를 해결하는 도중 적합성을 판단 할때, 최대 유연응력이 허용응력(항복응력)보다 한참 작음을 알 수 있었습니다. 설계를 함에 있어서 경제적인 요소도 중요하다고 배운바 있다 보니, 보를 8개로 지탱하지 않고 7개나, 6개 즉 보를 조금 줄여서 설계하더라도 적합한 설계가 될 것이라고 예상해 볼 수 있었습니다.6. MATLAB m 파일 설명(1) 하중선도TL = 6; % 주어진 L+a의 길이L = 5.5; % 주어진 L의 길이a = 0.5; % 주어진 a의 길이W = 10000; %주어진 트레일러의 무게 10tonw =UnitRamp, 제곱은 Unitquad, 세제곱은 UnitCubic를 사용하며,% ( x(i),0 )형태에서 뒤에 0은 경계값을 의미합니다. 즉 의 0을 의미합니다.% 이를 이용하여 위에 지정한 TL, L, w 등을 이용하여 하중선도 식을 나타내었습니다.% 또한 하중의 방향이 아랫방향이므로 -1을 곱하여 나타내었습니다.figure(1), plot(x,DistLoad,'linewidth',1.5) % 이를 그래프로 나타내기 위해 plot 명령어를 사용하였으며, 선의 굵기는 1.5로 지정했습니다.)title('Distributed Loading') % 제목은 Bending Momnet로 지정했습니다.xlabel('x [m]'); ylabel('Distributed Loading [N/m]') % x축은 단위가 m인 x축, y축은 N/m 단위의 분포하중을 의미합니다grid on % 그래프에서 격자를 사용하며,axis([-0.1 6.1 -3500 10]) % 그래프 범위는 x축으로 -0.1부터 6.1, y축은 -3500부터 10까지 나타내었습니다.(2) 전단력선도Ay=9154.513487; By=9759.166213; w=2229.5454; wb=1108.53;% Ay : A점에서의 y방향 반력, By : B점에서의 y방향 반력% w : 트레일러 하중(10t)에 의한 분포하중 (10(ton) / 8 / L * 중력가속도) = 2229.5454N/m% wb : W 460 * 113 보의 무게에 의한 분포하중 ( 113 * 중력가속도 ) = 1108.53 N/mTL = 6; % 주어진 L+a의 길이L = 5.5; % 주어진 L의 길이a = 0.5; % 주어진 a의 길이x=0:TL/1000:TL; % 간격을 TL/1000으로 나누므로, 총 데이터의 갯수는 1000 + 1 = 1001입니다.x=x';NumData=length(x); % 총 데이터의 갯수 1001입니다.DistLoad=[]; % DistLoad를 통해서 함수를 지정해줍니다.for i=1:NumData,% f정했습니다.xlabel('x [m]'); ylabel('ShearForce [N/m^2]') % x축은 단위가 m인 x축, y축은 N/m^2 단위의 전단력를 의미합니다grid on % 그래프에서 격자를 사용하며,axis([-0.1 6.1 -10000 10000]) % 그래프 범위는 x축으로 -0.1부터 6.1, y축은 -10000부터 10000까지 나타내었습니다.(3) 굽힘모멘트선도Ay=9154.513487; By=9759.166213; w=2229.5454; wb=1108.53;% Ay : A점에서의 y방향 반력, By : B점에서의 y방향 반력% w : 트레일러 하중(10t)에 의한 분포하중 (10(ton) / 8 / L * 중력가속도) = 2229.5454N/m% wb : W 460 * 113 보의 무게에 의한 분포하중 ( 113 * 중력가속도 ) = 1108.53 N/mTL = 6; % 주어진 L+a의 길이L = 5.5; % 주어진 L의 길이a = 0.5; % 주어진 a의 길이x=0:TL/1000:TL; % 간격을 TL/1000으로 나누므로, 총 데이터의 갯수는 1000 + 1 = 1001입니다.x=x';NumData=length(x); % 총 데이터의 갯수 1001입니다.DistLoad=[]; % DistLoad를 통해서 함수를 지정해줍니다.for i=1:NumData, % for 구문을 이용하여 첫번째 데이터 부터 지정된 1001개의 데이터까지 반복하는 작업을 합니다.DistLoad=[DistLoad; (Ay)*UnitRamp(x(i),0)-((wb+w)/2)*(UnitQuad(x(i),0)-UnitQuad(x(i),L)) + (By)*UnitRamp(x(i),L) - (wb/2)*(UnitQuad(x(i),L)-UnitRamp(x(i),TL))];end% 의 0승은 UnitStep로, 1승은 UnitRamp, 제곱은 Unitquad, 세제곱은 UnitCubic를 사용하며,% ( x(i),0 )형태에서 뒤에 0은 경계값을 의미합니다. 즉 의 0이 경반력