MSCS for coupled transients

The Minimal Serpent Coupling Script (MSCS) for transient simulations is a short (< 400 lines with comments) Python program intended to give a minimal working example of a wrapper program that can communicate with Serpent in the coupled transient calculation mode. MSCS provides a working example of externally coupled multi-physics simulations with Serpent and may be a good starting point for the users that are interested in running such simulations with Serpent.

Description

The coupling script communicates with Serpent using the file based communication mode.

Here, the coupling script calculates the fuel temperature solution itself. In many cases that part of the coupling script should be replaced by writing the input for an external solver, running the external solver and reading the results from the external solver output.

The temperature treatment of the interaction physics is done on-the-fly using the TMS temperature treatment technique for the base cross sections and interpolation of thermal scattering data for the thermal scattering libraries.

The problem solved here is a 200 cm long fuel rod in infinite lattice with axially black boundary conditions. The time dependent fuel temperature is solved by MSCS assuming no heat conduction out of fuel.

Files

MSCS.py

#############################################################
#                                                           #
#          Minimal Serpent Coupling Script v 0.3            #
#          For transient example                            #
#                                                           #
# Created by:  Ville Valtavirta              2016/05/06     #
# Modified by: Ville Valtavirta              2020/10/26     #
#                                                           #
#############################################################

import os
import signal
import math
import time

# Path to Serpent executable

sssexe = '/home/vvvillehe/Serpent2/2.1.31/sss2'

#######################################################
# Create the Serpent input-file for this run          #
# (process id or communication file must be appended) #
#######################################################

# Open original input for reading

file_in = open('./input','r')

# Open a new input file for writing

file_out = open('./coupledinput','w')

# Write original input to new file

for line in file_in:
    file_out.write(line)

# Close original input file

file_in.close()

# Append signalling mode

file_out.write('\n')
file_out.write('set comfile com.in com.out\n')

# Append interface names

file_out.write('\n')
file_out.write('ifc cool.ifc\n\n')
file_out.write('ifc fuel.ifc\n')

# Close new input file

file_out.close()

##############################################
# Write the initial coolant interface file   #
# (Coolant conditions will be held constant) #
##############################################

file_out = open('./cool.ifc','w')

# Write the header line (TYPE MAT OUT)

file_out.write('2 cool 1\n')

# Write the output line (OUTFILE NZ ZMIN ZMAX NR)

file_out.write('coolifc.out 10 -100 100 1\n')

# Write the mesh type

file_out.write('1\n')

# Write the mesh size (NX XMIN XMAX NY YMIN YMAX NZ ZMIN ZMAX)

file_out.write('1 -0.75 0.75 1 -0.75 0.75 10 -100 100\n')

# Write initial coolant temperatures and densities

for i in range(10):
    file_out.write('-0.9  520.0\n')

# Close interface file

file_out.close()

##############################################
# Write the initial fuel interface file      #
# (Fuel temperature will be updated)         #
##############################################

file_out = open('./fuel.ifc','w')

# Write the header line (TYPE MAT OUT)

file_out.write('2 fuel 0\n')

# Write the mesh type

file_out.write('1\n')

# Write the mesh size (NX XMIN XMAX NY YMIN YMAX NZ ZMIN ZMAX)

file_out.write('1 -0.75 0.75 1 -0.75 0.75 10 -100 100\n')

# Write initial fuel temperatures and densities

for i in range(10):
    file_out.write('-10.424 520.0\n')

# Close interface file

file_out.close()

# Archive the initial fuel interface

os.system('cp ./fuel.ifc ./fuel.ifc0')

################################
# Start the Serpent simulation #
################################

# Create a command string that will start the Serpent simulation

runcommand = sssexe+' -omp 3 ./coupledinput &'

# Execute the command string

os.system(runcommand)

############################################
# Initialize the fuel temperature solution #
############################################

TBOI = []
TEOI = []

