MODEL MATEMATYCZNY

w a T =a [x'“a (x*3T J+32 T ] , o<x<R , 0<z<L (5.2) z 2 2 x x 2 zz 2

gdzie: x i z oznaczają współrzędne przestrzenne, T^=T^(r,z) temperaturę odpowiednio fazy ciekłej e=l i fazy stałej e=2, parametr m określa geometrie wlewka (m=0 wlewek płaski, m=l wlewek okrągły), w jest poszukiwaną prędkością wyciągania, a# , e=l,2 są współczynnikami przewodzenia temperatury: a funkcja a=o(z) opisuje położenie granicy rozdz iału f a z .

Z uwagi na znikomą przewodność ciepła w kierunku wyciągania wlewka [5.11], równania (5.1)+(5.2) można uprościć opuszczając po prawej stronie drugi składnik, w rezultacie czego otrzymamy równania paraboliczne, w których współrzędna z (zgodna z kierunkiem wyciągania wlewka) spełnia role czasu.

59

-Rys.5.i. Modelowany obiekt Fig.5.1. Modelled object

Niech ę i %> oznaczają. odpowiednio podyktowane wymogami

n n

technologicznymi grubości naskórka i spadki temperatury na jego przekroju w ustalonych punktach kontrolnych £ €(0,L), n=l,N, a T* poszukiwaną

n

temperaturę zalewania. Załóżmy Jeszcze, Ze proces wymiany ciepła między wlewkiem a krystalizatorem określa warunek Fouriera i Ze znany jest współczynnik wymiany ciepła ct=a(z), ze<0,L).

Przyjęte założenia oraz wprowadzony do rozważań wektor parametrów {p}J={w,T*} pozwalaja sformułować matematyczny opis zagadnienia w postaci:

p a T = a X~m0 (K*iT ) , 0<x<CT , 0<z<L (5.3)

1 Z i 1 X X 1

p i T = a x ' ' i |*‘i T ) , o<x<R , 0<z<L (5.4) 1 z 2 2 x x 2

Tt=p2 , z=0 , 0<x<R (5.5)

£ T =0 . x=0 , 0<z£L

x 1 (5.6)

T =T , x=er , 0<z<L (S.7) masy odpowiednio fazy ciekłej oraz fazy stałej, Q jest ciepłem krzepnięcia,

5.3. METODA ROZWIĄZANIA

Postawione zadanie polegające na wyznaczeniu dla krzepnących wlewków prędkości wyciągania i temperatury zalewania zapewniających właściwe formowanie się naskórka rozwiązano wykorzystując metody optymalizacyjne.

Do opisu materraatycznego dołączono funkcję jakości:

granicy rozdziału faz i temperaturę powierzchni wlewka odpowiadające punktom kontrolnym, wyznaczane w wyniku rozwiązania zagadnienia brzegowego

oznaczają parametry wagowe.

Minimum funkcji (5.12) poszukuje się metodami gradientowymi [5.1,5.9,5.13]. W tym celu należy wyznaczyć pochodne danej funkcji względem poszukiwanych parametrów.

- 61

Mamy więc:

V F=;

k= 2 I K (V ' „ • * = 1 . 2 ( 5 . 1 3 )

n = 1

gdzie VkF=apkF({p}Z) , « ^ = a t« - {p} J) . a u £ "«p T2( R . {p>2) .

W wyrażeniach (5.13) występuje po N wartości »k n i , n=l7N odpowia­

dających punktom kontrolnym funkcji 9 = 9 T(z;{p}2) i U =d T (x,z:fDł2) k Pk i 21 pk 2 , l w | Do ich wyznaczenia wykorzystuje się układ równań i warunków (5.3)+(5.9).

Różniczkując równania (5.3)+(5.9) kolejno po p^ i p^ otrzymuje się dwa proste zagadnienia brzegowe postaci:

W P A W a x (x a xU l,l) ’ 0 < x < a ’ 0 < z < L ( 5 . 1 4 )

ia x^24= ^2X - W a x ^X m a x^2,l ^ ’ 0 < X < R ■ 0 < Z - L ( 5 . 1 5 )

U u = 0 , z = 0 , 0 < x < R ( 5 . 1 6 )

d xu u = ° , x = 0 , 0 < z < L ( 5 . 1 7 )

U . 4 + e 1a ltT ,>= 0 , x = a , 0 < z < L , e = l , 2 ( 5 . 1 8 )

~X1( a xU i,l+ a i<9L T i* + X2 1 a xU2,l+ a id L T2) ' + p tr 2Q B \> x = a ’ 0 < z < L ( 5 . 1 9 )

