Добавил:
Опубликованный материал нарушает ваши авторские права? Сообщите нам.
Вуз: Предмет: Файл:

книги / Прогнозирование теплового состояния изделий при эксплуатации в условиях воздействия солнечного излучения

..pdf
Скачиваний:
1
Добавлен:
12.11.2023
Размер:
6.84 Mб
Скачать

1.5.7. Учет граничных условий 2-го рода

Рассмотрим решение уравнения нестационарной теплопроводности при граничных условиях 2-го рода. Для этого найдем производную от 2-го интеграла уравнения (1.20). Используя

(1.24), (1.15) и (1.13а), получаем

 

 

 

 

J (e)

 

= 2π q

T

 

r

dL =

 

 

 

 

 

{T }

{

 

m}

 

 

 

 

 

T

 

 

 

 

 

m

S

2

 

 

 

 

 

 

q(α1

+ α2 + α3 )r dL = 2π

 

 

 

 

= 2π

 

 

q (C11r + C21z + C31 ) +

 

 

S2

 

 

 

 

 

 

 

 

S2

 

 

 

 

+ (C12r + C22 z + C32 ) + (C13r + C23 z + C33 ) r dL =

 

 

= 2π q{ (zi zk ) / (2S1 )r + (rk

 

 

 

rj ) / (2S1 )z +

(1.47)

 

 

S2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

+ (rj zk rk z j )

/ (2S1 )

+ (zk zi ) / (2S1 )r +

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

+ (ri

rk ) / (2S1 ) z + (rk zi

ri zk ) / (2S1 ) +

 

 

 

+(zi z j ) / (2S1 )r + (rj ri ) / (2S1 ) z +

+(ri z j rj zi ) / (2S1 ) }r dL.

1.5.8.Учет граничных условий 1-го рода

Задание граничных условий 1-го рода производится в процессе разработки вычислительной программы путем умножения вектора тепловой нагрузки на большое число порядка

1 1013. Возможен также учет граничных условий 1-го рода путем задания значения коэффициента теплообмена в граничных условиях третьего рода на 3 порядка больше по сравнению с реально возникающими.

1.5.9. Разрешающая система уравнений для элемента и ансамбля элементов

Суммируя все частные производные, получим для одного конечного элемента систему трех уравнений, которую представим в матричной форме. При явной схеме аппроксимации по времени (1.26) эта система уравнений принимает вид

31

1

[c]ek

{T

}m

=

1

[c]ek

[K ]ek {T

}m + {F

}m .

(1.48)

∆τ

∆τ

 

 

k +1

 

 

 

k

k

 

Если взамен (1.26) использовать неявную схему решения задачи по времени

T

 

Tτ Tτ−∆τ

 

Tk

Tk 1

 

 

 

 

 

=

=

=

T T

,

(1.49)

∂τ

 

 

∆τ

∆τ

 

∆τ

 

 

 

то система трех уравнений для одного конечного элемента представляется в форме

 

1

[c]e

+ [К]e {T

}m =

1

[c]m

 

 

}km1+ {F

}m .

(1.50)

{T

 

∆τ

∆τ

k

k

k

k

 

 

 

k

 

Полученные выражения (1.48) и (1.50) представляют собой разрешающие матричные системы трех уравнений для определения неизвестных температур в узловых точках одного конечного элемента в конце временного шага ∆τ .

Для всего ансамбля конечных элементов изделия системы уравнений (1.48) и (1.50) представляются в виде

1

[с]kэ {T

}M

=

1

[с]kэ [K ]kэ {T

}M + {F

}M ,

(1.51)

 

 

 

∆τ

 

 

k +1

∆τ

 

 

 

 

 

k

k

 

 

1

[c]э

+ [К

]э {T

}M =

1

[c]э

 

 

}kM1+ {F

}M .

(1.52)

