1. 서    론

아스팔트 콘크리트 포장에서 발생하는 피로 균열은 도로의 수명을 단축 시키는 주요 요인이다. 주로 기존 포장 도로에 존재하는 접합부 또는 균열에서 발생하는 응력 집중 및 이물질의 침투로 인하여 새로운 균열이 전파되며, 상호 연관적 균열(reflective crack)은 이의 대표적인 예이다. 아스팔트에서의 피로 균열에 영향을 미치는 주요 요인은 이동 하중의 특성, 온도 변화, 재료 특성 등으로 매우 복잡한 응답 특성을 보인다. 특히, 아스팔트 혼합물은 선형 탄성 이론을 이용하여 해석할 경우 한계가 있으며 점탄성적 재료 특성을 보이며 거동한다(Jo et al. 2006; Seo et al. 2008; Kwon 2014).

피로 수명을 예측하기 위하여 이용되는 현존 이론은 대부분 현장 또는 실험실에서의 관찰에 의존한 경험식에 바탕을 두고 있다. 특히, 피로와 관련하여 파괴 역학(fracture mechanics) 분야에서 널리 이용되는 물리량은 균열 끝단에서의 응력 확대 계수(stress intensity factor, SIF)(Sobczyk and Spencer 1992)와 에너지 해방율(energy release rate, ERR)(Kuai et al. 2009)이다. 따라서, 이들 파괴 역학적 물리량들에 대한 계산의 정확도를 높이는 것이 피로 균열 성장의 해석에 있어 매우 중요하다. 이들을 계산하기 위한 방법으로 유한 요소법(finite element method, FEM)이 널리 이용되고 있으며, 이때 국부적으로 초미세 유한 요소망(heavily refined FE mesh)을 사용하여 균열 끝단 주변에서 해에 대한 정확도를 높일 수 있다(Seo et al. 2014).

본 연구에서는 유한 요소법의 발전된 형태로서 최근 균열 성장을 모델링하기 위하여 널리 이용되고 있는 일반 유한 요소법(generalized FEM) (Babuska et al. 2004)을 이용하여 균열이 포함된 삼차원 경계값 문제를 해석하였다. 또한 균열 끝단 주변에서 상대적으로 작은 크기의 요소를 사용하여 국부적으로 해의 해상도를 높였다. 이때 각각의 균열 성장 단계에 따라 균열 끝단의 위치가 바뀌며, 따라서 상응하는 적응적 유한 요소망(adaptive FE mesh)도 지속적으로 달라진다. 그러나, 일반 유한 요소법에서 균열면은 반드시 체적 요소의 경계면에 일치할 필요가 없기 때문에 적응적 유한 요소망을 균열면의 기하학적 형상에 독립적으로 생성할 수 있다. 따라서 체적 요소 사이의 경계면에서 소위 응집 요소(cohesive element)를 이용하여 불연속을 표현하는 기존의 접근법(Turon et al. 2007)에 비하여 각각의 균열 성장 단계에서 매번 수행해야 하는 유한 요소망의 생성과정이 용이하다는 장점이 있다.

2. 점탄성 피로 균열 성장식

재료는 정적 파괴 강도보다 비교적 작은 크기의 반복 하중 하에 놓인 경우라 하더라도 점진적인 손상을 겪으며 이로 인하여 피로 파괴에 도달한다. 피로 균열 성장은 균열 끝단에서의 응력 확대 계수와 관계가 있음이 실험적으로 보고되었다. 이 중 가장 널리 쓰이는 식은 Paris-Erdogan 식(Anderson 2004)이며, 이 식은 각각의 균열 단계에서의 피로 균열 성장을 지수 법칙을 이용하여 표현한다.

 (1)

여기서, 은 하중 반복 횟수당 균열 성장량, 는 한개의 하중 사이클에 대한 응력 확대 계수의 변화량이다. 식 (1)은 로그-로그 스캐일에서 이 에 선형적으로 변화한다는 것을 의미하며, 이때 y 절편은 이고 기울기는 이다.

