meshIO.cpp

Go to the documentation of this file.

/* ***********************************************************************
//
//   Copyright (C) 2025 -- The 4D-STAR Collaboration
//   File Author: Emily Boudreaux
//   Last Modified: March 18, 2025
//
//   4DSSE is free software; you can use it and/or modify
//   it under the terms and restrictions the GNU General Library Public
//   License version 3 (GPLv3) as published by the Free Software Foundation.
//
//   4DSSE is distributed in the hope that it will be useful,
//   but WITHOUT ANY WARRANTY; without even the implied warranty of
//   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
//   See the GNU Library General Public License for more details.
//
//   You should have received a copy of the GNU Library General Public License
//   along with this software; if not, write to the Free Software
//   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
//
// *********************************************************************** */
#include "mfem.hpp"
#include <string>
#include <iostream>
#include <fstream>
#include <stdexcept>
 
#include "meshIO.h"
 
namespace serif::mesh {
    MeshIO::MeshIO(const std::string &mesh_file, double scale_factor)
    {
        mesh_file_ = mesh_file;
        std::ifstream mesh_stream(mesh_file);
        if (!mesh_stream)
        {
            throw std::runtime_error("Mesh file not found: " + mesh_file);
        }
        m_mesh = std::make_unique<mfem::Mesh>(mesh_stream, 1, 10);
 
        bool discont = false;
        int space_dim = m_mesh->SpaceDimension();
        int ordering = mfem::Ordering::byVDIM;
        // m_mesh->SetCurvature(3, discont, space_dim, ordering);
        m_mesh->Finalize(false, true);
        loaded_ = true;
        if (scale_factor != 1.0) { LinearRescale(scale_factor); }
    }
 
    MeshIO::~MeshIO() = default;
 
    void MeshIO::LinearRescale(const double scale_factor) const {
        if (!loaded_) {
            throw std::runtime_error("Mesh not loaded before rescaling.");
        }
 
        if (scale_factor <= 0.0) {
            throw std::invalid_argument("Scale factor must be positive.");
        }
 
        // Ensure there are nodes even for linear order meshes
        m_mesh->EnsureNodes();
        const mfem::GridFunction *nodes = m_mesh->GetNodes();
        double *node_data = nodes->GetData();
        const int data_size = nodes->Size();
 
        for (int i = 0; i < data_size; ++i) {
            node_data[i] *= scale_factor;
        }
 
        // nodes->Update(); /updates the fespace
        m_mesh->NodesUpdated();
    }
 
    bool MeshIO::IsLoaded() const
    {
        return loaded_;
    }
 
    mfem::Mesh &MeshIO::GetMesh() const noexcept {
        return *m_mesh;
    }
}