1 Background and Objectives

Zirconium-based alloys are extensively used in light water reactor for fuel cladding materials, which form the primary containment barrier. During service in reactor, long term neutron irradiation induced pronounced hardening: a significant increase in yield stress and a strong reduction in strain hardening capacity. Increase in strength is usually attributed to the high density of small irradiation induced loops which act as obstacles against the dislocation gliding. Reduction in strain hardening capacity is related with dis-location channeling where irradiation induced defect clusters are overcome and then annihilated or dragged by dislocations under sufficient stress, these defect-free channels will therefore constitute preferred areas for further dislocation gliding, leading to plastic strain localization at the grain scale and causing the premature failure [1]. Recent statistical TEM observations indicate that dislocation channels mostly occurred in the basal planes [2]. However for HCP zirconium, the principle slip system is the prism plane before neutron irradiation. The underlying mechanisms for such changes in principle slip systems before and after neutron irradiation remain unsolved.

The objectives of this study are to analyze interaction between gliding dislocations and irradiation induced (a) type dislocation loops in zirconium in the atomistic level by means of molecular dynamics method, to determine roles of different slip systems (basal and prism) on localized de-formation of irradiated zirconium, and then propose possible deformation mechanisms of irradiated zirconium under external shear deformation.

2 Simulation Methodologies

Fig.1 shows the simulation model. Two important slip systems are considered: prismatic and basal. For prismatic slip systems, the x axis is along the [1210] crystallographic direction, y along [1010] , and z along [0001]; for basal slip systems, the x axis is along the [1210] crystallographic direction, y along [0001] , and z along [1010]. The Periodic boundary conditions are applied to x and z, constrained boundary conditions are applied to y. In the inner region A are the mobile atoms which are free to move along 3 dimensions during integration process, whereas atoms in the constrained regions B and F are fixed during integration process and can only adjust their positions along the loading direction when the external shear strain is applied.

The model is 32.3 x 20.7 x 22.4nm, containing 640,000 atoms. Loop density is a function of the simulation cell size and calculated to be 6.7 x 1022 m(-3). Dimensional quantities are made unitless in the range of -0.5 to 0.5

One dislocation line and one dislocation loop are inserted into the model. The dislocation line is centered at (-0.25 0 0) with line direction along the z axis, Burgers vector b along the x axis and glide plane y=0. The irradiation induced dislocation loop is centered at (0.25 0 0) with varying Burgers vector, habit plane, nature (SIA or vacancy), position and size.

Four popular interatomic potentials for HCP zirconium metal are available in the open literature, Ackland1995[3] and Ackland2007 [4] are adopted in the molecular dynamics simulations.

Molecular dynamics simulations rely on the numerical integration of Newtron equation of motion for interesting atoms. The explicit time-reversible Velocity Verlet algorithm based on the Stormer-Verlet method[ is adopted. The NVT ensemble (constant number of atoms,volume and temperature) using a Nose-Hoover thermostat is adopted. Simulation temperature is kept at near OK;

During deformation a resolved shear strain increment of 1 x 10(-3) is applied to regions B and F in an opposite direction along x each step, and mobile atoms are relaxed until reaching equilibrium using the molecular dynamics method. The corresponding stress can be estimated from the internal force per unit area: Uyy = FINT/AXzi where FINT is defined as the total force from all the atoms in A and their images on all the atoms in B, Axz is the area of the x-z cross-section of the inner region A.

Information on microstructure and motion of defect zones (including: dislocation, loop, defect clusters and so on) during interaction process is very important to under- stand the deformation mechanisms. Local energy filtering and local geometry approach are used to identify the location of defect zone. The successive three atom planes are highlighted with different colors, and projected in some direction to identify the microstructures of the defect zone. To obtain moving information of the dislocation line and dislocation loop, the relative displacements of interested atoms for any two loading steps are calculated and plotted. In addition, an animation process is made to replay the interaction between gliding dislocation and dislocation loop.

3 Results and Discussion

3.1 Verification of the potential and code

One simple perfect model without either dislocation line or dislocation loop is built and relaxed to be equlilbrium. Based on the interatomic potential Ackland 1995, cohesive energy is calculated to to be 6.243 eV, close to experimental value of 6.32 eV. Then shear deformation simulations for basal and prismatic slip systems are conducted, respectively. The shear stress-strain curve within elastic deformation range is plotted in Figure 2, the gradient of the sheal stress-strain curve is equal to the corresponding shear modulus. The calculated shear modulus is in excellent agreement with the experimental values which further proves this package can calculate interatomic force correctly. i1 can be drawn that this simulation package is correct without bug.

3.2 Core structure of edge dislocation and dislocation loop

Three successive atom layers with red, white and dark from bottom to top respectively near y=0 are projected tc the gliding plane to illustrate the core structure of the edge dislocation. Fig.3 shows the core structure of (a) the prismatic dislocation and (b)the basal dislocation respectively the arrow shows the displacement direction of atoms in the top layer relative to those in the bottom layer. For prismatic dislocation, the displacements are along the [1210] which is parallel to b. Evidently, it does not dissociating in the prismatic plane. The extent of the core structure is about 14b. For the basal dislocation, the displacements are along the [1210] and [1010], the basal dislocation splits into two Shockley partials in the basal plane, both of them have the screw component. The plane sequences in the regions ] and III are dark...white...dark, corresponding to an perfect HCP structure, the plane sequences in the middle region II are dark...white...red, means there is one intrinsic 12 stacking fault with FCC structure, the width of the ribbon of the I2 stacking fault is nearly 9b. This dissociation process can be described by :

1/3[1210]->1/3[0110]+1/3[1100] (1)

The relaxed microstructures of two typical dislocation loops are shown as examples in fig.4 where the red atoms are of higher potential energy, Burger circuit is plotted. Fig.4(a) shows the (0001) atomic section of the cluster of 11 SIAs initially in 1/3/11201 configuration, after 250ps of relaxation, the single SIAs are of the basal crowdion(BC) configuration, and all SIAs reside in the {1010} plane, the distortion of the cluster perimeter is not significant. According to the Burgers circuit analysis, the Burgers vector is along 1/3[1120] with one 60ー rotation away from the normal of its habit plane. Fig.4(b) shows the (0001) atomic plane of the cluster of 49 vacancies initially in the {0110} plane, after 250ps of relaxation, these vacancies reside in the {0110} plane with the Burgers vector 1/3[1210]. These results are compatible with the experiment observations[5]

3.3 Interaction between (ai) edge prismatic dislocation and (a2) vacancy loop

The interaction between one prismatic dislocation and (a2) clusters of 49 vacancies is calculated. Figure 5 shows the 2-dimensional snapshots of different interaction stages, respectively. Red atoms represent they are near the defect zone with high potential energy. Figure 5(a) shows the initial structure before deformation. The dislocation is along the z in the gliding plane y=0, with Burgers vector (al) , and the dislocation loop habits in the {al} plane with Burgers vector (a2). The poles of the dislocation loop along the dislocation line z are interacting with the the basal plane, and dissociate into two Shockley partials. Figure 5(b) shows the dislocation line has start gliding. It will leave one row of vacancies at its original position which are sessile. Due to the attraction from the vacancy dislocation loop, dislocation line bend towards it. Under external shear stress, the dislocation loop also moves, the left parts are moving along the (預2) whereas the right parts are moving along the (a2), changing the loop shape. The left parts move a little slowly than the right parts. Figure 5(c) shows the dislocation line is interacting the dislocation loop. Be-cause the shape of the dislocation loop is changed so that most of the dislocation loop are not in the same plane as the gliding plane of the dislocation, then the dislocation does not interact with the dislocation loop, except of two poles of dislocation loop along the dislocation line z. Therefore the segments of the dislocation line far away from the dis-location loop can continue to move easily whereas the segments near the dislocation loop bows out strongly because they are pinned by that two poles of the dislocation loop. Figure 5(d) shows the dislocation line is breaking out of the dislocation loop and continues to glide.

In addition, the z projection of interaction between different dislocation segment and dislocation loop, it is clarified that the segments near the dislocation loop glide in the different plane whereas the segments far away from the dislocation loop always glide in the same plane, there must be one junction between these two segments which can not glide in the gliding plane. Such junctions make it difficult for the dislocation line to break out of the dislocation loop.

3.4 Interaction between (a1) edge dislocation and (al) SIA loop

