diff --git a/benchmarks/pu-met-fast-001/geometry.xml b/benchmarks/pu-met-fast-001/geometry.xml
deleted file mode 100644
index 7d1e489729..0000000000
--- a/benchmarks/pu-met-fast-001/geometry.xml
+++ /dev/null
@@ -1,25 +0,0 @@
-
-
-
-
-
- |
- 0
- 1
- -1
- |
-
-
- sphere
- 0. 0. 0. 6.3849
- vacuum
-
-
-
diff --git a/benchmarks/pu-met-fast-001/materials.xml b/benchmarks/pu-met-fast-001/materials.xml
deleted file mode 100644
index 4d5fce4cac..0000000000
--- a/benchmarks/pu-met-fast-001/materials.xml
+++ /dev/null
@@ -1,21 +0,0 @@
-
-
-
-
-
-
-
-
-
-
-
-
-
-
diff --git a/benchmarks/pu-met-fast-001/settings.xml b/benchmarks/pu-met-fast-001/settings.xml
deleted file mode 100644
index aa5b251cac..0000000000
--- a/benchmarks/pu-met-fast-001/settings.xml
+++ /dev/null
@@ -1,28 +0,0 @@
-
-
-
-
-
-
-
-
- 3000
- 20
- 10000
-
-
-
-
-
- box
- -1 -1 -1 1 1 1
-
-
-
diff --git a/benchmarks/pu-met-fast-001/tallies.xml b/benchmarks/pu-met-fast-001/tallies.xml
deleted file mode 100644
index 5f7f3264d3..0000000000
--- a/benchmarks/pu-met-fast-001/tallies.xml
+++ /dev/null
@@ -1,13 +0,0 @@
-
-
-
-
-
- 1 |
-
-
- total scatter absorption fission nu-fission
-
-
-
-
diff --git a/benchmarks/pu-met-fast-002/geometry.xml b/benchmarks/pu-met-fast-002/geometry.xml
deleted file mode 100644
index 6c2ceb58a9..0000000000
--- a/benchmarks/pu-met-fast-002/geometry.xml
+++ /dev/null
@@ -1,25 +0,0 @@
-
-
-
-
-
-
- 0
- 1
- -1
- |
-
-
- sphere
- 0. 0. 0. 6.6595
- vacuum
-
-
-
diff --git a/benchmarks/pu-met-fast-002/materials.xml b/benchmarks/pu-met-fast-002/materials.xml
deleted file mode 100644
index 779ec1eb44..0000000000
--- a/benchmarks/pu-met-fast-002/materials.xml
+++ /dev/null
@@ -1,22 +0,0 @@
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
diff --git a/benchmarks/pu-met-fast-002/settings.xml b/benchmarks/pu-met-fast-002/settings.xml
deleted file mode 100644
index 74e4e2b8e1..0000000000
--- a/benchmarks/pu-met-fast-002/settings.xml
+++ /dev/null
@@ -1,28 +0,0 @@
-
-
-
-
-
-
-
-
- 3000
- 20
- 10000
-
-
-
-
-
- box
- -1 -1 -1 1 1 1
-
-
-
diff --git a/benchmarks/pu-met-fast-005/geometry.xml b/benchmarks/pu-met-fast-005/geometry.xml
deleted file mode 100644
index 2b11068a42..0000000000
--- a/benchmarks/pu-met-fast-005/geometry.xml
+++ /dev/null
@@ -1,34 +0,0 @@
-
-
-
-
-
-
- 0
- 1
- -1
- |
-
- 0
- 2
- 1 -2
- |
-
-
- sphere
- 0. 0. 0. 5.0419
-
-
- sphere
- 0. 0. 0. 9.7409
- vacuum
-
-
-
diff --git a/benchmarks/pu-met-fast-005/materials.xml b/benchmarks/pu-met-fast-005/materials.xml
deleted file mode 100644
index 977733e54f..0000000000
--- a/benchmarks/pu-met-fast-005/materials.xml
+++ /dev/null
@@ -1,28 +0,0 @@
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
diff --git a/benchmarks/pu-met-fast-005/settings.xml b/benchmarks/pu-met-fast-005/settings.xml
deleted file mode 100644
index 831fbafbe7..0000000000
--- a/benchmarks/pu-met-fast-005/settings.xml
+++ /dev/null
@@ -1,28 +0,0 @@
-
-
-
-
-
-
-
-
- 3000
- 20
- 10000
-
-
-
-
-
- box
- -1 -1 -1 1 1 1
-
-
-
diff --git a/cross_sections_serpent.xml b/cross_sections_serpent.xml
new file mode 100644
index 0000000000..d61e95a161
--- /dev/null
+++ b/cross_sections_serpent.xml
@@ -0,0 +1,2705 @@
+
+
+
+ /opt/serpent/xsdata/endfb7/acedata
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/examples/basic/geometry.xml b/examples/basic/geometry.xml
index eef6f89060..21e85356ec 100644
--- a/examples/basic/geometry.xml
+++ b/examples/basic/geometry.xml
@@ -2,16 +2,14 @@
- |
- |
- |
- |
+ |
+ |
+ |
+ |
-
-
-
-
-
+
+
+
diff --git a/examples/basic/materials.xml b/examples/basic/materials.xml
index 450c988f19..88e20c69ef 100644
--- a/examples/basic/materials.xml
+++ b/examples/basic/materials.xml
@@ -2,12 +2,12 @@
-
+
-
+
diff --git a/examples/basic/settings.xml b/examples/basic/settings.xml
index b59a0de03e..a9a83dfe54 100644
--- a/examples/basic/settings.xml
+++ b/examples/basic/settings.xml
@@ -2,7 +2,9 @@
-
+
+ /home/paulromano/openmc/cross_sections_serpent.xml
+
@@ -11,9 +13,6 @@
10000
-
-
-
box
diff --git a/examples/lattice/geometry.xml b/examples/lattice/geometry.xml
index e5b2934595..8848e6810e 100644
--- a/examples/lattice/geometry.xml
+++ b/examples/lattice/geometry.xml
@@ -1,23 +1,23 @@
-
+ |
0
111
1 -2 3 -4
|
-
+ |
3
40
-5
|
-
+ |
3
41
5
|
-
+
rectangular
20 20
-10.0 -10.0
@@ -46,11 +46,11 @@
-
-
-
-
-
+
+
+
+
+
z-cylinder
0.0 0.0 0.4
diff --git a/examples/lattice/materials.xml b/examples/lattice/materials.xml
index 8070d64999..88e20c69ef 100644
--- a/examples/lattice/materials.xml
+++ b/examples/lattice/materials.xml
@@ -2,16 +2,16 @@
-
+
-
+
-
+
diff --git a/examples/lattice/settings.xml b/examples/lattice/settings.xml
index beb7193fad..adbec19c7d 100644
--- a/examples/lattice/settings.xml
+++ b/examples/lattice/settings.xml
@@ -2,7 +2,9 @@
-
+
+ /home/paulromano/openmc/cross_sections_serpent.xml
+
@@ -11,9 +13,6 @@
10000
-
-
-
box
diff --git a/examples/reflective/geometry.xml b/examples/reflective/geometry.xml
index 6f8c498fac..e080e49270 100644
--- a/examples/reflective/geometry.xml
+++ b/examples/reflective/geometry.xml
@@ -2,18 +2,18 @@
-
+ |
0
1
1 -2 3 -4 5 -6
|
-
-
-
-
-
-
+
+
+
+
+
+
|
diff --git a/examples/reflective/materials.xml b/examples/reflective/materials.xml
index 616065159a..eb55455f7a 100644
--- a/examples/reflective/materials.xml
+++ b/examples/reflective/materials.xml
@@ -2,7 +2,7 @@
-
+
diff --git a/examples/reflective/settings.xml b/examples/reflective/settings.xml
index e69a9ab675..3dd2f6e7d3 100644
--- a/examples/reflective/settings.xml
+++ b/examples/reflective/settings.xml
@@ -2,7 +2,9 @@
-
+
+ /home/paulromano/openmc/cross_sections_serpent.xml
+
@@ -11,9 +13,6 @@
10000
-
-
-
box
diff --git a/src/DEPENDENCIES b/src/DEPENDENCIES
index cb1119e040..432c5e1bc5 100644
--- a/src/DEPENDENCIES
+++ b/src/DEPENDENCIES
@@ -13,6 +13,8 @@ cross_section.o: material_header.o
cross_section.o: output.o
cross_section.o: string.o
+cross_section_header.o: constants.o
+
datatypes.o: datatypes_header.o
doppler.o: constants.o
@@ -35,6 +37,7 @@ fileio.o: string.o
fission.o: constants.o
fission.o: cross_section_header.o
fission.o: error.o
+fission.o: global.o
fission.o: interpolation.o
fission.o: search.o
@@ -46,6 +49,7 @@ geometry.o: global.o
geometry.o: output.o
geometry.o: particle_header.o
geometry.o: string.o
+geometry.o: tally.o
global.o: bank_header.o
global.o: constants.o
@@ -87,6 +91,7 @@ input_xml.o: mesh_header.o
input_xml.o: output.o
input_xml.o: string.o
input_xml.o: tally_header.o
+input_xml.o: xml-fortran/templates/cross_sections_t.o
input_xml.o: xml-fortran/templates/geometry_t.o
input_xml.o: xml-fortran/templates/materials_t.o
input_xml.o: xml-fortran/templates/settings_t.o
@@ -94,6 +99,7 @@ input_xml.o: xml-fortran/templates/tallies_t.o
interpolation.o: constants.o
interpolation.o: error.o
+interpolation.o: global.o
interpolation.o: search.o
logging.o: constants.o
@@ -160,6 +166,7 @@ plot.o: particle_header.o
search.o: constants.o
search.o: error.o
+search.o: global.o
source.o: bank_header.o
source.o: constants.o
@@ -172,9 +179,9 @@ source.o: physics.o
string.o: constants.o
string.o: error.o
+string.o: global.o
tally.o: constants.o
-tally.o: cross_section.o
tally.o: error.o
tally.o: global.o
tally.o: mesh.o
diff --git a/src/bank_header.f90 b/src/bank_header.f90
index b896382c98..0b42aa155a 100644
--- a/src/bank_header.f90
+++ b/src/bank_header.f90
@@ -9,7 +9,7 @@ module bank_header
!===============================================================================
type Bank
- integer(8) :: uid ! Unique ID
+ integer(8) :: id ! Unique ID
real(8) :: xyz(3) ! location of bank particle
real(8) :: uvw(3) ! diretional cosines
real(8) :: E ! energy
diff --git a/src/constants.f90 b/src/constants.f90
index ad0f0bc620..6ae68440dc 100644
--- a/src/constants.f90
+++ b/src/constants.f90
@@ -166,6 +166,12 @@ module constants
N_PT = 116, &
N_DA = 117
+ ! ACE table types
+ integer, parameter :: &
+ ACE_NEUTRON = 1, & ! continuous-energy neutron
+ ACE_THERMAL = 2, & ! thermal S(a,b) scattering data
+ ACE_DOSIMETRY = 3 ! dosimetry cross sections
+
! Tally macro reactions
integer, parameter :: N_MACRO_TYPES = 15
integer, parameter :: &
@@ -223,6 +229,7 @@ module constants
integer, parameter :: MAX_WORDS = 500
integer, parameter :: MAX_LINE_LEN = 250
integer, parameter :: MAX_WORD_LEN = 150
+ integer, parameter :: MAX_FILE_LEN = 255
! Unit numbers
integer, parameter :: UNIT_LOG = 11 ! unit # for writing log file
diff --git a/src/cross_section.f90 b/src/cross_section.f90
index 9043ac8b9a..885b5c9bac 100644
--- a/src/cross_section.f90
+++ b/src/cross_section.f90
@@ -1,7 +1,7 @@
module cross_section
use constants
- use cross_section_header, only: Nuclide, Reaction, SAB_Table, xsData
+ use cross_section_header, only: Nuclide, Reaction, SAB_Table, XsListing
use datatypes, only: dict_create, dict_add_key, dict_get_key, &
dict_has_key, dict_delete
use datatypes_header, only: DictionaryCI
@@ -10,7 +10,7 @@ module cross_section
use fileio, only: read_line, read_data, skip_lines
use global
use material_header, only: Material
- use output, only: message
+ use output, only: write_message
use string, only: split_string, str_to_int, str_to_real, &
lower_case, int_to_str
@@ -36,11 +36,10 @@ contains
integer :: i ! index in materials array
integer :: j ! index over nuclides in material
- integer :: index ! index in xsdatas array
+ integer :: index ! index in xs_listings array
integer :: index_nuclides ! index in nuclides
integer :: index_sab ! index in sab_tables
character(10) :: key ! name of isotope, e.g. 92235.03c
- character(MAX_LINE_LEN) :: msg ! output/error message
type(Material), pointer :: mat => null()
type(Nuclide), pointer :: nuc => null()
type(SAB_Table), pointer :: sab => null()
@@ -106,8 +105,8 @@ contains
do i = 1, n_materials
mat => materials(i)
do j = 1, mat % n_nuclides
- ! Get index in xsdatas array for this nuclide
- index = mat % xsdata(j)
+ ! Get index in xs_listings array for this nuclide
+ index = mat % xs_listing(j)
! Get name of nuclide
key = mat % names(j)
@@ -127,13 +126,13 @@ contains
key = mat % sab_name
if (.not. dict_has_key(temp_dict, key)) then
- ! Find the entry in xsdatas for this table
- if (dict_has_key(xsdata_dict, key)) then
- index = dict_get_key(xsdata_dict, key)
+ ! Find the entry in xs_listings for this table
+ if (dict_has_key(xs_listing_dict, key)) then
+ index = dict_get_key(xs_listing_dict, key)
else
- msg = "Cannot find cross-section " // trim(key) // " in specified &
- &xsdata file."
- call fatal_error(msg)
+ message = "Cannot find cross-section " // trim(key) // &
+ " in specified cross_sections.xml file."
+ call fatal_error()
end if
! Read the table and add entry to dictionary
@@ -156,9 +155,9 @@ contains
! Check to make sure S(a,b) table matched a nuclide
if (mat % sab_nuclide == 0) then
- msg = "S(a,b) table " // trim(mat % sab_name) // " did not match " &
- // "any nuclide on material " // trim(int_to_str(mat % uid))
- call fatal_error(msg)
+ message = "S(a,b) table " // trim(mat % sab_name) // " did not match " &
+ // "any nuclide on material " // trim(int_to_str(mat % id))
+ call fatal_error()
end if
end if
end do
@@ -177,7 +176,7 @@ contains
subroutine read_ACE_continuous(index_table, index)
integer, intent(in) :: index_table ! index in nuclides array
- integer, intent(in) :: index ! index in xsdatas array
+ integer, intent(in) :: index ! index in xs_listings array
integer :: in = 7 ! unit to read from
integer :: ioError ! error status for file access
@@ -188,56 +187,57 @@ contains
logical :: file_exists ! does ACE library exist?
logical :: found_xs ! did we find table in library?
character(7) :: readable ! is ACE library readable?
- character(MAX_LINE_LEN) :: msg ! output/error message
character(MAX_LINE_LEN) :: line ! single line to read
character(MAX_WORD_LEN) :: words(MAX_WORDS) ! words on a line
- character(MAX_WORD_LEN) :: filename ! name of ACE library file
+ character(MAX_FILE_LEN) :: filename ! name of ACE library file
character(10) :: tablename ! name of cross section table
type(Nuclide), pointer :: nuc => null()
- ! Check to make sure index in nuclides array and xsdata arrays are valid
+ ! Check to make sure index in nuclides array and xs_listings arrays are
+ ! valid
if (index_table > size(nuclides)) then
- msg = "Index of table to read is greater than length of nuclides."
- call fatal_error(msg)
- elseif (index > size(xsdatas)) then
- msg = "Index of xsdata entry is greater than length of xsdatas."
- call fatal_error(msg)
+ message = "Index of table to read is greater than length of nuclides."
+ call fatal_error()
+ elseif (index > size(xs_listings)) then
+ message = "Index of xs_listing entry is greater than length of " // &
+ "xs_listings."
+ call fatal_error()
end if
- filename = xsdatas(index)%path
- tablename = xsdatas(index)%id
+ filename = xs_listings(index) % path
+ tablename = xs_listings(index) % name
nuc => nuclides(index_table)
! Check if input file exists and is readable
inquire(FILE=filename, EXIST=file_exists, READ=readable)
if (.not. file_exists) then
- msg = "ACE library '" // trim(filename) // "' does not exist!"
- call fatal_error(msg)
+ message = "ACE library '" // trim(filename) // "' does not exist!"
+ call fatal_error()
elseif (readable(1:3) == 'NO') then
- msg = "ACE library '" // trim(filename) // "' is not readable! &
+ message = "ACE library '" // trim(filename) // "' is not readable! &
&Change file permissions with chmod command."
- call fatal_error(msg)
+ call fatal_error()
end if
! display message
- msg = "Loading ACE cross section table: " // tablename
- call message(msg, 6)
+ message = "Loading ACE cross section table: " // tablename
+ call write_message(6)
! open file
open(file=filename, unit=in, status='old', &
& action='read', iostat=ioError)
if (ioError /= 0) then
- msg = "Error while opening file: " // filename
- call fatal_error(msg)
+ message = "Error while opening file: " // filename
+ call fatal_error()
end if
found_xs = .false.
do while (.not. found_xs)
call read_line(in, line, ioError)
if (ioError < 0) then
- msg = "Could not find ACE table " // tablename // "."
- call fatal_error(msg)
+ message = "Could not find ACE table " // tablename // "."
+ call fatal_error()
end if
call split_string(line, words, n)
if (trim(words(1)) == trim(tablename)) then
@@ -1073,7 +1073,7 @@ contains
subroutine read_ACE_thermal(index_table, index)
integer, intent(in) :: index_table ! index in sab_tables array
- integer, intent(in) :: index ! index in xsdatas array
+ integer, intent(in) :: index ! index in xs_listings array
integer :: in = 7 ! unit to read from
integer :: ioError ! error status for file access
@@ -1084,47 +1084,46 @@ contains
logical :: file_exists ! does ACE library exist?
logical :: found_xs ! did we find table in library?
character(7) :: readable ! is ACE library readable?
- character(MAX_LINE_LEN) :: msg ! output/error message
character(MAX_LINE_LEN) :: line ! single line to read
character(MAX_WORD_LEN) :: words(MAX_WORDS) ! words on a line
character(MAX_WORD_LEN) :: filename ! name of ACE library file
character(10) :: tablename ! name of cross section table
type(SAB_Table), pointer :: table => null()
- filename = xsdatas(index)%path
- tablename = xsdatas(index)%id
+ filename = xs_listings(index) % path
+ tablename = xs_listings(index) % name
table => sab_tables(index_table)
! Check if input file exists and is readable
inquire(FILE=filename, EXIST=file_exists, READ=readable)
if (.not. file_exists) then
- msg = "ACE library '" // trim(filename) // "' does not exist!"
- call fatal_error(msg)
+ message = "ACE library '" // trim(filename) // "' does not exist!"
+ call fatal_error()
elseif (readable(1:3) == 'NO') then
- msg = "ACE library '" // trim(filename) // "' is not readable! &
+ message = "ACE library '" // trim(filename) // "' is not readable! &
&Change file permissions with chmod command."
- call fatal_error(msg)
+ call fatal_error()
end if
! display message
- msg = "Loading ACE cross section table: " // tablename
- call message(msg, 6)
+ message = "Loading ACE cross section table: " // tablename
+ call write_message(6)
! open file
open(file=filename, unit=in, status='old', &
& action='read', iostat=ioError)
if (ioError /= 0) then
- msg = "Error while opening file: " // filename
- call fatal_error(msg)
+ message = "Error while opening file: " // filename
+ call fatal_error()
end if
found_xs = .false.
do while (.not. found_xs)
call read_line(in, line, ioError)
if (ioError < 0) then
- msg = "Could not find ACE table " // tablename // "."
- call fatal_error(msg)
+ message = "Could not find ACE table " // tablename // "."
+ call fatal_error()
end if
call split_string(line, words, n)
if (trim(words(1)) == trim(tablename)) then
@@ -1226,6 +1225,8 @@ contains
! scattering
NE_out = NXS(4)
NMU = NXS(3) + 1
+ table % n_inelastic_e_out = NE_out
+ table % n_inelastic_mu = NMU
allocate(table % inelastic_e_out(NE_out, NE_in))
allocate(table % inelastic_mu(NMU, NE_out, NE_in))
@@ -1320,168 +1321,4 @@ contains
end function get_real
-!===============================================================================
-! GET_MACRO_XS
-!===============================================================================
-
- function get_macro_xs(p, mat, MT) result(xs)
-
- type(Particle), pointer :: p
- type(Material), pointer :: mat
- integer, intent(in) :: MT
- real(8) :: xs
-
- integer :: i, j
- integer :: n_nuclides
- integer :: IE
- real(8) :: density_i
- real(8) :: sigma_i
- real(8) :: f
- type(Nuclide), pointer :: nuc => null()
- type(Reaction), pointer :: rxn => null()
-
- ! initialize xs
- xs = ZERO
-
- ! loop over all nuclides in material
- n_nuclides = mat % n_nuclides
- do i = 1, n_nuclides
- nuc => nuclides(mat % nuclide(i))
-
- ! determine nuclide atom density
- density_i = mat % atom_density(i)
-
- ! search nuclide energy grid
- IE = nuc%grid_index(p % IE)
- f = (p%E - nuc%energy(IE))/(nuc%energy(IE+1) - nuc%energy(IE))
-
- ! handle special case of total cross section
- if (MT == 1) then
- xs = xs + mat % density * (ONE-f) * nuc%total(IE) + &
- & f * (nuc%total(IE+1))
- cycle
- end if
-
- ! loop over reactions in isotope
- do j = 1, nuc % n_reaction
- rxn => nuc % reactions(i)
-
- ! check for matching MT
- if (MT /= rxn % MT) cycle
-
- ! if energy is below threshold for this reaction, skip it
- if (IE < rxn % IE) cycle
-
- ! add to cumulative probability
- sigma_i = (ONE-f) * rxn%sigma(IE-rxn%IE+1) + &
- & f * (rxn%sigma(IE-rxn%IE+2))
- end do
-
- ! calculate nuclide macroscopic cross-section
- xs = xs + density_i * sigma_i
- end do
-
- end function get_macro_xs
-
-!===============================================================================
-! READ_XSDATA reads the data in a SERPENT xsdata file and builds a dictionary to
-! find cross-section information later on.
-!===============================================================================
-
- subroutine read_xsdata(path)
-
- character(*), intent(in) :: path
-
- type(xsData), pointer :: iso => null()
- character(MAX_LINE_LEN) :: line
- character(MAX_LINE_LEN) :: msg
- character(MAX_WORD_LEN) :: words(MAX_WORDS)
- character(MAX_WORD_LEN) :: filename
- integer :: n
- integer :: in = 7
- logical :: file_exists
- character(7) :: readable
- integer :: count
- integer :: index
- integer :: ioError
-
- msg = "Reading cross-section summary file..."
- call message(msg, 5)
-
- ! Construct filename
- filename = trim(path)
-
- ! Check if xsdata exists and is readable
- inquire(FILE=filename, EXIST=file_exists, READ=readable)
- if (.not. file_exists) then
- msg = "Cross section summary '" // trim(filename) // "' does not exist!"
- call fatal_error(msg)
- elseif (readable(1:3) == 'NO') then
- msg = "Cross section summary '" // trim(filename) // "' is not readable!" &
- & // "Change file permissions with chmod command."
- call fatal_error(msg)
- end if
-
- ! open xsdata file
- open(FILE=filename, UNIT=in, STATUS='old', &
- & ACTION='read', IOSTAT=ioError)
- if (ioError /= 0) then
- msg = "Error while opening file: " // filename
- call fatal_error(msg)
- end if
-
- ! determine how many lines
- count = 0
- do
- read(UNIT=in, FMT='(A)', IOSTAT=ioError) line
- if (ioError < 0) then
- ! reached end of file
- exit
- elseif (ioError > 0) then
- msg = "Unknown error while reading file: " // filename
- close(UNIT=in)
- call fatal_error(msg)
- end if
- count = count + 1
- end do
- allocate(xsdatas(count))
-
- ! read actual lines
- index = 0
- rewind(in)
- do
- read(UNIT=in, FMT='(A)', IOSTAT=ioError) line
- if (ioError < 0) exit
- index = index + 1
- call split_string(line, words, n)
- if (n == 0) cycle ! skip blank line
-
- ! Check to make sure there are enough arguments
- if (n < 9) then
- msg = "Not enough arguments on xsdata line: " // line
- close(UNIT=in)
- call fatal_error(msg)
- end if
-
- iso => xsdatas(index)
-
- ! store data
- iso%alias = words(1)
- iso%id = words(2)
- iso%type = str_to_int(words(3))
- iso%zaid = str_to_int(words(4))
- iso%isomeric = str_to_int(words(5))
- iso%awr = str_to_real(words(6))
- iso%temp = str_to_real(words(7))
- iso%binary = str_to_int(words(8))
- iso%path = words(9)
-
- ! create dictionary entry
- call dict_add_key(xsdata_dict, iso%alias, index)
- end do
-
- close(UNIT=in)
-
- end subroutine read_xsdata
-
end module cross_section
diff --git a/src/cross_section_header.f90 b/src/cross_section_header.f90
index ba6df05261..4698d1c2aa 100644
--- a/src/cross_section_header.f90
+++ b/src/cross_section_header.f90
@@ -1,5 +1,7 @@
module cross_section_header
+ use constants, only: MAX_FILE_LEN
+
implicit none
!===============================================================================
@@ -148,20 +150,20 @@ module cross_section_header
end type SAB_Table
!===============================================================================
-! XSDATA contains data read in from a SERPENT xsdata file
+! XSLISTING contains data read from a cross_sections.xml file
!===============================================================================
- type xsData
+ type XsListing
+ character(10) :: name
character(10) :: alias
- character(10) :: id
integer :: type
integer :: zaid
- integer :: isomeric
real(8) :: awr
real(8) :: temp
- integer :: binary
- character(150) :: path
- end type xsData
+ logical :: metastable
+ logical :: binary
+ character(MAX_FILE_LEN) :: path
+ end type XsListing
!===============================================================================
! NUCLIDEMICROXS contains cached microscopic cross sections for a
diff --git a/src/datatypes.f90 b/src/datatypes.f90
index 052e8d374e..855b91c141 100644
--- a/src/datatypes.f90
+++ b/src/datatypes.f90
@@ -16,9 +16,9 @@ module datatypes
implicit none
- integer, parameter :: HASH_SIZE = 4993
- integer, parameter :: HASH_MULTIPLIER = 31
- integer, parameter :: DICT_NULL = 0
+ integer, parameter, private :: HASH_SIZE = 4993
+ integer, parameter, private :: HASH_MULTIPLIER = 31
+ integer, parameter, private :: DICT_NULL = 0
!===============================================================================
! LIST Interfaces -- these allow one to use a single subroutine or function for
diff --git a/src/energy_grid.f90 b/src/energy_grid.f90
index c4b428774f..04c50bbb0e 100644
--- a/src/energy_grid.f90
+++ b/src/energy_grid.f90
@@ -4,7 +4,7 @@ module energy_grid
use datatypes, only: list_insert, list_size, list_delete
use datatypes_header, only: ListReal
use global
- use output, only: message
+ use output, only: write_message
contains
@@ -23,10 +23,9 @@ contains
type(Material), pointer :: mat => null()
type(Nuclide), pointer :: nuc => null()
integer :: i, j
- character(MAX_LINE_LEN) :: msg
- msg = "Creating unionized energy grid..."
- call message(msg, 5)
+ message = "Creating unionized energy grid..."
+ call write_message(5)
! loop over all materials
do i = 1, n_materials
diff --git a/src/error.f90 b/src/error.f90
index 83d95f695a..73e23c386b 100644
--- a/src/error.f90
+++ b/src/error.f90
@@ -2,7 +2,7 @@ module error
use ISO_FORTRAN_ENV
- use global, only: master, free_memory
+ use global, only: master, free_memory, message
implicit none
@@ -17,24 +17,22 @@ contains
! stream.
!===============================================================================
- subroutine warning(msg)
+ subroutine warning()
- character(*), intent(in) :: msg
-
- integer :: n_lines
- integer :: i
+ integer :: n_lines ! number of lines
+ integer :: i ! loop index for lines
! Only allow master to print to screen
if (.not. master) return
write(ou, fmt='(1X,A9)', advance='no') 'WARNING: '
- n_lines = (len_trim(msg)-1)/70 + 1
+ n_lines = (len_trim(message)-1)/70 + 1
do i = 1, n_lines
if (i == 1) then
- write(ou, fmt='(A70)') msg(70*(i-1)+1:70*i)
+ write(ou, fmt='(A70)') message(70*(i-1)+1:70*i)
else
- write(ou, fmt='(10X,A70)') msg(70*(i-1)+1:70*i)
+ write(ou, fmt='(10X,A70)') message(70*(i-1)+1:70*i)
end if
end do
@@ -46,23 +44,21 @@ contains
! the program is aborted.
!===============================================================================
- subroutine fatal_error(msg)
+ subroutine fatal_error()
- character(*), intent(in) :: msg
-
- integer :: n_lines
- integer :: i
+ integer :: n_lines ! number of lines
+ integer :: i ! loop index over lines
! Only allow master to print to screen
if (master) then
write(eu, fmt='(1X,A7)', advance='no') 'ERROR: '
- n_lines = (len_trim(msg)-1)/72 + 1
+ n_lines = (len_trim(message)-1)/72 + 1
do i = 1, n_lines
if (i == 1) then
- write(eu, fmt='(A72)') msg(72*(i-1)+1:72*i)
+ write(eu, fmt='(A72)') message(72*(i-1)+1:72*i)
else
- write(eu, fmt='(7X,A72)') msg(72*(i-1)+1:72*i)
+ write(eu, fmt='(7X,A72)') message(72*(i-1)+1:72*i)
end if
end do
write(eu,*)
diff --git a/src/fileio.f90 b/src/fileio.f90
index 001da2c177..96a7f83c8b 100644
--- a/src/fileio.f90
+++ b/src/fileio.f90
@@ -89,8 +89,8 @@ contains
integer, intent(in) :: n_lines ! number of lines to skip
integer, intent(out) :: ioError ! error status
- integer :: i ! index for number of lines
- character(MAX_LINE_LEN) :: tmp ! single line
+ integer :: i ! index for number of lines
+ character(MAX_LINE_LEN) :: tmp ! single line
do i = 1, n_lines
read(UNIT=unit, FMT='(A)', IOSTAT=ioError) tmp
@@ -110,8 +110,8 @@ contains
integer, intent(in) :: lines ! total number of lines
integer, intent(in) :: words_per_line ! number of words per line
- integer :: i ! line index
- integer :: loc ! locator for array
+ integer :: i ! line index
+ integer :: loc ! locator for array
loc = 0
do i = 1, lines
@@ -137,8 +137,8 @@ contains
integer, intent(in) :: lines ! total number of lines
integer, intent(in) :: words_per_line ! number of words per line
- integer :: i ! line index
- integer :: loc ! locator for array
+ integer :: i ! line index
+ integer :: loc ! locator for array
loc = 0
do i = 1, lines
diff --git a/src/fission.f90 b/src/fission.f90
index d26a1101e9..c544385361 100644
--- a/src/fission.f90
+++ b/src/fission.f90
@@ -3,6 +3,7 @@ module fission
use constants
use cross_section_header, only: Nuclide
use error, only: fatal_error
+ use global, only: message
use interpolation, only: interpolate_tab1
use search, only: binary_search
@@ -21,14 +22,13 @@ contains
real(8), intent(in) :: E ! energy of incoming neutron
real(8) :: nu ! number of total neutrons emitted per fission
- integer :: i ! loop index
- integer :: NC ! number of polynomial coefficients
- real(8) :: c ! polynomial coefficient
- character(MAX_LINE_LEN) :: msg ! error message
+ integer :: i ! loop index
+ integer :: NC ! number of polynomial coefficients
+ real(8) :: c ! polynomial coefficient
if (nuc % nu_t_type == NU_NONE) then
- msg = "No neutron emission data for table: " // nuc % name
- call fatal_error(msg)
+ message = "No neutron emission data for table: " // nuc % name
+ call fatal_error()
elseif (nuc % nu_t_type == NU_POLYNOMIAL) then
! determine number of coefficients
NC = int(nuc % nu_t_data(1))
@@ -57,9 +57,9 @@ contains
real(8), intent(in) :: E ! energy of incoming neutron
real(8) :: nu ! number of prompt neutrons emitted per fission
- integer :: i ! loop index
- integer :: NC ! number of polynomial coefficients
- real(8) :: c ! polynomial coefficient
+ integer :: i ! loop index
+ integer :: NC ! number of polynomial coefficients
+ real(8) :: c ! polynomial coefficient
if (nuc % nu_p_type == NU_NONE) then
! since no prompt or delayed data is present, this means all neutron
diff --git a/src/geometry.f90 b/src/geometry.f90
index 6b2254cf8f..9d5d4f7bd6 100644
--- a/src/geometry.f90
+++ b/src/geometry.f90
@@ -5,9 +5,10 @@ module geometry
use error, only: fatal_error
use geometry_header, only: Cell, Surface, Universe, Lattice
use global
- use output, only: message
+ use output, only: write_message
use particle_header, only: Particle
use string, only: int_to_str
+ use tally, only: score_surface_current
implicit none
@@ -95,14 +96,12 @@ contains
type(Particle), pointer :: p ! pointer to particle
logical, intent(inout) :: found ! particle found?
- character(MAX_LINE_LEN) :: msg ! error message
integer :: i ! index over cells
integer :: x, y
type(Cell), pointer :: c ! pointer to cell
type(Lattice), pointer :: lat ! pointer to lattice
- type(Universe), pointer :: lower_univ ! if particle is in lower
- ! universe, use this pointer
- ! to call recursively
+ type(Universe), pointer :: lower_univ ! if particle is in lower universe,
+ ! use this pointer to call recursively
found = .false.
@@ -118,7 +117,7 @@ contains
! set particle attributes
p % cell = univ % cells(i)
- p % universe = dict_get_key(universe_dict, univ % uid)
+ p % universe = dict_get_key(universe_dict, univ % id)
p % material = c % material
exit
elseif (c % type == CELL_FILL) then
@@ -127,8 +126,8 @@ contains
if (found) then
exit
else
- msg = "Could not locate particle in universe: "
- call fatal_error(msg)
+ message = "Could not locate particle in universe: "
+ call fatal_error()
end if
elseif (c % type == CELL_LATTICE) then
! Set current lattice
@@ -153,9 +152,9 @@ contains
if (found) then
exit
else
- msg = "Could not locate particle in lattice: " &
- & // int_to_str(lat % uid)
- call fatal_error(msg)
+ message = "Could not locate particle in lattice: " &
+ & // int_to_str(lat % id)
+ call fatal_error()
end if
end if
end if
@@ -191,7 +190,6 @@ contains
real(8) :: dot_prod ! dot product of direction and normal
real(8) :: norm ! "norm" of surface normal
logical :: found ! particle found in universe?
- character(MAX_LINE_LEN) :: msg ! output/error message?
type(Surface), pointer :: surf => null()
type(Cell), pointer :: c => null()
type(Lattice), pointer :: lat => null()
@@ -199,18 +197,31 @@ contains
surf => surfaces(abs(p % surface))
if (verbosity >= 10) then
- msg = " Crossing surface " // trim(int_to_str(surf % uid))
- call message(msg)
+ message = " Crossing surface " // trim(int_to_str(surf % id))
+ call write_message()
end if
if (surf % bc == BC_VACUUM .and. (.not. plotting)) then
! =======================================================================
! PARTICLE LEAKS OUT OF PROBLEM
+ ! Kill particle
p % alive = .false.
+
+ ! Score any surface current tallies -- note that the particle is moved
+ ! forward slightly so that if the mesh boundary is on the surface, it is
+ ! still processed
+
+ ! TODO: Find a better solution to score surface currents than physically
+ ! moving the particle forward slightly
+
+ p % xyz = p % xyz + 1e-6 * p % uvw
+ call score_surface_current(p)
+
+ ! Display message
if (verbosity >= 10) then
- msg = " Leaked out of surface " // trim(int_to_str(surf % uid))
- call message(msg)
+ message = " Leaked out of surface " // trim(int_to_str(surf % id))
+ call write_message()
end if
return
@@ -247,9 +258,9 @@ contains
p % uvw = (/ u, v, w /)
case (SURF_CYL_X)
! Find y-y0, z-z0 and dot product of direction and surface normal
- y = p % xyz(2) - surf % coeffs(2)
- z = p % xyz(3) - surf % coeffs(3)
- R = surf % coeffs(4)
+ y = p % xyz(2) - surf % coeffs(1)
+ z = p % xyz(3) - surf % coeffs(2)
+ R = surf % coeffs(3)
dot_prod = v*y + w*z
! Reflect direction according to normal
@@ -261,8 +272,8 @@ contains
case (SURF_CYL_Y)
! Find x-x0, z-z0 and dot product of direction and surface normal
x = p % xyz(1) - surf % coeffs(1)
- z = p % xyz(3) - surf % coeffs(3)
- R = surf % coeffs(4)
+ z = p % xyz(3) - surf % coeffs(2)
+ R = surf % coeffs(3)
dot_prod = u*x + w*z
! Reflect direction according to normal
@@ -275,7 +286,7 @@ contains
! Find x-x0, y-y0 and dot product of direction and surface normal
x = p % xyz(1) - surf % coeffs(1)
y = p % xyz(2) - surf % coeffs(2)
- R = surf % coeffs(4)
+ R = surf % coeffs(3)
dot_prod = u*x + v*y
! Reflect direction according to normal
@@ -301,19 +312,26 @@ contains
! Set vector
p % uvw = (/ u, v, w /)
case default
- msg = "Reflection not supported for surface " // &
- trim(int_to_str(surf % uid))
- call fatal_error(msg)
+ message = "Reflection not supported for surface " // &
+ trim(int_to_str(surf % id))
+ call fatal_error()
end select
! Reassign particle's cell and surface
p % cell = last_cell
p % surface = -p % surface
+ ! Score surface currents since reflection causes the direction of the
+ ! particle to change
+ call score_surface_current(p)
+
+ ! Set previous coordinate going slightly past surface crossing
+ p % last_xyz = p % xyz + 1e-6 * p % uvw
+
! Diagnostic message
if (verbosity >= 10) then
- msg = " Reflected from surface " // trim(int_to_str(surf%uid))
- call message(msg)
+ message = " Reflected from surface " // trim(int_to_str(surf%id))
+ call write_message()
end if
return
end if
@@ -332,8 +350,8 @@ contains
lower_univ => universes(c % fill)
call find_cell(lower_univ, p, found)
if (.not. found) then
- msg = "Could not locate particle in universe: "
- call fatal_error(msg)
+ message = "Could not locate particle in universe: "
+ call fatal_error()
end if
elseif (c % type == CELL_LATTICE) then
! Set current lattice
@@ -356,9 +374,9 @@ contains
call find_cell(lower_univ, p, found)
if (.not. found) then
- msg = "Could not locate particle in lattice: " &
- & // int_to_str(lat % uid)
- call fatal_error(msg)
+ message = "Could not locate particle in lattice: " // &
+ trim(int_to_str(lat % id))
+ call fatal_error()
end if
else
! set current pointers
@@ -379,8 +397,8 @@ contains
lower_univ => universes(c % fill)
call find_cell(lower_univ, p, found)
if (.not. found) then
- msg = "Could not locate particle in universe: "
- call fatal_error(msg)
+ message = "Could not locate particle in universe: "
+ call fatal_error()
end if
elseif (c % type == CELL_LATTICE) then
! Set current lattice
@@ -403,9 +421,9 @@ contains
call find_cell(lower_univ, p, found)
if (.not. found) then
- msg = "Could not locate particle in lattice: " &
- & // int_to_str(lat % uid)
- call fatal_error(msg)
+ message = "Could not locate particle in lattice: " // &
+ trim(int_to_str(lat % id))
+ call fatal_error()
end if
else
! set current pointers
@@ -430,9 +448,9 @@ contains
end do
! Couldn't find next cell anywhere!
- msg = "After particle crossed surface " // trim(int_to_str(p%surface)) // &
- & ", it could not be located in any cell and it did not leak."
- call fatal_error(msg)
+ message = "After particle crossed surface " // trim(int_to_str(p%surface)) &
+ // ", it could not be located in any cell and it did not leak."
+ call fatal_error()
end subroutine cross_surface
@@ -459,13 +477,12 @@ contains
real(8) :: x0 ! half the width of lattice element
real(8) :: y0 ! half the height of lattice element
logical :: found ! particle found in cell?
- character(MAX_LINE_LEN) :: msg ! output/error message
type(Lattice), pointer :: lat
type(Universe), pointer :: univ
if (verbosity >= 10) then
- msg = " Crossing lattice"
- call message(msg)
+ message = " Crossing lattice"
+ call write_message()
end if
lat => lattices(p % lattice)
@@ -530,11 +547,11 @@ contains
i_x = p % index_x
i_y = p % index_y
if (i_x < 1 .or. i_x > lat % n_x) then
- msg = "Reached edge of lattice."
- call fatal_error(msg)
+ message = "Reached edge of lattice."
+ call fatal_error()
elseif (i_y < 1 .or. i_y > lat % n_y) then
- msg = "Reached edge of lattice."
- call fatal_error(msg)
+ message = "Reached edge of lattice."
+ call fatal_error()
end if
! Find universe for next lattice element
@@ -543,8 +560,8 @@ contains
! Find cell in next lattice element
call find_cell(univ, p, found)
if (.not. found) then
- msg = "Could not locate particle in universe: "
- call fatal_error(msg)
+ message = "Could not locate particle in universe: "
+ call fatal_error()
end if
end subroutine cross_lattice
@@ -578,7 +595,6 @@ contains
real(8) :: a,b,c,k ! quadratic equation coefficients
real(8) :: quad ! discriminant of quadratic equation
logical :: on_surface ! is particle on surface?
- character(MAX_LINE_LEN) :: msg ! output/error message
type(Cell), pointer :: cell_p => null()
type(Cell), pointer :: parent_p => null()
type(Surface), pointer :: surf_p => null()
@@ -882,8 +898,8 @@ contains
end if
case (SURF_GQ)
- msg = "Surface distance not yet implement for general quadratic."
- call fatal_error(msg)
+ message = "Surface distance not yet implement for general quadratic."
+ call fatal_error()
end select
@@ -1114,12 +1130,11 @@ contains
integer, allocatable :: count_positive(:) ! # of cells on positive side
integer, allocatable :: count_negative(:) ! # of cells on negative side
logical :: positive ! positive side specified in surface list
- character(MAX_LINE_LEN) :: msg ! output/error message
type(Cell), pointer :: c
type(Surface), pointer :: surf
- msg = "Building neighboring cells lists for each surface..."
- call message(msg, 4)
+ message = "Building neighboring cells lists for each surface..."
+ call write_message(4)
allocate(count_positive(n_surfaces))
allocate(count_negative(n_surfaces))
diff --git a/src/geometry_header.f90 b/src/geometry_header.f90
index 7b44c0a570..5b3858755f 100644
--- a/src/geometry_header.f90
+++ b/src/geometry_header.f90
@@ -7,7 +7,7 @@ module geometry_header
!===============================================================================
type Universe
- integer :: uid ! Unique ID
+ integer :: id ! Unique ID
integer :: type ! Type
integer :: level ! Level of universe (0=base)
integer :: n_cells ! # of cells within
@@ -15,7 +15,6 @@ module geometry_header
real(8) :: x0 ! Translation in x-coordinate
real(8) :: y0 ! Translation in y-coordinate
real(8) :: z0 ! Translation in z-coordinate
- integer, allocatable :: tallies(:)
end type Universe
!===============================================================================
@@ -24,7 +23,7 @@ module geometry_header
!===============================================================================
type Lattice
- integer :: uid ! Universe number for lattice
+ integer :: id ! Universe number for lattice
integer :: type ! Type of lattice (rectangular, hex, etc)
integer :: level ! Level of lattice
integer :: n_x ! number of lattice cells in x-direction
@@ -34,7 +33,6 @@ module geometry_header
real(8) :: width_x ! width of lattice cell
real(8) :: width_y ! width of lattice cell
integer, allocatable :: element(:,:) ! specified universes
- integer, allocatable :: tallies(:)
end type Lattice
!===============================================================================
@@ -43,7 +41,7 @@ module geometry_header
!===============================================================================
type Surface
- integer :: uid ! Unique ID
+ integer :: id ! Unique ID
integer :: type ! Type of surface
real(8), allocatable :: coeffs(:) ! Definition of surface
integer, allocatable :: &
@@ -57,7 +55,7 @@ module geometry_header
!===============================================================================
type Cell
- integer :: uid ! Unique ID
+ integer :: id ! Unique ID
integer :: type ! Type of cell (normal, universe, lattice)
integer :: universe ! universe # this cell is in
integer :: fill ! universe # filling this cell
@@ -68,7 +66,6 @@ module geometry_header
& surfaces(:) ! List of surfaces bounding cell -- note that
! parentheses, union, etc operators will be listed
! here too
- integer, allocatable :: tallies(:)
end type Cell
! array index of universe 0
diff --git a/src/global.f90 b/src/global.f90
index 7a4f128414..9e601b24cc 100644
--- a/src/global.f90
+++ b/src/global.f90
@@ -2,8 +2,8 @@ module global
use bank_header, only: Bank
use constants
- use cross_section_header, only: Nuclide, SAB_Table, xsData, NuclideMicroXS, &
- MaterialMacroXS
+ use cross_section_header, only: Nuclide, SAB_Table, xsListing, &
+ NuclideMicroXS, MaterialMacroXS
use datatypes_header, only: DictionaryII, DictionaryCI
use geometry_header, only: Cell, Universe, Lattice, Surface
use material_header, only: Material
@@ -20,60 +20,79 @@ module global
implicit none
save
+ ! ============================================================================
+ ! GEOMETRY-RELATED VARIABLES
+
! Main arrays
- type(Cell), allocatable, target :: cells(:)
- type(Universe), allocatable, target :: universes(:)
- type(Lattice), allocatable, target :: lattices(:)
- type(Surface), allocatable, target :: surfaces(:)
- type(Material), allocatable, target :: materials(:)
- type(xsData), allocatable, target :: xsdatas(:)
- type(StructuredMesh), allocatable, target :: meshes(:)
- type(TallyObject), allocatable, target :: tallies(:)
+ type(Cell), allocatable, target :: cells(:)
+ type(Universe), allocatable, target :: universes(:)
+ type(Lattice), allocatable, target :: lattices(:)
+ type(Surface), allocatable, target :: surfaces(:)
+ type(Material), allocatable, target :: materials(:)
! Size of main arrays
- integer :: n_cells ! # of cells
- integer :: n_universes ! # of universes
- integer :: n_lattices ! # of lattices
- integer :: n_surfaces ! # of surfaces
- integer :: n_materials ! # of materials
- integer :: n_meshes ! # of structured meshes
- integer :: n_tallies ! # of tallies
+ integer :: n_cells ! # of cells
+ integer :: n_universes ! # of universes
+ integer :: n_lattices ! # of lattices
+ integer :: n_surfaces ! # of surfaces
+ integer :: n_materials ! # of materials
! These dictionaries provide a fast lookup mechanism -- the key is the
! user-specified identifier and the value is the index in the corresponding
! array
- type(DictionaryII), pointer :: cell_dict
- type(DictionaryII), pointer :: universe_dict
- type(DictionaryII), pointer :: lattice_dict
- type(DictionaryII), pointer :: surface_dict
- type(DictionaryII), pointer :: material_dict
- type(DictionaryII), pointer :: mesh_dict
- type(DictionaryII), pointer :: tally_dict
- type(DictionaryCI), pointer :: xsdata_dict
- type(DictionaryCI), pointer :: nuclide_dict
- type(DictionaryCI), pointer :: sab_dict
+ type(DictionaryII), pointer :: cell_dict => null()
+ type(DictionaryII), pointer :: universe_dict => null()
+ type(DictionaryII), pointer :: lattice_dict => null()
+ type(DictionaryII), pointer :: surface_dict => null()
+ type(DictionaryII), pointer :: material_dict => null()
+ type(DictionaryII), pointer :: mesh_dict => null()
+ type(DictionaryII), pointer :: tally_dict => null()
+
+ ! ============================================================================
+ ! CROSS SECTION RELATED VARIABLES
! Cross section arrays
- type(Nuclide), allocatable, target :: nuclides(:)
- type(SAB_Table), allocatable, target :: sab_tables(:)
- integer :: n_nuclides_total
- integer :: n_sab_tables
+ type(Nuclide), allocatable, target :: nuclides(:) ! Nuclide cross-sections
+ type(SAB_Table), allocatable, target :: sab_tables(:) ! S(a,b) tables
+ type(XsListing), allocatable, target :: xs_listings(:) ! cross_sections.xml listings
! Cross section caches
- type(NuclideMicroXS), allocatable :: micro_xs(:)
- type(MaterialMacroXS) :: material_xs
+ type(NuclideMicroXS), allocatable :: micro_xs(:) ! Cache for each nuclide
+ type(MaterialMacroXS) :: material_xs ! Cache for current material
- ! Tally map structure
- type(TallyMap), allocatable :: tally_maps(:)
+ integer :: n_nuclides_total ! Number of nuclide cross section tables
+ integer :: n_sab_tables ! Number of S(a,b) thermal scattering tables
+
+ ! Dictionaries to look up cross sections and listings
+ type(DictionaryCI), pointer :: nuclide_dict => null()
+ type(DictionaryCI), pointer :: sab_dict => null()
+ type(DictionaryCI), pointer :: xs_listing_dict => null()
! Unionized energy grid
integer :: n_grid ! number of points on unionized grid
real(8), allocatable :: e_grid(:) ! energies on unionized grid
- ! Histories/cycles/etc for both external source and criticality
- integer(8) :: n_particles ! # of particles (per cycle for criticality)
- integer :: n_cycles ! # of cycles
- integer :: n_inactive ! # of inactive cycles
+ ! ============================================================================
+ ! TALLY-RELATED VARIABLES
+
+ type(StructuredMesh), allocatable, target :: meshes(:)
+ type(TallyObject), allocatable, target :: tallies(:)
+
+ ! Tally map structure
+ type(TallyMap), allocatable :: tally_maps(:)
+
+ integer :: n_meshes ! # of structured meshes
+ integer :: n_tallies ! # of tallies
+
+ ! Flag for turning tallies on
+ logical :: tallies_on
+
+ ! ============================================================================
+ ! CRITICALITY SIMULATION VARIABLES
+
+ integer(8) :: n_particles = 10000 ! # of particles per cycle
+ integer :: n_cycles = 500 ! # of cycles
+ integer :: n_inactive = 50 ! # of inactive cycles
! External source
type(ExtSource), target :: external_source
@@ -88,71 +107,60 @@ module global
integer(8) :: source_index ! index for source particles
! cycle keff
- real(8) :: keff
+ real(8) :: keff = ONE
real(8) :: keff_std
- ! Flag for turning tallies on
- logical :: tallies_on
+ ! ============================================================================
+ ! PARALLEL PROCESSING VARIABLES
- ! Parallel processing variables
integer :: n_procs ! number of processes
integer :: rank ! rank of process
logical :: master ! master process?
logical :: mpi_enabled ! is MPI in use and initialized?
- ! Timing variables
+ ! ============================================================================
+ ! TIMING VARIABLES
+
type(Timer) :: time_total ! timer for total run
type(Timer) :: time_init ! timer for initialization
type(Timer) :: time_intercycle ! timer for intercycle synchronization
type(Timer) :: time_inactive ! timer for inactive cycles
type(Timer) :: time_compute ! timer for computation
- ! Paths to input file, cross section data, etc
- character(MAX_WORD_LEN) :: &
- & path_input, &
- & path_xsdata
+ ! ===========================================================================
+ ! VARIANCE REDUCTION VARIABLES
- ! Problem type
- integer :: problem_type
-
- ! The verbosity controls how much information will be printed to the
- ! screen and in logs
- integer :: verbosity
-
- ! Variance reduction options
logical :: survival_biasing = .false.
real(8) :: weight_cutoff = 0.25
real(8) :: weight_survive = 1.0
- ! Plotting options
+ ! ============================================================================
+ ! PLOTTING VARIABLES
+
logical :: plotting = .false.
real(8) :: plot_origin(3)
real(8) :: plot_width(2)
real(8) :: plot_basis(6)
real(8) :: pixel
+ ! ============================================================================
+ ! MISCELLANEOUS VARIABLES
+
+ character(MAX_WORD_LEN) :: path_input ! Path to input file
+ character(MAX_FILE_LEN) :: path_cross_sections ! Path to cross_sections.xml
+
+ ! Message used in message/warning/fatal_error
+ character(MAX_LINE_LEN) :: message
+
+ ! Problem type
+ integer :: problem_type = PROB_CRITICALITY
+
+ ! The verbosity controls how much information will be printed to the
+ ! screen and in logs
+ integer :: verbosity = 7
+
contains
-!===============================================================================
-! SET_DEFAULTS gives default values for many global parameters
-!===============================================================================
-
- subroutine set_defaults()
-
- ! Default problem type is external source
- problem_type = PROB_SOURCE
-
- ! Default number of particles
- n_particles = 10000
-
- ! Default verbosity
- verbosity = 5
-
- ! Defualt multiplication factor
- keff = ONE
-
- end subroutine set_defaults
-
!===============================================================================
! FREE_MEMORY deallocates all allocatable arrays in the program, namely the
! cells, surfaces, materials, and sources
@@ -170,10 +178,10 @@ contains
if (allocated(materials)) deallocate(materials)
if (allocated(lattices)) deallocate(lattices)
- ! Deallocate cross section data
- if (allocated(xsdatas)) deallocate(xsdatas)
+ ! Deallocate cross section data and listings
if (allocated(nuclides)) deallocate(nuclides)
if (allocated(sab_tables)) deallocate(sab_tables)
+ if (allocated(xs_listings)) deallocate(xs_listings)
! Deallocate energy grid
if (allocated(e_grid)) deallocate(e_grid)
diff --git a/src/initialize.f90 b/src/initialize.f90
index 82e0440cd4..15fa8bdc85 100644
--- a/src/initialize.f90
+++ b/src/initialize.f90
@@ -1,16 +1,17 @@
module initialize
use constants
- use cross_section, only: read_xs, read_xsdata
+ use cross_section, only: read_xs
use datatypes, only: dict_create, dict_add_key, dict_get_key, &
- dict_has_key, DICT_NULL, dict_keys
+ dict_has_key, dict_keys
use datatypes_header, only: ListKeyValueII, DictionaryII
use energy_grid, only: unionized_grid, original_indices
use error, only: fatal_error
use geometry, only: neighbor_lists
use geometry_header, only: Cell, Surface, Universe, Lattice, BASE_UNIVERSE
use global
- use input_xml, only: read_input_xml, cells_in_univ_dict
+ use input_xml, only: read_input_xml, read_cross_sections_xml, &
+ cells_in_univ_dict
use logging, only: create_log
use mcnp_random, only: RN_init_problem
use mpi_routines, only: setup_mpi
@@ -55,9 +56,6 @@ contains
! 63-bit RNGs.
call RN_init_problem(3, 0_8, 0_8, 0_8, 0)
- ! Set default values for settings
- call set_defaults()
-
! set up dictionaries
call create_dictionaries()
@@ -81,10 +79,10 @@ contains
if (.not. plotting) then
! Read cross section summary file to determine what files contain
! cross-sections
- call read_xsdata(path_xsdata)
+ call read_cross_sections_xml()
- ! With the AWRs from the xsdata, change all material specifications so that
- ! they contain atom percents summing to 1
+ ! With the AWRs from the xs_listings, change all material specifications
+ ! so that they contain atom percents summing to 1
call normalize_ao()
! Read ACE-format cross sections
@@ -185,7 +183,7 @@ contains
call dict_create(material_dict)
call dict_create(mesh_dict)
call dict_create(tally_dict)
- call dict_create(xsdata_dict)
+ call dict_create(xs_listing_dict)
call dict_create(nuclide_dict)
call dict_create(sab_dict)
@@ -214,18 +212,18 @@ contains
! We also need to allocate the cell count lists for each universe. The logic
! for this is a little more convoluted. In universe_dict, the (key,value)
- ! pairs are the uid of the universe and the index in the array. In
- ! ucount_dict, it's the uid of the universe and the number of cells.
+ ! pairs are the id of the universe and the index in the array. In
+ ! ucount_dict, it's the id of the universe and the number of cells.
key_list => dict_keys(universe_dict)
do while (associated(key_list))
! find index of universe in universes array
index = key_list%data%value
univ => universes(index)
- univ % uid = key_list%data%key
+ univ % id = key_list%data%key
! check for lowest level universe
- if (univ % uid == 0) BASE_UNIVERSE = index
+ if (univ % id == 0) BASE_UNIVERSE = index
! find cell count for this universe
count = dict_get_key(cells_in_univ_dict, key_list%data%key)
@@ -273,8 +271,7 @@ contains
integer :: j ! index over surface list
integer :: k
integer :: index ! index in surfaces/materials array
- integer :: uid ! user-specified uid
- character(MAX_LINE_LEN) :: msg ! output/error message
+ integer :: id ! user-specified id
type(Cell), pointer :: c => null()
type(Lattice), pointer :: l => null()
type(TallyObject), pointer :: t => null()
@@ -285,15 +282,15 @@ contains
c => cells(i)
do j = 1, c % n_surfaces
- uid = c % surfaces(j)
- if (uid < OP_DIFFERENCE) then
- if (dict_has_key(surface_dict, abs(uid))) then
- index = dict_get_key(surface_dict, abs(uid))
- c % surfaces(j) = sign(index, uid)
+ id = c % surfaces(j)
+ if (id < OP_DIFFERENCE) then
+ if (dict_has_key(surface_dict, abs(id))) then
+ index = dict_get_key(surface_dict, abs(id))
+ c % surfaces(j) = sign(index, id)
else
- msg = "Could not find surface " // trim(int_to_str(abs(uid))) // &
- & " specified on cell " // trim(int_to_str(c % uid))
- call fatal_error(msg)
+ message = "Could not find surface " // trim(int_to_str(abs(id))) // &
+ & " specified on cell " // trim(int_to_str(c % id))
+ call fatal_error()
end if
end if
end do
@@ -301,40 +298,40 @@ contains
! =======================================================================
! ADJUST UNIVERSE INDEX FOR EACH CELL
- uid = c % universe
- if (dict_has_key(universe_dict, uid)) then
- c % universe = dict_get_key(universe_dict, uid)
+ id = c % universe
+ if (dict_has_key(universe_dict, id)) then
+ c % universe = dict_get_key(universe_dict, id)
else
- msg = "Could not find universe " // trim(int_to_str(uid)) // &
- " specified on cell " // trim(int_to_str(c % uid))
- call fatal_error(msg)
+ message = "Could not find universe " // trim(int_to_str(id)) // &
+ " specified on cell " // trim(int_to_str(c % id))
+ call fatal_error()
end if
! =======================================================================
! ADJUST MATERIAL/FILL POINTERS FOR EACH CELL
- uid = c % material
- if (uid /= 0) then
- if (dict_has_key(material_dict, uid)) then
+ id = c % material
+ if (id /= 0) then
+ if (dict_has_key(material_dict, id)) then
c % type = CELL_NORMAL
- c % material = dict_get_key(material_dict, uid)
+ c % material = dict_get_key(material_dict, id)
else
- msg = "Could not find material " // trim(int_to_str(uid)) // &
- " specified on cell " // trim(int_to_str(c % uid))
- call fatal_error(msg)
+ message = "Could not find material " // trim(int_to_str(id)) // &
+ " specified on cell " // trim(int_to_str(c % id))
+ call fatal_error()
end if
else
- uid = c % fill
- if (dict_has_key(universe_dict, uid)) then
+ id = c % fill
+ if (dict_has_key(universe_dict, id)) then
c % type = CELL_FILL
- c % fill = dict_get_key(universe_dict, uid)
- elseif (dict_has_key(lattice_dict, uid)) then
+ c % fill = dict_get_key(universe_dict, id)
+ elseif (dict_has_key(lattice_dict, id)) then
c % type = CELL_LATTICE
- c % fill = dict_get_key(lattice_dict, uid)
+ c % fill = dict_get_key(lattice_dict, id)
else
- msg = "Specified fill " // trim(int_to_str(uid)) // " on cell " // &
- trim(int_to_str(c % uid)) // " is neither a universe nor a lattice."
- call fatal_error(msg)
+ message = "Specified fill " // trim(int_to_str(id)) // " on cell " // &
+ trim(int_to_str(c % id)) // " is neither a universe nor a lattice."
+ call fatal_error()
end if
end if
end do
@@ -346,13 +343,13 @@ contains
l => lattices(i)
do j = 1, l % n_x
do k = 1, l % n_y
- uid = l % element(j,k)
- if (dict_has_key(universe_dict, uid)) then
- l % element(j,k) = dict_get_key(universe_dict, uid)
+ id = l % element(j,k)
+ if (dict_has_key(universe_dict, id)) then
+ l % element(j,k) = dict_get_key(universe_dict, id)
else
- msg = "Invalid universe number " // trim(int_to_str(uid)) &
- // " specified on lattice " // trim(int_to_str(l % uid))
- call fatal_error(msg)
+ message = "Invalid universe number " // trim(int_to_str(id)) &
+ // " specified on lattice " // trim(int_to_str(l % id))
+ call fatal_error()
end if
end do
end do
@@ -366,13 +363,13 @@ contains
if (t % n_bins(T_CELL) > 0) then
do j = 1, t % n_bins(T_CELL)
- uid = t % cell_bins(j) % scalar
- if (dict_has_key(cell_dict, uid)) then
- t % cell_bins(j) % scalar = dict_get_key(cell_dict, uid)
+ id = t % cell_bins(j) % scalar
+ if (dict_has_key(cell_dict, id)) then
+ t % cell_bins(j) % scalar = dict_get_key(cell_dict, id)
else
- msg = "Could not find cell " // trim(int_to_str(uid)) // &
- & " specified on tally " // trim(int_to_str(t % uid))
- call fatal_error(msg)
+ message = "Could not find cell " // trim(int_to_str(id)) // &
+ & " specified on tally " // trim(int_to_str(t % id))
+ call fatal_error()
end if
end do
end if
@@ -382,13 +379,13 @@ contains
if (t % n_bins(T_SURFACE) > 0) then
do j = 1, t % n_bins(T_SURFACE)
- uid = t % surface_bins(j) % scalar
- if (dict_has_key(surface_dict, uid)) then
- t % surface_bins(j) % scalar = dict_get_key(surface_dict, uid)
+ id = t % surface_bins(j) % scalar
+ if (dict_has_key(surface_dict, id)) then
+ t % surface_bins(j) % scalar = dict_get_key(surface_dict, id)
else
- msg = "Could not find surface " // trim(int_to_str(uid)) // &
- & " specified on tally " // trim(int_to_str(t % uid))
- call fatal_error(msg)
+ message = "Could not find surface " // trim(int_to_str(id)) // &
+ & " specified on tally " // trim(int_to_str(t % id))
+ call fatal_error()
end if
end do
end if
@@ -398,13 +395,13 @@ contains
if (t % n_bins(T_UNIVERSE) > 0) then
do j = 1, t % n_bins(T_UNIVERSE)
- uid = t % universe_bins(j) % scalar
- if (dict_has_key(universe_dict, uid)) then
- t % universe_bins(j) % scalar = dict_get_key(universe_dict, uid)
+ id = t % universe_bins(j) % scalar
+ if (dict_has_key(universe_dict, id)) then
+ t % universe_bins(j) % scalar = dict_get_key(universe_dict, id)
else
- msg = "Could not find universe " // trim(int_to_str(uid)) // &
- & " specified on tally " // trim(int_to_str(t % uid))
- call fatal_error(msg)
+ message = "Could not find universe " // trim(int_to_str(id)) // &
+ & " specified on tally " // trim(int_to_str(t % id))
+ call fatal_error()
end if
end do
end if
@@ -414,13 +411,13 @@ contains
if (t % n_bins(T_MATERIAL) > 0) then
do j = 1, t % n_bins(T_MATERIAL)
- uid = t % material_bins(j) % scalar
- if (dict_has_key(material_dict, uid)) then
- t % material_bins(j) % scalar = dict_get_key(material_dict, uid)
+ id = t % material_bins(j) % scalar
+ if (dict_has_key(material_dict, id)) then
+ t % material_bins(j) % scalar = dict_get_key(material_dict, id)
else
- msg = "Could not find material " // trim(int_to_str(uid)) // &
- & " specified on tally " // trim(int_to_str(t % uid))
- call fatal_error(msg)
+ message = "Could not find material " // trim(int_to_str(id)) // &
+ & " specified on tally " // trim(int_to_str(t % id))
+ call fatal_error()
end if
end do
end if
@@ -430,13 +427,13 @@ contains
if (t % n_bins(T_CELLBORN) > 0) then
do j = 1, t % n_bins(T_CELLBORN)
- uid = t % cellborn_bins(j) % scalar
- if (dict_has_key(cell_dict, uid)) then
- t % cellborn_bins(j) % scalar = dict_get_key(cell_dict, uid)
+ id = t % cellborn_bins(j) % scalar
+ if (dict_has_key(cell_dict, id)) then
+ t % cellborn_bins(j) % scalar = dict_get_key(cell_dict, id)
else
- msg = "Could not find material " // trim(int_to_str(uid)) // &
- & " specified on tally " // trim(int_to_str(t % uid))
- call fatal_error(msg)
+ message = "Could not find material " // trim(int_to_str(id)) // &
+ & " specified on tally " // trim(int_to_str(t % id))
+ call fatal_error()
end if
end do
end if
@@ -445,13 +442,13 @@ contains
! ADJUST MESH INDICES FOR EACH TALLY
if (t % n_bins(T_MESH) > 0) then
- uid = t % mesh
- if (dict_has_key(mesh_dict, uid)) then
- t % mesh = dict_get_key(mesh_dict, uid)
+ id = t % mesh
+ if (dict_has_key(mesh_dict, id)) then
+ t % mesh = dict_get_key(mesh_dict, id)
else
- msg = "Could not find mesh " // trim(int_to_str(uid)) // &
- & " specified on tally " // trim(int_to_str(t % uid))
- call fatal_error(msg)
+ message = "Could not find mesh " // trim(int_to_str(id)) // &
+ & " specified on tally " // trim(int_to_str(t % id))
+ call fatal_error()
end if
end if
end do
@@ -534,10 +531,10 @@ contains
logical :: percent_in_atom ! nuclides specified in atom percent?
logical :: density_in_atom ! density specified in atom/b-cm?
character(10) :: key ! name of nuclide, e.g. 92235.03c
- character(MAX_LINE_LEN) :: msg ! output/error message
type(Material), pointer :: mat => null()
- ! first find the index in the xsdata array for each nuclide in each material
+ ! first find the index in the xs_listings array for each nuclide in each
+ ! material
do i = 1, n_materials
mat => materials(i)
@@ -545,9 +542,9 @@ contains
! given
if (.not. (all(mat%atom_percent > ZERO) .or. &
& all(mat%atom_percent < ZERO))) then
- msg = "Cannot mix atom and weight percents in material " // &
- & int_to_str(mat%uid)
- call fatal_error(msg)
+ message = "Cannot mix atom and weight percents in material " // &
+ & int_to_str(mat % id)
+ call fatal_error()
end if
percent_in_atom = (mat%atom_percent(1) > ZERO)
@@ -561,22 +558,22 @@ contains
! Check to make sure cross-section is continuous energy neutron table
n = len_trim(key)
if (key(n:n) /= 'c') then
- msg = "Cross-section table " // trim(key) // " is not a " // &
- "continuous-energy neutron table."
- call fatal_error(msg)
+ message = "Cross-section table " // trim(key) // &
+ " is not a continuous-energy neutron table."
+ call fatal_error()
end if
- if (dict_has_key(xsdata_dict, key)) then
- index = dict_get_key(xsdata_dict, key)
- mat % xsdata(j) = index
+ if (dict_has_key(xs_listing_dict, key)) then
+ index = dict_get_key(xs_listing_dict, key)
+ mat % xs_listing(j) = index
else
- msg = "Cannot find cross-section " // trim(key) // " in specified &
- &xsdata file."
- call fatal_error(msg)
+ message = "Cannot find cross-section " // trim(key) // &
+ " in specified cross_sections.xml file."
+ call fatal_error()
end if
! determine atomic weight ratio
- awr = xsdatas(index) % awr
+ awr = xs_listings(index) % awr
! if given weight percent, convert all values so that they are divided
! by awr. thus, when a sum is done over the values, it's actually
@@ -597,8 +594,8 @@ contains
if (.not. density_in_atom) then
sum_percent = ZERO
do j = 1, mat % n_nuclides
- index = mat % xsdata(j)
- awr = xsdatas(index) % awr
+ index = mat % xs_listing(j)
+ awr = xs_listings(index) % awr
x = mat % atom_percent(j)
sum_percent = sum_percent + x*awr
end do
diff --git a/src/input_xml.f90 b/src/input_xml.f90
index de2b252e18..ba881eb554 100644
--- a/src/input_xml.f90
+++ b/src/input_xml.f90
@@ -7,9 +7,9 @@ module input_xml
use geometry_header, only: Cell, Surface, Lattice
use global
use mesh_header, only: StructuredMesh
- use output, only: message
+ use output, only: write_message
use string, only: lower_case, int_to_str, str_to_int, str_to_real, &
- split_string
+ split_string, starts_with, ends_with
use tally_header, only: TallyObject
implicit none
@@ -46,27 +46,48 @@ contains
integer :: n
integer :: coeffs_reqd
logical :: file_exists
+ character(MAX_FILE_LEN) :: env_variable
character(MAX_WORD_LEN) :: type
- character(MAX_LINE_LEN) :: msg
character(MAX_LINE_LEN) :: filename
! Display output message
- msg = "Reading settings XML file..."
- call message(msg, 5)
+ message = "Reading settings XML file..."
+ call write_message(5)
! Check if settings.xml exists
filename = trim(path_input) // "settings.xml"
inquire(FILE=filename, EXIST=file_exists)
if (.not. file_exists) then
- msg = "Settings XML file '" // trim(filename) // "' does not exist!"
- call fatal_error(msg)
+ message = "Settings XML file '" // trim(filename) // "' does not exist!"
+ call fatal_error()
end if
+ ! Initialize XML scalar variables
+ cross_sections_ = ""
+ verbosity_ = 0
+
! Parse settings.xml file
call read_xml_file_settings_t(filename)
- ! Cross-section library path
- path_xsdata = trim(xslibrary % path)
+ ! Find cross_sections.xml file -- the first place to look is the
+ ! settings.xml file. If no file is found there, then we check the
+ ! CROSS_SECTIONS environment variable
+
+ if (len_trim(cross_sections_) == 0) then
+ ! No cross_sections.xml file specified in settings.xml, check environment
+ ! variable
+ call get_environment_variable("CROSS_SECTIONS", env_variable)
+ if (len_trim(env_variable) == 0) then
+ message = "No cross_sections.xml file was specified in " // &
+ "settings.xml or in the CROSS_SECTIONS environment " // &
+ "variable."
+ call fatal_error()
+ else
+ path_cross_sections = trim(env_variable)
+ end if
+ else
+ path_cross_sections = trim(cross_sections_)
+ end if
! Criticality information
if (criticality % cycles > 0) then
@@ -77,7 +98,7 @@ contains
end if
! Verbosity
- verbosity = verbosity_
+ if (verbosity_ > 0) verbosity = verbosity_
if (associated(source_ % coeffs)) then
! Determine external source type
@@ -88,19 +109,18 @@ contains
external_source % type = SRC_BOX
coeffs_reqd = 6
case default
- msg = "Invalid source type: " // trim(source_ % type)
- call fatal_error(msg)
+ message = "Invalid source type: " // trim(source_ % type)
+ call fatal_error()
end select
! Coefficients for external surface
n = size(source_ % coeffs)
if (n < coeffs_reqd) then
- msg = "Not enough coefficients specified for external source."
- print *, n, coeffs_reqd
- call fatal_error(msg)
+ message = "Not enough coefficients specified for external source."
+ call fatal_error()
elseif (n > coeffs_reqd) then
- msg = "Too many coefficients specified for external source."
- call fatal_error(msg)
+ message = "Too many coefficients specified for external source."
+ call fatal_error()
else
allocate(external_source % values(n))
external_source % values = source_ % coeffs
@@ -135,15 +155,14 @@ contains
integer :: coeffs_reqd
logical :: file_exists
character(MAX_LINE_LEN) :: filename
- character(MAX_LINE_LEN) :: msg
character(MAX_WORD_LEN) :: word
type(Cell), pointer :: c => null()
type(Surface), pointer :: s => null()
type(Lattice), pointer :: l => null()
! Display output message
- msg = "Reading geometry XML file..."
- call message(msg, 5)
+ message = "Reading geometry XML file..."
+ call write_message(5)
! ==========================================================================
! READ CELLS FROM GEOMETRY.XML
@@ -152,8 +171,8 @@ contains
filename = trim(path_input) // "geometry.xml"
inquire(FILE=filename, EXIST=file_exists)
if (.not. file_exists) then
- msg = "Geometry XML file '" // trim(filename) // "' does not exist!"
- call fatal_error(msg)
+ message = "Geometry XML file '" // trim(filename) // "' does not exist!"
+ call fatal_error()
end if
! Parse geometry.xml file
@@ -168,29 +187,30 @@ contains
c => cells(i)
! Copy data into cells
- c % uid = cell_(i) % uid
+ c % id = cell_(i) % id
c % universe = cell_(i) % universe
c % material = cell_(i) % material
c % fill = cell_(i) % fill
! Check to make sure that either material or fill was specified
if (c % material == 0 .and. c % fill == 0) then
- msg = "Neither material nor fill was specified for cell " // &
- trim(int_to_str(c % uid))
- call fatal_error(msg)
+ message = "Neither material nor fill was specified for cell " // &
+ trim(int_to_str(c % id))
+ call fatal_error()
end if
! Check to make sure that both material and fill haven't been
! specified simultaneously
if (c % material /= 0 .and. c % fill /= 0) then
- msg = "Cannot specify material and fill simultaneously"
- call fatal_error(msg)
+ message = "Cannot specify material and fill simultaneously"
+ call fatal_error()
end if
! Check to make sure that surfaces were specified
if (.not. associated(cell_(i) % surfaces)) then
- msg = "No surfaces specified for cell " // trim(int_to_str(c % uid))
- call fatal_error(msg)
+ message = "No surfaces specified for cell " // &
+ trim(int_to_str(c % id))
+ call fatal_error()
end if
! Allocate array for surfaces and copy
@@ -200,7 +220,7 @@ contains
c % surfaces = cell_(i) % surfaces
! Add cell to dictionary
- call dict_add_key(cell_dict, c % uid, i)
+ call dict_add_key(cell_dict, c % id, i)
! For cells, we also need to check if there's a new universe --
! also for every cell add 1 to the count of cells for the
@@ -228,7 +248,7 @@ contains
s => surfaces(i)
! Copy data into cells
- s % uid = surface_(i) % uid
+ s % id = surface_(i) % id
! Copy and interpret surface type
word = surface_(i) % type
@@ -274,8 +294,8 @@ contains
s % type = SURF_GQ
coeffs_reqd = 10
case default
- msg = "Invalid surface type: " // trim(surface_(i) % type)
- call fatal_error(msg)
+ message = "Invalid surface type: " // trim(surface_(i) % type)
+ call fatal_error()
end select
! Check to make sure that the proper number of coefficients
@@ -284,14 +304,14 @@ contains
n = size(surface_(i) % coeffs)
if (n < coeffs_reqd) then
- msg = "Not enough coefficients specified for surface: " // &
- trim(int_to_str(s % uid))
+ message = "Not enough coefficients specified for surface: " // &
+ trim(int_to_str(s % id))
print *, n, coeffs_reqd
- call fatal_error(msg)
+ call fatal_error()
elseif (n > coeffs_reqd) then
- msg = "Too many coefficients specified for surface: " // &
- trim(int_to_str(s % uid))
- call fatal_error(msg)
+ message = "Too many coefficients specified for surface: " // &
+ trim(int_to_str(s % id))
+ call fatal_error()
else
allocate(s % coeffs(n))
s % coeffs = surface_(i) % coeffs
@@ -310,13 +330,13 @@ contains
case ('periodic')
s % bc = BC_PERIODIC
case default
- msg = "Unknown boundary condition '" // trim(word) // &
- "' specified on surface " // trim(int_to_str(s % uid))
- call fatal_error(msg)
+ message = "Unknown boundary condition '" // trim(word) // &
+ "' specified on surface " // trim(int_to_str(s % id))
+ call fatal_error()
end select
! Add surface to dictionary
- call dict_add_key(surface_dict, s % uid, i)
+ call dict_add_key(surface_dict, s % id, i)
end do
@@ -330,8 +350,8 @@ contains
do i = 1, n_lattices
l => lattices(i)
- ! UID of lattice
- l % uid = lattice_(i) % uid
+ ! ID of lattice
+ l % id = lattice_(i) % id
! Read lattice type
word = lattice_(i) % type
@@ -342,15 +362,15 @@ contains
case ('hex', 'hexagon', 'hexagonal')
l % type = LATTICE_HEX
case default
- msg = "Invalid lattice type: " // trim(lattice_(i) % type)
- call fatal_error(msg)
+ message = "Invalid lattice type: " // trim(lattice_(i) % type)
+ call fatal_error()
end select
! Read number of lattice cells in each dimension
n = size(lattice_(i) % dimension)
if (n /= 2 .and. n /= 3) then
- msg = "Lattice must be two or three dimensions."
- call fatal_error(msg)
+ message = "Lattice must be two or three dimensions."
+ call fatal_error()
end if
n_x = lattice_(i) % dimension(1)
n_y = lattice_(i) % dimension(2)
@@ -359,18 +379,18 @@ contains
! Read lattice origin location
if (size(lattice_(i) % dimension) /= size(lattice_(i) % origin)) then
- msg = "Number of entries on must be the same as the " // &
- "number of entries on ."
- call fatal_error(msg)
+ message = "Number of entries on must be the same as " // &
+ "the number of entries on ."
+ call fatal_error()
end if
l % x0 = lattice_(i) % origin(1)
l % y0 = lattice_(i) % origin(2)
! Read lattice widths
if (size(lattice_(i) % width) /= size(lattice_(i) % origin)) then
- msg = "Number of entries on must be the same as the " // &
- "number of entries on ."
- call fatal_error(msg)
+ message = "Number of entries on must be the same as " // &
+ "the number of entries on ."
+ call fatal_error()
end if
l % width_x = lattice_(i) % width(1)
l % width_y = lattice_(i) % width(2)
@@ -384,7 +404,7 @@ contains
end do
! Add lattice to dictionary
- call dict_add_key(lattice_dict, l % uid, i)
+ call dict_add_key(lattice_dict, l % id, i)
end do
@@ -406,21 +426,20 @@ contains
character(MAX_WORD_LEN) :: units
character(MAX_WORD_LEN) :: name
character(MAX_LINE_LEN) :: filename
- character(MAX_LINE_LEN) :: msg
type(Material), pointer :: m => null()
type(nuclide_xml), pointer :: nuc => null()
type(sab_xml), pointer :: sab => null()
! Display output message
- msg = "Reading materials XML file..."
- call message(msg, 5)
+ message = "Reading materials XML file..."
+ call write_message(5)
! Check is materials.xml exists
filename = trim(path_input) // "materials.xml"
inquire(FILE=filename, EXIST=file_exists)
if (.not. file_exists) then
- msg = "Material XML file '" // trim(filename) // "' does not exist!"
- call fatal_error(msg)
+ message = "Material XML file '" // trim(filename) // "' does not exist!"
+ call fatal_error()
end if
! Parse materials.xml file
@@ -433,8 +452,8 @@ contains
do i = 1, n_materials
m => materials(i)
- ! Copy material uid
- m % uid = material_(i) % uid
+ ! Copy material id
+ m % id = material_(i) % id
! Copy density -- the default value for the units is given in the
! material_t.xml file and doesn't need to be specified here, hence case
@@ -452,15 +471,16 @@ contains
case ('atom/cm3', 'atom/cc')
m % density = 1.0e-24 * val
case default
- msg = "Unkwown units '" // trim(material_(i) % density % units) // &
- "' specified on material " // trim(int_to_str(m % uid))
- call fatal_error(msg)
+ message = "Unkwown units '" // trim(material_(i) % density % units) &
+ // "' specified on material " // trim(int_to_str(m % id))
+ call fatal_error()
end select
! Check to ensure material has at least one nuclide
if (.not. associated(material_(i) % nuclides)) then
- msg = "No nuclides specified on material " // trim(int_to_str(m % uid))
- call fatal_error(msg)
+ message = "No nuclides specified on material " // &
+ trim(int_to_str(m % id))
+ call fatal_error()
end if
! allocate arrays in Material object
@@ -468,7 +488,7 @@ contains
m % n_nuclides = n
allocate(m % names(n))
allocate(m % nuclide(n))
- allocate(m % xsdata(n))
+ allocate(m % xs_listing(n))
allocate(m % atom_density(n))
allocate(m % atom_percent(n))
@@ -481,13 +501,13 @@ contains
! Check if no atom/weight percents were specified or if both atom and
! weight percents were specified
if (nuc % ao == ZERO .and. nuc % wo == ZERO) then
- msg = "No atom or weight percent specified for nuclide " // &
+ message = "No atom or weight percent specified for nuclide " // &
trim(name)
- call fatal_error(msg)
+ call fatal_error()
elseif (nuc % ao /= ZERO .and. nuc % wo /= ZERO) then
- msg = "Cannot specify both atom and weight percents for a nuclide: " &
- // trim(name)
- call fatal_error(msg)
+ message = "Cannot specify both atom and weight percents for a " &
+ // "nuclide: " // trim(name)
+ call fatal_error()
end if
! Copy atom/weight percents
@@ -504,13 +524,13 @@ contains
m % sab_name = name
m % has_sab_table = .true.
elseif (size(material_(i) % sab) > 1) then
- msg = "Cannot have multiple S(a,b) tables on a single material."
- call fatal_error(msg)
+ message = "Cannot have multiple S(a,b) tables on a single material."
+ call fatal_error()
end if
end do
! Add material to dictionary
- call dict_add_key(material_dict, m % uid, i)
+ call dict_add_key(material_dict, m % id, i)
end do
@@ -527,13 +547,12 @@ contains
integer :: i ! loop over user-specified tallies
integer :: j ! loop over words
- integer :: uid ! user-specified identifier
+ integer :: id ! user-specified identifier
integer :: index ! index in meshes array
integer :: n ! size of arrays in mesh specification
integer :: n_words ! number of words read
logical :: file_exists ! does tallies.xml file exist?
character(MAX_LINE_LEN) :: filename
- character(MAX_LINE_LEN) :: msg
character(MAX_WORD_LEN) :: word
character(MAX_WORD_LEN) :: words(MAX_WORDS)
type(TallyObject), pointer :: t => null()
@@ -549,8 +568,8 @@ contains
end if
! Display output message
- msg = "Reading tallies XML file..."
- call message(msg, 5)
+ message = "Reading tallies XML file..."
+ call write_message(5)
! Parse tallies.xml file
call read_xml_file_tallies_t(filename)
@@ -569,8 +588,8 @@ contains
! Allocate tallies array
if (.not. associated(tally_)) then
n_tallies = 0
- msg = "No tallies present in tallies.xml file!"
- call warning(msg)
+ message = "No tallies present in tallies.xml file!"
+ call warning()
else
n_tallies = size(tally_)
allocate(tallies(n_tallies))
@@ -582,8 +601,8 @@ contains
do i = 1, n_meshes
m => meshes(i)
- ! copy mesh uid
- m % uid = mesh_(i) % id
+ ! copy mesh id
+ m % id = mesh_(i) % id
! Read mesh type
word = mesh_(i) % type
@@ -594,15 +613,15 @@ contains
case ('hex', 'hexagon', 'hexagonal')
m % type = LATTICE_HEX
case default
- msg = "Invalid mesh type: " // trim(mesh_(i) % type)
- call fatal_error(msg)
+ message = "Invalid mesh type: " // trim(mesh_(i) % type)
+ call fatal_error()
end select
! Determine number of dimensions for mesh
n = size(mesh_(i) % dimension)
if (n /= 2 .and. n /= 3) then
- msg = "Mesh must be two or three dimensions."
- call fatal_error(msg)
+ message = "Mesh must be two or three dimensions."
+ call fatal_error()
end if
m % n_dimension = n
@@ -616,22 +635,22 @@ contains
! Read mesh origin location
if (m % n_dimension /= size(mesh_(i) % origin)) then
- msg = "Number of entries on must be the same as the " // &
- "number of entries on ."
- call fatal_error(msg)
+ message = "Number of entries on must be the same as " // &
+ "the number of entries on ."
+ call fatal_error()
end if
m % origin = mesh_(i) % origin
! Read mesh widths
if (size(mesh_(i) % width) /= size(mesh_(i) % origin)) then
- msg = "Number of entries on must be the same as the " // &
- "number of entries on ."
- call fatal_error(msg)
+ message = "Number of entries on must be the same as " // &
+ "the number of entries on ."
+ call fatal_error()
end if
m % width = mesh_(i) % width
! Add mesh to dictionary
- call dict_add_key(mesh_dict, m % uid, i)
+ call dict_add_key(mesh_dict, m % id, i)
end do
! ==========================================================================
@@ -648,15 +667,15 @@ contains
t % n_bins = 0
t % stride = 0
- ! Copy material uid
- t % uid = tally_(i) % id
+ ! Copy material id
+ t % id = tally_(i) % id
! Check to make sure that both cells and surfaces were not specified
if (len_trim(tally_(i) % filters % cell) > 0 .and. &
len_trim(tally_(i) % filters % surface) > 0) then
- msg = "Cannot specify both cell and surface filters for tally " &
- // trim(int_to_str(t % uid))
- call fatal_error(msg)
+ message = "Cannot specify both cell and surface filters for tally " &
+ // trim(int_to_str(t % id))
+ call fatal_error()
end if
! TODO: Parse logical expressions instead of just each word
@@ -705,14 +724,14 @@ contains
t % mesh = tally_(i) % filters % mesh
if (t % mesh > 0) then
! Determine index in mesh array for this bin
- uid = t % mesh
- if (dict_has_key(mesh_dict, uid)) then
- index = dict_get_key(mesh_dict, uid)
+ id = t % mesh
+ if (dict_has_key(mesh_dict, id)) then
+ index = dict_get_key(mesh_dict, id)
m => meshes(index)
else
- msg = "Could not find mesh " // trim(int_to_str(uid)) // &
- " specified on tally " // trim(int_to_str(t % uid))
- call fatal_error(msg)
+ message = "Could not find mesh " // trim(int_to_str(id)) // &
+ " specified on tally " // trim(int_to_str(t % id))
+ call fatal_error()
end if
t % n_bins(T_MESH) = t % n_bins(T_MESH) + product(m % dimension)
@@ -759,15 +778,15 @@ contains
case ('flux')
t % macro_bins(j) % scalar = MACRO_FLUX
if (t % n_bins(T_ENERGYOUT) > 0) then
- msg = "Cannot tally flux with an outgoing energy filter."
- call fatal_error(msg)
+ message = "Cannot tally flux with an outgoing energy filter."
+ call fatal_error()
end if
case ('total')
t % macro_bins(j) % scalar = MACRO_TOTAL
if (t % n_bins(T_ENERGYOUT) > 0) then
- msg = "Cannot tally total reaction rate with an outgoing " &
- // "energy filter."
- call fatal_error(msg)
+ message = "Cannot tally total reaction rate with an " &
+ // "outgoing energy filter."
+ call fatal_error()
end if
case ('scatter')
t % macro_bins(j) % scalar = MACRO_SCATTER
@@ -790,16 +809,16 @@ contains
case ('absorption')
t % macro_bins(j) % scalar = MACRO_ABSORPTION
if (t % n_bins(T_ENERGYOUT) > 0) then
- msg = "Cannot tally absorption rate with an outgoing " &
+ message = "Cannot tally absorption rate with an outgoing " &
// "energy filter."
- call fatal_error(msg)
+ call fatal_error()
end if
case ('fission')
t % macro_bins(j) % scalar = MACRO_FISSION
if (t % n_bins(T_ENERGYOUT) > 0) then
- msg = "Cannot tally fission rate with an outgoing " &
+ message = "Cannot tally fission rate with an outgoing " &
// "energy filter."
- call fatal_error(msg)
+ call fatal_error()
end if
case ('nu-fission')
t % macro_bins(j) % scalar = MACRO_NU_FISSION
@@ -810,17 +829,17 @@ contains
! Check to make sure that current is the only desired response
! for this tally
if (n_words > 1) then
- msg = "Cannot tally other macro reactions in the same " &
+ message = "Cannot tally other macro reactions in the same " &
// "tally as surface currents. Separate other macro " &
// "reactions into a distinct tally."
- call fatal_error(msg)
+ call fatal_error()
end if
! Check to make sure that only the mesh filter was specified
if (t % mesh == 0 .or. t % n_bins(T_MESH) /= &
product(t % n_bins, t % n_bins > 0)) then
- msg = "Surface currents must be used with a mesh filter only."
- call fatal_error(msg)
+ message = "Surface currents must be used with a mesh filter only."
+ call fatal_error()
end if
! Since the number of bins for the mesh filter was already set
@@ -829,8 +848,8 @@ contains
t % n_bins(T_MESH) = t % n_bins(T_MESH) - product(m % dimension)
! Get pointer to mesh
- uid = t % mesh
- index = dict_get_key(mesh_dict, uid)
+ id = t % mesh
+ index = dict_get_key(mesh_dict, id)
m => meshes(index)
! We need to increase the dimension by one since we also need
@@ -844,8 +863,8 @@ contains
end if
case default
- msg = "Unknown macro reaction: " // trim(words(j))
- call fatal_error(msg)
+ message = "Unknown macro reaction: " // trim(words(j))
+ call fatal_error()
end select
end do
t % n_macro_bins = n_words
@@ -866,19 +885,18 @@ contains
logical :: file_exists ! does tallies.xml file exist?
character(MAX_LINE_LEN) :: filename
- character(MAX_LINE_LEN) :: msg
! Check if plot.xml exists
filename = trim(path_input) // "plot.xml"
inquire(FILE=filename, EXIST=file_exists)
if (.not. file_exists) then
- msg = "Plot XML file '" // trim(filename) // "' does not exist!"
- call fatal_error(msg)
+ message = "Plot XML file '" // trim(filename) // "' does not exist!"
+ call fatal_error()
end if
! Display output message
- msg = "Reading plot XML file..."
- call message(msg, 5)
+ message = "Reading plot XML file..."
+ call write_message(5)
! Parse plot.xml file
call read_xml_file_plot_t(filename)
@@ -908,4 +926,102 @@ contains
end subroutine read_plot_xml
+!===============================================================================
+! READ_CROSS_SECTIONS_XML reads information from a cross_sections.xml file. This
+! file contains a listing of the ACE cross sections that may be used.
+!===============================================================================
+
+ subroutine read_cross_sections_xml()
+
+ use xml_data_cross_sections_t
+
+ integer :: i ! loop index
+ integer :: n_listings ! number of listings in cross_sections.xml
+ logical :: file_exists ! does cross_sections.xml exist?
+ character(MAX_WORD_LEN) :: directory
+ type(XsListing), pointer :: listing => null()
+
+ ! Check if cross_sections.xml exists
+ inquire(FILE=path_cross_sections, EXIST=file_exists)
+ if (.not. file_exists) then
+ ! Could not find cross_sections.xml file
+ message = "Cross sections XML file '" // trim(path_cross_sections) // &
+ "' does not exist!"
+ call fatal_error()
+ end if
+
+ message = "Reading cross sections XML file..."
+ call write_message(5)
+
+ ! Initialize directory_ variable
+ directory_ = ""
+
+ ! Parse cross_sections.xml file
+ call read_xml_file_cross_sections_t(path_cross_sections)
+
+ ! Copy directory information if present
+ directory = trim(directory_)
+
+ ! Allocate xs_listings array
+ if (.not. associated(ace_tables_)) then
+ message = "No ACE table listings present in cross_sections.xml file!"
+ call fatal_error()
+ else
+ n_listings = size(ace_tables_)
+ allocate(xs_listings(n_listings))
+ end if
+
+
+ do i = 1, n_listings
+ listing => xs_listings(i)
+
+ ! copy a number of attributes
+ listing % name = trim(ace_tables_(i) % name)
+ listing % alias = trim(ace_tables_(i) % alias)
+ listing % zaid = ace_tables_(i) % zaid
+ listing % awr = ace_tables_(i) % awr
+ listing % temp = ace_tables_(i) % temperature
+
+ ! determine type of cross section
+ select case(ace_tables_(i) % type)
+ case ('neutron')
+ listing % type = ACE_NEUTRON
+ case ('thermal')
+ listing % type = ACE_THERMAL
+ case ('dosimetry')
+ listing % type = ACE_DOSIMETRY
+ end select
+
+ ! determine metastable state
+ if (ace_tables_(i) % metastable == 0) then
+ listing % metastable = .false.
+ else
+ listing % metastable = .true.
+ end if
+
+ ! determine whether binary/ascii
+ if (ace_tables_(i) % binary == 0) then
+ listing % binary = .false.
+ else
+ listing % binary = .true.
+ end if
+
+ ! determine path of cross section table
+ if (starts_with(ace_tables_(i) % path, '/')) then
+ listing % path = ace_tables_(i) % path
+ else
+ if (ends_with(directory,'/')) then
+ listing % path = trim(directory) // trim(ace_tables_(i) % path)
+ else
+ listing % path = trim(directory) // '/' // trim(ace_tables_(i) % path)
+ end if
+ end if
+
+ ! create dictionary entry for both name and alias
+ call dict_add_key(xs_listing_dict, listing % name, i)
+ call dict_add_key(xs_listing_dict, listing % alias, i)
+ end do
+
+ end subroutine read_cross_sections_xml
+
end module input_xml
diff --git a/src/interpolation.f90 b/src/interpolation.f90
index c0e6f631b6..231007e738 100644
--- a/src/interpolation.f90
+++ b/src/interpolation.f90
@@ -2,6 +2,7 @@ module interpolation
use constants
use error, only: fatal_error
+ use global, only: message
use search, only: binary_search
implicit none
@@ -34,7 +35,6 @@ contains
real(8) :: r ! interpolation factor
real(8) :: x0, x1 ! bounding x values
real(8) :: y0, y1 ! bounding y values
- character(MAX_LINE_LEN) :: msg
! determine starting location
if (present(loc_start)) then
@@ -76,8 +76,8 @@ contains
elseif (n_regions == 1) then
interp = data(loc_interp + 1)
elseif (n_regions > 1) then
- msg = "Multiple interpolation regions not yet supported."
- call fatal_error(msg)
+ message = "Multiple interpolation regions not yet supported."
+ call fatal_error()
end if
! handle special case of histogram interpolation
diff --git a/src/logging.f90 b/src/logging.f90
index 5c8279e228..377c1ab827 100644
--- a/src/logging.f90
+++ b/src/logging.f90
@@ -15,7 +15,7 @@ contains
character(MAX_WORD_LEN) :: path_log ! path of log file
logical :: file_exists ! does log file already exist?
- integer :: ioError ! error status for file access
+ ! integer :: ioError ! error status for file access
! Create filename for log file
path_log = trim(path_input) // ".log"
@@ -27,8 +27,8 @@ contains
end if
! Open log file for writing
- open(FILE=path_log, UNIT=UNIT_LOG, STATUS='replace', &
- & ACTION='write', IOSTAT=ioError)
+ ! open(FILE=path_log, UNIT=UNIT_LOG, STATUS='replace', &
+ ! & ACTION='write', IOSTAT=ioError)
end subroutine create_log
diff --git a/src/main.f90 b/src/main.f90
index d6caa6d0ef..320850efb6 100644
--- a/src/main.f90
+++ b/src/main.f90
@@ -6,7 +6,7 @@ program main
use initialize, only: initialize_run
use mcnp_random, only: RN_init_particle
use mpi_routines, only: synchronize_bank
- use output, only: message, header, print_runtime
+ use output, only: write_message, header, print_runtime
use particle_header, only: Particle
use plot, only: run_plot
use physics, only: transport
@@ -71,7 +71,7 @@ contains
call timer_start(time_compute)
msg = "Simulating cycle " // trim(int_to_str(i_cycle)) // "..."
- call message(msg, 8)
+ call write_message(8)
! Set all tallies to zero
n_bank = 0
@@ -88,7 +88,7 @@ contains
end if
! set random number seed
- i_particle = (i_cycle-1)*n_particles + p % uid
+ i_particle = (i_cycle-1)*n_particles + p % id
call RN_init_particle(i_particle)
! transport particle
diff --git a/src/material_header.f90 b/src/material_header.f90
index 5a4899653e..ac77394d76 100644
--- a/src/material_header.f90
+++ b/src/material_header.f90
@@ -7,15 +7,14 @@ module material_header
!===============================================================================
type Material
- integer :: uid ! unique identifier
+ integer :: id ! unique identifier
integer :: n_nuclides ! number of nuclides
character(10), allocatable :: names(:) ! isotope names
- integer, allocatable :: xsdata(:) ! index in xsdata list
+ integer, allocatable :: xs_listing(:) ! index in xs_listings list
integer, allocatable :: nuclide(:) ! index in nuclides array
real(8) :: density ! total atom density in atom/b-cm
real(8), allocatable :: atom_density(:) ! nuclide atom density in atom/b-cm
real(8), allocatable :: atom_percent(:) ! atom/weight percent (negative for weight)
- real(8), allocatable :: total_xs(:) ! macroscopic cross-section
! S(a,b) data references
logical :: has_sab_table = .false.
diff --git a/src/mesh_header.f90 b/src/mesh_header.f90
index f49a5c50ab..9ed686ec7e 100644
--- a/src/mesh_header.f90
+++ b/src/mesh_header.f90
@@ -8,7 +8,7 @@ module mesh_header
!===============================================================================
type StructuredMesh
- integer :: uid
+ integer :: id
integer :: type
integer :: n_dimension
integer, allocatable :: dimension(:)
diff --git a/src/mpi_routines.f90 b/src/mpi_routines.f90
index b19b3aa665..dc13f04f7f 100644
--- a/src/mpi_routines.f90
+++ b/src/mpi_routines.f90
@@ -4,7 +4,7 @@ module mpi_routines
use error, only: fatal_error
use global
use mcnp_random, only: rang, RN_init_particle, RN_skip
- use output, only: message
+ use output, only: write_message
use particle_header, only: Particle, initialize_particle
use tally_header, only: TallyObject
@@ -35,7 +35,6 @@ contains
integer :: bank_types(4) ! Datatypes
integer(MPI_ADDRESS_KIND) :: bank_disp(4) ! Displacements
integer(MPI_ADDRESS_KIND) :: base
- character(MAX_LINE_LEN) :: msg ! Error message
type(Bank) :: b
mpi_enabled = .true.
@@ -43,22 +42,22 @@ contains
! Initialize MPI
call MPI_INIT(ierr)
if (ierr /= MPI_SUCCESS) then
- msg = "Failed to initialize MPI."
- call fatal_error(msg)
+ message = "Failed to initialize MPI."
+ call fatal_error()
end if
! Determine number of processors
call MPI_COMM_SIZE(MPI_COMM_WORLD, n_procs, ierr)
if (ierr /= MPI_SUCCESS) then
- msg = "Could not determine number of processors."
- call fatal_error(msg)
+ message = "Could not determine number of processors."
+ call fatal_error()
end if
! Determine rank of each processor
call MPI_COMM_RANK(MPI_COMM_WORLD, rank, ierr)
if (ierr /= MPI_SUCCESS) then
- msg = "Could not determine MPI rank."
- call fatal_error(msg)
+ message = "Could not determine MPI rank."
+ call fatal_error()
end if
! Determine master
@@ -69,7 +68,7 @@ contains
end if
! Determine displacements for MPI_BANK type
- call MPI_GET_ADDRESS(b % uid, bank_disp(1), ierr)
+ call MPI_GET_ADDRESS(b % id, bank_disp(1), ierr)
call MPI_GET_ADDRESS(b % xyz, bank_disp(2), ierr)
call MPI_GET_ADDRESS(b % uvw, bank_disp(3), ierr)
call MPI_GET_ADDRESS(b % E, bank_disp(4), ierr)
@@ -113,7 +112,7 @@ contains
integer(8) :: finish ! ending index in local fission bank
integer(8) :: total ! total sites in global fission bank
integer(8) :: count ! index for source bank
- integer(8) :: index ! index for uid -- accounts for all nodes
+ integer(8) :: index ! index for id -- accounts for all nodes
integer :: send_to_left ! # of bank sites to send/recv to or from left
integer :: send_to_right ! # of bank sites to send/recv to or from right
integer(8) :: sites_needed ! # of sites to be sampled
@@ -122,7 +121,6 @@ contains
& temp_sites(:), & ! local array of extra sites on each node
& left_bank(:), & ! bank sites to send/recv to or from left node
& right_bank(:) ! bank sites to send/recv to or fram right node
- character(MAX_LINE_LEN) :: msg
#ifdef MPI
integer :: ierr
@@ -133,8 +131,8 @@ contains
real(8) :: t0, t1, t2, t3, t4
#endif
- msg = "Collecting number of fission sites..."
- call message(msg, 8)
+ message = "Collecting number of fission sites..."
+ call write_message(8)
#ifdef MPI
call MPI_BARRIER(MPI_COMM_WORLD, ierr)
@@ -159,8 +157,8 @@ contains
! Check if there are no fission sites
if (total == 0) then
- msg = "No fission sites banked!"
- call fatal_error(msg)
+ message = "No fission sites banked!"
+ call fatal_error()
end if
! Make sure all processors start at the same point for random sampling
@@ -171,7 +169,7 @@ contains
allocate(temp_sites(2*work))
count = 0_8 ! Index for local source_bank
- index = 0_8 ! Index for global source uid -- must account for all nodes
+ index = 0_8 ! Index for global source id -- must account for all nodes
if (total < n_particles) then
sites_needed = mod(n_particles,total)
@@ -180,8 +178,8 @@ contains
end if
p_sample = real(sites_needed,8)/real(total,8)
- msg = "Sampling fission sites..."
- call message(msg, 8)
+ message = "Sampling fission sites..."
+ call write_message(8)
! ==========================================================================
! SAMPLE N_PARTICLES FROM FISSION BANK AND PLACE IN TEMP_SITES
@@ -252,8 +250,8 @@ contains
t2 = MPI_WTIME()
t_sync(2) = t_sync(2) + (t2 - t1)
- msg = "Sending fission sites..."
- call message(msg, 8)
+ message = "Sending fission sites..."
+ call write_message(8)
! ==========================================================================
! SEND BANK SITES TO NEIGHBORS
@@ -281,8 +279,8 @@ contains
t_sync(3) = t_sync(3) + (t3 - t2)
#endif
- msg = "Constructing source bank..."
- call message(msg, 8)
+ message = "Constructing source bank..."
+ call write_message(8)
! ==========================================================================
! RECONSTRUCT SOURCE BANK
diff --git a/src/output.f90 b/src/output.f90
index 527463effe..4e7f0ccb0e 100644
--- a/src/output.f90
+++ b/src/output.f90
@@ -152,14 +152,13 @@ contains
end subroutine header
!===============================================================================
-! MESSAGE displays an informational message to the log file and the standard
-! output stream.
+! WRITE_MESSAGE displays an informational message to the log file and the
+! standard output stream.
!===============================================================================
- subroutine message(msg, level)
+ subroutine write_message(level)
- character(*), intent(in) :: msg
- integer, optional :: level
+ integer, optional :: level
integer :: n_lines
integer :: i
@@ -168,13 +167,13 @@ contains
if (.not. master .and. present(level)) return
if (.not. present(level) .or. level <= verbosity) then
- n_lines = (len_trim(msg)-1)/79 + 1
+ n_lines = (len_trim(message)-1)/79 + 1
do i = 1, n_lines
- write(ou, fmt='(1X,A)') trim(msg(79*(i-1)+1:79*i))
+ write(ou, fmt='(1X,A)') trim(message(79*(i-1)+1:79*i))
end do
end if
- end subroutine message
+ end subroutine write_message
!===============================================================================
! GET_TODAY determines the date and time at which the program began execution
@@ -232,13 +231,13 @@ contains
select case (p % type)
case (NEUTRON)
- write(ou,*) 'Neutron ' // int_to_str(p % uid)
+ write(ou,*) 'Neutron ' // int_to_str(p % id)
case (PHOTON)
- write(ou,*) 'Photon ' // int_to_str(p % uid)
+ write(ou,*) 'Photon ' // int_to_str(p % id)
case (ELECTRON)
- write(ou,*) 'Electron ' // int_to_str(p % uid)
+ write(ou,*) 'Electron ' // int_to_str(p % id)
case default
- write(ou,*) 'Unknown Particle ' // int_to_str(p % uid)
+ write(ou,*) 'Unknown Particle ' // int_to_str(p % id)
end select
write(ou,*) ' x = ' // real_to_str(p % xyz(1))
write(ou,*) ' y = ' // real_to_str(p % xyz(2))
@@ -258,20 +257,20 @@ contains
if (p % cell > 0) then
c => cells(p % cell)
- write(ou,*) ' Cell = ' // int_to_str(c % uid)
+ write(ou,*) ' Cell = ' // int_to_str(c % id)
else
write(ou,*) ' Cell not determined'
end if
if (p % surface > 0) then
s => surfaces(p % surface)
- write(ou,*) ' Surface = ' // int_to_str(s % uid)
+ write(ou,*) ' Surface = ' // int_to_str(s % id)
else
write(ou,*) ' Surface = None'
end if
u => universes(p % universe)
- write(ou,*) ' Universe = ' // int_to_str(u % uid)
+ write(ou,*) ' Universe = ' // int_to_str(u % id)
write(ou,*)
nullify(c)
@@ -315,26 +314,26 @@ contains
type(Lattice), pointer :: l => null()
type(Material), pointer :: m => null()
- write(ou,*) 'Cell ' // int_to_str(c % uid)
- temp = dict_get_key(cell_dict, c % uid)
+ write(ou,*) 'Cell ' // int_to_str(c % id)
+ temp = dict_get_key(cell_dict, c % id)
write(ou,*) ' Array Index = ' // int_to_str(temp)
u => universes(c % universe)
- write(ou,*) ' Universe = ' // int_to_str(u % uid)
+ write(ou,*) ' Universe = ' // int_to_str(u % id)
select case (c % type)
case (CELL_NORMAL)
write(ou,*) ' Fill = NONE'
case (CELL_FILL)
u => universes(c % fill)
- write(ou,*) ' Fill = Universe ' // int_to_str(u % uid)
+ write(ou,*) ' Fill = Universe ' // int_to_str(u % id)
case (CELL_LATTICE)
l => lattices(c % fill)
- write(ou,*) ' Fill = Lattice ' // int_to_str(l % uid)
+ write(ou,*) ' Fill = Lattice ' // int_to_str(l % id)
end select
if (c % material == 0) then
write(ou,*) ' Material = NONE'
else
m => materials(c % material)
- write(ou,*) ' Material = ' // int_to_str(m % uid)
+ write(ou,*) ' Material = ' // int_to_str(m % id)
end if
write(ou,*) ' Parent Cell = ' // int_to_str(c % parent)
string = ""
@@ -373,12 +372,12 @@ contains
character(MAX_LINE_LEN) :: string
type(Cell), pointer :: c => null()
- write(ou,*) 'Universe ' // int_to_str(univ % uid)
+ write(ou,*) 'Universe ' // int_to_str(univ % id)
write(ou,*) ' Level = ' // int_to_str(univ % level)
string = ""
do i = 1, univ % n_cells
c => cells(univ % cells(i))
- string = trim(string) // ' ' // int_to_str(c % uid)
+ string = trim(string) // ' ' // int_to_str(c % id)
end do
write(ou,*) ' Cells =' // trim(string)
write(ou,*)
@@ -395,7 +394,7 @@ contains
type(Lattice), pointer :: lat
- write(ou,*) 'Lattice ' // int_to_str(lat % uid)
+ write(ou,*) 'Lattice ' // int_to_str(lat % id)
write(ou,*) ' n_x = ' // int_to_str(lat % n_x)
write(ou,*) ' n_y = ' // int_to_str(lat % n_y)
write(ou,*) ' x0 = ' // real_to_str(lat % x0)
@@ -417,7 +416,7 @@ contains
integer :: i
character(MAX_LINE_LEN) :: string
- write(ou,*) 'Surface ' // int_to_str(surf % uid)
+ write(ou,*) 'Surface ' // int_to_str(surf % id)
select case (surf % type)
case (SURF_PX)
string = "X Plane"
@@ -493,7 +492,7 @@ contains
type(Nuclide), pointer :: nuc => null()
! Write identifier for material
- write(ou,*) 'Material ' // int_to_str(mat % uid)
+ write(ou,*) 'Material ' // int_to_str(mat % id)
! Write total atom density in atom/b-cm
write(ou,*) ' Atom Density = ' // trim(real_to_str(mat % density)) &
@@ -528,7 +527,7 @@ contains
type(TallyObject), pointer :: t
integer :: i
- integer :: uid
+ integer :: id
character(MAX_LINE_LEN) :: string
type(Cell), pointer :: c => null()
type(Surface), pointer :: s => null()
@@ -536,14 +535,14 @@ contains
type(Material), pointer :: m => null()
type(StructuredMesh), pointer :: sm => null()
- write(ou,*) 'Tally ' // int_to_str(t % uid)
+ write(ou,*) 'Tally ' // int_to_str(t % id)
if (t % n_bins(T_CELL) > 0) then
string = ""
do i = 1, t % n_bins(T_CELL)
- uid = t % cell_bins(i) % scalar
- c => cells(uid)
- string = trim(string) // ' ' // trim(int_to_str(c % uid))
+ id = t % cell_bins(i) % scalar
+ c => cells(id)
+ string = trim(string) // ' ' // trim(int_to_str(c % id))
end do
write(ou, *) ' Cell Bins:' // trim(string)
end if
@@ -551,9 +550,9 @@ contains
if (t % n_bins(T_SURFACE) > 0) then
string = ""
do i = 1, t % n_bins(T_SURFACE)
- uid = t % surface_bins(i) % scalar
- s => surfaces(uid)
- string = trim(string) // ' ' // trim(int_to_str(s % uid))
+ id = t % surface_bins(i) % scalar
+ s => surfaces(id)
+ string = trim(string) // ' ' // trim(int_to_str(s % id))
end do
write(ou, *) ' Surface Bins:' // trim(string)
end if
@@ -561,9 +560,9 @@ contains
if (t % n_bins(T_UNIVERSE) > 0) then
string = ""
do i = 1, t % n_bins(T_UNIVERSE)
- uid = t % universe_bins(i) % scalar
- u => universes(uid)
- string = trim(string) // ' ' // trim(int_to_str(u % uid))
+ id = t % universe_bins(i) % scalar
+ u => universes(id)
+ string = trim(string) // ' ' // trim(int_to_str(u % id))
end do
write(ou, *) ' Material Bins:' // trim(string)
end if
@@ -571,17 +570,17 @@ contains
if (t % n_bins(T_MATERIAL) > 0) then
string = ""
do i = 1, t % n_bins(T_MATERIAL)
- uid = t % material_bins(i) % scalar
- m => materials(uid)
- string = trim(string) // ' ' // trim(int_to_str(m % uid))
+ id = t % material_bins(i) % scalar
+ m => materials(id)
+ string = trim(string) // ' ' // trim(int_to_str(m % id))
end do
write(ou, *) ' Material Bins:' // trim(string)
end if
if (t % n_bins(T_MESH) > 0) then
string = ""
- uid = t % mesh
- sm => meshes(uid)
+ id = t % mesh
+ sm => meshes(id)
string = trim(string) // ' ' // trim(int_to_str(sm % dimension(1)))
do i = 2, sm % n_dimension
string = trim(string) // ' x ' // trim(int_to_str(sm % dimension(i)))
@@ -592,9 +591,9 @@ contains
if (t % n_bins(T_CELLBORN) > 0) then
string = ""
do i = 1, t % n_bins(T_CELLBORN)
- uid = t % cellborn_bins(i) % scalar
- c => cells(uid)
- string = trim(string) // ' ' // trim(int_to_str(c % uid))
+ id = t % cellborn_bins(i) % scalar
+ c => cells(id)
+ string = trim(string) // ' ' // trim(int_to_str(c % id))
end do
write(ou, *) ' Birth Region Bins:' // trim(string)
end if
diff --git a/src/particle_header.f90 b/src/particle_header.f90
index 735c4cf1b7..deec93b065 100644
--- a/src/particle_header.f90
+++ b/src/particle_header.f90
@@ -11,7 +11,7 @@ module particle_header
type Particle
! Basic data
- integer(8) :: uid ! Unique ID
+ integer(8) :: id ! Unique ID
integer :: type ! Particle type (n, p, e, etc)
! Physical data
@@ -68,6 +68,7 @@ contains
p % type = NEUTRON
p % wgt = ONE
+ p % last_wgt = ONE
p % alive = .true.
p % n_bank = 0
p % cell = 0
diff --git a/src/physics.f90 b/src/physics.f90
index f47f6d8f09..4beecb2cf2 100644
--- a/src/physics.f90
+++ b/src/physics.f90
@@ -11,9 +11,9 @@ module physics
use global
use interpolation, only: interpolate_tab1
use mcnp_random, only: rang
- use output, only: message, print_particle
+ use output, only: write_message, print_particle
use particle_header, only: Particle
- use tally, only: score_tally
+ use tally, only: score_tally, score_surface_current
use search, only: binary_search
use string, only: int_to_str
@@ -36,7 +36,6 @@ contains
real(8) :: distance ! distance particle travels
logical :: found_cell ! found cell which particle is in?
logical :: in_lattice ! is surface crossing in lattice?
- character(MAX_LINE_LEN) :: msg ! output/error message
type(Universe), pointer :: univ
if (p % cell == 0) then
@@ -45,9 +44,9 @@ contains
! if particle couldn't be located, print error
if (.not. found_cell) then
- write(msg, '(A,3ES11.3)') &
+ write(message, '(A,3ES11.3)') &
"Could not locate cell for particle at: ", p % xyz
- call fatal_error(msg)
+ call fatal_error()
end if
! set birth cell attribute
@@ -55,13 +54,13 @@ contains
end if
if (verbosity >= 9) then
- msg = "Simulating Particle " // trim(int_to_str(p % uid))
- call message(msg)
+ message = "Simulating Particle " // trim(int_to_str(p % id))
+ call write_message()
end if
if (verbosity >= 10) then
- msg = " Born in cell " // trim(int_to_str(cells(p%cell)%uid))
- call message(msg)
+ message = " Born in cell " // trim(int_to_str(cells(p%cell)%id))
+ call write_message()
end if
! find energy index, interpolation factor
@@ -79,9 +78,6 @@ contains
! Select smaller of the two distances
distance = min(d_to_boundary, d_to_collision)
- ! Save original coordinates of particle
- p % last_xyz = p % xyz
-
! Advance particle
p % xyz = p % xyz + distance * p % uvw
p % xyz_local = p % xyz_local + distance * p % uvw
@@ -100,6 +96,9 @@ contains
! collision
p % surface = 0
call collision(p)
+
+ ! Save coordinates at collision for tallying purposes
+ p % last_xyz = p % xyz
end if
end do
@@ -393,8 +392,8 @@ contains
! check for very low energy
if (p % E < 1.0e-100_8) then
p % alive = .false.
- ! msg = "Killing neutron with extremely low energy"
- ! call warning(msg)
+ message = "Killing neutron with extremely low energy"
+ call warning()
end if
! Score collision estimator tallies for any macro tallies -- this is done
@@ -413,6 +412,7 @@ contains
if (tallies_on) then
call score_tally(p, scattered, fissioned)
+ call score_surface_current(p)
end if
! Reset number of particles banked during collision
@@ -444,7 +444,6 @@ contains
real(8) :: prob ! cumulative probability
real(8) :: cutoff ! random number
real(8) :: atom_density ! atom density of nuclide in atom/b-cm
- character(MAX_LINE_LEN) :: msg ! output/error message
type(Material), pointer :: mat => null()
type(Nuclide), pointer :: nuc => null()
type(Reaction), pointer :: rxn => null()
@@ -464,8 +463,8 @@ contains
! Check to make sure that a nuclide was sampled
if (i > mat % n_nuclides) then
- msg = "Did not sample any nuclide during collision."
- call fatal_error(msg)
+ message = "Did not sample any nuclide during collision."
+ call fatal_error()
end if
! Find atom density and microscopic total cross section
@@ -528,9 +527,9 @@ contains
! Check to make sure partial fission reaction sampled
if (i > nuc % n_fission) then
- msg = "Did not sample any partial fission reaction for " // &
- "survival biasing in " // trim(nuc % name)
- call fatal_error(msg)
+ message = "Did not sample any partial fission reaction " // &
+ "for survival biasing in " // trim(nuc % name)
+ call fatal_error()
end if
rxn => nuc % reactions(nuc % index_fission(i))
@@ -636,10 +635,10 @@ contains
! Check to make sure inelastic scattering reaction sampled
if (i > nuc % n_reaction) then
- msg = "Did not sample any reaction for nuclide " // &
+ message = "Did not sample any reaction for nuclide " // &
trim(nuc % name) // " on material " // &
- trim(int_to_str(mat % uid))
- call fatal_error(msg)
+ trim(int_to_str(mat % id))
+ call fatal_error()
end if
rxn => nuc % reactions(i)
@@ -767,7 +766,6 @@ contains
real(8) :: mu_i1jk ! outgoing cosine k for E_in(i+1) and E_out(j)
real(8) :: prob ! probability for sampling Bragg edge
real(8) :: u, v, w ! directional cosines
- character(MAX_LINE_LEN) :: msg
type(SAB_Table), pointer :: sab => null()
! Get pointer to S(a,b) table
@@ -865,7 +863,7 @@ contains
j = 1
end if
else
- msg = "Invalid secondary energy mode on S(a,b) table " // &
+ message = "Invalid secondary energy mode on S(a,b) table " // &
trim(sab % name)
end if
@@ -1027,7 +1025,6 @@ contains
real(8) :: yield ! delayed neutron precursor yield
real(8) :: prob ! cumulative probability
logical :: actual_event ! did fission actually occur? (no survival biasing)
- character(MAX_LINE_LEN) :: msg ! error message
type(Nuclide), pointer :: nuc
! Get pointer to nuclide
@@ -1079,7 +1076,7 @@ contains
if (nu == 0 .or. n_bank == 3*n_particles) return
do i = n_bank + 1, min(n_bank + nu, 3*n_particles)
! Bank source neutrons by copying particle data
- fission_bank(i) % uid = p % uid
+ fission_bank(i) % id = p % id
fission_bank(i) % xyz = p % xyz
! sample cosine of angle
@@ -1099,9 +1096,9 @@ contains
NR = nuc % nu_d_precursor_data(loc + 1)
NE = nuc % nu_d_precursor_data(loc + 2 + 2*NR)
if (NR > 0) then
- msg = "Multiple interpolation regions not supported while &
+ message = "Multiple interpolation regions not supported while &
&sampling delayed neutron precursor yield."
- call fatal_error(msg)
+ call fatal_error()
end if
! interpolate on energy grid
@@ -1254,27 +1251,6 @@ contains
end subroutine inelastic_scatter
-!===============================================================================
-! N_ABSORPTION handles all absorbing reactions, i.e. (n,gamma), (n,p), (n,a),
-! etc.
-!===============================================================================
-
- subroutine n_absorption(p)
-
- type(Particle), pointer :: p
-
- integer :: cell_num ! user-specified cell number
- character(MAX_LINE_LEN) :: msg ! output/error message
-
- p % alive = .false.
- if (verbosity >= 10) then
- cell_num = cells(p % cell)%uid
- msg = " Absorbed in cell " // trim(int_to_str(cell_num))
- call message(msg)
- end if
-
- end subroutine n_absorption
-
!===============================================================================
! SAMPLE_ANGLE samples the cosine of the angle between incident and exiting
! particle directions either from 32 equiprobable bins or from a tabular
@@ -1302,7 +1278,6 @@ contains
real(8) :: c_k ! cumulative frequency at k
real(8) :: c_k1 ! cumulative frequency at k+1
real(8) :: p0,p1 ! probability distribution
- character(MAX_LINE_LEN) :: msg ! error message
! check if reaction has angular distribution -- if not, sample outgoing
! angle isotropically
@@ -1379,20 +1354,20 @@ contains
mu = mu0 + (sqrt(p0*p0 + 2*frac*(xi - c_k))-p0)/frac
end if
else
- msg = "Unknown interpolation type: " // trim(int_to_str(interp))
- call fatal_error(msg)
+ message = "Unknown interpolation type: " // trim(int_to_str(interp))
+ call fatal_error()
end if
if (abs(mu) > ONE) then
- msg = "Sampled cosine of angle outside [-1, 1)."
- call warning(msg)
+ message = "Sampled cosine of angle outside [-1, 1)."
+ call warning()
mu = sign(ONE,mu)
end if
else
- msg = "Unknown angular distribution type: " // trim(int_to_str(type))
- call fatal_error(msg)
+ message = "Unknown angular distribution type: " // trim(int_to_str(type))
+ call fatal_error()
end if
end function sample_angle
@@ -1501,15 +1476,14 @@ contains
real(8) :: E_max ! parameter for n-body dist
real(8) :: x, y, v ! intermediate variables for n-body dist
real(8) :: r1, r2, r3, r4, r5, r6
- character(MAX_LINE_LEN) :: msg ! error message
! TODO: If there are multiple scattering laws, sample scattering law
! Check for multiple interpolation regions
if (edist % n_interp > 0) then
- msg = "Multiple interpolation regions not supported while &
+ message = "Multiple interpolation regions not supported while &
&attempting to sample secondary energy distribution."
- call fatal_error(msg)
+ call fatal_error()
end if
! Determine which secondary energy distribution law to use
@@ -1524,9 +1498,9 @@ contains
NE = edist % data(2 + 2*NR)
NET = edist % data(3 + 2*NR + NE)
if (NR > 0) then
- msg = "Multiple interpolation regions not supported while &
+ message = "Multiple interpolation regions not supported while &
&attempting to sample equiprobable energy bins."
- call fatal_error(msg)
+ call fatal_error()
end if
! determine index on incoming energy grid and interpolation factor
@@ -1591,9 +1565,9 @@ contains
NR = edist % data(1)
NE = edist % data(2 + 2*NR)
if (NR > 0) then
- msg = "Multiple interpolation regions not supported while &
+ message = "Multiple interpolation regions not supported while &
&attempting to sample continuous tabular distribution."
- call fatal_error(msg)
+ call fatal_error()
end if
! find energy bin and calculate interpolation factor -- if the energy is
@@ -1650,9 +1624,9 @@ contains
if (ND > 0) then
! discrete lines present
- msg = "Discrete lines in continuous tabular distributed not &
+ message = "Discrete lines in continuous tabular distributed not &
&yet supported"
- call fatal_error(msg)
+ call fatal_error()
end if
! determine outgoing energy bin
@@ -1685,8 +1659,8 @@ contains
& p_l_k)/frac
end if
else
- msg = "Unknown interpolation type: " // trim(int_to_str(INTT))
- call fatal_error(msg)
+ message = "Unknown interpolation type: " // trim(int_to_str(INTT))
+ call fatal_error()
end if
! Now interpolate between incident energy bins i and i + 1
@@ -1762,17 +1736,17 @@ contains
! KALBACH-MANN CORRELATED SCATTERING
if (.not. present(mu_out)) then
- msg = "Law 44 called without giving mu_out as argument."
- call fatal_error(msg)
+ message = "Law 44 called without giving mu_out as argument."
+ call fatal_error()
end if
! read number of interpolation regions and incoming energies
NR = edist % data(1)
NE = edist % data(2 + 2*NR)
if (NR > 0) then
- msg = "Multiple interpolation regions not supported while &
+ message = "Multiple interpolation regions not supported while &
&attempting to sample Kalbach-Mann distribution."
- call fatal_error(msg)
+ call fatal_error()
end if
! find energy bin and calculate interpolation factor -- if the energy is
@@ -1830,9 +1804,9 @@ contains
if (ND > 0) then
! discrete lines present
- msg = "Discrete lines in continuous tabular distributed not &
+ message = "Discrete lines in continuous tabular distributed not &
&yet supported"
- call fatal_error(msg)
+ call fatal_error()
end if
! determine outgoing energy bin
@@ -1879,8 +1853,8 @@ contains
KM_R = R_k + (R_k1 - R_k)*(E_out - E_l_k)/(E_l_k1 - E_l_k)
KM_A = A_k + (A_k1 - A_k)*(E_out - E_l_k)/(E_l_k1 - E_l_k)
else
- msg = "Unknown interpolation type: " // trim(int_to_str(INTT))
- call fatal_error(msg)
+ message = "Unknown interpolation type: " // trim(int_to_str(INTT))
+ call fatal_error()
end if
! Now interpolate between incident energy bins i and i + 1
@@ -1905,17 +1879,17 @@ contains
! CORRELATED ENERGY AND ANGLE DISTRIBUTION
if (.not. present(mu_out)) then
- msg = "Law 44 called without giving mu_out as argument."
- call fatal_error(msg)
+ message = "Law 44 called without giving mu_out as argument."
+ call fatal_error()
end if
! read number of interpolation regions and incoming energies
NR = edist % data(1)
NE = edist % data(2 + 2*NR)
if (NR > 0) then
- msg = "Multiple interpolation regions not supported while &
+ message = "Multiple interpolation regions not supported while &
&attempting to sample correlated energy-angle distribution."
- call fatal_error(msg)
+ call fatal_error()
end if
! find energy bin and calculate interpolation factor -- if the energy is
@@ -1973,9 +1947,9 @@ contains
if (ND > 0) then
! discrete lines present
- msg = "Discrete lines in continuous tabular distributed not &
+ message = "Discrete lines in continuous tabular distributed not &
&yet supported"
- call fatal_error(msg)
+ call fatal_error()
end if
! determine outgoing energy bin
@@ -2009,8 +1983,8 @@ contains
& p_l_k)/frac
end if
else
- msg = "Unknown interpolation type: " // trim(int_to_str(INTT))
- call fatal_error(msg)
+ message = "Unknown interpolation type: " // trim(int_to_str(INTT))
+ call fatal_error()
end if
! Now interpolate between incident energy bins i and i + 1
@@ -2062,8 +2036,8 @@ contains
mu_out = mu_k + (sqrt(p_k*p_k + 2*frac*(r3 - c_k))-p_k)/frac
end if
else
- msg = "Unknown interpolation type: " // trim(int_to_str(JJ))
- call fatal_error(msg)
+ message = "Unknown interpolation type: " // trim(int_to_str(JJ))
+ call fatal_error()
end if
case (66)
diff --git a/src/search.f90 b/src/search.f90
index 445b72acdf..dd9b02c7ff 100644
--- a/src/search.f90
+++ b/src/search.f90
@@ -2,6 +2,7 @@ module search
use constants, only: ONE, MAX_LINE_LEN
use error, only: fatal_error
+ use global, only: message
contains
@@ -20,14 +21,13 @@ contains
integer :: L
integer :: R
real(8) :: testval
- character(MAX_LINE_LEN) :: msg
L = 1
R = n
if (val < array(L) .or. val > array(R)) then
- msg = "Value outside of array during binary search"
- call fatal_error(msg)
+ message = "Value outside of array during binary search"
+ call fatal_error()
end if
do while (R - L > 1)
diff --git a/src/source.f90 b/src/source.f90
index 8add3cb5ca..c448da0a1b 100644
--- a/src/source.f90
+++ b/src/source.f90
@@ -5,7 +5,7 @@ module source
use cross_section_header, only: Nuclide
use global
use mcnp_random, only: rang, RN_init_particle
- use output, only: message
+ use output, only: write_message
use particle_header, only: Particle, initialize_particle
use physics, only: watt_spectrum
@@ -33,7 +33,7 @@ contains
character(MAX_LINE_LEN) :: msg ! error message
msg = 'Initializing source particles...'
- call message(msg, 6)
+ call write_message(6)
! Determine maximum amount of particles to simulate on each processor
maxwork = ceiling(real(n_particles)/n_procs,8)
@@ -51,7 +51,7 @@ contains
! Initialize first cycle source bank
do i = 0, n_procs - 1
if (rank == i) then
- ! UID's of first and last source particles
+ ! ID's of first and last source particles
bank_first = i*maxwork + 1
bank_last = min((i+1)*maxwork, n_particles)
@@ -66,7 +66,7 @@ contains
! sample position
r = (/ (rang(), k = 1,3) /)
- p % uid = j
+ p % id = j
p % xyz = p_min + r*(p_max - p_min)
p % xyz_local = p % xyz
p % last_xyz = p % xyz
@@ -120,8 +120,8 @@ contains
! point to next source particle
p => source_bank(source_index)
- ! set uid
- p % uid = bank_first + source_index - 1
+ ! set id
+ p % id = bank_first + source_index - 1
end function get_source_particle
diff --git a/src/string.f90 b/src/string.f90
index 54453b2123..e3ed5f11ad 100644
--- a/src/string.f90
+++ b/src/string.f90
@@ -2,6 +2,7 @@ module string
use constants, only: MAX_WORDS, MAX_LINE_LEN, ERROR_INT, ERROR_REAL
use error, only: warning
+ use global, only: message
implicit none
@@ -31,7 +32,6 @@ contains
integer :: i ! current index
integer :: i_start ! starting index of word
integer :: i_end ! ending index of word
- character(MAX_LINE_LEN) :: msg
i_start = 0
i_end = 0
@@ -49,9 +49,9 @@ contains
if (i_end > 0) then
n = n + 1
if (i_end - i_start + 1 > len(words(n))) then
- msg = "The word '" // string(i_start:i_end) // "' is longer than " &
- & // "the space allocated for it."
- call warning(msg)
+ message = "The word '" // string(i_start:i_end) // &
+ "' is longer than the space allocated for it."
+ call warning()
end if
words(n) = string(i_start:i_end)
! reset indices
diff --git a/src/tally.f90 b/src/tally.f90
index fcd49e04f7..84e10e8791 100644
--- a/src/tally.f90
+++ b/src/tally.f90
@@ -1,12 +1,11 @@
module tally
use constants
- use cross_section, only: get_macro_xs
use error, only: fatal_error
use global
use mesh, only: get_mesh_bin, bin_to_mesh_indices, get_mesh_indices
use mesh_header, only: StructuredMesh
- use output, only: message, header
+ use output, only: write_message, header
use search, only: binary_search
use string, only: int_to_str, real_to_str
use tally_header, only: TallyScore, TallyMapItem, TallyMapElement
@@ -36,13 +35,12 @@ contains
real(8), save :: k1 = 0. ! accumulated keff
real(8), save :: k2 = 0. ! accumulated keff**2
real(8) :: std ! stdev of keff over active cycles
- character(MAX_LINE_LEN) :: msg ! output/error message
#ifdef MPI
integer :: ierr
#endif
- msg = "Calculate cycle keff..."
- call message(msg, 8)
+ message = "Calculate cycle keff..."
+ call write_message(8)
! set k1 and k2 at beginning of run
if (i_cycle == 1) then
@@ -103,7 +101,6 @@ contains
integer :: n ! number of bins
integer :: filter_bins ! running total of number of filter bins
integer :: score_bins ! number of scoring bins
- character(MAX_LINE_LEN) :: msg ! output/error message
type(TallyObject), pointer :: t => null()
type(StructuredMesh), pointer :: m => null()
@@ -243,8 +240,8 @@ contains
if (n > 0) then
score_bins = n
else
- msg = "Must have macro tally bins!"
- call fatal_error(msg)
+ message = "Must have macro tally bins!"
+ call fatal_error()
end if
! Allocate scores for tally
@@ -319,7 +316,6 @@ contains
logical :: in_mesh
logical :: has_energyout_bin
logical :: analog
- character(MAX_LINE_LEN) :: msg
type(TallyObject), pointer :: t
type(StructuredMesh), pointer :: m
@@ -334,11 +330,8 @@ contains
do i = 1, n_tallies
t => tallies(i)
- ! Handle surface current tallies separately
- if (t % surface_current) then
- call score_surface_current(p, t)
- cycle
- end if
+ ! Surface current tallies are treated separately
+ if (t % surface_current) cycle
! =======================================================================
! DETERMINE SCORING BIN COMBINATION
@@ -391,12 +384,7 @@ contains
call get_mesh_bin(m, p % xyz, mesh_bin, in_mesh)
if (.not. in_mesh) cycle
- if (t % surface_current) then
- msg = "Surface current mesh tally not yet implemented."
- call fatal_error(msg)
- else
- bins(T_MESH) = mesh_bin
- end if
+ bins(T_MESH) = mesh_bin
else
bins(T_MESH) = 1
end if
@@ -500,8 +488,8 @@ contains
case (MACRO_NU_FISSION)
cycle
case default
- msg = "Invalid macro reaction on analog tally."
- call fatal_error(msg)
+ message = "Invalid macro reaction on analog tally."
+ call fatal_error()
end select
end if
@@ -582,16 +570,17 @@ contains
end subroutine score_tally
!===============================================================================
-! SCORE_SURFACE_CURRENT
+! SCORE_SURFACE_CURRENT tallies surface crossings in a mesh tally by manually
+! determining which mesh surfaces were crossed
!===============================================================================
- subroutine score_surface_current(p, t)
+ subroutine score_surface_current(p)
type(Particle), pointer :: p
- type(TallyObject), pointer :: t
integer :: i ! loop indices
integer :: j ! loop indices
+ integer :: k ! loop indices
integer :: ijk0(3) ! indices of starting coordinates
integer :: ijk1(3) ! indices of ending coordinates
integer :: n_cross ! number of surface crossings
@@ -607,206 +596,211 @@ contains
logical :: x_same ! same starting/ending x index (i)
logical :: y_same ! same starting/ending y index (j)
logical :: z_same ! same starting/ending z index (k)
+ type(TallyObject), pointer :: t => null()
type(StructuredMesh), pointer :: m => null()
- ! Get pointer to mesh
- m => meshes(t % mesh)
+ do i = 1, n_tallies
+ ! Copy starting and ending location of particle
+ xyz0 = p % last_xyz
+ xyz1 = p % xyz
- ! Copy starting and ending location of particle
- xyz0 = p % last_xyz
- xyz1 = p % xyz
+ ! Get pointer to tally and mesh
+ t => tallies(i)
+ m => meshes(t % mesh)
- ! Determine indices for starting and ending location
- call get_mesh_indices(m, xyz0, ijk0, start_in_mesh)
- call get_mesh_indices(m, xyz1, ijk1, end_in_mesh)
+ ! Skip non-surface-current tallies
+ if (.not. t % surface_current) cycle
- ! Check to make sure start or end is in mesh
- if ((.not. start_in_mesh) .and. (.not. end_in_mesh)) return
+ ! Determine indices for starting and ending location
+ call get_mesh_indices(m, xyz0, ijk0, start_in_mesh)
+ call get_mesh_indices(m, xyz1, ijk1, end_in_mesh)
- ! Calculate number of surface crossings
- n_cross = sum(abs(ijk1 - ijk0))
- if (n_cross == 0) return
+ ! Check to make sure start or end is in mesh
+ if ((.not. start_in_mesh) .and. (.not. end_in_mesh)) cycle
- ! Copy particle's direction
- uvw = p % uvw
+ ! Calculate number of surface crossings
+ n_cross = sum(abs(ijk1 - ijk0))
+ if (n_cross == 0) cycle
- ! ==========================================================================
- ! SPECIAL CASES WHERE TWO INDICES ARE THE SAME
+ ! Copy particle's direction
+ uvw = p % uvw
- x_same = (ijk0(1) == ijk1(1))
- y_same = (ijk0(2) == ijk1(2))
- z_same = (ijk0(3) == ijk1(3))
+ ! =======================================================================
+ ! SPECIAL CASES WHERE TWO INDICES ARE THE SAME
- if (x_same .and. y_same) then
- ! Only z crossings
- if (uvw(3) > 0) then
- do i = ijk0(3), ijk1(3) - 1
- ijk0(3) = i
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + OUT_TOP
- call add_to_score(t % scores(score_index, 1), p % last_wgt)
- end if
- end do
- else
- do i = ijk0(3) - 1, ijk1(3), -1
- ijk0(3) = i
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + IN_TOP
- call add_to_score(t % scores(score_index, 1), p % last_wgt)
- end if
- end do
- end if
- return
- elseif (x_same .and. z_same) then
- ! Only y crossings
- if (uvw(2) > 0) then
- do i = ijk0(2), ijk1(2) - 1
- ijk0(2) = i
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + OUT_FRONT
- call add_to_score(t % scores(score_index, 1), p % last_wgt)
- end if
- end do
- else
- do i = ijk0(2) - 1, ijk1(2), -1
- ijk0(2) = i
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + IN_FRONT
- call add_to_score(t % scores(score_index, 1), p % last_wgt)
- end if
- end do
- end if
- return
- elseif (y_same .and. z_same) then
- ! Only x crossings
- if (uvw(1) > 0) then
- do i = ijk0(1), ijk1(1) - 1
- ijk0(1) = i
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + OUT_RIGHT
- call add_to_score(t % scores(score_index, 1), p % last_wgt)
- end if
- end do
- else
- do i = ijk0(1) - 1, ijk1(1), -1
- ijk0(1) = i
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + IN_RIGHT
- call add_to_score(t % scores(score_index, 1), p % last_wgt)
- end if
- end do
- end if
- return
- end if
+ x_same = (ijk0(1) == ijk1(1))
+ y_same = (ijk0(2) == ijk1(2))
+ z_same = (ijk0(3) == ijk1(3))
- ! ==========================================================================
- ! GENERIC CASE
-
- ! Bounding coordinates
- do i = 1, 3
- if (uvw(i) > 0) then
- xyz_cross(i) = m % origin(i) + ijk0(i) * m % width(i)
- else
- xyz_cross(i) = m % origin(i) + (ijk0(i) - 1) * m % width(i)
- end if
- end do
-
- do i = 1, n_cross
- ! Reset scoring bin index
- score_index = 0
-
- ! Calculate distance to each bounding surface. We need to treat special
- ! case where the cosine of the angle is zero since this would result in a
- ! divide-by-zero.
-
- do j = 1, 3
- if (uvw(j) == 0) then
- d(j) = INFINITY
+ if (x_same .and. y_same) then
+ ! Only z crossings
+ if (uvw(3) > 0) then
+ do j = ijk0(3), ijk1(3) - 1
+ ijk0(3) = j
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + OUT_TOP
+ call add_to_score(t % scores(score_index, 1), p % last_wgt)
+ end if
+ end do
else
- d(j) = (xyz_cross(j) - xyz0(j))/uvw(j)
+ do j = ijk0(3) - 1, ijk1(3), -1
+ ijk0(3) = j
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + IN_TOP
+ call add_to_score(t % scores(score_index, 1), p % last_wgt)
+ end if
+ end do
+ end if
+ cycle
+ elseif (x_same .and. z_same) then
+ ! Only y crossings
+ if (uvw(2) > 0) then
+ do j = ijk0(2), ijk1(2) - 1
+ ijk0(2) = j
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + OUT_FRONT
+ call add_to_score(t % scores(score_index, 1), p % last_wgt)
+ end if
+ end do
+ else
+ do j = ijk0(2) - 1, ijk1(2), -1
+ ijk0(2) = j
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + IN_FRONT
+ call add_to_score(t % scores(score_index, 1), p % last_wgt)
+ end if
+ end do
+ end if
+ cycle
+ elseif (y_same .and. z_same) then
+ ! Only x crossings
+ if (uvw(1) > 0) then
+ do j = ijk0(1), ijk1(1) - 1
+ ijk0(1) = j
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + OUT_RIGHT
+ call add_to_score(t % scores(score_index, 1), p % last_wgt)
+ end if
+ end do
+ else
+ do j = ijk0(1) - 1, ijk1(1), -1
+ ijk0(1) = j
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + IN_RIGHT
+ call add_to_score(t % scores(score_index, 1), p % last_wgt)
+ end if
+ end do
+ end if
+ cycle
+ end if
+
+ ! =======================================================================
+ ! GENERIC CASE
+
+ ! Bounding coordinates
+ do j = 1, 3
+ if (uvw(j) > 0) then
+ xyz_cross(j) = m % origin(j) + ijk0(j) * m % width(j)
+ else
+ xyz_cross(j) = m % origin(j) + (ijk0(j) - 1) * m % width(j)
end if
end do
- ! Determine the closest bounding surface of the mesh cell by calculating
- ! the minimum distance
+ do k = 1, n_cross
+ ! Reset scoring bin index
+ score_index = 0
- distance = minval(d)
-
- ! Now use the minimum distance and diretion of the particle to determine
- ! which surface was crossed
-
- if (distance == d(1)) then
- if (uvw(1) > 0) then
- ! Crossing into right mesh cell -- this is treated as outgoing
- ! current from (i,j,k)
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + OUT_RIGHT
+ ! Calculate distance to each bounding surface. We need to treat
+ ! special case where the cosine of the angle is zero since this would
+ ! result in a divide-by-zero.
+
+ do j = 1, 3
+ if (uvw(j) == 0) then
+ d(j) = INFINITY
+ else
+ d(j) = (xyz_cross(j) - xyz0(j))/uvw(j)
end if
- ijk0(1) = ijk0(1) + 1
- xyz_cross(1) = xyz_cross(1) + m % width(1)
- else
- ! Crossing into left mesh cell -- this is treated as incoming
- ! current in (i-1,j,k)
- ijk0(1) = ijk0(1) - 1
- xyz_cross(1) = xyz_cross(1) - m % width(1)
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + IN_RIGHT
+ end do
+
+ ! Determine the closest bounding surface of the mesh cell by
+ ! calculating the minimum distance
+
+ distance = minval(d)
+
+ ! Now use the minimum distance and diretion of the particle to
+ ! determine which surface was crossed
+
+ if (distance == d(1)) then
+ if (uvw(1) > 0) then
+ ! Crossing into right mesh cell -- this is treated as outgoing
+ ! current from (i,j,k)
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + OUT_RIGHT
+ end if
+ ijk0(1) = ijk0(1) + 1
+ xyz_cross(1) = xyz_cross(1) + m % width(1)
+ else
+ ! Crossing into left mesh cell -- this is treated as incoming
+ ! current in (i-1,j,k)
+ ijk0(1) = ijk0(1) - 1
+ xyz_cross(1) = xyz_cross(1) - m % width(1)
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + IN_RIGHT
+ end if
+ end if
+ elseif (distance == d(2)) then
+ if (uvw(2) > 0) then
+ ! Crossing into front mesh cell -- this is treated as outgoing
+ ! current in (i,j,k)
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + OUT_FRONT
+ end if
+ ijk0(2) = ijk0(2) + 1
+ xyz_cross(2) = xyz_cross(2) + m % width(2)
+ else
+ ! Crossing into back mesh cell -- this is treated as incoming
+ ! current in (i,j-1,k)
+ ijk0(2) = ijk0(2) - 1
+ xyz_cross(2) = xyz_cross(2) - m % width(2)
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + IN_FRONT
+ end if
+ end if
+ else if (distance == d(3)) then
+ if (uvw(3) > 0) then
+ ! Crossing into top mesh cell -- this is treated as outgoing
+ ! current in (i,j,k)
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + OUT_TOP
+ end if
+ ijk0(3) = ijk0(3) + 1
+ xyz_cross(3) = xyz_cross(3) + m % width(3)
+ else
+ ! Crossing into bottom mesh cell -- this is treated as incoming
+ ! current in (i,j,k-1)
+ ijk0(3) = ijk0(3) - 1
+ xyz_cross(3) = xyz_cross(3) - m % width(3)
+ if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
+ score_index = sum(t % stride(1:3) * ijk0) + IN_TOP
+ end if
end if
end if
- elseif (distance == d(2)) then
- if (uvw(2) > 0) then
- ! Crossing into front mesh cell -- this is treated as outgoing
- ! current in (i,j,k)
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + OUT_FRONT
- end if
- ijk0(2) = ijk0(2) + 1
- xyz_cross(2) = xyz_cross(2) + m % width(2)
- else
- ! Crossing into back mesh cell -- this is treated as incoming
- ! current in (i,j-1,k)
- ijk0(2) = ijk0(2) - 1
- xyz_cross(2) = xyz_cross(2) - m % width(2)
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + IN_FRONT
- end if
+
+ ! Check for errors
+ if (score_index < 0 .or. score_index > t % n_total_bins) then
+ message = "Score index outside range."
+ call fatal_error()
end if
- else if (distance == d(3)) then
- if (uvw(3) > 0) then
- ! Crossing into top mesh cell -- this is treated as outgoing
- ! current in (i,j,k)
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + OUT_TOP
- end if
- ijk0(3) = ijk0(3) + 1
- xyz_cross(3) = xyz_cross(3) + m % width(3)
- else
- ! Crossing into bottom mesh cell -- this is treated as incoming
- ! current in (i,j,k-1)
- ijk0(3) = ijk0(3) - 1
- xyz_cross(3) = xyz_cross(3) - m % width(3)
- if (all(ijk0 >= 0) .and. all(ijk0 <= m % dimension)) then
- score_index = sum(t % stride(1:3) * ijk0) + IN_TOP
- end if
+
+ ! Add to surface current tally
+ if (score_index > 0) then
+ call add_to_score(t % scores(score_index, 1), p % last_wgt)
end if
- end if
- ! Check for errors
- if (score_index < 0 .or. score_index > t % n_total_bins) then
- print *, ijk0
- print *, ijk1
- print *, score_index
- print *, t % stride(1:3), t % n_total_bins
- call fatal_error("Score_index outside range.")
- end if
+ ! Calculate new coordinates
+ xyz0 = xyz0 + distance * uvw
+ end do
- ! Add to surface current tally
- if (score_index > 0) then
- call add_to_score(t % scores(score_index, 1), p % last_wgt)
- end if
-
- ! Calculate new coordinates
- xyz0 = xyz0 + distance * uvw
end do
end subroutine score_surface_current
@@ -995,7 +989,7 @@ contains
t => tallies(i)
! Write header block
- call header("TALLY " // trim(int_to_str(t % uid)), 3, UNIT_TALLY)
+ call header("TALLY " // trim(int_to_str(t % id)), 3, UNIT_TALLY)
! Handle surface current tallies separately
if (t % surface_current) then
@@ -1227,19 +1221,19 @@ contains
select case(filter_type)
case (T_UNIVERSE)
index = t % universe_bins(bin) % scalar
- label = int_to_str(universes(index) % uid)
+ label = int_to_str(universes(index) % id)
case (T_MATERIAL)
index = t % material_bins(bin) % scalar
- label = int_to_str(materials(index) % uid)
+ label = int_to_str(materials(index) % id)
case (T_CELL)
index = t % cell_bins(bin) % scalar
- label = int_to_str(cells(index) % uid)
+ label = int_to_str(cells(index) % id)
case (T_CELLBORN)
index = t % cellborn_bins(bin) % scalar
- label = int_to_str(cells(index) % uid)
+ label = int_to_str(cells(index) % id)
case (T_SURFACE)
index = t % surface_bins(bin) % scalar
- label = int_to_str(surfaces(index) % uid)
+ label = int_to_str(surfaces(index) % id)
case (T_MESH)
m => meshes(t % mesh)
allocate(ijk(m % n_dimension))
diff --git a/src/tally_header.f90 b/src/tally_header.f90
index a525461795..832b12b70f 100644
--- a/src/tally_header.f90
+++ b/src/tally_header.f90
@@ -58,7 +58,7 @@ module tally_header
type TallyObject
! Basic data
- integer :: uid
+ integer :: id
integer :: type
real(8) :: volume
logical :: surface_current = .false.
diff --git a/src/utils/convert_xsdata.py b/src/utils/convert_xsdata.py
new file mode 100755
index 0000000000..c5f02f22ab
--- /dev/null
+++ b/src/utils/convert_xsdata.py
@@ -0,0 +1,147 @@
+#!/usr/bin/env python
+
+import os
+import sys
+from xml.dom.minidom import getDOMImplementation
+
+types = {1: "neutron", 2: "dosimetry", 3: "thermal"}
+
+class Xsdata(object):
+
+ def __init__(self, filename):
+ self._table_dict = {}
+ self.tables = []
+
+ for line in open(filename, 'r'):
+ words = line.split()
+
+ # If this listing is just an alias listing, only assign the alias
+ # attribute
+ name = words[1]
+ alias = words[0]
+ table = self.find_table(name)
+ if table:
+ if name not in table.alias:
+ table.alias.append(alias)
+ continue
+
+ table = XsdataTable()
+ table.name = name
+ table.type = types[int(words[2])]
+ table.zaid = int(words[3])
+ table.metastable = int(words[4])
+ table.awr = float(words[5])
+ table.temperature = 8.6173423e-11 * float(words[6])
+ table.binary = int(words[7])
+ table.path = words[8]
+
+ self.tables.append(table)
+ self._table_dict[name] = table
+
+ # Check for common directory
+ self.directory = os.path.dirname(self.tables[0].path)
+ for table in self.tables:
+ if not table.path.startswith(self.directory):
+ self.directory = None
+ break
+
+ def to_xml(self):
+ # Create XML document
+ impl = getDOMImplementation()
+ doc = impl.createDocument(None, "cross_sections", None)
+
+ # Get root element
+ root = doc.documentElement
+
+ # Add a directory node
+ if self.directory:
+ directoryNode = doc.createElement("directory")
+ text = doc.createTextNode(self.directory)
+ directoryNode.appendChild(text)
+ root.appendChild(directoryNode)
+
+ for table in self.tables:
+ table.path = os.path.basename(table.path)
+
+ # Add a node for each table
+ for table in self.tables:
+ node = table.to_xml_node(doc)
+ root.appendChild(node)
+
+ return doc
+
+ def find_table(self, name):
+ if name in self._table_dict:
+ return self._table_dict[name]
+ else:
+ return None
+
+
+class XsdataTable(object):
+
+ def __init__(self):
+ self.alias = []
+
+ def to_xml_node(self, doc):
+ node = doc.createElement("ace_table")
+ node.setAttribute("name", self.name)
+ for attribute in ["alias", "zaid", "type", "metastable",
+ "awr", "temperature", "binary", "path"]:
+ if hasattr(self, attribute):
+ # Join string for alias attribute
+ if attribute == "alias":
+ if not self.alias:
+ continue
+ string = " ".join(self.alias)
+ else:
+ string = "{0}".format(getattr(self,attribute))
+
+ # Skip metastable and binary if 0
+ if attribute == "metastable" and self.metastable == 0:
+ continue
+ if attribute == "binary" and self.binary == 0:
+ continue
+
+ # Create attribute node
+ # nodeAttr = doc.createElement(attribute)
+ # text = doc.createTextNode(string)
+ # nodeAttr.appendChild(text)
+ # node.appendChild(nodeAttr)
+ node.setAttribute(attribute, string)
+ return node
+
+
+if __name__ == '__main__':
+ # Read command line arguments
+ if len(sys.argv) < 3:
+ sys.exit("Usage: convert_xsdata.py xsdataFile xmlFile")
+ xsdataFile = sys.argv[1]
+ xmlFile = sys.argv[2]
+
+ # Read xsdata and create XML document object
+ xsdataObject = Xsdata(xsdataFile)
+ doc = xsdataObject.to_xml()
+
+ # Reduce number of lines
+ lines = doc.toprettyxml(indent=' ')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+ lines = lines.replace('\n ','')
+
+ # Write document in pretty XML to specified file
+ f = open(xmlFile, 'w')
+ f.write(lines)
+ f.close()
diff --git a/src/utils/convert_xsdir.py b/src/utils/convert_xsdir.py
new file mode 100755
index 0000000000..7ff1ea892a
--- /dev/null
+++ b/src/utils/convert_xsdir.py
@@ -0,0 +1,212 @@
+#!/usr/bin/env python
+
+import os
+import sys
+from xml.dom.minidom import getDOMImplementation
+
+class Xsdir(object):
+
+ def __init__(self, filename):
+ self.f = open(filename, 'r')
+ self.filename = os.path.abspath(filename)
+ self.directory = os.path.dirname(filename)
+ self.awr = {}
+ self.tables = []
+
+ # Read first section (DATAPATH)
+ line = self.f.readline()
+ words = line.split()
+ if words:
+ if words[0].lower().startswith('datapath'):
+ index = line.index('=')
+ self.datapath = line[index+1:].strip()
+
+ # Read second section
+ line = self.f.readline()
+ words = line.split()
+ assert len(words) == 3
+ assert words[0].lower() == 'atomic'
+ assert words[1].lower() == 'weight'
+ assert words[2].lower() == 'ratios'
+
+ while True:
+ line = self.f.readline()
+ words = line.split()
+
+ # Check for end of second section
+ if len(words) % 2 != 0 or words[0] == 'directory':
+ break
+
+ for zaid, awr in zip(words[::2], words[1::2]):
+ self.awr[zaid] = awr
+
+ # Read third section
+ while words[0] != 'directory':
+ words = self.f.readline().split()
+
+ while True:
+ words = self.f.readline().split()
+ if not words:
+ break
+
+ # Handle continuation lines
+ while words[-1] == '+':
+ extraWords = self.f.readline().split()
+ words = words + extraWords
+ assert len(words) >= 7
+
+ # Create XsdirTable object and add to line
+ table = XsdirTable(self.directory)
+ self.tables.append(table)
+
+ # All tables have at least 7 attributes
+ table.name = words[0]
+ table.awr = float(words[1])
+ table.filename = words[2]
+ table.access = words[3]
+ table.filetype = int(words[4])
+ table.address = int(words[5])
+ table.tablelength = int(words[6])
+
+ if len(words) > 7:
+ table.recordlength = int(words[7])
+ if len(words) > 8:
+ table.entries = int(words[8])
+ if len(words) > 9:
+ table.temperature = float(words[9])
+ if len(words) > 10:
+ table.ptable = (words[10] == 'ptable')
+
+ def to_xml(self):
+ # Create XML document
+ impl = getDOMImplementation()
+ doc = impl.createDocument(None, "cross_sections", None)
+
+ # Get root element
+ root = doc.documentElement
+
+ # Add a directory node
+ if self.directory:
+ directoryNode = doc.createElement("directory")
+ text = doc.createTextNode(self.directory)
+ directoryNode.appendChild(text)
+ root.appendChild(directoryNode)
+
+ for table in self.tables:
+ table.path = os.path.basename(table.path)
+
+ # Add a node for each table
+ for table in self.tables:
+ if table.name[-1] in ['e', 'p', 'u', 'h', 'g' ,'m', 'd']:
+ continue
+ node = table.to_xml_node(doc)
+ root.appendChild(node)
+
+ return doc
+
+
+class XsdirTable(object):
+
+ def __init__(self, directory=None):
+ self.directory = None
+ self.name = None
+ self.awr = None
+ self.filename = None
+ self.access = None
+ self.filetype = None
+ self.address = None
+ self.tablelength = None
+ self.recordlength = None
+ self.entries = None
+ self.temperature = None
+ self.ptable = False
+
+ @property
+ def path(self):
+ if self.directory:
+ return os.path.join(self.directory, self.filename)
+ else:
+ return self.filename
+
+ @path.setter
+ def path(self, value):
+ self.diretory = ''
+ self.filename = value
+
+ @property
+ def metastable(self):
+ # Only valid for neutron cross-sections
+ if not self.name.endswith('c'):
+ return
+
+ # Handle special case of Am-242 and Am-242m
+ if self.zaid == '95242':
+ return 1
+ elif self.zaid == '95642':
+ return 0
+
+ # All other cases
+ A = int(self.zaid) % 1000
+ if A > 600:
+ return 1
+ else:
+ return 0
+
+ @property
+ def zaid(self):
+ if self.name.endswith('c'):
+ return self.name[:self.name.find('.')]
+ else:
+ return 0
+
+ def to_xml_node(self, doc):
+ node = doc.createElement("ace_table")
+ node.setAttribute("name", self.name)
+ for attribute in ["alias", "zaid", "type", "metastable",
+ "awr", "temperature", "binary", "path"]:
+ if hasattr(self, attribute):
+ # Join string for alias attribute
+ if attribute == "alias":
+ if not self.alias:
+ continue
+ string = " ".join(self.alias)
+ else:
+ string = "{0}".format(getattr(self,attribute))
+
+ # Skip metastable and binary if 0
+ if attribute == "metastable" and self.metastable == 0:
+ continue
+ if attribute == "binary" and self.binary == 0:
+ continue
+
+ # Skip any attribute that is none
+ if getattr(self, attribute) is None:
+ continue
+
+ # Create attribute node
+ # nodeAttr = doc.createElement(attribute)
+ # text = doc.createTextNode(string)
+ # nodeAttr.appendChild(text)
+ # node.appendChild(nodeAttr)
+ node.setAttribute(attribute, string)
+ return node
+
+
+if __name__ == '__main__':
+ # Read command line arguments
+ if len(sys.argv) < 3:
+ sys.exit("Usage: convert_xsdir.py xsdirFile xmlFile")
+ xsdirFile = sys.argv[1]
+ xmlFile = sys.argv[2]
+
+ # Read xsdata and create XML document object
+ xsdirObject = Xsdir(xsdirFile)
+ doc = xsdirObject.to_xml()
+
+ # Reduce number of lines
+ lines = doc.toprettyxml(indent=' ')
+
+ # Write document in pretty XML to specified file
+ f = open(xmlFile, 'w')
+ f.write(lines)
+ f.close()
diff --git a/src/xml-fortran/templates/cross_sections_t.xml b/src/xml-fortran/templates/cross_sections_t.xml
new file mode 100644
index 0000000000..d49d798620
--- /dev/null
+++ b/src/xml-fortran/templates/cross_sections_t.xml
@@ -0,0 +1,23 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/src/xml-fortran/templates/geometry_t.xml b/src/xml-fortran/templates/geometry_t.xml
index 2752cb7cee..7aeec2a315 100644
--- a/src/xml-fortran/templates/geometry_t.xml
+++ b/src/xml-fortran/templates/geometry_t.xml
@@ -6,7 +6,7 @@
-
+
@@ -14,14 +14,14 @@
-
+
-
+
diff --git a/src/xml-fortran/templates/materials_t.xml b/src/xml-fortran/templates/materials_t.xml
index 41ca975ee2..dbb939ac12 100644
--- a/src/xml-fortran/templates/materials_t.xml
+++ b/src/xml-fortran/templates/materials_t.xml
@@ -31,7 +31,7 @@
-
+
diff --git a/src/xml-fortran/templates/settings_t.xml b/src/xml-fortran/templates/settings_t.xml
index 7309eb039c..167ff0d63d 100644
--- a/src/xml-fortran/templates/settings_t.xml
+++ b/src/xml-fortran/templates/settings_t.xml
@@ -3,10 +3,6 @@
-
-
-
-
@@ -23,11 +19,11 @@
-
+
| | |