ChemicalData["Aspirin","MoleculePlot"]

t=0.01889726*ChemicalData["Aspirin","AtomPositions"]

ChemicalData["Aspirin","BondTally"]

u=ChemicalData["Aspirin","VertexTypes"]

u1={0.18,0.18,0.18,0.18,0.26,0.26,0.26,0.26,0.26,0.26,0.26,0.26,0.26,0.04,0.04,0.04,0.04,0.04,0.04,0.04,0.04};

(*For[i=1,i<15,i++,If[u[[i]]="C",u1[[i]]=0.26,u1[[i]]=0.18]]
u1*)

ChemicalData["Aspirin","MolecularFormulaString"]

v=0.01889726*ChemicalData["Aspirin","VertexCoordinates"]

lp=ListPlot[v, Axes->False]

(*Plot1=Table[Graphics[{Circle[{v[[i,1]]/0.5292,v[[i,2]]/0.5292},If[i<5,0.18,0.26]]}],{i,1,13}];*)
Plot1=Table[Graphics[{Circle[{v[[i,1]],v[[i,2]]},u1[[i]]]}],{i,1,21}];
Show[Plot1,lp,{AspectRatio->Automatic,Axes -> False,Frame -> Automatic,PlotRange -> Automatic}]

ChemicalData["Aspirin","CHStructureDiagram"]

x=Table[t[[i,1]],{i,1,21}];
y=Table[t[[i,2]],{i,1,21}];
z=Table[t[[i,3]],{i,1,21}];

(* Distance matrix Bohr *)
Dij=Table[If[i≥j,If[sq=Sqrt[(x[[i]]-x[[j]])^2+(y[[i]]-y[[j]])^2+(z[[i]]-z[[j]])^2]; sq<3.0,sq,0]," "],{i,1,21},{j,1,21}];
(*Dij=Table[If[sq=Sqrt[(x[[i]]-x[[j]])^2+(y[[i]]-y[[j]])^2+(z[[i]]-z[[j]])^2]; sq<3.0,sq," "],{i,1,21},{j,1,21}];*)
MatrixForm[Dij]

(*Eij=Table[If[i≥j,If[sq=Sqrt[(x[[i]]-x[[j]])^2+(y[[i]]-y[[j]])^2+(z[[i]]-z[[j]])^2]; sq<3.0,0.529177*sq,0]," "],{i,1,21},{j,1,21}];
MatrixForm[Eij];*)

u2={0.18,0.18,0.18,0.18,0.26,0.26,0.26,0.26,0.26,0.26,0.26,0.26,0.26,0.56,0.56,0.56,0.56,0.56,0.56,0.56,0.50};

Fij=Dij-u2;
MatrixForm[Fij];
Gij=Table[If[i≥j,If[Fij[[i,j]]>0,Fij[[i,j]],0]," "],{i,1,21},{j,1,21}];
MatrixForm[Gij]

Hij=Transpose[Gij];
MatrixForm[Hij]

Kij=Hij-u2;
MatrixForm[Kij];
Lij=Table[If[i≤j,If[Kij[[i,j]]>0,Kij[[i,j]],0]," "],{i,1,21},{j,1,21}];
MatrixForm[Lij]

Lij=Table[If[i≤j,If[Lij[[i,j]]>0,If[j<14,0.5*Lij[[i,j]],Lij[[i,j]]],0]," "],{i,1,21},{j,1,21}]; Qij=Lij;
MatrixForm[Lij]

Nij=Table[If[i≤j,If[Lij[[i,j]]>0,2.25/Lij[[i,j]]^2,0]," "],{i,1,21},{j,1,21}];
MatrixForm[Nij]

Tsigma=Sum[If[i≤j,Nij[[i,j]],0],{i,1,21},{j,1,21}]+Sum[2.25/u2[[i]]^2,{i,1,13}]

Mij=Prepend[Qij,u];
Oij=Transpose[Mij];
Pij=Prepend[Oij,Prepend[u," "]];

MatrixForm[Pij]

Pij;

(* sigma bond cloud radii *)

For[i=2,i<23,i++,
For[j=2,j≤i,j++,
If[Pij[[i,j]]>0,Print[u[[j-1]],j-1,"-",u[[i-1]],i-1,"  ",Extract[Pij,{i,j}]]]]]