[cig-commits] commit: First efforts to generalize the interface code. Uses FTensor
Mercurial
hg at geodynamics.org
Fri Mar 9 13:30:07 PST 2012
changeset: 75:7b5ea67a2e69
tag: tip
user: Walter Landry <wlandry at caltech.edu>
date: Fri Mar 09 11:16:05 2012 -0800
files: compute_Cx.cxx wscript
description:
First efforts to generalize the interface code. Uses FTensor
diff -r 114d616b8587 -r 7b5ea67a2e69 compute_Cx.cxx
--- a/compute_Cx.cxx Mon Mar 05 17:40:26 2012 -0800
+++ b/compute_Cx.cxx Fri Mar 09 11:16:05 2012 -0800
@@ -2,6 +2,7 @@
#include "Model.hxx"
#include <iostream>
#include "compute_v_on_interface.hxx"
+#include "FTensor.hpp"
void compute_Cx(const Model &model, const double zx[Nx+1][Ny],
const double zy[Nx][Ny+1],
@@ -26,19 +27,34 @@ void compute_Cx(const Model &model, cons
vx,vy,dvx_y,dvy_y,
dvx_yy,dvy_yy);
+ const int d(0);
+ FTensor::Tensor1<double,2> v(vx,0), dv(0,dvy_y), ddv(dvx_yy,0), norm(1,0),
+ tangent(0,1);
+ FTensor::Tensor2_antisymmetric<double,2> eps(1);
+ const FTensor::Index<'a',2> a;
+
double dx=(x>middle) ? (middle-(x-h))
: (middle-(x+h));
- double zx_jump=eta_jump*vx;
- double dzx_x_jump=eta_jump*dvy_y;
- double p_jump=-2*dvy_y*eta_jump;
- double dp_x_jump=2*dvx_yy*eta_jump;
- double dzx_xx_jump=
- -dvx_yy*eta_jump + dp_x_jump;
- double dzx_xx_correction=
- -(zx_jump + dx*dzx_x_jump
- + dx*dx*dzx_xx_jump/2)/(h*h);
- const int dzx_xx_sign(x>middle ? -1 : 1);
- Cx=dzx_xx_sign*2*dzx_xx_correction;
+ double z_jump=eta_jump*v(d);
+
+ double dz_n_jump=eta_jump*(dv(0)*norm(a) - dv(1)*tangent(a))*eps(d,a);
+ double dz_t_jump=eta_jump*(dv(0)*norm(d) + dv(1)*tangent(d));
+ double dz_x_jump=dz_n_jump*norm(d) + dz_t_jump*tangent(d);
+
+ double p_jump=-2*eta_jump*dv(1);
+
+ double dp_t_jump=-2*eta_jump*ddv(1);
+ double dp_n_jump=2*eta_jump*ddv(0);
+ double dp_x_jump=dp_n_jump*norm(d) + dp_t_jump*tangent(d);
+
+ double dz_tt=eta_jump*(ddv(0)*norm(d) + ddv(1)*tangent(d));
+ double dz_xx_jump=-dz_tt + dp_x_jump;
+
+ double dz_xx_correction=
+ -(z_jump - dx*dz_x_jump
+ + dx*dx*dz_xx_jump/2)/(h*h);
+ const int dz_xx_sign(x>middle ? -1 : 1);
+ Cx=dz_xx_sign*2*dz_xx_correction;
double eta=exp(log_etax[i][j]);
double delta=(h-std::fabs(dx))/h;
@@ -56,7 +72,7 @@ void compute_Cx(const Model &model, cons
double dzx_xx_dv=-eta_jump*dvx_yy_dv + dp_x_dv;
double dp_dv=-2*eta_jump*dvy_y_dv;
- dCx_dvx=-dzx_xx_sign*2*(eta_jump*dvx_dv + dx*dzx_x_dv
+ dCx_dvx=-dz_xx_sign*2*(eta_jump*dvx_dv + dx*dzx_x_dv
+ dx*dx*dzx_xx_dv/2)/(h*h);
if(x+h/2>middle && x-h/2<middle)
diff -r 114d616b8587 -r 7b5ea67a2e69 wscript
--- a/wscript Mon Mar 05 17:40:26 2012 -0800
+++ b/wscript Fri Mar 09 11:16:05 2012 -0800
@@ -1,6 +1,3 @@
-# top = '.'
-# out = 'build'
-
def options(opt):
opt.load('compiler_cxx')
@@ -18,6 +15,7 @@ def build(bld):
target = 'Gamr_disc',
# cxxflags = ['-std=c++0x','-g','-D_GLIBCXX_DEBUG','-Wall'],
- cxxflags = ['-std=c++0x','-O3','-Wall'],
+ cxxflags = ['-std=c++0x','-O3','-Wall','-Drestrict='],
lib = ['boost_filesystem', 'boost_system', 'gsl', 'gslcblas'],
+ includes = ['/home/boo/programs/FTensor-bzr'],
)
More information about the CIG-COMMITS
mailing list