Comparison Plotting Scripts

A successful migration from the solid mechanics module to the tensor mechanics module system is completed when the simulation results are replicated. One the best ways to compare the results of the tensor mechanics simulation to the solid mechanic simulation is with a set of plots.

Example matplotlib Script

The Bison team is in the process of developing a plotting script, with matplotlib. An example of the plotting script used to compare results in the US PWR 16x16 TSQ002 Bison assessment case is shown below.

import numpy as np
from numpy import linalg
from matplotlib import pyplot, rc, patches
import matplotlib.ticker as mtick
import matplotlib

def readTimeSeries( filename, field1, field2 ):
  with open( filename, 'r' ) as f:
    header = f.readline().strip().split(',')
  field1_index = header.index(field1)
  #x_index = header.index()
  field2_index = header.index(field2)
  raw = np.genfromtxt( filename, delimiter=',' )[1:,:]
  field1 = raw[:,field1_index]
  field2 = raw[:,field2_index]
  return field1, field2

offset  = 0.
matplotlib.rcParams.update({'font.size': 16})

###### User must update these values:
# Filenames of the input files used in these simulations:
solid_mechanics_input_file = 'TSQ002_input'
tensor_mechanics_input_file = 'TSQ002_input_tm'
oneptfive_input_file = 'TSQ002_input_1pt5'

# Number of last time step (last four digits on the end of the radial displacement files)
solid_mechanics_last_time_step = '0184'
tensor_mechanics_last_time_step = '0183'
oneptfive_last_time_step = '0168'
##################

# Generate the names of the postprocessor csv output files
sm_csv_file = [solid_mechanics_input_file, '_out.csv']
tm_csv_file = [tensor_mechanics_input_file, '_out.csv']
oneptfive_csv_file = [oneptfive_input_file, '_out.csv']

sm_output_file = ''.join(sm_csv_file)
tm_output_file = ''.join(tm_csv_file)
oneptfive_output_file = ''.join(oneptfive_csv_file)

# Generate the names of the radial displacement files
sm_clad_radial_disp = [solid_mechanics_input_file, '_outfile_1_clad_dia_', solid_mechanics_last_time_step, '.csv']
tm_clad_radial_disp = [tensor_mechanics_input_file, '_outfile_1_clad_dia_', tensor_mechanics_last_time_step, '.csv']
oneptfive_clad_radial_disp = [oneptfive_input_file, '_outfile_1_clad_dia_', oneptfive_last_time_step, '.csv']

sm_fuel_radial_disp = [solid_mechanics_input_file, '_outfile_fuel_radial_displacement_pellet_dia_', solid_mechanics_last_time_step, '.csv']
tm_fuel_radial_disp = [tensor_mechanics_input_file, '_outfile_fuel_radial_displacement_pellet_dia_', tensor_mechanics_last_time_step, '.csv']
oneptfive_fuel_radial_disp = [oneptfive_input_file, '_outfile_fuel_radial_displacement_pellet_dia_', oneptfive_last_time_step, '.csv']

sm_clad_radial_disp_file = ''.join(sm_clad_radial_disp)
tm_clad_radial_disp_file = ''.join(tm_clad_radial_disp)
oneptfive_clad_radial_disp_file = ''.join(oneptfive_clad_radial_disp)

sm_fuel_radial_disp_file = ''.join(sm_fuel_radial_disp)
tm_fuel_radial_disp_file = ''.join(tm_fuel_radial_disp)
oneptfive_fuel_radial_disp_file = ''.join(oneptfive_fuel_radial_disp)

# Build the array of csv files to read the selected data
files = []
files.append(sm_output_file)
files.append(tm_output_file)
files.append(oneptfive_output_file)
print "Gathering values from these files: " + ' and '.join(files)

#Max, Average, and Min Temperatures in the Clad
columns = []
columns.append('max_clad_temp')
columns.append('avg_clad_temp')
columns.append('min_clad_temp')
print "Data collected: " + ' '.join(columns)

#set a figure number integer for the plots to make sure we get a clean figure for each new plot
fig_number = 1

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Set up a single plot on a single figure
pyplot.figure(fig_number)
ax1= pyplot.subplot(111)