식 (1)에서 와 은 재료 실험 결과로부터 결정하여야 할 재료의 물성치로서 하중의 크기에는 독립적이나, 일반적으로 재료, 환경, 온도, 하중의 진동수 등에 의존적이다. 또한 응력 집중 계수는 하중의 진동수의 변화에 영향을 받지 않기 때문에 탄성적인 거동을 보이는 재료에서의 피로 균열 성장과 관련이 있다. 따라서 아스팔트 콘크리트에서의 점탄성 효과를 반영하기 위해서 응력 확대 계수를 에너지 해방율로 대체하여 다음과 같은 수정된 식이 이용되고 있다(Kuai et al. 2009).

 (2)

여기서, 는 한 개의 하중 사이클에 대한 에너지 해방율의 변화량이다. 따라서 피로 수명을 정확하게 예측하기 위해서는 균열 끝단에서 를 정밀하게 계산하여야 하며, 는 분석적 또는 수치적으로 계산할 수 있다. 이 연구에서는 분석적인 해를 구하기 힘든 일반적인 삼차원 조건 하에서의 문제를 다루므로 일반 유한 요소법을 이용한 수치 해석적인 방법을 통하여 를 구하였다.

3. 수치 해석 기법

3.1 균열을 포함한 선형 탄성 문제의 순차 해석

점탄성 재료에서의 피로 균열 성장을 모델링하기 위해서는 일반적으로 무수히 많은 반복 하중에 대하여 비선형 시간 이력 해석을 수행하여야 하며, 이는 일반적인 삼차원 조건 하에서 많은 계산량이 요구된다(Zocher et al. 1997). 이러한 문제를 해결하기 위하여 본 연구에서는 주어진 경계값 문제를 여러 개의 균열 성장 단계로 분리한 후, 각각의 단계에서 균열 끝단 부근에서만 상대적으로 크기가 작은 요소를 가지는 적응적 유한 요소망을 사용하여 해석을 수행하였다.

각각의 균열 성장 단계는 균열의 진행 속도가 느리다는 가정 하에서 고정된 크기의 균열을 포함하는 준정적 문제로 정의할 수 있다. 또한, 이러한 준정적 문제는 선형 탄성 파괴 역학(linear elastic fracture mechanics)의 조건 하에 해석하였다. 따라서, 뉴턴-랩슨 반복법 등과 같은 반복 계산을 요하는 방법을 적용할 필요가 없으며, 또한, 적응적 유한 요소망을 사용할 시 서로 다른 유한 요소망 사이에 비선형 해를 사상(mapping)해줘야 하는 문제도 피할 수 있다. 결국 이를 통하여 계산량을 줄일 수 있을 뿐만 아니라, 사상을 위하여 복잡합 계산 알고리즘을 적용할 필요가 없다는 장점이 있다.

3.2 피로 수명 계산

3.2.1 에너지 해방율 계산: 탄성-점탄성 대응원리

점탄성적 응답 특성을 보이는 재료에서의 균열 성장량과 그에 상응하는 하중 반복 횟수를 계산하기 위해서 각각의 균열 성장 단계에서 식 (2)를 이용하였다. 이때 식 (2)를 이용하기 위하여 를 계산하여야 하며, 는 점탄성적인 물리량이다. 본 연구에서는 탄성-점탄성 대응원리(elastic-viscoelastic correspondence principle)(Schapery 1984)에 근거하여 라플라스 영역에서 계산한 탄성해를 시간 영역에서의 점탄성해로 변환하는 과정을 통하여 를 계산한다.

라플라스 영역에서의 에너지 해방율 는 다음과 같이 표현된다.

 (3)

여기서, 모든 변수들은 라플라스 영역에서 정의되며, 는 변환변수, 는 균열의 크기, 는 응력 확대 계수, , , 는 각각 전단 릴렉세이션함수, 체적 릴렉세이션함수, 작용하중의 크기를 반영하여 응력 확대 계수를 조정하는 함수이다. 라플라스 영역에서 정의된 식 (3)에 대응하는 시간 영역에서의 에너지 해방율 는 다음과 같이 표현된다.

 (4)