_ X a xU 2,l= a U 2d ' X = K ’ ° - z - L ( 5 . 2 0 )

P,a tU M = a ix"m a x ( x " a i|U u ) , 0 < x < o . 0 < z < L ( 5 . 2 1 )

p 3 U = a , x ~ " a ( x " a u ) , cr<x<R , 0 < z < L ( 5 . 2 2 )

1 Z M » X X • » »

Uw =l , z=0 , 0<x<R (5.23)

a U = 0 , x = 0 , 0 < z < L ( 5 . 2 4 ) x li

U + « a T = 0 , x=a , 0<z<L , e=l,2 (5.25) +2 2 x •

- N ta,UU +fi2al T .,+X2(d,U2 / V L T!)-p1''2t! ®2 ’ X=° * °-ZiL (5’26)

-X a U =cdj , x=R , 0<z<L (5.27) 2 x 2 2 22

gdzie <^=d#t(z)/dz , k=l,2.

Minimum funkcji F można poszukiwać jedną z metod opisanych w dodatku.

Korzystając np. z algorytmu gradientu sprzężonego (algorytm Fletchera - Reevsa) zakłada się początkowe przybliżenie dla składowych wektora parametrów {p}2, tzn. pk=p£°\ k=l,2, następnie pierwsze przybliżenie znajduje się metodą gradientu prostego lub metodą najszybszego spadku, a kolejne przybliżenia znajduje się zgodnie ze wzorem:

({p}2)<Ed>= ({ p }2) Cs >+i?.({ 0 )2)<,I> s= l , 2 f3, . . . (5.2B)

gdzie s jest numerem iteracji, ({0}j)<s> kierunkiem poszukiwania, a parametry (stałe dodatnie) są krokami iteracji. Parametry ft* wybiera się

2 (e )

w ten sposób, Ze funkcja F osiąga w kierunku ({©}k) minimum w punkcie ({p}2)<s**\ Kierunek ({0}Z)tB) liniową kombinacją wektorów (grad F)<s) i ({©Ij) . tzn.:

({0}2)U>=-(grad F ) <B)+ls« e> 2)l rU (5.29)

gdzie (grad F)<£* oznacza wartość gradientu funkcji F w danej iteracji, a parametr wyznacza się ze wzoru:

