ДЕРЕВООБРАБОТКА И ХИМИЧЕСКИЕ ТЕХНОЛОГИИ
МОДЕЛИРОВАНИЕ ПРОЦЕССА ПРОПИТКИ ДРЕВЕСИНЫ ЖИДКОСТЬЮ
Н.В. ГУБАНОВА, асп. каф. древесиноведения ВГЛТА
В настоящее время в области использования древесины на первый план выдвигаются задачи получения материалов из древесины с заданными свойствами, поскольку древесина является единственным природным возобновляемым материалом. Для широкого развития модифицирования древесины необходимо искать пути решения проблемы математического моделирования процессов модифицирования с целью совершенствования технологических процессов.
Одним из путей решения проблемы математического моделирования в гидромеханике является использование конечно-элементного подхода. В рамках данного подхода в последние 30-40 лет с применением компьютеров неоднократно успешно решался ряд важных теоретических задач гидравлики и гидромеханики. С использованием этого подхода можно создавать модели в высокой степени адекватные, универсальные, быстро реализуемые, легко корректируемые и надежные (позволяющие со стопроцентной вероятностью пройти от начала расчета до конца расчета без тупиков, частых для аналитических методов). В рамках метода конечных элементов сложная система разбивается на большое число однотипных подсистем, взаимодействующих между собой. При этом сложность системы сводится к просчету кооперативных явлений взаимодействия элементов между собой.
Целью данной работы являются:
1. Разработка модели структуры древесины хвойных пород (на примере сосны) и лиственных пород (на примере дуба).
2. Разработка модели пропитывания древесины жидкостью.
Применительно к моделированию процесса пропитки древесины мы предлагаем разбить весь объем движущейся жидкости на большое количество шаров (в трехмерном случае) или кругов (в двухмерном случае), представляющих собой элементы жидкости
и взаимодействующих между собой. Диаметр элементов может составлять порядка 1-10 мкм (в зависимости от решаемой задачи), при этом достаточно хорошо моделируется течение жидкости, и одновременно обеспечивается высокая скорость компьютерных расчетов. Шаровидная (или круговая) форма элементов принята, чтобы добиться изотропии свойств модельной жидкости. Основные свойства жидкости (плотность, модуль упругости, теплопроводность, и т.п.) пересчитываются на один элемент жидкости. Введение различных типов шаров (с соответствующими свойствами) позволяет рассматривать в модели одновременно различные среды: жидкости различных типов, газы, модификаторы, элементы древесины.
Для повышения универсальности модели древесина и жидкость состоят из кругов одинакового диаметра. Однако круги-элементы древесины выдерживаются неподвижными в процессе моделирования. Распределяя в пространстве элементы древесины можно добиться воспроизведения как структуры хвойных пород (сосны, рис. 1, а), так и структуры лиственных пород (дуба, рис. 1, б). Модель воспроизводит основные структурные элементы древесины: сосуды, перегородки, окаймленные и неокаймленные поры, лестничную перфорацию и др.
Для оценочных расчетов достаточно учесть минимум свойств жидкости, в частности, только ее механические свойства. В этом случае необходимо описать ньютоновское движение большого количества тел (элементов жидкости), взаимодействующих между собой вязкоупругими силами.
Принято решение использовать двухмерную модель, так как при этом значительно ускоряются расчеты, что важно на первом этапе.
Состояние каждого элемента-круга i определяется четырьмя переменными: декартовыми координатами его центра (x y) и декартовыми составляющими скорости (v
134
ЛЕСНОЙ ВЕСТНИК 3/2013
ДЕРЕВООБРАБОТКА И ХИМИЧЕСКИЕ ТЕХНОЛОГИИ
vyi). Взаимодействие элементов между собой будем считать вязкоупругим, что позволяет адекватно учитывать возникновение упругости при сжатии жидкости и потери энергии при течении жидкости.
Расчет сил, действующих на элементы, производится следующим образом. Некоторый элемент i испытывает силовое воздействие со стороны каждого из окружающих его элементов j
F^rf+Fj), (1)
7=1
где Fy.. и FB.. - силы упругого и вязкого взаимодействия элементов i и j;
N - общее количество элементов в модели жидкости.
При расчете сил для каждой пары элементов предварительно вычисляется расстояние г., между их центрами S(x, у) и S(x, у) (рис. 2)
rij =yj(xi-Xj)2+(yi-yJ)2. (2)
От способа задания силы между элементами F^r ) зависят свойства жидкости или газа. В простейшем случае можно считать взаимодействие упругим и подчиняющимся закону Гука (рис. 3). Этого достаточно для решения большинства задач о пропитке древесины
При внедрении элементов друг в друга возникает возвращающая сила, пропорциональная величине их внедрения (ветвь графика при 0 < r.. < d_ где d., - диаметр элемента). Для того, чтобы учесть склонность элементов жидкости соединяться между собой и обеспечивать неразрывность объему жидкости, то же самое пропорциональное взаимодействие сохраняется и в диапазоне d,. < r < dQ (где d0 - расстояние, при котором происходит отрыв элементов). При этом попытка удалить элементы друг от друга также вызовет возвращающую силу, стремящуюся вернуть элементы в состояние касания. Таким образом, упругое взаимодействие подчиняется следующему закону
Fv = {cidj-ryXXi-Xj)/r9, если гц < d0;(3) xiJ |0, если rtj > d0;
Fy. =lc(d3~rij)(yi-yj)lrij’ есш rij<do, (4) гуу |0, если Гу > dQ,
где Fy . . и Fy .. - декартовы составляющие "’силы F'1..;
У5
c - жесткость упругого взаимодействия элементов.
Для расчета F6.. выбрана общеприня-
.j
тая прямо пропорциональная зависимость вязкой силы от скорости движущегося в среде тела, при этом введен дополнительный коэффициент (r.. - dA представляющий собой величину взаимного проникновения элементов друг в друга
К=Щ-Уо)^-^У; (5)
Р^Щ-do )0v -vyJ). (6)
где k - коэффициент демпфирования.
Движение элементов в рамках классической механики описывается дифференциальными уравнениями, составляемыми на основе второго закона Ньютона. Для .-го элемента можно записать
Рис. 1. Представление в модели древесины хвойных пород (а) и лиственных пород (б)
ЛЕСНОЙ ВЕСТНИК 3/2013
135
ДЕРЕВООБРАБОТКА И ХИМИЧЕСКИЕ ТЕХНОЛОГИИ
V;
Ji
Рис. 2. Вязкоупругое взаимодействие двух элементов модели
Рис.3. Зависимость силы взаимодействия двух элементов i и j от расстояния между ними
т.
d2v ✓ \ Nn
Ш j=i j=i
где тЭ - масса элемента, t - время,
g - ускорение свободного падения,
Fn - сила действующая со стороны j-го
(8) элемента древесины на i-й элемент
жидкости.
Расписывая силы Fy, F и Fn получим следующие дифференциальные уравнения, описывающие движение i-го элемента:
((
d2 X:
т
э dt2
L=Z
м
j*i
w
(f
c(d3 - Гу) 0, ry>d0
(h~Xj)
\
rv
> ru<d0
\
+ k(rij-d0)(vxi-vXj)
+
NM
M
W
d3+ da (xt ху)
сд(—^—Гу) j
0, ra > d.
> ra < d<
rij
ij ''“ОД
ij ОД
(f
+ kM(ry-d0M)vxi
m
3 d?
Nr,
= -m3g + 2l
м
J*i
(f
\\
C(d3~rii)
0, ry>d0
bi-Vi) .
, r4<d0
Г-
V
A?
+2
j=i
W
,d3+dM 4(y,-yj) ^
сд(---v)--------------» rij<doM
0, r.. > d.
ry
ij ОД
+ k(rij-d0)(vyi-vyj)
+ k„{rirdoa)vyi
+
(9)
где Ад - количество элементов древесины, dд - диаметр элемента древесины, сд - коэффициент жесткости взаимодействия элемента жидкости с элементом древесины,
кд - коэффициент вязкости взаимодействия элемента жидкости с элементом древесины.
Совокупность большого количества уравнений последнего вида для всех Аэ элементов описывает эволюцию жидкости с течением времени. В целом окончательную систему уравнений можно записать следующим образом (см. формулу 10).
Необходимо отметить, что вращение элементов-кругов вокруг их центров в мо-
136
ЛЕСНОЙ ВЕСТНИК 3/2013
ДЕРЕВООБРАБОТКА И ХИМИЧЕСКИЕ ТЕХНОЛОГИИ
гг
Шг
т.
d2xt Д
' Л2 =h
7*1
7=1
Л2У1
W
{х. -*,)
и о
О, г1;. > d0
А+^д A-Xj)
Сд( rij) > r\j<d0g
Z ГЧ
o, rv > d0A
+ k(rlj-d0)(vxi -V .)
Л ^ Л
+ кд(Гц~ dQM ) V^,
ff
3 dt2
= ~Щ S + £
"д
+i
7=1
«Лс, *
7Й,
3 dt2
7=1
7*2
*д
+z
7=1
//
vv
7=1
7*1
vv
^ < d0
0, rv > </0
. ,d3+d„ .. _ з
сд( гц) > r\j< d0ji
z rij
°, »-ц*</од
сА~гг.)^2 Xj\ r2. < d0
+ k(rlJ-d0)(vl-v.)
+ kg(rird0JJ)vyl
0, r2. ></0
Л+</д ^ ,
Сд( - r2 j) ’ r2j<d0g
°. ry - d0g
ff
r2J
+ k{r1j-d0){vxl-vxj)
^ ^ л
+ МГ27 -^дК2
3 dt2
d2y2 £
m,—- = -»b£ + Z,
7=1 7*2 ff
N«
+i
7=1
VV
С(^Э-Г2/)--------. A<£/0
0, r2. > </0
Г27
+ Kr2J-d0)(vy 2-vyj)
VV
rd3+du 0’t-yj) ^ ,
Сд( 9 r2j) > r2j<d0g
Z r2y
0, r2y ^ <У0д
d0JI)Vy2
(10)
у/
'*<э л
dt2
d x
»3 V
—— = 2д 7=1
j*N3
ff
NД
+x
vv
(x„ -X.) ^
c(d3~rN3\j<d0
°. A7 - dO
A
+ k(rN.1j-d0 )(vxN3 ~vxj)
7=1
VV
c Md«
СД' - Ч77
0> rN3j - d0g ff
A
r < И
5 JV од
+ кд(гкэ] d0M)vxNi
m.
d У„3
3 dt2
*'Э
■-Щ g+X
"д
+i
7=1
ff
VV
7=1
7*^э у у
d^ л- d „
,, ,0’N3~yt) , Л
c(d3 rNj) , rNj < d0
°. \Ado
+ k(rN3j-d0)(vN3 -v.)
C r «,
стЛ ^ 'лгэу/ » rN„i v Ma
глгэ7 - d0M
'N3j ' “ОД
"э7
+ ^д(ГДэ7 d0JJ)VyN3
дели не учитывается, при этом практически не вносится погрешность, так как элементы образуют более крупные элементы жидкости (капли, слои). При движении таких крупных фрагментов учет их вращения происходит автоматически - движением элементов друг относительно друга.
При исследовании тех или иных специфических явлений разработанная модель легко позволяет вносить необходимые коррективы. Таким образом, разработанная модель, несмотря на простоту реализации, позволяет добиться высокой степени адекватности, легко корректируется в зависимос-
ЛЕСНОИ ВЕСТНИК 3/2013
137
ДЕРЕВООБРАБОТКА И ХИМИЧЕСКИЕ ТЕХНОЛОГИИ
ти от решаемых задач и гарантированно дает конечный результат при любой постановке задачи.
Разработанная модель в целом представляет собой систему из большого количества дифференциальных и алгебраических уравнений, а также условий включения тех или иных сил. Решение системы дифференциальных уравнений производится численно. Используется модифицированный метод Эйлера-Коши [5], который особенно эффективен при решении дифференциальных уравнений второго порядка. При этом координаты и скорости рассчитываются по формулам вида
x.+l = x. + v . At + a (At)2/2 (11)
v+i = v + a A^ (12)
где x, v, a - координата, скорость и ускорение элемента,
- номер шага интегрирования ( - текущий шаг, i +1 - последующий шаг),
At - шаг интегрирования.
Шаг численного интегрирования At системы дифференциальных уравнений определяется путем многократного проведения экспериментов с последовательно уменьшающимся в 2 раза шагом. Останавливаются на том шаге, после которого результаты моделирования практически не изменяются (изменение составляет не более 1-2 %). Определенный таким образом шаг составил At = 10-6 с и был использован во всех расчетах в данной работе.
Компьютерный эксперимент заключается в просчете проникновения жидкости в структуру древесины в течение некоторого промежутка времени (1000 шагов интегрирования). В процессе компьютерного эксперимента выводятся распределения концентрации жидкости внутри древесины C(x) и C(y), а также распределение давления в направлении пропитки P(x).
В результате проделанной работы можно сделать следующие выводы:
1. В рамках данной работы разработана математическая модель пропитки древесины жидкостью. Модель обладает высокой степенью универсальности и позволяет оценивать эффективность пропитки древесины
различных пород в зависимости от параметров процесса пропитки.
2. С увеличением диаметра сосуда dc скорость движения жидкости vx возрастает приблизительно по следующему закону уж = v0 + Рyjdc, где v0 - скорость движения жидкости в самом малом сосуде; в - коэффициент пропорциональности.
3. Чем меньше диаметр сосуда, тем больше давление жидкости в нем. Из-за этого, несмотря на то, что положение фронта жидкости в сосуде зависит от диаметра сосуда, масса жидкости, находящейся в сосудах, при движении практически не зависит от диаметра сосуда (за исключением очень малых диаметров сосудов поздней зоны).
4. Массовая концентрация жидкости в заполненной жидкостью древесине уменьшается в направлении пропитки (по слабому линейному закону).
Библиографический список
1. Гулд, Х. Компьютерное моделирование в физике / Х. Гулд, Я. Тобочник. - М.: Мир, 1990. - Т 2. -400 с.
2. Зенкевич, О.С. Метод конечных элементов в теории сооружений и в механике сплошной среды / О.С. Зенкевич, И. Чанг. - М.: Недра, 1974. -238 с.
3. Мельников, С.В. Планирование эксперимента в исследованиях сельскохозяйственных процессов. 2-е изд / С.В. Мельников, В.Р. Алешкин, П.М. Рощин. - Л.: Колос, 1980. - 168 с.
3. Советов, Б.Я. Моделирование систем : учеб. пособие для студентов вузов / Б.Я. Советов, С.А. Яковлев. - М.: Высш. шк., 1998. - 319 с.
4. Хеерман, Д.В. Методы компьютерного эксперимента в теоретической физике / Д.В. Хеерман.
- М.: Наука, 1990. - 176 с.
5. Шамаев, В.А. Модифицирование древесины: учеб. пособие для студентов вузов / В.А. Шамаев. - Воронеж: ВГЛТА, 2005. - 197 с.
6. Шамаев, В.А. Модифицирование древесины : учеб. пособие для студентов вузов / В.А. Шамаев.
- Воронеж: ВГЛТА, 2007. - 385 с.
4. Шамаев, В.А. Подшипники скольжения из модифицированной древесины / В.А. Шамаев // Вестник машиностроения. - 2010. - № 7. -С. 62-68.
9. Экштайн, В. Компьютерное моделирование взаимодействия частиц с поверхностью твердого тела / В. Экштайн. - М.: Мир, 1995. - 321 с.
138
ЛЕСНОЙ ВЕСТНИК 3/2013