[cig-commits] commit: Rearrange and cleanup

Mercurial hg at geodynamics.org
Sun Mar 25 20:02:37 PDT 2012


changeset:   118:dadba81bf86f
tag:         tip
user:        Walter Landry <wlandry at caltech.edu>
date:        Sun Mar 25 20:02:27 2012 -0700
files:       compute_Cxyz.cxx
description:
Rearrange and cleanup


diff -r e5469bbb9139 -r dadba81bf86f compute_Cxyz.cxx
--- a/compute_Cxyz.cxx	Sun Mar 25 19:43:20 2012 -0700
+++ b/compute_Cxyz.cxx	Sun Mar 25 20:02:27 2012 -0700
@@ -33,6 +33,8 @@ void compute_Cxyz(const double zx[Nx+1][
   const FTensor::Index<'c',2> c;
   pos(a)=interface.grid_pos(a);
 
+  FTensor::Tensor1<int,2> xyz(0,0);
+  xyz(d)=1;
   bool try_corners(true);
   for(int dd=0;dd<2;++dd)
     {
@@ -43,17 +45,13 @@ void compute_Cxyz(const double zx[Nx+1][
              dir2, dir3=components of mixed derivative we are correcting.
              So for zy,xy with boundary in x, dir2=y, dir3=y. */
           FTensor::Tensor1<double,2> v, dv, ddv;
-          FTensor::Tensor1<int,2> xyz(0,0), dir(0,0);
-          xyz(d)=1;
-          dir(dd)=1;
           compute_v_on_interface(zx,zy,log_etax,log_etay,distx,disty,
                                  interface.pos[dd],d,
                                  v,dv,ddv);
 
-          
           FTensor::Tensor1<double,2> dx(0,0);
-          dx(dd)=interface.pos[dd](dd) - pos(dd)
-            - h*sign(interface.pos[dd](dd)-pos(dd));
+          dx(dd)=pos(dd) - interface.pos[dd](dd)
+            - h*sign(pos(dd) - interface.pos[dd](dd));
           FTensor::Tensor1<double,2> z_jump;
           FTensor::Tensor2<double,2,2> dz_jump;
           FTensor::Tensor3_christof<double,2,2> ddz_jump;
@@ -63,7 +61,7 @@ void compute_Cxyz(const double zx[Nx+1][
           compute_jumps(v,dv,ddv,z_jump,dz_jump,ddz_jump,p_jump,dp_jump);
 
           double dz_dd_correction=
-            -xyz(a)*(z_jump(a) - dz_jump(a,b)*dx(b)
+            -xyz(a)*(z_jump(a) + dz_jump(a,b)*dx(b)
                      + ddz_jump(a,b,c)*dx(b)*dx(c)/2)/(h*h);
 
           const int dz_dd_factor((pos(dd)>interface.pos[dd](dd) ? -1 : 1)
@@ -73,10 +71,8 @@ void compute_Cxyz(const double zx[Nx+1][
 
           if(d==dd && interface.intersect_dp)
             {
-              dx(dd)=(pos(dd)>interface.pos[dd](dd))
-                ? ((pos(dd)-h/2)-interface.pos[dd](dd))
-                : ((pos(dd)+h/2)-interface.pos[dd](dd));
-
+              dx(dd)=pos(dd) - interface.pos[dd](dd)
+                - (h/2)*sign(pos(dd) - interface.pos[dd](dd));
               C+=(p_jump + dp_jump(a)*dx(a))/h;
             }
 
@@ -86,9 +82,8 @@ void compute_Cxyz(const double zx[Nx+1][
               FTensor::Tensor1<int,2> dir2(0,0), dir3(0,0);
               dir2((d+1)%2)=1;
               dir3((dd+1)%2)=1;
-              dx(dd)=(pos(dd)>interface.pos[dd](dd))
-                ? ((pos(dd)-h/2)-interface.pos[dd](dd))
-                : ((pos(dd)+h/2)-interface.pos[dd](dd));
+              dx(dd)=pos(dd) - interface.pos[dd](dd)
+                - (h/2)*sign(pos(dd) - interface.pos[dd](dd));
 
               C+=-(dz_jump(a,b)*dir2(a)*dir3(b)
                    + ddz_jump(a,b,c)*dir2(a)*dir3(b)*dx(c))/h;



More information about the CIG-COMMITS mailing list