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.2 # # For transient example # # # # Created by: Ville Valtavirta 2016/05/06 # # Modified by: Ville Valtavirta 2016/06/14 # # # ############################################################# import os import signal import math import time # Path to Serpent executable sssexe = '/home/vvvillehe/Serpent2/2.1.27/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: # Got the signal to resume sleeping = 0 elif int(line) == signal.SIGUSR2: # Got the signal to move to next time point iterating = 0 sleeping = 0 elif int(line) == signal.SIGTERM: # 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])) ########################### # 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)) 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)) file_out.close()