FEM: Port elmer transient analysis to multiframe results

This commit is contained in:
Stefan Tröger
2025-03-01 12:13:05 +01:00
committed by Benjamin Nauck
parent c815612dc6
commit 22c8d389d4
7 changed files with 168 additions and 67 deletions

View File

@@ -260,41 +260,84 @@ bool FemPostPipeline::canRead(Base::FileInfo File)
return File.hasExtension({"vtk", "vtp", "vts", "vtr", "vti", "vtu", "pvtu", "vtm"});
}
void FemPostPipeline::read(Base::FileInfo File)
vtkSmartPointer<vtkDataObject> FemPostPipeline::dataObjectFromFile(Base::FileInfo File)
{
// checking on the file
if (!File.isReadable()) {
throw Base::FileException("File to load not existing or not readable", File);
}
if (File.hasExtension("vtu")) {
readXMLFile<vtkXMLUnstructuredGridReader>(File.filePath());
return readXMLFile<vtkXMLUnstructuredGridReader>(File.filePath());
}
else if (File.hasExtension("pvtu")) {
readXMLFile<vtkXMLPUnstructuredGridReader>(File.filePath());
return readXMLFile<vtkXMLPUnstructuredGridReader>(File.filePath());
}
else if (File.hasExtension("vtp")) {
readXMLFile<vtkXMLPolyDataReader>(File.filePath());
return readXMLFile<vtkXMLPolyDataReader>(File.filePath());
}
else if (File.hasExtension("vts")) {
readXMLFile<vtkXMLStructuredGridReader>(File.filePath());
return readXMLFile<vtkXMLStructuredGridReader>(File.filePath());
}
else if (File.hasExtension("vtr")) {
readXMLFile<vtkXMLRectilinearGridReader>(File.filePath());
return readXMLFile<vtkXMLRectilinearGridReader>(File.filePath());
}
else if (File.hasExtension("vti")) {
readXMLFile<vtkXMLImageDataReader>(File.filePath());
return readXMLFile<vtkXMLImageDataReader>(File.filePath());
}
else if (File.hasExtension("vtk")) {
readXMLFile<vtkDataSetReader>(File.filePath());
return readXMLFile<vtkDataSetReader>(File.filePath());
}
else if (File.hasExtension("vtm")) {
readXMLFile<vtkXMLMultiBlockDataReader>(File.filePath());
return readXMLFile<vtkXMLMultiBlockDataReader>(File.filePath());
}
else {
throw Base::FileException("Unknown extension");
throw Base::FileException("Unknown extension");
}
void FemPostPipeline::read(Base::FileInfo File)
{
Data.setValue(dataObjectFromFile(File));
}
void FemPostPipeline::read(std::vector<Base::FileInfo>& files, std::vector<double>& values, Base::Unit unit, std::string& frame_type)
{
if (files.size() != values.size()) {
Base::Console().Error("Result files and frame values have different length.\n");
return;
}
// setup the time information for the multiblock
vtkStringArray* TimeInfo = vtkStringArray::New();
TimeInfo->SetName("TimeInfo");
TimeInfo->InsertNextValue(frame_type);
TimeInfo->InsertNextValue(unit.getString());
auto multiblock = vtkSmartPointer<vtkMultiBlockDataSet>::New();
for (ulong i = 0; i < files.size(); i++) {
// add time information
vtkFloatArray* TimeValue = vtkFloatArray::New();
TimeValue->SetNumberOfComponents(1);
TimeValue->SetName("TimeValue");
TimeValue->InsertNextValue(values[i]);
// checking on the file
auto File = files[i];
if (!File.isReadable()) {
throw Base::FileException("File to load not existing or not readable", File);
}
auto data = dataObjectFromFile(File);
data->GetFieldData()->AddArray(TimeValue);
data->GetFieldData()->AddArray(TimeInfo);
multiblock->SetBlock(i, data);
}
multiblock->GetFieldData()->AddArray(TimeInfo);
Data.setValue(multiblock);
}
void FemPostPipeline::scale(double s)

View File

@@ -87,9 +87,10 @@ public:
return "FemGui::ViewProviderFemPostPipeline";
}
// load data from files
// load data from files (single or as multiframe)
static bool canRead(Base::FileInfo file);
void read(Base::FileInfo file);
void read(std::vector<Base::FileInfo>& files, std::vector<double>& values, Base::Unit unit, std::string& frame_type);
void scale(double s);
// load from results
@@ -126,14 +127,15 @@ private:
template<class TReader>
void readXMLFile(std::string file)
vtkSmartPointer<vtkDataObject> readXMLFile(std::string file)
{
vtkSmartPointer<TReader> reader = vtkSmartPointer<TReader>::New();
reader->SetFileName(file.c_str());
reader->Update();
Data.setValue(reader->GetOutput());
return reader->GetOutput();
}
vtkSmartPointer<vtkDataObject> dataObjectFromFile(Base::FileInfo File);
};
} // namespace Fem

