The Maple command remove can be used to remove items from a list:

with(QuantumChemistry);
mol1 := MolecularGeometry("1,2-dichlorobenzene");
mol1 := [["Cl", 2.003, 1.6182, -0.0009], 

  ["Cl", 2.0039, -1.6175, 0.0017], ["C", 0.54, 0.6976, -0.0001], 

  ["C", 0.5404, -0.6972, -0.001], ["C", -0.6682, 1.3947, 0.0008], 

  ["C", -0.6674, -1.395, -0.0012], ["C", -1.876, 0.697, 0.0009], 

  ["C", -1.8756, -0.6978, -0.0002], 

  ["H", -0.6871, 2.4817, 0.0015], 

  ["H", -0.6859, -2.482, -0.0015], 

  ["H", -2.8169, 1.2396, 0.0016], ["H", -2.8161, -1.2409, -0.0001]

  ]

PlotMolecule(mol1);

mol2 := remove(has, mol1, "H");
mol2 := [["Cl", 2.003, 1.6182, -0.0009], 

  ["Cl", 2.0039, -1.6175, 0.0017], ["C", 0.54, 0.6976, -0.0001], 

  ["C", 0.5404, -0.6972, -0.001], ["C", -0.6682, 1.3947, 0.0008], 

  ["C", -0.6674, -1.395, -0.0012], ["C", -1.876, 0.697, 0.0009], 

  ["C", -1.8756, -0.6978, -0.0002]]

PlotMolecule(mol2);

The PlotMolecule command automatically draws chemical bonds when atoms are "close" relative to an expected bond distance between atoms.  In this molecule the C-S bond is long relative to standard C-S bonds such that PlotMolecule chooses not to represent the bond explicitly.  If you want to override the automatic bonding, you can pass bond connectivity to the PlotMolecule command.  Note that bond connectivity can be returned by MolecularData and MolecularGeometry.  Here is the molecule with the long C-S bond added to PlotMolecule:

with(QuantumChemistry);
file := "mol.xyz";
mol := ReadXYZ(file):
dist := BondDistances(mol);
bonding := indices(dist);
bonding := {seq(convert(b, set), b in bonding)};
bonding := bonding union {{1, 5}};
PlotMolecule(mol, bonds = bonding);

PlotMolecule-AdjustBonding.mw

Maple with the Quantum Chemistry Toolbox (QCT) with the QuantumChemistry:-MolecularGeometry command can retrieve 3D coordinates and connectivity from a database, using the molecule's name or cid number.  

Maple with the Quantum Chemistry Toolbox (QCT) can retrieve 3D coordinates and connectivity from a database, using the molecule's name or cid number, with the QuantumChemistry:-MolecularGeometry command.  The "Simplified Molecular-Input Line-Entry System" (SMILES) provides a shorthand notation for the connectivity of atoms in a molecule.  Some support for SMILES is likely to appear in the upcoming, next release of QCT.

Hi Oliveira,

       Here are some Maple books that I like:

1) Andre Heck's Introduction to Maple
2) Programming Guide https://www.maplesoft.com/documentation_center/maple2021/ProgrammingGuide.pdf

David

Here is the PES of ArF for bond lengths between 2.2 and 3.6 from HartreeFock and Parametric2RDM in the cc-pVDZ basis set. 

> with(QuantumChemistry):

> mols := [seq([["Ar",0,0,0],["F",0,0,2+0.2*i]], i=1..8)];
         mols := [[["Ar", 0, 0, 0], ["F", 0, 0, 2.2]],

           [["Ar", 0, 0, 0], ["F", 0, 0, 2.4]],

           [["Ar", 0, 0, 0], ["F", 0, 0, 2.6]],

           [["Ar", 0, 0, 0], ["F", 0, 0, 2.8]],

           [["Ar", 0, 0, 0], ["F", 0, 0, 3.0]],

           [["Ar", 0, 0, 0], ["F", 0, 0, 3.2]],

           [["Ar", 0, 0, 0], ["F", 0, 0, 3.4]],

           [["Ar", 0, 0, 0], ["F", 0, 0, 3.6]]]

> energies_hf := map(Energy, mols, basis="ccpvdz");
  energies_hf := [-626.151444023726, -626.162686642662,

    -626.167860375677, -626.170173135759, -626.171174175847,

    -626.171382087638, -626.171677653194, -626.171758395594]

> energies_p2rdm := map(Energy, mols, basis="ccpvdz", method=Parametric2RDM);
  energies_p2rdm := [-626.472462673700, -626.480383384324,

    -626.483830355819, -626.485335054846, -626.485825123186,

    -626.485799957140, -626.485973310631, -626.485977743599]

> x := [seq(i[2][4], i in mols)];
         x := [2.2, 2.4, 2.6, 2.8, 3.0, 3.2, 3.4, 3.6]

> xy_hf := [seq([x[i],energies_hf[i]], i=1..nops(x))];
 xy_hf := [[2.2, -626.151444023726], [2.4, -626.162686642662],

   [2.6, -626.167860375677], [2.8, -626.170173135759],

   [3.0, -626.171174175847], [3.2, -626.171382087638],

   [3.4, -626.171677653194], [3.6, -626.171758395594]]

> xy_p2rdm := [seq([x[i],energies_p2rdm[i]], i=1..nops(x))];
xy_p2rdm := [[2.2, -626.472462673700], [2.4, -626.480383384324],

  [2.6, -626.483830355819], [2.8, -626.485335054846],

  [3.0, -626.485825123186], [3.2, -626.485799957140],

  [3.4, -626.485973310631], [3.6, -626.485977743599]]

> plot(spline(x,energies_p2rdm,z), axes=boxed, view=[2.2..3.6,-626.47..-626.488]);

In the Parametric2RDM PES there is a very slight increase around 3.0 Angstroms.  From your note it appears that a similar increase is present in CoupledCluster.  The PES has flattened by 3.0 Angstroms because the atoms are breaking apart.  The increase is within the error of quantum chemistry methods.  Counterpoise correction would not likely matter.

David

Hi -

Use:

LinearAlgebra:-Eigenvalues(A);

David

Hi,

You should be able to do something like:

fsolve( diif( f(x), x) = 0, x);

David