diff --git a/_setup/outputOptions.py b/_setup/outputOptions.py
index d2556f637b81fceb3699e59a8578c1a222deaeb4..f3ced67c1f25b36cc017da06606f15c9069008f2 100644
--- a/_setup/outputOptions.py
+++ b/_setup/outputOptions.py
@@ -38,9 +38,9 @@ class OutputOptions:
         + "<m_y>,".ljust(self.columnWidth) + " " \
         + "<m_z>,".ljust(self.columnWidth)
       self.headerExtField = "columns are: t,".ljust(self.columnWidth-2) + " " \
-        + "H_x(0,0,0),".ljust(self.columnWidth) + " " \
-        + "H_y(0,0,0),".ljust(self.columnWidth) + " " \
-        + "H_z(0,0,0),".ljust(self.columnWidth)
+        + "H_x,".ljust(self.columnWidth) + " " \
+        + "H_y,".ljust(self.columnWidth) + " " \
+        + "H_z,".ljust(self.columnWidth)
       self.headerCurrent = "columns are: t,".ljust(self.columnWidth-2) + " " \
         + "Je(sl),".ljust(self.columnWidth) + " " \
         + "Je(zl)_x,".ljust(self.columnWidth) + " " \
@@ -99,7 +99,6 @@ class OutputOptions:
 
       self.traceUnprojectedExchangeEnergy = False
 
-      self._externalFieldEvaluationPoint = (0.0, 0.0, 0.0)
     else:
       for key, value in other.__dict__.items():
         self.__setattr__(key, value)
@@ -119,24 +118,6 @@ class OutputOptions:
 #    return string
 
 
-#------------------------------------------------------------------------------#
-
-
-  def SetExternalFieldEvaluationPoint(self, value):
-    self._externalFieldEvaluationPoint = value
-    self.headerExtField = "columns are: t,".ljust(self.columnWidth) \
-      + "H_x"+str(value)+",".ljust(self.columnWidth) \
-      + "H_y"+str(value)+",".ljust(self.columnWidth) \
-      + "H_z"+str(value)+",".ljust(self.columnWidth)
-
-
-#------------------------------------------------------------------------------#
-
-
-  def GetExternalFieldEvaluationPoint(self):
-    return self._externalFieldEvaluationPoint
-
-
 #------------------------------------------------------------------------------#
 
 
diff --git a/_setup/parameters.py b/_setup/parameters.py
index bfe21f426c612793189c60c31f2e6fd1b41aadd0..2e405199a87e29e930f3cbd191658a6bcd12af99 100644
--- a/_setup/parameters.py
+++ b/_setup/parameters.py
@@ -39,10 +39,11 @@ class Parameters:
       (magnetostaticField = strayfield & oersted field)
       .dmCoupling = ['none', 'bulk', 'interfacial'] (if interfacial, .dmi_axis (default = 'z'))
       .spintronicsCoupling = ['none', 'zhangLi', 'slonczewski']
