<commit_before>/*
   Copyright 2016 <NAME>

   Licensed under the Apache License, Version 2.0 (the "License");
   you may not use this file except in compliance with the License.
   You may obtain a copy of the License at

       http://www.apache.org/licenses/LICENSE-2.0

   Unless required by applicable law or agreed to in writing, software
   distributed under the License is distributed on an "AS IS" BASIS,
   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
   See the License for the specific language governing permissions and
   limitations under the License.
*/

#include "xs_mesh_homogenized.hpp"

#include <iomanip>
#include <iostream>
#include <numeric>
#include <sstream>
#include "pugixml.hpp"
#include "util/files.hpp"
#include "util/h5file.hpp"
#include "util/string_utils.hpp"

namespace mocc {
XSMeshHomogenized::XSMeshHomogenized(const CoreMesh &mesh)
    : XSMesh(mesh, MeshTreatment::PIN), mesh_(mesh), flux_(nullptr)
{
    // Set up the non-xs part of the xs mesh
    eubounds_ = mesh_.mat_lib().g_bounds();
    ng_       = eubounds_.size();

    int n_xsreg = std::accumulate(
        mesh_.macroplanes().begin(), mesh_.macroplanes().end(), 0,
        [](int n, const MacroPlane p) { return n + p.plane->n_pin(); });

    // Allocate space to store the cross sections
    this->allocate_xs(n_xsreg, ng_);

    // Set up the regions
    //
    // It is assumed that this type of XS mesh will be used to provide
    // cross sections to a MeshTreatment::PIN type sweeper, possibly using
    // flux-weighted homogenization with flux from a MeshTreatment::PLANE-type
    // sweeper. The mapping from XS mesh region to sweeper region will therefore
    // be from a homogenized pin at the macroplane level to all of the regions
    // in a MeshTreatment::PIN-type mesh within the corresponding macroplane.
    regions_.reserve(n_xsreg);
    n_reg_expanded_ = 0;
    int ixsreg      = 0;
    for (const auto &mplane : mesh_.macroplanes()) {
        VecI ireg;
        ireg.reserve(mplane.iz_max - mplane.iz_min + 1);
        for (unsigned ipin = 0; ipin < mplane.size(); ++ipin) {
            Position pos = mesh_.pin_position(ipin);
            for (int iz = mplane.iz_min; iz <= mplane.iz_max; ++iz) {
                pos.z = iz;
                ireg.push_back(mesh_.coarse_cell(pos));
            }
            regions_.emplace_back(ireg, &xstr_(ixsreg, 0), &xsnf_(ixsreg, 0),
                                  &xsch_(ixsreg, 0), &xsf_(ixsreg, 0),
                                  &xsrm_(ixsreg, 0), ScatteringMatrix());

            n_reg_expanded_ += ireg.size();
            ireg.clear();
            ixsreg++;
        }
    }

    // Homogenize initial cross sections
    ixsreg = 0;
    for (const auto &mplane : mesh_.macroplanes()) {
        for (const auto &mpin : mplane) {
            this->homogenize_region(ixsreg, *mpin, regions_[ixsreg]);
            ixsreg++;
        }
    }
}

/**
 * This creates a "homogenized" cross-section mesh using data prescribed by
 * one or more HDF5 files. The allowance of multiple files is so that an
 * \ref XSMeshHomogenized object can be composed from various sources. For
 * example, imagine a 3-D mesh, where the cross sections are obtained from
 * several different 2-D solutions, and the cross sections for each plane
 * may be drawn from any of the 2-D results.
 *
 * Since the primary use-case for this functionality is for using planar,
 * fine-mesh MoC solutions to provide cross sections to a coarse-mesh Sn
 * solver, the data is specified along with the upper bound macroplane index to
 * which the cross sections are to be applied. The data itself is assumed to
 * be coming from a call to \ref XSMesh::output().
 */
XSMeshHomogenized::XSMeshHomogenized(const CoreMesh &mesh,
                                     const pugi::xml_node &input)
    : XSMeshHomogenized(mesh)
{
    if (input.child("data").empty()) {
        throw EXCEPT("No data found in input tag.");
    }

    // First, validate the data tags. Make sure that they are the right size and
    // cover valid planes in the mesh. There are two ways to do this, either
    // with 2-D data that we apply to planes with ranges, or one 3-D file, which
    // should be 100% compatible, size-wise

    int n_data = 0;
    for (const auto &d : input.children("data")) {
        LogFile << "Found <data /> specification for "
                << d.attribute("file").value() << std::endl;
        n_data++;
    }

    bool single_compat = false;
    if (n_data == 1) {
        // see if the specified file seems compatible
        const auto &data = input.child("data");
        H5Node h5d(data.attribute("file").value(), H5Access::READ);
        auto dims = h5d.dimensions("/xsmesh/xstr/0");
        if ((dims[2] == mesh_.nx()) && (dims[1] == mesh_.ny()) &&
            (dims[0] == mesh_.macroplanes().size())) {
            single_compat = true;
        }
    }

    if (single_compat) {
        LogScreen << "Reading single data file for cross sections" << std::endl;
        this->read_data_single(input.child("data"));
    } else {
        this->read_data_multi(input);
    }

    return;
} // HDF5 constructor

XSMeshHomogenized::XSMeshHomogenized(XSMeshHomogenized &other,
                                     const std::vector<int> &regions):
    XSMeshHomogenized(other)
{
    assert(regions.size() == 0);
    int ixsreg = 0;
    for(auto &xsr: regions_) {
        xsr.reg_.resize(1);
        xsr.reg_[0] = ixsreg;
        ++ixsreg;
    }

    n_reg_expanded_ = regions_.size();

    return;
}

/**
 * Update the XS mesh, incorporating a new estimate of the scalar flux.<commit_msg>Fix bug in the aliased XSMeshHomogenized ctor<commit_after>/*
   Copyright 2016 <NAME>

   Licensed under the Apache License, Version 2.0 (the "License");
   you may not use this file except in compliance with the License.
   You may obtain a copy of the License at

       http://www.apache.org/licenses/LICENSE-2.0

   Unless required by applicable law or agreed to in writing, software
   distributed under the License is distributed on an "AS IS" BASIS,
   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
   See the License for the specific language governing permissions and
   limitations under the License.
*/

#include "xs_mesh_homogenized.hpp"

#include <iomanip>
#include <iostream>
#include <numeric>
#include <sstream>
#include "pugixml.hpp"
#include "util/files.hpp"
#include "util/h5file.hpp"
#include "util/string_utils.hpp"

namespace mocc {
XSMeshHomogenized::XSMeshHomogenized(const CoreMesh &mesh)
    : XSMesh(mesh, MeshTreatment::PIN), mesh_(mesh), flux_(nullptr)
{
    // Set up the non-xs part of the xs mesh
    eubounds_ = mesh_.mat_lib().g_bounds();
    ng_       = eubounds_.size();

    int n_xsreg = std::accumulate(
        mesh_.macroplanes().begin(), mesh_.macroplanes().end(), 0,
        [](int n, const MacroPlane p) { return n + p.plane->n_pin(); });

    // Allocate space to store the cross sections
    this->allocate_xs(n_xsreg, ng_);

    // Set up the regions
    //
    // It is assumed that this type of XS mesh will be used to provide
    // cross sections to a MeshTreatment::PIN type sweeper, possibly using
    // flux-weighted homogenization with flux from a MeshTreatment::PLANE-type
    // sweeper. The mapping from XS mesh region to sweeper region will therefore
    // be from a homogenized pin at the macroplane level to all of the regions
    // in a MeshTreatment::PIN-type mesh within the corresponding macroplane.
    regions_.reserve(n_xsreg);
    n_reg_expanded_ = 0;
    int ixsreg      = 0;
    for (const auto &mplane : mesh_.macroplanes()) {
        VecI ireg;
        ireg.reserve(mplane.iz_max - mplane.iz_min + 1);
        for (unsigned ipin = 0; ipin < mplane.size(); ++ipin) {
            Position pos = mesh_.pin_position(ipin);
            for (int iz = mplane.iz_min; iz <= mplane.iz_max; ++iz) {
                pos.z = iz;
                ireg.push_back(mesh_.coarse_cell(pos));
            }
            regions_.emplace_back(ireg, &xstr_(ixsreg, 0), &xsnf_(ixsreg, 0),
                                  &xsch_(ixsreg, 0), &xsf_(ixsreg, 0),
                                  &xsrm_(ixsreg, 0), ScatteringMatrix());

            n_reg_expanded_ += ireg.size();
            ireg.clear();
            ixsreg++;
        }
    }

    // Homogenize initial cross sections
    ixsreg = 0;
    for (const auto &mplane : mesh_.macroplanes()) {
        for (const auto &mpin : mplane) {
            this->homogenize_region(ixsreg, *mpin, regions_[ixsreg]);
            ixsreg++;
        }
    }
}

/**
 * This creates a "homogenized" cross-section mesh using data prescribed by
 * one or more HDF5 files. The allowance of multiple files is so that an
 * \ref XSMeshHomogenized object can be composed from various sources. For
 * example, imagine a 3-D mesh, where the cross sections are obtained from
 * several different 2-D solutions, and the cross sections for each plane
 * may be drawn from any of the 2-D results.
 *
 * Since the primary use-case for this functionality is for using planar,
 * fine-mesh MoC solutions to provide cross sections to a coarse-mesh Sn
 * solver, the data is specified along with the upper bound macroplane index to
 * which the cross sections are to be applied. The data itself is assumed to
 * be coming from a call to \ref XSMesh::output().
 */
XSMeshHomogenized::XSMeshHomogenized(const CoreMesh &mesh,
                                     const pugi::xml_node &input)
    : XSMeshHomogenized(mesh)
{
    if (input.child("data").empty()) {
        throw EXCEPT("No data found in input tag.");
    }

    // First, validate the data tags. Make sure that they are the right size and
    // cover valid planes in the mesh. There are two ways to do this, either
    // with 2-D data that we apply to planes with ranges, or one 3-D file, which
    // should be 100% compatible, size-wise

    int n_data = 0;
    for (const auto &d : input.children("data")) {
        LogFile << "Found <data /> specification for "
                << d.attribute("file").value() << std::endl;
        n_data++;
    }

    bool single_compat = false;
    if (n_data == 1) {
        // see if the specified file seems compatible
        const auto &data = input.child("data");
        H5Node h5d(data.attribute("file").value(), H5Access::READ);
        auto dims = h5d.dimensions("/xsmesh/xstr/0");
        if ((dims[2] == mesh_.nx()) && (dims[1] == mesh_.ny()) &&
            (dims[0] == mesh_.macroplanes().size())) {
            single_compat = true;
        }
    }

    if (single_compat) {
        LogScreen << "Reading single data file for cross sections" << std::endl;
        this->read_data_single(input.child("data"));
    } else {
        this->read_data_multi(input);
    }

    return;
} // HDF5 constructor

XSMeshHomogenized::XSMeshHomogenized(XSMeshHomogenized &other,
                                     const std::vector<int> &regions)
    : XSMeshHomogenized(other)
{
    assert(regions.size() == 0);

    int ixsreg = 0;
    int iz     = 0;
    for (const auto &mplane : mesh_.macroplanes()) {
        VecI ireg;
        ireg.reserve(mplane.iz_max - mplane.iz_min + 1);
        for (unsigned ipin = 0; ipin < mplane.size(); ++ipin) {
            Position pos = mesh_.pin_position(ipin);
            pos.z        = iz;

            regions_[ixsreg].reg_.resize(1);
            regions_[ixsreg].reg_.front() = mesh_.coarse_cell(pos);
            ixsreg++;
        }
        iz++;
    }

    for (auto &xsr : regions_) {
        int new_reg = xsr.reg_.front();
        xsr.reg_.resize(1);
        xsr.reg_[0] = new_reg;
        ++ixsreg;
    }

    n_reg_expanded_ = regions_.size();

    return;
}

/**
 * Update the XS mesh, incorporating a new estimate of the scalar flux.
{'fraction_non_alphanumeric': 0.11094237745718384, 'fraction_numerical': 0.007761748080654687, 'mean_word_length': 2.486470588235294, 'pattern_counts': {'":': 0, '<': 42, '<?xml version=': 0, '>': 27, 'https://': 0, 'lorem ipsum': 0, 'www.': 2, 'xml': 4}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 8, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '1189256', 'n_tokens_mistral': 3732, 'n_tokens_neox': 3569, 'n_words': 1179}