diff --git a/src/meshIO/private/meshIO.cpp b/src/meshIO/private/meshIO.cpp index ee8a6b3..3b340c3 100644 --- a/src/meshIO/private/meshIO.cpp +++ b/src/meshIO/private/meshIO.cpp @@ -2,11 +2,12 @@ #include #include #include +#include #include "meshIO.h" -MeshIO::MeshIO(const std::string &mesh_file) +MeshIO::MeshIO(const std::string &mesh_file, double scale_factor) { mesh_file_ = mesh_file; std::ifstream mesh_stream(mesh_file); @@ -19,6 +20,7 @@ MeshIO::MeshIO(const std::string &mesh_file) { mesh_ = mfem::Mesh(mesh_stream, 1, 10); loaded_ = true; + if (scale_factor != 1.0) { LinearRescale(scale_factor); } } } @@ -26,6 +28,32 @@ MeshIO::~MeshIO() { } +void MeshIO::LinearRescale(double scale_factor) { + 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 + mesh_.EnsureNodes(); + const mfem::GridFunction *nodes = mesh_.GetNodes(); + double *node_data = nodes->GetData(); // THIS IS KEY + + int data_size = nodes->Size(); + + for (int i = 0; i < data_size; ++i) { + node_data[i] *= scale_factor; + } + + // nodes->Update(); /updates the fespace + mesh_.NodesUpdated(); + + +} + bool MeshIO::IsLoaded() const { return loaded_; diff --git a/src/meshIO/public/meshIO.h b/src/meshIO/public/meshIO.h index 5975801..d7772b8 100644 --- a/src/meshIO/public/meshIO.h +++ b/src/meshIO/public/meshIO.h @@ -19,13 +19,19 @@ public: * @brief Constructor that initializes the MeshIO object with a mesh file. * @param mesh_file The name of the mesh file. */ - MeshIO(const std::string &mesh_file); + MeshIO(const std::string &mesh_file, double scale_factor = 1.0); /** * @brief Destructor for the MeshIO class. */ ~MeshIO(); + /** + * @brief Rescale the mesh by a linear factor. + * @param scale_factor The factor by which to scale the mesh. + */ + void LinearRescale(double scale_factor); + /** * @brief Get the mesh object. * @return Reference to the mesh object. diff --git a/src/poly/solver/private/polySolver.cpp b/src/poly/solver/private/polySolver.cpp index 897bbda..8d1be8d 100644 --- a/src/poly/solver/private/polySolver.cpp +++ b/src/poly/solver/private/polySolver.cpp @@ -18,7 +18,7 @@ const std::string SPHERICAL_MESH = std::string(getenv("MESON_SOURCE_ROOT")) + "/ PolySolver::PolySolver(double n, double order) : n(n), order(order), - meshIO(SPHERICAL_MESH), + meshIO(SPHERICAL_MESH, 5), mesh(meshIO.GetMesh()), feCollection(std::make_unique(order, mesh.SpaceDimension())), feSpace(std::make_unique(&mesh, feCollection.get())), @@ -33,6 +33,7 @@ PolySolver::PolySolver(double n, double order) }())), nonLinearSourceCoeff(std::make_unique(1.0)), gaussianCoeff(std::make_unique(0.1)) { + assembleNonlinearForm(); assembleConstraintForm(); }