이전에 유체가 어떻게 움직이는지 알아보았다. 이제는 유체가 왜 움직이는지, 유체가 흐르는 원인이 되는 힘이 뭔지, 그러한 힘에 반응해서 어떻게 흐르는지, 힘과 흐름에 어떤 관계가 있는지를 알아볼 차례이다. 유체의 흐름에 대한 지배 방정식은 나비에-스토크스 방정식인데, 이 글의 목표는 나비에-스토크스 방정식까지 도달하는 것이다. 나비에-스토크스 방정식까지 도달한 이후에 몇 가지 내용들을 더 살펴보고 글을 마무리할 것이다.
모든 유체역학은 질량보존법칙, 운동량보존법칙, 에너지보존법칙에 기반한다. 이러한 보존방정식은 한 점에 적용 가능한 미분의 형태일 수도 있고, 어떤 영역에 적용 가능한 적분의 형태일 수도 있다. 적분의 형태는 또 다시 어떤 부피가 공간에서 고정된 영역일 수도 있고, 동일한 유체 입자들을 수송하는 움직이는 영역일 수도 있다. 전자를 control volume이라고 하고 후자를 material volume이라고 한다. 이 글에서 두 가지를 모두 다룰 것이다.
본격적으로 보존방정식을 유도하기 전에 알아야 하는 몇 가지 정리들이 있다.
미분의 형태와 적분의 형태는 아래의 divergenve theorem of Gauss에 의해 서로 변형될 수 있다.
$$\int _V \frac{\partial F}{\partial x_i} dV = \int _A F dA_i $$
사실 이 divergence theorem도 물리적 의미를 생각해보면 당연한 것이다. 어떤 체적을 둘러싸는 면을 지나는 어떤 물리량의 flux의 총합은 그 체적 내부에서의 변화량의 총합과 같다는 것이다.
그리고 어떤 함수 $F$를 부피에 대해서 적분하는데, 해당 부피가 material volume이라서 위치가 변하는 경우를 생각해보자. 이 적분을 시간에 대해서 미분하려면 $F$라는 함수 자체의 시간적 변화와 적분 영역 자체의 시간적 변화를 모두 생각해야 한다. 이에 대한 theorem이 아래의 Leibnitz rule이다.
$$\frac{d}{dt}\int_{x=a(t)}^{b(t)} F(x,t)dx = \int_{a}^{b} \frac{\partial F}{\partial t} dx + \frac{db}{dt} F(b,t) - \frac{da}{dt} F(a,t)$$
고3때 수능을 공부할 때 대치동 학원가에서 만들어낸 매우 어려운 난이도의 문제를 선생님이 풀이해주시면서 이러한 라이프니츠 rule을 소개해주셨던 기억이 난다.
아무튼 위 방정식의 우변에서 첫 번째 항은 $F$라는 함수 자체의 변화에 대한 것이고, 두 번째와 세 번째 항은 시간이 지남에 따른 적분 구간의 변화에 대한 것이다. 이 식을 3차원으로 일반화하면 되는데, 2차원 선적분에서 적분 구간 변화에 대한 항은 3차원 부피 적분에서는 적분 영역을 둘러싼 곡면의 변화에 해당한다. 이를 생각하면 아래와 같은 일반화가 쉽게 이해가 된다.
$$\frac{d}{dt}\int_{V(t)} F(\mathbf{x}, t)dV = \int_{V(t)} \frac{\partial F}{\partial t} dV + \int_{A(t)} F \mathbf{u}\cdot \mathbf{dA}$$
여기서 $\mathbf{u}$는 boundary의 속도를 의미한다.
그런데 이 방정식은 다소 마음에 들지 않는다. 좌변과 우변의 첫 번째 항은 부피적분인데 우변의 두 번째 항은 면적분이다. 다 같은 부피적분으로 표현할 수 있으면 적분 자체를 날릴 수가 있다. 그래서 면적분을 부피적분으로 바꾸기 위해 여기서 divergence theorem을 사용한다. 그러면 우변의 두 번째 항이 $\int _V \nabla \cdot (F\mathbf{u})dV$가 돼서 이 식을 정리하면 아래의 Kinematic transport theorem (KTT)가 유도된다.
$$\frac{d}{dt} \int_V FdV = \int_V \left[ \frac{\partial F}{\partial t} + \frac{\partial}{\partial x_j}(Fu_j) \right] dV$$
KTT를 통해서 알 수 있는 것은, 내가 관심있는 변수를 움직이는 부피에 대해서 적분하고 시간에 대해서 미분한 결과이다.
이제 보존방정식을 유도할 수 있다. 먼저 질량보존을 생각해보자. 움직이는 부피에 대한 질량 보존을 어떻게 생각할 수 있을까? 움직이거나 변형하는 부피 $V$에 대해서 질량의 총합은 $\int \rho dV$이다. 이 질량의 총합은 시간이 지나도 변하면 안된다. 따라서 질량보존의 법칙은
$$\frac{d}{dt}\int_V \rho dV =0 $$
이다.
여기에 KTT를 적용하면
$$\frac{d}{dt}\int_V \rho dV = \int_V \left[ \frac{\partial \rho}{\partial t} + \frac{\partial}{\partial x_j}(\rho u_j) \right] dV = 0$$
어떤 material volume에 대해서라도 항상 0이어야 하므로 피적분함수 자체가 0임이 자명하다.
따라서 다음이 성립한다.
$$\frac{\partial \rho}{\partial t} + \frac{\partial}{\partial x_i} (\rho u_i) = \frac{\partial \rho}{\partial t} + u_i \frac{\partial \rho}{\partial x_i} + \rho \frac{\partial u_i}{\partial x_i} =\frac{D\rho}{Dt}+\rho \frac{\partial u_i}{\partial x_i} =0 $$
그리고 다양한 경우에 물 등의 유체는 비압축성이라고 가정을 한다. 비압축성이라는 것은 유체의 밀도가 변하지 않는다는 것이다. 그래서 위 식에서 가장 우변에서 $\frac{D\rho}{Dt}=0$이 되어서 비압축성 유체의 질량보존의 법칙은 아래와 같다.
$$\frac{u_i}{x_i} = 0$$
어디서 본 적이 있는 형태이다. Kinematics에서 linear strain rate의 합으로 volumetric strain rate을 정의했었는데, 그 volumetric strain rate가 0이라는 것이다. 어떻게 보면 당연한 것이, 밀도가 변하지 않는다면 각 방향으로의 변형의 합이 0이 되어야 하는 것이 질량보존이랑 의미가 통한다.
그 다음은 운동량보존이다. 운동량보존을 유도하기 전에 먼저 응력에 대해서 얘기를 해야 한다. 뉴턴 유체에서는 모멘텀의 경과와 동점성계수를 곱한 것이 전단응력인데, 따라서 전단응력이란 곧 모멘텀의 전달을 말한다.
위와 같은 유체 요소에서, 일반적으로 응력 $\tau$를 표시할 때 첫 번째 첨자는 그 면에 수직인 방향을 나타내고 두 번째 첨자는 응력의 방향을 나타낸다. 따라서 두 첨자가 같은 성분이 수직응력인 것이고, 두 첨자가 다른 성분들이 전단응력인 것이다.
그리고 이러한 성분들로 구성된 응력 텐서가 대칭이라는 것을 보이기 위해서, 아래의 그림을 고려하자.
이러한 평면에서 토크를 발생시키는 것은 오직 shear stress이다. centroid axis에서의 토크를 계산하면 다음과 같다.
$$T=\left[\tau_{12}+\frac{1}{2}\frac{\partial \tau _{12}}{\partial x_1} dx_1 \right]dx_2 \frac{dx_1}{2}+\left[\tau_{12}-\frac{1}{2}\frac{\partial \tau_{12}}{\partial x_1} dx_1 \right] dx_2 \frac{dx_1}{2}$$
$$-\left[\tau_{21}+\frac{1}{2}\frac{\partial \tau_{21}}{\partial x_2} dx_2 \right] dx_1 \frac{dx_2}{2}-\left[\tau_{21}-\frac{1}{2}\frac{\partial \tau_{21}}{\partial x_2} dx_2 \right] dx_1 \frac{dx_2}{2} $$
계산하여 정리하면
$$T=(\tau_{12}-\tau_{21})dx_1 dx_2 $$
인데, 원래 토크의 정의에 의하면 $T=I\omega _3$이고 여기서 $I$는 관성모멘트이고 $\omega _3$는 각가속도이다. 이 정의에 의해 토크를 계산하면 $I=dx_1 dx_2(dx_1^2 +dx_2^2)\rho /12$이므로
$$ (\tau_{12}-\tau_{21})dx_1 dx_2 = dx_1 dx_2(dx_1^2 +dx_2^2)\rho \omega _3 /12 $$
인데, 유체의 체적이 infinitesimal할 때, $dx_1$과 $dx_2$는 0이 되므로 결국 우변이 0이 되어서
$$\tau_{ij} = \tau_{ji}$$
이다. 즉 응력텐서는 대칭텐서이다.
운동량보존법칙은 사실 잘 알고 있는 뉴턴의 힘-운동량 법칙을 유체에 알맞게 작성한 것에 불과하다.
뉴턴의 법칙에 따르면 $F=ma$에서 $F=m du/dt$라고 할 수 있다. 즉 어떤 물체에 작용하는 모든 힘의 합은 질량에 속도를 미분한 것을 곱한다는 뜻인데, 여기서는 단위부피당 힘을 고려하기 위해 양변을 부피로 나눈다. 그러면 결과적으로 모든 힘을 합한 것은 밀도에 속도의 미분을 곱한 것과 같다는 결론을 얻을 수 있다.
그렇다면 유체에 작용하는 모든 힘을 찾아야 하는데, 일반적으로 힘은 body force와 surface force로 나눌 수 있다고 했다. 작용하는 body force는 중력만이 유일하다고 가정하고, surface force는 응력만이 유일하다.
위의 그림을 바탕으로, 모든 방향으로 일반화하기 전에 먼저 $x_1$방향에 대해 surface force를 계산하면 아래와 같다.
$$\left(\tau_{11}+\frac{\partial \tau_{11}}{\partial x_1} \frac{dx_1}{2} - \tau_{11}+\frac{\partial \tau_{11}}{\partial x_1}\frac{dx_1}{2} \right) dx_2 dx_3$$
$$+\left(\tau_{21}+\frac{\partial \tau_{21}}{\partial x_2} \frac{dx_2}{2} - \tau_{21}+\frac{\partial \tau_{21}}{\partial x_2}\frac{dx_2}{2} \right) dx_1 dx_3$$
$$+\left(\tau_{31}+\frac{\partial \tau_{31}}{\partial x_3} \frac{dx_3}{2} - \tau_{31}+\frac{\partial \tau_{31}}{\partial x_3}\frac{dx_3}{2} \right) dx_1 dx_2$$
이를 정리하면 아래와 같이 된다.
$$\left(\frac{\partial \tau_{11}}{\partial x_1} + \frac{\partial \tau_{21}}{\partial x_2} + \frac{\partial \tau_{31}}{\partial x_3} \right)dx_1 dx_2 dx_3 = \frac{\partial \tau_{j1}}{\partial x_j}dV$$
$x_2$와 $x_3$의 방향에 대해서도 동일한 결과를 얻을 수 있고 부피 $dV$로 나눠서 얻은 단위부피당 힘의 $i$번째 요소는 아래와 같다.
$$\frac{\partial \tau_{ji}}{\partial x_j} = \frac{\partial \tau_{ij}}{\partial x_j}$$
Surface force의 합을 구했으므로 단위부피당 $F=ma$와 같은 형태로 정리하면 아래와 같다.
$$\rho \frac{Du_i}{Dt} = \rho g_i + \frac{\partial \tau_{ij}}{\partial x_j}$$
이 방정식은 유체든 고체든 상관없이 적용되는 방정식이다. 이는 Cauchy's equation of motion이라는 이름으로도 불린다.
이제 위의 방정식의 마지막 항인 $\partial \tau_{ij} / \partial x_j$와 모멘텀이 어떤 관련이 있는지만 알아내면 모멘텀보존법칙의 유도가 끝난다. 지난 번의 Kinematics에서 보였듯이, deformation과 momentum은 관련이 큰다. 따라서 여기에서는 stress와 deformation의 관계를 나타내는 constitutive equation을 고려한다.
먼저, 유체가 정지해있을 때는 모멘텀이 없기때문에 당연히 모멘텀 전달도 없고 따라서 전단응력이 없다. 따라서 수직응력만이 존재하는데, 그 수직응력은 thermodynamic pressure, 즉 그냥 일반적으로 생각하는 압력이 유일하다. 이러한 압력은 항상 면에 수직인 방향으로 작용하므로 응력 텐서에서 대각원소만이 0이 아닌 값을 갖는다. 따라서 이 경우에 응력텐서는 크로네커 델타를 사용해 다음과 같이 쓸 수 있다.
$$\tau_{ij} = -p \delta_{ij}$$
그런데 유체가 정지해있지 않고 움직이는 경우에는 점성에 의한 응력, 즉 전단응력에 대한 항이 추가되어야 한다. 따라서 움직이는 유체에서의 응력 텐서는 아래와 같이 쓸 수 있다.
$$\tau_{ij}=-p\delta_{ij}+\sigma_{ij}$$
정지한 유체에서의 응력텐서는 등방텐서, 즉 모든 방향에 대해 그 값이 같았다. 하지만 움직이는 유체에서는 $\sigma_{ij}$에 의해 등방성을 잃게 되는데 여기서 비등방성에 해당하는 $\sigma$를 deviatoric stress tensor라고 부른다. 한국어로는 보통 편차응력텐서라고 한다.
이 편차응력텐서는 전단응력과 관련이 있기때문에 모멘텀 gradient와도 관련이 있을 것이다. 따라서 먼저 momentum gradient tensor를 살펴보면 지난 글의 끝에서 언급했듯이 아래와 같이 두 부분으로 분리할 수 있다.
$$\frac{\partial u_i}{\partial x_j} = \frac{1}{2}\left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i} \right) + \frac{1}{2} \left(\frac{\partial u_i}{\partial x_j} - \frac{\partial u_j}{\partial x_i} \right)$$
우변의 두 번째 항인 antisymmetric part는 변형이 없는 rotation을 나타내므로 이 part에서는 응력이 생성될 수 없다. 응력은 반드시 우변의 첫 번째 항인 strain rate tensor에서 생성된다.
$$e_{ij} = \frac{1}{2}\left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i} \right) $$
가 응력이 관련이 있는 부분인데, 이 strain rate tensor가 편차응력텐서와 선형적인 관계가 있다고 가정해보자. 그렇다면 어떠한 4차 텐서 $K_{ijmn}$에 대하여 다음의 선형 관계가 성립한다는 것이다.
$$\sigma_{ij}=K_{ijmn}e_{mn}$$
Constitutive equation이 선형 방정식이라는 것은 일반적으로 고체에서도 쓰이는 접근법이다. 그리고 곱해지는 텐서가 4차 텐서인 것은 편차응력텐서와 strain rate tensor가 모두 2차 텐서이기 때문에 어떤 텐서와 2차 텐서를 곱해서 2차 텐서가 나오기 위해서는 $K$가 4차 텐서여야 한다. 그리고 여기서는 3가지의 방향을 고려하기 때문에 텐서의 차수가 높아질수록 원소의 개수가 3배가 되고, 따라서 $K_{ijmn}$의 원소의 개수는 81개이다.
그리고 유체가 isotropic하다고 가정하자. 즉, 편차응력텐서와 strain rate tensor의 관계가 방향에 따라 다르지 않다고 가정하자. 여기서 편차응력텐서와 strain rate tensor의 관계라는 것은 다시 말해서 $K_{ijmn}$이므로 이러한 가정 하에서는 $K_{ijmn}$이 반드시 등방텐서이다.
그리고 tensor analysis에 대한 책을 보면, 짝수 차수의 등방 텐서는 모두 2차 크로네커 델타의 곱으로 구성되어있다는 것이 알려져있고, 특히 모든 4차 등방텐서는 다음과 같은 형태를 갖는다는 것이 알려져있다.
$$K_{ijmn} = \lambda \delta_{ij} \delta_{mn} + \mu \delta_{im} \delta_{jn} + \gamma \delta_{in} \delta_{jm}$$
그리고 편차응력텐서 $\sigma_{ij}$가 대칭텐서라고 가정하자. 이 가정을 말이 되는 것이, 위에서 응력텐서가 대칭텐서임을 보였고, 응력텐서는 압력텐서와 편차응력텐서의 합인데 압력텐서가 대칭이므로 편차응력텐서도 대칭일 것이다. 그렇다면 $K_{ijmn}$또한 $i$와 $j$에 대해서 대칭이어야 한다. 이 말은 위의 식에서 $i$와 $j$를 바꿔도 같은 결과가 나와야 한다는 것이므로 이를 위해서는 $\gamma = \mu$가 성립해야 한다.
기존에 $K$를 구성하던 요소는 81개라고 했는데 이제는 $\lambda$와 $\mu$ 두 개밖에 남지 않았다. 요소를 81개에서 2개로 줄이는데 필요한 가정은 유체의 등방성과 편차응력의 대칭성이다. 이렇게 얻은 $K_{ijmn}$을 위의 $\sigma_{ij}=K_{ijmn}e_{mn}$에 대입하면 아주 쉽게 다음의 관계를 얻을 수 있다.
$$\sigma_{ij} = 2\mu e_{ij} + \lambda e_{mm} \delta_{ij}$$
그렇다면 이제 편차응력텐서를 알았으니 응력텐서 $\tau_{ij}$로 돌아가면, 다음과 같이 쓸 수 있다.
$$\tau_{ij} = -p\delta_{ij} + 2\mu e_{ij} + \lambda e_{mm} \delta_{ij} $$
이제 두 상수 $\mu$와 $\lambda$의 관계를 얻기 위해 $i=j$일 때 겹치는 인덱스에 대한 합을 수행하면 다음을 얻는다.
$$\tau_{ii} = -3p+(2\mu+3\lambda ) e_{mm}$$
이를 $p$에 대해 정리하고 strain rate tensor $e_{mm}$의 정의를 대입하면 다음을 얻을 수 있다.
$$p=-\frac{1}{3} \tau_{ii} + \left(\frac{2}{3} \mu + \lambda \right) \nabla \cdot \mathbf{u}$$
그런데, 지금은 유체가 움직이고 있는 상황이기 때문의 $e_{ij}$의 diagonal terms는 같지 않다. 이 diagonal term이란 각 방향에 대한 linear strain rate을 말하는 것이므로 당연히 같지 않겠다. 그렇다면 위의 $\tau_{ii} = -3p+(2\mu+3\lambda ) e_{mm}$에서 알 수 있듯이, $e_{ij}$의 diagonal terms인 $ e_{mm}$에 곱해진 항이 있기 때문에 $\tau_{ij}$의 diagonal term인 $\tau_{ii}$가 같지 않다. 그래서 $\tau_{ij}$의 diagonal terms의 평균을 취해 mean pressure라는 값을 정의한다.
$$\bar{p} = -\frac{1}{3}\tau_{ii}$$
이를 위 식에 대입하면 다음을 얻는다.
$$p-\bar{p} = \left(\frac{2}{3}\mu + \lambda \right) \nabla \cdot \mathbf{u}$$
그리고 이 방정식의 우변에서, 많은 경우에 적용되는 Stokes 가정에 따르면 $\lambda + 2\mu / 3 =0$을 사용하는 것이 합리적이다.
따라서 이를 $\sigma_{ij} = 2\mu e_{ij} + \lambda e_{mm} \delta_{ij}$에 대입하면 다음과 같은표현을 얻는다.
$$ \sigma_{ij} = 2\mu \left(e_{ij} -\frac{1}{3}e_{mm} \delta_{ij}\right) $$
그리고 비압축성유체에서는 질량보존의 법칙에 의해 $e_{mm}=0$이므로 응력텐서는 다음과 같다.
$$\tau_{ij} = -p\delta_{ij}+2\mu e_{ij} = -p\delta_{ij} + \mu \left(\frac{\partial u_i}{\partial x_j} + \frac{\partial u_j}{\partial x_i} \right) $$
이제 응력텐서에 대한 정의를 마쳤다.
아까 위에서 정의한 Cauchy's equation of motion은 다음과 같았다.
$$\rho \frac{Du_i}{Dt} = \rho g_i + \frac{\partial \tau_{ij}}{\partial x_j}$$
응력텐서자리에 우리가 얻은 응력텐서를 대입하면 다음과 같다.
$$ \rho \frac{Du_i}{Dt} = -\rho g \delta _{i3} + \frac{\partial}{\partial x_j} \left[ -p \delta _{ij} +\mu \left( \frac{\partial u_i}{\partial x_j} +\frac{\partial u_j}{\partial x_i} \right) \right] = \rho g \delta _{i3} - \frac{\partial p}{\partial x_i} + \mu \left( \frac{\partial ^2 u_i}{\partial x_j ^2} + \frac{\partial ^2 u_j}{\partial x_i \partial x_j} \right)$$
그런데 여기서 가장 우변의 $\mu$에 곱해진 괄호에서
$$\frac{\partial ^2 u_j}{\partial x_i \partial x_j} = \frac{\partial}{\partial x_i} \frac{\partial u_j}{\partial x_j} $$
인데 질량보존법칙에 의해 $\partial u_j / \partial x_j = 0$이므로 0이 된다.
양변을 $\rho$로 나눠서 나타내면 최종적으로 다음을 얻는다.
$$\frac{\partial u_i}{\partial t} + u_j \frac{\partial u_i}{\partial x_j} = -\frac{1}{\rho} \frac{\partial p}{\partial x_i} -g \delta _{i3} + \nu \frac{\partial ^2 u_i}{\partial x_j ^2} $$
이게 바로 그 유명한 비압축성 나비에 스토크스 방정식의 momentum equation이다. 이 momentum equation과 질량보존법칙을 합쳐서 부르는 것이 나비에 스토크스 방정식이다.
이 글의 목적이었던 나비에 스토크스 방정식까지 도달하는 것은 성공했다. 이제 이와 관련하여 몇 가지 내용들을 더 살펴보고 글을 마무리한다.
먼저 rotating frame에 대한 내용인데, 위의 모든 내용들은 좌표계가 고정되어있을 때의 내용들이다. 하지만 지구 스케일에서 일어나는 유체역학 현상을 분석하는 경우처럼, 좌표계가 고정되어있지 않고 회전할 때는 좌표계의 회전을 고려해야 한다.
어떤 벡터 $\mathbf{P}$가 있을 때, 단위 방향 벡터 $\hat{i}$에 대해서 $\mathbf{P}=\hat{P}\hat{i}$와 같이 나타낼 수 있다. 이 때 좌표계가 회전한다면 $\mathbf{P}$의 미분은 좌표계의 변화까지 고려해서 다음과 같이 나타내야 한다.
$$\frac{d\mathbf{P}}{dt} = \frac{d\hat{P}}{dt} \hat{i} + \hat{P} \frac{d\hat{i}}{dt}$$
여기서 $\hat{i}$가 회전하는 중심축이 $\mathbf{\Omega}$라면 $d\hat{i}/dt=\mathbf{\Omega} \times \hat{i}$로 나타낼 수 있다.
위의 그림에서 해당 내용을 쉽게 생각해볼 수 있다. 그렇다면 속도와 가속도도 아래와 같이 위치를 미분해서 얻을 수 있다. 여기에서 아래첨자 $F$는 fixed observer가 관찰했다는 것이고, $R$은 rotating observer가 관찰했다는 것이다.
$$\mathbf{u}_F = \frac{d\mathbf{r}}{dt}_F = \mathbf{u}_R + \mathbf{\Omega} \times \mathbf{r}$$
$$\mathbf{a}_F = \frac{d\mathbf{u}_F}{dt}_F = \mathbf{a}_R + \mathbf{\Omega} \times \mathbf{\Omega} \times \mathbf{r} + 2 \mathbf{\Omega} \times \mathbf{u}_R $$
그리고 외적의 정의를 고려하면 위의 그림에서 $\mathbf{R}$를 도입했을 때, $\mathbf{\Omega}\times \mathbf{r} = \mathbf{\Omega} \times \mathbf{R}$임을 쉽게 알 수 있다. 그리고 vector identity를 이용하면 $\mathbf{\Omega} \times \mathbf{\Omega} \times \mathbf{R} = -\Omega ^2 \mathbf{R}$임을 알 수 있다. 따라서 다음이 성립한다.
$$\mathbf{a}_F = \mathbf{a} + 2\mathbf{\Omega} \times \mathbf{u} - \Omega ^2 \mathbf{R}$$
여기서 우변의 첫 번째 항은 Coriolis force라고 부르는 전향력의 가속도이고 세 번째 항은 회전반경 곱하기 각속도의 제곱이므로 구심가속도이다.
이제 나비에 스토크스 방정식에 이 가속도 표현을 대입하면 다음과 같다.
$$\frac{D\mathbf{u}}{Dt} = -\frac{1}{\rho} \nabla p + \nu \nabla ^2 \mathbf{u} + (\mathbf{g}_n + \Omega ^2 \mathbf{R} ) -2\mathbf{\Omega} \times \mathbf{u} $$
이 글을 처음 시작할 때 보존방정식은 질량, 모멘텀, 에너지 세 개가 있다고 했는데 앞의 두 개만 설명하고 에너지보존방정식은 설명하지 않았다. 여기서 간략하게 설명한다.
아래와 같은 equation of motion에
$$\rho \frac{Du_i}{Dt} = \rho g_i + \frac{\partial \tau _{ij}}{\partial x_j} $$
양변에 $u_i$를 곱하면 다음과 같은 결과를 얻는다.
$$ \rho \frac{D}{Dt} \left( \frac{1}{2} u_i^2 \right) = \rho u_i g_i + u_i \frac{\partial \tau_{ij}}{\partial x_j}$$
그리고 continuity equation에 의해 다음이 성립한다.
$$ \frac{1}{2} \rho u_i^2 \left[\frac{\partial \rho}{\partial t} + \frac{\partial}{\partial x_j} \left(\rho u_j \right) \right] = 0 $$
이 식의 값이 0이기때문에 위의 $u_i$가 곱해진 equation of motion에 더할 수 있다. 식을 더하고 물질미분을 국부미분과 대류미분으로 분리해서 식을 정리하면 다음이 성립한다.
$$ \frac{\partial}{\partial t} \left( \frac{1}{2} \rho u_i^2 \right) + \frac{\partial}{\partial x_j} \left[ u_j \frac{1}{2} \rho u_i^2 \right] = \rho u_i g_i + u_i \frac{\partial \tau_{ij}}{\partial x_j} $$
여기서 $E=\rho u_i^2 /2$라고 정의하면 다음과 같은 형태가 나타난다.
$$\frac{\partial E}{\partial t} + \nabla \cdot (\mathbf{u}E) = \rho \mathbf{u} \cdot \mathbf{g} + \mathbf{u} \cdot (\nabla \cdot \mathbf{\tau})$$
이게 에너지보존방정식이다.
유체역학을 고등학교 물리에서도 가르치는데, 그 때 가르치는 게 베르누이 방정식이다. 이는 사실 나비에 스토크스 방정식에서 점성항을 0이라고 가정하고 적분한 것이다. 둘이 전혀 관련이 없는 것이 아니다.
먼저 $\nu=0$인 비점성유동, inviscid flow를 설명하는 오일러 방정식은 다음과 같다.
$$\frac{\partial u_i}{\partial t} + u_j \frac{\partial u_i}{\partial x_j} = -\frac{\partial}{\partial x_i} (gz) - \frac{1}{\rho} \frac{\partial p}{\partial x_i} $$
여기서 대류가속도는 vorticity를 통해서 다음과 같이 표현될 수 있다.
$$u_j \frac{\partial u_i}{\partial x_j} = u_j \left( \frac{\partial u_i}{\partial x_j} - \frac{\partial u_j}{\partial x_i} \right) + u_j \frac{\partial u_j}{\partial x_i} = -(\mathbf{u} \times \mathbf{\omega} )_i + \frac{\partial}{\partial x_i} \left( \frac{1}{2} q^2 \right) $$
이제 오일러 방정식은 다음과 같이 나타내어진다.
$$ \frac{\partial u_i}{\partial t} + \frac{\partial}{\partial x_i} (\frac{1}{2} q^2) +\frac{1}{\rho} \frac{\partial p}{\partial x_i} + \frac{\partial}{\partial x_i} (gz) = (\mathbf{u} \times \omega )_i $$
여기서 barotropic flow라는 가정을 하는데, 이는 밀도가 오직 압력에 대한 함수라는 가정이다. 이 경우에는 pressure term을 다음과 같이 나타낼 수 있다.
$$\frac{1}{\rho} \frac{\partial p}{\partial x_i} = \frac{\partial}{\partial x_i} \int \frac{dp}{\rho} $$
그러면 오일러 방정식은 다음과 같이 쓸 수 있다.
$$\frac{\partial u_i}{\partial t} + \frac{\partial}{\partial x_i} \left[ \frac{1}{2} q^2 + \int \frac{dp}{\rho} + gz \right] = (\mathbf{u} \times \mathbf{\omega} )_i $$
그리고 베르누이 function을 다음과 같이 정의하면
$$B=\frac{1}{2} q^2 + \int \frac{dp}{\rho} + gz $$
오일러 방정식은 다음과 같이 정리된다.
$$\frac{\partial \mathbf{u}}{\partial t} + \nabla B = \mathbf{u} \times \mathbf{\omega} $$
여기서 steady flow의 경우를 먼저 살펴보면 오일러 방정식은 $\nabla B = \mathbf{u} \times \omega $가 된다. 좌변은 $B$의 gradient이기 때문에 $B=$ constant인 평면에 수직인 벡터이고, 우변은 $\mathbf{u}$와 $\omega$ 모두와 수직인 벡터이다. 따라서 constant $B$인 평면은 streamline과 vortex line을 반드시 포함해야 한다. 즉, streamline과 vortex line을 따라서 $B$는 constant하다. 더 나아가 만약 흐름이 irrotational하다면, $(\omega = 0)$ 우변이 0이 되므로 B는 모든 곳에서 constant하다.
베르누이 방정식의 unsteady form은 오직 flow가 irrotational할 때만 가능하다. 흐름이 irrotational하면 이전 글에서 언급했듯이, velocity potential의 존재가 보장된다.
$$\mathbf{u} = \nabla \phi$$
이를 unsteady Euler equation에 대입하면 오일러 방정식은 다음과 같이 된다.
$$\nabla \left[ \frac{\partial \phi}{\partial t} + \frac{1}{2} q^2 + \int \frac{dp}{\rho} + gz \right] = 0$$
적분하면 공간과 관련없는 함수 $F(t)$에 대해 다음과 같이 정리된다.
$$ \frac{\partial \phi}{\partial t} + \frac{1}{2} q^2 + \int \frac{dp}{\rho} + gz = F(t) $$
'유체역학' 카테고리의 다른 글
| [유체동역학] 5. 차원해석과 상사이론[Dimensional Analysis and Dynamic Similarity] (0) | 2025.03.05 |
|---|---|
| [유체동역학] 4. 비회전유동, 포텐셜유동[Irrotational Flow, Potential Flow] (0) | 2025.03.02 |
| [유체동역학] 2. 흐름의 운동학[Kinematics] (2) | 2025.03.01 |
| [유체동역학] 1. 정수역학[Fluid Statics] (0) | 2025.02.27 |
| [유체동역학] 0. 유체역학을 공부하기 위해 알아야 하는 것들 (0) | 2025.02.27 |