lt=j(grad F ) U>|/1 (grad F)(b_1>5

Dla s=l przyjmuje się, ±e ({0}2) <0>=-(gr ad F)<0>. Również w dalszych iteracjach,co pewien czas, jest stosowana odnowa algorytmu polegająca na

63

-zaprzestaniu budowy kierunków sprzężonych i użyciu kierunku największego spadku.

Proces iteracyjny kończymy, jeżeli spełniony jest Jeden z warunków:

i c3 są. to zadane dostatecznie małe liczby.

5.4. MODEL NUMERYCZNY

V celu stworzenia praktycznej możliwości realizacji obliczeń, podobnie Jak w rozdziale poprzednim, przekształcono układy równań (5.3)+(5.11), (5.14)+(5.20) oraz (5.21)+(5.27) wprowadzając nową bezwymiarową zmienną dla każdej z faz w następujący sposób:

W przyjętym nowym układzie współrzędnych obszary wlewka odpowiadające fazie ciekłej i fazie stałej są prostokątne. "Unieruchomiona" została również granica rozdziału faz, a funkcja a = a (z) opisująca jej położenie stała sie parametrem.

Otrzymuje sie kolejno:

(5.30)

| (grad F)(k>| < ^ (5.31)

(5.32)

gdzie F i F^,(s+l) oznaczają wartości funkcji w kolejnych iteracjach, a

x

~a O<x<o , 0<z<L

y = (5.33)

, o<x<R , 0<z<L

(5.34)

P , 1 3 T + ( y - 1 ) <R - a^{) '}V ³ T ^{) =}

i x js y *

= a3(R-o)'2[ y ( R - o ) + o ] ' m3 y{ [ y ( R - o ) + o ] " 3 yT a} , 0<y<l , 0<z<L ^{( 5 . 3 5}

T =p2 , z=0 . 0<y<l (5.36)

3yTt=0 , y=0 , 0<z<L (5.37)

T1(l,z)=T2(0,z)=T , 0<z<L (5.38)

—X T~ł3 T (1, z )+X (R—o) "*3 T (0,ż) =p y o' ■ 0<z<L (5.39)

l y ł 2 y2 i z

-X (R-of'a T =<x(T -T ), y = l , 0<z<L (5.40)

2 y 2 2 co

R-Ç Sí cr(C ) , n=l,N (5.41)

p S T-T (1.Ç ) . n=l,N (5.42)

n 2 n

dla układu równań i warunków (5.3)+(5.1) i odpowiednio:

p (3 U + y a ' la ' d U, )+a T + y a ~ io ’d T =

*1 z 1.1 y 1,1 y 1,1 z 1 z i y 1y i ( 5 4 3 )

y yI

= a lT ' V " 3 y (yBa yU i t) , 0<y<l , 0<z<L

P 1

1S

^{1 U 2}

,

¹

+ ( r l 1 (R' o) "1 a yu 2 ,1 ) + d iT 2+ t r l 1 ( * - * ) - V a V

a ₂(R-<y)‘2[ y ( R - p ) + a ] ‘" 3 { [y ( R - o ) + o ] " a U } , 0<y<l , 0<z<L

y y 2 , i

(5.44)

UJ4=0 , z=0 , 0<y<l (5.45)

dyUi(=0 , y=0 , 0<z<L (5.46)

U ^ + S a ' ^ T ^ O , y=l , 0<z<L (5.47)

U2)+6i(R-cj) ayT2=0 , y=0 , 0<z<L (5.48)

-X

i

la'lO U

yi » i

( l , z ) + *

i

o'2a 2 T (1,2) 1+X I

y y i 2

( R - a ) 1 3 U ^{y 2.i}( 0 , z ) + + a i ( R - a ) ' 2 a 2 yT 2 ( o , z ) ] = r 2 O s ’ + p 1r 2 e '1 . o<z<l

(5.49)

65

--X2(R-ct) 1ayU2 =ctU2j, y = l , 0<z<L (5.50)

Pl(azUU +yCJ^ a ’ayUl,2)=ai0 y ’ 0<y<3- • 0 < x < L (5.51)

P,[5 U +(y—1) (R—ct) ' » ' d U ] =

» z * (* y 2,2

=a (R-CT)'2(y(R-cr)+a]"'"a { [y(R-o)+<j]ma U } , 0<y<l , 0<z<L

-X1[o'1ayU i 2(l,2)+S2o'2a 2yT 1 (1,2) ]+X2[ (R-a)'‘a yU2^(0,2) + +a2(R-o)'2a yjfT 2 (0,2)]=ptr 2Q a ' , o<z<l

(5.52)

uu =1 » 2=0 - 0<y<l (5.53)

SyUi2=0 , y=0 , 0<z<L (5.54)

u1-2+®2CT lavT ,=0 , y=l , 0<z<L (5.55)

U2_2+e2(R-a)"'ayT2=° , y=0 . 0<z<L (5.56)

(5.57)

-X2(R-CT)'‘ayU2^=ciU2j , y=l , 0<z<L (5.58)

dla zagadnień brzegowych Stefana (5.14)+(5.20) i (5.21)+(5.27)•

Do rozwiązywania występujących w procesie iteracyjnym zagadnień brzegowych można wykorzystać metody różnicowe.

Dla układu równań (5.34)+(5.40) układ rozwiązujący wyprowadza sie wykorzystując równania (5.34)+(5.38) i (5.40), natomiast warunek (5.39) stanowi baze do budowy procesów iteracyjnych, w których wyznacza sie wartości funkcji a określających na każdym poziomie siatki różnicowej w starym układzie współrzędnych położenie frontu krzepnięcia* Podobnie postępujemy w przypadku równań wraź 1iwościowych (5.43) + (5.50) i

(5. 51 )-f(5.58) .

równomiernie siatki różnicowe flj i o krokach Az=L/J oraz Ay =1/1 2

i Ay2=l/I2 odpowiednio, tzn.:

n] j = {yi:yl=iAyl,i=0,I1}x{zj:zJ=jńz, j=0,J}

rf j = {y,:yi=iAy2ti=0,I2}x{zj:zj=jAz, j=0,J}

Równania różnicowe wyprowadza się na podstawie niejawnej gwiazdy czteropunktowej, zgodnie z zasadami podanymi w [5.6,5.7,5.12].

Układ rozwiązujący dla zagadnienia Stefana (5.34)+(5.40) tworzy J układów równań kolejno dla j=l,J o postaci:

+ b . T + > = i = 0 '1 r 1 ( 5 - 5 9 )

= T2,(* >= T ( 5 . 6 0 )

Vj

oj

2 X * ) j 2 . ^ ) b 2,(.) -jAU) 2,<s) ^.(e) _ d 2-(B) (5.61)

U i-J.J U U U i*l.J i.J 2

(.*> = T («> <«> fTi.<*> + 3 T *.(«)_ 4Tt.(.>i

B 2,J [ l.J 0,J 2.J J

(5.62)

gdzie T* j*>=T <!:>(iAy ,JAz), e=l,2, s jest numerem iteracji, a współczynniki odpowiednio wynoszą;

67

'j.j A y 3 (R-a‘* ł)

69

-przy czym T1>0=Pt> 1=0,1,-1, Tj T^^=T, i-0,I2, a^-atJAz), a

a ^ - a (j Az ).

Rozpoczynając proces iteracyjny, przyjmujemy, że a ^ =<7j-i‘ Proces iteracyjny kończymy, jeżeli:

. (s+l) _<®)i

_

l°j -°j I < c

^(5.63)

gdzie c jest zadaną dostatecznie małą liczbą, a jako wynik końcowy za a można przyjąć wynik uzyskany w ostatniej iteracji.

Podobne układy rozwiązujące uzyskamy również dla równań wrażliwoSciowych. Otrzymamy kolejno:

u!;to= °- uuo = O, 1=0,1,

(5.64)

J t-i.J 1.J 1J 1.J 1+l.J 1.J ‘

(5.65)

2J..O^,k*) W A \ , a W . ) i=57f

l.J l-i.J . i,J t,J 1J 1+1,i 1.J 2 ^(5.66)

(5.67)

gdzie k=l,2, j=l,J, U* iAy^, j A z ) ] , e-1,2, 6, 6, (JAz ), s jest numerem iteracji, a współczynniki odpowiednio wynoszą:

