임의 형상 고정단 평판의 고정밀도 고유치 해석을 위한 파동 함수 기반 무요소법

Meshless Method Based on Wave-type Function for Accurate Eigenvalue Analysis of Arbitrarily Shaped, Clamped Plates

Abstract

The paper proposes a practical meshless method for the free vibration analysis of clamped plates having arbitrary shapes by extending the non-dimensional dynamic influence function (NDIF) method, which was developed by the author in 1999. In the proposed method, the domain and boundary of the plate of interest are discretized using only nodes without elements unlike FEM and the system matrices are obtained by making domain nodes and boundary nodes satisfy the governing differential equation and boundary conditions, respectively. However, since the above system matrices are not square ones, the problem of free vibrations of clamped plates is not reduced to an algebraic eigenvalue problem. An additional theoretical treatment is considered to produce an algebraic eigenvalue problem. It is revealed from case studies that the proposed method is valid and accurate.

1. 서 론

이 논문의 저자는 임의 형상 멤브레인과 다양한 경계조건을 가진 임의 형상 평판의 고정밀도 고유치 해석을 위해 NDIF법을 처음 제안하여 지금까지 많은 논문을 발표하였다(1~7). 파동함수를 기저함수로 사용하는 NDIF법은 경계요소법(8)과 마찬가지로 해석 대상물의 경계만을 이산화하나, 요소(element)를 사용하지 않고 노드(node)만으로 이산화한다는 면에서 경계요소법과는 차별성을 가진다. 이러한 이유로 인해, NDIF법은 이론 정식이 간결하여 해의 정밀도가 유한요소법이나 경계요소법에 비해 우수하다(8, 9).

그러나 NDIF법은 유한요소법 및 경계요소법과는 달리 최종 시스템 행렬이 주파수 파라미터에 종속되는 단점을 가진다. 결과적으로 최종 시스템행렬식은 대수 고유치 문제(algebraic eigenvalue problem)의 형태를 가지지 못하게 되어, 시스템 행렬로부터 고유치와 고유모드를 추출하는 과정이 매우 번거롭게 된다(1~7).

이 논문에서는 이러한 NDIF법의 단점을 보완할 수 있는, 파동함수 기반 무요소법을 개발하여 임의 형상 고정단 평판의 고유치 해석에 적용하였다. 제안된 정식의 타당성을 검증하기 위해 직사각형 고정단 평판에 대한 예제 연구를 수행하였다.

 

2. 이론 정식화

2.1 지배방정식과 경계조건

평판의 진동을 조화 운동으로 가정하면 운동 지배 방정식은 다음과 같다(10).

여기서 W 는 평판의 진동 변위, 𝚲 는 주파수 파라미터이다. 그리고 고정단 경계조건은 다음과 같이 두 개의 식으로 주어진다.

여기서 n은 경계의 법선 방향을 의미한다.

2.2 근사해 가정 및 시스템 행렬식 추출

Fig. 1과 같이 해석 대상 평판의 내부와 경계를 총 N개의 노드로 이산화 한다. 평판 내부 노드(domain node)와 평판 경계 노드(boundary node)의 개수를 각각 NI 와 NB로 가정한다( N = NI + NB ). 평판 내부의 변위를 다음과 같이 가정하자.

Fig. 1Discretization of a plate not using elements but using nodes only

여기서 r 은 Fig. 2와 같이 평판의 임의 점 P 에 대한 위치 벡터이며, rj 는 j 번째 노드 점에 대한 위치 벡터이다. 기저함수 ϕ(| r - rj |)는 j 번째 노드에서 정의된 파동 함수로, 두 점 사이의 거리 | r - rj |만의 함수이다. 그리고 Aj 는 j번째 파동함수의 기여도 계수를 나타낸다.

Fig. 2Position vectors for arbitrary point P and the j-th node Pj in the plate

가정된 변위 해 식(4)는 NI 개의 내부 노드에서 지배방정식 식 (1)을 만족해야 하므로, 식 (4)를 식 (1)에 대입한다. 그러면 다음의 식을 얻을 수 있다.

다음으로, 가정된 변위 해 식 (4)는 NB개의 경계 노드에서 경계조건 식 (2), 식 (3)을 만족해야 하므로, 식 (4)를 식 (2), 식 (3)에 대입한다. 그러면 다음의 식들을 각각 얻을 수 있다.

여기서 ni는 i번째 경계 노드의 법선 방향을 의미한다.

다음으로 식 (5)에서 식 (7)은 하나의 시스템 행렬식으로 다음과 같이 나타낼 수 있다.

여기서

