diff --git a/dune/microstructure/prestrainedMaterial.hh b/dune/microstructure/prestrainedMaterial.hh index 3a274edcf111027ab41673537d10aeea589c6d9b..afd2dca574b83bf62a5e9fcb59882a711e922fbb 100644 --- a/dune/microstructure/prestrainedMaterial.hh +++ b/dune/microstructure/prestrainedMaterial.hh @@ -245,16 +245,6 @@ public: // } - - // ----------------------------------------------------------------- - // --- write material (grid)functions to VTK - void write_materialFunctions() - { - - - return; - - }; /////////////////////////////// // Getter /////////////////////////////// @@ -276,7 +266,120 @@ public: //getLameParameters() .. Throw Exception if isotropic_ = false! + /////////////////////////////// + // VTK - Writer + /////////////////////////////// + // ----------------------------------------------------------------- + // --- write material (grid)functions to VTK + void writeVTKMaterialFunctions(std::array<int, dim> nElements, const int level) + { + std::string outputPath = parameterSet_.get("outputPath", "../../outputs"); + using VTKGridType = YaspGrid<dim, EquidistantOffsetCoordinates<double, dim> >; + //VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},{80,80,80}); + VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},nElements); + + using GridViewVTK = VTKGridType::LeafGridView; + const GridViewVTK gridView_VTK = grid_VTK.leafGridView(); + + auto scalarP0FeBasis = makeBasis(gridView_VTK,lagrange<0>()); + auto scalarP1FeBasis = makeBasis(gridView_VTK,lagrange<1>()); + + std::vector<double> indicatorFunction_CoeffP0; + Functions::interpolate(scalarP0FeBasis, indicatorFunction_CoeffP0, indicatorFunction_); + auto indicatorFunction_P0 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP0FeBasis, indicatorFunction_CoeffP0); + + std::vector<double> indicatorFunction_CoeffP1; + Functions::interpolate(scalarP1FeBasis, indicatorFunction_CoeffP1, indicatorFunction_); + auto indicatorFunction_P1 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP1FeBasis, indicatorFunction_CoeffP1); + + VTKWriter<GridView> MaterialVtkWriter(gridView_VTK); + + MaterialVtkWriter.addVertexData( + indicatorFunction_P0, + VTK::FieldInfo("indicatorFunction_P0", VTK::FieldInfo::Type::scalar, 1)); + MaterialVtkWriter.addCellData( + indicatorFunction_P1, + VTK::FieldInfo("indicatorFunction_P1", VTK::FieldInfo::Type::scalar, 1)); + + MaterialVtkWriter.write(outputPath + "/MaterialFunctions-level"+ std::to_string(level) ); + std::cout << "wrote data to file:" + outputPath +"/MaterialFunctions-level" + std::to_string(level) << std::endl; + return; + }; + // ----------------------------------------------------------------- + // --- write prestrain (grid)functions to VTK + // void writeVTKPrestainFunctions(std::array<int, dim> nElements, const int level) + // { + // std::string outputPath = parameterSet_.get("outputPath", "../../outputs"); + // using VTKGridType = YaspGrid<dim, EquidistantOffsetCoordinates<double, dim> >; + // //VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},{80,80,80}); + // VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},nElements); + + // using GridViewVTK = VTKGridType::LeafGridView; + // const GridViewVTK gridView_VTK = grid_VTK.leafGridView(); + + // auto scalarP0FeBasis = makeBasis(gridView_VTK,lagrange<0>()); + // auto scalarP1FeBasis = makeBasis(gridView_VTK,lagrange<1>()); + + // std::vector<double> prestrainFunction_CoeffP0; + // Functions::interpolate(scalarP0FeBasis, prestrainFunction_CoeffP0, prestrain_); + // auto prestrainFunction_P0 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP0FeBasis, prestrainFunction_CoeffP0); + + // std::vector<double> prestrainFunction_CoeffP1; + // Functions::interpolate(scalarP1FeBasis, prestrainFunction_CoeffP1, prestrain_); + // auto prestrainFunction_P1 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP1FeBasis, prestrainFunction_CoeffP1); + + // VTKWriter<GridView> MaterialVtkWriter(gridView_VTK); + + // MaterialVtkWriter.addVertexData( + // prestrainFunction_P0, + // VTK::FieldInfo("prestrainFunction_P0", VTK::FieldInfo::Type::scalar, 1)); + // MaterialVtkWriter.addCellData( + // prestrainFunction_P1, + // VTK::FieldInfo("prestrainFunction_P1", VTK::FieldInfo::Type::scalar, 1)); + + // MaterialVtkWriter.write(outputPath + "/PrestrainFunctions-level"+ std::to_string(level) ); + // std::cout << "wrote data to file:" + outputPath + "/PrestrainFunctions-level" + std::to_string(level) << std::endl; + // return; + // }; + + // void writeVTKPrestainFunctions(std::array<int, dim> nElements, const int level) + // { +// using VTKGridType = YaspGrid<dim, EquidistantOffsetCoordinates<double, dim> >; +// // VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},{80,80,80}); +// // VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},{40,40,40}); +// VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},nElements); +// using GridViewVTK = VTKGridType::LeafGridView; +// const GridViewVTK gridView_VTK = grid_VTK.leafGridView(); + +// FTKfillerContainer<dim> VTKFiller; +// VTKFiller.vtkPrestrainNorm(gridView_VTK, B_Term, "PrestrainBNorm"); + +// // WORKS Too +// VTKFiller.vtkProblemCell(gridView_VTK, B_Term, muLocal,"VTKProblemCell");; + + +// // TEST +// auto scalarP0FeBasis = makeBasis(gridView_VTK,lagrange<0>()); +// auto scalarP1FeBasis = makeBasis(gridView_VTK,lagrange<1>()); + +// std::vector<double> B_CoeffP0; +// Functions::interpolate(scalarP0FeBasis, B_CoeffP0, B_Term); +// auto B_DGBF_P0 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP0FeBasis, B_CoeffP0); + + + +// VTKWriter<GridView> PrestrainVtkWriter(gridView_VTK); + +// PrestrainVtkWriter.addCellData( +// B_DGBF_P0, +// VTK::FieldInfo("B_P0", VTK::FieldInfo::Type::scalar, 1)); + +// PrestrainVtkWriter.write(outputPath + "/PrestrainFunctions-level"+ std::to_string(level) ); +// std::cout << "wrote data to file:" + outputPath +"/PrestrainFunctions-level" + std::to_string(level) << std::endl; + +// // }; + }; // end class diff --git a/src/Cell-Problem-New.cc b/src/Cell-Problem-New.cc index 22cd8d035574338753cf07e3b608db984a919e84..ad05b62ea2d77c6fe18c34994767d5aff60f2683 100644 --- a/src/Cell-Problem-New.cc +++ b/src/Cell-Problem-New.cc @@ -570,11 +570,6 @@ int main(int argc, char *argv[]) auto correctorComputer = CorrectorComputer(Basis_CE, material_, muTerm, lambdaTerm, gamma, log, parameterSet); correctorComputer.solve(); - - -////////////////////////////////////////////////// - - //--- check Correctors (options): if(parameterSet.get<bool>("write_L2Error", false)) correctorComputer.computeNorms(); @@ -595,6 +590,11 @@ int main(int argc, char *argv[]) auto effectiveQuantitiesComputer = EffectiveQuantitiesComputer(correctorComputer,material_); effectiveQuantitiesComputer.computeEffectiveQuantities(); + //--- write material indicator function to VTK + if (write_materialFunctions) + { + material_.writeVTKMaterialFunctions(nElements,level); + } @@ -652,95 +652,6 @@ int main(int argc, char *argv[]) - - if (write_materialFunctions) - { - using VTKGridType = YaspGrid<dim, EquidistantOffsetCoordinates<double, dim> >; -// VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},{80,80,80}); -// VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},{40,40,40}); - VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},nElements); - - using GridViewVTK = VTKGridType::LeafGridView; - const GridViewVTK gridView_VTK = grid_VTK.leafGridView(); - - auto scalarP0FeBasis = makeBasis(gridView_VTK,lagrange<0>()); - auto scalarP1FeBasis = makeBasis(gridView_VTK,lagrange<1>()); - - std::vector<double> mu_CoeffP0; - Functions::interpolate(scalarP0FeBasis, mu_CoeffP0, muTerm); - auto mu_DGBF_P0 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP0FeBasis, mu_CoeffP0); - - std::vector<double> mu_CoeffP1; - Functions::interpolate(scalarP1FeBasis, mu_CoeffP1, muTerm); - auto mu_DGBF_P1 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP1FeBasis, mu_CoeffP1); - - - std::vector<double> lambda_CoeffP0; - Functions::interpolate(scalarP0FeBasis, lambda_CoeffP0, lambdaTerm); - auto lambda_DGBF_P0 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP0FeBasis, lambda_CoeffP0); - - std::vector<double> lambda_CoeffP1; - Functions::interpolate(scalarP1FeBasis, lambda_CoeffP1, lambdaTerm); - auto lambda_DGBF_P1 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP1FeBasis, lambda_CoeffP1); - - VTKWriter<GridView> MaterialVtkWriter(gridView_VTK); - - MaterialVtkWriter.addVertexData( - mu_DGBF_P1, - VTK::FieldInfo("mu_P1", VTK::FieldInfo::Type::scalar, 1)); - MaterialVtkWriter.addCellData( - mu_DGBF_P0, - VTK::FieldInfo("mu_P0", VTK::FieldInfo::Type::scalar, 1)); - MaterialVtkWriter.addVertexData( - lambda_DGBF_P1, - VTK::FieldInfo("lambda_P1", VTK::FieldInfo::Type::scalar, 1)); - MaterialVtkWriter.addCellData( - lambda_DGBF_P0, - VTK::FieldInfo("lambda_P0", VTK::FieldInfo::Type::scalar, 1)); - - MaterialVtkWriter.write(outputPath + "/MaterialFunctions-level"+ std::to_string(level) ); - std::cout << "wrote data to file:" + outputPath +"/MaterialFunctions-level" + std::to_string(level) << std::endl; - - } - - -// if (write_prestrainFunctions) -// { -// using VTKGridType = YaspGrid<dim, EquidistantOffsetCoordinates<double, dim> >; -// // VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},{80,80,80}); -// // VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},{40,40,40}); -// VTKGridType grid_VTK({-1.0/2.0, -1.0/2.0, -1.0/2.0},{1.0/2.0, 1.0/2.0, 1.0/2.0},nElements); -// using GridViewVTK = VTKGridType::LeafGridView; -// const GridViewVTK gridView_VTK = grid_VTK.leafGridView(); - -// FTKfillerContainer<dim> VTKFiller; -// VTKFiller.vtkPrestrainNorm(gridView_VTK, B_Term, "PrestrainBNorm"); - -// // WORKS Too -// VTKFiller.vtkProblemCell(gridView_VTK, B_Term, muLocal,"VTKProblemCell");; - - -// // TEST -// auto scalarP0FeBasis = makeBasis(gridView_VTK,lagrange<0>()); -// auto scalarP1FeBasis = makeBasis(gridView_VTK,lagrange<1>()); - -// std::vector<double> B_CoeffP0; -// Functions::interpolate(scalarP0FeBasis, B_CoeffP0, B_Term); -// auto B_DGBF_P0 = Functions::makeDiscreteGlobalBasisFunction<double>(scalarP0FeBasis, B_CoeffP0); - - - -// VTKWriter<GridView> PrestrainVtkWriter(gridView_VTK); - -// PrestrainVtkWriter.addCellData( -// B_DGBF_P0, -// VTK::FieldInfo("B_P0", VTK::FieldInfo::Type::scalar, 1)); - -// PrestrainVtkWriter.write(outputPath + "/PrestrainFunctions-level"+ std::to_string(level) ); -// std::cout << "wrote data to file:" + outputPath +"/PrestrainFunctions-level" + std::to_string(level) << std::endl; - -// } -