# p1 = ax1.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Maximum Clad Temperature')
# p2 = ax1.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Maximum Clad Temperature')
# p3 = ax1.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Maximum Clad Temperature')
p4 = ax1.plot(xvalue[3], yvalue[3], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Average Clad Temperature')
p5 = ax1.plot(xvalue[4], yvalue[4], 'g', marker='>', markevery=10, label='2D-RZ TM Average Clad Temperature')
p6 = ax1.plot(xvalue[5], yvalue[5], 'b', marker='o', markevery=10, label='1.5D Average Clad Temperature')
# p1 = ax1.plot(xvalue[6], yvalue[6], 'goldenrod', linestyle='-.', label='2D-RZ SM Minimum Clad Temperature')
# p2 = ax1.plot(xvalue[7], yvalue[7], 'g', linestyle='-.', label='2D-RZ TM Minimum Clad Temperature')
# p2 = ax1.plot(xvalue[8], yvalue[8], 'b', linestyle='-.', label='1.5D Minimum Clad Temperature')

ax1.legend(loc='best', numpoints = 1, prop={'size':15}, framealpha=0.5)
#ax1.set_xlim([-10000000,141900000])
#ax1.set_ylim([575,650])

ax1.set_ylabel("Temperature (K)")
ax1.set_xlabel('Time (s)')

ax1.set_title('Temperature in the Clad')

ax1.yaxis.tick_left()
ax1.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)

pyplot.savefig('clad_temperature_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

############
#Max, Average, and Min Temperatures in the Fuel
columns = []
columns.append('max_fuel_temp')
columns.append('ave_temp_interior')
columns.append('min_fuel_temp')
print "Data collected: " + ' '.join(columns)

#set a figure number integer for the plots to make sure we get a clean figure for each new plot
fig_number = 1

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Set up a single plot on a single figure
pyplot.figure(fig_number)
ax1= pyplot.subplot(111)

# p1 = ax1.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Maximum Fuel Temperature')
# p2 = ax1.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Maximum Fuel Temperature')
# p3 = ax1.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Maximum Fuel Temperature')
p4 = ax1.plot(xvalue[3], yvalue[3], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Avg Interior Temperature')
p5 = ax1.plot(xvalue[4], yvalue[4], 'g', marker='>', markevery=10, label='2D-RZ TM Avg Interior Temperature')
p6 = ax1.plot(xvalue[5], yvalue[5], 'b', marker='o', markevery=10, label='1.5D Avg Interior Temperature')
# p1 = ax1.plot(xvalue[6], yvalue[6], 'goldenrod', linestyle='-.', label='2D-RZ SM Minimum Fuel Temperature')
# p2 = ax1.plot(xvalue[7], yvalue[7], 'g', linestyle='-.', label='2D-RZ TM Minimum Fuel Temperature')
# p2 = ax1.plot(xvalue[8], yvalue[8], 'b', linestyle='-.', label='1.5D Minimum Fuel Temperature')

ax1.legend(loc='lower center', numpoints = 1, prop={'size':15}, framealpha=0.5)
#ax1.set_xlim([-10000000,141900000])
ax1.set_ylim([200,800])

ax1.set_ylabel("Temperature (K)")
ax1.set_xlabel('Time (s)')

ax1.set_title('Temperature in the Fuel')

ax1.yaxis.tick_left()
ax1.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)

pyplot.savefig('fuel_temperature_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

############
#Fuel centerline temperature
columns = []
columns.append('FCT')
print "Data collected: " + ' '.join(columns)

#set a figure number integer for the plots to make sure we get a clean figure for each new plot
fig_number = 1

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Set up a single plot on a single figure
pyplot.figure(fig_number)
ax1= pyplot.subplot(111)

p1 = ax1.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Center Node Temperature')
p2 = ax1.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Center Node Temperature')
p3 = ax1.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Center Node Temperature')

ax1.legend(loc='best', numpoints = 1, prop={'size':15})
#ax1.set_xlim([-10000000,141900000])
ax1.set_ylim([0,1500])

ax1.set_ylabel("Temperature (K)")
ax1.set_xlabel('Time (s)')

ax1.set_title('Fuel Centerline Temperature')

ax1.yaxis.tick_left()
ax1.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)

pyplot.savefig('fuel_centerline_temperature_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

############
# Plenum Pressure
#headers of the columns which we want to plot
columns = []
columns.append('plenum_pressure')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Plenum Pressure')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g-', marker='>', markevery=10, label='2D-RZ TM Plenum Pressure')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b-', marker='o', markevery=10, label='1.5D Plenum Pressure')

ax2.legend(loc='lower center', numpoints = 1, prop={'size':15})
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([2000000,3400000])

ax2.set_ylabel("Plenum Pressure (Pa)")
ax2.set_xlabel('Time (s)')
ax2.set_title('Plenum Pressure ')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('plenum_pressure_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

#########
# Fuel vonMises Stress
columns = []
columns.append('vonmises_stress_fuel')

print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Fuel von Mises Stress')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Fuel von Mises Stress')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Fuel von Mises Stress')

ax2.legend(loc='lower right', numpoints = 1, prop={'size':15}, framealpha=0.5)
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([0.5e9,1.1e9])

ax2.set_ylabel("Average von Mises Stress (Pa)")
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))
ax2.set_xlabel('Time (s)')
ax2.set_title('Average Stress in the Fuel')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('fuel_vonMisesStress_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

#########
# Average Clad vonMises Stress
columns = []
columns.append('vonmises_stress_clad')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Clad von Mises Stress')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Clad von Mises Stress')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Clad von Mises Stress')

