From fbf2276a3a2fde4829df9dee057b2e4a50c7caa2 Mon Sep 17 00:00:00 2001
From: Stephan Hilb <stephan@ecshi.net>
Date: Fri, 23 Jul 2021 20:22:20 +0200
Subject: [PATCH] implement mesh area

---
 src/mesh.jl |  8 ++++++++
 src/run.jl  | 10 ++++++----
 2 files changed, 14 insertions(+), 4 deletions(-)

diff --git a/src/mesh.jl b/src/mesh.jl
index 10bcd1d..86a9f23 100644
--- a/src/mesh.jl
+++ b/src/mesh.jl
@@ -340,6 +340,14 @@ function save(filename::String, mesh::Mesh, fs...)
     vtk_save(vtk)
 end
 
+function area(mesh, cell)
+    A = SArray{Tuple{ndims_space(mesh), nvertices_cell(mesh)}}(
+	view(mesh.vertices, :, view(mesh.cells, :, cell)))
+    return det(A[:, 2:3] .- A[:, 1]) / 2
+end
+
+area(mesh) = mapreduce(cell -> area(mesh, cell), +, cells(mesh))
+
 function diam(mesh, cell)
     A = SArray{Tuple{ndims_space(mesh), nvertices_cell(mesh)}}(
 	view(mesh.vertices, :, view(mesh.cells, :, cell)))
diff --git a/src/run.jl b/src/run.jl
index 39a4a51..a921111 100644
--- a/src/run.jl
+++ b/src/run.jl
@@ -310,12 +310,14 @@ function denoise(img; name, params...)
 	    estimate!(ctx)
 	    pvd[k] = save_denoise(ctx, k)
 	    println()
+
+	    norm_step = sqrt((norm_l2(ctx.du)^2 + norm_l2(ctx.dp1)^2 + norm_l2(ctx.dp2)^2) / area(mesh))
+
 	    println("ndofs: $(ndofs(ctx.u.space)), est: $(norm_l2(ctx.est)))")
 	    println("primal energy: $(primal_energy(ctx))")
+	    println("norm_step: $(norm_step)")
 
-	    norm_step = sqrt(norm_l2(ctx.du)^2 + norm_l2(ctx.dp1)^2 + norm_l2(ctx.dp2))
-	    norm_step = sqrt(norm_l2(ctx.du)^2)
-            norm_step <= 1e-2 && break
+            norm_step <= 1e-1 && break
 	end
 	marked_cells = mark(ctx; theta = 0.5)
 	#println(marked_cells)
@@ -327,7 +329,7 @@ function denoise(img; name, params...)
 	ctx.g.data .= gnew.data
 	#interpolate!(ctx.g, x -> norm(x .- m) < norm(m .- 1) / 3)
 
-	k >= 50 && break
+	k >= 100 && break
     end
     vtk_save(pvd)
     return ctx
-- 
GitLab