‌

Making a Molecule

At base, the molecule should just be a list of atoms with other properties, but let’s write a nicer constructor than requiring an explicit list of atoms:

SetChemicalSystemUpValues@Molecule;

Molecule[atoms:{_Atom...}:{},
 bonds:{_Bond...}:{}]:=ChemicalSystemAdd[If[Length@atoms>0,
  bonds/.{Bond[i_Integer,j_Integer,t___]Bond[atoms[[i]],atoms[[j]],t]};
  With[{s=atoms[[1]]["System"]},ChemicalSystem@s],
  $DefaultChemicalSystem],
 Molecule@<|"Atoms"atoms|>
  ];

We’ll also code in a constructor for molecules in zmatrix form

zmConvert[{el_String}]:=AppendTo[$zmAtomStore,{el,{0,0,0} }];
zmConvert[{el_String,ref_Integer,dist_}]:=AppendTo[$zmAtomStore,{el,Last@$zmAtomStore[[ref]]+{dist,0,0} }];
zmConvert[{el_String,ref1_Integer,dist_,ref2_Integer,ang_}]:=AppendTo[$zmAtomStore,
 {el,With[{p1=Last@$zmAtomStore[[ref1]],p2=Last@$zmAtomStore[[ref2]]},
   p1+dist*Normalize[RotationMatrix[ang,
    Replace[(p2-p1){1,0,0},
     Except[_List]|{0,0,0}->{0,0,1}]].(p2-p1)
    ]
   ]}];
 zmConvert[{el_String,ref1_Integer,dist_,ref2_Integer,ang_,ref3_Integer,dihed_}]:=AppendTo[$zmAtomStore,
  {el,With[{
    p1=Last@$zmAtomStore[[ref1]],
    p2=Last@$zmAtomStore[[ref2]],
    p3=Last@$zmAtomStore[[ref3]]},
   With[{n=(p2-p1)(p3-p1)},
    p1+dist*Normalize@(
     If[dihed!=0,
      RotationMatrix[dihed,p2-p1],
      IdentityMatrix@3
      ].RotationMatrix[ang,n].(p2-p1)
     )
    ]]}];
   zmConvert[{a1_Integer,a2_Integer,t___}]:=AppendTo[$zmBondStore,
    Bond[
     a1-Length@Cases[$zmAtomStore[[;;a1]],{"X",___}],
     a2-Length@Cases[$zmAtomStore[[;;a2]],{"X",___}],
     t]];

ZMatrixToStandard[ZMTable:{_List..}]:=Block[{$zmAtomStore={},$zmBondStore={} },
 zmConvert/@ZMTable;
 {DeleteCases[$zmAtomStore,{"X",___}],$zmBondStore}
 ];

Molecule[ZMTable:{_List..}]:=With[{zmConf=ZMatrixToStandard@ZMTable},
 Molecule[Atom@@@First@zmConf,Last@zmConf]
 ];