Моделирование процесса сушки бревна

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

Рассмотрим процесс при котором бревно идеально обдувается с боков, торцы обмазаны и ничего не пропускают (чтобы не трещали). Условия сушки постоянные - температура и влажность заданы
tн, град;
Hн, %(отн.); соответственно.
сушку будем рассматривать влажного бревна, свежеспила к примеру, но не мокрого (полежавшего в воде). Т.е. начальная влажность на уровне 80-100%, а никак не предельные 200%.
Эти 200%, а также средняя влажность сосны 88% в свежесрубленном состоянии упоминаются тут - http://stanokwood.ru/articles/8/

Кто не читал -
Дрова и теплотворная способность. КПД и теплоутилизаторы. Все точки над И.
напомню, что под влажностью древесины обычно понимается отношение массы воды к массе сухого остатка, а точнее даже к массе лабораторно высушенной древесины, поэтому эта величина может быть более 100% (на 1 кг дерева приходится более килограмма воды).

Дерево влажное, а не мокрое. Это свидетельствует о том, что оно не водонасыщено и воды могло бы быть больше. Это говорит о том что в дереве присутствуют микропустоты и вода связана с деревом. Это соображение нам нужно чтобы ограничится диффузионным механизмом сушки. Т.е. выходит пар, а сама жидкая вода не движется. Вообще говоря предположение спорное. И был бы рад сведениям о подвижности жидкой воды в древесине.
При этом механизм связывания скорее физический (малые энергии взаимодействия), что подтверждается близкой к 2,3МДж/кг теплотой испарения воды из дерева при влажностях вплоть до 10% влажности.
До 10% влажности сушку рассматривать не будем, поэтому теплоту испарения положим постоянной и равной этим 2,3МДж/кг.

Далее, будем считать что бревно находится при постоянной температуре. Это довольно точное предположение так как времена установления температурного режима - часы, а сушки недели/месяцы.

Паропроницаемость древесины положим также постоянной. Предположение грубовато, так как луб, заболонь, ядро имеют разную плотность и очевидно по разному пропускают пар, но таких данных мне не найти, при этом основным испытуемым материалом обычно является заболонь, а именно она основная толщина сушки древесины. А избыточные барьерные свойства луба в какой-то мере компенсируются повышенной проницаемостью ядра при рассмотрении всего процесса.
Сама величина паропроницаемости для сосны 0,06 мг/(м*ч*Па) согласно http://www.dpva.info/Guide/GuideMatherials/BuildingMaterials/BuildingWood/

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

Еще нам потребуется граничное условие на наружной части бревна. Пренебрегая кинетикой фазового перехода (на самом деле она есть, но данных на нее не найти, поэтому расчет получится скорее оценкой снизу во времени сушки), а также принимая во внимание рассмотренный случай идеального обдува влажность приповерхностного слоя дерева должна быть равновесна с воздухом.
Для получения этой приповерхностной влажности воспользуемся подходом химпотенциалов из упомянутой мной статьи про дрова.
Для этого сначала аналитически аппрокисмируем зависимость равновесной влажности воздуха от влажности древесины при 20град. Добавим при том предельные точки, ну и точку 28 заменил на 26:
0%->0%(отн)
3% -> 5%(отн)
4% -> 10%(отн)
6% -> 24%(отн)
9% -> 44%(отн)
12% -> 63%(отн)
15% -> 73%(отн)
20% -> 85%(отн)
26% -> 93%(отн)
200%->100%(отн).
Решено подбирать функцию вида: H=100%*x^n/(a+x^n),
здесь n, a - свободные параметры, x - влажность древесины в %.
Итоговый результат давший довольно хорошее совпадение (перебор руками):
n=2,5
a=300.

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

Теперь что же с температурой. Воспользуемся тем же графиком зависимости равновесной влажности воздуха от влажности древесины. Пусть при 10%
10град =>48%(отн)
20град =>51%(отн)
30град =>54%(отн)
40град =>57%(отн)
50град =>60%(отн)
60град =>63%(отн)
70град =>67%(отн)
80град =>71%(отн)
90град =>75%(отн)
100град =>79%(отн)

Хорошо аппроксимируется экспонентой 45,596exp(0,0055x)
Если взять 4% влажность дерева то экспонента 7,7827exp(0,0134x)
При этом коэффициенты в показатели коррелируют между собой: 0,0055*10/4=0,0138, стало быть температурную зависимость можно записать в виде H=A*exp(0,055*T/x)
x тут влажность древесины в %-ах.

Итого, у нас есть функция раз
H=100%*x^2,5/(300+x^2,5),
и два
H=A*exp(0,055*T/x)
Определим A:
100%*x^2,5/(300+x^2,5)=A*exp(0,055*293/x)
Отсюда:
A=100%*x^2,5/(300+x^2,5)*exp(-0,055*293/x)
Итоговая формула:
H=100%*x^2,5/(300+x^2,5)*exp(0,055*(T-293)/x)=100%*x^2,5/(300+x^2,5)*exp(0,055*(t-20)/x)
T - абсолютная температура в Кельвинах, t - обычная в цельсиях.

