참조
카테고리 : 슬롯 머신 NOA 라이브러리
메인파일 : flowobs.gms
$onText
직사각형 영역에 있는 비압축성 유체의 정지 흐름
장애물의 존재.
다음에서 적응:
McKinney, D.C., Savitsky, A.G., 물과 물에 대한 기본 최적화 모델
에너지 관리. 1999년 6월(개정 6, 2003년 2월).
http://www.ce.utexas.edu/prof/mckynney/ce385d/papers/슬롯 머신-Tutorial.pdf
$offText
X /u1*u20/ 설정;
Y 설정 /u1*u20/;
* 물 이동 방정식의 영역 결정
vyside(X,Y)를 설정합니다.
vxside(X,Y)를 설정합니다.
vxside(X,Y) = 예;
vxside(X,Y)$(ord(X)=1) = 아니오;
vxside(X,Y)$(ord(X)=카드(X)) = 아니오;
vxside(X,Y)$(ord(Y)=1) = 아니오;
vxside(X,Y)$(ord(Y)=카드(Y)) = 아니오;
vxside('u10','u10') = 아니요;
vxside('u10','u11') = 아니요;
vyside(X,Y) = vxside(X,Y);
* 매개변수
x 방향의 스칼라 dx 스텝 공간 /1/;
y 방향의 스칼라 dy 스텝 공간 /1/;
스칼라 r 유체 밀도 /1000/;
매개변수 m(X,Y) 동점도 ;
m(X,Y) := 0.5;
변수
obj 목적변수
D(X,Y) 오류
P(X,Y) 압력
Vx(X,Y) x 방향 속도
Vy(X,Y) y 방향 속도 ;
* 변수 범위 및 초기화
Vx.1(X,Y) = 0.0;
Vy.1(X,Y) = 0.0;
Vy.l(x,y)$(vxside(x,y)) = 0.0;
D.lo(X,Y) = 0.0;
D.up(X,Y) = 10.0;
P.up(X,Y) = 20000;
P.lo(X,Y) = -20000;
P.1(X,Y) = 0.0;
*경계 조건
Vx.fx('u1',Y) = 0.5;
Vx.fx('u20',Y) = 0.5;
Vx.fx(X,'u1') = 0;
Vx.fx(X,'u20') = 0;
Vy.fx('u1',Y) = 0;
Vy.fx('u20',Y) = 0;
Vy.fx(X,'u1') = 0;
Vy.fx(X,'u20') = 0;
* 장애물 설명
vx.fx('u10','u10') = 0;
vx.fx('u10','u11') = 0;
방정식
For_Vx(X,Y)
For_Vy(X,Y)
Div_Vxy(X,Y)
eobj 목적 함수 ;
For_Vx(X,Y)$(vxside(X,Y))..
(P(X+1,Y)-P(X,Y))/(r*dx) =e=
m(x,Y)*((Vx(X+1,Y)-2*Vx(X,Y)+Vx(x-1,Y))/(dx*dy)
+ (Vx(X,Y+1)-2*Vx(X,Y)+Vx(X,Y-1))/(dx*dy));
For_Vy(X,Y)$(vxside(X,Y))..
(P(X,Y+1)-P(X,Y))/(r*dy) =e=
m(X,Y)*((Vy(X+1,Y)-2*Vy(X,Y)+Vy(X-1,Y))/(dx*dy)
+ (Vy(X,Y+1)-2*Vy(X,Y)+Vy(X,Y-1))/(dy*dx));
*--
*For_Vx(X,Y)$(Vxside(X,Y))..
* Vx(X,Y)*(Vx(X+1,Y)-Vx(X-1,Y))/(2*dx) +
* 0.25*(Vy(X+1,Y-1)+Vy(X+1,Y)+Vy(X,Y-1)+Vy(X,Y)) *
* (Vx(X,Y+1)-Vx(X,Y-1))/(2*dy) +
* (P(X+1,Y)-P(X,Y))/(r*dx)
* =e=
* m(X,Y)*((Vx(X+1,Y)-2*Vx(X,Y)+Vx(X-1,Y))/(dx*dx) +
* (Vx(X,Y+1)-2*Vx(X,Y)+Vx(X,Y-1))/(dy*dy));
*
*For_Vy(X,Y)$(Vyside(X,Y))..
* 0.25*(Vx(X-1,Y+1)+Vx(X-1,Y)+Vx(X,Y+1)+Vx(X,Y)) *
* (Vy(X+1,Y)-Vy(X-1,Y))/(2*dy) +
* Vy(X,Y)*(Vy(X,Y+1)-Vy(X,Y-1))/(2*dy) +
* (P(X,Y+1)-P(X,Y))/(r*dy)
* =e=
* m(X,Y)*((Vy(X+1,Y)-2*Vy(X,Y)+Vy(X-1,Y))/(dx*dx) +
* (Vy(X,Y+1)-2*Vy(X,Y)+Vy(X,Y-1))/(dy*dy));
*--
Div_Vxy(X,Y)$((ord(X) > 1) $ (ord(Y) > 1))..
(Vx(X,Y)-Vx(X-1,Y))/dx + (Vy(X,Y)-Vy(X,Y-1))/dy =e=
D(X,Y);
eobj..obj =e= SUM((X,Y),(D(X,Y)*D(X,Y)));
모델 흐름 /all/;
obj를 최소화하는 nlp를 사용하여 flowob을 해결합니다.
$ifThenI x%mode%==xbook
* 솔루션을 넣어
파일 res1 /압력.dat/
res1을 넣어;
loop(X, loop(Y, put P.l(x,y):10:4; ); put /;); 놓다 /;
파일 res2 / vx.dat/
res2를 넣어;
loop(X, loop(Y, put Vx.l(x,y):9:5; ); put /;); 놓다 /;
파일 res3 / vy.dat/
res3을 넣어;
loop(X, loop(Y, put Vy.l(x,y):5:1; ); put /;); 놓다 /;
파일 res4 /err.dat/
res4를 넣어
loop(X, loop(Y, put D.l(x,y):5:1; ); put /;); 놓다 /;
$endIf
* 흐름 종료
*************************************************** 2008년 11월 19일