1 files changed, 46 insertions, 0 deletions
diff --git a/src/derivs.F90 b/src/derivs.F90
new file mode 100644
index 0000000..1d425c8
--- /dev/null
+++ b/src/derivs.F90
@@ -0,0 +1,46 @@
+! $Header$
+
+#ifndef TGR_INCLUDED
+
+#include "cctk.h"
+#include "cctk_Parameters.h"
+
+module derivs
+  implicit none
+  private
+  public get_derivs
+  public get_derivs2
+contains
+#endif
+  subroutine get_derivs (a, f, pos, off, dx)
+    CCTK_REAL, intent(in)  :: a(*)
+    CCTK_REAL, intent(out) :: f(3)
+    integer,   intent(in)  :: pos, off(3)
+    CCTK_REAL, intent(in)  :: dx(3)
+    integer :: i
+    do i=1,3
+       f(i) = (a(pos+off(i)) - a(pos-off(i))) / (2*dx(i))
+    end do
+  end subroutine get_derivs
+  subroutine get_derivs2 (a, f, pos, off, dx)
+    CCTK_REAL, intent(in)  :: a(*)
+    CCTK_REAL, intent(out) :: f(3,3)
+    integer,   intent(in)  :: pos, off(3)
+    CCTK_REAL, intent(in)  :: dx(3)
+    integer :: i
+    do i=1,3
+       f(i,i) = (a(pos+off(i)) - 2*a(pos) + a(pos-off(i))) / dx(i)**2
+    end do
+    f(1,2) = (  a(pos-off(1)-off(2)) - a(pos+off(1)-off(2)) &
+         &    - a(pos-off(1)+off(2)) + a(pos+off(1)+off(2))) / (4*dx(1)*dx(2))
+    f(2,1) = f(1,2)
+    f(1,3) = (  a(pos-off(1)-off(3)) - a(pos+off(1)-off(3)) &
+         &    - a(pos-off(1)+off(3)) + a(pos+off(1)+off(3))) / (4*dx(1)*dx(3))
+    f(3,1) = f(1,3)
+    f(2,3) = (  a(pos-off(2)-off(3)) - a(pos+off(2)-off(3)) &
+         &    - a(pos-off(2)+off(3)) + a(pos+off(2)+off(3))) / (4*dx(2)*dx(3))
+    f(3,2) = f(2,3)
+  end subroutine get_derivs2
+#ifndef TGR_INCLUDED
+end module derivs
+#endif