[OpenFOAM 소스코드 파해치기] 열전도

일본어 원문링크

작성일 : 2009년 12월   5일
번역일 : 2016년  2월   10일

크롬 브라우저로 보시는 것을 권장해 드립니다.

 ------ OpenFOAM 소스코드 파해치기 시리즈 ------

OpenFOAM 소스코드 파해치기 목차로 이동

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

 ----- OpenFOAM 소스코드를 다루는 문법 기본 -----

OpenFOAM 소스코드를 다루는 문법 기본으로 이동

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

 

들어가기

열전도문제를 풀어보도록 하자.

사용 버전

OpenFOAM 1.6

 

파일

이번에 사용파는 파일이다.(08-heat.tar.gz).

• Make/files
• Make/options
• system/controlDict
• system/fvSchemes
• system/fvSolution
• constant/polyMesh/blockMeshDict
• 0/T
• heat.C

격자는 이전과 같다. 온도를 위한 경계조건 파일을 필요로 하며 방식은 이전의 스칼라필드의 경우와 동일하다. 똑한, 이번에 드디어 직접 생성한 방정식을 풀게 되므로, fvSchemes, fvSolution의 설정이 필요하다.

프로그램

 

이번에는 정상상태의 열전도문제를 계산하는 프로그램을 만들도록 한다. 프로그램의 설명부터 시작해보자. 본 프로그램은 laplacianFoam을 참고하여 작성되었다.

heat.C (일부)

...
	volScalarField T(
		IOobject("T", runTime.timeName(), mesh,
			IOobject::MUST_READ, IOobject::AUTO_WRITE),
		mesh
	);

	dimensionedScalar DT("DT", dimensionSet(0, 2, -1, 0, 0, 0, 0), 1.);

	solve(fvm::laplacian(DT, T));

	Info << "max T = " << max(T).value() << endl;
	Info << "min T = " << min(T).value() << endl;
	Info << "volume weithed average T = "
		<< T.weightedAverage(mesh.V()).value() << endl;

	runTime++;
	T.write();
	...

온도장변수를 T로 하고 스칼라필드로 정의한다. T의 정의로는 "T"라고 적혀있는 부분이 있으나 이 부분은 변수명을 의미하며 보통의 프로그램에서 정의하는 변수와 일치하도록 작성한다. 여기에서 이름은, 경계조건파일의 이름으로 참조되거나 방정식의 설정에서 참조된다.

DT는 확산계수이다. 방정식에서 사용될 때 단위가 요구되므로 단위를 지정하도록 한다. 이 부분이 맞지 않으면 계산에 실패하게 된다.

"solve..." 로 너무나 간단히 한줄로 풀고 있으나, 방정식을 시우고 그것을 푸는 코드에 해당한다. solve()는 방정식을 푸는 함수로 그 인수는 fvm:laplacian()가 라플라시안이고, DT와 T를 곱해 라플라시안을 씌워 0으로 놓은 방정식을 의미한다.

T.weightedAverage(mesh.V()) 는, mesh.V()가 격자의 셀 체적 스칼라필드, 즉, 셀체적리스트를 출력하며 그것을 하용하여 weightAverage를 함으로써 체적비를 고려한 온도의 평균을 얻어내게 된다.

정상상태의 문제이므로, 1 스텝을 진행하여 결과를 출력하여 종료하게 된다.

fvSchemes과 fvSolution

프로그램에서 푸는 방정식은 각각 이산화스킴과 대수방정식해법의 설정을 필요로 한다.

fvSchemes (일부)

...
laplacianSchemes
{
    default         none;
    laplacian(DT,T) Gauss linear corrected;
}
...

프로그램 안에 적혀있는 방정식에서 사용하는 필드의 연산자 (laplacian, div, grad 등) 에는, 각각의 항애 대해 이산화스킴을 지정해야한다. 위의 설정을 보면, 프로그램에 적힌 방정식 그대로를 적어 설정하고 있다.

fvSolution (일부)

...
solvers
{
    T
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-06;
        relTol          0;
    }
}
...

변수 T에 대한 방정식을 풀기위해 솔버를 설정하고 있다.

경계조건

변수 T의 방정식에 대한 경계조건의 설정으로, 0/T 파일이 필요하다.

0/T (일부)

...
solvers
{
    T
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-06;
        relTol          0;
    }
}
...
변수 T 에 대한 방정식을 풀기 위한 솔버 설정입니다.

경계조건

변수 T 의 방정식을 위한 경계조건설정으로、0/T 파일을 필요로 합니다.

T (일부)

...
dimensions      [0 0 0 1 0 0 0];

internalField   uniform 283;

boundaryField
{
    left
    {
        type            fixedValue;
        value           uniform 300;
    }

    right
    {
        type            fixedValue;
        value           uniform 500;
    }

    top
    {
        type            zeroGradient;
    }

    bottom
    {
        type            zeroGradient;
    }
}
...

dimensions는 단위에 대한 설정으로 숫자열은 각각 [kg, m, s, K, mol, A, Cd]의 지수를 의미한다. 온도장이므로 단위는 K, 즉 [0 0 0 1 0 0 0]이 된다.

internalField는 내부영역전체의 설정이다. 여기서는 초기온도로 283K 를 지정하고 있다.

각 경계에선, 사각형 영역의 좌측에 300K, 우측에 500K 로 설정하여 위아래는 단열조건으로 하고 있다.

실행

컴파일하여 실행해 본다.

$ wmake
$ blockMesh
$ ./heat

paraFoam에서 스칼라값을 표시해보면, 값이 셀의 순번에 따라 들어가 있는 것을 확인할 수 있다.

추가정보

 

비정상열전도문제를 푸는 프로그램으로 개조해보자.(08-heat2.tar.gz).

heat.C (일부)

...
	volScalarField T(
		IOobject("T", runTime.timeName(), mesh,
			IOobject::MUST_READ, IOobject::AUTO_WRITE),
		mesh
	);

	IOdictionary transportProperties(
		IOobject(
			"transportProperties",
			runTime.constant(), mesh,
			IOobject::MUST_READ, IOobject::NO_WRITE
		)
	);

	dimensionedScalar DT(transportProperties.lookup("DT"));

	while(runTime.loop()){
		Info << "Time = " << runTime.timeName() << endl;

		solve(fvm::ddt(T) - fvm::laplacian(DT, T));

		Info << "max T = " << max(T).value() << endl;
		Info << "min T = " << min(T).value() << endl;
		Info << "volume weithed average T = "
			<< T.weightedAverage(mesh.V()).value() << endl;

		runTime.write();
	}

	runTime.writeNow();

	Info << "execution time: " << runTime.elapsedCpuTime() << " s" << endl;
	Info << "clock time: " << runTime.elapsedClockTime() << " s" << endl;
	...

DT를 외부파일로부터 읽어들이는 것처럼 하고 있다. transportProperties의 정의가 그것에 해당하며 constant/transportProperties 로부터 DT를 읽어들이고 있다.

transportProperties (일부)

...
DT              DT [ 0 2 -1 0 0 0 0 ] 1.;
...

transportProperties.lookup("DT") 로 불러들인 DT를 참조하고 있다.

시간에 따라 값이 변화하므로 계산 루프를 작성하고 있다. 방정식은 시간미분항 "fvm::ddt(T)"를 추가하기면 하면 된다.

OpenFOAM 소스코드 파해치기 목차로 이동