Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile1D.cs =================================================================== diff -u -r5267 -r5290 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile1D.cs (.../SoilProfile1D.cs) (revision 5267) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile1D.cs (.../SoilProfile1D.cs) (revision 5290) @@ -159,23 +159,11 @@ get { IList sortedLayers = Layers.OrderBy(l => l.TopLevel).ToList(); - SoilLayer1D aquiferLayer = null; - int aquiferIndex = -1; + SoilLayer1D aquiferLayer = DeepestAquiferLayer; + int aquiferIndex = sortedLayers.IndexOf(DeepestAquiferLayer); - // Search deepest aquifer layer - for (var layerIndex = 0; layerIndex < sortedLayers.Count; layerIndex++) - { - SoilLayer1D layer = sortedLayers[layerIndex]; - if (IsAquiferLayer(layer)) - { - aquiferIndex = layerIndex; - aquiferLayer = layer; - break; - } - } - // aquifer may consist of more than 1 connected (aquifer) layers - // Search all layers above the first found aquifer to find top aquifer layer + // Search all layers above the deepest aquifer to find top aquifer layer if (aquiferIndex >= 0) { for (int layerIndex = aquiferIndex + 1; layerIndex < sortedLayers.Count; layerIndex++) @@ -195,7 +183,34 @@ return aquiferLayer; } } + + /// + /// Gets the deepest aquifer layer + /// + /// + /// The deepest aquifer layer + /// + public SoilLayer1D DeepestAquiferLayer + { + get + { + IList sortedLayers = Layers.OrderBy(l => l.TopLevel).ToList(); + SoilLayer1D aquiferLayer = null; + + for (var layerIndex = 0; layerIndex < sortedLayers.Count; layerIndex++) + { + SoilLayer1D layer = sortedLayers[layerIndex]; + if (IsAquiferLayer(layer)) + { + aquiferLayer = layer; + break; + } + } + return aquiferLayer; + } + } + /// /// Gets the list of all clusters of in-between aquifers: /// Item1 is the top aquifer and Item2 the bottom aquifer of the cluster. @@ -485,19 +500,38 @@ /// public SoilLayer1D GetHighestAquifer() { - SoilLayer1D topLevel; + SoilLayer1D layer; if (InBetweenAquiferLayer != null) { - topLevel = InBetweenAquiferLayer; + layer = InBetweenAquiferLayer; } else { - topLevel = BottomAquiferLayer; + layer = BottomAquiferLayer; } - return topLevel; + return layer; } + + /// + /// Gets the lowest layer of the highest cluster of aquifers. + /// + /// + public SoilLayer1D GetLowestLayerOfHighestAquiferCluster() + { + SoilLayer1D layer; + if (InBetweenAquiferLayer != null) + { + layer = GetInBetweenAquiferClusters.Count == 0 ? null : GetInBetweenAquiferClusters.First().Item2; + } + else + { + layer = DeepestAquiferLayer; + } + return layer; + } + /// /// Returns a that represents this instance. /// Index: DamEngine/trunk/src/Deltares.DamEngine.TestHelpers/Factories/FactoryForSoilProfiles.cs =================================================================== diff -u -r5275 -r5290 --- DamEngine/trunk/src/Deltares.DamEngine.TestHelpers/Factories/FactoryForSoilProfiles.cs (.../FactoryForSoilProfiles.cs) (revision 5275) +++ DamEngine/trunk/src/Deltares.DamEngine.TestHelpers/Factories/FactoryForSoilProfiles.cs (.../FactoryForSoilProfiles.cs) (revision 5290) @@ -1264,7 +1264,7 @@ point3, point4 } - ], soil); + ], soil, isAquifer); } public static SoilLayer2D CreatePentagonSoilLayer2D(Point2D point1, Point2D point2, Point2D point3, Point2D point4, @@ -1300,7 +1300,7 @@ point4, point5 } - ], soil); + ], soil, isAquifer); } public static SoilLayer2D CreateHexagonSoilLayer2D(Point2D point1, Point2D point2, Point2D point3, Point2D point4, Point2D point5, Point2D point6, Soil soil = null, bool isAquifer = false) @@ -1315,7 +1315,7 @@ point5, point6 } - ], soil); + ], soil, isAquifer); } public static SoilLayer2D CreateHeptagonSoilLayer2D(Point2D point1, Point2D point2, Point2D point3, Point2D point4, Point2D point5, Point2D point6, Point2D point7, Soil soil) @@ -1334,7 +1334,7 @@ ], soil); } - public static SoilLayer2D CreatePolygoneSoilLayer2D(List points, Soil soil) + public static SoilLayer2D CreatePolygoneSoilLayer2D(List points, Soil soil, bool isAquifer = false) { var soilLayer2D = new SoilLayer2D { @@ -1343,7 +1343,7 @@ OuterLoop = new GeometryLoop() }, Soil = soil, - IsAquifer = false + IsAquifer = isAquifer }; for (var i = 0; i < points.Count - 1; i++) { Index: DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Geotechnics/SoilProfile1DTests.cs =================================================================== diff -u -r4934 -r5290 --- DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Geotechnics/SoilProfile1DTests.cs (.../SoilProfile1DTests.cs) (revision 4934) +++ DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Geotechnics/SoilProfile1DTests.cs (.../SoilProfile1DTests.cs) (revision 5290) @@ -48,6 +48,8 @@ Assert.That(soilProfile1D.GetHighestAquifer().BottomLevel, Is.EqualTo(topLevelLayer4)); Assert.That(soilProfile1D.BottomAquiferLayer.TopLevel, Is.EqualTo(topLevelLayer3)); Assert.That(soilProfile1D.BottomAquiferLayer.BottomLevel, Is.EqualTo(topLevelLayer4)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().TopLevel, Is.EqualTo(topLevelLayer3)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().BottomLevel, Is.EqualTo(topLevelLayer4)); }); // 3 clusters with aquifer, aquitard, aquifer => no in-between aquifer and bottom aquifer is the highest layer of cluster 5-6 @@ -60,6 +62,8 @@ Assert.That(soilProfile1D.GetHighestAquifer().BottomLevel, Is.EqualTo(topLevelLayer6)); Assert.That(soilProfile1D.BottomAquiferLayer.TopLevel, Is.EqualTo(topLevelLayer5)); Assert.That(soilProfile1D.BottomAquiferLayer.BottomLevel, Is.EqualTo(topLevelLayer6)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().TopLevel, Is.EqualTo(topLevelLayer6)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().BottomLevel, Is.EqualTo(bottomLevelLayer6)); }); // All layers are aquifers => no in-between aquifer and bottom aquifer is the highest aquifer (layer 1) @@ -72,6 +76,8 @@ Assert.That(soilProfile1D.GetHighestAquifer().BottomLevel, Is.EqualTo(topLevelLayer2)); Assert.That(soilProfile1D.BottomAquiferLayer.TopLevel, Is.EqualTo(topLevelLayer1)); Assert.That(soilProfile1D.BottomAquiferLayer.BottomLevel, Is.EqualTo(topLevelLayer2)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().TopLevel, Is.EqualTo(topLevelLayer6)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().BottomLevel, Is.EqualTo(bottomLevelLayer6)); }); // All layers are aquitard => no in-between aquifer and no bottom aquifer @@ -82,6 +88,7 @@ Assert.That(soilProfile1D.InBetweenAquiferLayer, Is.Null); Assert.That(soilProfile1D.GetHighestAquifer(), Is.Null); Assert.That(soilProfile1D.BottomAquiferLayer, Is.Null); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster(), Is.Null); }); // Aquifers are at top and bottom of the profile => no in-between aquifer and bottom aquifer is the highest layer of cluster 3-6 @@ -94,6 +101,8 @@ Assert.That(soilProfile1D.GetHighestAquifer().BottomLevel, Is.EqualTo(topLevelLayer4)); Assert.That(soilProfile1D.BottomAquiferLayer.TopLevel, Is.EqualTo(topLevelLayer3)); Assert.That(soilProfile1D.BottomAquiferLayer.BottomLevel, Is.EqualTo(topLevelLayer4)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().TopLevel, Is.EqualTo(topLevelLayer6)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().BottomLevel, Is.EqualTo(bottomLevelLayer6)); }); // Top aquitard cluster and bottom aquifer cluster => no in-between aquifer and bottom aquifer is the highest layer of cluster 4-6 @@ -106,6 +115,8 @@ Assert.That(soilProfile1D.GetHighestAquifer().BottomLevel, Is.EqualTo(topLevelLayer5)); Assert.That(soilProfile1D.BottomAquiferLayer.TopLevel, Is.EqualTo(topLevelLayer4)); Assert.That(soilProfile1D.BottomAquiferLayer.BottomLevel, Is.EqualTo(topLevelLayer5)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().TopLevel, Is.EqualTo(topLevelLayer6)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().BottomLevel, Is.EqualTo(bottomLevelLayer6)); }); // 2 in-between aquifers present (layers 2 and 4) and bottom aquifer is the lowest layer @@ -123,6 +134,8 @@ Assert.That(soilProfile1D.GetHighestAquifer().BottomLevel, Is.EqualTo(topLevelLayer3)); Assert.That(soilProfile1D.BottomAquiferLayer.TopLevel, Is.EqualTo(topLevelLayer6)); Assert.That(soilProfile1D.BottomAquiferLayer.BottomLevel, Is.EqualTo(bottomLevelLayer6)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().TopLevel, Is.EqualTo(topLevelLayer2)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().BottomLevel, Is.EqualTo(topLevelLayer3)); }); // 1 in-between aquifers cluster present (layers 2-4) and bottom aquifer is the lowest layer @@ -138,6 +151,8 @@ Assert.That(soilProfile1D.GetHighestAquifer().BottomLevel, Is.EqualTo(topLevelLayer3)); Assert.That(soilProfile1D.BottomAquiferLayer.TopLevel, Is.EqualTo(topLevelLayer6)); Assert.That(soilProfile1D.BottomAquiferLayer.BottomLevel, Is.EqualTo(bottomLevelLayer6)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().TopLevel, Is.EqualTo(topLevelLayer4)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().BottomLevel, Is.EqualTo(topLevelLayer5)); }); // 1 in-between aquifer present (layer 2) and bottom aquifer is the highest layer of cluster 4-6 @@ -153,6 +168,8 @@ Assert.That(soilProfile1D.GetHighestAquifer().BottomLevel, Is.EqualTo(topLevelLayer3)); Assert.That(soilProfile1D.BottomAquiferLayer.TopLevel, Is.EqualTo(topLevelLayer4)); Assert.That(soilProfile1D.BottomAquiferLayer.BottomLevel, Is.EqualTo(topLevelLayer5)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().TopLevel, Is.EqualTo(topLevelLayer2)); + Assert.That(soilProfile1D.GetLowestLayerOfHighestAquiferCluster().BottomLevel, Is.EqualTo(topLevelLayer3)); }); } Index: DamEngine/trunk/src/Deltares.DamEngine.Calculators/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverter.cs =================================================================== diff -u -r5284 -r5290 --- DamEngine/trunk/src/Deltares.DamEngine.Calculators/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverter.cs (.../PlLinesToWaternetConverter.cs) (revision 5284) +++ DamEngine/trunk/src/Deltares.DamEngine.Calculators/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverter.cs (.../PlLinesToWaternetConverter.cs) (revision 5290) @@ -20,6 +20,7 @@ // All rights reserved. using System; +using System.Collections; using System.Collections.Generic; using System.Data; using System.Linq; @@ -29,7 +30,6 @@ using Deltares.DamEngine.Data.Geometry; using Deltares.DamEngine.Data.Geotechnics; using Deltares.DamEngine.Data.Standard; -using Deltares.MacroStability.CSharpWrapper.Input; using CharacteristicPointType = Deltares.DamEngine.Data.Geotechnics.CharacteristicPointType; using SoilProfile = Deltares.DamEngine.Data.Geotechnics.SoilProfile; @@ -39,13 +39,18 @@ { public enum LayerType { - BottomAquifer, - TopLayerInBetweenAquiferCluster, - BottomLayerInBetweenAquiferCluster, - HighestAquifer, + BottomAquiferCluster, + InBetweenAquiferCluster, + HighestAquiferCluster, LowestLayer } + private enum BoundaryType + { + Top, + Bottom + } + private const double deviationX = 1e-06; private const double toleranceAlmostEqual = 1e-09; @@ -147,8 +152,8 @@ ThrowWhenPlLinesIsNull(plLines); ThrowWhenSoilProfileIsNull(soilProfile); - Point2D[] bottomAquiferCoordinates = DetermineLayerBoundaryPoints(LayerType.BottomAquifer, soilProfile).ToArray(); - if (!bottomAquiferCoordinates.Any()) + DetermineLayerBoundaryPoints(LayerType.BottomAquiferCluster, soilProfile, out Point2D[] bottomAquiferCoordinates, out _); + if (bottomAquiferCoordinates == null) { throw new NoNullAllowedException(string.Format(Resources.NoBottomAquiferLayer, soilProfile.Name)); } @@ -270,9 +275,8 @@ waternet.HeadLineList.Add(headLine); for (var i = 0; i < inBetweenAquiferCount; i++) { - Point2D[] inBetweenAquiferUpperCoordinates = DetermineLayerBoundaryPoints(LayerType.TopLayerInBetweenAquiferCluster, soilProfile, i).ToArray(); - Point2D[] inBetweenAquiferLowerCoordinates = DetermineLayerBoundaryPoints(LayerType.BottomLayerInBetweenAquiferCluster, soilProfile, i).ToArray(); - + DetermineLayerBoundaryPoints(LayerType.InBetweenAquiferCluster, soilProfile, out Point2D[] inBetweenAquiferUpperCoordinates, out Point2D[] inBetweenAquiferLowerCoordinates, i); + if (inBetweenAquiferUpperCoordinates.Any() && inBetweenAquiferLowerCoordinates.Any()) { WaternetLine waternetLineTop = CreateWaternetLine(inBetweenAquiferUpperCoordinates); @@ -348,20 +352,21 @@ /// /// Determines all the points of the layer boundary for the specified layer type: - /// - For LayerType = BottomAquifer, it will return the top boundary of the bottom aquifer (= waternet line for PL3) - /// - For TopLayerInBetweenAquiferCluster, it will return the top boundary of the in-between aquifer (= waternet line for PL4) + /// - For LayerType = BottomAquifer, layerBoundaryTop is the top boundary of the bottom aquifer (= waternet line for PL3) + /// - For LayerType = InBetweenAquiferCluster, it will return the top boundary of the in-between aquifer (= waternet line for PL4) /// - For BottomLayerInBetweenAquiferCluster, it will return the bottom boundary of the in-between aquifer (= waternet line for PL4) - /// - For HighestAquifer, it will return the top boundary of the highest aquifer (= waternet line for PL1 for "Hydrostatic" waternet) - /// - For LowestLayer, it will return the bottom boundary of the lowest layer (= waternet line for PL1 for "Full hydrostatic" waternet) + /// - For LayerType = HighestAquifer, it will return the top boundary of the highest aquifer (= waternet line for PL1 for "Hydrostatic" waternet) + /// - For LayerType = LowestLayer, it will return the bottom boundary of the lowest layer (= waternet line for PL1 for "Full hydrostatic" waternet) /// /// The type of layer boundary. /// The soil profile 2D. + /// All the points of the layer top boundary. + /// All the points of the layer bottom boundary. /// In case of several in-between aquifers, the index of the in-between aquifer must be specified. /// All the points of the layer boundary. - internal static IEnumerable DetermineLayerBoundaryPoints(LayerType layerType, SoilProfile2D soilProfile, int indexInBetweenAquifer = -1) + internal static void DetermineLayerBoundaryPoints(LayerType layerType, SoilProfile2D soilProfile, out Point2D[] layerBoundaryTop, out Point2D[] layerBoundaryBottom, int indexInBetweenAquifer = -1) { double[] xCoordinates = DetermineAllXCoordinatesOfSoilProfile(soilProfile); - SoilLayer1D previousAquiferLayer = null; var xCoordinatesAll = new List(); foreach (double xCoordinate in xCoordinates) { @@ -376,55 +381,79 @@ } } - var layerBoundaryPoints = new List(); + var layerBoundaryTopPoints = new List(); + var layerBoundaryBottomPoints = new List(); + double previousAquiferTop = double.NaN; + double previousAquiferBottom = double.NaN; foreach (double xCoordinate in xCoordinatesAll) { SoilProfile1D crossSection = soilProfile.GetSoilProfile1D(xCoordinate); - // Determine if the layer is in range of the previous layer - // If not, return empty coordinates, because the layer is interrupted - SoilLayer1D currentAquifer = GetSoilLayer1D(layerType, crossSection, indexInBetweenAquifer); - if (previousAquiferLayer != null && currentAquifer != null - && !AreHorizontallyConnected(previousAquiferLayer, currentAquifer)) + // Determine if the cluster of layers is in range of the previous cluster of layers. + // If not, return empty coordinates, because the cluster of layer is interrupted. + double currentAquiferTop = GetLevel(layerType, BoundaryType.Top, crossSection, indexInBetweenAquifer); + double currentAquiferBottom = GetLevel(layerType, BoundaryType.Bottom, crossSection, indexInBetweenAquifer); + + if (!double.IsNaN(previousAquiferTop) && !double.IsNaN(previousAquiferBottom) + && !double.IsNaN(currentAquiferTop) + && !double.IsNaN(currentAquiferBottom) + && !AreHorizontallyConnected(previousAquiferTop, previousAquiferBottom, currentAquiferTop, currentAquiferBottom)) { - return Enumerable.Empty(); + layerBoundaryTop = null; + layerBoundaryBottom = null; + return; } - - if (currentAquifer != null) + + if (!double.IsNaN(currentAquiferTop) && !double.IsNaN(currentAquiferBottom)) { - previousAquiferLayer = currentAquifer; - layerBoundaryPoints.Add(layerType is LayerType.BottomLayerInBetweenAquiferCluster or LayerType.LowestLayer ? new Point2D(xCoordinate, currentAquifer.BottomLevel) : new Point2D(xCoordinate, currentAquifer.TopLevel)); + previousAquiferTop = currentAquiferTop; + previousAquiferBottom = currentAquiferBottom; + layerBoundaryTopPoints.Add(new Point2D(xCoordinate, currentAquiferTop)); + layerBoundaryBottomPoints.Add(new Point2D(xCoordinate, currentAquiferBottom)); } } // Perform a short validation that the coordinates are fully defined from the beginning to the end // of the profile. If not, the layer is not continuous. - if (!IsLayerBoundaryContinuous(layerBoundaryPoints, xCoordinates.First(), xCoordinates.Last())) + if (!IsLayerBoundaryContinuous(layerBoundaryTopPoints, xCoordinates.First(), xCoordinates.Last()) + || !IsLayerBoundaryContinuous(layerBoundaryBottomPoints, xCoordinates.First(), xCoordinates.Last())) { - return Enumerable.Empty(); + layerBoundaryTop = null; + layerBoundaryBottom = null; + return; } - - return layerBoundaryPoints; + layerBoundaryTop = layerBoundaryTopPoints.ToArray(); + layerBoundaryBottom = layerBoundaryBottomPoints.ToArray(); } - private static bool IsLayerBoundaryContinuous(List boundaryPoints, double leftGeometryBoundary, double rightGeometryBoundary) + private static bool IsLayerBoundaryContinuous(IEnumerable boundaryPoints, double leftGeometryBoundary, double rightGeometryBoundary) { - return boundaryPoints.Any() && (Math.Abs(boundaryPoints.First().X - leftGeometryBoundary) < toleranceAlmostEqual) - && Math.Abs(boundaryPoints.Last().X - rightGeometryBoundary) < toleranceAlmostEqual; + IEnumerable point2Ds = boundaryPoints.ToList(); + return point2Ds.Any() && (Math.Abs(point2Ds.First().X - leftGeometryBoundary) < toleranceAlmostEqual) + && Math.Abs(point2Ds.Last().X - rightGeometryBoundary) < toleranceAlmostEqual; } - private static SoilLayer1D GetSoilLayer1D(LayerType layerType, SoilProfile1D soilProfile1D, int indexInBetweenAquifer) + private static double GetLevel(LayerType layerType, BoundaryType boundaryType, SoilProfile1D soilProfile1D, int indexInBetweenAquifer) { + SoilLayer1D layer = GetLayer(layerType, boundaryType, soilProfile1D, indexInBetweenAquifer); + if (layer != null) + { + return boundaryType == BoundaryType.Top ? layer.TopLevel : layer.BottomLevel; + } + + return double.NaN; + } + + private static SoilLayer1D GetLayer(LayerType layerType, BoundaryType boundaryType, SoilProfile1D soilProfile1D, int indexInBetweenAquifer) + { switch (layerType) { - case LayerType.BottomAquifer: - return soilProfile1D.BottomAquiferLayer; - case LayerType.TopLayerInBetweenAquiferCluster: - return soilProfile1D.GetInBetweenAquiferClusters[indexInBetweenAquifer].Item1; - case LayerType.BottomLayerInBetweenAquiferCluster: - return soilProfile1D.GetInBetweenAquiferClusters[indexInBetweenAquifer].Item2; - case LayerType.HighestAquifer: - return soilProfile1D.GetHighestAquifer(); + case LayerType.BottomAquiferCluster: + return boundaryType == BoundaryType.Top ? soilProfile1D.BottomAquiferLayer : soilProfile1D.DeepestAquiferLayer; + case LayerType.InBetweenAquiferCluster: + return boundaryType == BoundaryType.Top ? soilProfile1D.GetInBetweenAquiferClusters[indexInBetweenAquifer].Item1 : soilProfile1D.GetInBetweenAquiferClusters[indexInBetweenAquifer].Item2; + case LayerType.HighestAquiferCluster: + return boundaryType == BoundaryType.Top ? soilProfile1D.GetHighestAquifer() : soilProfile1D.GetLowestLayerOfHighestAquiferCluster(); case LayerType.LowestLayer: return soilProfile1D.Layers.Last(); } @@ -511,7 +540,7 @@ } case IntrusionVerticalWaterPressureType.HydroStatic: { - Point2D[] highestAquiferCoordinates = DetermineLayerBoundaryPoints(LayerType.HighestAquifer, soilProfile2D).ToArray(); + DetermineLayerBoundaryPoints(LayerType.HighestAquiferCluster, soilProfile2D, out Point2D[] highestAquiferCoordinates, out _); return CreateWaternetLine(highestAquiferCoordinates); } case IntrusionVerticalWaterPressureType.Linear: @@ -521,7 +550,7 @@ } case IntrusionVerticalWaterPressureType.FullHydroStatic: { - Point2D[] lowestBoundaryCoordinates = DetermineLayerBoundaryPoints(LayerType.LowestLayer, soilProfile2D).ToArray(); + DetermineLayerBoundaryPoints(LayerType.LowestLayer, soilProfile2D, out _, out Point2D[] lowestBoundaryCoordinates); return CreateWaternetLine(lowestBoundaryCoordinates); } case null: @@ -613,48 +642,39 @@ { SoilProfile1D crossSectionLeft = soilProfile.GetSoilProfile1D(xCoordinate - deviationX); SoilProfile1D crossSectionRight = soilProfile.GetSoilProfile1D(xCoordinate + deviationX); - if (crossSectionLeft != null && crossSectionRight != null) + if (crossSectionLeft == null || crossSectionRight == null) { - SoilLayer1D currentAquiferLeft = GetSoilLayer1D(layerType, crossSectionLeft, indexInBetweenAquifer); - SoilLayer1D currentAquiferRight = GetSoilLayer1D(layerType, crossSectionRight, indexInBetweenAquifer); - if (currentAquiferLeft != null && currentAquiferRight != null) - { - if (layerType is LayerType.BottomAquifer or LayerType.TopLayerInBetweenAquiferCluster) - { - return Math.Abs(currentAquiferLeft.TopLevel - currentAquiferRight.TopLevel) > GeometryConstants.Accuracy; - } + return false; + } - return Math.Abs(currentAquiferLeft.BottomLevel - currentAquiferRight.BottomLevel) > GeometryConstants.Accuracy; - } + double aquiferLeftTop = GetLevel(layerType, BoundaryType.Top, crossSectionLeft, indexInBetweenAquifer); + double aquiferLefBottom = GetLevel(layerType, BoundaryType.Bottom, crossSectionLeft, indexInBetweenAquifer); + double aquiferRightTop = GetLevel(layerType, BoundaryType.Top, crossSectionRight, indexInBetweenAquifer); + double aquiferRightBottom = GetLevel(layerType, BoundaryType.Bottom, crossSectionRight, indexInBetweenAquifer); + if (!double.IsNaN(aquiferLeftTop) && !double.IsNaN(aquiferLefBottom) && !double.IsNaN(aquiferRightTop) && !double.IsNaN(aquiferRightBottom)) + { + return Math.Abs(aquiferLeftTop - aquiferRightTop) > GeometryConstants.Accuracy || Math.Abs(aquiferLefBottom - aquiferRightBottom) > GeometryConstants.Accuracy; } return false; } - private static bool AreHorizontallyConnected(SoilLayer1D leftSoilLayer, SoilLayer1D rightSoilLayer) + private static bool AreHorizontallyConnected(double leftLayerTop, double leftLayerBottom, double rightLayerTop, double rightLayerBottom) { - // The layers are identical - if (ReferenceEquals(leftSoilLayer, rightSoilLayer)) - { - return true; - } - // Left soil layer envelopes whole right soil layer - if (leftSoilLayer.BottomLevel <= rightSoilLayer.BottomLevel - && leftSoilLayer.TopLevel >= rightSoilLayer.TopLevel) + if (leftLayerBottom <= rightLayerBottom && leftLayerTop >= rightLayerTop) { return true; } // Right soil layer envelopes whole left soil layer - if (rightSoilLayer.BottomLevel <= leftSoilLayer.BottomLevel - && rightSoilLayer.TopLevel >= leftSoilLayer.TopLevel) + if (rightLayerBottom <= leftLayerBottom && rightLayerTop >= leftLayerTop) { return true; } - return (rightSoilLayer.TopLevel <= leftSoilLayer.TopLevel && rightSoilLayer.TopLevel >= leftSoilLayer.BottomLevel) // Top level lies inbetween the left soil layer - || (rightSoilLayer.BottomLevel >= leftSoilLayer.BottomLevel && rightSoilLayer.BottomLevel <= leftSoilLayer.TopLevel); // Bottom level lies inbetween the left soil layer + return (rightLayerTop <= leftLayerTop && rightLayerTop >= leftLayerBottom) // Top level lies in-between the left soil layer + || (rightLayerBottom >= leftLayerBottom && rightLayerBottom <= leftLayerTop); // Bottom level lies in-between the left soil layer } private static bool IsBelowSoilProfile(SoilProfile1D soilProfile, PlLine line) Index: DamEngine/trunk/src/Deltares.DamEngine.Calculators.Tests/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverterTests.cs =================================================================== diff -u -r5281 -r5290 --- DamEngine/trunk/src/Deltares.DamEngine.Calculators.Tests/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverterTests.cs (.../PlLinesToWaternetConverterTests.cs) (revision 5281) +++ DamEngine/trunk/src/Deltares.DamEngine.Calculators.Tests/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverterTests.cs (.../PlLinesToWaternetConverterTests.cs) (revision 5290) @@ -40,7 +40,7 @@ public class PlLinesToWaternetConverterTests { private const string mapTestFiles = @"KernelWrappers\MacroStabilityCommon\TestFiles\"; - + private const double precision3Decimals = 0.00051; private const double precision5Decimals = 0.0000051; private const double penetrationLength = 0.5; @@ -56,16 +56,16 @@ // Setup var line = new Line(new Point2D(0, 0), new Point2D(90, 9)); var xCoordinatesToBeAdded = new List() - + { 30, 60 }; List splitLines = []; - + // Call PlLinesToWaternetConverter.SplitLineAtXCoordinate(xCoordinatesToBeAdded, line, splitLines); - + // Assert Assert.That(splitLines, Has.Count.EqualTo(3)); Assert.Multiple(() => @@ -86,7 +86,7 @@ Assert.That(splitLines[2].EndPoint.Z, Is.EqualTo(9)); }); } - + [Test] [TestCase(0)] [TestCase(0.01)] @@ -141,7 +141,7 @@ { const double pL1Level = 0; PlLines plLines = CreateAllPlLines(pL1Level); - + Assert.Multiple(() => { Assert.That(soilProfile1D.GetInBetweenAquiferClusters, Has.Count.EqualTo(2).Within(precision3Decimals)); @@ -279,7 +279,7 @@ }); } } - + [TestCase(IntrusionVerticalWaterPressureType.Standard, -6.110)] [TestCase(IntrusionVerticalWaterPressureType.Linear, 1.212)] [TestCase(IntrusionVerticalWaterPressureType.HydroStatic, -2.110)] @@ -376,9 +376,9 @@ // check that an exception is raised when no aquifer present soilProfile1D = CreateSoilProfile1DWithoutAquifersForTest(); - Assert.That(() => PlLinesToWaternetConverter.CreateWaternetBasedOnPlLines(plLines, soilProfile1D, surfaceLine, penetrationLength, IntrusionVerticalWaterPressureType.Standard), + Assert.That(() => PlLinesToWaternetConverter.CreateWaternetBasedOnPlLines(plLines, soilProfile1D, surfaceLine, penetrationLength, IntrusionVerticalWaterPressureType.Standard), Throws.InstanceOf().With.Message.EqualTo("The deepest layer is not an aquifer in soilprofile 'DefaultNameSoilProfile1D'.")); - + } public class Given2DSoilProfileWithTwoClustersOfInBetweenAquifersAndVerticalLayerSeparations @@ -610,15 +610,105 @@ SoilProfile2D soilProfile, IEnumerable expectedBottomAquiferCoordinates) { // Call - IEnumerable bottomAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.BottomAquifer, soilProfile); - IEnumerable highestAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.HighestAquifer, soilProfile); + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.BottomAquiferCluster, soilProfile, out Point2D[] bottomAquifer, out _); + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.HighestAquiferCluster, soilProfile, out Point2D[] highestAquifer, out _); // Assert AssertGeometry(bottomAquifer, expectedBottomAquiferCoordinates); AssertGeometry(highestAquifer, expectedBottomAquiferCoordinates); } /// + /// |-------------------------------------------------------------| Level 0 m + /// | Aquitard 1 | + /// |-------------------------------------------------------------| Level -10 m + /// | Aquitard 2 B /------------------------\ C | Level -12 m + /// | / Aquifer 1 \ | + /// | / point X \ | + /// |---------------/--------------*---------------\--------------| Level -20 m + /// | A Aquifer 2 D | + /// |------------------------------*------------------------------| Level -30 m + /// point Y + /// + [Test] + [TestCase(-20,-30)] + [TestCase(-20.1,-30)] + [TestCase(-19.9,-30)] + [TestCase(-20,-30.1)] + [TestCase(-20,-29.9)] + public void Given2DSoilProfileWithNotContinuousAquiferConnectedToContinuousAquifer_WhenDeterminingLayerBoundaryPointsOfAquifer_ThenExpectedCoordinatesReturned(double zPointX, double zPointY) + { + var pointA = new Point2D(leftCoordinate + (rightCoordinate - leftCoordinate) / 3, -20); + var pointB = new Point2D(pointA.X + 0.01, -12); + var pointC = new Point2D(leftCoordinate + 2 * (rightCoordinate - leftCoordinate) / 3, -12); + var pointD = new Point2D(pointC.X + 0.01, -20); + var pointX = new Point2D(middleXCoordinate, zPointX); + var pointY = new Point2D(middleXCoordinate, zPointY); + + SoilLayer2D soilUpperLayer = CreateRectangularSoilLayer2D(0, -10, leftCoordinate, rightCoordinate, false); + SoilLayer2D soilLayerAquitard2 = FactoryForSoilProfiles.CreatePolygoneSoilLayer2D([ + ..new[] + { + new Point2D(leftCoordinate, -10), + new Point2D(rightCoordinate, -10), + new Point2D(rightCoordinate, -20), + pointD, + pointC, + pointB, + pointA, + new Point2D(leftCoordinate, -20) + } + ], null); + SoilLayer2D soilLayerAquifer1 = FactoryForSoilProfiles.CreatePentagonSoilLayer2D(pointA, pointB, pointC, pointD, pointX, null, null, true); + SoilLayer2D soilLayerAquifer2 = FactoryForSoilProfiles.CreatePolygoneSoilLayer2D([ + ..new[] + { + new Point2D(leftCoordinate, -20), + pointA, + pointX, + pointD, + new Point2D(rightCoordinate, -20), + new Point2D(rightCoordinate, -30), + pointY, + new Point2D(leftCoordinate, -30) + } + ], null, true); + + var soilProfile = new SoilProfile2D + { + Geometry = new GeometryData + { + Left = leftCoordinate, + Right = rightCoordinate, + Bottom = -30 + } + }; + soilProfile.Surfaces.Add(soilUpperLayer); + soilProfile.Surfaces.Add(soilLayerAquitard2); + soilProfile.Surfaces.Add(soilLayerAquifer1); + soilProfile.Surfaces.Add(soilLayerAquifer2); + soilProfile.Geometry.Surfaces.Add(new GeometrySurface()); + + // Call + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.BottomAquiferCluster, soilProfile, out Point2D[] bottomAquifer, out _); + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.HighestAquiferCluster, soilProfile, out Point2D[] highestAquifer, out _); + + // Assert + var expectedBottomAquiferCoordinates = new[] + { + new GeometryPoint(leftCoordinate, -20), + new GeometryPoint(pointA.X, pointA.Z), + new GeometryPoint(pointB.X, pointB.Z), + new GeometryPoint(middleXCoordinate, -12), + new GeometryPoint(pointC.X, pointC.Z), + new GeometryPoint(pointD.X, pointD.Z), + new GeometryPoint(rightCoordinate, -20) + }; + AssertGeometry(bottomAquifer, expectedBottomAquiferCoordinates); + AssertGeometry(highestAquifer, expectedBottomAquiferCoordinates); + } + + /// /// --------------------------------------------------------------- Level 0 m /// top layer /// --------------------------------------------------------------- Level -10 m @@ -657,9 +747,8 @@ // Call bool isInBetweenLayerPresent = PlLinesToWaternetConverter.AreInBetweenAquiferClustersPresent(soilProfile, out int count); - IEnumerable topLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.TopLayerInBetweenAquiferCluster, soilProfile, 0); - IEnumerable bottomLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.BottomLayerInBetweenAquiferCluster, soilProfile, 0); - + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.InBetweenAquiferCluster, soilProfile, out Point2D[] topLevelInBetweenAquifer, out Point2D[] bottomLevelInBetweenAquifer, 0); + // Assert Assert.Multiple(() => { @@ -719,8 +808,7 @@ // Call bool isInBetweenLayerPresent = PlLinesToWaternetConverter.AreInBetweenAquiferClustersPresent(soilProfile, out int count); - IEnumerable topLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.TopLayerInBetweenAquiferCluster, soilProfile, 0); - IEnumerable bottomLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.BottomLayerInBetweenAquiferCluster, soilProfile, 0); + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.InBetweenAquiferCluster, soilProfile, out Point2D[] topLevelInBetweenAquifer, out Point2D[] bottomLevelInBetweenAquifer, 0); // Assert Assert.Multiple(() => @@ -812,8 +900,7 @@ // Call bool isInBetweenLayerPresent = PlLinesToWaternetConverter.AreInBetweenAquiferClustersPresent(soilProfile, out int count); - IEnumerable topLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.TopLayerInBetweenAquiferCluster, soilProfile, 0); - IEnumerable bottomLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.BottomLayerInBetweenAquiferCluster, soilProfile, 0); + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.InBetweenAquiferCluster, soilProfile, out Point2D[] topLevelInBetweenAquifer, out Point2D[] bottomLevelInBetweenAquifer, 0); // Assert Assert.Multiple(() => @@ -905,8 +992,7 @@ // Call bool isInBetweenLayerPresent = PlLinesToWaternetConverter.AreInBetweenAquiferClustersPresent(soilProfile, out int count); - IEnumerable topLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.TopLayerInBetweenAquiferCluster, soilProfile, 0); - IEnumerable bottomLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.BottomLayerInBetweenAquiferCluster, soilProfile, 0); + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.InBetweenAquiferCluster, soilProfile, out Point2D[] topLevelInBetweenAquifer, out Point2D[] bottomLevelInBetweenAquifer, 0); // Assert Assert.Multiple(() => @@ -925,6 +1011,7 @@ [ new(leftCoordinate, -20), new(middleXCoordinate, -20), + new(middleXCoordinate, -20), new(rightCoordinate, -20) ]; AssertGeometry(topLevelInBetweenAquifer, expectedTop); @@ -1023,8 +1110,7 @@ // Call bool isInBetweenLayerPresent = PlLinesToWaternetConverter.AreInBetweenAquiferClustersPresent(soilProfile, out int count); - IEnumerable topLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.TopLayerInBetweenAquiferCluster, soilProfile, 0); - IEnumerable bottomLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.BottomLayerInBetweenAquiferCluster, soilProfile, 0); + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.InBetweenAquiferCluster, soilProfile, out Point2D[] topLevelInBetweenAquifer, out Point2D[] bottomLevelInBetweenAquifer, 0); // Assert Assert.Multiple(() => @@ -1142,8 +1228,7 @@ // Call bool isInBetweenLayerPresent = PlLinesToWaternetConverter.AreInBetweenAquiferClustersPresent(soilProfile, out int count); - IEnumerable topLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.TopLayerInBetweenAquiferCluster, soilProfile, 0); - IEnumerable bottomLevelInBetweenAquifer = PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.BottomLayerInBetweenAquiferCluster, soilProfile, 0); + PlLinesToWaternetConverter.DetermineLayerBoundaryPoints(PlLinesToWaternetConverter.LayerType.InBetweenAquiferCluster, soilProfile, out Point2D[] topLevelInBetweenAquifer, out Point2D[] bottomLevelInBetweenAquifer, 0); // Assert Assert.Multiple(() =>