Guidelines for performing PIEDA

Performing PIEDA for some system consists of at most 4 steps.

1. getting the BDA energies (only for covalently bound fragments).

Use the provided ala10.bda.inp as is, except you need to replace:

 n0bda=   10 r0bda(1)= ...

with the data for your own system. This data is now generated during any FMO
run and is punched into the .dat file. An easy way to generate it is to
take any FMO input, and set nbody=0. Then n0bda and r0bda can be found in
the punch file (.dat).

The necessary input options are:

 $fmoprp modorb=1 ipieda=1 $end

The basis set/HMOs have to be adjusted, too.

This run cannot use GDDI with more than 1 group. Parallel runs using 1
group are permitted.
This step is not needed for molecular clusters.

2. free state (FMO0) ala10.free.inp

Take a typical FMO input file and set the following options:

 $contrl runtyp=fmo0 $end
 $fmoprp maxaoc=N modorb=1 $end

maxaoc has to be set to the maximum number of orbitals per fragment.
This value is printed during all usual FMO runs (except runtyp=fmo0!):

 Max AOs per frg:         7

Here N is 7 so you would have maxaoc=7.

Finally, you must save the file F30 during this run and use it in step 3.
The file contains the free state orbitals.

3. free state polarisation (PL0) ala10.pl0.inp

Rename file F30 from step 2 into F40 for the input in this step and copy
it to the job scratch directory on the master node where you will run
the job. Do not copy it to slave nodes.

Input file:

For covalently bound fragments, copy r0bda(1) from step 1. 
Make sure n0bda is NOT set. Copy e0bda from the .dat file of step 1.

Copy efmo0 from the .dat file of step 2.

 $fmoprp irest=2 maxit=1 modorb=1 ipieda=2 modpar=77 $end

4. PIEDA ala10.pieda-full.inp  

For covalently bound fragments, copy r0bda(1) from step 1.
Make sure n0bda is NOT set. Copy e0bda from the .dat file of step 1.

Extract efmo0 from the .dat file of step 2 and add it to your input.

Similarly, extract epl0ds and eint0 from the .dat file of step 3.

 $fmoprp modorb=1 ipieda=2 $end

--------------

We used respap=0 resppc=4 in steps 1-4 in the original PIEDA publication to
improve the accuracy. However, since then we came up with an approximation
improvement and now one can reduce these thresholds toi, say, respap=0
resppc=3. In the special PIEDA case described below, the default values may be used.

It is possible to greatly simplify the calculations by doing the analysis
based on the PL state (accomplished with ipieda=1, see ala10.pieda-full.inp). 
In this case the free state is not used and only pair interactions are 
decomposed. For this run, steps 2 and 3 are not needed, and, if you wish, you 
can also omit step 1 (in which case connected pairs would have -14 a.u. 
interaction energy).
Note that you will be required to provide r0bda, however, e0bda may be
omitted.


PIEDA program limitations.

It is possible to use modorb=3 in all files, replacing modorb=1.

It is permitted to do RHF instead of MP2 runs (loosing the dispersion).
All jobs in the 4 steps must be of the same type though.

Multilayer runs are not supported. DFT and MCSCF cannot be used. 
Any CC (supported by FMO) can be used instead of MP2.  

At least some approximation threshold must be set for covalently bound
systems (any one of RESPAP/RESPPC/RESDIM/RITRIM/RCORSD). 

Gradient runs are not supported. nbody=2 must be used. Restarts are
supported during step 4 only (irest=2 only). 

PIEDA/PCM is not supported.