Full Text

Введение Технологии аддитивного производства предоставляют широкий спектр для научно-технических исследований [1-4]. Аддитивное производство позволяет создавать цельные трехмерные объекты практически любой геометрической формы. Большинство существующих в настоящее время аддитивных технологий изготовления металлических изделий основаны на их использовании в качестве исходного материала для формирования изделий порошковых систем [4-8]. При этом существенно сужается круг конструкционных материалов, из которых возможно получение готовых изделий высокого качества. В частности, формирование изделий из большинства алюминиевых сплавов, а также сплавов активных металлов, таких как титановые и магниевые сплавы, приводит к повышенной пористости материала полученного изделия с существенным снижением его механических характеристик [9-12]. Производительность формирования деталей из порошковых материалов в существующих аддитивных технологиях является крайне низкой, что практически исключает перспективы применения этих технологий для производства крупногабаритных изделий. Кроме того, несмотря на использование в этих технологиях оборудования, которое по точности в большинстве случаев превосходит традиционное оборудование для обработки металла, при порошковом 3D-синтезе точность изготовления изделия оказывается невысокой в связи с усадкой материала, составляющей, по данным исследователей, от 8 до 30 %. Технология гибридного изготовления изделий объединяет лучшие характеристики аддитивного формирования заготовки и последующего механического удаления материала в процессе создания металлических изделий [13-15]. Процесс может реализовываться на одной платформе с гибридным послойным нанесением проволочного материала и обработкой с помощью станка с ЧПУ и является оптимальным для изготовления крупногабаритных деталей низкой и средней форм сложности. Формирование крупногабаритных металлических изделий с использованием оплавления проволочного материала позволяет избавиться от проблем, связанных с низкой производительностью существующих методов, высокой стоимостью применяемого оборудования, ограниченностью типов применяемых материалов, обусловленной традиционным использованием в качестве исходного материала для аддитивного формирования изделий порошковых систем, оплавляемых мощным тепловым источником [13-14]. Математическое описание процессов, сопровождающих аддитивное производство, представляет собой масштабную мультидисциплинарную проблему, решение которой требует привлечения понятий, законов и методов теории тепломассообмена, материаловедения, металлургии, термомеханики и других наук. Как элемент современного интеллектуального цифрового производства, аддитивные технологии неразрывно связаны с широким использованием конечно-элементного компьютерного моделирования средствами специализированных программных комплексов, реализующих математические модели процессов. Процесс изготовления детали методом аддитивного производства сопровождается сложными термомеханическими явлениями, приводящими к формированию полей технологических остаточных напряжений и возможному короблению детали [16-23]. Причина появления внутренних напряжений в выращиваемом объекте связана с существенным пространственно-временным неоднородным распределением температурно-конверсионных полей. Механизм возникновения остаточных напряжений вызван несовместностью неупругих деформаций [24], в первую очередь деформаций температурной усадки при охлаждении, структурной усадки вследствие протекания фазовых превращений (кристаллизации расплава), отличающейся деформационной историей различных точек материала из-за неоднородности температур, температурных градиентов и скоростей температур. Зачастую для устранения возникающих при изготовлении современных деталей геометрических дефектов проводят предварительную натурную отработку технологии, выбор определенного расположения детали в установке, технологического режима с целью минимизации дефектов, что является серьезным препятствием для организации полноценного цифрового производства. Подход, связанный с построением и использованием математических и компьютерных моделей термомеханического поведения детали в процессе аддитивного производства, позволяет без проведения натурных экспериментов на этапе проектирования параметров технологического процесса и параметров конструкции прогнозировать качественные и количественные характеристики напряженного состояния и коробления будущей детали, а также исследовать эффективность возможных мер по снижению остаточных напряжений и деформаций (термообработка). Известны исследования, направленные на моделирование процессов послойного синтеза [25-28] и оптимизацию термического цикла [29], оценку влияния параметров процессов на изменение формы готового изделия [30]. Возможность моделирования технологического процесса сплавления представляет значительный интерес и является резервом для оптимизации технологических режимов изготовления деталей, разработки управляющих программ, минимизации дефектов и повышения качества изготовления сложных деталей [31-32]. Численное моделирование процесса оплавления проволочного материала рассмотрено в работах [33-34], проводится сравнение с процессом производства методом оплавления для реального образца. Расчеты ведутся методом конечных элементов, осаждение металла моделируется выделением и активацией элементов, для вычислений используется коммерческий продукт Marc и Mentat. Модель учитывает термические и деформационные процессы, наблюдаемые в образце, и данные вычислительного эксперимента хорошо согласуются с натурным экспериментом. В представленной статье рассмотрены особенности численного моделирования процесса аддитивного наращивания изделий методом сплавления с использованием пакета ANSYS, проведен анализ возможностей этой широко применяемой в ПНИПУ программы для заявленной цели. 1. Объект исследования Целью данной работы является верификация созданного численного алгоритма решения трехмерной задачи производства металлических изделий с использованием наплавки проволочных материалов. В качестве эталона взяты результаты эксперимента по наплавлению проволоки на металлическую подложку с использованием дуговой наплавки в защитных газах вольфрамовым (неплавящимся) электродом, изложенные в работе [33]. На рис. 1, а показана рассматриваемая пластина с размерами 275×100×12 мм. На верхнюю поверхность наплавляется металлическая стенка длиной 165 мм, поперечное сечение которой показано на рис. 1, б. Рис. 1. Общий вид изделия после эксперимента с системой креплений и датчиков (а); размеры наплавляемой зоны в сечении, мм (б) [33] Fig. 1. General view of the product after the experiment with the fixing system and sensors (a); the dimensions of the weld area in the cross section (b) [33] Процесс дуговой наплавки производится в два этапа. На первом для предварительного разогрева пластины-основания осуществляется проход горелкой туда-обратно со скоростью 5 мм/c и охватом площади, в 1,5 раза превышающей площадь основания наплавляемой стенки. Потребляемая мощность горелки 2,7 кВт, диаметр пятна 6 мм. Последующая наплавка производится в 10 слоев по непрерывной траектории (на обратном ходе процесс не прерывается). Материал - никелевый сплав Inconel 718 [35]. 2. Математическая постановка Принимая во внимание малость деформаций и пренебрежимо малое диссипативное тепловыделение, можно разделить краевую задачу нестационарной теплопроводности и краевую задачу термомеханики о напряженно-деформированном состоянии (НДС), которые в такой постановке являются несвязанными. Для их решения применяется технология «умерщвления» и последующего «оживления» (Elements Birth and Death в ANSYS) части материала, первоначально в модели отсутствующего, а затем возникающего в процессе наплавки. При этом в качестве расчетной рассматривается область, занимаемая уже готовым изделием. Непрерывное наращивание материала производится дискретно, на каждом подэтапе расчета, соответствующем «оживлению» очередной подобласти из «мертвых» элементов, решается краевая задача теплопроводности и термомеханики, причем результат решения предыдущего подэтапа служит начальными условиями для последующего. На k-м подэтапе решения постановка краевой задачи нестационарной теплопроводности по определению полей температур в области с границей Sk с учетом принятых гипотез включает в себя [36]: уравнение теплопроводности , , (1) где , , - соответственно теплоемкость, теплопроводность и плотность неоднородно легированного материала; - удельная мощность внешнего источника тепла; граничные условия , (2) где первое слагаемое правой части описывает конвективный теплоперенос, а второе - излучение (закон Стефана-Больцмана); - коэффициент черноты; - постоянная Стефана-Больцмана; - коэффициент теплопередачи; - температура окружающей среды; - внешняя единичная нормаль к границе охлаждаемого тела; начальные условия , , (3) где - начальное распределение температуры для -го подэтапа; - температура, определенная на конец предыдущего. В данных соотношениях учитывается, что на протяжении подэтапа область исследования остается неизменной. Здесь через и обозначены зоны, занимаемые «живыми» и «мертвыми» элементами соответственно. При этом теплофизические свойства материала в зоне «мертвых» элементов подвержены деградации: , , Несвязанная квазистатическая краевая задача механики деформируемого твердого тела (МДТТ) с учетом несущественности вклада массовых сил на k-м подэтапе включает [37-38]: уравнения равновесия , (4) где - тензор напряжений; геометрические соотношения Коши . (5) где - вектор перемещений; - тензор полных деформаций; граничные условия в перемещениях (6) и напряжениях (7) где - части границы с заданными перемещениями и нагрузками соответственно. Термомеханические параметры материала в зоне «мертвых» элементов исключают физическую нелинейность, являются идеально упругими с деградированными значениями: где - тензор четвертого ранга упругих констант материала. Общая система уравнений краевой задачи МДТТ включает также определяющие соотношения. Для описания вязкоупругопластического поведения исследуемого сплава использована модель «Ананд» [41], включенная в перечень физических моделей пакета ANSYS. Модель имеет следующий вид: (8) (9) где - тензор четвертого ранга упругих констант материала; - тензор упругих деформаций; - тензор полных деформаций; - тензор температурных деформаций; - тензор вязких деформаций; - коэффициент температурного расширения материала; - температура начала отсчета температурной деформации; - единичный тензор второго ранга. Скорость вязкой деформации вычисляется по формуле , где - девиатор тензора напряжений; - среднее напряжение; - эквивалентное напряжение, (символ «:» означает скалярное произведение, ранее введенный символ «··» означает тензорное произведение); - эквивалентная скорость вязкой деформации, связанная с остальными параметрами соотношением (10) Модель включает также эволюционное уравнение (11) предполагающее возможность деформационного упрочнения и разупрочнения. Здесь ; ; - интенсивность напряжений; - сопротивление деформации (скалярная функция упрочнения); - величина насыщения функции упрочнения; - универсальная газовая постоянная; - энергия активации; - абсолютная температура; , , , , , , - эмпирические коэффициенты. С учетом особенностей поведения элементов, «оживляемых» по технологии, используемой в ANSYS, соотношения (8) преобразуются к виду (12) где - полная деформация, вычисленная к концу -го подэтапа (актуальное состояние на момент «оживления» элемента является для него естественным). 3. Алгоритм численной реализации 3.1. Решение задачи теплопроводности Алгоритм, используемый для расчета в конечно-элементном пакете ANSYS температурных полей при численном моделировании процесса электронно-лучевого наплавления, предполагает выполнение следующих расчетных процедур. 1. Создание конечно-элементной модели, включающей разделенные на отдельные горизонтальные слои объемы, занятые будущим изделием, а также платформу-основание с соответствующими теплофизическими свойствами. 2. «Умерщвление» (команда EKILL) части элементов изделия, отсутствующих в реальном процессе наращивания до его начала. 3. В цикле по наращиваемым слоям, начиная с нижнего: 3.1. Задание на верхней границе слоя условий конвективного теплообмена по формуле (2); 3.2. В цикле по подэтапам прохода зон наплавления очередного расчетного слоя: 3.2.1. Удаление на нижней границе части слоя под зоной ранее заданных условий конвективного теплообмена, если таковые имеются; 3.2.2. «Оживление» (команда EALIVE) всех элементов k-й зоны; 3.2.3. Нагрев части элементов зоны распределенным по объему источником тепловой энергии (см. формулу (1)) в течение времени воздействия луча ; 3.2.4. Удаление источника тепла и ожидание в течение временного интервала, равного времени ожидания . 4. Ожидание до полного (частичного) охлаждения системы. Время воздействия луча с учетом возможного перекрытия траекторий движения горелки вычисляется по формуле (13) где - площадь пятна; - скорость горелки; - расстояние между соседними треками. Время ожидания для k-й зоны наплавления определяется выражением (14) где - площадь зоны. Удельная мощность тепловыделения на поверхности пятна нагрева наплавляемого материала определяется по формуле (15) где - мощность электронного луча; - доля мощности, расходуемая на нагрев. Удельная объемная мощность теплового источника, используемая в уравнении (1), принимает вид (16) где - объем наплавляемой зоны; - толщина слоя. Толщина наплавляемого слоя задается исходя из наблюдаемой в эксперименте общей высоты наплавленной части изделия и фактического количества слоев. Теплофизические свойства материала взяты из работы [39]. При этом в соответствии с представленными в статье данными плотность ρ = 8170 кг/м3 и теплопроводность λ = = 30 Вт/(м·К) полагались постоянными. Температурная зависимость энтальпии представлена на рис. 2. Теплоемкость материала, входящая в уравнение теплопроводности (1), вычислялась по формуле . Следует отметить, что перепад на кривой в диапазоне 1500-1600 К на рис. 2 обусловлен выделением тепла при кристаллизации сплава. Доля мощности луча, расходуемая на нагрев K, в расчете принималась равной 0,2. Рис. 2. Зависимость энтальпии H(T), кДж/кг от температуры T, К для сплава Inconel 718 Fig. 2. The dependence of enthalpy H(T), kJ/kg on temperature T, K for alloy Inconel 718 На рис. 3. представлена конечно-элементная модель задачи. При решении использованы 20-узловые изопараметрические элементы Solid279. Рис. 3. Дискретная модель задачи о наплавлении Fig. 3. Discrete model of the problem of deposition 3.2. Решение задачи определения напряженно-деформированного состояния конструкции Дискретная модель идентична изображенной на рис. 3. Перед началом расчета температурный элемент Solid279 заменяется на Solid286, использующий перемещения в качестве степеней свободы, и добавляются термомеханические свойства. Алгоритм аналогичен приведенному в п. 3.1. 1. «Умерщвление» (команда EKILL) части элементов изделия, отсутствующих в реальном процессе наращивания до его начала. Задание граничных условий в перемещениях (симметричное закрепление по плоскости разреза, закрепление по вертикали в зонах захватов (рис. 1, 5) и в центральном нижнем узле). 2. В цикле по наращиваемым слоям, начиная с нижнего: 3.2. В цикле по подэтапам прохода зон наплавления очередного расчетного слоя: 3.2.1. «Оживление» (команда EALIVE) всех элементов k-й зоны; 3.2.2. Приложение ранее вычисленных для данного момента времени температур к узлам модели в течение времени воздействия луча ; 3.2.3. Приложение ранее вычисленных для момента конца времени ожидания температур и выдержка секунд. 4. Чтение температурного поля для времени частичного охлаждения системы, расчет НДС. 5. Приложение температуры окружающей среды, высвобождение от захватов, расчет НДС. Идентификация выбранной модели «Ананд» (8)-(11) для исследуемого материала проводилась по данным эксперимента на растяжение с заданной скоростью при различных температурах, представленного в работе [40]. Искомые константы модели определялись в несколько этапов методом Нелдера-Мида в пакете Matlab. При этом минимизировалась относительная невязка между расчетными и экспериментальными значениями напряжений на диаграмме (рис. 4). Рис. 4. Сравнение экспериментальных (точки) и расчетных значений напряжений в одноосных испытаниях на растяжение при температурах 720, 850, 900, 1150 ºС (сверху-вниз) Fig. 4. The comparison of experimental (points) and calculated values of stresses in uniaxial tensile tests at temperatures of 720, 850, 900, 1150 ºС (from the top down) В результате найдены следующие значения материальных констант: = 2,0447e+05 Па, = 2,0864e+03, = 0,3335, = 2,5008e+04, = 0,4363, = 4,0352e+06, = 3,2148e+07 Па, = 0,2273, = 0,4461. Как видно из рис. 4, относительная погрешность расчета не превышает 10 %. Материал полагался изотропным, модуль Юнга и коэффициент Пуассона 153 ГПа и 0,32 соответственно не зависят от температуры. Среднее значание КЛТР на температурном интервале до температуры солидуса = 6,5 10-6 К-1. Линейная усадка на интервале между температурами солидуса и ликвидуса составляет 0,3 % и введена в температурную зависимость . Эффекты упрочнения-разупрочнения, сопровождающие процесс наращивания, заложены в модель «Ананд». 4. Результаты расчета и обсуждение На рис. 5 показано расположение датчиков температуры в эксперименте из работы [33]. При этом датчики 2, 3, 5, 6 размещались на верхней поверхности пластины-основания, а 1 и 4 - на нижней. Начальное направление движения горелки - справа-налево. На рис. 6-8 представлены графики эволюции температуры в данных точках. Как видно из рисунков, качественно полученные данные хорошо согласуются с экспериментом. Наихудшее количественное совпадение наблюдается в точках на верхней плоскости основания, расположенных близко к зоне нагрева (точки 3 и 5, рис. 8). Кроме того, кривые содержат осцилляции, при этом температура в некоторых точках опускается ниже температуры окружающей среды. В целом экспериментальные графики лежат выше расчетных. Наличие скачков можно объяснить близостью зоны нагрева и связанных с этим больших температурных градиентов, существенно искажающих результат на достаточно грубой (см. рис. 3) сетке. Более высокий уровень нагрева в эксперименте объясняется в том числе неучетом тепловой радиации, исходящей от горелки, распространяющейся далеко за пределы расчетного размера пятна. В других точках измерения (см. рис. 6, 7) погрешность расчетных данных существенно ниже. Максимальная абсолютная невязка не превышает 15 °С. На рис. 9 в качестве иллюстрации показаны температурные поля на момент завершения предварительного прогрева пластины-основания (см. рис. 9, а) и на момент достижения конца четвертого слоя (рис. 9, б). Рис. 5. Схема размещения и нумерация датчиков температуры [33] Fig. 5. The layout and numbering of temperature sensors [33] На рис. 10 показаны результаты решения задачи механики деформируемого твердого тела (МДТТ), изображены остаточные перемещения в центральных сечениях пластины. Из рисунка видно, что остаточное коробление конструкции прогнозируется с приемлемой точностью. Предсказана форма пластины в виде чаши с высотой подъема краев около 1 мм. Наименьшая точность расчета - в поперечном сечении, где относительная невязка подъема достигает 20 %. Рис. 6. Температура, °С, в точках 1 (красные линии) и 4 (черные линии). Тонкие кривые - эксперимент, жирные - расчет Fig. 6. The temperature (°С) in points 1 (red lines) and 4 (black lines). Thin curves show the experiment, thick curves show the calculation Рис. 7. Температура, °С, в точках 2 (красные линии) и 6 (черные). Тонкие кривые - эксперимент, жирные - расчет Fig. 7. The temperature (°С) in points 2 (red lines) and 6 (black lines). Thin curves show the experiment, thick curves show the calculation Рис. 8. Температура, °С, в точках 3 (красные линии) и 5 (черные). Тонкие кривые - эксперимент, жирные - расчет Fig. 8. The temperature (°С) in points 3 (red lines) and 5 (black lines). Thin curves show the experiment, thick curves show the calculation а б Рис. 9. Распределение температуры, К: (а) 291 с; (б) 506 с Fig. 9. The distribution of temperature, K: (а) 291 s, (b) 506 s На рис. 11 представлена картина распределения остаточных характеристик в конструкции. Рис. 10. Относительные остаточные перемещения, мм, в продольном (синие линии) и поперечном (черные) сечениях плиты на нижней поверхности. Тонкие кривые - эксперимент, жирные - расчет Fig. 10. Relative residual displacements, mm in the longitudinal (blue lines) and transverse (black lines) sections of a plate on the lower surface. Thin curves show the experiment, thick curves show the calculation а б Рис. 11. Вертикальные перемещения в конструкции после остывания и освобождения, м (а); интенсивность остаточных напряжений, Па (б) Fig. 11. Vertical displacements in the construction after cooling and release, m (a); the intensity of residual stresses, Pa (b) Из рис. 11, а видно, что доминирующим является продольный прогиб, поперечный примерно на порядок ниже. Наибольшая интенсивность остаточных напряжений наблюдается в зоне стыка наплавляемого материала с платформой. Выводы 1. Разработан и реализован в виде программы на APDL пакета ANSYS алгоритм расчета нестационарных температурных полей и напряженно-деформированного состояния конструкции в процессе ее создания дуговой 3D-наплавкой проволочных материалов. 2. По экспериментальным данным других авторов проведена идентификация теплофизических и термомеханических параметров материала - никелевого сплава Inconel 718. 3. Верификация разработанной модели показала приемлемую точность полученных результатов. 4. Для повышения точности определения температуры имеет смысл в дальнейшем усложнить математическую модель процесса, учтя такие факторы, как лучистый перенос тепловой энергии сварочной дуги на поверхность изделия и термоконвекционные явления в зоне расплава.

About the authors

O Yu Smetannikov

Perm National Research Polytechnic University

D N Trushnikov

Perm National Research Polytechnic University

P V Maksimov

Perm National Research Polytechnic University

M L Bartolomey

Perm National Research Polytechnic University

A V Kovyazin

Perm National Research Polytechnic University

References

Statistics

Views

Abstract - 403

PDF (Russian) - 166