Annals of Nuclear Energy 141 (2020) 107353
Contents lists available at ScienceDirect
Annals of Nuclear Energy journal homepage: www.elsevier.com/locate/anucene
Multi-scale coupling of CFD code and sub-channel code based on a generic coupling architecture Xilin Zhang a,⇑, Kanglong Zhang b, V.H. Sanchez-Espinoza b, Hongli Chen a a b
School of Physical Sciences, University of Science and Technology of China, Hefei, Anhui 230027, China Institute for Neutron Physics and Reactor, Karlsruhe Institute of Technology, Eggenstein-Leopoldshafen, 76344, Germany
a r t i c l e
i n f o
Article history: Received 26 September 2019 Received in revised form 30 December 2019 Accepted 20 January 2020
Keywords: Coupling Sub-channel code CFD code ICoCo MEDCoupling
a b s t r a c t This paper describes a multi-scale thermal–hydraulic coupling system, which combines the capabilities of the open-source TrioCFD code and of the sub-channel code SubChanFlow (SCF) aiming to improve the description and prediction of the multi-scale thermal–hydraulic physical phenomenon inside the reactor vessel. This is a parallel coupling system based on the ICoCo (Interface for Code Coupling) concept, where each code is first wrapped by the ICoCo interface and then it is compiled to a shared library. A parallel C++ script was developed as the supervisor, which utilizes the two shared libraries and supervises the calculation of the coupled system. The MEDCoupling libraries developed by CEA provide a generic way to handle the mesh interpolation and field mapping between different domains of TrioCFD and SCF. Moreover, the explicit temporal coupling method and the domain-decomposition approach are adopted and presented in this paper. An academic problem is developed to test the multi-scale coupling code. The results indicate that the inter-code data exchange works well and that the coupled code TrioCFD/SCF can deal with various time-dependent boundary conditions. Ó 2020 Elsevier Ltd. All rights reserved.
1. Introduction Current thermal hydraulic simulations of nuclear power plants are based on different numerical simulation methods that historically were developed parallel to the progress of the reactor technology and regulatory requirements. There are three main kinds of thermal–hydraulic codes commonly used by the nuclear community: 1) system codes; 2) sub-channel codes; 3) CFD codes. They correspond to three scales in sequence: The macro scale where the simulation could cover the full nuclear power system, and give reasonable results with accessible consumption of the computation resources. Some of the well-known system codes are RELAP5 family, TRACE family, CATHARE family, ATHLET family, etc. The component scale where the results’ resolution is refined in the core and some crucial safety parameters e.g. the CHF could be sufficiently inspected. Some of the well-known sub-channel codes are COBRA-IV (Wheeler et al., 1976), COBRA-EN (Basile et al., 1999), FLICA (Toumi et al., 2000), SUBCHANFLOW (Imke and Sanchez, 2012), KMC-SUB (Li et al., 2017), etc.
⇑ Corresponding author. E-mail address:
[email protected] (X. Zhang). https://doi.org/10.1016/j.anucene.2020.107353 0306-4549/Ó 2020 Elsevier Ltd. All rights reserved.
The micro scale where the CFD code could give very fine descriptions of the local thermal–hydraulic phenomenon e.g. the flow distribution in the lower plenum. One branch is the commercial codes, like ANSYS CFX (ANSYS CFX, 2009), ANSYS FLUENT (ANSYS FLUENT, 2009), Star-CCM+ (CD-adapco STAR CCM+, 2017), etc. The other branch is open-source codes, representatives are OpenFOAM (OpenFOAM Foundation Ltd, 2015), TrioCFD (Angeli et al., 2015), PyFR (Witherden et al., 2014), etc. Traditional numerical tools e.g. system thermal hydraulic codes and sub-channel codes use simplified models with a large number of empirical correlations validated with test data. These simplified models could produce reasonable well the main physical phenomena of nuclear reactors. A large number of verification and validation work has been performed for each reactor design e.g. OECD validation matrix summarizes the validation work performed so far in the frame of international activities. In last years, coarse mesh 3D system thermal hydraulic codes with six equations models in the two-fluid approach were developed und validation work was done specially for 3D phenomena inside the RPV e.g. using coolant mixing, and boron dilution data. Meanwhile, new versions of subchannel codes with four, six, and nine conservations equations considering two and three fluids were developed to simulate not only single fuel assemblies but also full core at both subchannel or fuel assembly level.
2
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353
However, the system codes and sub-channel codes cannot precisely handle the phenomena where significant 3D effects dominate e.g. in the downcomer, lower and upper plena. This drawback mainly comes from the simplification during the codes modeling. However, sometimes we indeed have great interests in the 3D dominating flow conditions. For instance, the asymmetrical temperature distribution during a Main Steam Line Break and the boron concentration during a boron dilution scenario, in Pressurized Water Reactors (PWR) (Bertolotto, 2011). For pool-type Generation IV reactors, these phenomena include thermal mixing and thermal stratification due to the large amount of primary coolant and the existence of hot/cold plenum (Roelofs et al., 2019; Tenchine, 2010). This requirement and the huge computer power available nowadays at low cost has also fostered the massive use of CFD for design optimization and safety investigations of nuclear reactors. But the current status of the two-phase flow models and the impossibility of the development of highly resolved full cores of a nuclear power plant at micro-level (i.e. each subchannel and fuel rod) requires to look for alternative solutions to improve the simulation of safety-relevant phenomena within the nuclear power plant circuits and core. One promising alternative is the development of multiscale computational methods by combining different types of thermal hydraulic codes using innovative coupling schemes to overcome the limitations of these codes when applied as single applications. During the last two decades, multi-scale thermal–hydraulic coupling codes have been developed worldwide e.g. the coupling of system code with sub-channel code (Liu and Cheng, 2015; Zhang and Sanchez-espinoza, 2019; Zhang, 2019), the coupling of system code with CFD code (Bertolotto, 2011; Pialla et al., 2015). However, the traditional coupling is based on the master–slave scheme which is not so generic, and the data transfer method adopted is too simple to meet the requirements of high-fidelity multi-scale multi-physics simulations.This paper discusses the coupling of a sub-channel code – Subchanflow (SCF) and an open-source CFD code – TrioCFD based on a generic architecture. This is a coupling system based on the innovative ICoCo interface and coordinated by a newly-developed C++ supervisor. The codes involved (SCF, TrioCFD, and the C++ supervisor code) are shortly described in Chapter 2. The coupling scheme including the coupling architecture, the coupling interface, and the coupling algorithm (the data exchange methods, the code execution mode, the temporal coupling scheme, and the computation domain division) are presented in Chapter 3. The testing of the coupled code is performed in Chapter 4. Finally, the summary and outlook complete this paper. 2. Description of the codes involved The research goal of this paper is to couple Sub-channel code and CFD code. Other codes (system codes, mechanical analysis codes, 3D neutron diffusion codes and so on) will be considered in the future work. As shown in Fig. 1, to achieve the supervised scheme, not only the codes to be coupled (Sub-channel code and CFD code) but also the supervisor code (C++ script/Python script/Salome platform (Babur et al., 2015; CEA, 2019)) are needed.
Supervisor (C++, Python, Salome) Method call
ICoCo (C++)
ICoCo (C++)
Code 1 Resolution of Problem 1 (any language)
Code 2 Resolution of Problem 2 (any language)
Fig. 1. The architecture of the ICoCo-based SCF/TrioCFD coupled code (Deville and Perdu, 2012).
ment methods have been improved and the coolant properties in the code have been updated. SCF is part of the NURESIM simulation platform where it is coupled with DYN3D, CRONOS, and COBAYA3. Recently, it was coupled with TRACE via a socket-based interface named ECI (the Exterior Communication Interface) (Zhang and Sanchez-espinoza, 2019; Zhang, 2019) and further equipped with an ICoCo API interface (Deville and Perdu, 2012) at KIT (Zhang and Sanchez, 2019). For now, the ICoCo interface allows SCF to write in and extract diverse fields on various kinds of meshes. These fields include the thermal power, the coolant temperature, the boron concentration, the coolant velocity, the system pressure and so on. The meshes are built according to ICoCo standards, which will be illustrated in section 2.3. For the time being, the fuel mesh and coolant mesh could be built. All of the information defining the meshes come from the SCF input card. No additional data is required. Moreover, the coolant mesh types further include Cartesian/Hybrid mesh for square fuel assembly, while unstructured mesh with hexagon cells or triangle cells for hexagon fuel assembly. The fuel mesh always consists of discrete square or cuboid cells. 2.2. The CFD code - TrioCFD TrioCFD is an open-source CFD code based on the TrioCFD platform (Trio_U Software for Thermo-hydraulics) being developed by the Thermo-hydraulics Service and Fluid Mechanics (STMF) of the Department of Nuclear Energy at the CEA. Two main models in TrioCFD are the Reynolds-Averaged Navier-Stokes (RANS) and the Large Eddy Simulation (LES). The Direct Numerical Simulation (DNS) is also available but mostly for academic purposes. The V&V is now a major process providing an evaluation of the reliability level of the computed solutions, as well as the correct implementation of the desired models ICoCo is already a build-in coupling-oriented interface in TrioCFD as an inherent module (Pialla et al., 2015). Through it, TrioCFD (the former Trio_U) has been coupled with the system code CATHARE (Bavière et al., 2014). For the sake of coupling issues, TrioCFD is developed based on MED and MEDCoupling libraries. This is also to meet the requirement of ICoCo standards. To ensure correct storage and exchange of MED format data (especially under the parallel environment), the MED libraries in TrioCFD are used as the basic libraries for the whole supervised scheme.
2.1. The sub-channel code - SCF 2.3. The C++ supervisor SCF is a thermal–hydraulic sub-channel code developed at KIT for the simulation of fuel rod bundles and cores of light water and innovative reactor systems. It solves the mixture equations for mass, momentum, and energy and an additional one for lateral flow. It is completely written in Fortran 2003. The data manage-
As shown in Fig. 1, there are three forms of supervisor code: C++ script, Python script and Salome platform. Among the three options, the SALOME platform is an opensource platform for pre- and post-processing of simulations.
3
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353
Besides, it also allows integrating different kinds of physical solvers into one single simulation platform. These physical solvers may involve material science, solid mechanics, thermal hydraulics, nuclear physics and so on. Especially, users could organize and monitor the simulation workflow in the SALOME-YACS module (CEA, 2019) conveniently through a GUI. Moreover, its textual user interface (TUI) can also parameterize and automate the workflow. To couple codes via Salome platform, we need to integrate the codes (any language) separately into Salome via YACSGEN (CEA, 2019) module. During the integration, ICoCo wrapping (C++), SWIG (Simplified Wrapper and Interface Generator) (Beazley, 2003) wrapping (Python) and CORBA (Vinoski, 1997) wrapping (IDL, Interface Description Language) should be performed successively. Moreover, SALOME does not support MPI (Message Passing Interface) functions for the time being. The MPI capability is indispensable in TrioCFD, which must be considered during the supervisor selection. Compared with Salome platform as the supervisor, C++ script as the supervisor is a more direct way and can avoid the SWIG wrapping and CORBA wrapping. The supervisor written by Python will be considered in the future or if necessary. The flowchart of the supervisor code is shown in Fig. 4. Examples of C++ supervisor codes can be found in the annex of the ICoCo report (Deville and Perdu, 2012) and the installation directory of TrioCFD.
3. The selected multi-scale coupling scheme The coupling scheme can be divided into 4 parts: the coupling architecture, the coupling interface, the coupling algorithm and the workflow of coupling computation.
3.1. The coupling architecture There are two main options for the coupling architecture: the master-slave scheme; the supervised scheme. The master-slave scheme means that one code will control itself and the other codes. Examples include: OpenFOAM (master)/ATHLET (slave) via zeromq/msgpack (Pialla et al., 2015), ANSYS CFX (master)/ATHLET (slave) via CFX User Fortran (Papukchiev et al., 2015), ANSYS Fluent (master)/PKM (slave) via Fluent UDF (Chen et al., 2015) and so on. However, the code interfaces usually just fit one single pair of codes. For this reason, the master–slave scheme is not as generic as the supervised scheme. In the supervised scheme, as shown in Fig. 1, a supervisor program controls the coordination of the codes and the inter-codes data exchange. The supervisor could be C++ script, Python script or SALOME platform (for GUI-assisted coupling). As explained in Section 2.3, the C++ script is selected. To enable the supervision, each code coupled should provide efficient interfaces. The code wrapping method (especially about the data exchange) will be discussed in the section below. 3.2. The code interface The Interface for Code Coupling (ICoCo) concept, proposed by CEA, is selected as the code wrapping concept. The ICoCo concept requires the codes to implement a list of API (Application Programming Interface) functions to perform the following functions: initialize the code, compute the time step, initialize with the time step, solve/iterate with the time step, extract the output field, set the input field and terminate the code. Fig. 2 illustrates the calling sequence of ICoCo API functions in the supervisor program.
Problem
0..*
setDataFile setMPIComm 0..*
initialize
getInputFieldsNames getOutputFieldsNames getInputFieldTemplate getOutputField isStationary
0..* computeTimeStep 0..*
save restore forget
initTimeStep
presentTime can be called at any time in the green box
Resolution and field exchange on the computation timestep
validateTimeStep
1..* setInputField 0..* 1..*
solveTimeStep iterateTimeStep
abortTimeStep
terminate
~Problem Fig. 2. Execution flow chart using ICoCo interface (Deville and Perdu, 2012).
solveTimeStep iterateTimeStep
4
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353
The data exchange between the codes is necessary for the coupling. However, different codes have different forms of input/outputs. Especially, for CFD codes, the input/output data can be large, complex and heterogeneous. A generic and efficient data format to store the data into memories/files is needed. The MED (modèle d’échange de données, Data Exchange Model) libraries (CEA, 2019), developed by CEA and EDF, is a good solution. MED libraries can allow the codes to store mesh/field data into memory/files, and perform interpolations/projections of field data defined on different meshes. Fig. 3 shows how the MED libraries are developed. As shown in Fig. 3, the MED libraries (gray blocks) mainly consist of three parts: Data format in memory: MEDCoupling Data format in files: MEDLoader (Read/write MED files), where the MED file format is developed based on the famous hdf5 file format. Data interpolation in memory: Remapper (sequential), ParaMEDMEM (parallel). The data dealt with by MED libraries is mesh/field data, which is usual for CFD codes. When implementing the ICoCo API functions about data exchange into other codes (e.g. sub-channel codes), we must ensure that each code to be coupled can input/output the same data type defined in MEDCoupling. 3.3. The coupling algorithm In the supervised scheme, the coupling algorithm (code execution mode, temporal coupling scheme, and computation domain division) is written only in the supervisor code rather than in the codes coupled. 3.3.1. Code execution mode For the supervised scheme, there are two code execution modes: 1) sequential mode; 2) parallel mode. In the sequential mode, the codes coupled are executed one by one. In the parallel mode, the codes are executed at the same time. Considering that parallel computation is necessary for CFD codes, the parallel mode is selected. The parallel mode is achieved just based on the Message Passing Interface (MPI) (Kranzlmüller et al., 2004) technology. The data communication between different processors is achieved by MPI_Allreduce function of MPI libraries (for simple data types, i.e. int/double/bool) and TrioDEC class (for field data, inheriting from InterpKernelDEC class defined in the MEDCoupling namespace). 3.3.2. Temporal coupling scheme There are three alternative temporal coupling schemes:
1) Explicit. The codes will exchange data just after each timestep iteration. This is supported by the ‘‘solveTimeStep” function of the ICoCo concept. 2) Semi-implicit. The codes will exchange data after several sub-iterations during each time-step iteration. This is supported by the ‘‘iterateTimeStep” function of the ICoCo concept. 3) Implicit. The governing equations of the codes are discretized and solved together in the same way. This does not fit the research goal of this paper to couple two separate codes. In this paper, the simplest scheme ‘‘Explicit” is firstly selected. 3.3.3. Computation domain division As discussed in the Introduction, different codes are good at simulating different parts of the computation domain. There are two methods to assign the computation domain for each code: domain decomposition method; domain overlapping method. For the domain decomposition method, different codes will simulate different parts of the computation domain. The data exchange will be conducted at the boundaries of sub-domains. For the domain overlapping method, the results of one code for one part will be corrected by another code. The data exchange will be conducted through virtual distributed source terms (somehow similar to the porous medium method). In this paper, the simpler domain decomposition method is adopted. The domain overlapping method will be considered in future work. 3.4. The workflow of coupling computation The flowchart of the supervisor code is presented in Fig. 4. In the beginning, each sub-code should be initialized with its input files. After the calculation timesteps for the sub-codes are obtained, the minimum value is set as the real timestep for the calculations afterward. After finishing solving each timestep, the sub-codes will exchange data via Data Exchange Channels (DEC) in the supervisor code. If each timestep is solved correctly, the physical time will proceed with the calculation timestep. Otherwise, the physical time will not change and the calculation time steps of the subcodes will be computed again. The supervisor code will terminate the running of all the sub-codes when any sub-code stops the calculation. Notably, before the transient calculation, the steady-state calculation must be performed once by SUBCHANFLOW (Referred to as SCF in Fig. 4.). 4. Testing of the multi-scale coupling code In this paper, an academic problem is developed to test the capability of the multi-scale T/H coupling code. The test is used to verify the data transfer between the codes and the stability of the coupled code under different transient conditions. It is worth to mention that these tests are of academic nature. 4.1. Specification of the academic problem
Fig. 3. The MED libraries architecture (CEA, 2019).
There are 3 computation domains in this academic problem, see Fig. 5. Two TrioCFD solvers are launched to simulate the pipe region (red part in Fig. 5) and the plenum region (purple part in Fig. 5), respectively. The core region (blue part in Fig. 5) is simulated by the SCF solver. In Table 1, the geometry models and meshes of these computational domains are described. As shown in Fig. 5, the data transfer between solvers are conducted at the core inlet and the core outlet. At these 2 coupling interfaces, the boundary meshes for different solvers are shown in Fig. 6.
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353
5
Create Data Exchange Channels(DECs) between codes Set input files and MPI settings for each code Initialize each code Compute timestep and initialize with timestep Yes Terminate all the codes
Check if any code can stop now No Get template of input field; Get output field data;
Send/Receive field data via DEC
Set the input of each code
Check if the first coupling is Yes
No
Only SCF perform steady-state calculation
Solve current timestep
No
Check if the solving is successful Yes
Abort current timestep
Validate current timestep
Physical time doesn't change
Physical time proceeds with t=t+dt
Fig. 4. The flowchart and the coupling algorithm of supervisor code.
In the multi-scale coupling code, each solver uses its unique physical models to simulate specific phenomenon at a certain scale. Table 2 explains all the important aspects of the physical models adopted in the multi-scale coupling code corresponding to this academic problem. In TrioCFD solver, considering that the Prandtl number of lead is low, the turbulent Prandtl number Prt is calculated from the following analytical function of Kays (Bieder et al., 2018), rather than taken as constant (the default value is 0.9 for most fluids like air or water).
Prt ¼ 0:85 þ 0:7
k
tt qc p
ð4:1Þ
In SCF solver, for the rod bundle with lead as coolant, the following correlations (Equations 4.2–4.6) are applied to determine the pressure drop and heat transfer. The friction factor for flow in rod bundles f b is calculated by the equation below, recommended by the CEFR designers (Xu, 2011).
h i 4=3 f b ¼ f 1 þ aðP=dÞ 0:58 þ 0:42eb
ð4:2Þ
where f is the friction factor for flow in round tubes, maybe obtained from the Blasius correlation; P is the rod pitch of rod bundles; d is the rod diameter; a; b are the constants dependent on the bundle configuration. For the square array,
6
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353 Table 2 The physical models adopted in the multi-scale coupling code.
Coolant Solver Flow and heat transfer models
Thermophysical properties of coolant
Fig. 5. The computation domains.
a ¼ 0:1166
b¼
h i3=4 3 2 740:3 106 ðP=dÞ 1:273ðP=dÞ 1 ½1:122ðP=dÞ 1
b¼
h
2
28:45 10 ðP=dÞ 1:102ðP=dÞ 1
Lead TrioCFD CFD model Std. k -e turbulence model Boussinesq assumption
Lead SCF Sub-channel model
Lead TrioCFD CFD model Std. k-e turbulence model Boussinesq assumption
Pressure drop correlation (Eqs. (4.2)– (4.4)) Heat transfer correlation (Equations (4.5)–(4.6))
Turbulent Prandtl Turbulent Prandtl number is number is calculated from calculated from Kays’ function Kays’ function (Eq. (Equation (4.1)) (4.1)) The dynamic viscosity, the thermal expansion coefficient, and the thermal conductivity are functions of temperature (Fazio, 2015). Density as constant Density as Density is (440 °C) constant (440 °C) function of temperature (Fazio, 2015) The specific isobaric heat capacity is taken as constant (440 °C)
Nu ¼ 0:047 1 e3:8ðP=d1Þ Pe0:77 þ 250
i3=4 ð4:4Þ
4:5
½1:044ðP=dÞ 1
The calculation method of Nusselt number is also dependent on the bundle configuration (Fazio, 2015). For the square array (1:2 P=d 1:5; 10 Pe 2500),
Nu ¼ 7:55ðP=dÞ 14ðP=dÞ
Plenum
For the triangular array (1:1 P=d 1:95; 30 Pe 5000),
a ¼ 0:1066 3
Core
ð4:3Þ
4:5
For the triangular array,
6
Pipe
5
þ 0:007Pe0:64þ0:246ðP=dÞ
ð4:5Þ
ð4:6Þ
where Pe is the Peclet number, Pe ¼ RePr. To perform a whole numerical simulation, it is necessary to make clear the geometry model, the discretization method (e.g. mesh), the physical models (e.g. governing equations, simplified methods, constitutive equations), and the initial & boundary & sources conditions. Table 3 lists the design parameters under normal conditions.
Table 1 The geometry models and meshes of computational domains. Pipe
Core
Plenum
Geometry model Dimension
short square tube 0.48 0.48 m 0.48 m
3 3 square fuel assemblies from SNCLFR-100 0.48 m 0.48 m 3.40 m
Mesh
Mesh type 3D cells 2D cells Refinement
Cartesian mesh 20 20 20 hexahedrons 20 20 quadrangles (Inlet & Outlet) refined near the wall boundary
Cylinder with 1 square inlet and 3 outlet nozzles Cylinder height 1.5 m Cylinder diameter 1.44 m Nozzle centerline height 1.1 m Nozzle diameter 0.1 m Inlet dimension 0.48 m 0.48 m Inlet height 0.4 m Cartesian mesh Unstructured mesh 3 3 20 hexahedrons 134,325 tetrahedrons 3 3 quadrangles (Inlet & Outlet) 200 triangles (Inlet) —— refined near the wall boundary
(a) The pipe outlet
(b) The core inlet/outlet Fig. 6. The boundary meshes at 2 coupling interfaces.
(c) The plenum inlet
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353
7
h
Table 3 The design parameters under normal condition. Parameter
Value
Coolant temperature at pipe inlet Coolant velocity at pipe inlet Power generated in core region Average coolant temperature at core outlet Power factors of fuel assemblies
400 °C 1.5 m/s 4.2485 107 W 480 °C 0.75, 0.80, 1.10, 1.20, 1.30, 1.20, 1.10, 0.80, 0.75
4.2. Data transfer capability In this problem, the boundary conditions of 3 computation domains are listed in Table 4. However, only the boundary conditions in blue or red are involved in the coupling interfaces. The detailed data transfer at the 2 coupling interfaces is shown in Fig. 7. Via ICoCo interface, all the output fields of solvers are all fields represented in MEDCouplingFieldDouble format. Since that the solvers are executed in parallel, the data transfer of MEDCoupling fields is based on the ParaMEDMEM module. For thermal–hydraulic codes, there are 3 kinds of physical fields: velocity, temperature, and pressure. The data transfer of these fields should follow the corresponding physical laws: mass conservation law (velocity), energy conservation law (temperature), and maximum principle (pressure and backflow temperature). Before data transfer, the following field conversions need to be performed: 1) velocity field v ? mass velocity field qv , to ensure the conR versation of X ðPv Þds; 2) temperature field T ? enthalpy flux field qv h (equivalent to qv T in this academic problem), to ensure the conversation of R X ðP v T Þds. The conservation laws and the maximum principle will be carefully checked in this case. Different spatially varying forms of fields (velocity, temperature, and pressure) are considered. In this case, the field plane is perpendicular to the z-axis. The fields to be transferred as shown in Fig. 7 are set manually with the general cosine distribution function below:
x x p y y o p zi zo p o i U ¼ Acos 2p þ cos 2p i þ cos 2p þ k 4 k 4 k 4 i þ A þ d Direction ð4:7Þ where U is the pressure value or the velocity component in the zdirection; A is the amplitude; xi , yi and zi are the coordinates of the ith cell; xo , yo and zo are the coordinates of the field center; k is the wavelength of the distribution; d is the offset value to ensure the overall value not too small, the default value is 0.0; Direction has 2 alternative values: 1 or 1. The data transfer is checked through the comparison of the source field and the transferred target field. The different distribution forms with different wavelengths (1.0 m, 0.1, 0.01 m) and the compound distribution of the former three distribution forms are considered in this test. The results of the compound distribution (1 m + 0.1 m + 0.01 m) are shown in Table 5. As shown in Table 5, because the mesh of sub-channel code is usually coarser than that of CFD code, the data transferred from CFD code to sub-channel code will lose certain details but the fidelity is enough for the calculation of sub-channel code. Thanks to the finer mesh and the MED libraries, the data transfer from sub-channel code to CFD code can keep the details of the original distribution very well. Moreover, according to the checking of the field profile and the important parameters (e.g. mass, energy, and force), in all the scenarios, there is no conflicts with the existing physical laws. However, according to the results about data transfer from subchannel code to CFD code, to keep all the details of sub-channel code, it’s better to ensure the meshes of CFD code can well represent the boundary mesh cells of sub-channel codes. 4.3. Transient calculation stability The checking of data transfer capability just can be used to check whether the MED libraries are used correctly. To check whether the coupled code can work normally, some tests need to be performed. In this section, different temporally varying forms of spatially uniform fields (velocity and temperature) at pipe inlet are considered. The velocity field and the temperature field at pipe inlet are set manually with the general cosine evolution function below:
t to U ¼ A sin 2p þ A þ d Direction T
Table 4 Boundary conditions of computation domains. Pipe
Core
Plenum
Inlet
velocity, temperature
velocity, temperature
Outlet
backflow temperature, pressure
velocity, temperature pressure
backflow temperature, pressure
ð4:8Þ
Where: U is the value of the physical parameter (like velocity component in z-direction); A is the amplitude; t is the present physical time;
Fig. 7. The data transfer flow in the multi-scale coupling code.
8
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353
Table 5 Data transfer between solvers. Field description
Source field
Target field
Physical laws
Pressure field Core ? Pipe SCF ? TrioCFD (1st coupling interface)
Maximum principle R X Pds
Temperature field Pipe ? Core TrioCFD ? SCF (1st coupling interface)
Energy conservation R X P v Tds
Velocity field Pipe ? Core TrioCFD ? SCF (1st coupling interface)
Mass conservation R X P v ds
Pressure field Plenum ? Core TrioCFD ? SCF (2nd coupling interface)
Maximum principle R X Pds
Temperature field Core ? Plenum SCF ? TrioCFD (2nd coupling interface)
Energy conservation R X P v Tds
Velocity field Core ? Plenum SCF ? TrioCFD (2nd coupling interface)
Mass conservation R X P v ds
9
2.0
2.0
1.5
1.5 Flowrate factor
Flowrate factor
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353
1.0 Pipe_in Pipe_out Core_in Core_out Plenum_in Plenum_out
0.5
0.0
3
1.0 Pipe_in Pipe_out Core_in Core_out Plenum_in Plenum_out
0.5
5
7
Time (s)
0.0 6.00
9
6.04
6.08
6.12
6.16
6.20
6.12 Time (s)
6.16
6.20
Time (s)
3-9s
6.0-6.2s Fig. 8. The evolution of the mass flowrate.
490
490 Pipe_in Pipe_out Core_in Core_out Plenum_in Plenum_out
430
460 o
Temperature( C)
o
Temperature( C)
460
400
400
370
370
340
430
3
5
Time (s)
7
9
340 6.00
Pipe_in Pipe_out Core_in Core_out Plenum_in Plenum_out
6.04
3-9s
6.08
6.0-6.2s Fig. 9. The evolution of the average temperature.
to is the time when the transient condition starts; T is the period of the evolution of physical parameter; d is the offset value to ensure the overall value not too small, the default value is 0.0; Direction has 2 alternative values: 1 or 1. In this academic problem, three evolution periods (1.0 s, 0.1 s, 0.01 s) and the compound evolution (1.0 s + 0.1 s + 0.01 s) are considered. The results of the compound evolution (1.0 s + 0.1 s + 0. 01 s) are given in Figs. 8–12. Fig. 8 shows the evolutions of the mass flowrates at pipe inlet, pipe outlet, core inlet, core outlet, plenum inlet, and plenum outlet. The mass flowrates at pipe outlet and core inlet agree with each other, and so do the mass flowrates at core outlet and plenum inlet. Because the speed of pressure wave is so high, all the mass flowrates follow the evolution of the mass flowrate at pipe inlet. However, the evolution of temperature depends on heat conduction and heat convection. As shown in Fig. 9, although the change of temperature field at the pipe is highly frequent, the
average temperature at pipe outlet changes more gradually. During the whole transient, the heating power of core region doesn’t change. In this problem, due to the thermal inertias of computation domains and the combination of evolutions of physical fields (velocity and temperature) at pipe inlet, the average temperature at core outlet changes more frequently than that at the plenum outlet. The evolution of velocity distribution at plenum inlet is shown in Fig. 10. Due to the strong turbulent mixing effects, the outlet velocity distribution of pipe is usually uniform in this case, and thus the velocity distribution at plenum inlet is nearly uniform during all the transient processes. However, the absolute value of velocity changes periodically just like the inlet velocity of pipe. For pool-type reactors (e.g. Lead-cooled fast reactor), the thermal stratification, the recirculation flow and the mixing of coolant jets with different temperatures are important thermal hydraulic phenomena with significant 3D effects. The non-uniform power distribution of core region is set to enhance the jet mixing. As
10
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353
Fig. 10. The evolution of velocity distribution at plenum inlet.
Fig. 11. The evolution of temperature distribution at plenum inlet.
shown in Fig. 11, corresponding to the power distribution of core, the inlet temperature distribution of plenum is non-uniform. The strong recirulation flow and jet mixing are found in Fig. 12. 5. Conclusions The open-source TrioCFD code and the KIT SubChanFlow (SCF) code were successfully coupled through a standardized interface
named ICoCo (Interface for Code Coupling) coordinated by a parallel C++ supervisor. The domain-decomposition approach was applied, where SCF and TrioCFD simulate different thermal–hydraulic domains sharing one boundary interface. The MEDCoupling library took charge of the mesh interpolation and field mapping between the different meshes of SCF and TrioCFD. The explicit temporal scheme was implemented for data synchronization during time advancement.
X. Zhang et al. / Annals of Nuclear Energy 141 (2020) 107353
Fig. 12. The streamline in the plenum (7.0 s).
An academic problem was developed to test the coupled code. The data transfer capability at the coupling interfaces and the transient calculation stability are reasonable and good. For the former one, various spatially varying forms of vector/scalar field data at the coupling interface are considered. For the latter one, various temporally varying forms of spatially uniform velocity field at pipe inlet are considered. In the near future, the validation of the coupled code will be performed according to (Oberkampf and Roy, 2010; Versluis, 2013) by using appropriate experimental data. Its application into the Light Water Reactors and GEN IV reactors will also be conducted in the future. CRediT authorship contribution statement Xilin Zhang: Conceptualization, Methodology, Software, Formal analysis, Investigation, Visualization, Data curation, Writing - original draft. Kanglong Zhang: Conceptualization, Methodology, Software, Writing - review & editing. V.H. Sanchez-Espinoza: Conceptualization, Methodology, Writing - review & editing, Supervision. Hongli Chen: Writing - review & editing, Supervision. Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Acknowledgements The authors appreciate the technical support and fruitful discussions of CEA (Commissariat à l’énergie Atomique, France) and of KIT (Karlsruhe Institute of Technology, Germany). This work is supported by the scholarship from China Scholarship Council (CSC, China). Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi.org/10.1016/j.anucene.2020.107353. References Angeli, Pierre Emmanuel Bieder, U., Fauchet, G., 2015. Overview of the TrioCFD Code : Main Features , V & V Procedures and Typical Applications To Nuclear Engineering, in: In Proceedings of 16th International Topical Meeting on NUclear REactor Thermal Hydraulics (NURETH-16). Chicago(USA), pp. 252–265.
11
ANSYS CFX, 2009. ANSYS CFX 12.0 User’s Guide. New York (USA). ANSYS FLUENT, 2009. ANSYS FLUENT 12.0 User’s guide. New York (USA). Babur, Ö., Smilauer, V., Verhoeff, T., Van Den Brand, M., 2015. A survey of open source multiphysics frameworks in engineering. Procedia Comput. Sci. 51, 1088–1097. Basile, D., Beghi, M., Chierici, R., Salina, E., Brega, E., 1999. COBRA-EN: an upgraded version of the COBRA-3C/MIT code for thermal hydraulic transient analysis of light water reactor fuel assemblies and cores. ENEL-CRTN, Milano. Bavière, R., Tauveron, N., Perdu, F., Garré, E., Li, S., 2014. A first system / CFD coupled simulation of a complete nuclear reactor transient using CATHARE2 and TRIO U . Preliminary validation on the Phénix Reactor Natural Circulation Test 277, 124–137. Beazley, D.M., 2003. Automated scientific software scripting with SWIG. Futur. Gener. Comput. Syst. 19, 599–609. Bertolotto, D., 2011. Coupling a System Code with Computational Fluid Dynamics for the Simulation of Complex Coolant Reactivity Effects. École Polytechnique Fédérale de Lausanne. Bieder, U., Ziskind, G., Rashkovan, A., 2018. CFD analysis and experimental validation of steady state mixed convection sodium flow. Nucl. Eng. Des. 326, 333–343. CD-adapco STAR CCM+, 2017. STAR CCM+ User’s Guide (Version 12.04). New York (USA). CEA, 2019. SALOME version 9.4.0 Release. Notes. CEA. Chen, Z., Chen, X.N., Rineiski, A., Zhao, P., Chen, H., 2015. Coupling a CFD code with neutron kinetics and pin thermal models for nuclear reactor safety analyses. Ann. Nucl. Energy 83, 41–49. Deville, E., Perdu, F., 2012. Documentation of the Interface for Code Coupling : ICoCo. CEA. Fazio, C. et al., 2015. Handbook on lead-bismuth eutectic alloy and lead properties, materials compatibility, thermal-hydraulics and technologies-2015 edition. Organisation for Economic Co-Operation and Development. Imke, U., Sanchez, V.H., 2012. Validation of the Subchannel Code SUBCHANFLOW Using the NUPEC PWR Tests (PSBT). Sci. Technol. Nucl. Install. 2012, 1–12. Kranzlmüller, D., Kacsuk, P., Dongarra, J., Gabriel, E., Fagg, G., Bosilca, G., Angskun, T., Squyres, J., Sahay, V., Kambadur, P., Barrett, B., Lumsdaine, A., Castain, R., Daniel, D., Graham, R., Woodall, T., 2004. Open MPI: Goals, Concept, and Design of a Next Generation MPI Implementation. Lecture Notes in Computer Science. Springer, Berlin, Heidelberg. Li, S., Cao, L., Khan, M.S., Chen, H., 2017. Development of a sub-channel thermal hydraulic analysis code and its application to lead cooled fast reactor. Appl. Therm. Eng. 117, 443–451. Liu, X.J., Cheng, X., 2015. Sub-channel/system coupled code development and its application to SCWR-FQT loop. Nucl. Eng. Des. 285, 39–47. Oberkampf, W.L., Roy, C.J., 2010. Verification and Validation in Scientific Computing. Cambridge University Press. OpenFOAM Foundation Ltd, 2015. OpenFOAM User’s guide. London (UK). Papukchiev, A., Jeltsov, M., Kööp, K., Kudinov, P., Lerchl, G., 2015. Comparison of different coupling CFD-STH approaches for pre-test analysis of a TALL-3D experiment. Nucl. Eng. Des. 290, 135–143. Pialla, D., Tenchine, D., Li, S., Gauthe, P., Vasile, A., Baviere, R., Tauveron, N., Perdu, F., Maas, L., Cocheme, F., Huber, K., Cheng, X., 2015. Overview of the system alone and system/CFD coupled calculations of the PHENIX Natural Circulation Test within the THINS project. Nucl. Eng. Des. 290, 78–86. Roelofs, F., Gerschenfeld, A., Tarantino, M., Tichelen, K. Van, Pointer, W.D., 2019. Thermal-hydraulic challenges in liquid-metal-cooled reactors. In: Roelofs, Ferry (Ed.), Thermal Hydraulics Aspects of Liquid Metal Cooled Nuclear Reactors. Woodhead Publishing, pp. 17–43. Tenchine, D., 2010. Some thermal hydraulic challenges in sodium cooled fast reactors. Nucl. Eng. Des. 240, 1195–1217. Toumi, I., Bergeron, A., Gallo, D., Royer, E., Caruge, D., 2000. FLICA-4: A threedimensional two-phase flow computer code with advanced numerical methods for nuclear applications. Nucl. Eng. Des. 200, 139–155. Versluis, R.M., 2013. NEAMS Software Verification and Validation Plan Requirements (Version 0). U.S Department of Energy. Vinoski, S., 1997. CORBA: integrating diverse applications within distributed heterogeneous environments. IEEE Commun. Mag. 35, 46–55. Wheeler, C.L., Stewart, C.W., Cena, R.J., Rowe, D.S., Sutey, A.M., 1976. COBRA-IV-I: an interim version of COBRA for thermal-hydraulic analysis of rod bundle nuclear fuel elements and cores. Witherden, F.D., Farrington, A.M., Vincent, P.E., 2014. PyFR: An open source framework for solving advection-diffusion type problems on streaming architectures using the flux reconstruction approach. Comput. Phys. Commun. 185, 3028–3040. Zhang, K.L., Sanchez-espinoza, V.H., 2019. Optimization and verification of the coupled code TRACE/SubChanFlow using the VVER-1000 coolant mixing benchmark data. Nucl. Eng. Des. 353, 110238. Zhang, K.L., Sanchez, V., 2019. In: Instructions on Setting of SUBCHANFLOW and TRACE with ICoCo in SALOME Platform on Debian and Mandriva in VirtualBox. Karlsruhe Institute of Technology. Xu, M. et al., 2011. Thermo-hydraulics of fast reactor (in Chinese). China Atomic Energy Press. Zhang, K.L., 2019. Coupling of TRACE and SubChanFlow based on the Exterior Communication Interface. Prog. Nucl. Energy 103040.