Support ADMBase's dtlapse and dtshift variables.

Parallelise the main loops with OpenMP.

Clean up the parameter file syntax.


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

8 files changed, 175 insertions, 221 deletions

diff --git a/param.ccl b/param.ccl
index 90772f6..1c03d3f 100644
--- a/param.ccl
+++ b/param.ccl
@@ -24,14 +24,22 @@ EXTENDS KEYWORD shift_evolution_method
   "exact" :: "Use shift from exact solution"
 }
 
-EXTENDS KEYWORD initial_lapse 
+EXTENDS KEYWORD initial_lapse
 {
   "exact" :: "Initial lapse from exact solution"
 }
-EXTENDS KEYWORD initial_shift 
+EXTENDS KEYWORD initial_shift
 {
   "exact" :: "Initial shift from exact solution"
 }
+EXTENDS KEYWORD initial_dtlapse
+{
+  "exact" :: "Initial dtlapse from exact solution"
+}
+EXTENDS KEYWORD initial_dtshift
+{
+  "exact" :: "Initial dtshift from exact solution"
+}
 
 EXTENDS KEYWORD evolution_method
 {
@@ -63,12 +71,12 @@ BOOLEAN exblend_gauge "Blend the lapse and shift with the exact solution?"
 
 REAL exblend_rout "Outer boundary of blending region"
 {
-: :: "Positive means radial value, negative means use outer bound of grid"
+*:* :: "Positive means radial value, negative means use outer bound of grid"
 } -1.0
 
 REAL exblend_width "Width of blending zone"
 {
-: :: "Positive means width in radius, negative means width = exbeldn_width*dx"
+*:* :: "Positive means width in radius, negative means width = exbeldn_width*dx"
 } -3.0
 
 ################################################################################
@@ -79,12 +87,12 @@ private:
 
 REAL slice_Gauss_ampl "Amplitude of Gauss slice in exact"
 {
-0.0: :: "Positive please"
+0.0:* :: "Positive please"
 } 0.0
 
 REAL slice_Gauss_width "Width of Gauss slice in exact"
 {
-0.0: :: "Positive please"
+0.0:* :: "Positive please"
 } 1.0
 
 # this parameter is used directly in our schedule.ccl file
@@ -603,7 +611,7 @@ REAL multi_BH__z3 "multi-BH: z coord of black hole number 3"
 # black hole #4
 REAL multi_BH__mass4 "multi-BH: mass of black hole number 4"
 {
-0.0: :: "any real number >= 0"
+0.0:* :: "any real number >= 0"
 } 0.0
 REAL multi_BH__x4 "multi-BH: x coord of black hole number 4"
 {
@@ -626,12 +634,12 @@ REAL multi_BH__z4 "multi-BH: z coord of black hole number 4"
 
 REAL Alvi__mass1 "Alvi: mass of BH number 1"
 {
-0.0: :: "any real number >= 0"
+0.0:* :: "any real number >= 0"
 } 1.0
 
 REAL Alvi__mass2 "Alvi: mass of BH number 2"
 {
-0.0: :: "any real number >= 0"
+0.0:* :: "any real number >= 0"
 } 1.0
 
 REAL Alvi__separation "Alvi: spatial separation of the black holes"
@@ -872,19 +880,19 @@ REAL Gowdy_wave__amplitude "Gowdy-wave: amplitude parameter"
 
 REAL boost_rotation_symmetric__scale "boost-rotation symmetric: length scale"
 {
-0.0: :: "Positive please"
+0.0:* :: "Positive please"
 } 1.0
 
 REAL boost_rotation_symmetric__amp \
   "boost-rotation symmetric: dimensionless amplitude"
 {
-0.0: :: "Positive please"
+0.0:* :: "Positive please"
 } 0.1
 
 REAL boost_rotation_symmetric__min_d \
   "boost-rotation symmetric: dimensionless safety distance"
 {
-(0.0: :: "any positive real number"
+(0.0:* :: "any positive real number"
 } 0.01
 
 ################################################################################
diff --git a/schedule.ccl b/schedule.ccl
index 59b4a50..7146bf4 100644
--- a/schedule.ccl
+++ b/schedule.ccl
@@ -64,7 +64,7 @@ if ( (CCTK_Equals(initial_data,"slice")) && ! (CCTK_Equals(evolution_method,"sli
 # Lapse and/or shift from a trivial slice through an exact solution.
 # At the initial time:
 
-if ( (CCTK_Equals(initial_lapse,"exact")) || (CCTK_Equals(initial_shift,"exact")) )
+if ( (CCTK_Equals(initial_lapse,"exact")) || (CCTK_Equals(initial_shift,"exact")) || (CCTK_Equals(initial_dtlapse,"exact")) || (CCTK_Equals(initial_dtshift,"exact")) )
 {
   schedule Exact__gauge in ADMBase_InitialGauge
   {
diff --git a/src/Bona_Masso_data.F77 b/src/Bona_Masso_data.F77
index 02ea72e..e79db6f 100644
--- a/src/Bona_Masso_data.F77
+++ b/src/Bona_Masso_data.F77
@@ -18,7 +18,8 @@ c $Header$
      $     dxgxx, dxgyy, dxgzz, dxgxy, dxgyz, dxgxz,
      $     dygxx, dygyy, dygzz, dygxy, dygyz, dygxz,
      $     dzgxx, dzgyy, dzgzz, dzgxy, dzgyz, dzgxz,
-     $     alp, ax, ay, az, betax, betay, betaz,
+     $     alp, dtalp, ax, ay, az,
+     $     betax, betay, betaz, dtbetax, dtbetay, dtbetaz,
      $     bxx, bxy, bxz, byx, byy, byz, bzx, bzy, bzz)
 
       implicit none
@@ -32,7 +33,8 @@ c $Header$
      $       dxgxx, dxgyy, dxgzz, dxgxy, dxgyz, dxgxz,
      $       dygxx, dygyy, dygzz, dygxy, dygyz, dygxz,
      $       dzgxx, dzgyy, dzgzz, dzgxy, dzgyz, dzgxz,
-     $       alp, ax, ay, az, betax, betay, betaz,
+     $       alp, dtalp, ax, ay, az,
+     $       betax, betay, betaz, dtbetax, dtbetay, dtbetaz,
      $       bxx, bxy, bxz, byx, byy, byz, bzx, bzy, bzz
 
 C     gxx is g_xx etc.
@@ -253,6 +255,21 @@ C     Calculate t-derivatives, and extrinsic curvature.
      $     + bxx * gxz + bxy * gyz + bxz * gzz
      $     + bzx * gxx + bzy * gxy + bzz * gxz) / alp
 
+C     Calculate time derivatives of  lapse and shift
+
+C     alp = 1.d0 / sqrt(-gutt)
+      dtalp = 0.5d0 / sqrt(-gutt)**3 * (gutt_p - gutt_m) / 2.d0
+
+C     betax = - gutx / gutt
+C     betay = - guty / gutt
+C     betaz = - gutz / gutt
+      dtbetax = - (gutx_p - gutx_m) / gutt / 2.d0
+     $          + gutx * (gutt_p - gutt_m) / gutt**2 / 2.d0
+      dtbetay = - (guty_p - guty_m) / gutt / 2.d0
+     $          + guty * (gutt_p - gutt_m) / gutt**2 / 2.d0
+      dtbetaz = - (gutz_p - gutz_m) / gutt / 2.d0
+     $          + gutz * (gutt_p - gutt_m) / gutt**2 / 2.d0
+
 C     Calculate second derivatives of the conformal factor
 
       if (psi_on) then
diff --git a/src/blended_boundary.F77 b/src/blended_boundary.F77
index dfc3185..752ce52 100644
--- a/src/blended_boundary.F77
+++ b/src/blended_boundary.F77
@@ -31,8 +31,8 @@ C $Header$
      $     exact_psixx, exact_psiyy, exact_psizz,
      $     exact_psixy, exact_psiyz, exact_psixz
 
-      CCTK_REAL alpe, axe, aye, aze
-      CCTK_REAL betaxe,betaye,betaze
+      CCTK_REAL alpe, dtalpe, axe, aye, aze
+      CCTK_REAL betaxe,betaye,betaze, dtbetaxe,dtbetaye,dtbetaze
       CCTK_REAL bxxe,bxye,bxze,byxe,byye,byze,bzxe,bzye,bzze
 
       CCTK_REAL det, uxx, uxy, uxz, uyy, uyz, uzz
@@ -119,7 +119,8 @@ C     (also to tell the models about the conformal_state)
      $                 dxgxxe, dxgyye, dxgzze, dxgxye, dxgyze, dxgxze,
      $                 dygxxe, dygyye, dygzze, dygxye, dygyze, dygxze,
      $                 dzgxxe, dzgyye, dzgzze, dzgxye, dzgyze, dzgxze,
-     $                 alpe, axe, aye, aze, betaxe, betaye, betaze,
+     $                 alpe, dtalpe, axe, aye, aze,
+     $                 betaxe, betaye, betaze, dtbetaxe, dtbetaye, dtbetaze,
      $                 bxxe, bxye, bxze, byxe, 
      $                 byye, byze, bzxe, bzye, bzze)
 
diff --git a/src/boundary.F77 b/src/boundary.F77
index 323de80..31ee55c 100644
--- a/src/boundary.F77
+++ b/src/boundary.F77
@@ -70,28 +70,30 @@ C     time t+1. This is to provide boundary conditions for testing
 C     elliptic gauge conditions. If they are not used, they will be
 C     overwritten by Exact__gauge.
 
-#define EXACTDATAPOINT                                          \
-               call Exact__Bona_Masso_data(                     \
-                   decoded_exact_model,                         \
-                   x(i,j,k), y(i,j,k), z(i,j,k), tplusone,      \
-                   gxx(i,j,k), gyy(i,j,k), gzz(i,j,k),          \
-                   gxy(i,j,k), gyz(i,j,k), gxz(i,j,k),          \
-                   kxx(i,j,k), kyy(i,j,k), kzz(i,j,k),          \
-                   kxy(i,j,k), kyz(i,j,k), kxz(i,j,k),          \
-                   exact_psi,                                   \
-                   exact_psix, exact_psiy, exact_psiz,          \
-                   exact_psixx, exact_psiyy, exact_psizz,       \
-                   exact_psixy, exact_psiyz, exact_psixz,       \
-                   dxgxxjunk, dxgyyjunk, dxgzzjunk,             \
-                   dxgxyjunk, dxgyzjunk, dxgxzjunk,             \
-                   dygxxjunk, dygyyjunk, dygzzjunk,             \
-                   dygxyjunk, dygyzjunk, dygxzjunk,             \
-                   dzgxxjunk, dzgyyjunk, dzgzzjunk,             \
-                   dzgxyjunk, dzgyzjunk, dzgxzjunk,             \
-                   alp(i,j,k), axjunk, ayjunk, azjunk,          \
-                   betax(i,j,k), betay(i,j,k), betaz(i,j,k),    \
-                   bxxjunk, bxyjunk, bxzjunk,                   \
-                   byxjunk, byyjunk, byzjunk,                   \
+#define EXACTDATAPOINT                                                  \
+               call Exact__Bona_Masso_data(                             \
+                   decoded_exact_model,                                 \
+                   x(i,j,k), y(i,j,k), z(i,j,k), tplusone,              \
+                   gxx(i,j,k), gyy(i,j,k), gzz(i,j,k),                  \
+                   gxy(i,j,k), gyz(i,j,k), gxz(i,j,k),                  \
+                   kxx(i,j,k), kyy(i,j,k), kzz(i,j,k),                  \
+                   kxy(i,j,k), kyz(i,j,k), kxz(i,j,k),                  \
+                   exact_psi,                                           \
+                   exact_psix, exact_psiy, exact_psiz,                  \
+                   exact_psixx, exact_psiyy, exact_psizz,               \
+                   exact_psixy, exact_psiyz, exact_psixz,               \
+                   dxgxxjunk, dxgyyjunk, dxgzzjunk,                     \
+                   dxgxyjunk, dxgyzjunk, dxgxzjunk,                     \
+                   dygxxjunk, dygyyjunk, dygzzjunk,                     \
+                   dygxyjunk, dygyzjunk, dygxzjunk,                     \
+                   dzgxxjunk, dzgyyjunk, dzgzzjunk,                     \
+                   dzgxyjunk, dzgyzjunk, dzgxzjunk,                     \
+                   alp(i,j,k), dtalp(i,j,k),                            \
+                   axjunk, ayjunk, azjunk,                              \
+                   betax(i,j,k), betay(i,j,k), betaz(i,j,k),            \
+                   dtbetax(i,j,k), dtbetay(i,j,k), dtbetaz(i,j,k),      \
+                   bxxjunk, bxyjunk, bxzjunk,                           \
+                   byxjunk, byyjunk, byzjunk,                           \
                    bzxjunk, bzyjunk, bzzjunk)
 
       if (cctk_bbox(1) .eq. 1) then 
diff --git a/src/gauge.F77 b/src/gauge.F77
index b10bfa1..f8b35e2 100644
--- a/src/gauge.F77
+++ b/src/gauge.F77
@@ -18,23 +18,25 @@ C $Header$
 
       integer i,j,k
       integer nx,ny,nz
+      logical set_lapse, set_dtlapse, set_shift, set_dtshift
 
       CCTK_REAL tt, xx, yy, zz
-      CCTK_REAL 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,
-     $     alpjunk, axjunk, ayjunk, azjunk, 
-     $     betaxjunk, betayjunk, betazjunk,
-     $     bxxjunk, bxyjunk, bxzjunk, 
-     $     byxjunk, byyjunk, byzjunk, 
-     $     bzxjunk, bzyjunk, bzzjunk
+      CCTK_REAL gxxtmp, gyytmp, gzztmp, 
+     $     gxytmp, gyztmp, gxztmp,
+     $     hxxtmp, hyytmp, hzztmp, 
+     $     hxytmp, hyztmp, hxztmp,
+     $     dxgxxtmp, dxgyytmp, dxgzztmp, 
+     $     dxgxytmp, dxgyztmp, dxgxztmp,
+     $     dygxxtmp, dygyytmp, dygzztmp, 
+     $     dygxytmp, dygyztmp, dygxztmp,
+     $     dzgxxtmp, dzgyytmp, dzgzztmp, 
+     $     dzgxytmp, dzgyztmp, dzgxztmp,
+     $     alptmp, dtalptmp, axtmp, aytmp, aztmp, 
+     $     betaxtmp, betaytmp, betaztmp,
+     $     dtbetaxtmp, dtbetaytmp, dtbetaztmp,
+     $     bxxtmp, bxytmp, bxztmp, 
+     $     byxtmp, byytmp, byztmp, 
+     $     bzxtmp, bzytmp, bzztmp
       CCTK_REAL
      $     exact_psi,
      $     exact_psix, exact_psiy, exact_psiz,
@@ -69,18 +71,41 @@ C     is why this remains a bug for now... :( :(
 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
 
 C
-C     Set both lapse and shift?
+C     Set lapse and/or shift?
 C
-      if (
-     $     ( is_initial_slice
-     $       .and. (CCTK_Equals(initial_lapse,         "exact").ne.0)
-     $       .and. (CCTK_Equals(initial_shift,         "exact").ne.0) )
-     $     .or.
-     $     ( is_later_slice
-     $       .and. (CCTK_Equals(lapse_evolution_method,"exact").ne.0)
-     $       .and. (CCTK_Equals(shift_evolution_method,"exact").ne.0) )
-     $   ) then
+      if ( is_initial_slice ) then
+         set_lapse   = CCTK_Equals(initial_lapse,   "exact").ne.0
+         set_shift   = CCTK_Equals(initial_shift,   "exact").ne.0
+         set_dtlapse = CCTK_Equals(initial_dtlapse, "exact").ne.0
+         set_dtshift = CCTK_Equals(initial_dtshift, "exact").ne.0
+      end if
+      if ( is_later_slice ) then
+         set_lapse   = CCTK_Equals(lapse_evolution_method, "exact").ne.0
+         set_shift   = CCTK_Equals(shift_evolution_method, "exact").ne.0
+         set_dtlapse = .false.
+         set_dtshift = .false.
+      end if
 
+      if ( set_lapse .or. set_shift .or. set_dtlapse .or. set_dtshift) then
+
+C$omp    parallel do private (
+C$omp$     tt, xx, yy, zz,
+C$omp$     alptmp, dtalptmp, axtmp, aytmp, aztmp, 
+C$omp$     betaxtmp, betaytmp, betaztmp,
+C$omp$     dtbetaxtmp, dtbetaytmp, dtbetaztmp,
+C$omp$     bxxtmp, bxytmp, bxztmp, 
+C$omp$     byxtmp, byytmp, byztmp, 
+C$omp$     bzxtmp, bzytmp, bzztmp,
+C$omp$     dxgxxtmp, dxgyytmp, dxgzztmp, 
+C$omp$     dxgxytmp, dxgyztmp, dxgxztmp,
+C$omp$     dygxxtmp, dygyytmp, dygzztmp, 
+C$omp$     dygxytmp, dygyztmp, dygxztmp,
+C$omp$     dzgxxtmp, dzgyytmp, dzgzztmp, 
+C$omp$     dzgxytmp, dzgyztmp, dzgxztmp,
+C$omp$     exact_psi,
+C$omp$     exact_psix, exact_psiy, exact_psiz,
+C$omp$     exact_psixx, exact_psiyy, exact_psizz,
+C$omp$     exact_psixy, exact_psiyz, exact_psixz)
          do k=1,nz
             do j=1,ny
                do i=1,nx
@@ -110,163 +135,43 @@ C     (also to tell the models about the conformal_state)
                   call Exact__Bona_Masso_data(
      $                 decoded_exact_model,
      $                 xx, yy, zz, tt,
-     $                 gxxjunk, gyyjunk, gzzjunk, 
-     $                 gxyjunk, gyzjunk, gxzjunk,
-     $                 hxxjunk, hyyjunk, hzzjunk, 
-     $                 hxyjunk, hyzjunk, hxzjunk,
+     $                 gxxtmp, gyytmp, gzztmp, 
+     $                 gxytmp, gyztmp, gxztmp,
+     $                 hxxtmp, hyytmp, hzztmp, 
+     $                 hxytmp, hyztmp, hxztmp,
      $                 exact_psi,
      $                 exact_psix, exact_psiy, exact_psiz,
      $                 exact_psixx, exact_psiyy, exact_psizz,
      $                 exact_psixy, exact_psiyz, exact_psixz,
-     $                 dxgxxjunk, dxgyyjunk, dxgzzjunk, 
-     $                 dxgxyjunk, dxgyzjunk, dxgxzjunk,
-     $                 dygxxjunk, dygyyjunk, dygzzjunk, 
-     $                 dygxyjunk, dygyzjunk, dygxzjunk,
-     $                 dzgxxjunk, dzgyyjunk, dzgzzjunk, 
-     $                 dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $                 alp(i,j,k), axjunk, ayjunk, azjunk, 
-     $                 betax(i,j,k), betay(i,j,k), betaz(i,j,k),
-     $                 bxxjunk, bxyjunk, bxzjunk, 
-     $                 byxjunk, byyjunk, byzjunk, 
-     $                 bzxjunk, bzyjunk, bzzjunk)
-                  betax(i,j,k) = betax(i,j,k) + shift_add_x
-                  betay(i,j,k) = betay(i,j,k) + shift_add_y
-                  betaz(i,j,k) = betaz(i,j,k) + shift_add_z
-               end do
-            end do
-         end do
-
-CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
-
-C
-C     Set lapse only?
-C
-      elseif (
-     $         ( is_initial_slice
-     $           .and. (CCTK_Equals(initial_lapse,         "exact").ne.0) )
-     $         .or.
-     $         ( is_later_slice
-     $           .and. (CCTK_Equals(lapse_evolution_method,"exact").ne.0) )
-     $       ) then
-
-         do k=1,nz
-            do j=1,ny
-               do i=1,nx
-
-                  tt = cctk_time
-                  xx = x(i,j,k) - cctk_time * shift_add_x
-                  yy = y(i,j,k) - cctk_time * shift_add_y
-                  zz = z(i,j,k) - cctk_time * shift_add_z
-
-C     Initialize the psi of exact
-C     (also to tell the models about the conformal_state)
-                  if (conformal_state .ne. 0) then
-                    exact_psi = 1.0D0
-                  else
-                    exact_psi = 0.0D0
+     $                 dxgxxtmp, dxgyytmp, dxgzztmp, 
+     $                 dxgxytmp, dxgyztmp, dxgxztmp,
+     $                 dygxxtmp, dygyytmp, dygzztmp, 
+     $                 dygxytmp, dygyztmp, dygxztmp,
+     $                 dzgxxtmp, dzgyytmp, dzgzztmp, 
+     $                 dzgxytmp, dzgyztmp, dzgxztmp,
+     $                 alptmp, dtalptmp, axtmp, aytmp, aztmp, 
+     $                 betaxtmp, betaytmp, betaztmp,
+     $                 dtbetaxtmp, dtbetaytmp, dtbetaztmp,
+     $                 bxxtmp, bxytmp, bxztmp, 
+     $                 byxtmp, byytmp, byztmp, 
+     $                 bzxtmp, bzytmp, bzztmp)
+
+                  if ( set_lapse ) then
+                     alp(i,j,k) = alptmp
                   end if
-                  exact_psix = 0.0D0
-                  exact_psiy = 0.0D0
-                  exact_psiz = 0.0D0
-                  exact_psixx = 0.0D0
-                  exact_psixx = 0.0D0
-                  exact_psizz = 0.0D0
-                  exact_psixy = 0.0D0
-                  exact_psiyz = 0.0D0
-                  exact_psixz = 0.0D0
-
-                  call Exact__Bona_Masso_data(
-     $                 decoded_exact_model,
-     $                 xx, yy, zz, tt,
-     $                 gxxjunk, gyyjunk, gzzjunk, 
-     $                 gxyjunk, gyzjunk, gxzjunk,
-     $                 hxxjunk, hyyjunk, hzzjunk, 
-     $                 hxyjunk, hyzjunk, hxzjunk,
-     $                 exact_psi,
-     $                 exact_psix, exact_psiy, exact_psiz,
-     $                 exact_psixx, exact_psiyy, exact_psizz,
-     $                 exact_psixy, exact_psiyz, exact_psixz,
-     $                 dxgxxjunk, dxgyyjunk, dxgzzjunk, 
-     $                 dxgxyjunk, dxgyzjunk, dxgxzjunk,
-     $                 dygxxjunk, dygyyjunk, dygzzjunk, 
-     $                 dygxyjunk, dygyzjunk, dygxzjunk,
-     $                 dzgxxjunk, dzgyyjunk, dzgzzjunk, 
-     $                 dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $                 alp(i,j,k), axjunk, ayjunk, azjunk, 
-     $                 betaxjunk, betayjunk, betazjunk,
-     $                 bxxjunk, bxyjunk, bxzjunk, 
-     $                 byxjunk, byyjunk, byzjunk, 
-     $                 bzxjunk, bzyjunk, bzzjunk)
-
-               end do
-            end do
-         end do
-
-CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
-
-C
-C     Set shift only?
-C
-      elseif (
-     $         ( is_initial_slice
-     $           .and. (CCTK_Equals(initial_shift,"exact").ne.0) )
-     $         .or.
-     $         ( is_later_slice
-     $           .and. (CCTK_Equals(shift_evolution_method,"exact").ne.0) )
-     $       ) then
-
-         do k=1,nz
-            do j=1,ny
-               do i=1,nx
-
-                  tt = cctk_time
-                  xx = x(i,j,k) - cctk_time * shift_add_x
-                  yy = y(i,j,k) - cctk_time * shift_add_y
-                  zz = z(i,j,k) - cctk_time * shift_add_z
-
-C     Initialize the psi of exact
-C     (also to tell the models about the conformal_state)
-                  if (conformal_state .ne. 0) then
-                    exact_psi = 1.0D0
-                  else
-                    exact_psi = 0.0D0
+                  if ( set_shift ) then
+                     betax(i,j,k) = betaxtmp + shift_add_x
+                     betay(i,j,k) = betaytmp + shift_add_y
+                     betaz(i,j,k) = betaztmp + shift_add_z
+                  end if
+                  if ( set_dtlapse ) then
+                     dtalp(i,j,k) = dtalptmp
+                  end if
+                  if ( set_dtshift ) then
+                     dtbetax(i,j,k) = dtbetaxtmp
+                     dtbetay(i,j,k) = dtbetaytmp
+                     dtbetaz(i,j,k) = dtbetaztmp
                   end if
-                  exact_psix = 0.0D0
-                  exact_psiy = 0.0D0
-                  exact_psiz = 0.0D0
-                  exact_psixx = 0.0D0
-                  exact_psixx = 0.0D0
-                  exact_psizz = 0.0D0
-                  exact_psixy = 0.0D0
-                  exact_psiyz = 0.0D0
-                  exact_psixz = 0.0D0
-
-                  call Exact__Bona_Masso_data(
-     $                 decoded_exact_model,
-     $                 xx, yy, zz, tt,
-     $                 gxxjunk, gyyjunk, gzzjunk, 
-     $                 gxyjunk, gyzjunk, gxzjunk,
-     $                 hxxjunk, hyyjunk, hzzjunk, 
-     $                 hxyjunk, hyzjunk, hxzjunk,
-     $                 exact_psi,
-     $                 exact_psix, exact_psiy, exact_psiz,
-     $                 exact_psixx, exact_psiyy, exact_psizz,
-     $                 exact_psixy, exact_psiyz, exact_psixz,
-     $                 dxgxxjunk, dxgyyjunk, dxgzzjunk, 
-     $                 dxgxyjunk, dxgyzjunk, dxgxzjunk,
-     $                 dygxxjunk, dygyyjunk, dygzzjunk, 
-     $                 dygxyjunk, dygyzjunk, dygxzjunk,
-     $                 dzgxxjunk, dzgyyjunk, dzgzzjunk, 
-     $                 dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $                 alpjunk, axjunk, ayjunk, azjunk, 
-     $                 betax(i,j,k), betay(i,j,k), betaz(i,j,k),
-     $                 bxxjunk, bxyjunk, bxzjunk, 
-     $                 byxjunk, byyjunk, byzjunk, 
-     $                 bzxjunk, bzyjunk, bzzjunk)
-                  betax(i,j,k) = betax(i,j,k) + shift_add_x
-                  betay(i,j,k) = betay(i,j,k) + shift_add_y
-                  betaz(i,j,k) = betaz(i,j,k) + shift_add_z
-
                end do
             end do
          end do
diff --git a/src/initialize.F77 b/src/initialize.F77
index 94c728a..c1ef459 100644
--- a/src/initialize.F77
+++ b/src/initialize.F77
@@ -20,8 +20,9 @@ C $Header$
       integer nx,ny,nz
 
       CCTK_REAL tt, xx, yy, zz
-      CCTK_REAL alpjunk, axjunk, ayjunk, azjunk, 
+      CCTK_REAL alpjunk, dtalpjunk, axjunk, ayjunk, azjunk, 
      $     betaxjunk, betayjunk, betazjunk,
+     $     dtbetaxjunk, dtbetayjunk, dtbetazjunk,
      $     bxxjunk, bxyjunk, bxzjunk, 
      $     byxjunk, byyjunk, byzjunk, 
      $     bzxjunk, bzyjunk, bzzjunk
@@ -59,6 +60,24 @@ C     Set data pointwise.
       ny = cctk_lsh(2)
       nz = cctk_lsh(3)
 
+C$omp parallel do private(
+C$omp$     tt, xx, yy, zz,
+C$omp$     alpjunk, dtalpjunk, axjunk, ayjunk, azjunk, 
+C$omp$     betaxjunk, betayjunk, betazjunk,
+C$omp$     dtbetaxjunk, dtbetayjunk, dtbetazjunk,
+C$omp$     bxxjunk, bxyjunk, bxzjunk, 
+C$omp$     byxjunk, byyjunk, byzjunk, 
+C$omp$     bzxjunk, bzyjunk, bzzjunk,
+C$omp$     dxgxxjunk, dxgyyjunk, dxgzzjunk, 
+C$omp$     dxgxyjunk, dxgyzjunk, dxgxzjunk,
+C$omp$     dygxxjunk, dygyyjunk, dygzzjunk, 
+C$omp$     dygxyjunk, dygyzjunk, dygxzjunk,
+C$omp$     dzgxxjunk, dzgyyjunk, dzgzzjunk, 
+C$omp$     dzgxyjunk, dzgyzjunk, dzgxzjunk,
+C$omp$     exact_psi,
+C$omp$     exact_psix, exact_psiy, exact_psiz,
+C$omp$     exact_psixx, exact_psiyy, exact_psizz,
+C$omp$     exact_psixy, exact_psiyz, exact_psixz)
       do k=1,nz
          do j=1,ny
             do i=1,nx
@@ -102,8 +121,9 @@ C     (also to tell the models about the conformal_state)
      $              dygxyjunk, dygyzjunk, dygxzjunk,
      $              dzgxxjunk, dzgyyjunk, dzgzzjunk, 
      $              dzgxyjunk, dzgyzjunk, dzgxzjunk,
-     $              alpjunk, axjunk, ayjunk, azjunk, 
+     $              alpjunk, dtalpjunk, axjunk, ayjunk, azjunk, 
      $              betaxjunk, betayjunk, betazjunk,
+     $              dtbetaxjunk, dtbetayjunk, dtbetazjunk,
      $              bxxjunk, bxyjunk, bxzjunk,
      $              byxjunk, byyjunk, byzjunk,
      $              bzxjunk, bzyjunk, bzzjunk)
diff --git a/src/xyz_blended_boundary.F77 b/src/xyz_blended_boundary.F77
index d8357bb..8106361 100644
--- a/src/xyz_blended_boundary.F77
+++ b/src/xyz_blended_boundary.F77
@@ -29,8 +29,8 @@ C $Header$
       CCTK_REAL dxgxxe, dxgyye, dxgzze, dxgxye, dxgyze, dxgxze
       CCTK_REAL dygxxe, dygyye, dygzze, dygxye, dygyze, dygxze
       CCTK_REAL dzgxxe, dzgyye, dzgzze, dzgxye, dzgyze, dzgxze
-      CCTK_REAL alpe, axe, aye, aze
-      CCTK_REAL betaxe,betaye,betaze
+      CCTK_REAL alpe, dtalpe, axe, aye, aze
+      CCTK_REAL betaxe,betaye,betaze, dtbetaxe,dtbetaye,dtbetaze
       CCTK_REAL bxxe,bxye,bxze,byxe,byye,byze,bzxe,bzye,bzze
       CCTK_REAL det, uxx, uxy, uxz, uyy, uyz, uzz
       CCTK_REAL
@@ -120,7 +120,8 @@ C     (also to tell the models about the conformal_state)
      $                 dxgxxe, dxgyye, dxgzze, dxgxye, dxgyze, dxgxze,
      $                 dygxxe, dygyye, dygzze, dygxye, dygyze, dygxze,
      $                 dzgxxe, dzgyye, dzgzze, dzgxye, dzgyze, dzgxze,
-     $                 alpe, axe, aye, aze, betaxe, betaye, betaze,
+     $                 alpe, dtalpe, axe, aye, aze,
+     $                 betaxe, betaye, betaze, dtbetaxe, dtbetaye, dtbetaze,
      $                 bxxe, bxye, bxze, byxe, 
      $                 byye, byze, bzxe, bzye, bzze)