{T

 

∆τ

∆τ

k

 

 

k

k

k

 

 

 

k

 

Полученные выражения (1.51) и (1.52) представляют собой разрешающие матричные системы уравнений всего ансамбля элементов для определения неизвестных температур в узловых точках дискретной схемы изделия в конце временного шага ∆τ .

В ряде работ по вычислительной математике и методу конечных элементов показано, что численное решение задачи, полученное с помощью системы уравнений (1.51) на основе явной конечно-разностной схемы (1.26), является устойчивым лишь при относительно небольших шагах интегрирования по времени. А численное решение задачи, полученное с помощью систе-

32

мы уравнений (1.52) на основе неявной конечно-разностной схемы (1.49), является безусловно устойчивым при любых шагах интегрирования по времени, т.е. в этом случае шаг по времени ∆τ может быть взят достаточно большим. Поэтому на практике преимущественно используются схема аппроксимации (1.49) и разрешающая система уравнений (1.52). Полученное матричное уравнение (1.52) непосредственно используется для составления вычислительной программы.

1.6.Метод взвешенных невязок получения разрешающей системы уравнений

Разрешающая система уравнений (1.18) и (1.19) для одного конечного элемента выше была получена на основе вариационной формулировки нестационарной задачи теплопроводности путем минимизации функционала, связанного с дифференциальным уравнением теплопроводности. Существуют прямые способы получения решения системы дифференциальных уравнений в частных производных [22], например, метод взвешенных невязок [23, 24]. Преимуществом этого способа является то, что отправной точкой для него служит непосредственно само дифференциальное уравнение и он не нуждается в вариационной формулировке физической задачи. Это очень важно, так как не всякое дифференциальное уравнение имеет вариационный аналог. Метод взвешенных невязок может быть реализован

вдвух вариантах [12]:

1)по методу Галеркина;

2)по методу наименьших квадратов.

Метод Галеркина был предложен в 1915 г. в статье [25] как приближенный метод решения краевых задач математической физики. Заметим, что этот метод является частным случаем метода взвешенных невязок. Применение метода Галеркина в сочетании с конечно-элементной моделью рассматривалось в ряде работ [5, 12, 26–30]. Согласно методу Галеркина приближенное численное решение уравнения нестационарной теплопроводности находится из выражения

33

Wβ Rv dv + Wβ Rs ds = 0,

(1.53)

V

s

 

где wβ – взвешивающие функции; Rv – невязка уравнения нестационарнойтеплопроводности; Rs – невязка граничных условий.

Следуя методу Галеркина, представим неизвестную искомую функцию, т.е. поле температур, в виде произведения двух функций

T (r, z,τ) = Nβ (r, z) T (τ),

(1.54)

где Nβ (r, z) – функция формы, зависящая только от координат; T (τ) – функция, зависящая только от времени.

Отметим, что в качестве взвешивающих функций в методе Галеркина применяются функции формы конечного элемента, т.е.

 

 

Nβ (r, z) = [Ni , N j , Nk ] ,

(1.55)

где Ni

= 0,5(ai

+ bi r + ci z) / S1 ;

 

N j

= 0,5(aj

+ bj r + cj z) / S1 ;

 

Nk

= 0,5(ak

+ bk r + ck z) / S1 ;

 

ai = rj zk rk z j ; bi

= z j

zk ; ci

= rk rj ;

aj = rk zi ri zk ; bj

= zk

zi ; cj

= ri rk ;

ak = ri z j rj zi ; bk = zi z j ; ck = rj ri ;

S1 – площадь элемента.

Функция Т(τ) представляется в виде вектора узловых тем-

ператур

 

T (τ) = {T} = {Ti ,Tj ,Tk }T .

(1.56)

Сучетом (1.31) и (1.32) выражение (1.30) представляется

ввиде

T (r, z,τ) = [Ni , N j , Nk ] {Ti ,Tj ,Tk }T .

(1.57)

34

Производная по времени аппроксимируется выражением

T (r, z,τ)

 

T (τ)

 

∂τ

= T (r, z,τ) = Nβ (r, z)

∂τ

= Nβ (r, z)T (τ) . (1.58)

 

 

Подставляя (1.54) и (1.58) в уравнение нестационарной теплопроводности (1.2) и применяя операцию взвешивания с весовыми функциями формы Nβ , получим первое слагаемое выра-

жения (1.53):

Wβ Rv dV = Nβ [Nβ (cρT ) T (λ Nβ ) Q]dV = 0 . (1.59)