for i in range(10):
    TBOI.append(520.0)
    TEOI.append(520.0)

# Reset time step

curtime = 0

########################
# Loop over time steps #
########################

simulating = 1

while simulating == 1:

    #########################
    # Picard iteration loop #
    #########################

    iterating = 1

    while iterating == 1:
        ###################
        # Wait for signal #
        ###################

        sleeping = 1

        while sleeping == 1:

            # Sleep for two seconds

            time.sleep(2)

            # Open file to check if we got a signal

            fin = open('./com.out','r')

            # Read line

            line = fin.readline()

            # Close file

            fin.close()

            # Check signal

            if int(line) != -1:
                if int(line) == signal.SIGUSR1.value:
                    # Got the signal to resume

                    sleeping = 0

                elif int(line) == signal.SIGUSR2.value:
                    # Got the signal to move to next time point

                    iterating = 0
                    sleeping = 0
                elif int(line) == signal.SIGTERM.value:
                    # Got the signal to end the calculation

                    iterating = 0
                    sleeping = 0
                    simulating = 0
                else:
                    # Unknown signal

                    print("\nUnknown signal read from file, exiting\n")

                    # Exit

                    quit()

                # Reset the signal in the file

                file_out = open('./com.out','w')

                file_out.write('-1')

                file_out.close()

        # Check if simulation has finished and break out of iterating
        # loop

        if (simulating == 0):
            break

        ###########################
        # Read power distribution #
        ###########################

        # Reset power distribution

        P = []

        # Open output file

        file_in = open('./coolifc.out','r')

        # Loop over output file to read power distribution

        for line in file_in:
            # Split line to values

            strtuple = line.split()

            # Store power

            P.append(float(strtuple[8]))

        file_in.close()

        ###########################
        # Calculate TH-solution   #
        ###########################

        # Fuel specific heat capacity

        cp = 300 # J/(kg*K)

        # Calculate EOI temperatures at (nz) axial nodes
        # No heat transfer, just deposition

        for i in range(10):

            # Calculate EOI temperature based on BOI temperature
            # and energy deposition during current time interval

            # Calculate mass of this node (in kg)

            m = (math.pi*0.4335**2*20.0)*10.424*1e-3

            # Calculate initial heat in this axial node

            Q = TBOI[i]*(cp*m)

            # The interface output is joules in case of time dependent
            # simulation, no need to multiply with time step

            dQ = P[i]

            # Calculate new temperature based on new amount of heat

            TEOI[i] = (Q + dQ)/(cp*m)

        ###########################
        # Update interface        #
        ###########################

        file_out = open('./fuel.ifc','w')

        # Write the header line (TYPE MAT OUT)

        file_out.write('2 fuel 0\n')

        # Write the mesh type

        file_out.write('1\n')

        # Write the mesh size (NX XMIN XMAX NY YMIN YMAX NZ ZMIN ZMAX)

        file_out.write('1 -0.75 0.75 1 -0.75 0.75 10 -100 100\n')

        # Write updated fuel temperatures

        for i in range(10):
            # Use the base density throughout the simulation
            # Write density and temperature at this layer

            file_out.write('-10.424 {}\n'.format(TEOI[i]))

        file_out.close()

        ############################
        # Signal Serpent (SIGUSR1) #
        ############################

        ##########################################################
        # We could actually check if the temperature changed by  #
        # a significant margin and send a SIGUSR2 to indicate    #
        # that Serpent should not iterate the transport solution #
        ##########################################################

        file_out = open('./com.in','w')

        file_out.write(str(signal.SIGUSR1.value))

        file_out.close()

    # Increment time step

    curtime += 1

    ########################
    # Do some archiving... #
    ########################

    # Copy the interface and interface output to
    # archive files for later examination
    # Note: These are EOI fields

    os.system('cp ./fuel.ifc ./fuel.ifc{:d}'.format(curtime))
    os.system('cp ./coolifc.out ./coolifc.out{:d}'.format(curtime))

    ####################################
    # Check if simulation has finished #
    ####################################

    if (simulating == 0):
        break

    ############################
    # Moving to next time step #
    ############################

    # Copy EOI temperatures to BOI vector

    for i in range(10):
        TBOI[i] = TEOI[i]

    ##################################################
    # Here we could calculate an initial guess for   #
    # the the EOI temperatures of the next time      #
    # interval and update the interface.             #
    # But to keep the script minimal, we won't do it #
    ##################################################

    ############################
    # Signal Serpent (SIGUSR1) #
    ############################

    file_out = open('./com.in','w')

    file_out.write(str(signal.SIGUSR1.value))

    file_out.close()