위의 식들의 행렬 요소 값들은 다음 식들에 의해 주어진다.

식(8)을 살펴보면, SM1과 SM2의 크기가 N×(N+NB)로 미지수의 개수 보다 식의 개수가 더 많은 행렬식임을 알 수 있다. 결과적으로 식 (8)은 일반적인 대수 고유치 문제와 다른 형태를 가지기 때문에, 수치해석 상용프로그램에서 제공하는 고유치 솔버(solver)로 풀리지 않게 된다. 이러한 문제점을 해결하기 위한 방안이 2.3절에서 고려된다.

2.3 대수 고유치 문제로의 정식화

식 (4)에 NB개의 기여도 계수를 추가하는 방법에 의해, 진동 변위 해를 식 (17)과 같이 새로이 가정한다.

식(17)에서 ψ (| r - rj |)를 어떤 식으로 가정할지에 대한 여러 가지 방안이 있으나, 이 논문에서는 다음과 같이 원래 주어진 기저함수 ϕ(| r - rj |)를 경계 노드 법선 방향으로 미분한 형태로 가정한다.

여기서 nj는 j번째 경계 노드의 법선 방향을 의미한다.

식 (17)을 지배방정식 식 (1)과 경계조건 식 (2), 식 (3)에 적용하면, 다음과 같은 시스템행렬식을 얻을 수 있다.

여기서

위의 식들에서 새로이 정의된 행렬들의 요소 값들은 다음 식에 의해 주어진다.

식(19)에 포함된 시스템 행렬 SM3 과 SM4의 크기는 (N + NB ) × (N + NB )로 정사각 행렬이며, 식 (19)는 식 (8)과는 달리 대수 고유치 문제의 형태임을 알 수 있다. 수치해석 상용 프로그램의 대수 고유치 해석 솔버를 이용하기 위하여 식 (19)는 다음과 같이 변형한다.

여기서 I 는 단위 행렬이며,

식 (27)로부터 고유치와 고유벡터를 구할 수 있으며, 고유벡터를 식 (17)에 대입하면 모드 형상을 그릴 수 있다.

 

3. 예제 연구

이 논문에서 제안된 정식을 검증하기 위해 고정단 직사각형 평판(가로 1.2 m, 세로 0.9 m)에 대한 예제 연구를 수행하였다. 직사각형 평판을 사용한 이유는 노드 이산화를 하기가 용이하며, 저자가 기수행한 NDIF법 연구결과와 비교가 가능하기 때문이다. 참고로, 보다 다양한 형상을 가진 평판에 대한 예제 연구는 임의 형상에 대한 노드 자동 이산 프로그램이 개발되는 후속 연구에서 수행될 예정이다.

이 예제 연구에서는 식 (17)에 사용될 기저함수로, 파동함수 중 식 형태가 상대적으로 간결한 ϕ(| r - rj |) = e-cr2 이 사용되었으며(11,12), 여기서 r = |r - rj |이고 c는 기저 함수의 형상 파라미터이다. Fig. 3은 63개의 노드(가로 9개, 세로 7개)를 사용하여 등간격으로 영역을 이산화한 직사각형 평판을 보여준다. 이 평판에 대한 고유치 해석 결과는 Table 1에 요약되었다. Table 1에서 이 논문에서 제안된 방법(proposed method) 에 의해 구해진 고유치를 ANSYS(FEM) 해와 비교해보면, 제안된 방법이 단지 63개의 노드만을 사용했음에도 불구하고, 961개의 노드를 사용한 ANSYS 결과와 1.5 %이내의 오차를 가짐을 확인할 수 있다. 참고로 저자가 과거에 제안했던 NDIF법(NDIF method)은 ANSYS 결과와 0.1 % 이내 오차에 들 정도로 아주 정확한 결과를 제공함을 확인할 수 있다. 하지만, 서론부에서 설명한 바와 같이 NDIF법은 고유치 추출 방법의 번거로움이 존재하는 단점을 지닌다.

Fig. 3Equal-space discretization of the rectangular plate using 63 nodes

Table 1Eigenvalues of the clamped plate by the proposed method using 63 nodes (equal space; the parentheses values: percentage errors for the ANSYS results)

Fig. 4는 제안된 방법에 의해 구해진 1차에서 6차 모드 형상을 보여준다. Fig. 4의 모드 형상을 ANSYS에 의해 구해진 모드 형상과 비교해보면, 이 논문에서 제안된 방법이 성공적으로 모드를 추출할 수 있음을 확인할 수 있다.

Fig. 4Mode shapes of the clamped rectangular plate (63 nodes, c = 7)

