Homogenní lineární pásek s buzením a izolací
Obecná rovnice vedení tepla v materiálu hustoty ρ[kg.m−3
] o tepelné kapacitě cp [J.kg−1
.K−1
] a tepelné vodivosti
λ[W.m−1
.K−1
] má tvar
ρcp
∂T
∂t
= ∇(λ∇T) + q, (1)
kde q[W.m−3
] je objemová hustota pohlcovaného výkonu.
Pravá strana rovnice (1) vyjadřuje tepelnou bilanci elementárního objemu (včetně vnějších vlivů), kdežto levá strana
zahrnuje vyrovnání případného hromadění energie.
Nás bude v rámci budování testovacího řešení zajímat zjednodušený případ jednorozměrného homogenního (λ=konst)
materiálu v ustáleném stavu, potom
∂T
∂t
:= 0 ⇒ λ∆T + q = 0 (2)
představuje mírné zobecnění Laplaceovy rovnice, rovnici Poissonovu.
Za okrajové podmínky zvolíme jeden konec pásku v kontaktu s termostatem (Dirichletova podmmínka, T[N]≡TN )
a druhý pod vlivem např. excitačního laseru (předpokládáme osvit pouze prvního elementu, tedy pouze q[0]=0).
Podstatné je, že přítomnost q na konci pásku nepředstavuje okrajovou podmínku, tu musíme explicitně doplnit –
zvolíme např., že osvětlený konec pásku je tepelně izolován od okolí (homogenní Neumannova podmínka, T[−1]=T[1]).
Budeme úlohu pro jednoduchost řešit v 1D, s rovnoměrnou diskretizací o kroku h. Pro všechny uzly kromě krajních
tedy máme Laplaceovu rovnici ve tvaru
1 ≤ i ≤ N − 1 :
T[i−1] − 2T[i] + T[i+1]
h2
= 0. (3)
Tuto soustavu si můžeme přepsat do tvaru 2T[i]=T[i−1]+T[i+1], a vypsat explicitně v prních dvou vnitřních uzlech:
i = 1 : 2T[1] = T[0] + T[2]
i = 2 : 2T[2] = T[1] + T[3].
Rovnici s nižším indexem můžem použít k vyloučení T[1] z rovnice pro i=2:
i = 2 : 3T[2] = T[0] + 2T[3],
čímž se nám podařilo v rovnici pro dané i nahradit posloupnost třech po sobě jdoucích uzlů uzlem aktuálním, o jedno
vyšším a nultým. Tímto způsobem můžeme postupovat uzel za uzlem, až T[0] vložíme do všech rovnic:
1 ≤ i ≤ N − 1 : T[0] = (i + 1)T[i] − iT[i+1] (4)
Vypíšeme-li nyní získané vztahy (4) v nejvyšších dvou vnitřních uzlech,
i = N − 1 : T[0] = NT[N −1] − (N − 1)T[N]
i = N − 2 : T[0] = (N − 1)T[N −2] − (N − 2)T[N −1],
můžeme první z rovnic využít k elimiaci T[N −1] z rovnice druhé:
i = N − 2 : 2T[0] = NT[N −2] − (N − 2)T[N],
čímž jsme v dané rovnici dosáhli výskytu pouze prvního, posledního a aktuálního uzlu. Takto můžeme postupovat
pro všechny další nižší rovnice, až obecně obdržíme
1 ≤ i ≤ N − 1 : T[i] =
N − i
N
T[0] +
i
N
T[N] (5),
1
takže řešení Laplaceovy rovnice nutně spojuje okrajové hodnoty lineárně.
Hodnotu na ozářeném konci pásku neznéme, víme ale, že ozářený kousek pásku je tepelně izolován. Můžeme tedy
psát
i = 0 :
T[−1] − 2T[0] + T[1]
h2
=
−2T[0] + 2T[1]
h2
= −
q[0]
λ
odkud již přímo s využitím (5) pro i=1 dostáváme
T[0] − T[N] =
q[0]h2
N
2λ
kde na pravé straně stojí přímo hledená tepelná oprava získaná působením laseru. Získaný vztah je však na první
poheld závadný: hledaná oprava podle něj explicitně závisí na počtu elementů diskretizace - takové chování je nepřijatelné.
Jistou pomocí vnést do získaného vztahu fyzikální veličiny, v našem případě délku vzorku l=hN:
T[0] − T[N] =
q[0]hl
2λ
Nyní je ovšem jasné, že krok sítě se nám proti počtu elementů nikdy nepodaří zcela vyvážit a teplený exces stále
nepřijatelně závisí na parametrech simulace. Zde je potřebný fyzikální vhled: výkon laseru při fokusaci na vzorek
zpravidla popisujeme plošným výkonem P soustředěným na povrchu, P [W.m−2
], přičemž samozřejmě pro osvit pouze
prvního elementu P =q[0]h. Nyní již snadno
T[0] − T[N] =
Pl
2λ
(6)
Přítomnost h při použití plného výkonu laseru tak byla vlastně v pořádku: jelikož se u q jedná o objemový výkon,
zmenšování diksretizačního kroku vedlo nutně k celkově menšímu množství pohlcené energie a tepelný exces tedy
musel na h záviset.
Po provedení plné simulace můžeme v rámci kontroly řešení porovnat teplotní exces v místě buzení s hodnotou
získanou analyticky ve vztahu (6). Ukázka reálné simulace (metodou konečných prvků) vedení tepla podle předchozí
úlohy, ve 2D:
Ohřátí vzorku GaAs tloušťky 0.2 mm při Ramanské spektroskopii s laserem He-Ne 50 mW fokusovaného do
stopy cca 4 µm. Krychlová oblast o hraně 120 dílků s výběrem omezeným symetrií. Doba výpočtu cca 1
hod. Barevná škála odpovídá změně teploty oproti údaji termostatu (připojen k dolní části vzorku). Vpravo
výřez v okolí místa excitace.
2