ax2.legend(loc='best', numpoints = 1, prop={'size':15})
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([0.5e9,1.1e9])

ax2.set_ylabel("Average von Mises Stress (Pa)")
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))
ax2.set_xlabel('Time (s)')
ax2.set_title('Average Stress in the Clad')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('clad_vonMisesStress_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

#########
# Fuel Contact Penetration
columns = []
# columns.append('max_penetration')
columns.append('gap')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

# p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Maximum Contact Penetration')
# p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Maximum Contact Penetration')
# p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Maximum Contact Penetration')
p4 = ax2.plot(xvalue[0], -1*yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Midpoint Gap')
p5 = ax2.plot(xvalue[1], -1*yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Midpoint Gap')
p5 = ax2.plot(xvalue[2], -1*yvalue[2], 'b', marker='o', markevery=10, label='1.5D Midpoint Gap')

ax2.legend(loc='best', numpoints = 1, prop={'size':15})
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([0.5e9,1.1e9])

ax2.set_ylabel("Pellet-Clad Gap (m)")
ax2.set_xlabel('Time (s)')
ax2.set_title('Pellet-Clad Gap Evolution')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('fuel_gap_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

#########
# Fuel Nodal Contact Pressure
columns = []
columns.append('contact_pressure')
# columns.append('min_contact_pressure')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D')
# p4 = ax2.plot(xvalue[3], yvalue[3], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Minimum Contact Pressure')
# p5 = ax2.plot(xvalue[4], yvalue[4], 'g--', label='2D-RZ TM Minimum Contact Pressure')
# p6 = ax2.plot(xvalue[5], yvalue[5], 'b--', label='1.5D Minimum Contact Pressure')

ax2.legend(loc='upper left', numpoints = 1, prop={'size':15})
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([0.5e9,1.1e9])

ax2.set_ylabel("Contact Pressure (Pa)")
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))
ax2.set_xlabel('Time (s)')
ax2.set_title('Contact Pressure at Midpoint')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('fuel_contactPressure_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

############
# Fission Gas Release
columns = []
columns.append('fis_gas_released')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Gas Released')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Gas Released')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Gas Released')

ax2.legend(loc='upper left', numpoints = 1, prop={'size':15}, framealpha=0.5)
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([0,4])

ax2.set_ylabel("Fission Gas Release Percent")
ax2.set_xlabel('Time (s)')
ax2.set_title('Percent Fission Gas Released')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('fission_gas_release_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

############
# Fission Gas Generated
columns = []
columns.append('fis_gas_generated')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Generated Gas')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Generated Gas')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Generated Gas')

ax2.legend(loc='upper left', numpoints = 1, prop={'size':15}, framealpha=0.5)
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([0,4])

ax2.set_ylabel("Fission Gas Generated")
ax2.set_xlabel('Time (s)')
ax2.set_title('Fission Gas Generated ')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('fission_gas_generated_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

#########
# Fuel Pellet volume
columns = []
columns.append('pellet_volume')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Fuel pellet volume')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ Fuel pellet volume')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Fuel pellet volume')

ax2.legend(loc='best', numpoints = 1, prop={'size':15})
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([-6.2e-6,-6.6e-6])

ax2.set_ylabel("Fuel Pellet Volume ($m^3$)")
ax2.set_xlabel('Time (s)')
ax2.set_title('Fuel Pellet Volume Change')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('swelling_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

#########
# Clad internal volume
columns = []
columns.append('clad_inner_vol')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Clad Internal Volume')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Clad Internal Volume')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Clad Internal Volume')

ax2.legend(loc='lower right', numpoints = 1, prop={'size':15})
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([-6.2e-6,-6.6e-6])

ax2.set_ylabel("Clad Internal Volume ($m^3$)")
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))
ax2.set_xlabel('Time (s)')
ax2.set_title('Clad Volume Change')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('clad_internalVolume_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

#########
# Gas internal volume
columns = []
columns.append('gas_volume')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Gas Volume')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Gas Volume')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Gas Volume')

