[[This is a redo of my "cvs import" of 2002/06/11, this time using proper

cvs operations (commit/delete/add) to preserve the full CVS history of this
thorn.]]

This is a major cleanup/revision of AEIThorns/Exact.

Major user-visible changes:
* major expansion of doc/documentation.tex
* major expansion of documentation in param.ccl file
* rename all parameters, systematize spacetime/coordinate/parameter names
  (there is a perl script in par/convert-pars.pl to convert old parameter
  files to the new names)
* [from Mitica Vulcanov] many additions and fixes to
  cosmological solutions and Schwarzschild-Lemaitre
* fix stress-energy tensor computations so they work -- before they were
  all disabled in CVS (INCLUDES lines were commented out in interface.ccl)
  due to requiring excessive friendship with evolution thorns
  and/or public parameters; new code copies parameters to restricted
  grid scalars, which Cactus automagically "pushes" to friends
* added some more tests to testsuite, though these don't yet work fully

Additional internal changes:
* rename many Fortran subroutines (and a few C ones too)
  so their names start with the thorn name
  to reduce the chances of name collisions with other thorns
* move all metrics to subdirectory so the main source directory isn't
  so cluttered
* move two files containing subroutines which were never called
  (they didn't belong in this thorn, but somehow got into cvs by accident)
  into new archive/ directory
* some (small) improvements in efficiency -- the exact_model parameter
  is now decoded from a keyword (string) to an integer once at INITIAL,
  and that integer tested by the stress-energy tensor code,
  rather than requiring a separate series of string tests at each grid
  point (!) like the old stress-energy tensor code did

Removed Files:
  Comparison.c			moved to  ../archive/Comparison.c
  ComparisonSolutions.F 	moved to  ../archive/ComparisonSolutions.F


git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinInitialData/Exact/trunk@98 e296648e-0e4f-0410-bd07-d597d9acff87

2 files changed, 0 insertions, 514 deletions

diff --git a/src/Comparison.c b/src/Comparison.c
deleted file mode 100644
index bc5554a..0000000
--- a/src/Comparison.c
+++ /dev/null
@@ -1,331 +0,0 @@
-/* Please god comment me ... */
-
-#include "cctk.h"
-
-static int comp_initialized = 0;
-static int comp_active = -1;
-
-static int comp_metric = 0;
-static int comp_curv   = 0;
-static int comp_gauge  = 0;
-
-#define f_comparisonmetric FORTRAN_NAME(comparisonmetric_,COMPARISONMETRIC,comparisonmetric)
-void FMODIFIER f_comparisonmetric(FORT_ARGS_PROTO, PROTO_FORT_FIELDS,
-				  double *,
-				  double *,
-				  double *,
-				  double *,
-				  double *,
-				  double *);
-#define f_comparisoncurvature FORTRAN_NAME(comparisoncurvature_,COMPARISONCURVATURE,comparisoncurvature)
-void FMODIFIER f_comparisoncurvature(FORT_ARGS_PROTO, PROTO_FORT_FIELDS,
-				  double *,
-				  double *,
-				  double *,
-				  double *,
-				  double *,
-				  double *);
-#define f_comparisongauge FORTRAN_NAME(comparisongauge_,COMPARISONGAUGE,comparisongauge)
-void FMODIFIER f_comparisongauge(FORT_ARGS_PROTO, PROTO_FORT_FIELDS,
-				  double *,
-				  double *,
-				  double *,
-				  double *);
-
-static int comptrips=1;
-
-void Comparison(pGH *in) {
-  void DoComparison(pGH *high, pGH *med, pGH *low,
-		    int gfi, double *he, double *me, double *le);
-
-  int i,j,k;
-  double *h[6];
-  double *m[6];
-  double *l[6];
-
-  int compevery = 1;
-  
-  pGH *high, *med, *low;
-
-  /* Grab the approriate grid heirarchies */
-  if (in->convlevel != 0) return;
-
-  high = in;
-  med = GetGHbyLevel(1,in->level, in->mglevel);
-  low = GetGHbyLevel(2,in->level, in->mglevel);
-
-  /* Only compare when grids are aligned. */
-  if (med) compevery = compevery * 2;
-  if (low) compevery = compevery * 2;
-
-  if (iteration % compevery != 0) return;
-
-  /* See if I'm active */
-  if (comp_active < 0) {
-    comp_active = Contains("comparison","yes");
-  }
-  if (!comp_active) return;
-
-  if (!comp_initialized) {
-    pGF *tGF;
-    comp_initialized =1;
-    for (i=0; i<high->ngridFuncs;i++) {
-      tGF = high->gridFuncs[i];
-      if (Contains_bounded("compfields",tGF->name)) {
-	tGF->do_comparison = 1;
-	if (tGF->gfno >= GFI_GXX &&
-	    tGF->gfno <= GFI_GZZ)
-	  comp_metric = 1;
-
-	if (tGF->gfno >= GFI_HXX &&
-	    tGF->gfno <= GFI_HZZ)
-	  comp_curv = 1;
-
-	if (tGF->gfno == GFI_ALP   ||
-	    tGF->gfno == GFI_BETAX ||
-	    tGF->gfno == GFI_BETAY ||
-	    tGF->gfno == GFI_BETAZ)
-	  comp_gauge = 1;
-      }
-    }
-  }
-  if (!(comp_metric || comp_gauge || comp_curv)) {
-    comp_active = 0;
-    printf ("You must compare gauge, curvature, or metric");
-    return;
-  }
-
-  /* OK at this point we know we are going to do something, so
-     allocate the memory for the comparisons */
-  for (i=0;i<6;i++) {
-    h[i] = (double *)malloc(high->npoints*sizeof(double));
-    if (med)
-      m[i] = (double *)malloc( med->npoints*sizeof(double));
-    else
-      m[i] = NULL;
-    if (low) 
-      l[i] = (double *)malloc( low->npoints*sizeof(double));
-    else
-      l[i] = NULL;
-  }
-
-  EnableGFDataStorage(high, high->gridFuncs[high->ngridFuncs-1]);
-  if (med) 
-    EnableGFDataStorage( med,  med->gridFuncs[med->ngridFuncs-1]);
-  if (low) 
-    EnableGFDataStorage( low,  low->gridFuncs[low->ngridFuncs-1]);
-
-  if (comp_metric) {
-    SetupFortranArrays(high);
-    f_comparisonmetric(FORT_ARGS(high), PASS_FORT_FIELDS(high),
-		       h[0],h[1],h[2],h[3],h[4],h[5]);
-
-    if (med) {
-      SetupFortranArrays(med);
-      f_comparisonmetric(FORT_ARGS(med), PASS_FORT_FIELDS(med),
-			 m[0],m[1],m[2],m[3],m[4],m[5]);
-    }
-    
-    if (low) {
-      SetupFortranArrays(low);
-      f_comparisonmetric(FORT_ARGS(low), PASS_FORT_FIELDS(low),
-			 l[0],l[1],l[2],l[3],l[4],l[5]);
-    }
-
-    for (i=GFI_GXX; i<=GFI_GZZ;i++) {
-      if (high->gridFuncs[i]->do_comparison) {
-	DoComparison(high,med,low,i,
-		     h[i-GFI_GXX],m[i-GFI_GXX],l[i-GFI_GXX]);
-      }
-    }
-  }
-
-  if (comp_curv) {
-    SetupFortranArrays(high);
-    f_comparisoncurvature(FORT_ARGS(high), PASS_FORT_FIELDS(high),
-		       h[0],h[1],h[2],h[3],h[4],h[5]);
-
-    if (med) {
-      SetupFortranArrays(med);
-      f_comparisoncurvature(FORT_ARGS(med), PASS_FORT_FIELDS(med),
-			    m[0],m[1],m[2],m[3],m[4],m[5]);
-    } 
-
-    if (low) {
-      SetupFortranArrays(low);
-      f_comparisoncurvature(FORT_ARGS(low), PASS_FORT_FIELDS(low),
-			    l[0],l[1],l[2],l[3],l[4],l[5]);
-    }
-
-    for (i=GFI_HXX; i<=GFI_HZZ;i++) {
-      if (high->gridFuncs[i]->do_comparison) {
-	DoComparison(high,med,low,i,
-		     h[i-GFI_HXX],m[i-GFI_HXX],l[i-GFI_HXX]);
-      }
-    }
-  }
-
-  if (comp_gauge) {
-    SetupFortranArrays(high);
-    f_comparisongauge(FORT_ARGS(high), PASS_FORT_FIELDS(high),
-		       h[0],h[1],h[2],h[3]);
-
-    if (med) {
-      SetupFortranArrays(med);
-      f_comparisongauge(FORT_ARGS(med), PASS_FORT_FIELDS(med),
-			m[0],m[1],m[2],m[3]);
-    }
-
-    if (low) {
-      SetupFortranArrays(low);
-      f_comparisongauge(FORT_ARGS(low), PASS_FORT_FIELDS(low),
-			l[0],l[1],l[2],l[3]);
-    }
-
-    if (high->gridFuncs[GFI_ALP]->do_comparison) {
-      DoComparison(high,med,low,GFI_ALP,h[0],m[0],l[0]);
-    }
-
-    if (high->gridFuncs[GFI_BETAX]->do_comparison) {
-      DoComparison(high,med,low,GFI_BETAX,h[1],m[1],l[1]);
-    }
-
-    if (high->gridFuncs[GFI_BETAY]->do_comparison) {
-      DoComparison(high,med,low,GFI_BETAY,h[2],m[2],l[2]);
-    }
-
-    if (high->gridFuncs[GFI_BETAZ]->do_comparison) {
-      DoComparison(high,med,low,GFI_BETAZ,h[3],m[3],l[3]);
-    }
-  }
-
-  for (i=0;i<6;i++) {
-    free(h[i]); if (m[i]) free(m[i]); if (l[i]) free(l[i]);
-  }
-  DisableGFDataStorage(high, high->gridFuncs[high->ngridFuncs-1]);
-  if (med)
-    DisableGFDataStorage( med,  med->gridFuncs[med->ngridFuncs-1]);
-  if (low)
-    DisableGFDataStorage( low,  low->gridFuncs[low->ngridFuncs-1]);
-  comptrips ++ ;
-}
-
-void DoComparison(pGH *high, pGH *med, pGH *low, int gfi,
-		  double *he, double *me, double *le) {
-
-  pGF *ws, *cur;
-  double max[3], nm1[3], nm2[3];
-  double sigtop, sigbot, sig3w, sig10, sig21;
-
-  int i;
-
-  printf ("Comparing %s\n",high->gridFuncs[gfi]->name);
-
-  /* Output the analytic solution if we need it */
-  cur = high->gridFuncs[gfi]; ws = high->gridFuncs[high->ngridFuncs-1];
-  if (cur->do_1dio) {
-    ws->do_1dio = comptrips;
-    for (i=0;i<high->npoints;i++) 
-      ws->data[i] = he[i];
-    sprintf(ws->name,"%s_exact",cur->name);
-    IO_Write1D(high,ws);
-    ws->do_1dio = 0;
-    ws->lastio_it[1]--;
-  }
-
-  /* Diffs go into workspace */
-  for (i=0;i<high->npoints;i++) 
-    ws->data[i] = cur->data[i] - he[i];
-
-  if (cur->do_1dio) {
-    ws->do_1dio = comptrips;
-    sprintf(ws->name,"%s_diff",cur->name);
-    IO_Write1D(high,ws);
-    ws->do_1dio = 0;
-    ws->lastio_it[1]--;
-  }
-
-  max[0] = pGF_MaxVal(high,ws);
-  nm1[0] = pGF_Norm1 (high,ws);
-  nm2[0] = pGF_Norm2 (high,ws);
-
-  if (med) {
-    cur = med->gridFuncs[gfi]; ws = med->gridFuncs[med->ngridFuncs-1];
-    for (i=0;i<med->npoints;i++)
-      ws->data[i] = cur->data[i] - me[i];
-    if (cur->do_1dio) {
-      ws->do_1dio = comptrips;
-      sprintf(ws->name,"%s_diff",cur->name);
-      IO_Write1D(med,ws);
-      ws->do_1dio = 0;    
-      ws->lastio_it[1]--;
-    }
-
-    max[1] = pGF_MaxVal(med,ws);
-    nm1[1] = pGF_Norm1 (med,ws);
-    nm2[1] = pGF_Norm2 (med,ws);
-  }
-
-  if (low) {
-    cur = low->gridFuncs[gfi]; ws = low->gridFuncs[low->ngridFuncs-1];
-    for (i=0;i<low->npoints;i++)
-      ws->data[i] = cur->data[i] - le[i];
-    if (cur->do_1dio) {
-      ws->do_1dio = comptrips;
-      sprintf(ws->name,"%s_diff",cur->name);
-      IO_Write1D(low,ws);
-      ws->do_1dio = 0;
-      ws->lastio_it[1]--;
-    }
-    
-    max[2] = pGF_MaxVal(low,ws);
-    nm1[2] = pGF_Norm1 (low,ws);
-    nm2[2] = pGF_Norm2 (low,ws);
-  } else {
-    max[2] = 0.0; nm1[2] = 0.0; nm2[2] = 0.0;
-  }
-
-  if (low) {
-    printf (" ------ :   High     Med       Low\n");
-    printf ("Max Diff: %lf %lf %lf\n",
-	    max[0],max[1],max[2]);
-    printf ("NM1 Diff: %lf %lf %lf\n",
-	  nm1[0],nm1[1],nm1[2]);
-    printf ("NM2 Diff: %lf %lf %lf\n",
-	    nm2[0],nm2[1],nm2[2]);
-  } else {
-    if (med) {
-      printf (" ------ :   High     Med      \n");
-      printf ("Max Diff: %lf %lf\n",
-	      max[0],max[1]);
-      printf ("NM1 Diff: %lf %lf\n",
-	      nm1[0],nm1[1]);
-      printf ("NM2 Diff: %lf %lf\n",
-	      nm2[0],nm2[1]);
-    }
-  }
-
-
-  if (low) {
-    sigbot = nm2[0] - nm2[1];
-    if (sigbot == 0) sig3w = 0.0;
-    else sig3w = log(fabs((nm2[1]-nm2[2])/(nm2[0]-nm2[1])))/log(2.0);
-
-    if (nm2[1] == 0) sig21 = 0.0;
-    else sig21 = log(fabs(nm2[2]/nm2[1]))/log(2.0);
-  } else {
-    sig3w = 0.0; sig21 = 0;
-  }
-    
-  if (nm2[0] == 0) sig10 = 0.0;
-  else sig10 = log(fabs(nm2[1]/nm2[0]))/log(2.0);
-
-  if (low)
-    printf ("sigma: 3way %lf  hm %lf  ml %lf\n",
-	    sig3w, sig10, sig21);
-  else
-    if (med)
-      printf ("sigma: hm %lf\n",sig10);
-  
-}
- 
diff --git a/src/ComparisonSolutions.F b/src/ComparisonSolutions.F
deleted file mode 100644
index b602d9b..0000000
--- a/src/ComparisonSolutions.F
+++ /dev/null
@@ -1,183 +0,0 @@
-#include "cctk.h"
-#include "cctk_Arguments.h"
-
-      subroutine ComparisonMetric(CCTK_ARGUMENTS, 
-     $     gxx_ex, gxy_ex, gxz_ex,
-     $     gyy_ex, gyz_ex, gzz_ex)
-
-      implicit none
-
-      DECLARE_CCTK_ARGUMENTS
-
-      CCTK_REAL gxx_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL gxy_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL gxz_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL gyy_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL gyz_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL gzz_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      
-      integer i,j,k
-
-C     Dummy arguments of subroutine boostrotdata: these are calculated
-C     (at a point) but not used.
-
-      CCTK_REAL hxxjunk, hyyjunk, hzzjunk, 
-     $     hxyjunk, hyzjunk, hxzjunk,
-     $     dxgxxjunk, dxgyyjunk, dxgzzjunk, 
-     $     dxgxyjunk, dxgyzjunk, dxgxzjunk,
-     $     dygxxjunk, dygyyjunk, dygzzjunk, 
-     $     dygxyjunk, dygyzjunk, dygxzjunk,
-     $     dzgxxjunk, dzgyyjunk, dzgzzjunk, 
-     $     dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $     alpjunk, axjunk, ayjunk, azjunk, 
-     $     betaxjunk, betayjunk, betazjunk,
-     $     bxxjunk, bxyjunk, bxzjunk, 
-     $     byxjunk, byyjunk, byzjunk, 
-     $     bzxjunk, bzyjunk, bzzjunk
-
-C     Call boostrotdata pointwise. Most of what it calculates is
-C     thrown away, variables ending in ...junk.
-      do k=1,cctk_lsh(3)
-         do j=1,cctk_lsh(2)
-            do i=1,cctk_lsh(1)
-               call exactdata(x(i,j,k), y(i,j,k), z(i,j,k), cctk_time,
-     $              gxx_ex(i,j,k), gyy_ex(i,j,k), gzz_ex(i,j,k), 
-     $              gxy_ex(i,j,k), gyz_ex(i,j,k), gxz_ex(i,j,k),
-     $              hxxjunk, hyyjunk, hzzjunk, 
-     $              hxyjunk, hyzjunk, hxzjunk,
-     $              dxgxxjunk, dxgyyjunk, dxgzzjunk, 
-     $              dxgxyjunk, dxgyzjunk, dxgxzjunk,
-     $              dygxxjunk, dygyyjunk, dygzzjunk, 
-     $              dygxyjunk, dygyzjunk, dygxzjunk,
-     $              dzgxxjunk, dzgyyjunk, dzgzzjunk, 
-     $              dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $              alpjunk, axjunk, ayjunk, azjunk, 
-     $              betaxjunk, betayjunk, betazjunk,
-     $              bxxjunk, bxyjunk, bxzjunk, 
-     $              byxjunk, byyjunk, byzjunk, 
-     $              bzxjunk, bzyjunk, bzzjunk)
-            end do
-         end do
-      end do
-
-      return
-      end
-
-cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
-
-      subroutine ComparisonCurvature(CCTK_ARGUMENTS,
-     $     hxx_ex, hxy_ex, hxz_ex,
-     $     hyy_ex, hyz_ex, hzz_ex)
-
-      implicit none
-
-      DECLARE_CCTK_ARGUMENTS
-
-      CCTK_REAL hxx_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL hxy_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL hxz_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL hyy_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL hyz_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL hzz_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-
-      integer i,j,k
-
-      CCTK_REAL gxxjunk, gyyjunk, gzzjunk, 
-     $     gxyjunk, gyzjunk, gxzjunk,
-     $     dxgxxjunk, dxgyyjunk, dxgzzjunk, 
-     $     dxgxyjunk, dxgyzjunk, dxgxzjunk,
-     $     dygxxjunk, dygyyjunk, dygzzjunk, 
-     $     dygxyjunk, dygyzjunk, dygxzjunk,
-     $     dzgxxjunk, dzgyyjunk, dzgzzjunk, 
-     $     dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $     alpjunk, axjunk, ayjunk, azjunk, 
-     $     betaxjunk, betayjunk, betazjunk,
-     $     bxxjunk, bxyjunk, bxzjunk, 
-     $     byxjunk, byyjunk, byzjunk, 
-     $     bzxjunk, bzyjunk, bzzjunk
-
-      do k=1,cctk_lsh(3)
-         do j=1,cctk_lsh(2)
-            do i=1,cctk_lsh(1)
-               call exactdata(x(i,j,k), y(i,j,k), z(i,j,k), cctk_time,
-     $              gxxjunk, gyyjunk, gzzjunk, 
-     $              gxyjunk, gyzjunk, gxzjunk,
-     $              hxx_ex(i,j,k), hyy_ex(i,j,k), hzz_ex(i,j,k), 
-     $              hxy_ex(i,j,k), hyz_ex(i,j,k), hxz_ex(i,j,k),
-     $              dxgxxjunk, dxgyyjunk, dxgzzjunk, 
-     $              dxgxyjunk, dxgyzjunk, dxgxzjunk,
-     $              dygxxjunk, dygyyjunk, dygzzjunk, 
-     $              dygxyjunk, dygyzjunk, dygxzjunk,
-     $              dzgxxjunk, dzgyyjunk, dzgzzjunk, 
-     $              dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $              alpjunk, axjunk, ayjunk, azjunk, 
-     $              betaxjunk, betayjunk, betazjunk,
-     $              bxxjunk, bxyjunk, bxzjunk, 
-     $              byxjunk, byyjunk, byzjunk, 
-     $              bzxjunk, bzyjunk, bzzjunk)
-            end do
-         end do
-      end do
-      return
-      end
-
-cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
-
-c     Note the exact shift comes in even if the shift is not
-c     allocated (eg if shift is "none"). In this case just
-c     don't use it, since it won't be compared against acctk_lsh(2)thing.
-
-      subroutine ComparisonGauge(CCTK_ARGUMENTS,
-     $     alp_ex, betax_ex, betay_ex, betaz_ex)
-
-      implicit none
-
-      DECLARE_CCTK_ARGUMENTS
-
-      CCTK_REAL alp_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL betax_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL betay_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-      CCTK_REAL betaz_ex(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))
-
-      integer i,j,k
-
-      CCTK_REAL gxxjunk, gyyjunk, gzzjunk, 
-     $     gxyjunk, gyzjunk, gxzjunk,
-     $     dxgxxjunk, dxgyyjunk, dxgzzjunk, 
-     $     dxgxyjunk, dxgyzjunk, dxgxzjunk,
-     $     dygxxjunk, dygyyjunk, dygzzjunk, 
-     $     dygxyjunk, dygyzjunk, dygxzjunk,
-     $     dzgxxjunk, dzgyyjunk, dzgzzjunk, 
-     $     dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $     hxxjunk, hyyjunk, hzzjunk, 
-     $     hxyjunk, hyzjunk, hxzjunk,
-     $     axjunk, ayjunk, azjunk, 
-     $     bxxjunk, bxyjunk, bxzjunk, 
-     $     byxjunk, byyjunk, byzjunk, 
-     $     bzxjunk, bzyjunk, bzzjunk
-
-      do k=1,cctk_lsh(3)
-         do j=1,cctk_lsh(2)
-            do i=1,cctk_lsh(1)
-               call exactdata(x(i,j,k), y(i,j,k), z(i,j,k), cctk_time,
-     $              gxxjunk, gyyjunk, gzzjunk, 
-     $              gxyjunk, gyzjunk, gxzjunk,
-     $              hxxjunk, hyyjunk, hzzjunk, 
-     $              hxyjunk, hyzjunk, hxzjunk,
-     $              dxgxxjunk, dxgyyjunk, dxgzzjunk, 
-     $              dxgxyjunk, dxgyzjunk, dxgxzjunk,
-     $              dygxxjunk, dygyyjunk, dygzzjunk, 
-     $              dygxyjunk, dygyzjunk, dygxzjunk,
-     $              dzgxxjunk, dzgyyjunk, dzgzzjunk, 
-     $              dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $              alp_ex(i,j,k), axjunk, ayjunk, azjunk, 
-     $              betax_ex(i,j,k), betay_ex(i,j,k), betaz_ex(i,j,k),
-     $              bxxjunk, bxyjunk, bxzjunk, 
-     $              byxjunk, byyjunk, byzjunk, 
-     $              bzxjunk, bzyjunk, bzzjunk)
-            end do
-         end do
-      end do
-      
-      return
-      end