V V

Применив аналогичную операцию к граничным условиям, например, к граничному условию 3-го рода по выражению (1.6), получим второе слагаемое выражения (1.53):

Wβ Rs ds = Nβ λ

Nβ

Tds + Nβ α(Tc NβT )ds = 0. (1.60)

n

s

s

s

Преобразуя второй интеграл в выражении (1.59) по формуле Остроградского и объединяя уравнения (1.59) и (1.60) согласно (1.53), получим интегральное уравнение относительно неизвестного поля температур рассматриваемой конструкции. Затем это интегральное уравнение после ряда преобразований приводится к диф- ференциально-матричномувформе

[C]{

T

} + [K ]{

T

} = {

F

},

(1.61)

которое является разрешающим уравнением и используется для составления вычислительной программы. Подробный вывод разрешающего уравнения (1.61) с помощью метода Галеркина имеется в ряде работ [12, 26–30], и поэтому здесь ограничимся описанием принципиального подхода к решению задачи. Отметим, что уравнение (1.61) аналогично уравнениям (1.48) и (1.50), полученным на основе вариационной формулировки нестационарной задачи теплопроводности.

35

ГЛАВА 2. МЕТОДИКА РАСЧЕТА ТЕПЛОВЫХ НАГРУЗОК ПРИ ХРАНЕНИИ ИЗДЕЛИЙ В ЕСТЕСТВЕННЫХ КЛИМАТИЧЕСКИХ УСЛОВИЯХ

2.1. Виды тепловых нагрузок

Рассмотрим процесс теплового нагружения изделия, находящегося на открытой площадке, под воздействием окружающей среды. При этом наружная поверхность изделия подвергается следующим тепловым нагрузкам [1, 6]:

теплообмен с окружающей средой (воздухом);

воздействие солнечного излучения;

воздействие атмосферных осадков и других атмосферных явлений.

Задача расчета состоит в определении температурных режимов эксплуатации и прогноза сроков сохраняемости характеристик изделий. Наибольшие температуры в изделии возникают

вслучае сочетания максимума температуры воздуха, максимума интенсивности солнечного излучения, минимальной скорости ветра и отсутствия атмосферных осадков. Наименьшие температуры в изделиях возникают при минимуме температуры воздуха, минимуме солнечного излучения, максимальной скорости движения воздуха [1, 31]. Наибольшие температуры возникают при хранении их в условиях тропического сухого климата [3], а наименьшие – в условиях холодного и очень холодного климата [2]. Краткие сведения о температуре воздуха в тропическом климате согласно ГОСТ 24482-80 [3] приведены в табл. 2.1–2.6.

Так как изделия из полимерных материалов обладают значительной тепловой инерцией, даже для их прогрева при расчете тепловых нагрузок необходимо рассматривать несколько суточных циклов, например, не менее трех суток, по истечении которых изменение температуры в изделии выходит на квазистационарный режим. Для определения эквивалентной температуры и прогноза сроков сохраняемости изделий необходимо изучение их теплового состояния не менее чем за год.

Рассмотрим каждую из тепловых нагрузок.

36

2.2. Солнечное излучение

Солнечное излучение характеризуется уровнем суммарного солнечного излучения и состоит из двух компонент: прямого и рассеянного. Прямое солнечное излучение поступает на изделие непосредственно от солнца в виде пучка параллельных лучей. Рассеянное (диффузное) излучение вызвано рассеиванием облаками и воздушной пылью солнечного излучения. Интенсивность солнечного излучения зависит от географической широты, времени года, а также от климатических факторов (наличия облачности, влажности и запыленности атмосферы). Поток солнечной радиации, действующий на наружный слой земной

атмосферы, в среднем составляет q

 

= 1380 Вт/м2, а у поверх-

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

с.н

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ности земли q

 

= 1050 Вт/м2 [32].

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

с.з

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Таблица

2 . 1

Средняя температура воздуха по месяцам и за год, °С,

 

 

 

 

 

в тропическом климатическом районе

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Пункт

1

 

2

 

 

3

 

4

 

 

5

 

6

 

 

7

 

8

 

9

 

10

11

12

 

Год

