Update file paths (5abbf5a8) · Commits · SCALE / Analysis / 2025 gFHR

_analysis/1gxs/one-group-xs-vs-burnup.ipynb

+2 −603

Original line number	Diff line number	Diff line
		@@ -106,610 +106,9 @@
		},
		{
		"cell_type": "code",
		"execution_count": 25,
		"execution_count": null,
		"metadata": {},
		"outputs": [
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"0\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"0\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"0\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"0\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"1\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"1\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"1\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"1\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"2\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"2\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"2\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"2\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"3\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"3\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"3\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"3\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"4\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"4\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"4\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"4\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"5\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"5\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"5\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"5\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"6\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"6\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"6\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"6\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"7\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"7\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"7\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"7\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"8\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"8\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"8\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"8\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"9\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"9\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"9\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		},
		{
		"name": "stdout",
		"output_type": "stream",
		"text": [
		"9\n"
		]
		},
		{
		"name": "stderr",
		"output_type": "stream",
		"text": [
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n",
		"Loaded Scale kappa library /Users/f78/scale-data/kappa.h5\n"
		]
		}
		],
		"outputs": [],
		"source": [
		"losxs = {}\n",
		"\n",

_analysis/nuclide-densities/zone-wise-dens-comp-with-benchmark-perturbed-violin.ipynb

+20 −19

Original line number	Diff line number	Diff line
		%% Cell type:markdown id: tags:

		# Visualize zone-wise nuclide comparison to benchmark via violin plot

		%% Cell type:code id: tags:

		``` python
		import seaborn as sns
		import matplotlib.pyplot as plt
		from matplotlib import rcParams
		import pandas as pd
		import json
		import os
		import re
		import glob
		import numpy as np
		from lib_utils import get_densities_from_f71

		acts = ["u235", "u238", "pu239", "pu240", "pu241", "pu242", "am241", "am242", \
		"am242m", "cm244", "cm245"]
		fps = ["kr85",
		"sr90",
		"ag110m",
		"i131",
		"xe135",
		"cs133", "cs134", "cs137",
		"nd143", "nd145", "nd148",
		"sm147", "sm149", "sm150", "sm151", "sm152",
		"eu153", "eu154", "eu155",
		"gd154", "gd155", "gd157"
		]

		stdcomp_fname = "/Users/f78/scale_dev/packages/InMemoryStdCompLib/InMemoryStdCompLibrary.json"
		# Download file from: https://code.ornl.gov/scale/code/scale-public/-/blob/main/packages/InMemoryStdCompLib/InMemoryStdCompLibrary.json
		stdcomp_fname = "InMemoryStdCompLibrary.json"
		```

		%% Cell type:markdown id: tags:

		### Benchmark composition

		%% Cell type:code id: tags:

		``` python
		def get_nuc_name_id_converter():
		"""
		Get converter between SCALE nuclide id and SCALE nuclide name

		Returns
		=======
		converter : dict
		key: SCALE id
		value: SCALE name
		"""

		assert os.path.isfile(stdcomp_fname)

		f = open(stdcomp_fname, "r")
		stdcmp = json.load(f)
		f.close()

		converter = {}
		for nuc in stdcmp['compositions']:
		if len(nuc['elements']) > 1:
		continue
		for el in nuc['elements']:
		converter[el['elementid']] = nuc['name']

		return converter
		```

		%% Cell type:code id: tags:

		``` python
		# Read Titan input
		titan_input_fname = "/Users/f78/neams/physor2024-pbr-performance/titan-gfhr/gfhr-tally-on.json"
		assert os.path.isfile(titan_input_fname)
		f = open(titan_input_fname, "r")
		titan_input = json.load(f)
		# Read KP compositions from file
		kp_material_input_fname = "../../_input_preparation/gfhr_benchmark_comps.json"
		assert os.path.isfile(kp_material_input_fname)
		f = open(kp_material_input_fname, "r")
		kp_material_input = json.load(f)
		f.close()

		# Get SCALE id/name converter
		nuc_id_to_name = get_nuc_name_id_converter()

		# get list of nuclide ids that are included in relevant nuclide list
		titan_nuc_list = []
		for mat_name, items in titan_input['materials'].items():
		kp_nuc_list = []
		for mat_name, items in kp_material_input['materials'].items():
		regExp = re.compile(r"FZ(?P<z>\d+)R(?P<r>\d+)")
		if not re.match(regExp, mat_name):
		continue
		assert len(items) == 2, f"Expected material vector to have length of 2, got {len(items)}"
		assert items[0] == 'CONC', f"Expected first value in vector to be 'CONC', got {items[0]}"
		titan_nuc_dens = items[1][1]
		for id, dens in titan_nuc_dens:
		kp_nuc_dens = items[1][1]
		for id, dens in kp_nuc_dens:
		assert str(id)[-1] == '0', f"Expected nuclide id to end in 0, got {id}"
		nuc_name = nuc_id_to_name[str(int(id/10))]
		nuc_name = f"{nuc_name.split('-')[0]}{nuc_name.split('-')[1]}"
		if nuc_name not in titan_nuc_list:
		titan_nuc_list.append(nuc_name)
		if nuc_name not in kp_nuc_list:
		kp_nuc_list.append(nuc_name)

		gfhr_nuc_list = []
		for nuc in acts+fps:
		if nuc == "nd148":
		print("found it")
		if nuc in titan_nuc_list and nuc not in gfhr_nuc_list:
		if nuc in kp_nuc_list and nuc not in gfhr_nuc_list:
		gfhr_nuc_list.append(nuc)
		if nuc == "nd148":
		print("added it")

		print(len(titan_nuc_list))
		print(len(kp_nuc_list))
		print(len(gfhr_nuc_list))

		print(gfhr_nuc_list)

		# Get list of SCALE materials, with nuclide ids correctly
		# converted from Titan to SCALE
		# converted from JSON file to SCALE
		benchmark_comps = {}
		assert 'materials' in titan_input
		for mat_name, items in titan_input['materials'].items():
		assert 'materials' in kp_material_input
		for mat_name, items in kp_material_input['materials'].items():

		regExp = re.compile(r"FZ(?P<z>\d+)R(?P<r>\d+)")
		if not re.match(regExp, mat_name):
		continue
		z = int(re.match(regExp, mat_name).group('z'))
		r = int(re.match(regExp, mat_name).group('r'))

		assert len(items) == 2, f"Expected material vector to have length of 2, got {len(items)}"
		assert items[0] == 'CONC', f"Expected first value in vector to be 'CONC', got {items[0]}"

		# get all the nuclide densities
		titan_nuc_dens = items[1][1]
		kp_nuc_dens = items[1][1]
		scale_nuc_dens = {}
		for id, dens in titan_nuc_dens:
		for id, dens in kp_nuc_dens:
		assert str(id)[-1] == '0', f"Expected nuclide id to end in 0, got {id}"
		nuc_name = nuc_id_to_name[str(int(id/10))]
		nuc_name = f"{nuc_name.split('-')[0]}{nuc_name.split('-')[1]}"
		scale_nuc_dens[nuc_name] = dens

		# filter for nuclides we want to compare
		target_comps = {}
		for nuc in gfhr_nuc_list:
		if nuc in scale_nuc_dens:
		target_comps[nuc] = scale_nuc_dens[nuc]
		else:
		print(f"Setting zero density for {nuc} in r/z {r}/{z}")
		target_comps[nuc] = 0

		benchmark_comps[(r,z)] = target_comps

		```

		%% Output

		found it
		added it
		89
		24
		['u235', 'u238', 'pu239', 'pu240', 'pu241', 'pu242', 'am241', 'sr90', 'xe135', 'cs133', 'cs134', 'cs137', 'nd143', 'nd145', 'nd148', 'sm147', 'sm149', 'sm150', 'sm151', 'sm152', 'eu153', 'eu154', 'eu155', 'gd157']
		Setting zero density for gd157 in r/z 1/1
		Setting zero density for gd157 in r/z 2/1
		Setting zero density for gd157 in r/z 3/1
		Setting zero density for gd157 in r/z 4/1
		Setting zero density for gd157 in r/z 4/2
		Setting zero density for gd157 in r/z 4/3
		Setting zero density for gd157 in r/z 4/10

		%% Cell type:markdown id: tags:

		### SLICE final composition

		%% Cell type:code id: tags:

		``` python
		slice_zone_f71_list = glob.glob('../../tmp/eq_core_generation_r4z10_new/outer-it-5/arp_origen_output_it5_res-5pct/Z*.f71')
		slice_comps = {}

		for file in sorted(slice_zone_f71_list):
		df = get_densities_from_f71(file)
		dir = file.split("/")[4]
		subdir = file.split("/")[5]
		zone = file.split("/")[6]
		it_outer = re.search(r"it-(?P<it>\d+)",dir).group("it")
		it_inner = subdir[-1]
		Z = re.search(re.compile(r"Z(?P<Z>\d+)_R(?P<R>\d+)"),zone).group("Z")
		R = re.search(re.compile(r"Z(?P<Z>\d+)_R(?P<R>\d+)"),zone).group("R")
		Z = int(Z)
		R = int(R)
		target_comps = {}
		for nuc in gfhr_nuc_list:
		target_comps[nuc] = df.loc[nuc].values[0]
		slice_comps[(R,Z+1)] = target_comps
		```

		%% Cell type:markdown id: tags:

		### Determine relative difference

		%% Cell type:code id: tags:

		``` python
		df_slice = pd.DataFrame.from_dict(slice_comps)
		print(df_slice.loc["u235"])

		df_benchmark = pd.DataFrame.from_dict(benchmark_comps)
		print(df_benchmark.loc["u235"])
		```

		%% Output

		1 1 0.002561
		2 1 0.002563
		3 1 0.002567
		4 1 0.002568
		1 2 0.002496
		2 2 0.002502
		3 2 0.002513
		4 2 0.002514
		1 3 0.002427
		2 3 0.002437
		3 3 0.002455
		4 3 0.002456
		1 4 0.002359
		2 4 0.002372
		3 4 0.002396
		4 4 0.002397
		1 5 0.002295
		2 5 0.002311
		3 5 0.002341
		4 5 0.002340
		1 6 0.002237
		2 6 0.002255
		3 6 0.002290
		4 6 0.002289
		1 7 0.002186
		2 7 0.002207
		3 7 0.002246
		4 7 0.002244
		1 8 0.002144
		2 8 0.002166
		3 8 0.002209
		4 8 0.002205
		1 9 0.002109
		2 9 0.002133
		3 9 0.002178
		4 9 0.002174
		1 10 0.002081
		2 10 0.002106
		3 10 0.002154
		4 10 0.002150
		Name: u235, dtype: float64
		1 1 0.002543
		2 1 0.002548
		3 1 0.002556
		4 1 0.002561
		1 2 0.002477
		2 2 0.002485
		3 2 0.002500
		4 2 0.002507
		1 3 0.002409
		2 3 0.002421
		3 3 0.002443
		4 3 0.002450
		1 4 0.002346
		2 4 0.002360
		3 4 0.002388
		4 4 0.002394
		1 5 0.002287
		2 5 0.002304
		3 5 0.002337
		4 5 0.002343
		1 6 0.002236
		2 6 0.002254
		3 6 0.002292
		4 6 0.002296
		1 7 0.002192
		2 7 0.002211
		3 7 0.002253
		4 7 0.002256
		1 8 0.002153
		2 8 0.002175
		3 8 0.002220
		4 8 0.002222
		1 9 0.002123
		2 9 0.002146
		3 9 0.002193
		4 9 0.002195
		1 10 0.002092
		2 10 0.002119
		3 10 0.002169
		4 10 0.002172
		Name: u235, dtype: float64

		%% Cell type:code id: tags:

		``` python
		# df_diff = (df_slice - df_benchmark) / df_benchmark
		df_diff = df_slice.subtract(df_benchmark)
		df_diff = df_diff.div(df_benchmark).replace(np.inf, 0)
		df_diff.reset_index(level=0, inplace=True)
		df_diff.index = [x for x in range(1, len(gfhr_nuc_list)+1)]
		df_diff.columns = ['Nuclide', '1', '2', '3', '4', '1.1', '2.1', '3.1', '4.1', '1.2',
		'2.2', '3.2', '4.2', '1.3', '2.3', '3.3', '4.3', '1.4', '2.4',
		'3.4', '4.4', '1.5', '2.5', '3.5', '4.5', '1.6', '2.6', '3.6',
		'4.6', '1.7', '2.7', '3.7', '4.7', '1.8', '2.8', '3.8', '4.8',
		'1.9', '2.9', '3.9', '4.9']

		df_diff
		```

		%% Output

		Nuclide 1 2 3 4 1.1 2.1 \
		1 u235 0.007314 0.006119 0.004265 0.002478 0.007972 0.007061
		2 u238 -0.001621 -0.001677 -0.001767 -0.001926 -0.001487 -0.001532
		3 pu239 -0.028862 -0.029419 -0.031867 -0.028674 -0.038408 -0.038228
		4 pu240 -0.028790 -0.027683 -0.024462 -0.024701 -0.023861 -0.023654
		5 pu241 -0.016502 -0.015191 -0.015134 -0.010259 -0.022800 -0.021889
		6 pu242 0.015287 0.016794 0.021025 0.024167 0.015152 0.016930
		7 am241 -0.148491 -0.149957 -0.153853 -0.157433 -0.149039 -0.149971
		8 sr90 -0.004477 -0.003038 -0.000704 0.001489 -0.004737 -0.003698
		9 xe135 0.055855 0.054677 0.049790 0.017035 -0.013263 -0.014925
		10 cs133 0.008133 0.008578 0.008995 0.009068 0.001650 0.002341
		11 cs134 -0.054744 -0.053347 -0.050502 -0.045891 -0.057346 -0.056020
		12 cs137 -0.007202 -0.005751 -0.003311 -0.000995 -0.007571 -0.006483
		13 nd143 0.001737 0.001710 0.000995 0.000516 -0.001503 -0.001545
		14 nd145 0.003226 0.004626 0.006951 0.008870 0.002877 0.003888
		15 nd148 -0.003714 -0.002271 0.000170 0.002451 -0.004323 -0.003214
		16 sm147 -0.033658 -0.034072 -0.034984 -0.036220 -0.031460 -0.032261
		17 sm149 0.317865 0.341225 0.381823 0.344456 0.252896 0.252597
		18 sm150 0.315456 0.316392 0.317931 0.319789 0.313459 0.314419
		19 sm151 0.066107 0.068004 0.066917 0.066577 0.052638 0.054132
		20 sm152 0.044915 0.046291 0.049105 0.049937 0.047957 0.049021
		21 eu153 -0.008004 -0.007182 -0.005640 -0.003668 -0.012284 -0.011438
		22 eu154 0.022697 0.025089 0.025391 0.032285 0.016945 0.019314
		23 eu155 -0.037738 -0.038002 -0.032867 -0.030582 -0.018227 -0.014919
		24 gd157 0.000000 0.000000 0.000000 0.000000 -0.076106 -0.076066

		3.1 4.1 1.2 ... 3.7 4.7 1.8 2.8 \
		1 0.005032 0.003047 0.007433 ... -0.004950 -0.007585 -0.006483 -0.005998
		2 -0.001649 -0.001800 -0.001426 ... -0.001898 -0.001842 -0.001822 -0.001810
		3 -0.037891 -0.035114 -0.044018 ... -0.046778 -0.062000 -0.052535 -0.052166
		4 -0.021288 -0.020445 -0.017824 ... 0.004653 -0.001365 0.009576 0.007441
		5 -0.020654 -0.017337 -0.027704 ... -0.028936 -0.032860 -0.030650 -0.030427
		6 0.020494 0.023966 0.017424 ... 0.042535 0.044473 0.046877 0.046037
		7 -0.152744 -0.156133 -0.154296 ... -0.177787 -0.182077 -0.184719 -0.182357
		8 -0.001338 0.000983 -0.003603 ... 0.009721 0.012023 0.010738 0.010448
		9 -0.012483 -0.014284 -0.021371 ... -0.024871 -0.026683 -0.030137 -0.025804
		10 0.004003 0.004976 -0.001403 ... 0.008991 0.011212 0.010310 0.010309
		11 -0.052750 -0.048514 -0.057585 ... -0.042661 -0.041786 -0.041393 -0.043577
		12 -0.004034 -0.001546 -0.006518 ... 0.007340 0.009609 0.008369 0.008062
		13 -0.001208 -0.001281 -0.005113 ... -0.005806 -0.004105 -0.006369 -0.006124
		14 0.006220 0.008365 0.004280 ... 0.018442 0.020386 0.019829 0.019521
		15 -0.000698 0.001784 -0.003443 ... 0.010636 0.012962 0.011760 0.011406
		16 -0.033210 -0.034796 -0.030598 ... -0.032940 -0.033236 -0.031630 -0.032029
		17 0.254169 0.205171 0.382545 ... 0.483115 0.404257 0.500269 0.493262
		18 0.316856 0.318552 0.317501 ... 0.355252 0.341616 0.368979 0.365307
		19 0.055650 0.054520 0.054963 ... 0.068204 0.056807 0.070196 0.069638
		20 0.051078 0.052464 0.051840 ... 0.076949 0.077862 0.083319 0.082721
		21 -0.008905 -0.006673 -0.012294 ... 0.001938 0.000548 0.003763 0.002545
		22 0.020848 0.027530 0.012234 ... 0.015263 0.019272 0.011351 0.012078
		23 -0.014960 -0.007728 -0.002184 ... 0.032012 0.035200 0.029888 0.032359
		24 -0.055949 0.000000 -0.097376 ... -0.052913 -0.018420 -0.060879 -0.056388

		3.8 4.8 1.9 2.9 3.9 4.9
		1 -0.006848 -0.009239 -0.005551 -0.005995 -0.006916 -0.009919
		2 -0.001983 -0.001876 -0.001937 -0.001969 -0.002138 -0.001987
		3 -0.046755 -0.062798 -0.041831 -0.041586 -0.037162 -0.056397
		4 0.006824 -0.001109 -0.003947 -0.004298 -0.003176 -0.007195
		5 -0.027830 -0.032439 -0.019933 -0.019464 -0.018145 -0.025156
		6 0.046054 0.047000 0.042718 0.043001 0.043844 0.046456
		7 -0.178393 -0.183375 -0.181885 -0.179026 -0.175311 -0.181214
		8 0.011467 0.013419 0.009631 0.010158 0.011266 0.013787
		9 -0.021191 -0.026274 0.019581 0.017793 0.015010 -0.012873
		10 0.011253 0.013060 0.011025 0.011721 0.012618 0.014601
		11 -0.039786 -0.039865 -0.037945 -0.039090 -0.035590 -0.036271
		12 0.009216 0.011036 0.007117 0.007669 0.008969 0.011399
		13 -0.004397 -0.002958 -0.004515 -0.003669 -0.002195 -0.000896
		14 0.020255 0.021805 0.018528 0.019058 0.019891 0.022088
		15 0.012541 0.014421 0.010582 0.011076 0.012342 0.014817
		16 -0.033001 -0.033383 -0.033185 -0.032667 -0.033463 -0.033307
		17 0.482014 0.407616 0.703074 0.688084 0.657956 0.512249
		18 0.360890 0.345144 0.369349 0.366606 0.361944 0.346489
		19 0.070038 0.057798 0.087581 0.086397 0.084695 0.067279
		20 0.080149 0.080290 0.076072 0.076861 0.075092 0.077897
		21 0.004471 0.002334 0.005445 0.005095 0.006814 0.004578
		22 0.016735 0.020068 0.023140 0.024411 0.028266 0.028553
		23 0.034115 0.035013 0.001135 0.003677 0.008348 0.020418
		24 -0.041459 -0.010051 0.190972 0.189690 0.188008 0.000000

		[24 rows x 41 columns]

		%% Cell type:code id: tags:

		``` python
		rcParams['figure.figsize'] = 12,3.5
		rcParams['font.size'] = 12.0
		rcParams['font.family'] = 'serif'
		rcParams['font.serif'] = 'DejaVu Serif'


		fontdict = {
		'fontsize': 12,
		'fontweight' : "bold",
		'font' : 'DejaVu Serif'
		}

		my_pal = {nuclide: "lightblue" if nuclide in fps else "salmon" for nuclide in df_diff["Nuclide"]}
		# my_pal = {species: "r" if species ==
		# "versicolor" else "b" for species in df["species"].unique()}


		sns.set_theme(style="ticks")

		df_flat = df_diff.melt(id_vars="Nuclide", var_name="Zone", value_name="Difference (%)")
		df_flat["Difference (%)"] = df_flat["Difference (%)"].multiply(100)
		vp = sns.violinplot(df_flat, x="Nuclide", y="Difference (%)", hue="Nuclide", order=gfhr_nuc_list, width=0.8, dodge=False, palette=my_pal, linewidth=0.5)

		new_labels = [v.get_text().capitalize() for v in vp.get_xticklabels()]
		vp.set_xticklabels(new_labels, rotation=45, font='DejaVu Serif')
		vp.set_xlabel('')
		vp.set_ylabel('Difference [%]', fontdict=fontdict)
		# vp.set_yticks([x/100 for x in range(-50, 60, 25)])
		vp.set_ylim(-20, 20)
		# vp.set_title(f'{key}', weight="bold", fontdict=fontdict)
		plt.grid(axis = 'y')
		plt.legend([],[], frameon=False)

		plt.savefig(f"zone-wise-dens-slice-vs-benchmark-res-5pct-violin.pdf", format='pdf', bbox_inches='tight', transparent=False)
		plt.show()
		```

		%% Output



		%% Cell type:code id: tags:

		``` python
		print(f"{df_diff.iloc[0].values[0]} {100min(df_diff.iloc[0].values[1:]):.2f}% {100max(df_diff.iloc[0].values[1:]):.2f}%")
		print(f"{df_diff.iloc[14].values[0]} {100min(df_diff.iloc[14].values[1:]):.2f}% {100max(df_diff.iloc[14].values[1:]):.2f}%")
		```

		%% Output

		u235 -0.99% 0.80%
		nd148 -0.43% 1.48%

_analysis/nuclide-densities/zone-wise-dens-comp-with-benchmark-violin.ipynb

+20 −19

Original line number	Diff line number	Diff line
		%% Cell type:markdown id: tags:

		# Visualize zone-wise nuclide comparison to benchmark via violin plot

		%% Cell type:code id: tags:

		``` python
		import seaborn as sns
		import matplotlib.pyplot as plt
		from matplotlib import rcParams
		import pandas as pd
		import json
		import os
		import re
		import glob
		import numpy as np
		from lib_utils import get_densities_from_f71

		acts = ["u235", "u238", "pu239", "pu240", "pu241", "pu242", "am241", "am242", \
		"am242m", "cm244", "cm245"]
		fps = ["kr85",
		"sr90",
		"ag110m",
		"i131",
		"xe135",
		"cs133", "cs134", "cs137",
		"nd143", "nd145", "nd148",
		"sm147", "sm149", "sm150", "sm151", "sm152",
		"eu153", "eu154", "eu155",
		"gd154", "gd155", "gd157"
		]

		orig_nuc_list = acts+fps

		stdcomp_fname = "/Users/f78/scale_dev/packages/InMemoryStdCompLib/InMemoryStdCompLibrary.json"
		# Download file from: https://code.ornl.gov/scale/code/scale-public/-/blob/main/packages/InMemoryStdCompLib/InMemoryStdCompLibrary.json
		stdcomp_fname = "InMemoryStdCompLibrary.json"
		```

		%% Cell type:markdown id: tags:

		### Benchmark composition

		%% Cell type:code id: tags:

		``` python
		def get_nuc_name_id_converter():
		"""
		Get converter between SCALE nuclide id and SCALE nuclide name

		Returns
		=======
		converter : dict
		key: SCALE id
		value: SCALE name
		"""

		assert os.path.isfile(stdcomp_fname)

		f = open(stdcomp_fname, "r")
		stdcmp = json.load(f)
		f.close()

		converter = {}
		for nuc in stdcmp['compositions']:
		if len(nuc['elements']) > 1:
		continue
		for el in nuc['elements']:
		converter[el['elementid']] = nuc['name']

		return converter
		```

		%% Cell type:code id: tags:

		``` python
		# Read Titan input
		titan_input_fname = "/Users/f78/neams/physor2024-pbr-performance/titan-gfhr/gfhr-tally-on.json"
		assert os.path.isfile(titan_input_fname)
		f = open(titan_input_fname, "r")
		titan_input = json.load(f)
		# Read KP compositions from file
		kp_material_input_fname = "../../_input_preparation/gfhr_benchmark_comps.json"
		assert os.path.isfile(kp_material_input_fname)
		f = open(kp_material_input_fname, "r")
		kp_material_input = json.load(f)
		f.close()

		# Get SCALE id/name converter
		nuc_id_to_name = get_nuc_name_id_converter()

		# get list of nuclide ids that are included in relevant nuclide list
		titan_nuc_list = []
		for mat_name, items in titan_input['materials'].items():
		kp_nuc_list = []
		for mat_name, items in kp_material_input['materials'].items():
		regExp = re.compile(r"FZ(?P<z>\d+)R(?P<r>\d+)")
		if not re.match(regExp, mat_name):
		continue
		assert len(items) == 2, f"Expected material vector to have length of 2, got {len(items)}"
		assert items[0] == 'CONC', f"Expected first value in vector to be 'CONC', got {items[0]}"
		titan_nuc_dens = items[1][1]
		for id, dens in titan_nuc_dens:
		kp_nuc_dens = items[1][1]
		for id, dens in kp_nuc_dens:
		assert str(id)[-1] == '0', f"Expected nuclide id to end in 0, got {id}"
		nuc_name = nuc_id_to_name[str(int(id/10))]
		nuc_name = f"{nuc_name.split('-')[0]}{nuc_name.split('-')[1]}"
		if nuc_name not in titan_nuc_list:
		titan_nuc_list.append(nuc_name)
		if nuc_name not in kp_nuc_list:
		kp_nuc_list.append(nuc_name)
		gfhr_nuc_list = []
		for nuc in acts+fps:
		if nuc in titan_nuc_list and nuc not in gfhr_nuc_list:
		if nuc in kp_nuc_list and nuc not in gfhr_nuc_list:
		gfhr_nuc_list.append(nuc)

		print(len(titan_nuc_list))
		print(len(kp_nuc_list))
		print(len(gfhr_nuc_list))

		print(gfhr_nuc_list)

		# Get list of SCALE materials, with nuclide ids correctly
		# converted from Titan to SCALE
		# converted from JSON to SCALE
		benchmark_comps = {}
		assert 'materials' in titan_input
		for mat_name, items in titan_input['materials'].items():
		assert 'materials' in kp_material_input
		for mat_name, items in kp_material_input['materials'].items():

		regExp = re.compile(r"FZ(?P<z>\d+)R(?P<r>\d+)")
		if not re.match(regExp, mat_name):
		continue
		z = int(re.match(regExp, mat_name).group('z'))
		r = int(re.match(regExp, mat_name).group('r'))

		assert len(items) == 2, f"Expected material vector to have length of 2, got {len(items)}"
		assert items[0] == 'CONC', f"Expected first value in vector to be 'CONC', got {items[0]}"

		# get all the nuclide densities
		titan_nuc_dens = items[1][1]
		kp_nuc_dens = items[1][1]
		scale_nuc_dens = {}
		for id, dens in titan_nuc_dens:
		for id, dens in kp_nuc_dens:
		assert str(id)[-1] == '0', f"Expected nuclide id to end in 0, got {id}"
		nuc_name = nuc_id_to_name[str(int(id/10))]
		nuc_name = f"{nuc_name.split('-')[0]}{nuc_name.split('-')[1]}"
		scale_nuc_dens[nuc_name] = dens

		# filter for nuclides we want to compare
		target_comps = {}
		for nuc in gfhr_nuc_list:
		if nuc in scale_nuc_dens:
		target_comps[nuc] = scale_nuc_dens[nuc]
		else:
		print(f"Setting zero density for {nuc} in r/z {r}/{z}")
		target_comps[nuc] = 0

		benchmark_comps[(r,z)] = target_comps

		```

		%% Output

		89
		24
		['u235', 'u238', 'pu239', 'pu240', 'pu241', 'pu242', 'am241', 'sr90', 'xe135', 'cs133', 'cs134', 'cs137', 'nd143', 'nd145', 'nd148', 'sm147', 'sm149', 'sm150', 'sm151', 'sm152', 'eu153', 'eu154', 'eu155', 'gd157']
		Setting zero density for gd157 in r/z 1/1
		Setting zero density for gd157 in r/z 2/1
		Setting zero density for gd157 in r/z 3/1
		Setting zero density for gd157 in r/z 4/1
		Setting zero density for gd157 in r/z 4/2
		Setting zero density for gd157 in r/z 4/3
		Setting zero density for gd157 in r/z 4/10

		%% Cell type:markdown id: tags:

		### SLICE final composition

		%% Cell type:code id: tags:

		``` python
		slice_zone_f71_list = glob.glob('../../tmp/eq_core_generation_r4z10_new/outer-it-5/arp_origen_output_it5/Z*.f71')
		slice_comps = {}

		for file in sorted(slice_zone_f71_list):
		df = get_densities_from_f71(file)
		dir = file.split("/")[4]
		subdir = file.split("/")[5]
		zone = file.split("/")[6]
		it_outer = re.search(r"it-(?P<it>\d+)",dir).group("it")
		it_inner = subdir[-1]
		Z = re.search(re.compile(r"Z(?P<Z>\d+)_R(?P<R>\d+)"),zone).group("Z")
		R = re.search(re.compile(r"Z(?P<Z>\d+)_R(?P<R>\d+)"),zone).group("R")
		Z = int(Z)
		R = int(R)
		target_comps = {}
		for nuc in gfhr_nuc_list:
		target_comps[nuc] = df.loc[nuc].values[0]
		slice_comps[(R,Z+1)] = target_comps
		```

		%% Cell type:markdown id: tags:

		### Determine relative difference

		%% Cell type:code id: tags:

		``` python
		df_slice = pd.DataFrame.from_dict(slice_comps)
		print(df_slice.loc["u235"])

		df_benchmark = pd.DataFrame.from_dict(benchmark_comps)
		print(df_benchmark.loc["u235"])
		```

		%% Output

		1 1 0.002494
		2 1 0.002496
		3 1 0.002500
		4 1 0.002501
		1 2 0.002428
		2 2 0.002434
		3 2 0.002445
		4 2 0.002446
		1 3 0.002357
		2 3 0.002367
		3 3 0.002385
		4 3 0.002386
		1 4 0.002287
		2 4 0.002301
		3 4 0.002325
		4 4 0.002326
		1 5 0.002222
		2 5 0.002238
		3 5 0.002269
		4 5 0.002268
		1 6 0.002163
		2 6 0.002182
		3 6 0.002217
		4 6 0.002215
		1 7 0.002111
		2 7 0.002133
		3 7 0.002172
		4 7 0.002169
		1 8 0.002068
		2 8 0.002091
		3 8 0.002134
		4 8 0.002131
		1 9 0.002033
		2 9 0.002057
		3 9 0.002103
		4 9 0.002099
		1 10 0.002004
		2 10 0.002030
		3 10 0.002078
		4 10 0.002075
		Name: u235, dtype: float64
		1 1 0.002543
		2 1 0.002548
		3 1 0.002556
		4 1 0.002561
		1 2 0.002477
		2 2 0.002485
		3 2 0.002500
		4 2 0.002507
		1 3 0.002409
		2 3 0.002421
		3 3 0.002443
		4 3 0.002450
		1 4 0.002346
		2 4 0.002360
		3 4 0.002388
		4 4 0.002394
		1 5 0.002287
		2 5 0.002304
		3 5 0.002337
		4 5 0.002343
		1 6 0.002236
		2 6 0.002254
		3 6 0.002292
		4 6 0.002296
		1 7 0.002192
		2 7 0.002211
		3 7 0.002253
		4 7 0.002256
		1 8 0.002153
		2 8 0.002175
		3 8 0.002220
		4 8 0.002222
		1 9 0.002123
		2 9 0.002146
		3 9 0.002193
		4 9 0.002195
		1 10 0.002092
		2 10 0.002119
		3 10 0.002169
		4 10 0.002172
		Name: u235, dtype: float64

		%% Cell type:code id: tags:

		``` python
		# df_diff = (df_slice - df_benchmark) / df_benchmark
		df_diff = df_slice.subtract(df_benchmark)
		df_diff = df_diff.div(df_benchmark).replace(np.inf, 0)
		df_diff.reset_index(level=0, inplace=True)
		df_diff.index = [x for x in range(1, len(gfhr_nuc_list)+1)]
		df_diff.columns = ['Nuclide', '1', '2', '3', '4', '1.1', '2.1', '3.1', '4.1', '1.2',
		'2.2', '3.2', '4.2', '1.3', '2.3', '3.3', '4.3', '1.4', '2.4',
		'3.4', '4.4', '1.5', '2.5', '3.5', '4.5', '1.6', '2.6', '3.6',
		'4.6', '1.7', '2.7', '3.7', '4.7', '1.8', '2.8', '3.8', '4.8',
		'1.9', '2.9', '3.9', '4.9']

		df_diff
		```

		%% Output

		Nuclide 1 2 3 4 1.1 2.1 \
		1 u235 -0.019066 -0.020194 -0.021926 -0.023652 -0.019758 -0.020528
		2 u238 -0.003995 -0.004047 -0.004132 -0.004278 -0.003964 -0.003999
		3 pu239 -0.022193 -0.022735 -0.025136 -0.022438 -0.030313 -0.030257
		4 pu240 -0.006297 -0.005176 -0.001794 -0.001483 -0.002200 -0.002043
		5 pu241 0.024870 0.026310 0.026382 0.031002 0.017698 0.018866
		6 pu242 0.120405 0.122094 0.126847 0.130427 0.118987 0.121042
		7 am241 -0.088078 -0.089549 -0.093497 -0.097211 -0.090711 -0.091404
		8 sr90 0.030134 0.031629 0.034052 0.036324 0.029681 0.030773
		9 xe135 0.040191 0.039151 0.034346 0.001903 -0.028075 -0.029589
		10 cs133 0.044609 0.045075 0.045512 0.045623 0.037863 0.038563
		11 cs134 0.026584 0.028095 0.031187 0.036169 0.023429 0.024851
		12 cs137 0.033437 0.034952 0.037500 0.039912 0.032909 0.034056
		13 nd143 0.032469 0.032456 0.031747 0.031251 0.028555 0.028534
		14 nd145 0.038128 0.039578 0.041985 0.043993 0.037536 0.038602
		15 nd148 0.037284 0.038791 0.041342 0.043722 0.036498 0.037667
		16 sm147 0.032948 0.032526 0.031596 0.030355 0.034468 0.033679
		17 sm149 0.310977 0.333116 0.371074 0.330242 0.262224 0.261302
		18 sm150 0.376884 0.377878 0.379535 0.381472 0.374847 0.375812
		19 sm151 0.080285 0.082194 0.080936 0.079726 0.068359 0.069866
		20 sm152 0.078038 0.079478 0.082436 0.083608 0.080523 0.081675
		21 eu153 0.050344 0.051216 0.052852 0.054921 0.045777 0.046669
		22 eu154 0.100816 0.103417 0.103737 0.110608 0.094345 0.096970
		23 eu155 0.025576 0.025265 0.030917 0.034586 0.045754 0.049140
		24 gd157 0.000000 0.000000 0.000000 0.000000 -0.028032 -0.028046

		3.1 4.1 1.2 ... 3.7 4.7 1.8 2.8 \
		1 -0.022277 -0.024181 -0.021737 ... -0.038580 -0.041232 -0.042411 -0.041370
		2 -0.004099 -0.004211 -0.004015 ... -0.004862 -0.004623 -0.004972 -0.004915
		3 -0.030072 -0.029552 -0.034864 ... -0.038104 -0.059200 -0.044434 -0.044035
		4 0.000588 0.002860 0.003322 ... 0.025148 0.019891 0.029995 0.027616
		5 0.020229 0.022468 0.011807 ... 0.008162 0.001242 0.004940 0.005596
		6 0.125215 0.129348 0.120134 ... 0.143497 0.146535 0.146004 0.145675
		7 -0.093732 -0.097128 -0.098511 ... -0.126599 -0.132206 -0.136270 -0.132918
		8 0.033244 0.035646 0.030639 ... 0.043591 0.045819 0.044205 0.044012
		9 -0.027183 -0.029833 -0.036183 ... -0.040716 -0.046053 -0.046547 -0.042050
		10 0.040240 0.041350 0.034766 ... 0.044318 0.046801 0.045095 0.045192
		11 0.028407 0.032933 0.022855 ... 0.037439 0.038125 0.038606 0.036169
		12 0.036628 0.039210 0.033855 ... 0.047756 0.049784 0.048558 0.048308
		13 0.028918 0.028821 0.024492 ... 0.022628 0.024053 0.021725 0.021990
		14 0.041046 0.043323 0.038718 ... 0.052367 0.054354 0.053260 0.053090
		15 0.040313 0.042896 0.037264 ... 0.051396 0.053528 0.052284 0.051988
		16 0.032804 0.031361 0.034488 ... 0.029090 0.030047 0.029184 0.029069
		17 0.260729 0.205611 0.384855 ... 0.471447 0.388043 0.486372 0.479651
		18 0.378320 0.379985 0.379153 ... 0.417078 0.402032 0.430781 0.427017
		19 0.071023 0.067605 0.071452 ... 0.084462 0.066969 0.086463 0.085969
		20 0.083936 0.086228 0.083948 ... 0.108444 0.111311 0.114314 0.113790
		21 0.049321 0.051536 0.045442 ... 0.058444 0.056972 0.059412 0.058417
		22 0.098639 0.104270 0.088868 ... 0.089973 0.090284 0.084633 0.085708
		23 0.049381 0.059221 0.062988 ... 0.099849 0.099577 0.097291 0.099913
		24 -0.008451 0.000000 -0.049468 ... -0.005795 0.022848 -0.015070 -0.009636

		3.8 4.8 1.9 2.9 3.9 4.9
		1 -0.041123 -0.043551 -0.042194 -0.042011 -0.041760 -0.044777
		2 -0.005003 -0.004697 -0.005141 -0.005124 -0.005203 -0.004839
		3 -0.038403 -0.060322 -0.034015 -0.033708 -0.029025 -0.054175
		4 0.027238 0.019760 0.016061 0.015517 0.016929 0.013231
		5 0.008941 0.001227 0.015702 0.016623 0.018689 0.008404
		6 0.146734 0.148731 0.140874 0.141800 0.143817 0.147678
		7 -0.127060 -0.133667 -0.132760 -0.128895 -0.123316 -0.131123
		8 0.045256 0.047099 0.042920 0.043579 0.044929 0.047345
		9 -0.037304 -0.046018 0.002098 0.000601 -0.001878 -0.033174
		10 0.046290 0.048313 0.045413 0.046256 0.047343 0.049517
		11 0.040235 0.039872 0.042017 0.040703 0.044459 0.043432
		12 0.049617 0.051134 0.047168 0.047814 0.049284 0.051399
		13 0.023758 0.024844 0.023213 0.024117 0.025661 0.026539
		14 0.054079 0.055638 0.051755 0.052461 0.053570 0.055782
		15 0.053289 0.054915 0.050969 0.051561 0.053008 0.055222
		16 0.028614 0.029593 0.027156 0.028030 0.027788 0.029435
		17 0.468561 0.389677 0.685890 0.671183 0.641225 0.491214
		18 0.422647 0.405352 0.430828 0.428053 0.423464 0.406443
		19 0.086113 0.067605 0.103890 0.102753 0.100786 0.076867
		20 0.111540 0.113595 0.106678 0.107583 0.106177 0.110947
		21 0.060759 0.058545 0.060828 0.060773 0.062923 0.060646
		22 0.091174 0.090728 0.096875 0.098560 0.103225 0.099479
		23 0.101949 0.098870 0.066434 0.069165 0.074348 0.082987
		24 0.004891 0.030074 0.247070 0.246650 0.243508 0.000000

		[24 rows x 41 columns]

		%% Cell type:code id: tags:

		``` python
		rcParams['figure.figsize'] = 12,3.5
		rcParams['font.size'] = 12.0
		rcParams['font.family'] = 'serif'
		rcParams['font.serif'] = 'DejaVu Serif'


		fontdict = {
		'fontsize': 12,
		'fontweight' : "bold",
		'font' : 'DejaVu Serif'
		}

		my_pal = {nuclide: "lightblue" if nuclide in fps else "salmon" for nuclide in df_diff["Nuclide"]}
		# my_pal = {species: "r" if species ==
		# "versicolor" else "b" for species in df["species"].unique()}


		sns.set_theme(style="ticks")

		df_flat = df_diff.melt(id_vars="Nuclide", var_name="Zone", value_name="Difference (%)")
		df_flat["Difference (%)"] = df_flat["Difference (%)"].multiply(100)
		vp = sns.violinplot(df_flat, x="Nuclide", y="Difference (%)", hue="Nuclide", order=gfhr_nuc_list, width=0.8, dodge=False, palette=my_pal, linewidth=0.5)

		new_labels = [v.get_text().capitalize() for v in vp.get_xticklabels()]
		vp.set_xticklabels(new_labels, rotation=45, font='DejaVu Serif')
		vp.set_xlabel('')
		vp.set_ylabel('Difference [%]', fontdict=fontdict)
		# vp.set_yticks([x/100 for x in range(-50, 60, 25)])
		vp.set_ylim(-20, 20)
		# vp.set_title(f'{key}', weight="bold", fontdict=fontdict)
		plt.grid(axis = 'y')
		plt.legend([],[], frameon=False)

		plt.savefig(f"zone-wise-dens-slice-vs-benchmark-violin.pdf", format='pdf', bbox_inches='tight', transparent=False)
		plt.show()
		```

		%% Output



		%% Cell type:code id: tags:

		``` python
		print(f"{df_diff.iloc[0].values[0]} {100min(df_diff.iloc[0].values[1:]):.2f}% {100max(df_diff.iloc[0].values[1:]):.2f}%")
		print(f"{df_diff.iloc[14].values[0]} {100min(df_diff.iloc[14].values[1:]):.2f}% {100max(df_diff.iloc[14].values[1:]):.2f}%")
		```

		%% Output

		u235 -4.48% -1.91%
		nd148 3.65% 5.52%

_analysis/power/core_profiles.ipynb

+1 −1

Original line number	Diff line number	Diff line
		%% Cell type:markdown id: tags:

		# Full core power and flux profile

		%% Cell type:code id: tags:

		``` python
		import os
		import sys
		import json
		import numpy as np
		sys.path.insert(0, os.path.abspath('/Users/f78/slice-method-for-pbr-eq-core-generation/scripts'))
		sys.path.insert(0, os.path.abspath('~/slice-method-for-pbr-eq-core-generation/scripts'))
		import matplotlib
		import matplotlib.ticker as mtick
		import matplotlib.pyplot as plt
		from collections import OrderedDict

		%load_ext autoreload
		%autoreload 2

		%matplotlib inline
		%load_ext wurlitzer
		```

		%% Cell type:code id: tags:

		``` python
		power_files = {}
		power_files["fresh"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-0.power_normalized.json'
		power_files["outer it. 1"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-1.power_normalized.json'
		power_files["outer it. 2"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-2.power_normalized.json'
		power_files["outer it. 3"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-3.power_normalized.json'
		power_files["outer it. 4"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-4.power_normalized.json'
		power_files["outer it. 5"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-5.power_normalized.json'
		power_files["outer it. 6"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-6.power_normalized.json'

		flux_files = {}
		flux_files["fresh"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-0.flux_normalized.json'
		flux_files["outer it. 1"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-1.flux_normalized.json'
		flux_files["outer it. 2"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-2.flux_normalized.json'
		flux_files["outer it. 3"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-3.flux_normalized.json'
		flux_files["outer it. 4"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-4.flux_normalized.json'
		flux_files["outer it. 5"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-5.flux_normalized.json'
		flux_files["outer it. 6"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-6.flux_normalized.json'

		for file in power_files.values():
		assert os.path.isfile(file), file
		for file in flux_files.values():
		assert os.path.isfile(file), file

		radial_bounds = [0, 30, 60, 90, 120]
		axial_bounds = [x for x in range(0, 310+31, 31)]
		```

		%% Cell type:code id: tags:

		``` python
		power_axial = OrderedDict()
		flux_axial = OrderedDict()
		power_radial = OrderedDict()
		flux_radial = OrderedDict()
		for key, file in power_files.items():
		with open(file, "r", encoding="utf-8") as f:
		data = json.load(f)
		axial = []
		for p in (data["axial_profile"]):
		axial.append(p)
		axial.append(axial[-1])
		power_axial[key] = axial
		radial = []
		for pow in data["radial_profile"]:
		rad_power = []
		for p in pow:
		rad_power.append(p)
		rad_power.append(rad_power[-1])
		radial.append(rad_power)
		power_radial[key] = radial
		for key, file in flux_files.items():
		with open(file, "r", encoding="utf-8") as f:
		data = json.load(f)
		axial = []
		for p in (data["axial_profile"]):
		axial.append(p)
		axial.append(axial[-1])
		flux_axial[key] = axial
		radial = []
		for pow in data["radial_profile"]:
		rad_power = []
		for p in pow:
		rad_power.append(p)
		rad_power.append(rad_power[-1])
		radial.append(rad_power)
		flux_radial[key] = radial
		```

		%% Cell type:code id: tags:

		``` python
		matplotlib.rcParams['font.size'] = 14.0
		matplotlib.rcParams['font.family'] = 'serif'

		fig, ax = plt.subplots(1, 2, sharey=True, figsize=(10,6))


		ax[0].set_xlim(0, 1.5)
		ax[1].set_xlim(-0.4, 0.4)
		ax[0].set_ylim(min(axial_bounds)-10, max(axial_bounds)+10)
		ax[1].set_ylim(min(axial_bounds)-10, max(axial_bounds)+10)
		ax[0].set_ylabel('Axial height [cm]', fontweight="bold")
		ax[0].set_xlabel('Normalized power', fontweight="bold")
		ax[1].set_xlabel('Difference', fontweight="bold")
		ax[1].xaxis.set_major_formatter(mtick.PercentFormatter(1.0))

		ax[0].xaxis.set_major_locator(mtick.MultipleLocator(0.5))
		ax[0].xaxis.set_minor_locator(mtick.MultipleLocator(0.1))
		ax[1].xaxis.set_major_locator(mtick.MultipleLocator(0.2))
		ax[1].xaxis.set_minor_locator(mtick.MultipleLocator(0.05))
		ax[0].yaxis.set_major_locator(mtick.MultipleLocator(50))
		ax[0].yaxis.set_minor_locator(mtick.MultipleLocator(10))

		colors = ["black", "blue", "red", "gold", "purple", "green", "gray", "cyan", "pink"]

		cnt = 0
		for key, p in power_axial.items():
		ax[0].step(p, [x for x in axial_bounds], where='pre',\
		linestyle="solid", linewidth=1, color=colors[cnt], label=f"{key}")
		cnt += 1

		diffs = OrderedDict()
		keys = list(power_axial.keys())
		for i, key in enumerate(keys):
		if i == 0:
		continue
		diffs[key] = np.array(power_axial[keys[i]]) - np.array(power_axial[keys[i-1]])

		cnt = 1
		for key, p in diffs.items():
		rms = np.sum([v**2 for v in p])/len(p)
		ax[1].step(p, [x for x in axial_bounds], where='pre',\
		linestyle="-.", linewidth=1, color=colors[cnt], \
		# label=f"{key}, rms={100*rms:.4f}%")
		label=f"rms={100*rms:.3f}%")
		cnt += 1

		for x in ax:
		x.legend(ncol=1, shadow=True, fancybox=True, handlelength=1.0, handletextpad=0.3, loc='best', fontsize=10)

		plt.savefig(f"core1_axial_power_profile.pdf", format='pdf', bbox_inches='tight', transparent=False)
		plt.show()
		```

		%% Output



		%% Cell type:code id: tags:

		``` python
		matplotlib.rcParams['font.size'] = 14.0
		matplotlib.rcParams['font.family'] = 'serif'

		fig, ax = plt.subplots(1, 2, sharey=True, figsize=(10,6))


		ax[0].set_xlim(0, 1.6)
		ax[1].set_xlim(-0.4, 0.4)
		ax[0].set_ylim(min(axial_bounds)-10, max(axial_bounds)+10)
		ax[1].set_ylim(min(axial_bounds)-10, max(axial_bounds)+10)
		ax[0].set_ylabel('Axial height [cm]', fontweight="bold")
		ax[0].set_xlabel('Normalized flux', fontweight="bold")
		ax[1].set_xlabel('Difference', fontweight="bold")
		ax[1].xaxis.set_major_formatter(mtick.PercentFormatter(1.0))

		ax[0].xaxis.set_major_locator(mtick.MultipleLocator(0.5))
		ax[0].xaxis.set_minor_locator(mtick.MultipleLocator(0.1))
		ax[1].xaxis.set_major_locator(mtick.MultipleLocator(0.2))
		ax[1].xaxis.set_minor_locator(mtick.MultipleLocator(0.05))
		ax[0].yaxis.set_major_locator(mtick.MultipleLocator(50))
		ax[0].yaxis.set_minor_locator(mtick.MultipleLocator(10))

		colors = ["black", "blue", "red", "gold", "purple", "green", "gray", "cyan", "pink"]

		cnt = 0
		for key, p in flux_axial.items():
		ax[0].step(p, [x for x in axial_bounds], where='pre',\
		linestyle="solid", linewidth=1, color=colors[cnt], label=f"{key}")
		cnt += 1

		diffs = OrderedDict()
		keys = list(flux_axial.keys())
		for i, key in enumerate(keys):
		if i == 0:
		continue
		diffs[key] = np.array(flux_axial[key]) - np.array(flux_axial[keys[i-1]])


		cnt = 1
		for key, p in diffs.items():
		rms = np.sum([v**2 for v in p])/len(p)
		ax[1].step(p, [x for x in axial_bounds], where='pre',\
		linestyle="-.", linewidth=1, color=colors[cnt], \
		# label=f"{key}, rms={100*rms:.4f}%")
		label=f"rms={100*rms:.4f}%")
		cnt += 1

		for x in ax:
		x.legend(ncol=1, shadow=True, fancybox=True, handlelength=1.0, handletextpad=0.3, loc='best', fontsize=10)

		plt.savefig(f"core1_axial_flux_profile.pdf", format='pdf', bbox_inches='tight', transparent=False)
		plt.show()
		```

		%% Output



		%% Cell type:code id: tags:

		``` python
		```

_analysis/power/core_profiles_comp_with_benchmark.ipynb

+1 −1

Original line number	Diff line number	Diff line
		%% Cell type:markdown id: tags:

		# Full core power and flux profile - comparison with KP benchmark model

		%% Cell type:code id: tags:

		``` python
		import os
		import sys
		import json
		import numpy as np
		sys.path.insert(0, os.path.abspath('/Users/f78/slice-method-for-pbr-eq-core-generation/scripts'))
		sys.path.insert(0, os.path.abspath('~/slice-method-for-pbr-eq-core-generation/scripts'))
		import matplotlib
		import matplotlib.ticker as mtick
		import matplotlib.pyplot as plt
		from collections import OrderedDict

		%load_ext autoreload
		%autoreload 2

		%matplotlib inline
		%load_ext wurlitzer
		```

		%% Cell type:code id: tags:

		``` python
		power_files = {}
		power_files["Benchmark"] = '../../full_core_kp_inventory/core-dodec-array.t6-depl-shift.v7.1-252.power_normalized.json'
		power_files["SLICE"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-5.power_normalized.json'
		power_files["SLICE (95% $T_{pass}$)"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-6.res-5pct.power_normalized.json'
		power_files["SLICE (90% $T_{pass}$)"] = '../../tmp/eq_core_generation_r4z10_new/core-array.csas6-shift.v7.1-252.oi-6.res-10pct.power_normalized.json'


		for file in power_files.values():
		assert os.path.isfile(file), file

		core_height = 310
		radial_bounds = [0, 30, 60, 90, 120]
		axial_bounds = [x for x in range(0, 310+31, 31)]
		```

		%% Cell type:code id: tags:

		``` python
		power_axial = OrderedDict()
		power_radial = OrderedDict()
		for key, file in power_files.items():
		with open(file, "r", encoding="utf-8") as f:
		data = json.load(f)
		axial = []
		for p in (data["axial_profile"]):
		axial.append(p)
		axial.append(axial[-1])
		power_axial[key] = axial
		radial = []
		for pow in data["radial_profile"]:
		rad_power = []
		for p in pow:
		rad_power.append(p)
		rad_power.append(rad_power[-1])
		radial.append(rad_power)
		power_radial[key] = radial
		```

		%% Cell type:code id: tags:

		``` python
		matplotlib.rcParams['font.size'] = 14.0
		matplotlib.rcParams['font.family'] = 'serif'

		fig, ax = plt.subplots(1, 2, sharey=True, figsize=(10,6))


		ax[0].set_xlim(0, 1.5)
		ax[1].set_xlim(-0.1, 0.1)
		ax[0].set_ylim(-10, core_height+10)
		ax[1].set_ylim(-10, core_height+10)
		ax[0].set_ylabel('Axial height [cm]', fontweight="bold")
		ax[0].set_xlabel('Normalized power', fontweight="bold")
		ax[1].set_xlabel('Difference', fontweight="bold")
		ax[1].xaxis.set_major_formatter(mtick.PercentFormatter(1.0))

		ax[0].xaxis.set_major_locator(mtick.MultipleLocator(0.5))
		ax[0].xaxis.set_minor_locator(mtick.MultipleLocator(0.1))
		# ax[1].xaxis.set_major_locator(mtick.MultipleLocator(0.02))
		ax[1].xaxis.set_minor_locator(mtick.MultipleLocator(0.01))
		ax[0].yaxis.set_major_locator(mtick.MultipleLocator(50))
		ax[0].yaxis.set_minor_locator(mtick.MultipleLocator(10))

		colors = ["black", "blue", "red", "gold", "green", "gray", "cyan", "pink"]

		cnt = 0
		for key, p in power_axial.items():
		ax[0].step(p, [x for x in axial_bounds], where='pre',\
		linestyle="solid", linewidth=1, color=colors[cnt], label=f"{key}")
		cnt += 1

		diffs = OrderedDict()
		keys = list(power_axial.keys())
		for i, key in enumerate(keys):
		if i == 0:
		continue
		diffs[key] = np.array(power_axial[keys[i]]) - np.array(power_axial[keys[0]])

		cnt = 1
		for key, p in diffs.items():
		rms = np.sum([v**2 for v in p])/len(p)
		print(f"Min: {min(p)}")
		print(f"Max: {max(p)}")
		ax[1].step(p, [x for x in axial_bounds], where='pre',\
		linestyle="-.", linewidth=1, color=colors[cnt], \
		# label=f"{key}, rms={100*rms:.4f}%")
		label=f"rms={100*rms:.2f}%")
		cnt += 1

		for x in ax:
		x.legend(ncol=1, shadow=True, fancybox=True, handlelength=1.0, handletextpad=0.3, loc='best', fontsize=10)

		plt.savefig(f"axial_power_profile_slice_vs_benchmark.pdf", format='pdf', bbox_inches='tight', transparent=False)
		plt.show()
		```

		%% Output

		Min: -0.06805000560370256
		Max: 0.08895411564163425
		Min: -0.049521229725611904
		Max: 0.07429574097809777
		Min: -0.03856242860987191
		Max: 0.05007867676043909