Commit a050260c authored by GILLES Sebastien's avatar GILLES Sebastien

#1444 Ensight output is now rewritten and works properly for all the models inside the library.

parent 72d8d011
......@@ -30,7 +30,6 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
namespace // anonymous
{
using DofWithValuesType = std::pair<Data::DofInformation::const_shared_ptr, double>;
void CreateCaseFile(const std::string& output_directory,
const Data::UnknownInformation::vector_const_shared_ptr& unknown_list,
......@@ -57,6 +56,24 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
};
class IndexMatcher
{
public:
IndexMatcher(std::size_t processor_wise_index,
std::size_t ensight_position);
std::size_t GetProcessorWiseIndex() const noexcept;
std::size_t GetEnsightPosition() const noexcept;
private:
std::size_t processor_wise_index_;
std::size_t ensight_position_;
};
void CreateUnknownFile(const Data::TimeIteration& time_iteration,
const std::string& output_directory,
const std::string& unknown_name,
......@@ -74,8 +91,10 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
{
assert(!(!time_iteration_ptr));
if (time_iteration_ptr->GetNumberingSubsetId() == numbering_subset_id)
{
std::cout << "ADD " << time_iteration_ptr->GetSolutionFilename() << std::endl;
ret.push_back(time_iteration_ptr.get());
}
}
return ret;
......@@ -234,13 +253,16 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
for (auto unknown_index = 0ul; unknown_index < Nunknown; ++unknown_index)
{
std::cout << "UNKNOWN INDEX = " << unknown_index << std::endl;
const auto numbering_subset_id = numbering_subset_id_list[unknown_index];
const auto& unknown_ptr = selected_unknown_list[unknown_index];
assert(!(!unknown_ptr));
const auto& unknown = *unknown_ptr;
// Ensight expects 3D dimension no matter what - even if the mesh is actually 2D the third component
// must be specified with 0..
const auto unknown_dimension =
(unknown.GetNature() == UnknownNS::Nature::vectorial ? mesh.GetDimension() : 1ul);
(unknown.GetNature() == UnknownNS::Nature::vectorial ? 3ul : 1ul);
// Extract the relevant time iteration files
const auto relevant_time_iteration_file_list =
......@@ -250,16 +272,18 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
if (relevant_time_iteration_file_list.empty())
continue;
const auto expected_size = mesh.NprocessorWiseCoord() // all \a Coords really as we run from the full mesh in sequential
* unknown_dimension;
const auto NprogramWiseCoord = mesh.NprocessorWiseCoord() * unknown_dimension; // what we actually want here
// is the number of program-wise coords that are also vertices (#248) but as EnsightOutput is run sequentially
// by reading the full mesh the method called here yields the expected value.
const std::size_t Nprocessor = post_processing.Nprocessor();
std::vector<double> ensight_file_content(expected_size, 0.);
std::vector<double> ensight_file_content(NprogramWiseCoord, 0.);
std::vector<std::vector<double>> ensight_file_content_per_time_step(relevant_time_iteration_file_list.size(),
ensight_file_content);
for (auto rank = 0u; rank < Nprocessor; ++rank)
{
const auto dof_list = ComputeDofListOnVertex(post_processing,
......@@ -271,8 +295,9 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
// Value is the matching position in the Ensight vector; it is the \a Coords position directly
// if the unknown is scalar and a bit more if we consider a vectorial one (as there are two or three dofs
// to write per \a Coord).
std::vector<std::size_t> coords_index_per_dof; // \todo Change name!
coords_index_per_dof.reserve(dof_list.size());
std::vector<IndexMatcher> index_matcher_list; // \todo Change name!
const auto Nprocessor_wise_dof = dof_list.size();
index_matcher_list.reserve(Nprocessor_wise_dof);
// \\ \todo THis doesn't exactly work: consider the case of several unknowns in same NS; I couldn't use thread_local
// vector index in this case!
// On the top of my head, we should here use a vector of a struct which also holds the program-wise index
......@@ -281,8 +306,6 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
// - Program-wise index in MoReFEM
// - Computed Ensight index (e.g. coords_index * unknown_dimension + dof.GetUnknownComponent())
std::cout << "Ndofs = " << dof_list.size() << std::endl;
for (const auto& dof_ptr : dof_list)
{
assert(!(!dof_ptr));
......@@ -297,9 +320,23 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
assert(coords_list.size() == 1ul);
const auto coords_index = coords_list.back();
coords_index_per_dof.push_back(coords_index * unknown_dimension + dof.GetUnknownComponent());
IndexMatcher index_matcher(dof.GetProcessorWiseIndex(),
coords_index * unknown_dimension + dof.GetUnknownComponent());
index_matcher_list.emplace_back(index_matcher);
}
assert(Nprocessor_wise_dof == index_matcher_list.size());
// static bool t = false;
//
// if (!t)
// {
// t = true;
// Utilities::PrintContainer<>::Do(index_matcher_list);
// }
auto iterator_ensight_sol = ensight_file_content_per_time_step.begin();
for (const auto& time_iteration_ptr : time_iteration_list)
......@@ -310,26 +347,30 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
const auto solution = LoadVectorForTimeIteration(time_iteration, rank);
auto& ensight_solution = *iterator_ensight_sol;
for (const auto& item : solution)
{
if (item < 0.)
std::cout << "NEGATIVE VALUE READ: " << item << std::endl;
}
assert(expected_size == ensight_solution.size());
auto& ensight_solution = *iterator_ensight_sol;
const auto Nelt_on_proc = coords_index_per_dof.size();
assert(NprogramWiseCoord == ensight_solution.size());
std::cout << "CHECK -> " << solution.size() << "\t" << Nelt_on_proc << std::endl;
assert(solution.size() >= Nelt_on_proc); // solution might include dofs not considered here
// (not on vertex, related to a different unknown, etc...)
std::cout << "CHECK -> " << solution.size() << "\t" << Nprocessor_wise_dof << std::endl;
assert(solution.size() >= Nprocessor_wise_dof); // solution might include dofs not considered here
// (not on vertex, related to a different unknown, etc...)
for (auto i = 0ul; i < Nelt_on_proc; ++i)
for (auto i = 0ul; i < Nprocessor_wise_dof; ++i)
{
const auto ensight_sol_index = coords_index_per_dof[i];
const auto& index_matcher = index_matcher_list[i];
const auto ensight_sol_index = index_matcher.GetEnsightPosition();
// Solution: dof index should be used!
// Ensight solution: the one computed in coords_index_per_dof.
// Ensight solution: the one computed in index_matcher_list.
assert(ensight_sol_index < ensight_solution.size());
ensight_solution[ensight_sol_index] = solution[i];
ensight_solution[ensight_sol_index] = solution[index_matcher.GetProcessorWiseIndex()];
}
CreateUnknownFile(time_iteration,
......@@ -342,7 +383,6 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
}
}
}
......@@ -480,6 +520,25 @@ namespace MoReFEM::PostProcessingNS::OutputFormat
}
IndexMatcher::IndexMatcher(std::size_t processor_wise_index,
std::size_t ensight_position)
: processor_wise_index_(processor_wise_index),
ensight_position_(ensight_position)
{ }
std::size_t IndexMatcher::GetProcessorWiseIndex() const noexcept
{
return processor_wise_index_;
}
std::size_t IndexMatcher::GetEnsightPosition() const noexcept
{
return ensight_position_;
}
} // namespace anonymous
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment