AP05_5.vp 1 Introduction A necessary condition for extensive development of nu- clear energy is to maintain and increase the safety of nuclear fuel. At the same time, there is a need to operate NPPs (Nu- clear Power Plants) economically. The importance of nuclear fuel as a basic part of each NPP is more and more apparent, so our work aims at a more correct description of the behaviour of cylindrical fuel rods under more stringent operational conditions, including high burn up. This implies that modern fuel rod computer models must be able to predict with greater reliability the margin for cladding tube integrity loss under normal operation and operational transients. One of the critical moments that occur during fuel duty is the pellet – cladding contact (and also modified pellet – pellet contact after a thermo-mechanical interaction). The difficulty of modelling the contact event can be shown by a very simple example: one core in one loading of one reactor contains ~104 fuel rods, and each fuel rod contains ~102 fuel pellets. Where, when, how and in what detail do we have to predict and model it? 2 Mathematical description of the problem – unilateral contact Mathematical modelling of unilateral contact problems is a new area in applied mathematics. The finite element method (FEM) is used. The method is based on the fact that the zone of unilateral contact of two solid bodies need not be defined a priori; the correct definition is one of the results of the solution of the problem. The classical analysis of this problem, started by Hertz in 1896, was limited to simple geometries. The age of high- -speed computers has brought a qualitative change into the analysis of contact problems. On the basis of suitable discre- tization, by means of finite differences or FEM, the problems can be solved approximately even for complex geometrical situations and boundary conditions. A. Signorini originally formulated it as a technical problem in 1933 for the case of unilateral contact of an elastic body with an ideal rigid and smooth base (known as Signorini’s problem). From the mathematical point of view the unilateral contact problem performs the application of variational inequalities, which to- day form the basis of convex analysis and make a connection between mathematical physics problems and optimisation problems. Unlike classical variational equations, variational inequalities are solved not on the whole function space but on some closed part of it. The problem tends to seek the condi- tional extreme of the functional of the deformation energy. These problems are nonlinear, due to the existence of mate- rial nonlinearities, and also because the set on which the solution is being found can alter during the solution process. The mathematical details can be found, e.g., in [1, 2, 3]. 3 Systems used for calculation 3.1 COSMOS/M The problem has been solved by means of the COS- MOS /M [4] software system. It is a complete, modular finite element system, which includes for instance modules to solve linear and nonlinear static problems in addition to problems of heat transfer. Developing a reliable model capable of pre- dicting the behaviour of structural systems represents one of the most difficult problems to face the analyst. FEM provides an effective vehicle for performing these problems due to its versatility and the great advancement in its adaptation to computer use. The success of a finite element analysis de- pends largely on how accurately the geometry, the material behaviour and the boundary conditions of the actual problem are defined. One has to take into consideration that all real structures are nonlinear in some way. In addition, the unilat- eral contact problems are nonlinear due to their optimization substance as was mentioned above. The unilateral contact forms the basis for adequate pellet-cladding modelling and it was used in all calculated problems. © Czech Technical University Publishing House http://ctn.cvut.cz/ap/ 5 Czech Technical University in Prague Acta Polytechnica Vol. 45 No. 5/2005 Physical and Numerical Difficulties in Computer Modelling of Pellet-Cladding Contact Problems for Burned-Up Fuel M. Dostál, J. Zymák, M. Valach The importance of fuel reliability is growing due to the deregulated electricity market and the demands on operability and availability to the electricity grid of nuclear units. Under these conditions of fuel exploitation, the problems of PCMI (Pellet-Cladding Mechanical Interaction) are very important from the point of view of fuel rod integrity and reliability. Severe loading is thermophysically and mechanically expressed as a greater probability of cladding failure especially during power maneuvering. We have to be able to make a realistic prediction of safety margins, which is very difficult by using computer simulation methods. NRI (Nuclear Research Institute) has recently been engaged in developing 2D and 3D FEM (Finite Element Method) based models dealing with this problem. The latest effort in this field has been to validate 2D r-z models developed in the COSMOS/M system against calculations using the FEMAXI-V code. This paper presents a preliminary comparison between classical FEM based integral code calculations and new models that are still under development. The problem has not been definitely solved. The presented data is of a preliminary nature, and several difficult problems remain to be solved. Keywords: PCMI, FEM, contact, pellet, cladding. After the first thermoelastic cases, all further problems of pellet-cladding interaction were solved as nonlinear problems in these ways: geometrical nonlinearity and the thermal, mostly nonlinear, dependencies of all used material proper- ties. Heat transfer with inner sources included was always solved as the initial solution for the subsequent mechanical problem. 3.2 FEMAXI-V The second used system was FEMAXI-V [5] – a Japanese light water reactor fuel analysis code. It deals with a single fuel rod and predicts thermal and mechanical response of a fuel rod to irradiation, including FP (fission product) gas release. The thermal analysis predicts rod temperature distribution on the basis of pellet heat generation, changes in pellet thermal conductivity and gap thermal conductance, and (transient) change in surface heat transfer to coolant, using radial one-dimensional geometry. The mechanical analysis performs elastic/plastic, creep and PCMI calculations by FEM. The FP gas release model calculates diffusion of FP gas atoms and accumulation in bubbles, release and increase of internal pressure in the rod. FEMAXI-V consists of two main parts: one for analyzing the temperature distribution, thermally induced deforma- tion, FP gas release, etc., and the other for analyzing the mechanical behaviour of the fuel rod. In the thermal analysis part, the temperature distribution is calculated as a one-dimensional axisymmetrical problem in the radial direction, and with this temperature, such tempera- ture-dependent values as FP gas release, gap gas flow in the axial direction, and their feedback effects on gap thermal conduction are also calculated. In the mechanical analysis part outside the thermal analy- sis part, users can select either analysis of the entire length of a fuel rod, or analysis of one pellet length. In the former, the axisymmetrical finite element method (FEM) is applied to the entire length of the rod; in the latter, the axisymmetrical FEM is applied to half the pellet length (for a symmetrical reason), and the mechanical interaction between pellet and cladding, i.e. local PCMI, is analysed. In the mechanical analysis, the magnitude of the pellet strain caused by thermal expansion, densification, swelling and relocation is calculated first, and a stiffness equation is formulated with consideration given to cracking, elastic- ity/plasticity and creep of the pellet. Then, the stress and strain of pellet and cladding are calculated by solving of the stiffness equation with bound- ary conditions corresponding to the pellet-cladding contact mode. When PCMI occurs and the states of the pellet-clad- ding contact change, calculation is re-started with the new boundary conditions of contact from the time when the change occurs. 4 Validation case postulated in NRI The validity of the new models in COSMOS /M is exam- ined by a comparison with FEMAXI-V. The power history is shown in Fig. 1 and covers 6 days: linear increase of the linear heat rate from 0 to 300 W/cm, 300 W/cm, linear decrease to 150 W/cm, 150 W/cm, linear increase to 300 W/cm and 300 W/cm. One hour time steps were used in both codes. 5 2D COSMOS/M results against FEMAXI-V results 5.1 Input data models synchronization The dimensions of a pellet for this validation case are selected to represent a VVER-440 fuel rod: inner diameter 6 © Czech Technical University Publishing House http://ctn.cvut.cz/ap/ Acta Polytechnica Vol. 45 No. 5/2005 Czech Technical University in Prague Fig. 1: Linear heat rate vs. time 1.4 mm, outer diameter 7.68 mm, height 10.0 mm, and clad- ding thickness 0.69 mm. A radial gap of 20 �m at room temperature is assumed, considering pellet swelling at a high burn up. The non-linear, thermal dependent material properties of the pellet and cladding were applied from FEMAXI-V into COSMOS /M. 2D meshing of the pellet and cladding is given in FEMAXI-V (see Fig. 2), so it was rearranged in the same way in COSMOS /M. A detail of an 8-node element with four integration points is shown in Fig. 3. 5.2 Checking the thermal field to achieve comparability Temperature is an important factor in predicting PCMI. The reason is that almost all processes and effects (thermal expansion, creep, etc.) are temperature dependent. Temperature in time follows the linear heat rate. The central temperatures are quite well comparable (see Fig. 4 and Fig. 5). The peak values are almost the same: 994 °C FEMAXI-V and cca 1245 K (972 °C) COSMOS /M. The pellet outer surface temperature and the cladding inner and outer surface temperatures are also consistent. This is a good base point for further comparison. © Czech Technical University Publishing House http://ctn.cvut.cz/ap/ 7 Czech Technical University in Prague Acta Polytechnica Vol. 45 No. 5/2005 Fig. 2: Mesh of the pellet and cladding in FEMAXI-V Fig. 3: Integration points for a rectangular 8-node element Fig. 4: Temperature (°C) vs. time (days) from FEMAXI-V. 8 © Czech Technical University Publishing House http://ctn.cvut.cz/ap/ Acta Polytechnica Vol. 45 No. 5/2005 Czech Technical University in Prague Fig. 5: Temperature (K) vs. time (seconds) from COSMOS/M a) b) c) d) Fig. 6: a) Eq. stress at the outer clad. surface, b) Eq. stress at the inner clad. surface, c) Circ. stress at the outer clad. surface, d) Circ. stress at the inner clad. surface 5.3 Presentation of the stress-strain field At a given outer pressure (coolant pressure) acting on the fuel rod, the contact pressure increases the stress in the clad- ding. Reaching the yield stress value can be the primary cause of cladding failure. The critical value is also dependent, e.g., on the fission rate (fluence). As comparable stresses, we chose the circumferential and equivalent (effective) stresses in the most exposed places, where the highest stress due to the “hour glassing” effect of the pellet is expected. In the FEMAXI-V mesh, there is a node with coordination axial 6 (IZ � 6 in Fig. 2) and radial 1 (IR � 1 in Fig. 2), i.e., the inner surface of the cladding and of course also node axial 6 (IZ � 6 in Fig. 2) and radial 2 (IR � 2 in Fig. 2), i.e., the outer surface of the cladding. These values vs. time can be seen in Figs. 6a–6d. The stresses correspond well with the gap size, the gap is closed after about 12–18 hours, from where the stresses show a high increase and reach the peak value (for equivalent stress at the outer surface of the cladding about 65 MPa) and then they relax during the steady power. A gap opens with the power decrease, and the stresses fall almost to their initial values. When there is low power (fourth day), the stresses remain also low. The second power increase induces the same behaviour, only with higher peak values. The COSMOS /M results are shown in Fig. 7. The top line corresponds with the line in Fig. 6b, second top Fig. 6a, third Fig. 6c and fourth Fig. 6d. The difference can be seen on the first day – the initial val- ues are zero, in contrast with FEMAXI-V. From the contact point COSMOS /M induces the same behaviour as FEMAXI-V with higher peak values (85 MPa (effective stress) for the outer surface). On the second and sixth days relaxation takes place due to steady power and established contact conditions. During day four, the stresses remain low. Slight numerical in- stability can also be seen at the beginning of the second day (24th–25th hour). All calculations were made with contact coefficient of fric- tion 0,1 and the following questions remain open: � What is the influence of bonding or sliding with friction on the transient behaviour? � What are the real (realistic) local friction coefficients for medium and highly burned-up fuel? 6 Conclusion The paper summarizes the preliminary results from a lo- cal axisymmetric model developed in the COSMOS/M sys- tem. For validation we used the results from the FEMAXI-V integral code. The absolute values of the calculated stresses are higher for COSMOS /M, and the stress progress in time differs on the first day of the power history. We were not fully able to define a totally compatible initial modelling state in both codes. We regard this as a weakness, and we will work on a more detailed analysis. © Czech Technical University Publishing House http://ctn.cvut.cz/ap/ 9 Czech Technical University in Prague Acta Polytechnica Vol. 45 No. 5/2005 Fig. 7: Stresses (MPa) from COSMOS/M vs. time (seconds) References [1] Timošenko, Š.: Strength of Materials I., II. Praha: Tech- nicko-vědecké vydavatelství, 1951 (in Czech). [2] Hlaváček, I., Haslinger, J., Nečas, J., Lovíšek, J.: Solution of Variational Inequalities in Mechanics. Bratislava: Alfa, 1982 (in Slovak). [3] Hinton, E., Owen, D. R. J.: Finite Element Programming. London: AP, 1977. [4] COSMOS/M 2.8 for Windows – online manuals, SRAC 2003. [5] Suzuki, M.: FEMAXI-V manual, Japan, Jaeri 2000. Ing. Martin Dostál e-mail: m.dostal@post.cz Department of Nuclear Reactors Czech Technical University in Prague Faculty of Nuclear Sciences and Physical Engineering V Holešovičkách 2 180 00 Praha 8, Czech Republic Nuclear Research Institute Řež plc. Reactor Technology Department 250 68 Husinec-Řež 130, Czech Republic RNDr. Jiří Zymák, CSc. e-mail: zym@ujv.cz Ing. Mojmír Valach, CSc. e-mail: val@ujv.cz Nuclear Research Institute Řež plc. Reactor Technology Department 250 68 Husinec-Řež 130, Czech Republic 10 © Czech Technical University Publishing House http://ctn.cvut.cz/ap/ Acta Polytechnica Vol. 45 No. 5/2005 Czech Technical University in Prague