여기서, 모든 변수들은 시간() 영역에서 정의되며, 는 라플라스 영역에서의 에 상응한다. 또한, , 는 각각 라플

라스 영역에서 와

의 관계를 가지는 컴플라이언스 함수들이다. 이때, 식

(3)에서 식 (4)로의 라플라스 역변환은 일반적인 삼차원 조건 하에서의 복잡한 경계조건과 재료 함수식을 고려하여 푸리에 급수 방법(Dubner and Abate 1968)을 이용하여 수치적으로 수행하였다.

3.2.2 에너지 해방율 변화량의 안정화

식 (2)를 이용하기 위하여 주어진 균열 크기에서 에너지 해방율의 변화량은 고정된 값이어야 한다. 일반적으로 반복 하중 하에서 에너지 해방율은 시간에 따라 변화하지만, 에너지 해방율의 변화량은 급속히 안정화되어 고정된 값을 가진다는 것이 관찰되었다(Kuai et al. 2009). 따라서 이러한 안정화 된 값을 식 (2)에서 사용할 수 있다. 이때 각각의 균열 성장 단계에서의 해석 시간은 에너지 해방율이 안정화되는 시간보다 길어야 하며, 이러한 조건은 균열의 성장 속도가 느린 경우에 만족된다.

3.2.3 균열 성장의 기하학적 표현

이 연구에서는 성장하는 균열면을 삼각형으로 구성된 요소망을 이용하여 명시적으로 표현하였다(Fig. 2 참조). 이러한 면 요소망은 체적 요소로 구성된 유한 요소망에 독립적이기 때문에 균열면의 기하학적 표현에 대한 정확도는 체적 요소의 크기에 영향을 받지 않는다. 따라서, 일반적인 삼차원 조건 하에서 균열면의 기하학적 형상을 정교하게 표현할 수 있다는 장점이 있다.

3.2.4 균열 성장량에 따른 하중 반복 횟수 계산

에너지 해방율은 일반적인 삼차원 조건 하에서 파괴 모드 I, II, III 대하여 각각 정의되므로, 모든 모드에 대한 효과를 반영하기 위하여 한 개의 등가 에너지 해방율 을 정의할 수 있다. 에너지 해방율의 변화량을 이용하여 균열의 성장량 뿐만 아니라 진행 방향도 계산할 수 있으며, 이를 위해 를 정의하고 이용하는 방법은 다양하다. 이 연구에서는 Schollmann’s criterion (Schollmann et al. 2002)을 이용하였다.

식 (2)를 이용하여 균열 끝단의 번재 절점에서 균열의 성장량을 다음과 같이 계산할 수 있다.

 (5)

여기서, 는 사용자에 의해 정의되는 값인 최대 균열 성장량, 는 번재 절점에서의 등가 에너지 해방율의 변화량, 은 균열 끝단의 절점을 따라서 의 최대값이다. 식 (5)에 대응하는 하중 반복 횟수는 다음과 같이 표현된다.

 (6)

여기서, 는 균열 성장 단계에 대한 색인이다.

3.3 일반 유한 요소법: 균열을 포함한 삼차원 선형탄성 문제의 해석

일반 유한 요소법은 유한 요소법의 발전된 형태로서, 파괴 역학 이론으로부터 알려진 분석적 해 또는 수치적으로 계산한 해에 대한 정보를 이용하여 기존의 유한 요소법의 해 공간을 확장하는 기법이다. 일반 유한 요소법의 형상 함수(shape function) 는 단위 오목 분할 함수(partition of unity) 에 확장 함수(enrichment function) 를 곱하여 다음과 같이 정의한다.

 (no summation on )        (7)

