[cig-commits] commit: Correctly switch the direction for the velocity of the corner correction.

[Mercurial]

changeset:   134:546f3776347a
tag:         tip
user:        Walter Landry <wlandry at caltech.edu>
date:        Thu Mar 29 11:50:48 2012 -0700
files:       compute_Cxyz.cxx
description:
Correctly switch the direction for the velocity of the corner correction.


diff -r 550bf5907a58 -r 546f3776347a compute_Cxyz.cxx
--- a/compute_Cxyz.cxx	Thu Mar 29 11:48:01 2012 -0700
+++ b/compute_Cxyz.cxx	Thu Mar 29 11:50:48 2012 -0700
@@ -146,43 +146,48 @@ void compute_Cxyz(const double zx[Nx+1][
         }
   /* Handle corner cutting */
   if(try_corners)
-    for(int dd=0;dd<2;++dd)
-      for(int ee=0;ee<2;++ee)
-        if(interface.intersect_corner[dd][ee])
-          {
-            FTensor::Tensor1<double,2> v, dv, ddv;
-            FTensor::Tensor2<double,2,2> nt_to_xy;
-            compute_v_on_interface(zx,zy,log_etax,log_etay,distx,disty,
-                                   interface.corner_pos[dd][ee],d,
-                                   v,dv,ddv,nt_to_xy);
+    {
+      FTensor::Tensor1<double,2> zyx(0,0);
+      zyx((d+1)%2)=1;
 
-            FTensor::Tensor1<double,2> z_jump;
-            FTensor::Tensor2<double,2,2> dz_jump;
-            FTensor::Tensor3_christof<double,2,2> ddz_jump;
+      for(int dd=0;dd<2;++dd)
+        for(int ee=0;ee<2;++ee)
+          if(interface.intersect_corner[dd][ee])
+            {
+              FTensor::Tensor1<double,2> v, dv, ddv;
+              FTensor::Tensor2<double,2,2> nt_to_xy;
+              compute_v_on_interface(zx,zy,log_etax,log_etay,distx,disty,
+                                     interface.corner_pos[dd][ee],d,
+                                     v,dv,ddv,nt_to_xy);
 
-            double eta_jump;
-            switch(d)
-              {
-              case 0:
-                eta_jump=exp(log_etax[i][j])-exp(log_etay[i-1+dd][j+ee]);
-                break;
+              FTensor::Tensor1<double,2> z_jump;
+              FTensor::Tensor2<double,2,2> dz_jump;
+              FTensor::Tensor3_christof<double,2,2> ddz_jump;
 
-              case 1:
-                eta_jump=exp(log_etay[i][j])-exp(log_etax[i+dd][j-1+ee]);
-                break;
-              default:
-                abort();
-                break;
-              }
-            compute_jumps(eta_jump,v,dv,ddv,nt_to_xy,z_jump,dz_jump,ddz_jump);
+              double eta_jump;
+              switch(d)
+                {
+                case 0:
+                  eta_jump=exp(log_etax[i][j])-exp(log_etay[i-1+dd][j+ee]);
+                  break;
 
-            FTensor::Tensor1<double,2> dx;
-            for(int aa=0;aa<2;++aa)
-              dx(aa)=pos(aa) - interface.corner_pos[dd][ee](aa)
-                - (h/2)*sign(pos(aa) - interface.corner_pos[dd][ee](aa));
+                case 1:
+                  eta_jump=exp(log_etay[i][j])-exp(log_etax[i+dd][j-1+ee]);
+                  break;
+                default:
+                  abort();
+                  break;
+                }
+              compute_jumps(eta_jump,v,dv,ddv,nt_to_xy,z_jump,dz_jump,ddz_jump);
 
-            C+=-sign(dx(0))*sign(dx(1))*xyz(a)
-              *(z_jump(a) + dz_jump(a,b)*dx(b)
-                + ddz_jump(a,b,c)*dx(b)*dx(c)/2)/(h*h);
-          }
+              FTensor::Tensor1<double,2> dx;
+              for(int aa=0;aa<2;++aa)
+                dx(aa)=pos(aa) - interface.corner_pos[dd][ee](aa)
+                  - (h/2)*sign(pos(aa) - interface.corner_pos[dd][ee](aa));
+
+              C+=sign(dx(0))*sign(dx(1))*zyx(a)
+                *(z_jump(a) + dz_jump(a,b)*dx(b)
+                  + ddz_jump(a,b,c)*dx(b)*dx(c)/2)/(h*h);
+            }
+    }
 }