input

% --- Input for MSCS testing

set title "Serpent-MSCS externally coupled transient calculation"

% --- Fuel Pin definition:

pin 1
fuel   0.4335
gas    0.442000
clad   0.502500
cool

% --- Lattice (type = 1, pin pitch = 1.5):

lat 10  1  0.0 0.0 1 1 1.5
1

% --- Boundary of the geometry (200 cm active length)

surf 2 cuboid -0.75 0.75 -0.75 0.75 -100 100

% --- Cell definitions:

cell  3  0  fill 10  -2     % Pin-cell
cell 99  0  outside   2     % Outside world

% --- Fuel material:

mat fuel    -10.424 tft 520.0 3000.0 rgb 100 160 140
 92235.03c   -0.015867
 92238.03c   -0.86563
  8016.03c   -0.1185

% --- Cladding material:

mat clad     -6.55 rgb 200 200 200
 40090.03c   -0.98135
 24052.03c   -0.00100
 26056.03c   -0.00135
 28058.03c   -0.00055
 50120.03c   -0.01450
  8016.03c   -0.00125

% --- Gas gap:

mat gas -1E-4 rgb 255 255 255
 2004.06c     1.0

% --- Coolant (Temperature won't change it will be given via interface):

mat cool     -0.90 moder lwtr 1001 rgb 150 160 240
 1001.03c     0.66667
 8016.03c     0.33333
 5010.03c     0.0001

% --- Thermal scattering data for light water:
%     On-the-fly treatment for SAB-data between 474 K -- 624 K

therm lwtr 520 lwj3.07t lwj3.09t

% --- Reflective boundary conditions in XY

set bc 3 3 1

% --- Initial neutron source
%     A transient source generated with a separate criticality
%     source simulation should typically be used

% --- A point source of 1 MeV neutrons at (0,0,0)

src init sp 0.0 0.0 0.0 se 1.0

% --- Simulation time structure (100 intervals between 0 and 5e-3 s)

tme tsim 2 100 0 5e-3

% --- Total neutron population and number of batches and time structure

set nps 20000 10 tsim

% --- Maximum number of coupled calculation iterations set to 2

set ccmaxiter 2

% --- Initial power for normalization (600 kW/cm, unrealistically high):
%     When separately generated transient source is used, the
%     normalization will correspond to that of the source generation
%     simulation. Here the power will be normalized to the power generation
%     during the first time interval

set power 120000000.00

% --- Fission heating detector

det HEAT dr -8 void

% --- Disable group constant generation

set gcu -1

Setup

MSCS has been tested with Python 3.8.13 and Serpent 2.2.0.

1. Save the contents of the MSCS.py file (above) to a file of the same name.
2. Replace the absolute path to the Serpent executable on line 18 of MSCS.py.
3. Save the contents of the input file to a file called input in the same directory where the MSCS.py file is located.
4. Add some cross section libraries to the input-file using the set acelib input option.
5. The test case is now ready to run.

Note: Before running the MSCS you should familiarize yourself with finding and killing background processes in your operating system. Since MSCS runs Serpent as a background process, killing MSCS will not kill the Serpent process automatically.

Run the simulation from the folder containing both files using the command python MSCS.py

The output of the Serpent run will be printed to the terminal.