2015년 5월 17일 일요일

Julia Code : COMBI REDUCER with Planetary & Harmonic Drive v05


MIT에서 만든 Julia Language에 대해서 관심을 계속 가지고 있었는데, 이제 v0.3 이상으로 올라가면서 꽤 쓸만해 진 모양이다.  특히 아마존 웹 서비스에다가 이걸 올려서 클라우드 서비스까지 해 주기 때문에 어렵게 PC에 개발환경 구축하느라 진땀 빼지 않아도 되어서 더 좋다.

나에게 필요한 것은, 컴퓨터공학적인 최적화라던가 대규모 데이타를 효과적으로 다루는 것이 아니고, 높은 수준의 공학용 계산기 내지는 간단한 스크립트 엔진이다.

여러가지 많이 검토해 봤는데 Julia가 역시 최신기술의 집약체이다 보니 모든 조건에 다 부합하는 것 같다.

아래는 어떤 문제를 풀기 위해 Julia를 실무적으로 사용해 본 예이다.
원래는 Matlab m-code로 작성해서 FreeMat, Octave 따위로 돌려보던 것인데, Julia로 테스트해 보기 위해 컨버젼도 해 보고 완전히 새로 작성도 해 보았다.

m-code를 Julia code로 자동 컨버젼해주는 패키지도 Julia 패키지 라이브러리에 있으나, 일단 직접 수동으로 컨버젼 해 보았더니 함수 이름 바꿔준다던가 정도의 문법 수정으로 20분만에 끝났다.

그리고 Julia에 좀 더 부합하도록 몇가지 기능을 추가하고 수식도 새로 적용하고 해서 새로 만들었다.


---------------------------------------------------------------------------------------------------------------


본 Julia Code는 기존의 버전을 Julia 문법에 좀 더 부합하게 고치고, Interactive 기능을 추가한 버전이다.

최외곽 서큘러 스플라인 잇수(Number of Teeth)와 그 안쪽으로 들어가면서 조립되는 각 기어의 잇수의 차이를 기준으로 사양이 자동으로 잡히도록 하였고, Sun Gear의 반지름이 음수로 나오면 색상이 붉은 색으로 변하면서 경고 문구가 표시되도록 하였다.

기타 파라미터들은 관계를 조정하여 기하학적인 모순 상황이 발생하지 않도록 하였다.

타원 변형된 Flex Spline의 원주길이는 다음 공식을 기본으로 하였다.


* Annotations

$a$ : 긴 축의 지름 (Known)

$x$ : 짧은 축의 지름 (Unknown)

$\varepsilon$ : 이심률 (Known)

$L$ : 타원의 원주길이 (Known)


* Notations

$\varepsilon = \sqrt{ 1 - \frac{x^2}{a^2} }$

$L = 2aE(\varepsilon) = \pi a \left[  1 - \left( \frac{1}{2} \right)^2 \varepsilon^2 -  \left( \frac{1 \cdot 3}{2 \cdot 4} \right)^2 \frac{\varepsilon^2}{3} - \left( \frac{1 \cdot 3 \cdot 5}{2 \cdot 4 \cdot 6} \right)^2 \frac{\varepsilon^2}{5} - ... \right]$


* Remarks

$E(\varepsilon)$은 제2종 타원정적분 연산기호이다.  타원정적분은 엄밀하게 풀이하기가 어렵기 때문에, 이를 풀어쓴 우변의 무한급수 식을 이용하기로 한다.  급수의 항이 늘어날 수록 해(Solution)의 개수가 늘어나기 때문에 적합한 해를 판별하는 판단이 더 필요하게 되므로, 간소하게 하기 위해 다음과 같이 첫 번째 항만 살린 근사식을 만들었다.

$L = \pi a \left[  1 - \left( \frac{1}{2} \right)^2 \varepsilon^2 \right]$

해를 신속하게 구하기 위해 SageMath를 이용하였다.
(SageMath 대신 WolframAlpha를 이용해도 되나, 해가 복소수로 나오는 등 형태가 다를 수도 있다.)


a,L = var('a,L')
eq = L == pi*a* (  1-((1/2)^2)*(( sqrt(1-(x^2)/(a^2)) )^2)   )
eq
print solve(eq,x)


SageMath를 이용하여 얻어진 해는 다음과 같다.


L == 1/4*pi*a*(x^2/a^2 + 3)
[
x == -sqrt(-3*a^2 + 4*L*a/pi),
x == sqrt(-3*a^2 + 4*L*a/pi)
]



2개의 해가 얻어졌는데, 이 중에서 양수에 해당하는 두 번째 것을 취한다.  이 식을 Julia Code에 집어넣고 사용한다.


* Julia Coding

PyPlot 패키지는 Python 라이브러리 중에서, Matlab과 거의 비슷하게 Plot을 해 주는 MatPlotLib 기능을 Julia에서 사용할 수 있도록 해 주는 라이브러리이다.  따라서 이것은 내부적으로는 Python MatPlotLib가 그대로 동작하며, Julia Code에서 좀 더 사용하기 쉽도록 Matlab 문법과 거의 유사하게 간소화되어 있다.  기계공학도에게 Julia가 Python보다 접근하기 쉬운 점은, 굳이 객체지향 방법으로 코딩하지 않아도 된다는 점이다.  Matlab을 다루듯이 편하게 쓰면 된다.

