From a927ca27c2116ecf1c19ff91d8f72b3469d93a2a Mon Sep 17 00:00:00 2001 From: Benedikt Zoennchen Date: Tue, 22 Oct 2019 16:47:30 +0200 Subject: [PATCH] add tests to solve the eikonal equation on the unstructured mesh for 2 scenarios (bridge and filled chicken). --- .../resources/poly/filled_chicken.poly | 28 ++++++++ .../improver/distmesh/Parameters.java | 3 +- .../improver/eikmesh/gen/GenEikMesh.java | 10 ++- .../mesh/EikonalSolverFMMTriangulation.java | 31 ++++++++ .../potential/solver/TestFMMEikMesh.java | 72 +++++++++++++++++++ 5 files changed, 140 insertions(+), 4 deletions(-) create mode 100644 VadereMeshing/resources/poly/filled_chicken.poly diff --git a/VadereMeshing/resources/poly/filled_chicken.poly b/VadereMeshing/resources/poly/filled_chicken.poly new file mode 100644 index 000000000..f8dc8fbe8 --- /dev/null +++ b/VadereMeshing/resources/poly/filled_chicken.poly @@ -0,0 +1,28 @@ +# nVertices dimension nAttributes boundaryMarker +8 2 0 0 +# vertexId x y +1 0.500000 0.500000 +2 34.500000 0.500000 +3 26.000000 21.000000 +4 26.000000 41.000000 +5 9.000000 41.000000 +6 9.000000 21.000000 +7 0.500000 59.500000 +8 34.500000 59.500000 +# +# nSegments boundaryMarker +8 0 +# lineId vertexId1 vertexId2 +1 2 8 +2 8 7 +3 7 1 +4 1 2 +5 3 4 +6 4 5 +7 5 6 +8 6 3 +# +# nHoles +1 +# vertexId x y (of a vertex which lies inside the hole) +1 17.500000 31.000000 \ No newline at end of file diff --git a/VadereMeshing/src/org/vadere/meshing/mesh/triangulation/improver/distmesh/Parameters.java b/VadereMeshing/src/org/vadere/meshing/mesh/triangulation/improver/distmesh/Parameters.java index 2b2571df3..3da4ffd10 100755 --- a/VadereMeshing/src/org/vadere/meshing/mesh/triangulation/improver/distmesh/Parameters.java +++ b/VadereMeshing/src/org/vadere/meshing/mesh/triangulation/improver/distmesh/Parameters.java @@ -15,11 +15,12 @@ public class Parameters { public final static boolean uniform = false; public final static String method = "Distmesh"; // "Distmesh" or "Density" public final static double qualityMeasurement = 0.85; + public final static double qualityConvergence = 0.001; public final static double MINIMUM = 0.25; public final static double DENSITYWEIGHT = 2; public final static int NPOINTS = 100000; public final static int SAMPLENUMBER = 10; public final static int SAMPLEDIVISION = 10; public final static int SEGMENTDIVISION = 0; - public final static int MAX_NUMBER_OF_STEPS = 100; + public final static int MAX_NUMBER_OF_STEPS = 250; } diff --git a/VadereMeshing/src/org/vadere/meshing/mesh/triangulation/improver/eikmesh/gen/GenEikMesh.java b/VadereMeshing/src/org/vadere/meshing/mesh/triangulation/improver/eikmesh/gen/GenEikMesh.java index fcf2d6ba9..06b698a12 100644 --- a/VadereMeshing/src/org/vadere/meshing/mesh/triangulation/improver/eikmesh/gen/GenEikMesh.java +++ b/VadereMeshing/src/org/vadere/meshing/mesh/triangulation/improver/eikmesh/gen/GenEikMesh.java @@ -57,7 +57,7 @@ public class GenEikMesh private Collection fixPoints; private double scalingFactor; private double deps; - private static final int MAX_STEPS = 250; + private static final int MAX_STEPS = Parameters.MAX_NUMBER_OF_STEPS; private int nSteps; private double initialEdgeLen; private double maxMovement; @@ -301,6 +301,10 @@ public class GenEikMesh } + public Collection getFixVertices() { + return refiner.getFixPoints(); + } + public void initialize() { initializeStep(); } @@ -343,8 +347,8 @@ public class GenEikMesh public boolean isFinished() { synchronized (getMesh()) { - boolean converged = dQuality < 0.000000001; - return initializationFinished() && quality > Parameters.qualityMeasurement && converged || (maxMovement > 0 && maxMovement / initialEdgeLen < Parameters.DPTOL) || nSteps >= MAX_STEPS; + boolean converged = dQuality < Parameters.qualityConvergence; + return initializationFinished() && quality >= Parameters.qualityMeasurement && converged || (maxMovement > 0 && maxMovement / initialEdgeLen < Parameters.DPTOL) || nSteps >= MAX_STEPS; } } diff --git a/VadereSimulator/src/org/vadere/simulator/models/potential/solver/calculators/mesh/EikonalSolverFMMTriangulation.java b/VadereSimulator/src/org/vadere/simulator/models/potential/solver/calculators/mesh/EikonalSolverFMMTriangulation.java index dff69fe65..1f9632dd0 100644 --- a/VadereSimulator/src/org/vadere/simulator/models/potential/solver/calculators/mesh/EikonalSolverFMMTriangulation.java +++ b/VadereSimulator/src/org/vadere/simulator/models/potential/solver/calculators/mesh/EikonalSolverFMMTriangulation.java @@ -116,11 +116,42 @@ public class EikonalSolverFMMTriangulation point.distance(p)); } } + /** + * Constructor for certain target points. + * + * @param timeCostFunction the time cost function t(x). Note F(x) = 1 / t(x). + * @param triangulation the triangulation the propagating wave moves on. + * @param targetVertices Points where the propagating wave starts i.e. points that are part of the target area. + */ + public EikonalSolverFMMTriangulation(@NotNull final ITimeCostFunction timeCostFunction, + @NotNull final IIncrementalTriangulation triangulation, + @NotNull final Collection targetVertices + ) { + this.triangulation = triangulation; + this.calculationFinished = false; + this.timeCostFunction = timeCostFunction; + this.narrowBand = new PriorityQueue<>(pointComparator); + this.targetVertices = new HashSet<>(); + this.distFunc = p -> IDistanceFunction.createToTargetPoints(targetVertices).apply(p); + + for(V vertex : targetVertices) { + this.targetVertices.add(vertex); + for(V neighbouringVertices : getMesh().getAdjacentVertexIt(vertex)) { + this.targetVertices.add(neighbouringVertices); + } + } + } + /** * Constructor for certain target shapes. * diff --git a/VadereSimulator/tests/org/vadere/simulator/models/potential/solver/TestFMMEikMesh.java b/VadereSimulator/tests/org/vadere/simulator/models/potential/solver/TestFMMEikMesh.java index b05d53c70..91b6e50dd 100644 --- a/VadereSimulator/tests/org/vadere/simulator/models/potential/solver/TestFMMEikMesh.java +++ b/VadereSimulator/tests/org/vadere/simulator/models/potential/solver/TestFMMEikMesh.java @@ -1,5 +1,6 @@ package org.vadere.simulator.models.potential.solver; +import org.jetbrains.annotations.NotNull; import org.junit.Ignore; import org.junit.Test; import org.vadere.meshing.examples.MeshExamples; @@ -8,19 +9,31 @@ import org.vadere.meshing.mesh.gen.PFace; import org.vadere.meshing.mesh.gen.PHalfEdge; import org.vadere.meshing.mesh.gen.PMesh; import org.vadere.meshing.mesh.gen.PVertex; +import org.vadere.meshing.mesh.impl.PMeshPanel; +import org.vadere.meshing.mesh.impl.PSLG; import org.vadere.meshing.mesh.inter.IIncrementalTriangulation; +import org.vadere.meshing.mesh.triangulation.IEdgeLengthFunction; +import org.vadere.meshing.mesh.triangulation.improver.eikmesh.gen.GenEikMesh; +import org.vadere.meshing.mesh.triangulation.improver.eikmesh.impl.PEikMesh; +import org.vadere.meshing.utils.io.poly.MeshPSLGWriter; import org.vadere.meshing.utils.io.poly.MeshPolyReader; import org.vadere.meshing.utils.io.poly.MeshPolyWriter; +import org.vadere.meshing.utils.io.poly.PSLGGenerator; import org.vadere.simulator.models.potential.solver.calculators.EikonalSolver; import org.vadere.simulator.models.potential.solver.calculators.mesh.EikonalSolverFMMTriangulation; import org.vadere.simulator.models.potential.solver.timecost.UnitTimeCostFunction; +import org.vadere.util.geometry.GeometryUtils; import org.vadere.util.geometry.shapes.VPoint; import org.vadere.util.geometry.shapes.VRectangle; import org.vadere.util.logging.Logger; +import org.vadere.util.math.IDistanceFunction; import java.io.IOException; import java.io.InputStream; +import java.util.ArrayList; +import java.util.Collection; import java.util.Collections; +import java.util.stream.Collectors; public class TestFMMEikMesh { private static Logger log = Logger.getLogger(TestFMMEikMesh.class); @@ -57,4 +70,63 @@ public class TestFMMEikMesh { System.out.println(mesh.toPythonTriangulation(v -> triangulation.getMesh().getDoubleData(v, "potential"))); } + + @Ignore + @Test + public void testFilledChickenFMM() throws IOException { + testTriangulationFMM("/poly/filled_chicken.poly", new VPoint(2,2), 1.0); + } + + @Ignore + @Test + public void testBridge() throws IOException { + testTriangulationFMM("/poly/bridge.poly", new VPoint(42.0, 46.0), 2.0); + } + + private void testTriangulationFMM(@NotNull final String file, @NotNull final VPoint targetPoint, final double h0) throws IOException { + // 1. read the base PSLG-file + final InputStream inputStream = MeshExamples.class.getResourceAsStream(file); + PSLG pslg = PSLGGenerator.toPSLGtoVShapes(inputStream); + + // 2. generate the unstructured mesh using EikMesh + Collection fixPoints = new ArrayList<>(); + fixPoints.add(targetPoint); + + var eikMesh = new GenEikMesh<>( + IDistanceFunction.create(pslg.getSegmentBound(), pslg.getHoles()), + p -> h0, + fixPoints, + h0, + pslg.getBoundingBox(), + pslg.getAllPolygons(), + () -> new PMesh()); + + PMeshPanel panel = new PMeshPanel(eikMesh.getMesh(), 600, 800); + panel.display(); + + // long for: eikMesh.generate(); + while (!eikMesh.isFinished()) { + eikMesh.improve(); + panel.repaint(); + } + //eikMesh.generate(); + + // 3. solve the eikonal equation on the given mehsh + EikonalSolver solver = new EikonalSolverFMMTriangulation( + new UnitTimeCostFunction(), + eikMesh.getTriangulation(), + eikMesh.getFixVertices().stream().filter(v -> v.distance(targetPoint) < GeometryUtils.DOUBLE_EPS).collect(Collectors.toList())); + long ms = System.currentTimeMillis(); + log.info("start FFM"); + solver.initialize(); + log.info("FFM finished"); + log.info("time: " + (System.currentTimeMillis() - ms)); + + + // 4. print the result to the console i.e. standard out + MeshPolyWriter meshPolyWriter = new MeshPolyWriter<>(); + System.out.println(meshPolyWriter.to2DPoly(eikMesh.getMesh(), 1, i -> "potential", v -> false)); + + //System.out.println(eikMesh.getMesh().toPythonTriangulation(v -> eikMesh.getMesh().getDoubleData(v, "potential"))); + } } -- GitLab