% r = °

l.k,(s)

71

2Jc,<«> _

" U '

a 2Az . (i+J ) A y2( R - o j )+Oj ^

p , A y2( R - O j) 2 - i A y2 (R-Oj )+Oj -1

2 (R—O j )Ay^

. .2JlM._

d o.j

-D ,T ^

R — a

2X . U ) , 2* 2

i.j i,j-i

W y ż e j 6kjI k = l , 2 j e s t s y m b o l e m C r o n e c k e r a © j = & t(jAz) n i e z d e f i n i o w a n e d o t y c h c z a s p a r a m e t r y k o l e j n o o z n a c z a j a :

D*T ’ = S T (iAy ,jAz) , 0=1 . 2

t.J y • •

D h ^ a j M i A y ^ , J A z ) , e=l,2

(iAy , j A z ) , e = l , 2

i.j y y • •

W «

0* =

f6 Y

+J.. D V 1

“ ■i P, l. U J Oj i,jj

[ D y ^ j A z - D gj - : ł - i

‘.J P4 t A z ( R - o ) -JJ

T ak Jak i p o p r z e d n i o r o z p o c z y n a j ą c p r o c e s y i t e r a c y j n e ,

= 0 ^ , k = l , 2 . P r o c e s y i t e r a c y j n e k o ń c z y m y , j e ż e l i :

a p o z o s t a ł e

p r z y j m u j e m y , ż e

gdzie s ą to zadane dostatecznie małe liczby, a jako wyniki końcowe za k

0j można przyjąć wyniki uzyskane w ostatnich iteracjach.

5.5. PRZYKŁADY OBLICZEŃ

Prezentowany w pracy algorytm przetestowano na licznych przykładach.

Cz q ś ć obliczeń testujących prowadzono w ten sposób, że najpierw zakładano

znajomość wszystkich parametrów opisujących proces wymiany ciepła miedzy wlewkiem a krystalizatorera, w tyra temperatury zalewania i prędkości wyciągania wlewka, następnie wyznaczano położenia frontu krzepnięcia i spadki temperatury na przekroju naskórka, rozwiązując odpowiednie proste zagadnienie brzegowe. Wyznaczone tak dyskretne funkcje wykorzystywano przy konstruowaniu postulowanych grubości naskórka i spadków temperatury na Jego przekroju w przyjętych punktach kontrolnych. Polega to na zaburzaniu procesem błędu otrzymanych rozwiązań, które wprowadzano następnie jako dane wejściowe do obliczeń dla zagadnienia odwrotnego. Następnie rozwiązywano proponowaną w pracy metodą zagadnienie odwrotne, a wyznaczone w wyniku Ts«-temperatura solidusu). Rozwiązując zagadnienie proste założono, że prędkość wyciągania wynosi w=0.02G [ra/s], temperatura zalewania T =1540 [ C], a krystalizator ma długość L=0+8[ra]. Parametry terraofizyczne materiału oraz współczynnik wymiany ciepła a. miedzy wlewkiem a krystalizatorera przyjęto jak w pracy [5.15], tzn. , że = 29 [W/raK], c^ = c^ = 695 [J/kgKj, Q = 272142 [J/kg], r = r = 7400 [kg/m3], T = 1500 [*C], T„ = 1440 [*C],

1 2 U '

a a e 1700+600 [W/m^].

Rozwiązując zagadnienie proste na modelowany przekształcony zgodnie z (5.33) obiekt nałożono siatki różnicowe o takich samych liczbach węzłów dla fazy ciekłej i fazy stałej odpowiednio po 20x40 węzłów.

Przystępując do rozwiązywania zagadnienia odwrotnego założono, źe są cztery punkty kontrolne oraz że interesuje nas tylko spadek temperatury na przekroju wlewka. Celem obliczeń jest dobór właściwej w sensie obranego kryterium prędkości wyciągania.

Przy powyższych założeniach, przyjmując, że p^=w*°*=0.01, po 17 iteracjach otrzymano, że najlepsza w sensie obranego kryterium prędkość wyciągania wlewka powinna wynosić W=0.0207 [ra/s]. Wyniki obliczeń zamieszczono na rysunku 5.2.

— - rzeczywisty spodek temperatury

'Rys.5.2. Wyniki obliczeń Fig.5.2. The results of computations

Przykład 3 . 2

Rozważano wytwarzany z miedzi wlewek w kształcie walca o średnicy 2R = 0 .18[m]. Rozwiązując zagadnienie proste założono, że prędkość wyciągania

75

-wynosi w-0.003 [m/s], temperatura zalewania T*=1180 [ *C J, a krystalizator ma długość L=0.18 [m]. Parametry termofizyczne materiału dobrano na podstawie literatury (patrz [5.2,5.3,5.10] w szczególności 0^5 4 0, c -420 [J/kgK], Xj=349, X2=360 [ W/raK ], ^=8300, y2=8700[kg/ra3], T=1083[*C] a Q=205[kJ/kg].

Rozwiązując zagadnienie proste na modelowany przekształcony zgodnie z

(5.33) obiekt nałożono siatki różnicowe o takich samych liczbach węzłów dla fazy ciekłej i fazy stałej odpowiednio po 20x80 węzłów.

Przystępując do rozwiązywania zagadnienia odwrotnego założono, że śą trzy punkty kontrolne i że interesuje nas tylko grubość naskórka.

Poszukiwano obu parametrów, tj. temperatury zalewania T § i prędkości

wyciągania optymalnych w sensie obranego kryterium.

Przy powyższych założeniach, przyjmując, że p^0> = w(0> = 0.001, ą

P2°>=T,Wl>=1100, po 32 iteracjach otrzymano, że poszukiwane parametry powinny wynosić odpowiednio W=0.0027 [m/s] a T*=1146 [*C].

Wyniki obliczeń zamieszczono na rysunku 5.3.

Rys.5.3. Wyniki obliczeń Fig.5.3. The results of computations

5.6. PODSUMOWAŃIE

W rozdziale tym przedstawiono algorytm obliczeń prowadzących do wyznaczenia dla wlewków o kształcie płyty lub walca odlewanych sposobem ciągłym temperatury zalewania i prędkości wyciągania zapewniających właściwe ich chłodzenie. W modelu wykorzystuje się informację o przybliżonych grubościach naskórka i spadkach temperatury na jego przekroju w ustalonych punktach kontrolnych.

Prezentowany algorytm oparto na rozwiązaniu odwrotnego zagadnienia mieszanego dla równania przewodnictwa ciepła, polegającego na wyznaczeniu cech różnych typów (T -warunek początkowy, w -współczynnik * równania

różniczkowego), uzyskanego wraźliwościową metodą optymalizacyjną.

Zaproponowana metoda rozwiązania wymagała dołączenia do opisu matematycznego funkcji jakości, kryterium wyboru nieznanych parametrów, którą minimalizowano. Minimum funkcji poszukiwano metodami gradientowymi.

Przeprowadzone obliczenia testujące potwierdziły zadowalającą skuteczność algorytmu.

6. DOBÓR GRUBOŚCI FORMY ODLEWNICZEJ PRZY ZACHOWANIU

W dokumencie Zagadnienia odwrotne w termodynamice procesów odlewniczych (Stron 58-77)