<br><font size=2 face="sans-serif">Hi Berk,</font>
<br>
<br><font size=2 face="sans-serif">I don't really understand what you are
getting at with your &quot;as well&quot; statement below.</font>
<br>
<br><font size=2 face="sans-serif">With cubic splines you simply can't
have </font>
<br><font size=2 face="sans-serif">&nbsp; - the function values</font>
<br><font size=2 face="sans-serif">&nbsp; - the first order derivates</font>
<br><font size=2 face="sans-serif">&nbsp; - and the second order derivatives</font>
<br><font size=2 face="sans-serif">right with the same set of coefficients.
if you want all at the same time with the same set of coefficients,</font>
<br><font size=2 face="sans-serif">you need quintic splines, which are
more expensive (it wouldn't be a catastrophy for the</font>
<br><font size=2 face="sans-serif">computation time, but would mean a tremendous
amount of code hacking).</font>
<br>
<br><font size=2 face="sans-serif">True, getting the 2nd order derivates
right minimizes the error in the potential calculation to o(h**3), h=x_i-x_i+1,</font>
<br><font size=2 face="sans-serif">but the point is, it does *not* minimize
the error in the forces at the same time.</font>
<br>
<br><font size=2 face="sans-serif">So it essentially boils down two these
four choices:</font>
<br>
<br><font size=2 face="sans-serif">&nbsp;1. get better force accuracy and
keep the current accuracy for the potential, but you have to pay for it
with quintic splines</font>
<br><font size=2 face="sans-serif">&nbsp;2. have either accurate potentials
and inaccurate forces (the way it is right now)</font>
<br><font size=2 face="sans-serif">&nbsp;3. improve force accuracy at the
expense of some potential accuracy (what I proposed)</font>
<br><font size=2 face="sans-serif">&nbsp;4. approximate potentials and
forces using separate cubic spline polynomials</font>
<br>
<br><font size=2 face="sans-serif">I don't consider the 4th choice a good
one, but I listed it for completeness anyway.</font>
<br><font size=2 face="sans-serif">It means you would have two different
sets of coefficient tables, one for the potentials and</font>
<br><font size=2 face="sans-serif">one for the forces, which, in turn,
would lead to a change in the inner loops of the force kernels.</font>
<br><font size=2 face="sans-serif">Such an approach would also introduce
a small systematical bias, since the simulated forces</font>
<br><font size=2 face="sans-serif">are no longer exact derivatives of the
measured potential function, which I, personally, would</font>
<br><font size=2 face="sans-serif">stay away from.</font>
<br>
<br><font size=2 face="sans-serif">So overall I would go for the 3rd choice,
because</font>
<br>
<br><font size=2 face="sans-serif">&nbsp;- I get better energy conservation
during the integration</font>
<br><font size=2 face="sans-serif">&nbsp;- I keep force and potential consistent
both in absolute values and in accuracy (both to o(h**2))</font>
<br>
<br><font size=2 face="sans-serif">the only thing I loose is some accuracy
in the potential, which</font>
<br><font size=2 face="sans-serif">&nbsp;- is approximate anyway</font>
<br><font size=2 face="sans-serif">&nbsp;- can always be improved by increasing
the number of table sites</font>
<br>
<br><font size=2 face="sans-serif">I consider the first two arguments more
important, but then again, I am not certain, that this really is the best
way to go,</font>
<br><font size=2 face="sans-serif">since I don't really have strong experience
using tabulated potentials and have not performed any studies on this.</font>
<br>
<br><tt><font size=3>Mathias PUETZ wrote:<br>
<br>
&gt;<i><br>
</i>&gt;<i> Hi,<br>
</i>&gt;<i><br>
</i>&gt;<i> the cubic spline function<br>
</i>&gt;<i><br>
</i>&gt;<i> spline_i (x) = a+b(x - x_i) + c (x - x_i)**2 + d (x - x_i)**3<br>
</i>&gt;<i><br>
</i>&gt;<i> &nbsp;has four variable coefficients: a,b,c,d.<br>
</i>&gt;<i><br>
</i>&gt;<i> &nbsp;For each table segment i you have four linear equations
which <br>
</i>&gt;<i> determine the coefficients for the spline function spline_i(x)<br>
</i>&gt;<i><br>
</i>&gt;<i> v(x_i) = spline_i(x_i)<br>
</i>&gt;<i> v(x_i+1) = spline_i(x_i+1)<br>
</i>&gt;<i> v''(x_i) = spline_i''(x_i)<br>
</i>&gt;<i> v''(x_i+1) = spline_i''(x_i+1)<br>
</i>&gt;<i><br>
</i>&gt;<i> solving these equations for a,b,c,d gives exactly the spline
formula <br>
</i>&gt;<i> you are using.<br>
</i>&gt;<i><br>
</i>&gt;<i> Since you force v and v'' to be continous, v' (the forces)
cannot be <br>
</i>&gt;<i> continuous at x_i<br>
</i>&gt;<i> unless v(x) itself is a polynomial function with a degree of
less than <br>
</i>&gt;<i> 4, which it typically isn't.<br>
</i>&gt;<i> The error will be small in most cases, but it is noticable
if you <br>
</i>&gt;<i> print out the forces<br>
</i>&gt;<i> calculated by your spline formula and the exact force function
in the <br>
</i>&gt;<i> vicinity of the x_i<br>
</i>&gt;<i><br>
</i>&gt;<i> I believe it would be more appropriate to enforce<br>
</i>&gt;<i><br>
</i>&gt;<i> v(x_i) = spline_i(x_i)<br>
</i>&gt;<i> v(x_i+1) = spline_i(x_i+1)<br>
</i>&gt;<i> v'(x_i) = spline_i'(x_i)<br>
</i>&gt;<i> v'(x_i+1) = spline_i'(x_i+1)<br>
</i>&gt;<i><br>
</i>&gt;<i> which will give you slightly different spline coefficients,
but <br>
</i>&gt;<i> guarantees that both v(x) and v'(x)<br>
</i>&gt;<i> will be continuous at x_i.<br>
</i>&gt;<i><br>
</i>You are right.<br>
<br>
But I think we would still prefer to have continuous second derivatives
<br>
as well.<br>
So the problem is then not in the interloops, but in the generation of
<br>
the tables.<br>
There we should not use the second derivative of the potential,<br>
but rather the second derivative of the cubic spline fit of the potential.<br>
<br>
Berk.</font></tt><font size=2 face="sans-serif"><br>
Viele Grüsse / Best regards,<br>
Dr. Mathias Pütz<br>
<br>
IT Specialist for Application Perfomance<br>
<br>
Deep Computing - Strategic Growth BusinessDeep <br>
IBM Systems &amp; Technology Group<br>
<br>
e-mail: &nbsp;mpuetz@de.ibm.com<br>
mobile: + 49-(0)160-7120602<br>
fax: &nbsp; &nbsp; &nbsp; &nbsp; + 49-(0)6131-84-6660<br>
<br>
snailmail:<br>
&nbsp;IBM Deutschland GmbH<br>
&nbsp;Department B458<br>
&nbsp;Hechtsheimer Str. 2 / Building 12<br>
&nbsp;55131 Mainz<br>
&nbsp;Germany<br>
</font>