Interact 패키지는 매우 간단하게 기본적인 GUI를 만들 수 있도록 해 준다.  여기 포함된 @manipulate 매크로를 이용하여, 자유롭게 조절하고 싶은 파라미터들을 설정해 주기만 하면 된다.  내부적으로는 무한루프를 돌리는 간단한 형태다.



#######################################################
# COMBI REDUCER with Planetary & Harmonic Drive
# V05
# 2015.05.18
# Dymaxion.kim
#######################################################
 
using PyPlot
using Interact
 
# Basic Parameter
alpha=[0:0.02:2*pi]
max_n_cs = 300
max_dn_fs = 20
max_dn_ig = 20
max_dn_sg = max_n_cs - max_dn_fs - max_dn_ig
# Making New Figure
f=figure(figsize=(8,8))
# Interact Macro
@manipulate for m = 0.1:0.1:2, # Module of The Gear
    n_cs = 50:1:max_n_cs, # Number of Teeth for Circular Spline
    dn_fs = 1:1:max_dn_fs, # Different Number of Teeth for Flex Spline
    dn_ig = 1:1:max_dn_ig, # Different Number of Teeth for Internal Gear
    dn_sg = 1:1:max_dn_sg # Different Number of Teeth for Sun Gear
    
    withfig(f) do
        grid("on")
        axis("auto")
        title("COMBI REDUCER v05")
        
        # Draw Circular Spline
        r_cs = m*n_cs/2
        plot( r_cs*sin(alpha), r_cs*cos(alpha), linewidth=2, color="Magenta" )
        annotate(["n_cs",int(n_cs)],xy=(0,r_cs),ha="center",va="bottom", color="magenta")
        
        # Draw Flex Spline in Idle Status
        r_fs = m*(n_cs-dn_fs)/2
        plot( r_fs*sin(alpha), r_fs*cos(alpha), linestyle="--", color="cyan" )

        # Draw Internal Gear in Idle Status
        r_ig = m*(n_cs-dn_fs-dn_ig)/2
        plot( r_ig*sin(alpha), r_ig*cos(alpha), linestyle="--", color="cyan" )
        
        # Draw Sun Gear
        r_sg = m*(n_cs-dn_fs-dn_ig-dn_sg)/2
        n_sg = 2*r_sg/m
        if r_sg > 0
            plot( r_sg*sin(alpha), r_sg*cos(alpha), linewidth=2, color="blue" )
            annotate(["n_sg",int(n_sg)],xy=(0,0),ha="center",va="center", color="blue")
        else
            plot( r_sg*sin(alpha), r_sg*cos(alpha), linewidth=4, linestyle="--", color="red" )
            annotate(["n_sg is in minus",int(n_sg)],xy=(0,0),ha="center",va="center", color="red")
        end
        
        # Draw Planet Gears
        t_fs = r_fs-r_ig # Thickness of Flex Spline's Plastic
        r_pg = (r_cs-t_fs-r_sg)/2
        plot( r_pg*sin(alpha), r_pg*cos(alpha)+(r_sg+r_pg), linewidth=2, color="magenta" )
        plot( r_pg*sin(alpha), r_pg*cos(alpha)-(r_sg+r_pg), linewidth=2, color="magenta" )
        n_pg = 2*r_pg/m
        annotate(["n_pg",int(n_pg)],xy=(0,r_sg+r_pg),ha="center",va="center", color="magenta")
        annotate(["Thickness",t_fs],xy=(0,r_sg+2*r_pg),ha="center",va="top", color="blue")
        
        # Draw Flex Spline in assembled
        L_fs = 2*pi*r_fs
        Dl_fs = 2*r_cs
        Ds_fs = sqrt(-3*Dl_fs^2 + 4*L_fs*Dl_fs/pi)
        plot( (Ds_fs/2)*sin(alpha), (Dl_fs/2)*cos(alpha), linewidth=2, color="blue" )
        n_fs = 2*r_fs/m
        annotate(["n_fs",int(n_fs)],xy=(-0.9*r_fs,0),ha="left",va="bottom", color="blue")

        # Draw Internal Gear in assembled
        L_ig = 2*pi*r_ig
        Dl_ig = Dl_fs - 2*t_fs
        #Ds_ig = sqrt(-3*Dl_ig^2 + 4*L_ig*Dl_ig/pi)
        Ds_ig = Ds_fs - 2*t_fs
        plot( (Ds_ig/2)*sin(alpha), (Dl_ig/2)*cos(alpha), linewidth=2, color="blue" )
        n_ig = 2*r_ig/m
        annotate(["n_ig",int(n_ig)],xy=(-0.9*r_ig,0),ha="left",va="top",color="blue")
        
        # Reduction Ratio
        RRp = (n_sg+n_ig)/n_sg
        RRh = n_fs/(n_fs-n_cs)
        RRt = RRp * RRh
        annotate(["Total Reduction Ratio",RRt],xy=(0,-r_cs),ha="center",va="top", color="black")

    end
end




아래 그림은 http://juliabox.org에서 구동한 결과이다.

검증결과 동일한 조건의 파라미터를 입력했을 때 기존 버전과 동일한 결과를 얻을 수 있었다.





댓글 없음:

댓글 쓰기