View File

@@ -15,7 +15,16 @@
</Documentation>
<Methode Name="read">
<Documentation>
<UserDocu>Read in vtk file</UserDocu>
<UserDocu>
read(filepath)
read([filepaths], [values], unit, frame_type)
Reads in a single vtk file or creates a multiframe result by reading in multiple result files. If multiframe is wanted, 4 argumenhts are needed:
1. List of result files each being one frame,
2. List of values valid for each frame (e.g. [s] if time data),
3. the unit of the value as FreeCAD.Units.Unit,
4. the Description of the frame type
</UserDocu>
</Documentation>
</Methode>
<Methode Name="scale">
@@ -25,7 +34,16 @@
</Methode>
<Methode Name="load">
<Documentation>
<UserDocu>Load a single result object or create a multiframe result by loading multiple result frames. If multistep is wanted, 4 argumenhts are needed: 1. List of result object each being one frame, 2. List of values valid for each frame (e.g. [s] if time data), 3. the unit of the value, 4. the Description of the frames</UserDocu>
<UserDocu>
load(result_object)
load([result_objects], [values], unit, frame_type)
Load a single result object or create a multiframe result by loading multiple result frames. If multiframe is wanted, 4 argumenhts are needed:
1. List of result files each being one frame,
2. List of values valid for each frame (e.g. [s] if time data),
3. the unit of the value as FreeCAD.Units.Unit,
4. the Description of the frame type
</UserDocu>
</Documentation>
</Methode>
<Methode Name="getFilter">

View File

@@ -51,6 +51,76 @@ PyObject* FemPostPipelinePy::read(PyObject* args)
PyMem_Free(Name);
Py_Return;
}
PyObject *files;
PyObject *values = nullptr;
PyObject *unitobj = nullptr;
const char* value_type;
if (PyArg_ParseTuple(args, "O|OO!s", &files, &values, &(Base::UnitPy::Type), &unitobj, &value_type)) {
if (values == nullptr) {
// single argument version was called!
if (!PyUnicode_Check(files)) {
PyErr_SetString(PyExc_TypeError, "argument must be file path");
return nullptr;
}
const char* path = PyUnicode_AsUTF8(files);
getFemPostPipelinePtr()->read(Base::FileInfo(path));
}
else if (values != nullptr && unitobj != nullptr) {
//multistep version!
if ( !(PyTuple_Check(files) || PyList_Check(files)) ||
!(PyTuple_Check(values) || PyList_Check(values)) ) {
std::string error = std::string("Files and values must be list of strings and number respectively.");
throw Base::TypeError(error);
}
// extract the result objects
Py::Sequence file_list(files);
Py::Sequence::size_type size = file_list.size();
std::vector<Base::FileInfo> file_result;
file_result.resize(size);
for (Py::Sequence::size_type i = 0; i < size; i++) {
auto path = Py::Object(file_list[i]);
if (!path.isString()) {
throw Base::TypeError("File path must be string");
}
file_result[i] = Base::FileInfo(path.as_string());
}
//extract the values
Py::Sequence values_list(values);
size = values_list.size();
std::vector<double> value_result;
value_result.resize(size);
for (Py::Sequence::size_type i = 0; i < size; i++) {
auto value = Py::Object(values_list[i]);
if (!value.isNumeric()) {
std::string error = std::string("Values must be numbers");
throw Base::TypeError(error);
}
value_result[i] = Py::Float(value).as_double();
}
// extract the unit
Base::Unit unit = *(static_cast<Base::UnitPy*>(unitobj)->getUnitPtr());
// extract the value type
std::string step_type = std::string(value_type);
// Finally call the c++ function!
getFemPostPipelinePtr()->read(file_result, value_result, unit, step_type);
Py_Return;
}
}
return nullptr;
}

View File