Асуан

15,7

17,3

21,4

26,4

 

31,3

33,1

 

33,9

33,5

31,3

 

28,6

22,9

17,7

 

26,1

(Египет)

 

 

 

 

Биканер

13,9

17,2

23,3

29,3

 

34,7

35,7

 

33,4

31,5

30,9

 

26,8

20,4

15,2

 

26,0

(Индия)

 

 

 

 

Ханой

16,7

17,2

19,9

23,6

 

27,2

28,8

 

28,6

28,2

27,2

 

24,6

21,2

18,1

 

23,4

(Вьетнам)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Таблица

2 . 2

Суточный перепад температуры воздуха по месяцам

 

 

и за год, °С, в тропическом климатическом районе

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Пункт

 

1

 

 

2

 

3

 

4

 

 

5

 

6

 

 

7

 

 

8

 

9

 

 

10

 

11

 

12

 

Год

Асуан

 

13,3

14,5

16,2

16,7

16,1

16,1

15,0

 

15,0

15,5

 

15,0

 

13,9

13,3

15,0

(Египет)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Калькутта

 

13,8

13,5

13,3

11,8

10,0

7,7

 

6,0

 

 

6,0

 

6,6

 

8,6

 

11,4

13,5

10,2

(Индия)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Ханой

 

6,7

6,2

6,1

7,2

 

8,9

 

7,7

 

7,2

 

 

6,6

 

6,7

 

7,2

 

7,8

 

7,2

 

7,2

(Вьетнам)

 

 

 

 

 

 

 

 

 

 

 

37

Таблица 2 . 3

Абсолютный максимум температуры воздуха по месяцам и за год, °С, в тропическом климатическом районе

Пункт

1

2

3

4

5

6

7

8

9

10

11

12

Год

Асуан

37,8

38,9

43,3

46,6

48,0

50,6

51,1

49,0

47,2

44,4

41,7

37,2

51,1

(Египет)

 

 

 

 

 

 

 

 

 

 

 

 

 

Биканер

31,1

37,2

42,8

47,2

49,4

48,9

46,7

43,3

43,9

42,2

37,2

32,2

49,4

(Индия)

 

 

 

 

 

 

 

 

 

 

 

 

 

Ханой

33,3

34,4

36,7

39,4

42,8

40,0

40,0

38,3

37,3

35,6

36,1

36,7

42,8

(Вьетнам)

 

 

 

 

 

 

 

 

 

 

 

 

Таблица

2 . 4

Абсолютный минимум температуры воздуха по месяцам

и за год, °С, в тропическом климатическом районе

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Пункт

1

2

 

3

4

5

6

 

7

8

9

10

11

 

12

 

Год

Асуан

3,0

1,7

6,0

9,3

11,0

19,0

 

20,2

19,4

17,0

12,2

3,0

 

2,5

 

1,7

(Египет)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Биканер

–2,2

–0,6

–0,6

8,9

16,7

17,8

 

20,6

21,1

13,0

8,3

0,6

 

0,0

–2,2

(Индия)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Ханой

5,6

6,1

11,7

10,0

15,6

20,6

 

21,7

21,1

17,2

13,9

6,7

 

6,7

 

5,6

(Вьетнам)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Таблица

2 . 5

 

Суточный ход температуры воздуха, °С,

 

 

 

 

 

 

в тропическом климатическом районе

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Пункт

 

 

Ночь

 

 

Утро

 

 

День

 

Вечер

Асуан (Египет)

 

 

20,7

 

 

22,6

 

 

33,0

 

 

27,1

 

Биканер (Индия)

 

 

23,5

 

 

19,5

 

 

29,5

 

 

31,7

 

Калькутта (Индия)

 

23,4

 

 

22,2

 

 

29,4

 

 

28,0

 

Ханой (Вьетнам)

 

 

24,5

 

 

21,9

 

 

25,8

 

 

23,8

 

 

 

 

 

 

 

 

 

 

 

 

 

Таблица

2 . 6

Статистические характеристики распределения температуры воздуха в тропическом климатическом районе

 

Средняя

Стандартное

Коэффициент

Коэффициент

Пункт

годовая темпе-