ax2.legend(loc='center right', numpoints = 1, prop={'size':15})
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([-6.2e-6,-6.6e-6])

ax2.set_ylabel("Gas Volume ($m^3$)")
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))
ax2.set_xlabel('Time (s)')
ax2.set_title('Gas Volume Change')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('gas_internalVolume_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

#########
# Flux Measures
#headers of the columns which we want to plot
columns = []
columns.append('flux_from_fuel')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Flux from Fuel')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ TM Flux from Fuel')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Flux from Fuel')

ax2.legend(loc='lower center', numpoints = 1, prop={'size':15})
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([2000,3000])

ax2.set_ylabel("Integrated Temperature Flux (K/m^2)")
ax2.set_xlabel('Time (s)')
ax2.set_title('Temperature Flux ')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('fuel_temperatureFlux_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

#########
# Power history
columns = []
columns.append('rod_input_power')
columns.append('rod_total_power')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

# p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod',  linestyle = '--', label='2D-RZ SM Input Power History')
# p2 = ax2.plot(xvalue[1], yvalue[1], 'g--', label='2D-RZ TM Input Power History')
# p3 = ax2.plot(xvalue[2], yvalue[2], 'b--', label='1.5D Input Power History')
p4 = ax2.plot(xvalue[3], yvalue[3], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Calculated Power History')
p5 = ax2.plot(xvalue[4], yvalue[4], 'g', marker='>', markevery=10, label='2D-RZ TM Calculated Power History')
p5 = ax2.plot(xvalue[5], yvalue[5], 'b', marker='o', markevery=10, label='1.5D Calculated Power History')

ax2.legend(loc='best', numpoints = 1, prop={'size':15}, framealpha=0.5)
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([0.5e9,1.1e9])

ax2.set_ylabel("Total Power (W)")
ax2.set_xlabel('Time (s)')
ax2.set_title('Power History')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('powerHistory_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()

############
# Run time for simulation
#headers of the columns which we want to plot
columns = []
columns.append('alive_time')
print "Data collected: " + ' '.join(columns)

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'time', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', marker='D', markevery=10, label='2D-RZ SM Runtime (24 processors)')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', marker='>', markevery=10, label='2D-RZ Runtime (24 processors)')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', markevery=10, label='1.5D Runtime (1 processor)')

ax2.legend(loc='upper left', numpoints = 1, prop={'size':15}, framealpha=0.5)
#ax2.set_xlim([-10000000,141900000])
#ax2.set_ylim([0,2000])## for discrete use 1500

ax2.set_ylabel('Simulation Runtime (s)')
ax2.set_xlabel('Time (s)')
ax2.set_title('Simulation Total Runtime ')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('runtime_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Clear the figure of the current plot to make room for a new plot
pyplot.clf()
#########
# Radial Displacements
#Filenames for the csv files to use for the comparison
files = []
files.append(sm_clad_radial_disp_file)
files.append(tm_clad_radial_disp_file)
files.append(oneptfive_clad_radial_disp_file)
files.append(sm_fuel_radial_disp_file)
files.append(tm_fuel_radial_disp_file)
files.append(oneptfive_fuel_radial_disp_file)

print "Now gathering radial displacement values from these files: " + ' and '.join(files)

#headers of the columns which we want to plot
columns = []
columns.append('disp_x')

#set up x and y lists to match the files
xvalue = [None] * len(files) * len(columns)
yvalue = [None] * len(files) * len(columns)

#Loop over the filenames to print out the same values in different folders
for i in xrange(len(columns)):
  for j in xrange(len(files)):
    print j+len(files)*i
    xvalue[j+len(files)*i], yvalue[j+len(files)*i] = readTimeSeries(files[j], 'y', columns[i])

# Using the matplotlib tutorial now
pyplot.figure(fig_number)
ax2= pyplot.subplot(111)

p1 = ax2.plot(xvalue[0], yvalue[0], 'goldenrod', label='2D-RZ SM Clad Displacement')
p2 = ax2.plot(xvalue[1], yvalue[1], 'g', label='2D-RZ TM Clad Displacement')
p3 = ax2.plot(xvalue[2], yvalue[2], 'b', marker='o', linestyle='none', label='1.5D Clad Displacement')
p4 = ax2.plot(xvalue[3], yvalue[3], 'goldenrod', linestyle = '--', label='2D-RZ SM Fuel Displacement')
p5 = ax2.plot(xvalue[4], yvalue[4], 'g--', label='2D-RZ TM Fuel Displacement')
p5 = ax2.plot(xvalue[5], yvalue[5], 'b', marker='s', linestyle='none', label='1.5D Fuel Displacement')

ax2.legend(loc='lower left', numpoints = 1, prop={'size':15}, framealpha=0.5)
ax2.set_xlim([-0.2,4.5])
ax2.set_ylim([-2.0e-4,1.0e-4])

ax2.set_ylabel("Radial Displacement (m)")
ax2.yaxis.set_major_formatter(mtick.FormatStrFormatter('%.1e'))
ax2.set_xlabel('Axial Position along fuel rod (m)')
ax2.set_title('Simulation End Radial Displacement')
ax2.yaxis.tick_left()
ax2.tick_params(labeltop='off')

pyplot.subplots_adjust(wspace=0.15)
pyplot.subplots_adjust(bottom=0.12)
pyplot.subplots_adjust(top=0.9)
pyplot.subplots_adjust(left=0.12)
pyplot.subplots_adjust(right=0.88)
pyplot.savefig('radial_displacement_TSQ002_comparison.pdf',bbox_inches='tight',pad_inches=0.1)

# Close all of the open figures
pyplot.close()
(assessment/US_PWR_16_x_16/analysis/TSQ002/plot_TSQ002_mechanics_comparison.py)

Note that this script assumes the use of Postprocessor names which have been added to the TSQ002 input files.


[Postprocessors]
  [./ave_temp_interior]
    type = SideAverageValue
    boundary = 9
    variable = temp
    execute_on = 'initial linear'
  [../]
  [./clad_inner_vol]
    type = InternalVolume
    boundary = 7
  [../]
  [./pellet_volume]
    type = InternalVolume
    boundary = 8
  [../]
  [./avg_clad_temp]
    type = SideAverageValue
    boundary = 7
    variable = temp
    execute_on = 'initial timestep_end'
  [../]
  [./max_fuel_temp]
    type = NodalExtremeValue
    block = pellet_type_1
    value_type = max
    variable = temp
    execute_on = 'initial timestep_end'
  [../]
  [./min_fuel_temp]
    type = NodalExtremeValue
    block = pellet_type_1
    value_type = min
    variable = temp
    execute_on = 'initial timestep_end'
  [../]
  [./max_clad_temp]
    type = NodalExtremeValue
    block = clad
    value_type = max
    variable = temp
    execute_on = 'initial timestep_end'
  [../]
  [./min_clad_temp]
    type = NodalExtremeValue
    block = clad
    value_type = min
    variable = temp
    execute_on = 'initial timestep_end'
  [../]
  [./fis_gas_generated]
    type = ElementIntegralFisGasGeneratedSifgrs
    variable = temp
    block = pellet_type_1
  [../]
  [./fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    variable = temp
    block = pellet_type_1
  [../]
  [./fis_gas_grain]
    type = ElementIntegralFisGasGrainSifgrs
    variable = temp
    block = pellet_type_1
    outputs = exodus
  [../]
  [./fis_gas_boundary]
    type = ElementIntegralFisGasBoundarySifgrs
    variable = temp
    block = pellet_type_1
    outputs = exodus
  [../]
  [./gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  [../]
  [./flux_from_clad]
    type = SideFluxIntegral
    variable = temp
    boundary = 5
    diffusivity = thermal_conductivity
  [../]
  [./flux_from_fuel]
    type = SideFluxIntegral
    variable = temp
    boundary = 10
    diffusivity = thermal_conductivity
  [../]
  [./_dt]
    type = TimestepSize
  [../]
  [./num_lin_it]
    type = NumLinearIterations
  [../]
  [./num_nonlin_it]
    type = NumNonlinearIterations
  [../]
  [./tot_lin_it]
    type = CumulativeValuePostprocessor
    postprocessor = num_lin_it
  [../]
  [./tot_nonlin_it]
    type = CumulativeValuePostprocessor
    postprocessor = num_nonlin_it
  [../]
  [./alive_time]
    type = PerformanceData
    event = ALIVE
  [../]
  [./rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet_type_1
  [../]
  [./rod_input_power]
    type = FunctionValuePostprocessor
    function = power_history
    scale_factor = 3.81381 # rod height
  [../]
  [./average_fission_rate]
    type = ElementAverageValue
    variable = fission_rate
    block = pellet_type_1
  [../]
  [./average_burnup]
    type = ElementAverageValue
    variable = burnup
    block = pellet_type_1
  [../]
  [./FCT]
    type = NodalVariableValue
    nodeid = 30330 #coords (0.0, 2.10133)
    variable = temp
    execute_on = 'initial timestep_end'
  [../]
  [./FCT_slice4]
    type = NodalVariableValue
    nodeid = 37085 #coords (0.0, 1.71896)
    variable = temp
    execute_on = 'initial timestep_end'
  [../]
  [./fis_gas_percent]
    type = FGRPercent
    value1 = fis_gas_released
    value2 = fis_gas_generated
  [../]
  [./vonmises_stress_fuel]
    type = ElementAverageValue
    block = pellet_type_1
    variable = vonmises_stress
  [../]
  [./vonmises_stress_clad]
    type = ElementAverageValue
    block = clad
    variable = vonmises_stress
  [../]

  ## Nodal comparison values
  [./gap_slice6]
    type = NodalVariableValue
    variable = penetration
    nodeid = 23579 #coords (0.0041275, 2.48172)
  [../]
  [./gap]
    type = NodalVariableValue
    variable = penetration
    nodeid = 30299 #coords (0.0041275, 2.10133)
  [../]
  [./gap_slice4]
    type = NodalVariableValue
    variable = penetration
    nodeid = 37054 #coords (0.0041275, 1.71896)
  [../]
  [./contact_pressure_slice6]
    type = NodalVariableValue
    variable = contact_pressure
    nodeid = 23579 #coords (0.0041275, 2.48172)
  [../]
  [./contact_pressure]
    type = NodalVariableValue
    variable = contact_pressure
    nodeid = 30299 #coords (0.0041275, 2.10133)
  [../]
  [./contact_pressure_slice4]
    type = NodalVariableValue
    variable = contact_pressure
    nodeid = 37054 #coords (0.0041275, 1.71896)
  [../]
[]
(assessment/US_PWR_16_x_16/analysis/TSQ002/TSQ002_tm.i)

Users may need to change the names of their existing postprocessors and add additional postprocessors to use this matplotlib script.