[cig-commits] commit: Fix the way corners cutting is detected.

Mercurial hg at geodynamics.org
Tue Apr 3 11:23:18 PDT 2012


changeset:   146:0b2d53a5d42f
user:        Walter Landry <wlandry at caltech.edu>
date:        Mon Apr 02 20:04:00 2012 -0700
files:       compute_coefficients/Interface.hxx
description:
Fix the way corners cutting is detected.


diff -r 67eb40d1a628 -r 0b2d53a5d42f compute_coefficients/Interface.hxx
--- a/compute_coefficients/Interface.hxx	Mon Apr 02 20:03:01 2012 -0700
+++ b/compute_coefficients/Interface.hxx	Mon Apr 02 20:04:00 2012 -0700
@@ -32,11 +32,11 @@ public:
                                     const double dist[][ny])
   {
     intersect_sides[0]=(sign(dist[i][j])!=sign(dist[i+1][j])
-                     && sign(dist[i][j+1])!=sign(dist[i+1][j+1])
-                     && sign(dist[i][j])==sign(dist[i][j+1]));
+                        && sign(dist[i][j+1])!=sign(dist[i+1][j+1])
+                        && sign(dist[i][j])==sign(dist[i][j+1]));
     intersect_sides[1]=(sign(dist[i][j])!=sign(dist[i][j+1])
-                     && sign(dist[i+1][j])!=sign(dist[i+1][j+1])
-                     && sign(dist[i][j])==sign(dist[i+1][j]));
+                        && sign(dist[i+1][j])!=sign(dist[i+1][j+1])
+                        && sign(dist[i][j])==sign(dist[i+1][j]));
     int center(sign(center_dist));
     for(int n0=0;n0<2;++n0)
       for(int n1=0;n1<2;++n1)
@@ -44,11 +44,34 @@ public:
           intersect_corner[n0][n1]=(center!=sign(dist[i+n0][j+n1]));
           if(intersect_corner[n0][n1])
             {
-              double delta(center_dist-dist[i+n0][j+n1]);
-              corner_pos[n0][n1](0)=
-                grid_pos(0)+sign(n0-.5)*center_dist*h/(2*delta);
-              corner_pos[n0][n1](1)=
-                grid_pos(1)+sign(n1-.5)*center_dist*h/(2*delta);
+              FTensor::Tensor1<double,2> norm, tangent;
+              norm(0)=dist[i+1][j]-dist[i][j];
+              norm(1)=dist[i][j+1]-dist[i][j];
+              tangent(0)=-norm(1);
+              tangent(1)=norm(0);
+
+              /* TODO: this assumes straight boundaries */
+              if(std::fabs(tangent(0))>std::fabs(tangent(1)))
+                {
+                  /* Top/Bottom */
+                  corner_pos[n0][n1](0)=grid_pos(0)
+                    + sign(n0-0.5)*h*(0.5
+                                      - std::fabs(dist[i+n0][j+n1]
+                                                  /(dist[i+(n0+1)%2][j+n1]
+                                                    - dist[i+n0][j+n1])));
+                  corner_pos[n0][n1](1)=grid_pos(1) + sign(n1-0.5)*h/2;
+                }
+              else
+                {
+                  /* Left/Right */
+                  corner_pos[n0][n1](0)=grid_pos(0) + sign(n0-0.5)*h/2;
+
+                  corner_pos[n0][n1](1)=grid_pos(1)
+                    + sign(n1-0.5)*h*(0.5
+                                      - std::fabs(dist[i+n0][j+n1]
+                                                  /(dist[i+n0][j+(n1+1)%2]
+                                                    - dist[i+n0][j+n1])));
+                }
             }
         }
     for(int n0=0;n0<2;++n0)



More information about the CIG-COMMITS mailing list