여기서, 는 유한 요소망에서의 절점에 대한 색인, 는 각각의 절점에서 사용하는 확장 함수에 대한 색인이 이다. 확장 함수는 주어진 문제의 해의 특성을 반영한다. 이 연구에서는 균열에서 먼 영역에서는 고차 다항식, 균열면과 균열 끝단에서는 각각 불연속 및 특이 함수를 사용하였다(Belytschko and Black 1999). 즉, 확장 함수는 각각의 절점 에서 균열의 기하학적 형상에 의존하여 활성화된다. 따라서 균열은 유한 요소망이 아닌 확장 함수에 의해 표현되기 때문에 균열의 전파는 유한 요소망에 독립적이다.

4. 수치 해석 예제

4.1 원판 시편 조밀 인장 시험

수치 해석 예제로서 원형 단순 인장 시험편 시험(disc-shaped compact tension test)을 고려하였다. 이 실험의 절차는 ASTM D7313-07에 규정되어 있으며, 구체적인 실험 세팅은 Fig. 1에 나타내었다. Kuai et al.(2009)를 참조하여 시편의 치수 및 하중 조건을 결정하였으며, 하중 함수는 다음과 같다.

 (8)

여기서, 와 는 각각 하중의 크기 및 진동수이다.

본 연구에서의 수치 해석 결과와 실험 결과를 비교하기 위하여 Kuai et al.(2009)의 실험에서 사용한 재료의 물성치와 동일한 값을 사용하였다. 다음과 같이 5개의 계수를 가지는 프로니 급수를 이용하여 표현하였다.

 (9)

여기서, 는 지연 시간, 와 는 프로니 급수의 계수이다. Kuai et al.(2009: 75)의 Fig. 3에는 점탄성 거동에 대한 재료의 물성치에 대하여 실험 결과로부터 얻은 크립 컴플라이언스 마스터 곡선이 보고되어 있으며, 이 곡선을 재현하도록 식 (9)에서의 계수들을 결정하였다. 또한, 피로 균열 성장을 지배하는 식 (2)의 계수 와 은 Kuai et al.(2009: 80)에서 적용한 값과 동일한 0.0015 및 1.3225를 사용하였다. 이들은 온도에 따라 변하는 값이며 이 연구에서는 온도 20 °C에 상응하는 값을 사용하였다. 그러나, 하중의 크기 및 진동수에 따라 변하지 않기 때문에 점탄성 재료의 물성치로 적합하다.

실제 실험에서는 뭉툭한 노치로부터 균열이 최초 성장하기 위하여 일정수의 반복 하중이 필요하다. 이러한 초기 반복 하중 횟수는 노치 끝단 및 시편의 기하학적 형상에 영향을 받으며, 이것을 전체 반복 하중 횟수에서 빼 줌으로써 시편 마다 서로 다른 초기 반복 하중 횟수의 영향을 배제하였다.

4.2 수치 해석 결과

수치 해석은 4절점 사면체 요소를 사용하여 수행 하였다. Fig. 2에서 보듯이 균열 끝단 주변에서 상대적으로 작은 크기의 요소를 사용하였으며 이에 상응하는 전체 자유도의 수는 371,382∼498,732개이다. 이와 같이 균열의 성장에 따라 계속 변하는 적응적 유한 요소망의 사용은 균열의 진전이 예상 되는 모든 경로 주변에서 상대적으로 작은 크기의 요소를 가지는 단일 유한 요소망을 사용하는 방법에 비하여 전체 자유도 수를 줄일 수 있다. 특히, 균열의 진전 방향을 예상할 수 없는 경우에 적응적 유한 요소망의 사용을 통한 자유도의 감소는 더욱 크다. 또한, 적응적 유한 요소망은 각각의 균열 성장 단계에서 자동적으로 생성하였다. 이는 3.2.3장에서 설명한 유한 요소망에 독립적으로 균열면을 표현하는 기법을 통하여 균열 끝단의 위치를 명시적으로 알 수 있기 때문에 가능하다. 각각의 균열 성장 단계에서 선형 해석을 수행하였으며, 따라서 서로 다른 유한 요소망 사이에 비선형 해를 사상해 줄 필요가 없다.