отклонение

 

ратура, °С

за год, °С

асимметрии

эксцесса

Асуан (Египет)

26,1

8,8

–0,12

–0,77

Биканер (Индия)

26,0

9,5

–0,50

–0,42

Калькутта (Индия)

26,4

5,6

–0,64

0,19

38

 

 

 

 

Сведения по солнечному излучению для тропического сухого климата содержатся в ГОСТ 24482-80 [3]. В нем приведены средняя энергетическая экспозиция суммарного солнечного излучения за сутки, полученная для действительных условий облачности (табл. 2.7), и средние значения максимума и минимума энергетической экспозиции суммарного солнечного излучения за сутки (табл. 2.8). Так как для определения полей температур в изделиях требуются текущие значения потока суммарного солнечного излучения, то эти данные необходимо пересчитать на почасовые значения, т.е. получить распределение плотности потока во времени. Методика такого пересчета приведена в работе [6]. На основе статистических исследований временного распределения суммарной солнечной радиации на горизонтальных поверхностях за сутки построены обобщенные графики, связывающие отношения величин часовой суммарной и суточной суммарной радиации с продолжительностью дня или с рассматриваемым интервалом времени. Такой график приведен на рис. 2.1 (сплошные линии) [6]. Цифры над кривыми соответствуют времени от полудня, которому принадлежит данная кривая. Коэффициент К = Qчас / Qсут равен отношению часовой (те-

кущей) суммарной радиации на горизонтальной поверхности к суточной суммарной радиации на горизонтальной поверхности.

Долгота дня определяется по зависимости [6]

 

τдня = (2 /15)arccos(tgϕ tg δ) ,

(2.1)

где ϕ – географическая широта местности (положительная для

северного полушария); δ – склонение, т.е. угловое положение Солнца в солнечный полдень относительно плоскости экватора (положительное для северного полушария).

Склонение Солнца определяется по формуле Купера [6]:

δ = 23,45sin (360(284 + n) / 365),

(2.2)

где n – порядковый номер дня года.

Таким образом, зная долготу дня в функции широты ϕ

и склонения δ и суммарное излучение за сутки, можно определить излучение за один час.

39

40

Таблица 2 . 7

Средняя энергетическая экспозиция суммарного солнечного излучения за сутки Q, МДж/м2 , и отношение ее к возможной Q 0 при условии безоблачного неба по месяцам и за год, %

Пункт

Вид из-

1

2

3

4

5

6

7

8

9

10

11

12

Год

 

лучения

 

 

 

 

 

 

 

 

 

 

 

 

 

Каир

Q

12,4

15,8

19,9

23,7

26,6

28,3

28,1

25,8

22,2

17,7

13,3

11,3

20,4

(Египет)

Q/Q0

0,76

0,76

0,75

0,76

0,81

0,85

0,85

0,83

0,81

0,80

0,78

0,76

0,80

Калькутта

Q

14,9

26,1

29,1

22,4

23,1

17,7

15,9

16,0

15,7

16,2

15,7

14,7

17,6

(Индия)

Q/Q0

0,84

0,85

0,86

0,83

0,80

0,60

0,55

0,60

0,66

0,76

0,85

0,86

0,74

Таблица 2 . 8 Средние из абсолютных максимумов и минимумов энергетической экспозиции

суммарного солнечного излучения за сутки по месяцам, МДж/м2

Пункт

 

1

2

3

4

5

6

7

8

9

10

11

12

Каир

max

13,6

18,4

22,5

26,6

29,6

30,1

29,9

27,7

23,5

19,4

14,3

12,6

(Египет)

 

 

 

 

 

 

 

 

 

 

 

 

 

min

10,9

12,8

18,1

22,0

25,1

26,0

26,2

23,9

20,7

16,0

11,8

10,0

Калькутта

max

15,7

20,3

23,0

25,1

26,8

21,3

19,0

20,2

18,4

17,5

16,8

16,3

(Индия)

 

 

 

 

 

 

 

 

 

 

 

 

 

min

13,4

15,9

18,8

20,6

20,7

13,2

11,4

12,1

12,6

12,7

13,1

13,0

Соседние файлы в папке книги