@@ -152,8 +152,6 @@ class Proxy(solverbase.Proxy):
obj.addProperty("App::PropertyLink", "ElmerResult", "Base", "", 4 | 8)
obj.addProperty("App::PropertyLinkList", "ElmerTimeResults", "Base", "", 4 | 8)
obj.addProperty("App::PropertyLink", "ElmerOutput", "Base", "", 4 | 8)
obj.addProperty(

View File

@@ -279,64 +279,36 @@ class Results(run.Results):
# a line has the form like this:
# <DataSet timestep=" 5.000E-02" group="" part="0" file="FreeCAD_t0001.vtu"/>
# so .split("\"") gives as 2nd the time and as 7th the filename
files = []
values = []
for i in range(0, len(pvdContent) - 2):
# get time
lineArray = pvdContent[i + 1].split('"')
time = float(lineArray[1])
filename = os.path.join(self.directory, lineArray[7])
if os.path.isfile(filename):
self._createTimeResults(time, i + 1)
self.solver.ElmerTimeResults[i].read(filename)
# for eigen analyses the resulting values are by a factor 1000 to high
# therefore scale all *EigenMode results
self.solver.ElmerTimeResults[i].ViewObject.transformField(
"displacement EigenMode1", 0.001
)
self.solver.ElmerTimeResults[i].recomputeChildren()
# recompute() will update the result mesh data
# but not the shape and bar coloring
self.solver.ElmerTimeResults[i].ViewObject.updateColorBars()
values.append(time)
files.append(filename)
else:
self.pushStatus(f"\nResult file for time {time} is missing.\n")
self.fail()
return
if self.solver.ElmerResult is None:
self._createResults()
self.solver.ElmerResult.read(files, values, FreeCAD.Units.TimeSpan, "Time")
# for eigen analyses the resulting values are by a factor 1000 to high
# therefore scale all *EigenMode results
self.solver.ElmerResult.ViewObject.transformField("displacement EigenMode1", 0.001)
self.solver.ElmerResult.recomputeChildren()
self.solver.Document.recompute()
# recompute() updated the result mesh data
# but not the shape and bar coloring
self.solver.ElmerResult.ViewObject.updateColorBars()
def _createTimeResults(self, time, counter):
# if self.solver.ElmerTimeResults[counter] exists, but time is different
# recreate, other wise append
# FreeCAD would replaces dots in object names with underscores, thus do the same
newName = self.solver.Name + "_" + str(time).replace(".", "_") + "_" + "Result"
if counter > len(self.solver.ElmerTimeResults):
pipeline = self.analysis.Document.addObject("Fem::FemPostPipeline", newName)
# App::PropertyLinkList does not support append
# thus we have to use a temporary list to append
tmplist = self.solver.ElmerTimeResults
tmplist.append(pipeline)
self.solver.ElmerTimeResults = tmplist
self._finishTimeResults(time, counter - 1)
else:
# recreate if time is not equal
if self.solver.ElmerTimeResults[counter - 1].Name != newName:
# store current list before removing object since object removal will automatically
# remove entry from self.solver.ElmerTimeResults
tmplist = self.solver.ElmerTimeResults
self.analysis.Document.removeObject(self.solver.ElmerTimeResults[counter - 1].Name)
tmplist[counter - 1] = self.analysis.Document.addObject(
"Fem::FemPostPipeline", newName
)
self.solver.ElmerTimeResults = tmplist
self._finishTimeResults(time, counter - 1)
def _finishTimeResults(self, time, counter):
# we purposely use the decimal dot in the label
self.solver.ElmerTimeResults[counter].Label = f"{self.solver.Name}_{time}_Result"
self.solver.ElmerTimeResults[counter].ViewObject.OnTopWhenSelected = True
self.analysis.addObject(self.solver.ElmerTimeResults[counter])
# to assure the user sees something, set the default to Surface
self.solver.ElmerTimeResults[counter].ViewObject.DisplayMode = "Surface"
def _getResultFile(self):
postPath = None

View File

@@ -362,8 +362,6 @@ class Writer:
"Order of time stepping method 'BDF'",
)
solver.BDFOrder = (2, 1, 5, 1)
if not hasattr(self.solver, "ElmerTimeResults"):
solver.addProperty("App::PropertyLinkList", "ElmerTimeResults", "Base", "", 4 | 8)
if not hasattr(self.solver, "OutputIntervals"):
solver.addProperty(
"App::PropertyIntegerList",