diff --git a/applications/FluidDynamicsHydraulicsApplication/custom_python/add_custom_utilities_to_python.cpp b/applications/FluidDynamicsHydraulicsApplication/custom_python/add_custom_utilities_to_python.cpp index 3ae187d9b58a..ce04859882be 100644 --- a/applications/FluidDynamicsHydraulicsApplication/custom_python/add_custom_utilities_to_python.cpp +++ b/applications/FluidDynamicsHydraulicsApplication/custom_python/add_custom_utilities_to_python.cpp @@ -26,12 +26,26 @@ void AddCustomUtilitiesToPython(pybind11::module& m) { namespace py = pybind11; py::class_(m, "HydraulicFluidAuxiliaryUtilities") - .def_static("CalculateWettedPetimeter", [](ModelPart &rModelPart, const Flags &rSkinFlag, const Variable &rDistanceVariable, const bool IsHistorical){return HydraulicFluidAuxiliaryUtilities::CalculateWettedPetimeter(rModelPart, rSkinFlag, rDistanceVariable, IsHistorical);}) - .def_static("CalculateWettedArea", [](ModelPart &rModelPart, const Flags &rSkinFlag, const Variable &rDistanceVariable, const bool IsHistorical){return HydraulicFluidAuxiliaryUtilities::CalculateWettedArea(rModelPart, rSkinFlag, rDistanceVariable, IsHistorical);}) - .def_static("InitialWaterDepth", [](ModelPart &rModelPart){return HydraulicFluidAuxiliaryUtilities::InitialWaterDepth(rModelPart);}) - .def_static("SetInletVelocity", [](ModelPart &rModelPart, double InletVelocity, const Variable &rDistanceVariable){return HydraulicFluidAuxiliaryUtilities::SetInletVelocity(rModelPart, InletVelocity, rDistanceVariable);}) - .def_static("FreeInlet", [](ModelPart &rModelPart){return HydraulicFluidAuxiliaryUtilities::FreeInlet(rModelPart);}) - .def_static("SetInletFreeSurface", [](ModelPart &rModelPart, const Flags &rSkinFlag, const Variable &rDistanceVariable){return HydraulicFluidAuxiliaryUtilities::SetInletFreeSurface(rModelPart, rSkinFlag, rDistanceVariable);}); + .def_static("CalculateWettedPetimeter", [](ModelPart &rModelPart, const Flags &rSkinFlag, const Variable &rDistanceVariable, const bool IsHistorical) + { return HydraulicFluidAuxiliaryUtilities::CalculateWettedPetimeter(rModelPart, rSkinFlag, rDistanceVariable, IsHistorical); }) + .def_static("CalculateWettedArea", [](ModelPart &rModelPart, const Flags &rSkinFlag, const Variable &rDistanceVariable, const bool IsHistorical) + { return HydraulicFluidAuxiliaryUtilities::CalculateWettedArea(rModelPart, rSkinFlag, rDistanceVariable, IsHistorical); }) + .def_static("InitialWaterDepth", [](ModelPart &rModelPart) + { return HydraulicFluidAuxiliaryUtilities::InitialWaterDepth(rModelPart); }) + .def_static("SetInletVelocity", [](ModelPart &rModelPart, double InletVelocity, const Variable &rDistanceVariable) + { return HydraulicFluidAuxiliaryUtilities::SetInletVelocity(rModelPart, InletVelocity, rDistanceVariable); }) + .def_static("FreeInlet", [](ModelPart &rModelPart) + { return HydraulicFluidAuxiliaryUtilities::FreeInlet(rModelPart); }) + .def_static("SetInletFreeSurface", [](ModelPart &rModelPart, const Flags &rSkinFlag, const Variable &rDistanceVariable) + { return HydraulicFluidAuxiliaryUtilities::SetInletFreeSurface(rModelPart, rSkinFlag, rDistanceVariable); }) + .def_static("FixCornerNodeVelocity", [](ModelPart &rModelPart, double MaximumAngle) + { return HydraulicFluidAuxiliaryUtilities::FixCornerNodeVelocity(rModelPart, MaximumAngle); }) + .def_static("TurnOffGravityOnAirElements", [](ModelPart &rModelPart) + { return HydraulicFluidAuxiliaryUtilities::TurnOffGravityOnAirElements(rModelPart); }) + .def_static("MaximumWaterDepthChange", [](ModelPart &rModelPart) + { return HydraulicFluidAuxiliaryUtilities::MaximumWaterDepthChange(rModelPart); }) + .def_static("CalculateArtificialViscosity", [](ModelPart &rModelPart, double LimiterCoefficient) + { return HydraulicFluidAuxiliaryUtilities::CalculateArtificialViscosity(rModelPart, LimiterCoefficient); }); } } // namespace Kratos::Python diff --git a/applications/FluidDynamicsHydraulicsApplication/custom_utilities/hydraulic_fluid_auxiliary_utilities.cpp b/applications/FluidDynamicsHydraulicsApplication/custom_utilities/hydraulic_fluid_auxiliary_utilities.cpp index e108a434f851..229bd5416285 100644 --- a/applications/FluidDynamicsHydraulicsApplication/custom_utilities/hydraulic_fluid_auxiliary_utilities.cpp +++ b/applications/FluidDynamicsHydraulicsApplication/custom_utilities/hydraulic_fluid_auxiliary_utilities.cpp @@ -302,9 +302,143 @@ void HydraulicFluidAuxiliaryUtilities::SetInletFreeSurface(ModelPart &rModelPart if (rNode.Is(rSkinFlag)){ double inlet_dist = rNode.GetValue(rDistanceVariable); rNode.FastGetSolutionStepValue(DISTANCE) = inlet_dist; - rNode.Fix(DISTANCE); + // The distance is fixed in the water nodes since the inlet has velocity. + if (inlet_dist<0.0){ + rNode.Fix(DISTANCE); + } + } + }); +} + +void HydraulicFluidAuxiliaryUtilities::TurnOffGravityOnAirElements(ModelPart &rModelPart){ + block_for_each(rModelPart.Elements(), [](Element &rElement){ + auto &r_geom = rElement.GetGeometry(); + Vector distances = ZeroVector(r_geom.PointsNumber()); + const Vector gravity = ZeroVector(3); + for (unsigned int i_nodes = 0; i_nodes < r_geom.PointsNumber(); i_nodes++){ + distances[i_nodes] = r_geom[i_nodes].FastGetSolutionStepValue(DISTANCE); + } + if (FluidAuxiliaryUtilities::IsPositive(distances)){ + for (unsigned int i_nodes = 0; i_nodes < r_geom.PointsNumber(); i_nodes++) + { + r_geom[i_nodes].FastGetSolutionStepValue(BODY_FORCE) = gravity; + } + } + }); +} + +void HydraulicFluidAuxiliaryUtilities::FixCornerNodeVelocity( + ModelPart &rModelPart, + const double MaximumAngle) +{ + // Obtain for each condition each neighbor condition + block_for_each(rModelPart.Nodes(), [&](NodeType &rNode){ rNode.SetValue(AUX_INDEX, 0); }); + auto &r_cond = rModelPart.Conditions(); + const double angle_corner_rad = MaximumAngle / 180.0 * Globals::Pi; // + double acceptable_cos = cos(angle_corner_rad); + double max_height = block_for_each>(rModelPart.Nodes(), [&](NodeType &rNode){return rNode.Z();}); + + block_for_each(rModelPart.Nodes(),[&](NodeType &rNode){ + if(rNode.Z() >=max_height){ + rNode.Set(MARKER,true); + } + }); + for (auto cond_it = r_cond.begin(); cond_it != r_cond.end(); cond_it++) + { + // reference for area normal of the face + const auto &r_face_normal = cond_it->GetValue(NORMAL); + const double An = norm_2(face_normal); + const auto &r_neighb = cond_it->GetValue(NEIGHBOUR_CONDITIONS); + const auto &r_geom = cond_it->GetGeometry(); + auto edgelist = r_geom.GenerateEdges(); + for (unsigned int c_itr = 0; c_itr < neighb.size(); c_itr++){ + auto rCondition = neighb[c_itr]; + const auto& r_face_neig_normal = rCondition.GetValue(NORMAL); + const double cos_angle = 1 / (An * norm_2(face_neig_normal)) * inner_prod(face_normal, face_neig_normal); + const auto &r_geom_neig = rCondition.GetGeometry(); + auto edgelist_neig = r_geom_neig.GenerateEdges(); + array_1d ids_neig; + array_1d ids_cond; + if (cos_angle < acceptable_cos){ + for (IndexType edge = 0; edge < edgelist.size(); edge++){ + for (IndexType i = 0; i < edgelist[edge].size(); i++){ + ids_cond[i] = edgelist[edge][i].Id(); + } + std::sort(ids_cond.begin(), ids_cond.end()); + for (IndexType edge_neg = 0; edge_neg < edgelist_neig.size(); edge_neg++){ + for (IndexType i = 0; i < edgelist_neig[edge_neg].size(); i++){ + ids_neig[i] = edgelist_neig[edge_neg][i].Id(); + } + std::sort(ids_neig.begin(), ids_neig.end()); + if (ids_cond == ids_neig){ + for (IndexType i = 0; i < edgelist[edge].size(); i++){ + double repeated = edgelist[edge][i].GetValue(AUX_INDEX); + repeated += 1; + edgelist[edge][i].SetValue(AUX_INDEX, repeated); + } + } + } + } + } + } + } + block_for_each(rModelPart.Nodes(),[&](NodeType &rNode){ + if (rNode.GetValue(AUX_INDEX)*0.5>=3 && rNode.Is(MARKER)){ + Vector veloctiy = ZeroVector(3); + rNode.FastGetSolutionStepValue(VELOCITY) = veloctiy; + rNode.Fix(VELOCITY_X); + rNode.Fix(VELOCITY_Y); + rNode.Fix(VELOCITY_Z); + double value= rNode.GetValue(AUX_INDEX); + rNode.SetValue(AUX_INDEX,value*0.5); } + else{ + rNode.SetValue(AUX_INDEX,0.0); + } + }); +} + +bool HydraulicFluidAuxiliaryUtilities::MaximumWaterDepthChange(const ModelPart &rModelPart){ + // make sure the Ids of the Nodes start with one + const double &min_Z = block_for_each>(rModelPart.Nodes(), [&](const auto &rNode){ return rNode.Z(); }); + //TODO: Try to create a min reduction with a pointer to avoid this second loop + bool maximum_water_depth_change = false; + block_for_each(rModelPart.Nodes(), [&](const Node &rNode){ + if (std::abs(rNode.Z() - min_Z) < 1.0e-10 && !maximum_water_depth_change){ + const double inlet_phi = rNode.GetValue(AUX_DISTANCE); + const double domain_phi = rNode.FastGetSolutionStepValue(DISTANCE); + if (std::abs(domain_phi) > std::abs(inlet_phi)){ + maximum_water_depth_change=true; + } + }}); + return maximum_water_depth_change; +} + +void HydraulicFluidAuxiliaryUtilities::CalculateArtificialViscosity(ModelPart &rModelPart, double LimiterCoefficient) +{ + const auto &r_process_info = rModelPart.GetProcessInfo(); + block_for_each(rModelPart.Elements(), [&](Element &rElement){ + double artificial_viscosity; + rElement.Calculate(ARTIFICIAL_DYNAMIC_VISCOSITY,artificial_viscosity ,r_process_info); + if (artificial_viscosity > LimiterCoefficient){ + artificial_viscosity = LimiterCoefficient; + } + double neg_nodes = 0.0; + double pos_nodes = 0.0; + for (const auto &r_node : rElement.GetGeometry()){ + const double distance = r_node.FastGetSolutionStepValue(DISTANCE); + if (distance > 0){ + pos_nodes += 1; + } + else{ + neg_nodes += 1; + } + } + if (neg_nodes > 0 && pos_nodes > 0){ + artificial_viscosity = 0.0; + } + rElement.SetValue(ARTIFICIAL_DYNAMIC_VISCOSITY, artificial_viscosity); }); } -} // namespace Kratos +} // namespace Kratos \ No newline at end of file diff --git a/applications/FluidDynamicsHydraulicsApplication/custom_utilities/hydraulic_fluid_auxiliary_utilities.h b/applications/FluidDynamicsHydraulicsApplication/custom_utilities/hydraulic_fluid_auxiliary_utilities.h index ad176c74bccb..779ddb6ea747 100644 --- a/applications/FluidDynamicsHydraulicsApplication/custom_utilities/hydraulic_fluid_auxiliary_utilities.h +++ b/applications/FluidDynamicsHydraulicsApplication/custom_utilities/hydraulic_fluid_auxiliary_utilities.h @@ -52,6 +52,11 @@ class KRATOS_API(FLUID_DYNAMICS_HYDRAULICS_APPLICATION) HydraulicFluidAuxiliaryU using PointsArrayType = typename GeometryType::PointsArrayType; using ModifiedShapeFunctionsFactoryType = std::function; + struct NodeConditionDataStruct + { + NodeType Node; // Condition normal + double CornerNode; // Gauss point shape functions values + }; ///@} ///@name Static Operations @@ -116,18 +121,49 @@ class KRATOS_API(FLUID_DYNAMICS_HYDRAULICS_APPLICATION) HydraulicFluidAuxiliaryU * @param rDistancesVariable Variable name of the inlet distance. */ static void SetInletFreeSurface(ModelPart &rModelPart, const Flags &rSkinFlag, const Variable &rDistanceVariable); - - ///@} + /** + * @brief For all elements that are completly air the gravity is turned off. + * @param rModelPart FluidModelPart + */ + static void TurnOffGravityOnAirElements(ModelPart &rModelPart); + + /** + * @brief Fixing the velocity to corner nodes in the skin modelpart. + * @param rModelPart FluidModelPart + * @param MaximumAngle Minimum angle for considering corner. + */ + + static void FixCornerNodeVelocity( + ModelPart &rModelPart, + const double MaximumAngle); -private : + /** + * @brief The maximum water depth at the inlet is adjusted based on the maximum depth within the entire domain. + * If the water depth at the inlet is less than the maximum depth across the domain, the inlet water depth is adjusted to match its maximum depth. + * @param rModelPart Inlet Model Part + + */ + static bool MaximumWaterDepthChange(const ModelPart &rModelPart); + /** + * @brief Calculate the artificial viscosity which is proportional to the maximum velocity gradient + * @param rModelPart Fluid Model Part + * @param LimiterCoefficient a coefficient to limit the maximum artificial viscosity + * + */ + static void CalculateArtificialViscosity( + ModelPart &rModelPart, + const double LimiterCoefficient); + ///@} + +private: struct EdgeDataContainer { NodeType::Pointer pNodeI = nullptr; NodeType::Pointer pNodeJ = nullptr; SizeType NumberOfRepetitions = 0; - }; + }; -}; // namespace Kratos + }; // namespace Kratos } diff --git a/applications/FluidDynamicsHydraulicsApplication/tests/channel_test/channel3D.mdpa b/applications/FluidDynamicsHydraulicsApplication/tests/channel_test/channel3D.mdpa new file mode 100644 index 000000000000..59d8e6be4622 --- /dev/null +++ b/applications/FluidDynamicsHydraulicsApplication/tests/channel_test/channel3D.mdpa @@ -0,0 +1,855 @@ +Begin ModelPartData +// VARIABLE_NAME value +End ModelPartData + +Begin Properties 0 +End Properties + +Begin Nodes + 1 1.7070000000 0.7070000000 0.5000000000 + 2 1.5349366881 0.5349366881 0.5000000000 + 3 1.7070000000 0.7070000000 0.2500000000 + 4 1.5302233047 0.8837766953 0.5000000000 + 5 1.5302233047 0.8837766953 0.2500000000 + 6 1.2811487999 0.6085142866 0.5000000000 + 7 1.3125573443 0.7206386397 0.2718671235 + 8 1.4448315433 0.4448315433 0.2314823555 + 9 1.3547263984 0.3547263984 0.5000000000 + 10 1.7070000000 0.7070000000 0.0000000000 + 11 1.3534466094 1.0605533906 0.5000000000 + 12 1.5349366881 0.5349366881 0.0000000000 + 13 1.5302233047 0.8837766953 0.0000000000 + 14 1.3534466094 1.0605533906 0.2500000000 + 15 1.1672288491 0.8743356303 0.5000000000 + 16 1.2200622166 0.5158682548 0.2068847863 + 17 1.2811487999 0.6085142866 0.0000000000 + 18 1.0981855456 0.4291021300 0.5000000000 + 19 1.2636438378 0.2636438378 0.2339933608 + 20 1.0731147078 0.7802214890 0.2417815670 + 21 1.3547263984 0.3547263984 0.0000000000 + 22 1.3534466094 1.0605533906 0.0000000000 + 23 0.9790005665 0.6861073477 0.5000000000 + 24 1.1672288491 0.8743356303 0.0000000000 + 25 1.1725612773 0.1725612773 0.5000000000 + 26 1.0022756219 0.4198280899 0.2230270649 + 27 1.0981855456 0.4291021300 0.0000000000 + 28 0.9790005665 0.6861073477 0.0000000000 + 29 1.1725612773 0.1725612773 0.0000000000 + 30 0.7928932188 0.5000000000 0.5000000000 + 31 0.8783159734 0.2210481469 0.5000000000 + 32 0.7928932188 0.5000000000 0.2500000000 + 33 1.0000000000 0.0000000000 0.5000000000 + 34 1.0000000000 0.0000000000 0.2500000000 + 35 0.7928932188 0.5000000000 0.0000000000 + 36 0.7434364093 0.2743859813 0.2893260368 + 37 0.8783159734 0.2210481469 0.0000000000 + 38 1.0000000000 0.0000000000 0.0000000000 + 39 0.7566319468 0.0000000000 0.5000000000 + 40 0.6143779241 0.2223257049 0.5000000000 + 41 0.5294945785 0.5000000000 0.5000000000 + 42 0.5612593978 0.3500973463 0.3021092157 + 43 0.7566319468 0.0000000000 0.0000000000 + 44 0.6143779241 0.2223257049 0.0000000000 + 45 0.5294945785 0.5000000000 0.0000000000 + 46 0.6291810510 0.0000000000 0.2315240966 + 47 0.3963679286 0.5000000000 0.2418345289 + 48 0.4961022788 0.1829756047 0.1912173316 + 49 0.5017301552 0.0000000000 0.5000000000 + 50 0.3729003501 0.2340289887 0.5000000000 + 51 0.2632412788 0.5000000000 0.5000000000 + 52 0.5017301552 0.0000000000 0.0000000000 + 53 0.3729003501 0.2340289887 0.0000000000 + 54 0.3729003501 0.0000000000 0.2340289887 + 55 0.2632412788 0.5000000000 0.0000000000 + 56 0.2397293112 0.2666791196 0.2795525567 + 57 0.2440705450 0.0000000000 0.5000000000 + 58 0.2440705450 0.0000000000 0.0000000000 + 59 0.0000000000 0.5000000000 0.5000000000 + 60 0.0000000000 0.5000000000 0.2500000000 + 61 0.0000000000 0.2500000000 0.5000000000 + 62 0.0000000000 0.2500000000 0.2500000000 + 63 0.0000000000 0.5000000000 0.0000000000 + 64 0.0000000000 0.2500000000 0.0000000000 + 65 0.0000000000 0.0000000000 0.5000000000 + 66 0.0000000000 0.0000000000 0.2500000000 + 67 0.0000000000 0.0000000000 0.0000000000 +End Nodes + + +Begin Elements Element3D4N// GUI group identifier: Fluid + 1 0 43 38 37 34 + 2 0 39 31 33 34 + 3 0 38 29 37 34 + 4 0 58 54 52 53 + 5 0 57 50 49 54 + 6 0 52 44 53 48 + 7 0 44 53 48 45 + 8 0 52 44 48 46 + 9 0 48 52 46 54 + 10 0 46 48 54 49 + 11 0 48 52 54 53 + 12 0 44 48 46 36 + 13 0 48 46 36 40 + 14 0 54 48 53 56 + 15 0 54 48 56 50 + 16 0 21 8 12 17 + 17 0 9 6 2 8 + 18 0 49 50 40 48 + 19 0 44 48 36 45 + 20 0 33 31 25 34 + 21 0 43 37 44 36 + 22 0 52 46 43 44 + 23 0 49 40 39 46 + 24 0 29 21 27 19 + 25 0 21 19 8 16 + 26 0 21 19 16 27 + 27 0 16 21 27 17 + 28 0 27 16 17 28 + 29 0 16 21 17 8 + 30 0 19 8 16 9 + 31 0 19 16 27 26 + 32 0 16 27 26 28 + 33 0 17 16 8 7 + 34 0 19 16 26 18 + 35 0 17 16 7 20 + 36 0 9 18 6 16 + 37 0 25 18 9 19 + 38 0 39 40 31 36 + 39 0 67 58 64 66 + 40 0 55 45 47 53 + 41 0 48 53 56 47 + 42 0 12 3 10 13 + 43 0 16 7 20 6 + 44 0 57 66 65 61 + 45 0 48 54 49 50 + 46 0 40 31 36 30 + 47 0 2 4 1 3 + 48 0 16 8 7 6 + 49 0 48 56 50 47 + 50 0 46 48 49 40 + 51 0 28 24 20 17 + 52 0 50 41 40 42 + 53 0 40 50 42 48 + 54 0 50 42 48 47 + 55 0 42 40 48 36 + 56 0 48 42 36 45 + 57 0 42 40 36 41 + 58 0 42 48 47 45 + 59 0 47 42 45 32 + 60 0 13 22 24 14 + 61 0 41 51 47 50 + 62 0 29 27 37 26 + 63 0 14 11 15 4 + 64 0 62 61 66 56 + 65 0 62 61 56 60 + 66 0 8 12 17 7 + 67 0 44 37 35 36 + 68 0 5 13 3 7 + 69 0 46 36 40 39 + 70 0 36 42 41 32 + 71 0 36 42 32 45 + 72 0 35 32 45 36 + 73 0 35 32 36 26 + 74 0 32 36 26 30 + 75 0 26 32 30 23 + 76 0 35 32 26 28 + 77 0 32 36 30 41 + 78 0 30 40 41 36 + 79 0 46 43 44 36 + 80 0 28 35 27 26 + 81 0 63 60 64 55 + 82 0 58 54 53 56 + 83 0 59 51 61 60 + 84 0 5 4 14 7 + 85 0 5 4 7 3 + 86 0 14 5 7 13 + 87 0 62 64 60 56 + 88 0 56 54 50 57 + 89 0 50 56 57 61 + 90 0 18 23 6 16 + 91 0 17 16 20 28 + 92 0 25 31 18 26 + 93 0 31 18 26 30 + 94 0 50 42 47 41 + 95 0 42 47 41 32 + 96 0 27 29 19 26 + 97 0 66 64 62 56 + 98 0 56 50 47 51 + 99 0 56 50 51 61 + 100 0 47 56 51 60 + 101 0 35 45 44 36 + 102 0 15 23 20 6 + 103 0 23 18 30 26 + 104 0 16 26 18 23 + 105 0 7 20 6 15 + 106 0 7 20 15 14 + 107 0 6 7 15 4 + 108 0 35 37 27 26 + 109 0 25 18 19 26 + 110 0 53 56 47 55 + 111 0 56 47 55 60 + 112 0 31 39 36 34 + 113 0 46 36 39 34 + 114 0 9 6 8 16 + 115 0 8 7 6 2 + 116 0 7 6 2 4 + 117 0 18 9 19 16 + 118 0 24 20 17 7 + 119 0 24 20 7 14 + 120 0 17 24 7 13 + 121 0 56 54 57 66 + 122 0 16 20 28 26 + 123 0 46 43 36 34 + 124 0 17 12 13 7 + 125 0 37 43 34 36 + 126 0 32 26 28 20 + 127 0 8 3 12 7 + 128 0 8 7 2 3 + 129 0 19 25 26 34 + 130 0 32 26 20 23 + 131 0 31 36 30 26 + 132 0 31 36 26 34 + 133 0 29 37 34 26 + 134 0 15 7 14 4 + 135 0 58 54 56 66 + 136 0 31 25 34 26 + 137 0 45 47 53 48 + 138 0 35 37 26 36 + 139 0 37 26 36 34 + 140 0 58 53 64 56 + 141 0 16 20 26 23 + 142 0 29 34 19 26 + 143 0 23 20 6 16 + 144 0 53 55 64 56 + 145 0 7 14 13 24 + 146 0 12 3 13 7 + 147 0 60 64 55 56 + 148 0 56 51 60 61 + 149 0 61 66 56 57 + 150 0 2 7 4 3 + 151 0 58 64 66 56 +End Elements + +Begin Conditions SurfaceCondition3D3N// GUI group identifier: _HIDDEN__SKIN_ + 1 0 67 58 64 + 2 0 43 38 37 + 3 0 38 29 37 + 4 0 12 10 13 + 5 0 13 22 24 + 6 0 63 64 55 + 7 0 43 37 44 + 8 0 58 52 53 + 9 0 52 43 44 + 10 0 52 44 53 + 11 0 53 44 45 + 12 0 29 21 27 + 13 0 21 12 17 + 14 0 28 35 27 + 15 0 35 45 44 + 16 0 44 37 35 + 17 0 17 12 13 + 18 0 21 17 27 + 19 0 27 17 28 + 20 0 24 28 17 + 21 0 45 55 53 + 22 0 58 53 64 + 23 0 29 27 37 + 24 0 35 37 27 + 25 0 24 17 13 + 26 0 53 55 64 + 27 0 67 58 66 + 28 0 43 38 34 + 29 0 33 39 34 + 30 0 57 65 66 + 31 0 58 52 54 + 32 0 52 43 46 + 33 0 39 49 46 + 34 0 49 57 54 + 35 0 46 43 34 + 36 0 52 46 54 + 37 0 54 46 49 + 38 0 58 54 66 + 39 0 39 46 34 + 40 0 54 57 66 + 41 0 38 29 34 + 42 0 12 10 3 + 43 0 1 2 3 + 44 0 25 33 34 + 45 0 29 21 19 + 46 0 21 12 8 + 47 0 2 9 8 + 48 0 9 25 19 + 49 0 8 12 3 + 50 0 21 8 19 + 51 0 19 8 9 + 52 0 29 19 34 + 53 0 2 8 3 + 54 0 19 25 34 + 55 0 35 32 45 + 56 0 32 30 41 + 57 0 59 60 51 + 58 0 60 63 55 + 59 0 55 45 47 + 60 0 41 51 47 + 61 0 47 45 32 + 62 0 55 47 60 + 63 0 47 51 60 + 64 0 41 47 32 + 65 0 22 14 24 + 66 0 14 11 15 + 67 0 30 32 23 + 68 0 32 35 28 + 69 0 28 24 20 + 70 0 15 23 20 + 71 0 20 24 14 + 72 0 28 20 32 + 73 0 20 23 32 + 74 0 15 20 14 + 75 0 67 64 66 + 76 0 63 60 64 + 77 0 59 61 60 + 78 0 65 66 61 + 79 0 66 64 62 + 80 0 62 64 60 + 81 0 62 60 61 + 82 0 62 61 66 + 83 0 22 13 14 + 84 0 10 3 13 + 85 0 1 4 3 + 86 0 11 14 4 + 87 0 14 13 5 + 88 0 5 13 3 + 89 0 5 3 4 + 90 0 5 4 14 + 91 0 65 57 61 + 92 0 39 33 31 + 93 0 33 25 31 + 94 0 2 1 4 + 95 0 4 11 15 + 96 0 59 61 51 + 97 0 39 31 40 + 98 0 57 49 50 + 99 0 49 39 40 + 100 0 49 40 50 + 101 0 50 40 41 + 102 0 25 9 18 + 103 0 9 2 6 + 104 0 23 30 18 + 105 0 30 41 40 + 106 0 40 31 30 + 107 0 6 2 4 + 108 0 9 6 18 + 109 0 18 6 23 + 110 0 15 23 6 + 111 0 41 51 50 + 112 0 57 50 61 + 113 0 25 18 31 + 114 0 30 31 18 + 115 0 15 6 4 + 116 0 50 51 61 +End Conditions + +Begin SubModelPart FluidParts_Fluid // Group Fluid // Subtree FluidParts + Begin SubModelPartNodes + 1 + 2 + 3 + 4 + 5 + 6 + 7 + 8 + 9 + 10 + 11 + 12 + 13 + 14 + 15 + 16 + 17 + 18 + 19 + 20 + 21 + 22 + 23 + 24 + 25 + 26 + 27 + 28 + 29 + 30 + 31 + 32 + 33 + 34 + 35 + 36 + 37 + 38 + 39 + 40 + 41 + 42 + 43 + 44 + 45 + 46 + 47 + 48 + 49 + 50 + 51 + 52 + 53 + 54 + 55 + 56 + 57 + 58 + 59 + 60 + 61 + 62 + 63 + 64 + 65 + 66 + 67 + End SubModelPartNodes + Begin SubModelPartElements + 1 + 2 + 3 + 4 + 5 + 6 + 7 + 8 + 9 + 10 + 11 + 12 + 13 + 14 + 15 + 16 + 17 + 18 + 19 + 20 + 21 + 22 + 23 + 24 + 25 + 26 + 27 + 28 + 29 + 30 + 31 + 32 + 33 + 34 + 35 + 36 + 37 + 38 + 39 + 40 + 41 + 42 + 43 + 44 + 45 + 46 + 47 + 48 + 49 + 50 + 51 + 52 + 53 + 54 + 55 + 56 + 57 + 58 + 59 + 60 + 61 + 62 + 63 + 64 + 65 + 66 + 67 + 68 + 69 + 70 + 71 + 72 + 73 + 74 + 75 + 76 + 77 + 78 + 79 + 80 + 81 + 82 + 83 + 84 + 85 + 86 + 87 + 88 + 89 + 90 + 91 + 92 + 93 + 94 + 95 + 96 + 97 + 98 + 99 + 100 + 101 + 102 + 103 + 104 + 105 + 106 + 107 + 108 + 109 + 110 + 111 + 112 + 113 + 114 + 115 + 116 + 117 + 118 + 119 + 120 + 121 + 122 + 123 + 124 + 125 + 126 + 127 + 128 + 129 + 130 + 131 + 132 + 133 + 134 + 135 + 136 + 137 + 138 + 139 + 140 + 141 + 142 + 143 + 144 + 145 + 146 + 147 + 148 + 149 + 150 + 151 + End SubModelPartElements +End SubModelPart +Begin SubModelPart AutomaticInlet3D_Inlet // Group Inlet // Subtree AutomaticInlet3D + Begin SubModelPartNodes + 59 + 60 + 61 + 62 + 63 + 64 + 65 + 66 + 67 + End SubModelPartNodes + Begin SubModelPartElements + End SubModelPartElements + Begin SubModelPartConditions + 75 + 76 + 77 + 78 + 79 + 80 + 81 + 82 + End SubModelPartConditions +End SubModelPart +Begin SubModelPart Outlet3D_Outlet // Group Outlet // Subtree Outlet3D + Begin SubModelPartNodes + 1 + 3 + 4 + 5 + 10 + 11 + 13 + 14 + 22 + End SubModelPartNodes + Begin SubModelPartElements + End SubModelPartElements + Begin SubModelPartConditions + 83 + 84 + 85 + 86 + 87 + 88 + 89 + 90 + End SubModelPartConditions +End SubModelPart +Begin SubModelPart Slip3D // Condition Slip3D + Begin SubModelPartNodes + End SubModelPartNodes + Begin SubModelPartElements + End SubModelPartElements + Begin SubModelPartConditions + End SubModelPartConditions + Begin SubModelPart Slip3D_Wall1 // Group Wall1 // Subtree Slip3D + Begin SubModelPartNodes + 1 + 2 + 3 + 8 + 9 + 10 + 12 + 19 + 21 + 25 + 29 + 33 + 34 + 38 + 39 + 43 + 46 + 49 + 52 + 54 + 57 + 58 + 65 + 66 + 67 + End SubModelPartNodes + Begin SubModelPartElements + End SubModelPartElements + Begin SubModelPartConditions + 27 + 28 + 29 + 30 + 31 + 32 + 33 + 34 + 35 + 36 + 37 + 38 + 39 + 40 + 41 + 42 + 43 + 44 + 45 + 46 + 47 + 48 + 49 + 50 + 51 + 52 + 53 + 54 + End SubModelPartConditions + End SubModelPart + Begin SubModelPart Slip3D_Wall2 // Group Wall2 // Subtree Slip3D + Begin SubModelPartNodes + 11 + 14 + 15 + 20 + 22 + 23 + 24 + 28 + 30 + 32 + 35 + 41 + 45 + 47 + 51 + 55 + 59 + 60 + 63 + End SubModelPartNodes + Begin SubModelPartElements + End SubModelPartElements + Begin SubModelPartConditions + 55 + 56 + 57 + 58 + 59 + 60 + 61 + 62 + 63 + 64 + 65 + 66 + 67 + 68 + 69 + 70 + 71 + 72 + 73 + 74 + End SubModelPartConditions + End SubModelPart + Begin SubModelPart Slip3D_Top // Group Top // Subtree Slip3D + Begin SubModelPartNodes + 1 + 2 + 4 + 6 + 9 + 11 + 15 + 18 + 23 + 25 + 30 + 31 + 33 + 39 + 40 + 41 + 49 + 50 + 51 + 57 + 59 + 61 + 65 + End SubModelPartNodes + Begin SubModelPartElements + End SubModelPartElements + Begin SubModelPartConditions + 91 + 92 + 93 + 94 + 95 + 96 + 97 + 98 + 99 + 100 + 101 + 102 + 103 + 104 + 105 + 106 + 107 + 108 + 109 + 110 + 111 + 112 + 113 + 114 + 115 + 116 + End SubModelPartConditions + End SubModelPart + Begin SubModelPart Slip3D_Floor // Group Floor // Subtree Slip3D + Begin SubModelPartNodes + 10 + 12 + 13 + 17 + 21 + 22 + 24 + 27 + 28 + 29 + 35 + 37 + 38 + 43 + 44 + 45 + 52 + 53 + 55 + 58 + 63 + 64 + 67 + End SubModelPartNodes + Begin SubModelPartElements + End SubModelPartElements + Begin SubModelPartConditions + 1 + 2 + 3 + 4 + 5 + 6 + 7 + 8 + 9 + 10 + 11 + 12 + 13 + 14 + 15 + 16 + 17 + 18 + 19 + 20 + 21 + 22 + 23 + 24 + 25 + 26 + End SubModelPartConditions + End SubModelPart +End SubModelPart diff --git a/applications/FluidDynamicsHydraulicsApplication/tests/test_hydraulic_fluid_auxiliary_utilities.py b/applications/FluidDynamicsHydraulicsApplication/tests/test_hydraulic_fluid_auxiliary_utilities.py index 8cabf41a71c9..fbc39ffeee62 100644 --- a/applications/FluidDynamicsHydraulicsApplication/tests/test_hydraulic_fluid_auxiliary_utilities.py +++ b/applications/FluidDynamicsHydraulicsApplication/tests/test_hydraulic_fluid_auxiliary_utilities.py @@ -1,56 +1,25 @@ import KratosMultiphysics as Kratos import KratosMultiphysics.FluidDynamicsHydraulicsApplication as KratosFluidHydraulics + import KratosMultiphysics.KratosUnittest as UnitTest import KratosMultiphysics.kratos_utilities as KratosUtils - +import KratosMultiphysics.FluidDynamicsApplication as KratosCFD class HydraulicFluidAuxiliaryUtilitiesTest(UnitTest.TestCase): def setUp(self): - """Generates the model part to work on.""" - - # Set up the domain model parts + # Read test Model Part self.model = Kratos.Model() - test_model_part = self.model.CreateModelPart("TestModelPart") - inlet_model_part = test_model_part.CreateSubModelPart("Inlet") - test_model_part.ProcessInfo.SetValue(Kratos.DOMAIN_SIZE, 3) - - # Add VELOCITY variable to be used as DOF and DISTANCE - test_model_part.AddNodalSolutionStepVariable(Kratos.VELOCITY) - test_model_part.AddNodalSolutionStepVariable(Kratos.DISTANCE) - - # Modelpart: The inlet condition - # 7 - 8 - 9 - # | \ | \ | - # 4 - 5 - 6 - # | \ | \ | - # 1 - 2 - 3 - - # Populate model part nodes - test_model_part.CreateNewNode(1, 0.0, 0.0, 0.0) - test_model_part.CreateNewNode(2, 0.0, 1.0, 0.0) - test_model_part.CreateNewNode(3, 0.0, 2.0, 0.0) - test_model_part.CreateNewNode(4, 0.0, 0.0, 1.0) - test_model_part.CreateNewNode(5, 0.0, 1.0, 1.0) - test_model_part.CreateNewNode(6, 0.0, 2.0, 1.0) - test_model_part.CreateNewNode(7, 0.0, 0.0, 2.0) - test_model_part.CreateNewNode(8, 0.0, 1.0, 2.0) - test_model_part.CreateNewNode(9, 0.0, 2.0, 2.0) - inlet_model_part.AddNodes([1, 2, 3, 4, 5, 6, 7, 8, 9]) - - # Create inlet conditions - prop_0 = test_model_part.Properties[0] - inlet_model_part.CreateNewCondition("SurfaceCondition3D3N", 1, [1, 2, 4], prop_0) - inlet_model_part.CreateNewCondition("SurfaceCondition3D3N", 2, [2, 5, 4], prop_0) - inlet_model_part.CreateNewCondition("SurfaceCondition3D3N", 3, [2, 3, 5], prop_0) - inlet_model_part.CreateNewCondition("SurfaceCondition3D3N", 4, [3, 6, 5], prop_0) - inlet_model_part.CreateNewCondition("SurfaceCondition3D3N", 5, [4, 5, 7], prop_0) - inlet_model_part.CreateNewCondition("SurfaceCondition3D3N", 6, [5, 8, 7], prop_0) - inlet_model_part.CreateNewCondition("SurfaceCondition3D3N", 7, [5, 6, 8], prop_0) - inlet_model_part.CreateNewCondition("SurfaceCondition3D3N", 8, [6, 9, 8], prop_0) + self.model_part = self.model.CreateModelPart("MainModelPart") + self.model_part.ProcessInfo[Kratos.DOMAIN_SIZE] = 3 + self.model_part.AddNodalSolutionStepVariable(Kratos.DISTANCE) + self.model_part.AddNodalSolutionStepVariable(Kratos.BODY_FORCE) + self.model_part.AddNodalSolutionStepVariable(Kratos.VELOCITY) + self.model_part.AddNodalSolutionStepVariable(Kratos.NORMAL) + input_filename= "channel_test/channel3D" + Kratos.ModelPartIO(input_filename).ReadModelPart(self.model_part) def testCalculateWettedPerimeter(self): - - inlet_skin_model_part = self.model.GetModelPart("TestModelPart.Inlet") + inlet_skin_model_part = self.model.GetModelPart("MainModelPart.AutomaticInlet3D_Inlet") level_set_z= 0.3 for node in inlet_skin_model_part.Nodes: node.Set(Kratos.INLET, True) @@ -59,12 +28,11 @@ def testCalculateWettedPerimeter(self): for node in inlet_skin_model_part.Nodes: node.SetValue(Kratos.AUX_DISTANCE, node.Z - level_set_z) WettedPerimeter= KratosFluidHydraulics.HydraulicFluidAuxiliaryUtilities.CalculateWettedPetimeter(inlet_skin_model_part,Kratos.INLET, Kratos.AUX_DISTANCE, False) - theoretical_wetted_perimeter= 2.0 +2*(level_set_z) + theoretical_wetted_perimeter= 0.5 +2*(level_set_z) self.assertAlmostEqual(WettedPerimeter, theoretical_wetted_perimeter, 12) def testCalculateWettedArea(self): - - inlet_skin_model_part = self.model.GetModelPart("TestModelPart.Inlet") + inlet_skin_model_part = self.model.GetModelPart("MainModelPart.AutomaticInlet3D_Inlet") level_set_z= 0.3 for node in inlet_skin_model_part.Nodes: node.Set(Kratos.INLET, True) @@ -73,8 +41,38 @@ def testCalculateWettedArea(self): for node in inlet_skin_model_part.Nodes: node.SetValue(Kratos.AUX_DISTANCE, node.Z - level_set_z) WettedArea= KratosFluidHydraulics.HydraulicFluidAuxiliaryUtilities.CalculateWettedArea(inlet_skin_model_part,Kratos.INLET,Kratos.AUX_DISTANCE, False) - theoretical_wetted_area= 2.0*(level_set_z) + theoretical_wetted_area= 0.5*(level_set_z) self.assertAlmostEqual(WettedArea, theoretical_wetted_area, 12) + + def testFixCornerVelocity(self): + corner_angle=10 + Kratos.NormalCalculationUtils().CalculateOnSimplex(self.model_part,3) + Kratos.FindConditionsNeighboursProcess(self.model_part,3,3).Execute() + Kratos.VariableUtils().AddDof(Kratos.VELOCITY_X,self.model_part) + Kratos.VariableUtils().AddDof(Kratos.VELOCITY_Y,self.model_part) + Kratos.VariableUtils().AddDof(Kratos.VELOCITY_Z,self.model_part) + KratosFluidHydraulics.HydraulicFluidAuxiliaryUtilities.FixCornerNodeVelocity(self.model_part , corner_angle) + tested_node = self.model_part.GetNode(33) + self.assertAlmostEqual(tested_node.GetSolutionStepValue(Kratos.BODY_FORCE_Z),0.0,12) + + def testNonAirGravity(self): + for node in self.model_part.Nodes: + distance= node.Z-0.2 + node.SetSolutionStepValue(Kratos.DISTANCE,distance) + node.SetSolutionStepValue(Kratos.BODY_FORCE, [0.0, 0.0, -9.81]) + KratosFluidHydraulics.HydraulicFluidAuxiliaryUtilities.TurnOffGravityOnAirElements(self.model_part) + tested_node = self.model_part.GetNode(33) + self.assertAlmostEqual(tested_node.GetSolutionStepValue(Kratos.BODY_FORCE_Z),0.0,12) + + + def testMaximumWaterDepthChange(self): + inlet_skin_model_part = self.model.GetModelPart("MainModelPart.AutomaticInlet3D_Inlet") + for node in inlet_skin_model_part.Nodes: + node.SetSolutionStepValue(Kratos.DISTANCE, node.Z-0.3) + node.SetValue(Kratos.AUX_DISTANCE, node.Z-0.1) + max_depth_change= KratosFluidHydraulics.HydraulicFluidAuxiliaryUtilities.MaximumWaterDepthChange(self.model_part) + self.assertTrue(max_depth_change) + if __name__ == '__main__': UnitTest.main() \ No newline at end of file