<div dir="ltr">I escaped a cluttered Fortran/C/Pascal/Basic/Perl jumble by "living" in stackoverflow. Mostly by reading questions and answers, and learning how not to get my questions closed in milliseconds. Break your program into small pieces, and look for how each piece might be done differently. <div><br></div><div>If you try to find ways to do things without ever manually checking the sizes of arrays, you'll naturally be guided towards "pythonic" methods and techniques. There is almost always a better way.</div><div><br></div><div>Computers are so fast these days, don't worry about ultimate speed until you really hit a brick wall. Letting python, numpy, or scipy methods handle housekeeping for you means you are usually defaulting to very well written compiled routines.</div><div><br></div><div>So for example, let python do the thinking...</div><div><br></div><div><div><font face="monospace, monospace">r = np.sqrt((diff ** 2).sum(2))</font></div><div><font face="monospace, monospace">r[r < 0.0001] = 0.0001 </font></div><div><font face="monospace, monospace"># or just</font></div><div><font face="monospace, monospace">r = np.maximum(np.sqrt((diff ** 2).sum(2)), 0.0001)</font></div><div><font face="monospace, monospace"><br></font></div><div><font face="monospace, monospace"># and</font></div><div><font face="monospace, monospace"><br></font></div><div><font face="monospace, monospace">Zij = Z * Z[:, None]</font></div><div><font face="monospace, monospace">triangle = np.triu(np.ones_like(Zij), k=1)</font></div><div><font face="monospace, monospace">I = Zij*trangle.sum()</font></div><div><font face="monospace, monospace"><br></font></div><div><font face="monospace, monospace"># or just</font></div><div><font face="monospace, monospace"><br></font></div><div><font face="monospace, monospace">I = np.triu(Z * Z[:, None], k=1).sum()</font></div></div><div><br></div></div><div class="gmail_extra"><br><div class="gmail_quote">On Sun, Apr 16, 2017 at 8:25 PM, Stephen P. Molnar <span dir="ltr"><<a href="mailto:s.molnar@sbcglobal.net" target="_blank">s.molnar@sbcglobal.net</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="HOEnZb"><div class="h5">On 04/15/2017 04:07 PM, Hjalmar Turesson wrote:<br>
</div></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div><div class="h5">
Hi Stephen,<br>
<br>
<br>
I see that N = 10 (ie. the length of Z), and that r is indexed by i and<br>
j, which run up to N. However, r is only a 5x5 array. Is this correct?<br>
<br>
Otherwise, something like this might work:<br>
<br>
import numpy as np<br>
<br>
start=1<br>
finish=31<br>
points=300<br>
s = np.linspace(start, finish, points)<br>
<br>
MASS = np.array([12.011, 1.008, 79.9, 35.453, 18.998])<br>
<br>
r = np.array([[0., 2.059801, 3.60937686, 3.32591826, 2.81569212],<br>
[2.059801, 0., 4.71452879, 4.45776445, 4.00467382],<br>
[3.60937686, 4.71452879, 0., 5.66500917, 5.26602175],<br>
[3.32591826, 4.45776445, 5.66500917, 0., 5.02324896],<br>
[2.81569212, 4.00467382, 5.26602175, 5.02324896, 0.]])<br>
<br>
r[r == 0.] = 0.0001<br>
<br>
Z = np.array([6.0, 210.0, 102.0, 54.0, 35.0, 17.0, 9.0, 595.0, 315.0,<br>
153.0])<br>
<br>
N = Z.size<br>
<br>
I = np.zeros(s.shape)<br>
<br>
for i in range(1, N):<br>
for j in range(j):<br>
I += Z[i] * Z[j] * np.sin(s * r[i, j])/(s * r[i, j])<br>
<br>
<br>
Do you have an example of the correct output?<br>
<br>
Best,<br>
Hjalmar<br>
<br>
On Sat, Apr 15, 2017 at 6:29 AM, Stephen P. Molnar<br></div></div><div><div class="h5">
<<a href="mailto:s.molnar@sbcglobal.net" target="_blank">s.molnar@sbcglobal.net</a> <mailto:<a href="mailto:s.molnar@sbcglobal.net" target="_blank">s.molnar@sbcglobal.net</a><wbr>>> wrote:<br>
<br>
I am a newcomer to Python and am attempting to translate a FORTRAN<br>
program that I wrote about 20 years ago into Python, specifically v-3.6.<br>
<br>
Let me say that I am not asking someone to write the program for me,<br>
but only to point me is the right direction.<br>
<br>
So far, I have bumbled my way to the point that I can get all of the<br>
input data resulting from a quantum mechanical calculation of a very<br>
simple organic molecule in to a Python program, but am encountering<br>
problems with processing the data.<br>
<br>
The equation that I a want to evaluate (the one that I programmed in<br>
FORTRAN) is equation (7) in the attached scan of one of the pages of<br>
the following document.<br>
<br>
Stephen P. Molnar and James W. King, Theory and Applications of the<br>
Integrated Molecular Transform and the Normalized Molecular Moment<br>
Structure Descriptors: QSAR and QSPR Paradigms, Int. J Quantum<br>
Chem., 85, 662 (2001).<br>
<br>
the equation is<br>
<br>
I(s) = sum from i=2toN sum from j=1toi-1 Z_sub_i*Z_sub_j *<br>
<br>
sin s*r_sub_ij<br>
----------------- Equation (7)<br>
s*r_sub_ij<br>
<br>
What I have managed to do so far is for a simple organic molecule<br>
containing 5 atoms (CHFClBr):<br>
<br>
import numpy as np<br>
import numpy, pandas, scipy.spatial.distance as dist<br>
<br>
r=[] #distances between pairs of atoms<br>
z = [] #atomic numbers of atoms<br>
Z_dist = [] #products of successive atomic numbers<br>
I = [] #vector calculated in equation (7)<br>
<br>
start=1<br>
finish=31<br>
points=300<br>
s = np.linspace(start,finish,point<wbr>s)<br>
<br>
<br>
name = input("Enter Molecule ID: ")<br>
name = str(name)<br>
name_in = name+'.dat'<br>
<br>
df = pandas.read_table(name_in, skiprows=2, sep=" ",<br>
skipinitialspace=True)<br>
z = numpy.array(df["ZA"])<br>
N = numpy.ma.size(z) ## of atoms in molecule<br>
a = numpy.array(df[["X", "Y", "Z"]])<br>
dist.squareform(dist.pdist(a, "euclidean"))<br>
anrows, ancols = np.shape(a)<br>
a_new = a.reshape(anrows, 1, ancols)<br>
diff = a_new - a<br>
<br>
r = (diff ** 2).sum(2)<br>
r = np.sqrt(r)<br>
<br>
for j in range(0,N-1):<br>
for k in range(j+1,N):<br>
Z_diff = (z[j]*z[k])<br>
<br>
<br>
For the test molecule I get the following:<br>
<br>
MASS = [ 12.011 1.008 79.9 35.453 18.998] (Atomic weights<br>
of atoms)<br>
<br>
r: [ 0. 2.059801 3.60937686 3.32591826 2.81569212]<br>
[ 2.059801 0. 4.71452879 4.45776445 4.00467382]<br>
[ 3.60937686 4.71452879 0. 5.66500917 5.26602175]<br>
[ 3.32591826 4.45776445 5.66500917 0. 5.02324896]<br>
[ 2.81569212 4.00467382 5.26602175 5.02324896 0<br></div></div>
<tel:02324896%20%200>. ]<div><div class="h5"><br>
<br>
z:<br>
<br>
6.0<br>
210.0<br>
102.0<br>
54.0<br>
35.0<br>
17.0<br>
9.0<br>
595.0<br>
315.0<br>
153.0<br>
<br>
I have checked these calculations with a spreadsheet calculation,<br>
consequently I'm correct as far as I have gotten.<br>
<br>
However, my problem is calculating the value of I for the 300<br>
x-coordinate values required for the molecular index calculation.<br>
<br>
The problem that I am currently having is calculating the product of<br>
the pre-sine function and the sine term for the 300 values of s.<br>
<br>
Unfortunately, the wealth of function in Python makes it difficult<br>
for me to ascertain just how to proceed. It seems that I have the<br>
problem of seeing the trees because of the size of the forest. At<br>
this point Goggling only makes my confusion deeper. Any pointers in<br>
the correct direction will be greatly appreciated.<br>
<br>
Thanks in advance.<br>
<br>
--<br>
Stephen P. Molnar, Ph.D. Life is a fuzzy set<br></div></div>
<a href="http://www.molecular-modeling.net" rel="noreferrer" target="_blank">www.molecular-modeling.net</a> <<a href="http://www.molecular-modeling.net" rel="noreferrer" target="_blank">http://www.molecular-modeling<wbr>.net</a>><span class=""><br>
Stochastic and multivariate<br>
<a href="tel:%28614%29312-7528" value="+16143127528" target="_blank">(614)312-7528</a> (c)<br>
Skype: smolnar1<br>
<br>
<br>
______________________________<wbr>_________________<br>
SciPy-User mailing list<br></span>
<a href="mailto:SciPy-User@python.org" target="_blank">SciPy-User@python.org</a> <mailto:<a href="mailto:SciPy-User@python.org" target="_blank">SciPy-User@python.org</a>><br>
<a href="https://mail.python.org/mailman/listinfo/scipy-user" rel="noreferrer" target="_blank">https://mail.python.org/mailma<wbr>n/listinfo/scipy-user</a><br>
<<a href="https://mail.python.org/mailman/listinfo/scipy-user" rel="noreferrer" target="_blank">https://mail.python.org/mailm<wbr>an/listinfo/scipy-user</a>><span class=""><br>
<br>
<br>
<br>
<br>
______________________________<wbr>_________________<br>
SciPy-User mailing list<br>
<a href="mailto:SciPy-User@python.org" target="_blank">SciPy-User@python.org</a><br>
<a href="https://mail.python.org/mailman/listinfo/scipy-user" rel="noreferrer" target="_blank">https://mail.python.org/mailma<wbr>n/listinfo/scipy-user</a><br>
<br>
</span></blockquote>
Hjalmar<br>
<br>
Thanks for your response.<br>
<br>
I'm sorry that my original request wasn't a bit clearer. I presume that the paper I attached to the request help was deleted before the email was posted to the list.<br>
<br>
Here is how I incorporated your suggestion in the code.<span class=""><br>
<br>
import numpy as np<br>
import numpy, pandas, scipy.spatial.distance as dist<br>
<br>
r = []<br></span>
s=[]<br>
z = []<br>
Z = []<br>
I = []<span class=""><br>
<br>
start=1<br>
finish=31<br>
points=300<br>
s = np.linspace(start,finish,point<wbr>s)<br></span>
np.savetxt('Number of Points',s,delimiter=' ')<span class=""><br>
<br>
name = input("Enter Molecule ID: ")<br>
name = str(name)<br>
name_in = name+'.dat'<br>
<br>
df = pandas.read_table(name_in, skiprows=2, sep=" ", skipinitialspace=True)<br></span>
Z = numpy.array(df["ZA"])<br>
print('z = ',z)<br>
N = numpy.ma.size(z)<span class=""><br>
a = numpy.array(df[["X", "Y", "Z"]])<br>
dist.squareform(dist.pdist(a, "euclidean"))<br>
anrows, ancols = np.shape(a)<br>
a_new = a.reshape(anrows, 1, ancols)<br>
diff = a_new - a<br>
<br>
r = (diff ** 2).sum(2)<br></span>
r = np.sqrt(z)<span class=""><br>
<br>
r[r == 0.] = 0.0001<br>
<br></span><span class="">
N = Z.size<br>
<br>
I = np.zeros(s.shape)<br>
<br>
for i in range(1, N):<br>
for j in range(j):<br>
I += Z[i] * Z[j] * np.sin(s * r[i, j])/(s * r[i, j])<br></span>
print('I: ',I)<br>
<br>
np.savetxt('I',I,delimiter=' ')<br>
<br>
The for j in range(j): has a marginal note (I'm using the Spyder IDE) "undefined name 'y' " and I is a vector with three hundred 0.0 terms.<br>
<br>
also, please note that I removed 'Z = np.array([6.0, 210.0, 102.0, 54.0, 35.0, 17.0, 9.0, 595.0, 315.0, 153.0])' from your suggestion as that is the result of the per-sine term 'Z[i] * Z[j]'<br>
<br>
Regards,<br>
<br>
Steve<div class="HOEnZb"><div class="h5"><br>
<br>
-- <br>
Stephen P. Molnar, Ph.D. Life is a fuzzy set<br>
<a href="http://www.molecular-modeling.net" rel="noreferrer" target="_blank">www.molecular-modeling.net</a> Stochastic and multivariate<br>
<a href="tel:%28614%29312-7528" value="+16143127528" target="_blank">(614)312-7528</a> (c)<br>
Skype: smolnar1<br>
______________________________<wbr>_________________<br>
SciPy-User mailing list<br>
<a href="mailto:SciPy-User@python.org" target="_blank">SciPy-User@python.org</a><br>
<a href="https://mail.python.org/mailman/listinfo/scipy-user" rel="noreferrer" target="_blank">https://mail.python.org/mailma<wbr>n/listinfo/scipy-user</a><br>
</div></div></blockquote></div><br></div>