Fig. 2는 서로 다른 균열 성장 단계에서 하중이 작 용하는 방향으로의 법선(normal) 응력을 나타낸다. 적응적 유한 요소망을 사용하여 균열 끝단 근처에서 발 생하는 응력 집중을 잘 포착하는 것이 관찰된다. 특히, 이와 같이 균열 끝단 근처에서의 국부적인 거동에 대한 해상도를 높임으로써 파괴 역학적 물리량인 에너지 해방율을 정확하게 계산할 수 있으며, 이를 통하여 정확한 피로 수명의 예측이 가능하다.

Fig. 3은 Fig. 2에 상응하는 각각의 균열 성장 단계에서 진전된 균열면의 형상을 나타낸다. 균열 끝단은 시편의 위에서 봤을 때 균열 진행 방향으로 볼록한 형태를 띠는 것이 확인된다. 이 현상은 삼차원 해석을 통해서만 관찰할 수 있으며, 식 (5)의 적용 시 균열 끝단을 따라 절점에서의 값들이 서로 다르기 때문에 발생하는 자연스러운 결과이다.

Table 1은 각각의 균열 성장 단계에서 계산한 최대 등가 에너지 해방율 변화량 의 값을 나타낸다. 이들은 안정화된 후의 값으로서 식 (5)에서 사용되며, 각각의 균열 끝단 절점에서 균열의 성장량 및 방향을 결정하는 매우 중요한 물리량이다. 균열이 진전함에 따라 에너지 해방율의 값이 증가하며, 따라서 균열의 성장 속도도 빨라진다는 것을 알 수 있다.

Kuai et al.(2009: 76)의 Fig. 4에서는 피로 수명에 대한 실험 결과로서 반복 하중 횟수에 따른 균열 성장량을 나타내었으며, 이에 대한 스냅샷을 Fig. 4에 포함시켰다. 하중의 진동수와 온도를 각각 10 Hz와 20 °C로 고정한 상태에서, 하중의 크기를 200N, 250N, 300N로 변화시키면서 관찰한 결과이다. Fig. 5는 수치 해석 결과를 실험 결과와 비교한 그래프이다. 두 개의 서로 다른 접근법에 의한 결과가 비교적 일치하는 것을 알 수 있다. 특히, Table 2에서 보듯이 최종 피로 수명에 대한 상대 오차는 불과 0.3 % 이하이다. 이러한 결과는 서로 다른 하중 크기에 대하여 동일한 및 을 이용하여 얻은 결과이며, 따라서 점탄성 재료에 대한 피로 균열 성장 식 (2)와 재료 물성치로서의 및 의 사용을 정당화한다.

5. 결    론

아스팔트 콘크리트를 점탄성 재료로 모델링 한 후, 일반적인 삼차원 조건 하에서 피로 균열 성장을 수치 모사하였다. 효율적인 해석을 위하여 주어진 문제를 균열 성장 단계에 따라 여러 개의 준정적 문제로 분리한 후 적응적 유한 요소망을 이용하여 풀었다. 이때 각각의 삼차원 경계값 문제를 풀기 위하여 일반 유한 요소법을 이용하여 선형 탄성 해석을 수행하였다. 일반 유한 요소법에서 균열면은 체적 요소를 관통할 수 있다는 장점이 있으며, 특히 본 연구에서는 균열면의 기하학적 형상을 유한 요소망에 독립적이고 명시적으로 표현하였다. 점탄성 재료의 피로 균열의 진행 속도 및 방향은 에너지 해방율의 변화량에 의해 결정된다. 점탄성적인 물리량인 에너지 해방율을 계산하기 위해 라플라스 영역에서의 탄성해를 시간 영역에서의 점탄성해로 변환시키는 탄성-점탄성 대응 원리를 이용하였다. 이러한 전략을 통하여 비선형 시간 이력 해석을 수행하는 경우에 비하여 계산량을 줄일 수 있을 뿐만 아니라, 적응적 유한 요소망의 사용 시 필요한 과정인 비선형 해의 사상을 피할 수 있었다. 수치 해석 결과를 기존의 문헌에 보고된 실험 결과와 비교함으로써 기술한 수치 해석 기법의 적용 및 우수성을 예증하였다.