Эта формула позволит нам определять влажность в порах древесины в процессе сушки, а ее решение вида x=f(H,T) влажность древесины на границе с воздухом
Найдем эту функцию.
выразим x из полученной ранее формулы
x=(300H/(100-H))^(2/5)

Формула с температурой тоже легко переворачивается:
x=0,055*T/ln(H/A)
подставляя 20 град находим A:
0,055*293/ln(H/A)=(300H/(100-H))^(2/5)
ln(H/A)=0,055*293*(300H/(100-H))^(-2/5)
A=H/EXP(0,055*293*(300H/(100-H))^(-2/5))
Подставляем и получаем итоговый результат:
x=(T/293)*(300H/(100-H))^(2/5))

Кстати, полученные формулы H(x,T) и x(H,T) универсальны и подходят для любого типа древесины, как и говорилось в предыдущей статье про дрова.

Жидкая вода на три порядка более плотная чем пар. Этот факт приводит к тому что количество влаги в паровой фазе внутри древесины незначительно по сравнению с жидкой(связанной) фазой. На основании этого будем считать что перенос пара квазистационарный, жидкость же при том иммобилизована и может перемещаться только в виде пара.

В результате возникают следующие уравнения.
Тут долго разъяснять, поэтому пока просто формулы.
2Pi*r*dr*L*Ro*dx/dt=d(dm/dt)
-d(dm/dt)=2Pi*(r+dr)*L*(J+dJ)-2Pi*r*L*J=2Pi*r*L*dJ+2Pi*dr*L*J

-r*dr*Ro*dx/dt=r*dJ+dr*J
-Ro*dx/dt=dJ/dr+J/r
J=-Mu*dP/dr=-Mu*d(H*Pнас)/dr=-Mu*Pнас*dH/dr

Ro*dx/dt=Mu*Pнас*(d2H/dr^2+(1/r)*dH/dr)

Итоговая система
dx/dt=(Mu*Pнас/R0)*(d2H/dr^2+(1/r)*dH/dr)
H=(x*100)^2,5/(300+(x*100)^2,5)*exp(0,00055*(t-20)/x)
(dH/dr)[r=0]=0 => (dx/dr)[r=0]=0
x[r=R]=x(Hн,Tн)
Задача замкнута, нужно разобраться с приведением в единую размерность и можно решать.

Cоставил конечно-разностную сетку. Сделал алгоритм в excel-e.
Вот некоторые результаты.
Исходные данные, бревно ~сосна, диаметр 30см, исходная влажность 80% сушим воздухом 50% влажность при разных температурах.
Динамика сушки на картинках.
20 град -
50 град -
90 град -
Влажность конечно не в %, а в долях, проценты в 100 раз больше.

Сушится дни и дни. На улице (в составе дома в т.ч.) так годы фактически.

Сам расчет тут. https://drive.google.com/file/d/0B4FXn2v9OuYXdmMxREZpZm9NMHc/view?usp=sharing

В качестьве PS. Любая модель это ряд допущений. Одними из сильных допущениями здесь на мой взгляд являются постоянство коэффициента паропроницаемости и "нещеляние" бревна. Второе ладно, просто рассматриваются только мягкие режимы сушки с получением идеального результата.
Что же касается паропроницаемости, поясню следующим образом. В нормативах указано что паропроницаемость сосны/ели поперек волокон, что в условиях А, что Б 0,06. Речь при том о влажностях 15-25%. Вода сама по себе абсолютно паропроницаема, т.е. доля объема в дереве занятая водой паропрозрачна. Например если бревно 80%-й влажности, то оно содержит ~40% от объема воды и паропроницаемость такого бревна должна вырасти на сообразную величину.
0,06 соответсвует средней влажности в 20%, т.е. воды там около 10%. Оценку паропроницаемости для остальных влажностей исходя из случайного (не кластерного) расположения воды в объеме можно оценить следующим образом:
0,055/(1-Объемная_доля_воды).
Т.е. при сушке бревна с 80% до 20 мы имеем дело с меняющейся паропроницаемостью около 0,09. Т.е. самый начальный период сушки будет проходить несколько эффективнее, однако эта же эффективность будет нивелирована неабсолютным обдувом бревна в процессе сушки.
В среднем же за период сушки паропроницаемость будет не сильно отличаться от выбранного значения и составит около 0,065-0,07. Т.е. тут имеет место быть источник погрешности на уровне 10%. В совокупности с влиянием усушки, неоднородностью древесины, грубостью сетки и т.д. общую точность расчета можно оценить на уровне 20%.