From 06903d0f2e0b9ebe34efbeec1dda3aaaef6dd26d Mon Sep 17 00:00:00 2001
From: schnetter <schnetter@f2ea251b-07d6-4a45-8f2c-a162b3fa7596>
Date: Tue, 24 Nov 2009 21:52:23 +0000
Subject: Add three thorns ID_Bin_BH, ID_Bin_NS, and ID_Mag_NS that import
 LORENE initial data.

git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinInitialData/Meudon_Bin_NS/trunk@2 f2ea251b-07d6-4a45-8f2c-a162b3fa7596
---
 src/Bin_NS.cc           | 158 ++++++++++++++++++++++++++++++++++++++++++++++++
 src/check_parameters.cc |  20 ++++++
 src/make.code.defn      |   7 +++
 3 files changed, 185 insertions(+)
 create mode 100644 src/Bin_NS.cc
 create mode 100644 src/check_parameters.cc
 create mode 100644 src/make.code.defn

(limited to 'src')

diff --git a/src/Bin_NS.cc b/src/Bin_NS.cc
new file mode 100644
index 0000000..4e2696f
--- /dev/null
+++ b/src/Bin_NS.cc
@@ -0,0 +1,158 @@
+#include <cassert>
+#include <cstdio>
+#include <vector>
+
+#include <cctk.h>
+#include <cctk_Arguments.h>
+#include <cctk_Parameters.h>
+
+#include <bin_ns.h>
+
+using namespace std;
+
+
+
+extern "C"
+void ID_Bin_NS_initialise (CCTK_ARGUMENTS)
+{
+  DECLARE_CCTK_ARGUMENTS;
+  DECLARE_CCTK_PARAMETERS;
+  
+  CCTK_INFO ("Setting up LORENE Bin_NS initial data");
+  
+  // Meudon data are distributed in SI units (MKSA).  Here are some
+  // conversion factors.
+  
+  // Defined constants
+  CCTK_REAL const c_light = 299792458.0; // speed of light [m/s]
+  
+  // Constants of nature (IAU, CODATA):
+  CCTK_REAL const G_grav = 6.67428e-11; // gravitational constant [m^3/kg/s^2]
+  CCTK_REAL const M_sun  = 1.98892e+30; // solar mass [kg]
+  
+  // Cactus units in terms of SI units:
+  // (These are derived from M = M_sun, c = G = 1, and using 1/M_sun
+  // for the magnetic field)
+  CCTK_REAL const cactusM = M_sun;
+  CCTK_REAL const cactusL = cactusM * G_grav / pow(c_light,2);
+  CCTK_REAL const cactusT = cactusL / c_light;
+  
+  // Other quantities in terms of Cactus units
+  CCTK_REAL const coord_unit = cactusL / 1.0e+3;         // from km
+  CCTK_REAL const rho_unit   = cactusM / pow(cactusL,3); // from kg/m^3
+  CCTK_REAL const ener_unit  = pow(cactusL,2);           // from c^2
+  CCTK_REAL const vel_unit   = cactusL / cactusT / c_light; // from c
+  
+  
+  
+  CCTK_INFO ("Setting up coordinates");
+  
+  int const npoints = cctk_lsh[0] * cctk_lsh[1] * cctk_lsh[2];
+  
+  vector<double> xx(npoints), yy(npoints), zz(npoints);
+  
+#pragma omp parallel for
+  for (int i=0; i<npoints; ++i) {
+    xx[i] = x[i] * coord_unit;
+    yy[i] = y[i] * coord_unit;
+    zz[i] = z[i] * coord_unit;
+  }
+  
+  
+  
+  CCTK_VInfo (CCTK_THORNSTRING, "Reading from file \"%s\"", filename);
+  
+  Bin_NS bin_ns (npoints, &xx[0], &yy[0], &zz[0], filename);
+  
+  CCTK_VInfo (CCTK_THORNSTRING, "omega [rad/s]:       %g", bin_ns.omega);
+  CCTK_VInfo (CCTK_THORNSTRING, "dist [km]:           %g", bin_ns.dist);
+  CCTK_VInfo (CCTK_THORNSTRING, "dist_mass [km]:      %g", bin_ns.dist_mass);
+  CCTK_VInfo (CCTK_THORNSTRING, "mass1_b [M_sun]:     %g", bin_ns.mass1_b);
+  CCTK_VInfo (CCTK_THORNSTRING, "mass2_b [M_sun]:     %g", bin_ns.mass2_b);
+  CCTK_VInfo (CCTK_THORNSTRING, "mass_ADM [M_sun]:    %g", bin_ns.mass_adm);
+  CCTK_VInfo (CCTK_THORNSTRING, "L_tot [G M_sun^2/c]: %g", bin_ns.angu_mom);
+  CCTK_VInfo (CCTK_THORNSTRING, "rad1_x_comp [km]:    %g", bin_ns.rad1_x_comp);
+  CCTK_VInfo (CCTK_THORNSTRING, "rad1_y [km]:         %g", bin_ns.rad1_y);
+  CCTK_VInfo (CCTK_THORNSTRING, "rad1_z [km]:         %g", bin_ns.rad1_z);
+  CCTK_VInfo (CCTK_THORNSTRING, "rad1_x_opp [km]:     %g", bin_ns.rad1_x_opp);
+  CCTK_VInfo (CCTK_THORNSTRING, "rad2_x_comp [km]:    %g", bin_ns.rad2_x_comp);
+  CCTK_VInfo (CCTK_THORNSTRING, "rad2_y [km]:         %g", bin_ns.rad2_y);
+  CCTK_VInfo (CCTK_THORNSTRING, "rad2_z [km]:         %g", bin_ns.rad2_z);
+  CCTK_VInfo (CCTK_THORNSTRING, "rad2_x_opp [km]:     %g", bin_ns.rad2_x_opp);
+  assert (bin_ns.np == npoints);
+  
+  
+  
+  CCTK_INFO ("Filling in Cactus grid points");
+  
+#pragma omp parallel for
+  for (int i=0; i<npoints; ++i) {
+    
+    alp[i] = bin_ns.nnn[i];
+    
+    betax[i] = bin_ns.beta_x[i];
+    betay[i] = bin_ns.beta_y[i];
+    betaz[i] = bin_ns.beta_z[i];
+    
+    CCTK_REAL g[3][3];
+    g[0][0] = bin_ns.g_xx[i];
+    g[0][1] = bin_ns.g_xy[i];
+    g[0][2] = bin_ns.g_xz[i];
+    g[1][1] = bin_ns.g_yy[i];
+    g[1][2] = bin_ns.g_yz[i];
+    g[2][2] = bin_ns.g_zz[i];
+    g[1][0] = g[0][1];
+    g[2][0] = g[0][2];
+    g[2][1] = g[1][2];
+    
+    CCTK_REAL ku[3][3];
+    ku[0][0] = bin_ns.k_xx[i];
+    ku[0][1] = bin_ns.k_xy[i];
+    ku[0][2] = bin_ns.k_xz[i];
+    ku[1][1] = bin_ns.k_yy[i];
+    ku[1][2] = bin_ns.k_yz[i];
+    ku[2][2] = bin_ns.k_zz[i];
+    ku[1][0] = ku[0][1];
+    ku[2][0] = ku[0][2];
+    ku[2][1] = ku[1][2];
+    
+    CCTK_REAL k[3][3];
+    for (int a=0; a<3; ++a) {
+      for (int b=0; b<3; ++b) {
+        k[a][b] = 0.0;
+        for (int c=0; c<3; ++c) {
+          for (int d=0; d<3; ++d) {
+            k[a][b] += g[a][c] * g[b][d] * ku[c][d];
+          }
+        }
+      }
+    }
+    
+    gxx[i] = g[0][0];
+    gxy[i] = g[0][1];
+    gxz[i] = g[0][2];
+    gyy[i] = g[1][1];
+    gyz[i] = g[1][2];
+    gzz[i] = g[2][2];
+    
+    kxx[i] = k[0][0];
+    kxy[i] = k[0][1];
+    kxz[i] = k[0][2];
+    kyy[i] = k[1][1];
+    kyz[i] = k[1][2];
+    kzz[i] = k[2][2];
+    
+    rho[i] = bin_ns.nbar[i] / rho_unit;
+    
+    eps[i] = rho[i] * bin_ns.ener_spec[i] / ener_unit;
+    
+    vel[i          ] = bin_ns.u_euler_x[i] / vel_unit;
+    vel[i+  npoints] = bin_ns.u_euler_y[i] / vel_unit;
+    vel[i+2*npoints] = bin_ns.u_euler_z[i] / vel_unit;
+    
+  } // for i
+  
+  
+  
+  CCTK_INFO ("Done.");
+}
diff --git a/src/check_parameters.cc b/src/check_parameters.cc
new file mode 100644
index 0000000..bef3136
--- /dev/null
+++ b/src/check_parameters.cc
@@ -0,0 +1,20 @@
+#include <cctk.h>
+#include <cctk_Arguments.h>
+#include <cctk_Parameters.h>
+
+
+
+extern "C"
+void ID_Bin_NS_check_parameters (CCTK_ARGUMENTS)
+{
+  DECLARE_CCTK_ARGUMENTS;
+  DECLARE_CCTK_PARAMETERS;
+  
+  if (not CCTK_EQUALS (initial_data , "ID_Bin_NS") or
+      not CCTK_EQUALS (initial_lapse, "ID_Bin_NS") or
+      not CCTK_EQUALS (initial_shift, "ID_Bin_NS") or
+      not CCTK_EQUALS (initial_hydro, "ID_Bin_NS"))
+  {
+    CCTK_PARAMWARN ("The parameters ADMBase::initial_data, ADMBase::initial_lapse, ADMBase::initial_shift, and HydroBase::initial_hydro must all be set to the value \"ID_Bin_NS\"");
+  }
+}
diff --git a/src/make.code.defn b/src/make.code.defn
new file mode 100644
index 0000000..bc238db
--- /dev/null
+++ b/src/make.code.defn
@@ -0,0 +1,7 @@
+# Main make.code.defn file for thorn ID_Bin_NS
+
+# Source files in this directory
+SRCS = Bin_NS.cc check_parameters.cc
+
+# Subdirectories containing source files
+SUBDIRS = 
-- 
cgit v1.2.3