(10)The interaction of one basal edge dislocation line and one (al) cluster of 11 SIAs is calculated. Figure 6 shows the Y projection of atoms in the x-z planes near y=0 at different deformation stages, under such view, the section is the basal plane, atoms are plotted in red, black and white depending on their potential energy form high to cohesive energy. Figure 6(a) shows the structure of dislocation line and dislocation loop at the initial stage when there is no external stress. The edge dislocation line dissociates into two Shockley partials in the basal plane between them the stacking fault forms with FCC structure. The width of the stacking fault is about 15b. Two segments of dislocation loop intersecting with basal plane are calculated by drawing Burgers circuit. The SIAs is in the crowdion position alone Burgers vector (a1) and habit in the (a3) plane. Fie(12)ure 6(b) shows the microstructure of dislocation line anc dislocation loop under shear strain 2.0%. Under externa: shear stress, partials leave from original position and art moving toward the dislocation loop. The leading partia bows out strongly due to the repulsion from SIA dislocation loop, the trailing partial moves a little faster than the leading one which makes the stacking fault zone narrowed., for example, the segment near the dislocation loop are narrowed to be lb, whereas the segments far away from the dislocation loop to be 9b. Figure 6(c) shows the structure of dislocation line and dislocation loop under shear strain 8.0%. Segments of partials far away from the dislocation loop may continue to gliding, the dislocation loop has little effect on the gliding of them. Segments of partials neat the dislocation loop are cutting the dislocation loop, due to the strong repulsion from the SIA dislocation loop and the moving forward trailing partials, the stacking fault between these two partials becomes narrowed and is lock by the dislocation loop. Figure 6(d) shows the structure of dislocation line and dislocation loop under shear strain 10.0%. With increasing shear strain, locked part of partials are re-leased from the SIA dislocation loop, the leading partial can continue to glide as before interaction. The trailing partial bows out by the attraction from SIA dislocation loop, but will continue to glide as before interaction.

According to Saada and Wshburn process[6], when a al dislocation gliding in the basal plane interacts with a al loop, kinks are produced during interaction processes. Be-cause the kinks are in the same plane as the dislocation line, they have no special effects on the dislocation gliding.

4 Conclusions

This study is the first molecular dynamics simulation concerning dislocation deformation mechanism in irradiated HCP zirconium. The primary simulation findings are summarized as following:

・Edge dislocation dissociates into two Shockley partials in the basal plane, however it does not dissociate in the prismatic plane;

・The dislocation slipping between planes with wide distance may contribute to plastic deformation and is considered to be physically reasonable.

・Two versions of interatomic potential Ackland1995 and Ackland2007 are compared in the MD simulation, only the latter can give correct stacking fault in the basal and prism plane.

・For interactions between (a1) prismatic edge dislocation and (al) dislocation loop, kinks are produced during interaction processes. The dislocation line glides in the same plane as before interaction. The junctions are glissile.

・For interactions between (a1) prismatic edge dislocation and (a~13) dislocation loop, jogs are produced during interaction processes. Segments of the dislocation line cutting through the center of the dislocation loop glide in the different plane with those of far away from the dislocation loop. The junctions are sessile.

・For interactions between (al) basal edge dislocation and (a11213) dislocation loop, kinks always are produced during interaction processes. The dislocation line glides in the same plane as before interaction. The junctions are glissile.

Based on the simulation results, dislocation breaking ability is proposed to determine how easily it is for the dislocation to break out of the dislocation loop and subsequently form the dislocation channeling. Two cases over three in the prismatic plane are of lower dislocation breaking ability, in contrast all cases in the basal plane are of hinger dislocation breaking ability. Such differences are attributed to the productions of interactions between dislocation line and dislocation loop, and lead to low possibility of dislocation channelings in the prismatic planes. These results are in good agreement with TEM observations of dislocation channels in the irradiated Zircaloy.

Figure 1: Schematic visualization of molecular dynamics simulation model

Figure 2: Strain-stress curve for the two crystals oriented for basal or prism glide

Figure 3: Relaxed core structure of the edge dislocationare

Figure 4: Relaxed configuration of (a) the 11 SIAs cluster initially in 1/3 [1120] configuration and (b) the 49 vacancies cluster initially in 1/3[1210] configuration

Figure 5: Interaction between (a1) edge dislocation and (a2) vacancy loop: Y projection of atoms in x-z planes near y=0

Figure 6: Interaction between (a1) edge dislocation and (a1) SLA loop: Y projection of atoms in x-z planes near y=0