+      .externalFieldRecordingPoint = (x1, x2, x3) in Omega, default "first mesh point"
   '''
   
 
-  def __init__(self, A=None, Ms=None, K=0.0, D=0.0, dmi_axis='z', gamma0=2.212761569e+05, mu0=1.2566370614359173e-06, alpha=1.0, T_start=None, T_end=None, timeStepSize=None, H_ext=None, easy_axis=(1,0,0), g=2.0023193043617, muB=9.27400968e-24, e=1.6021766208e-19, P=None, xi=None, dmCoupling="none", dmImplicitTreatment=False, maxwellCoupling="strayField", spintronicsCoupling="none", slonczewskiOutOfPlane_axis='z', hbar=1.054571800e-34, Je=None, flag_Je_stationary=False, d=None, p=(1,0,0), TOL=1.0e-8, preconditionerType="jacobi", alphaPreconditioner=None, QmatType="householder", QmatAxis="e3", guessType="old", theta=0.5, m0=None):
+  def __init__(self, A=None, Ms=None, K=0.0, D=0.0, dmi_axis='z', gamma0=2.212761569e+05, mu0=1.2566370614359173e-06, alpha=1.0, T_start=None, T_end=None, timeStepSize=None, H_ext=None, easy_axis=(1,0,0), g=2.0023193043617, muB=9.27400968e-24, e=1.6021766208e-19, P=None, xi=None, dmCoupling="none", dmImplicitTreatment=False, maxwellCoupling="strayField", spintronicsCoupling="none", slonczewskiOutOfPlane_axis='z', hbar=1.054571800e-34, Je=None, flag_Je_stationary=False, d=None, p=(1,0,0), TOL=1.0e-8, preconditionerType="jacobi", alphaPreconditioner=None, QmatType="householder", QmatAxis="e3", guessType="old", theta=0.5, m0=None, externalFieldRecordingPoint="first mesh point"):
     self.A = A
     self.Ms = Ms
     self.K = K
@@ -94,6 +95,7 @@ class Parameters:
 
     self.m0 = m0
 
+    self.externalFieldRecordingPoint = externalFieldRecordingPoint
 
 #------------------------------------------------------------------------------#
 
diff --git a/integrators/_details/interfaces/communicationInterface.py b/integrators/_details/interfaces/communicationInterface.py
index f934c207793c6bb01bd249e12f2e05dd102d2813..32dad90c6689953bb52a31c7e2c10b9cfe1a9b82 100644
--- a/integrators/_details/interfaces/communicationInterface.py
+++ b/integrators/_details/interfaces/communicationInterface.py
@@ -62,7 +62,7 @@ class CommunicationInterface(_interface._Interface):
 
   def _EvaluateExternalFieldAtPoint(self):
     import types
-    from ngsolve import x, y, z
+    from ngsolve import x, y, z, CoefficientFunction
     if self.parameters.H_ext == None:
       return [0.0, 0.0, 0.0]
 
@@ -71,12 +71,9 @@ class CommunicationInterface(_interface._Interface):
     for j in range(3):
       if isinstance(self.parameters.H_ext[j], types.FunctionType):
         time = self._time_cf / self.parameters.Ms / self.parameters.gamma0
-        cf = self.parameters.H_ext[j](time, x, y, z)
-        pnt = self.output.GetExternalFieldEvaluationPoint()
-        mip = self._geometry.GetMesh()( \
-          pnt[0] / self._geometry.scaling, \
-          pnt[1] / self._geometry.scaling, \
-          pnt[2] / self._geometry.scaling)
+        cf = CoefficientFunction(self.parameters.H_ext[j](time, x, y, z))
+        pnt = self._externalFieldRecordingPoint
+        mip = self._geometry.GetMesh()(pnt[0], pnt[1], pnt[2])
         ext[j] = cf(mip)
       else:
         ext[j] = self.parameters.H_ext[j]
diff --git a/integrators/_details/interfaces/computeInterface.py b/integrators/_details/interfaces/computeInterface.py
index 542dedebe75908177c16ed28c05165344d758a0f..11b27a199b5302591abed00fa08bace9e336d015 100644
--- a/integrators/_details/interfaces/computeInterface.py
+++ b/integrators/_details/interfaces/computeInterface.py
@@ -233,7 +233,7 @@ class ComputeInterface(_interface._Interface):
     for step in range(self._steps):
       self.Step()
     
-    print("Results have been written to: /", self.output.foldername)
+    print("Results have been written to:", self.output.foldername)
   
   
 #------------------------------------------------------------------------------#
diff --git a/integrators/_details/interfaces/parameterInterface.py b/integrators/_details/interfaces/parameterInterface.py
index 5f8f8833e656c8570bedc26937889b1b07de3b96..1638b05025df3c1d9df2d4e0b94ab07250044ed4 100644
--- a/integrators/_details/interfaces/parameterInterface.py
+++ b/integrators/_details/interfaces/parameterInterface.py
@@ -220,14 +220,14 @@ class ParameterInterface(_interface._Interface):
 
 
   def SetTOL(self, value):
-      self.parameters.TOL = value
+    self.parameters.TOL = value
 
 
 #------------------------------------------------------------------------------#
 
 
   def RecordExternalFieldAt(self, value):
-    self.output.SetExternalFieldEvaluationPoint(value)
+    self.parameters.externalFieldRecordingPoint = value
 
 
 #------------------------------------------------------------------------------#
diff --git a/integrators/_integrator.py b/integrators/_integrator.py
index 0d6676c1d090bad32ee18813e9d1223b063c2b02..6b1dd0e383f5402debcb56dd96d515f4d5f24ab2 100644
--- a/integrators/_integrator.py
+++ b/integrators/_integrator.py
@@ -94,6 +94,8 @@ class _Integrator( \
 
     self._maxwellCouplingChanged = False
 
+    self._externalFieldRecordingPoint = None
+
 
 #------------------------------------------------------------------------------#
 
@@ -153,6 +155,16 @@ class _Integrator( \
 
 
   def _SetUpOutput(self):
+    pnt = self.parameters.externalFieldRecordingPoint
+    if pnt == "first mesh point":
+      self._externalFieldRecordingPoint = \
+        list(self._geometry.GetMesh().ngmesh.Points()[1])
+    else:
+      self._externalFieldRecordingPoint = ( \
+        pnt[0] / self._geometry.scaling, \
+        pnt[1] / self._geometry.scaling, \
+        pnt[2] / self._geometry.scaling)
+
     if(self._saveMagEvery == None): return
     if(self._saveMagEvery == 0): return
 
diff --git a/integrators/integrator.py b/integrators/integrator.py
index dab3d59a508ea662a7818306ee72b9badf5e119f..643601d18e2569540ea57b08798ce131b1b88bfb 100644
--- a/integrators/integrator.py
+++ b/integrators/integrator.py
@@ -293,6 +293,7 @@ class Integrator:
     '''
     value [(float, float, float)]
     - external field will be recorded at given point in R^3
+    - make sure that point is part of the physical domain
     '''
     self._integrator.RecordExternalFieldAt(value)