CC BY
42
ТЕХНОЛОГИЧЕСКИЕ ПРОЦЕССЫ В СЕРВИСЕ
УДК 669: 538.975:517.958
ВОПРОСЫ ОЦЕНКИ НЕКОТОРЫХ ПАРАМЕТРОВ ПРИ МОДЕЛИРОВАНИИ ПРОЦЕССОВ ПРОИЗВОДСТВА И ВОССТАНОВЛЕНИЯ ИЗДЕЛИЙ ИЗ ЖАРОПРОЧНЫХ СПЛАВОВ
Лукашев Евгений Алексеевич, доктор технических наук, профессор, [email protected], ФГОУВПО «Российский государственный университет туризма и сервиса», г. Москва,
Яковлев Николай Николаевич, кандидат физико-математических наук, генеральный директор,
Тураевское машиностроительное КБ «Союз», Московская обл,
Радкевич Евгений Владимирович, доктор физико-математических наук, профессор,
ФГОУВПО «Московский государственный университет им. М. В. Ломоносова», г. Москва,
Коровин Илья Александрович, аспирант,
ФГОУВПО «Российский государственный университет туризма и сервиса», г. Москва,
The authors set up and solve the problem of estimating interference parameters of the crystallization process: the porosity of solid-liquid region, added mass and shrinkage of the alloys, as part of the technological problem of recovering the resource of high-tech products from high-temperature alloys.
В рамках технологической задачи восстановления ресурса высокотехнологических изделий из жаропрочных сплавов поставлена и решена задача оценки взаимовлияния параметров процесса кристаллизации: пористости твердожидкой области, «присоединенной массы» и усадки сплава.
Key words: heat-resistant alloy, gas turbine blades, process of directional crystallization
Ключевые слова: жаропрочный сплав, лопатки газовых турбин, процесс направленной кристаллизации
Задачи сервиса, связанные с восстановлением работоспособности машин, агрегатов и отдельных деталей особенно актуальны для отраслей, в которых предъявляются высокие требования к надежности техники. Особо сложными являются задачи восстановления ресурса и надежности высокотехнологических изделий, к которым, например, относятся литые охлаждаемые лопатки газотурбинных двигателей. Решение таких задач требует основательной теоретической и экспериментальной проработки многих вопросов материаловедения в части получения материала с заданными свойствами при производстве изделия, снижения характеристик материала при эксплуатации из-за накопления дефектов, восстановления характеристик материала, ресурса и надежности изделия при производстве ремонтных работ.
Эффективным способом восстановления структуры и свойств монокристаллических лопаток газотурбинных двигателей из жаро-
прочного сплава ЖС6У, как отмечают авторы [1], после длительной наработки является горячее изостатическое прессование в сочетании с восстановительной термической обработкой. Этот процесс требует особого внимания к соблюдению технологического режима, поскольку из-за пластических деформаций, возникающих в сплаве при горячем изостатическом прессовании, может начаться рекристаллизация, которая приведет к нарушению монокристаллической структуры сплава.
Эксплуатационные характеристики изделия снижаются из-за накопления дефектов и других повреждений в материале. Основной причиной накопления дефектов в литых монокристаллических лопатках является ползучесть. Накопление дефектов из-за ползучести происходит в у-матрице жаропрочного сплава, и эти дефекты затем развиваются в микротрещины и поры. Процессы накопления дефектов приводят к исчерпанию ресурса изделия. При ползучести решеточные дисло-
31
ТЕХНОЛОГИЧЕСКИЕ ПРОЦЕССЫ В СЕРВИСЕ
кации взаимодействуют с границами раздела фаз, образуя сетку на межфазных границах у и у' фаз. При этом энергия межфазных слоев растет, за счет этого микроструктура сплава становится нестабильной. Формирование неравновесных границ раздела при взаимодействии с решеточными дислокациями может привести к появлению микротрещин на межфазных и малоугловых границах зерен в местах скопления дислокаций [2]. К повреждениям также относят и диффузионное огрубление частиц упрочняющей у'-фазы по механизму Лифшица-Слезова [3], а это, в свою очередь, сопровождается снижением долговечности жаропрочного сплава. Вследствие постепенной деградации микроструктуры в ходе накопления повреждений происходит ухудшение механических свойств сплава и эксплуатационных характеристик лопаток. Например, после наработки в 700—800 часов начинается монотонное снижение сопротивления усталости жаропрочного сплава.
Способ высокотемпературной термовакуумной обработки используется для устранения нежелательных изменений структуры сплава лопаток, возникающих после определенного времени эксплуатации газотурбинных двигателей. Целью этой обработки является повышение ресурса лопаток. Рассматриваемый способ представляет собою фактически повторную термическую обработку изделия в вакууме (как защитной среде), которая включает гомогенизацию, т. е. образование твердого раствора с полным растворением вторичных огрубленных выделений у'-фазы, и старение, которое обеспечивает формирование оптимального структурного состояния у'-фазы в матрице. При гомогенизации происходит также отжиг и релаксация дислокаций, накопленных из-за ползучести. Метод высокотемпературной термовакуумной обработки позволяет почти полностью восстановить размеры и морфологию у'-фазы, а также устранить дефекты, образующиеся из-за ползучести, путем отжига. Свойства материала лопаток при этом частично восстанавливаются, но начальных эксплуатационных характеристик достичь не удается.
Это происходит из-за того, что метод высокотемпературной термовакуумной обработки не позволяет повышать плотность материала и снижать его пористость [1]. Кроме того, изменение размера микропор может привести
к крайне нежелательному результату, а именно к увеличению среднего размера пор. В то же время ремонтно-восстановительная обработка должна быть нацелена на устранение таких эксплуатационных повреждений, как поры, рыхлость, внутренние трещины, скопления дислокаций. Эту задачу, по мнению авторов [4], решает операция горячего изостатического прессования.
В [5] рассмотрено влияние давления на область гомогенности промежуточных фаз. Было показано, что равновесные атомные объемы компонентов сплава уменьшаются при увеличении давления. Для некоторых систем наблюдается следующая закономерность: компонент с большим атомным объемом сжимается сильнее, чем с меньшим. При повышенном давлении энтальпия сплавления оказывается меньше, что означает улучшение взаимной растворимости компонентов твердого раствора. Кроме того, всестороннее сжатие приводит к увеличению критических температур фазовых превращений «порядок — беспорядок». Всестороннее сжатие приводит также к увеличению области гомогенности промежуточных фаз, если существование их в данной системе возможно. В [5] отмечается также, что при всестороннем сжатии может возникнуть и новая фаза, которая при отсутствии сжатия является метастабильной. Этот эффект наблюдается при исследовании фазовых диаграмм некоторых двойных сплавов под давлением. При анализе поведения сплавов под давлением с использованием методов псевдопотенциала и статических концентрационных волн, представленном в [5], не учитывалась дефектная структура твердого тела. Однако именно удаление дефектов, дислокаций, микротрещин и пор является целью изостатического прессования в сочетании в восстановительной термической обработкой [1, 4]. Это требует в теоретическом описании учета исходной и приобретенной пористости материала.
Таким образом, поставленные авторами [1, 4] технологические задачи показывают, что существенно понимание и теоретическое описание процессов кристаллизации, перекристаллизации, рекристаллизации с учетом пористой структуры материала. Попытка подобного описания была предпринята в [6]. В этой работе при принятии упрощающих предположений было показано, что пористость и усадка сплава при кристаллизации
32 научный журнал ВЕСТНИК АССОЦИАЦИИ ВУЗОВ ТУРИЗМА И СЕРВИСА 2010 / № 4
Вопросы оценки некоторых параметров при моделировании процессов производства...
являются взаимосвязанными параметрами, которые во многом определяют ход процесса формирования структуры сплава при направленной кристаллизации. Это требует дальнейшего теоретического изучения процессов рекристаллизации, перекристаллизации и, прежде всего, кристаллизации с учетом формирования пористой структуры.
В [6] в качестве одного из блоков математической модели направленной кристаллизации использована модель пористой среды [7, 8], которая была адаптирована к поставленной задаче. Другим блоком модели является модель спинодального распада Кана-Хилларда, позволяющая описывать внутреннюю неустойчивость фазы расплава при направленной кристаллизации. Было показано, что подобное описание позволяет на модели изучать процесс формирования дендритной структуры сплава при направленной кристаллизации. Аналогичное описание внутренней неустойчивости фазы твердого раствора при гомогенизации в рамках технологических задач, поставленных авторами [1,4], требует также учета пористости материала в дополнение к оценкам, сделанным в [5]. Предполагается, что дальнейшее уточнение модели позволит оценить влияние технологических параметров на процесс ликвации легирующих компонентов в междендритном пространстве. Это, в свою очередь, позволит наметить пути повышения жаропрочности дисперсионно упрочненных сплавов.
В модели Био [7, 8], использованной в [6], одним из определяющих параметров является «присоединенная масса», задающая поведение пористого тела при механическом воздействии на него на низких и высоких частотах. Этот параметр характеризует такое явление, как увлечение части жидкости твердым телом при движении в ней, а также торможение жидкости при ее движении относительно твердого тела. Вклад указанных процессов особенно существенен для пористых тел, пропитанных жидкостью.
Оценку «присоединенной массы» можно сделать исходя из экспериментальных данных и достаточно простых геометрических соображений. Подобного рода экспериментальные данные имеются в коллоидной химии и мембранной технологии. По этим данным поверхность скольжения отстоит от геометрической поверхности твердого тела для случая водных растворов на расстоянии порядка 15 нм
(1,5• 10 -8м). Для круглого канала эта часть жидкости будет занимать кольцевой слой. Тогда объем жидкости в поровом канале составит
к, = ^Д
1 4
а часть движущейся жидкости —
г, П -s)h
2 4 .
Объем кольцевого слоя неподвижной жидкости составит п ,
V2 - V, = П h(2d5 - ё2)
4
Относительный объем будет равен
^ , - Ц =ё( 2-ё)=к
V, V, d ^ d) .
Величина пористости твердого тела является характеристикой, которая не несет информации о поровой структуре, так как является отношением всего объема жидкости к полному объему пористого тела. Поровую структуру наиболее полно отражает величина — распределение пор по размерам. Однако эта характеристика измеряется с трудом и всегда с погрешностью, которая зависит от метода измерений и вводимых допущений, необходимых при построении аппарата обработки экспериментальных данных.
Интегральной характеристикой поровой структуры является либо дисперсность, либо удельная поверхность пор пористого тела (чем меньше размер пор, тем дисперсность материала и удельная поверхность пор больше). Для каналов круглой формы она определяется в виде
. _ ndh _ 4
УД_ (nd2/4)h _d
Из этого определения видно, что АУД может
I
быть оценена, как r, где r — характеристический размер поровой структуры. Кроме того, учет формы пор позволяет r трактовать как минимальный характерный размер. Могут быть введены и другие уточняющие положения, например, с использованием весовых коэффициентов, значения которых будут зависеть от классификации, принятой для описания структуры пористости.
С уменьшением размера пор и поровых каналов доля жидкости, принадлежащей рассмотренному выше кольцевому слою, будет увеличиваться. Эту долю неподвижной жидко-
33
ТЕХНОЛОГИЧЕСКИЕ ПРОЦЕССЫ В СЕРВИСЕ
сти можно оценить исходя из массового и объемного балансов. Общая масса пористого тела составляет m=m+m , где m — полная масса пористого тела, m и m — массы твердой и жидкой частей. Средняя плотность пористого тела, пропитанного жидкостью, зависит от плотностей вещества твердого скелета (р) и жидкости (pL) с учетом пористости (в):
Рср
m
V
Vs Ps + VlPl
V V
(1 - в )Ps + вРL ,
где V, Vs и VL — полный объем, объемы твердого скелета и жидкости. Полную массу пористого тела m можно представить также в виде m=m +m + m , где m и m .. — массы подвижной и неподвижной частей жидкости. Из этого следует
Р _ m — VSPs i VL.aPL i VL.bPL .
'cp V V V v .
здесь
VLa + VL.b = Ъ = в V VP.
Доля неподвижной жидкости составляет
5
к = — d
2 -
5
d
и, кроме того,
V V
= к; = !-к.
VL ’ Vl
Из определения пористости в следует
в
V±
V
Vr b + Vl
V
Умножение этого равенства на PL дает
PPl
VLbPL + VL.aPL — kVlPl , (1 - K)VLpL =
V V V V
— k@Pl + (1 - k)PPl .
Таким образом, масса (плотность) жидкости, принадлежащей твердому скелету в виде присоединенной массы, составит
Padd =KPPl (1),
т.е. эта величина определяется пористостью в и структурой пористости, характеризуемой величиной к.
Пористость сплава определяется, прежде всего, режимом кристаллизации и последующих стадий обработки путем гомогенизации и отжига. Одним из определяющих пористость параметров является величина усадки, т.е. различие в удельных объемах жидкой и твер-
дой фаз при кристаллизации. Таким образом, должна быть поставлена задача связать усадку с другими параметрами, характеризующими ход процесса кристаллизации путем описания твердо-жидкой области.
Для получения необходимых оценок следует рассмотреть простую систему, содержащую два компонента A и B, в которой компонент A при кристаллизации переходит в твердую фазу, а B остается в жидкой фазе. Для характеристики этой системы можно ввести некоторые соотношения:
q = -
m
в
тл — отношение масс компонента B (mB) и компонента A (mA ) в исходном распла-
Q = —
ве; ^ ma — отношение атомных масс компонентов В и А (М, М).
Выразим число молей компонентов А и В в жидкой фазе (nA, nB), состав которой при кристаллизации будет меняться (жидкая фаза в ходе процесса обедняется компонентом А и обогащается компонентом B),
Ma
L
т
L
A
L
П , =
n
A
B
m
L
в
Mb
тогда концентрацию компонента А в жидкой фазе можно выразить как мольную долю
L
n
A
С' =
A
L ’
L
П . + n
A
B
соответственно, для компонента B:
Св
n
L
в
L
n
A
+n
L
B
В этом случае СА+ СВ=1. Опустим индексы у концентраций и обозначим с мольную долю компонента А, тогда мольная доля компонента В будет равна 1—с .
Представим концентрацию с как отношение
L
n
A
С =
L
n
где nL — число молей компонентов A и B в жидкой фазе; соответственно,
Ml
L
т
L
n =
здесь mL = mA + mB и ML — условная атомная масса вещества жидкой фазы. Выражение для ML может быть получено из следующего
34 научный журнал ВЕСТНИК АССОЦИАЦИИ ВУЗОВ ТУРИЗМА И СЕРВИСА 2010 / № 4
Вопросы оценки некоторых параметров при моделировании процессов производства...
L
nA
'"A 1L
m
M) ( + mB )/ML) ,■
ML
mL + mi
A B i
,, mL. mL
M.—— +MB—B AM BM„
zcMA + cMB—j-
n л
(2.а);
ML = cM„+cM,
nj in, + nB
B nLJinA + nA
L
n
c
1-c
cMA + cMB-—— = cMA + (1-c)MB
(2.6).
Введем величину y _ m, которая представляет собою отношение массы компонента А в жидкой фазе к общей массе компонента А в системе. Поскольку
mLjMA
mL / ML
L
m
то с можно выразить через величины q и Q , т. е.
= mL ML = mL cMa + (1-c)Mb ; mL MA mL + mB M. ’
c
L
m
A
L
B
[c + {l-c)Q
Выразим mA через y: mL = ymA, тогда предыдущее соотношение будет иметь вид
(последнее выражение представляет собой формулу пересчета текущей концентрации компонента А в фазе расплава из относительной массовой концентрации в концентрацию, выраженную в мольных долях). Так как Q и q — константы (q задается составом исходного расплава), то y = y(c). Таким образом, в ходе процесса кристаллизации у изменяется от 1 до 0, т. е. по предположению в начале процесса весь компонент А находится в жидкой фазе и mA _ mA (у = 1), а в конце процесса, если в жидкой фазе компонента А не остается, то у = О; если же наступает пересыщение по компоненту B, т. е. кристаллизуется эвтектика, то у = у, (О < у, < 1), а значение c соответствует эвтектическому составу.
Поскольку по исходному упрощающему предположению в системе кристаллизуется только один компонент А, а плотность твердой фазы является величиной постоянной, то для последующих вычислений удобно ввести параметр R:
R _
py
(3),
где р — плотность твердой фазы (т. е. плотность компонента А в твердом состоянии, р— const), V — исходный объем расплава (системы в жидком состоянии). Параметр R при заданном составе исходной смеси компонентов А и B является константой. Плотность жидкой фазы — величина переменная, определяется массой и объемом этой фазы
c =
утА
ymA + mLB
(1-c )
Так как компонент B остается в жидкой фазе,
то mB _ mB. Отсюда следует
pL
mL
VL
Далее введем величину, которая представляет собою относительную усадку сплава при кристаллизации
c
утА
ymA +m
'c + (1-c)Q_
у
у + AmB /mA )
c + (l-c)Q
или
c
у
у+q
[c + (1-c)Q
Из этого выражения получаем
cy + qc = cy + {l-c)Qy и у
qc
(l-c)Q
или
с
Qy
q + Qy
V - Vs - VL z _---------
V
(4),
где Vs и VL — объемы твердой и жидкой фаз соответственно.
Определим плотность жидкой фазы (р), поскольку эта величина фигурирует в уравнениях гидродинамики, входящих в модель [6]. Используя определение усадки (4), можем записать VL = (1-Z)V-VS. Далее необходимо получить выражение для Vs. Для этого используем определение плотности твердой фазы ps
S
s m
(p _уу). Поскольку кристаллизуется один
35
ТЕХНОЛОГИЧЕСКИЕ ПРОЦЕССЫ В СЕРВИСЕ
a S L
компонент A, то m = mA - mA, то, используя
L
=m_
определение величины y (y = m), получаем mS = mA - mAy = mA (l - y). Из этих соотношений следует
mS =pSvS = ma(l - y) и Vs = p(l -y).
Используя последнее равенство, найдем отношение
vp: A=
V ' V ~ psV
(1 - y).
Согласно (3) —p = R, поэтому
v
p v
= R(1 - y) и VS=VR(1-y).
Определив объем твердой фазы, для объема жидкой фазы получаем выражение VL = V(1— z—R(1—y)). Для нахождения pL необходимо выразить массу жидкой фазы mL через введенные параметры: mL = mA + mB = mA + mB.
m AL
Из определения y (y = —)
mA
mB
и определения q (q = —) получаем
mA
mL =ymA+qmA = (y + q )mA Отсюда следует
(y + q)mA
L
Lm
p = Vl~V(1 y z - R(1 - y))
В правой части этого выражения присутствует
mA
отношение V , поэтому, используя определение параметра R (r = mA ), можем записать
R(y + q)
Vps
LS
илиР = Р
R(y + q)
1 - z - R(1 - y)
(5).
P =
pS 1 - z - R(1 - y)
Таким образом, плотность жидкой фазы p определяется параметрами R, q (постоянные) и переменными y и z или, поскольку концентрация компонента А в фазе расплава может быть выражена и в виде y, и в виде с, то p определяется переменными с и z. Это, в свою очередь, позволяет согласно (1) и (5) выразить padd через характерные параметры и переменные задачи. Для этого следует выразить соответствующим образом пористость Д По определению в = VL и, так как согласно (4)
Vs + Vl
Vs + Vl = (1 - z)V, а Vl =(1 - z)V - Vs и
vs=vr i-y), то в=1 - R1m (6).
Другой величиной, определяющей наряду с пористостью в значение «присоединенной массы» p , является к. Выше этот параметр был оценен для кольцевого слоя неподвижной жидкости вблизи поверхности твердого скелета пористого тела
к=4{2-4]).
d 1 d
При относительно малой кривизне поверхности пор, т. е. при 8<d , этот параметр может быть представлен в обобщающем виде 5
к=а—, d
где ф — коэффициент формы пор, d — характерный линейный размер (диаметр пор). Тогда для padd, исходя из (1), (5) и (6), можем получить выражение
Padd _ KPPl _ ф~7
d
1-R
i-y
1-z
R (y+q)
1-z-R (1—y)
,R (y+q)
(7).
d 1-z
Если исходить из того, что параметр к может быть представлен как отношение неподвижной части жидкости к общему объему жид-
V
кости в пористом теле, т. е. к = vl и, в свою очередь, объем неподвижного слоя жидкости — в виде VLb =w(L У13 5 (здесь у — число пор в единице объема; у — величина, отражающая дисперсность, пропорциональна 1/d), то для padd получим выражение
Padd _ KPPl _ фб
i i}1-y S 1 -R—---p
[V ((1-z)-R(1-y)
V ((1-z)-R(1-y)
R (y+q) ,
12/3
'фб
Padd тrU3 P
1-z\r (1-z)-R(1-y) ’
s 1
113 ■
V13 \{1-z)-R(y + q)
1-z)-R (y + q) R (y + q)
1-z
_ фб S Padd тл 13 P
(1-z)+R(1-y) ’ R (y+q)
13
V13 {1-z)[{1-z)-R(y + q'i
(8).
Сравнение выражений (7) и (8) показывает, что в (8) имеется множитель
36
научный журнал ВЕСТНИК АССОЦИАЦИИ ВУЗОВ ТУРИЗМА И СЕРВИСА 2010 / № 4
Вопросы оценки некоторых параметров при моделировании процессов производства...
([ - z)- R(y + q) , который может интерпретироваться как поправка к оценкам, принятым в (7). В то же время выражение (7) обладает большей общностью, поскольку параметр d, характеризующий развитость поровой структуры — ее дисперсность, может меняться в процессе кристаллизации. Если в начале процесса плоская граница затравки означает, что значение dотносительно велико, р .. — соответственно, мало, то по мере развития процесса и формирования поровой структуры твердо-жидкой области значение d уменьшается, что дает рост padd. Таким образом, соотношение (7) является рекуррентным соотношением, т. е. значение d должно определяться по ходу процесса кристаллизации. Однако варьированием параметра ц в выражении (8) можно установить, насколько необходимым является использование формы (7).
Приведенный выше алгоритм расчета оценок позволяет установить зависимость усадки z от концентрации в форме у или c по экспериментальным данным для плотности расплава при различных составах (для систем, которые удовлетворяют принятым допущениям: кристаллизация одного компонента и т. д.). Это следует из определения величины z v - Vs - VL , и других соотношений:
z _---------
V
vl L m L L
pl , m _ mA + mB m B 4
1 t L 1L
'Vs \yA - mA У -p yA + mB )
S L
~ lA - ПА
V
далее, используя определения
т. е. z _
_ ml =
y _ m и q m , для z получаем
mA 1ПА
v —1ymAy - у )—1-jmAy+q)
z _---P----------P---------
V .
Представление m в виде m = RpV позволяет для z записать
z _1 - Ry - у)-руК{у + q) (9).
Из (9) следует, что экспериментальная зависимость р = pL(y) дает расчетную зависимость z = z(y).
Таким образом, поставленная в работе задача уточнения оценок, использованных в [6] для описания процесса направленной кристаллизации, решена. Следующим шагом исследования должна быть визуализация развития процесса с учетом полученных оценок в вычислительном эксперименте по типу выполненного авторами [9].
Литература
1. Логунов А. В., Маринин С. Ф., Поклад В. А., Разумовский И. М. Высокотемпературное газостатическое уплотнение монокристаллов жаропрочных никелевых сплавов // Технология легких сплавов, 2005. № 1—4. С. 71-77.
2. Бокштейн С. З., Гинзбург С. С., Кишкин С. Т., Разумовский И. М., Строганов Г. Б. Авторадиография поверхностей раздела и структурная стабильность сплавов. М.: Металлургия, 1987. 272 с.
3. Лифшиц И. М., Слезов В. В. О кинетике диффузионного распада пересыщенных твердых растворов // ЖЭТФ, 1958. Т. 35. № 2 (8). С. 479-492.
4. Логунов А. В., Голованов В. И., Разумовский И. М., Маринин С. Ф., Тихонов А. А., Поклад В. А., Фролов В. С., Быков Ю. Г. Газостатическое прессование — эффективный способ устранения неплотностей в лопатках ГТД // Перспективные материалы, 2004. № 3. С. 12-16.
5. Портной К.И., Богданов В.И., Фукс Д.Л. Расчет взаимодействия и стабильности фаз. М.: Металлургия, 1981. 248 с.
6. Яковлев Н. Н., Лукашев Е. А., Радкевич Е. В. Проблемы реконструкции процесса направленной кристаллизации. // Доклады РАН, 2008. Т. 421. № 5. С. 625-629.
7. Biot M.A. Theory of propagation of elastic waves in a fluid-saturated porous solid. I. Low-frequency range. // J. Acoustical Society ofAmerica, 1956. V. 28. N 2. P. 168-178.
8. Biot M. A. Theory of propagation of elastic waves in a fluid-saturated porous solid. II. Higher-frequency range. // J. Acoustical Society ofAmerica, 1956. V. 28. N 2. P. 179-191.
9. Зайцев Н.А., Рыков Ю. Г. Численный расчет одной модели, описывающей кристаллизацию металлов. 1. Одномерный случай // Препринт ИПМ им. М. В. Келдыша № 72. М., 2007.
37