ジルコニウム合金は、軽水炉の燃料被覆管材料として利用されており、その使用条件下での変形機構を解明することは、原子力発電の安全性確保のための基幹をなすものである。本論文では、ジルコニウム合金での変形挙動が中性子照射を受けることによって変化する機構を、分子動力学法を活用した計算機シミュレーションによって、転位と欠陥集合体の相互作用に基づくミクロな立場から検討し、燃料被覆管の破損限界評価のための基礎的な研究を行っている。

本論文は7章で構成されている。

第1章では、本論文における研究の背景並びに目的を述べている。燃料被覆管の健全性を確保する工学的要請、及び変形と破損に至るミクロからマクロまでの物理的現象を明らかにすることの意義をまとめ、金属ジルコニウム中の転位とその欠陥集合体との相互作用に起因する変形挙動の機構解明を計算機シミュレーションによって行うことの重要性を議論している。これに基づいて、すべり運動をする転位と転位ループの相互作用を明らかにし、中性子照射を受けたジルコニウム合金において主すべり面が変化する機構を解明することを本論文の目的と設定したことを述べている。

第2章は、開発した計算機シミュレーション手法について、詳述している。異方性を持つ稠密六方晶であるジルコニウムの変形を取り扱うために、分子動力学法に基づいて一般的な結晶構造と転位及び欠陥集合体を対象とした変形を取り扱いうる新たなコードを開発した成果、及びシミュレーション結果の可視化法について論じている。

第3章は、開発した分子動力学シミュレーション手法の検証を行っている。特に原子間ポテンシャルの有効性について多数のポテンシャルを比較して議論を行い、Ackland教授との継続的な議論に基づいた新しいポテンシャルがジルコニウムの転位挙動と変形を扱う上で最適であることを導いている。また、六方晶金属における転位芯の構造を議論し、柱面にすべり面を持つ2種類の刃状転位の安定性と可動性について、積層欠陥エネルギーの評価等に基づいて検討するとともに、底面がすべり面となる刃状転位の構造を明らかにしている。

第4章は、六方晶の柱面上にすべり面を持つ転位と各種の格子間原子ループ及び空孔型ループとの相互作用を系統的にシミュレーションした結果から明らかにしている。

また第5章は、底面上のすべり系における転位とループとの相互作用について、同様に系統的に計算機シミュレーションを行った結果について、論じている。

これらの結果、柱面をすべり面とする刃状転位とループはバーガースベクトルの向きが同じ場合、キンク及び空孔の形成はあるものの、転位のすべり運動に対してループは大きな障害とはならないことを明らかにしている。さらに、バーガースベクトルの向きが異なる柱面上の刃状転位とループの相互作用は、転位のすべり運動の障害となる不動のジョグを形成することを明らかにした。また、底面がすべり面である転位がループと相互作用する場合は、ループのバーガースベクトルにかかわらず常にキンクが形成されるが、キンクはすべり運動が可能であって、転位運動の障害とはならないことを明らかにしている。

第6章は、ループのサイズが変形挙動に与える影響について、サイズを数種類変化させたシミュレーション結果に基づいて、議論している。軽水炉の使用温度で想定されるサイズより小さなループであっても、新しいAcklandポテンシャルに基づいたシミュレーションでは吸収や引きずり現象は確認されず、ループのサイズ効果は小さいことを見いだしている。

第7章は結論であり、これまでの計算機シミュレーションの結果をまとめて、総合的な議論を行っている。ジルコニウム合金の変形挙動は、転位チャネリングの形成過程と転位のループセグメントの吸収容易性とともに、転位がループを突破する困難さによって評価できることを示し、この困難さを評価する指標を提示することに成功している。この指標は、転位とループの相互作用の結果、可動な欠陥集合体が形成するか否かを示すものであって、中性子照射によってループが形成したジルコニウム合金では、六方晶底面でのすべりが卓越するとの実験事実を説明するとしている。

以上を要するに、本論文は、原子炉の安全性評価の上で枢要な燃料被覆管材料としての使用環境におけるジルコニウムの変形機構を最新の原子間ポテンシャルに基づいた分子動力学シミュレーションによって明らかにすることに成功しており、システム量子工学、特に原子炉燃料システム工学の発展に寄与することが少なくない。よって本論文は、博士(工学)の学位請求論文として合格と認められる。