한편, 제안된 방법에 사용된 기저 함수의 형상 파라미터 c값이 고유치에 어떠한 영향을 미치는지 확인하기 위하여 Table 2에서와 같이 c값의 변화에 따른 고유치 해석을 수행하였다. Table 2에서 제안된 방법과 ANSYS 결과를 비교해보면, 첫 번째 고유치는 c = 5.5일 때, 두 번째 고유치는 c = 5.0일 때, 세번째 고유치는 c = 5.5일 때 ANSYS 결과와 오차가 제일 적음을 확인할 수 있다. 이들 세 고유치에 대한 평균 오차(ave. error) 값을 계산하면, c = 5.5일때 평균 오차가 0.48 %로 제일 적음을 확인할 수 있다.

Table 2Eigenvalues of the clamped plate using 63 nodes as c varies (equal space; the parentheses values: percentage errors for the ANSYS results)

다음으로, 노드 이산화 방법이 고유치에 미치는 영향을 확인하기 위하여, 노드 이산화를 Fig. 5와 같이 임의 간격으로 하였다. 이 경우에 대한 고유치 해석 결과는 Table 3과 같으며, c = 4.0일 때 고유치 평균 오차가 0.46 %로 오차가 조금 줄어들었음을 확인할 수 있다. 참고로, c값을 너무 작게 설정할 경우 고유치가 수렴하지 못하는 경우가 발생하므로, 수렴되지 않는 결과를 제공하는 c값의 최댓값을 찾은 다음에, c값을 단계적으로 크게 하면서 고유치 해석을 수행하는 방법이 이 논문에서 이용되었다. 상기와 같이 c값에 따라 고유치가 추출이 안되는 경우가 발생하는데, 이러한 특성은 제안된 방법이 유한요소법에 비해 가지는 단점이라고 말할 수 있으며 이를 개선하기 위한 추가 연구가 필요하다.

Fig. 5Random-space discretization of the rectangular plate using 63 nodes

Table 3Eigenvalues of the clamped plate using 63 nodes as c varies (random space; the parentheses values: percentage errors for the ANSYS results)

한편, 노드의 개수가 고유치 해석 오차에 미치는 영향을 확인하기 위해, 노드의 개수를 130개(가로 13개, 세로 10개)로 증가시켰다. 노드는 모두 가로 세로 0.1 m의 등간격으로 이산화 되었다. 이 경우에 대한 해석 결과는 Table 4와 같다. Table 4에서 알 수 있듯이, c = 12.5일 때 ANSYS 결과와 가장 작은 평균 오차를 보이며, 평균 오차가 0.24 %로 63개의 노드를 사용한 Table 2와 Table 3의 결과에 비해 오차가 확연히 줄어들었음을 확인할 수 있다.

Table 4Eigenvalues of the clamped plate using 130 nodes as c varies (equal space; the parentheses values: percentage errors for the ANSYS results)

마지막으로 130개의 노드를 임의 간격으로 Fig. 6와 같이 이산화 하였다. 이 경우에 대한 고유치 해석 결과는 Table 5와 같다. Table 5를 살펴보면, c = 9.5일때 평균 오차가 0.20 %로 최소가 됨을 알 수 있으며, 130개의 노드를 균일하게 이산화 한 Table 4의 결과와 비교해볼 때 오차가 많이 감소했음을 확인할 수 있다.

Fig. 6Random-space discretization of the rectangular plate using 130 nodes

Table 5Eigenvalues of the clamped plate using 130 nodes as c varies (random space; the parentheses values: percentage errors for the ANSYS results)

 

4. 결 론

이 논문에서는 유한요소법이나 경계요소법과는 달리, 요소를 이용하지 않고 단지 노드로만 고정단 평판을 이산화하여 정확한 고유치 해석을 수행할 수 있는 무요소법에 대한 이론을 정식화하였다. 예제 연구를 통해, 기저 함수의 형상 파라미터가 고유치의 정확도에 민감한 영향을 미침을 확인하였으며, 노드의 개수를 증가 시키면서 임의 간격으로 노드 이산화를 하였을 때 가장 정확한 고유치가 얻어짐을 검증하였다. 향후 연구에서는 보다 다양한 임의 형상에 대한 예제 연구가 수행될 예정이며, 정확한 고유치를 제공할 수 있는 최적의 형상 파라미터 값을 찾는 방법에 대한 연구도 병행될 것이다. 또한, 제안된 방법의 고유치 추출 계산 시간을 최소화하기 위한 연구도 향후 연구에서 진행될 예정이다.