OPTIMIZE of MCSCF crashes for > 255 basis functions

[hansh]

Dear Dalton experts,

I do not know whether this is a more general problem with the code or perhaps specific to my installation(s):
As in the previous version 2011 of DALTON, also in the current one, all geometry optimization jobs for
MCSCF energies crash upon starting the Static Property Section if more than 255 basis functions are
used (for smaller basis sets they run without any trouble). In case of DALTON 2013, the following error
messages appear:

DALTON-2013.1-Source/build_gfortran$ ./dalton: line 475: 21960 Bus error (core dumped)
Error in /home/hh/Dalton13/DALTON-2013.1-Source/build_gfortran/dalton.x, exit code 135

WARNING for ERROR (Dalton program stopped with exit code 135)

Below I am attaching the relevant parts of the in/output. Many thanks in advance for all hints
and assistance for overcoming this problem, Hans

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

Compilation information
-----------------------

Who compiled | hh
Host | CoolerMaster2012
System | Linux-3.5.0-46-generic
CMake generator | Unix Makefiles
Processor | x86_64
64-bit integers | ON
MPI | OFF
Fortran compiler | /usr/bin/gfortran
Fortran compiler version | GNU Fortran (Ubuntu/Linaro 4.6.3-1ubuntu5) 4.6.3
C compiler | /usr/bin/gcc
C compiler version | gcc (Ubuntu/Linaro 4.6.3-1ubuntu5) 4.6.3
C++ compiler | /usr/bin/g++
C++ compiler version | g++ (Ubuntu/Linaro 4.6.3-1ubuntu5) 4.6.3
Static linking | ON
Last Git revision | 653a3c9fdcde6b463c7e208ddf10abd66f7c54f6
Configuration time | 2014-03-02 20:49:50.253875


Content of the .dal input file
----------------------------------

**DALTON INPUT
.RUN WAVEFUNCTION
.OPTIMIZE
**WAVE FUNCTIONS
.HF
.MP2
.MCSCF
.NSYM
8
*SCF IN
.DOUBLY OCCUPIED
2 1 1 0 2 0 0 0
*CONFIGURATION INPUT
.SYMMETRY
2
.SPIN MUL
2
.INACTIVE
1 0 0 0 1 0 0 0
.ELECTRONS
7
.CAS SPACE
4 3 3 1 3 2 2 1
*OPTIMIZATION
.MAX CI
5
.MAX MACRO ITERATIONS
80
.MAX MICRO ITERATIONS
80
.MAXABS
40
.MAXAPM
25
!*POPULATION ANALYSIS
!.MULLIKEN
**PROPERTIES
.VIBANA
**END OF DALTON INPUT


Content of the .mol file
----------------------------

BASIS
cc-pV6Z
Wave Functions: RHF, MP2, and MCSCF (CAS) / cc-pVQ/5Z/6Z
Test Purpose: N = 11 Z var
Atomtypes=1
Charge=6. Atoms=2
N 0.0 0.0 -1.2
N 0.0 0.0 1.2


*******************************************************************
*********** Output from DALTON general input processing ***********
*******************************************************************

--------------------------------------------------------------------------------
Overall default print level: 0
Print level for DALTON.STAT: 1

HERMIT 1- and 2-electron integral sections will be executed
"Old" integral transformation used (limited to max 255 basis functions)
Wave function sections will be executed (SIRIUS module)
--------------------------------------------------------------------------------


Chosen parameters for *OPTIMI :
-------------------------------

Default 1st order method will be used: BFGS update.
Optimization will be performed in redundant internal coordinates.
Model Hessian will be used as initial Hessian.
The model Hessian parameters of Roland Lindh will be used.


Trust region method will be used to control step (default).

Convergence threshold for gradient set to : 1.00D-04
Convergence threshold for energy set to : 1.00D-06
Convergence threshold for step set to : 1.00D-04
Number of convergence criteria set to : 2
-------------------------------------------------------------------------------

last lines of output:

.------------------------------------------------.
| Starting in Static Property Section (ABACUS) - |
`------------------------------------------------'



Date and time (Linux) : Tue Mar 4 16:14:25 2014
Host name : CoolerMaster2012

[hansh]

To update my previous post, let me add that obviously the problem is not related to invoking
the 64-bit integer installation, but it also occurs in 32-bit integer version (compiled by ifort) and
then flags the following error:

/Dalton13/DALTON-2013.2-Source/G_In1$ forrtl: severe (174): SIGSEGV, segmentation fault occurred
Image PC Routine Line Source
dalton_copy.x 000000000115500B Unknown Unknown Unknown
dalton_copy.x 0000000001153853 Unknown Unknown Unknown
dalton_copy.x 0000000000FC3B56 Unknown Unknown Unknown
dalton_copy.x 0000000000FCAED9 Unknown Unknown Unknown
dalton_copy.x 000000000041377D Unknown Unknown Unknown
dalton_copy.x 00000000010F77BA Unknown Unknown Unknown
dalton_copy.x 00000000010F719F Unknown Unknown Unknown
dalton_copy.x 0000000000413DC3 Unknown Unknown Unknown
dalton_copy.x 000000000040663B Unknown Unknown Unknown
dalton_copy.x 00000000004057EC Unknown Unknown Unknown
libc.so.6 00007F0C70FC176D Unknown Unknown Unknown
dalton_copy.x 00000000004056E9 Unknown Unknown Unknown

Error in /home/hh/Dalton13/DALTON-2013.2-Source/b_iftn/dalton.x, exit code 174

Again the job crashes upon entering the
"Starting in Static Property Section (ABACUS)".

Details of the installation employed:

Who compiled | hh
Host | CoolerMaster2012
System | Linux-3.5.0-46-generic
CMake generator | Unix Makefiles
Processor | x86_64
64-bit integers | OFF
MPI | OFF
Fortran compiler | /opt/intel/composer_xe_2011_sp1.9.293/bin/intel64/
| ifort
Fortran compiler version | ifort (IFORT) 12.1.3 20120212
C compiler | /usr/bin/gcc
C compiler version | gcc (Ubuntu/Linaro 4.6.3-1ubuntu5) 4.6.3
C++ compiler | /usr/bin/g++
C++ compiler version | g++ (Ubuntu/Linaro 4.6.3-1ubuntu5) 4.6.3
BLAS | /usr/lib/libblas.so
LAPACK | /usr/lib/liblapack.so
Static linking | ON
Last Git revision | 0be1cb1e4b1461124063216a63270ec99bfebef1
Configuration time | 2014-03-06 15:30:55.456738

All suggestions for fixing this problem are greatly appreciated.
Many thanks in advance, Hans

[erik]

Dear Hans

What was your setup command? (you can see this in your "build" directory in the file "setup_command")

I don't know which compilers you are pointing to at your system but I have run such calculations (above 255 basis functions) daily using 32-bit intel (ifort) compilers, but something obviously goes wrong in your case. You cna try to post the mol file as well so we can try to reproduce the error

With best regards
Erik

[hansh]

Dear Erik,

many thanks for your kind response. I invoked

./setup --fc=ifort --cc=gcc --cxx=g++ --static b_iftn

for setting up my installation. The versions of the compilers are those displayed in my second post above,
and the .dal and .mol input is the same as in my first post above.

In the meantime, I also tried an installation based on gfortran (without 64-bin integers), i.e., the following:

Compilation information
-----------------------

Who compiled | hh
Host | CoolerMaster2012
System | Linux-3.5.0-46-generic
CMake generator | Unix Makefiles
Processor | x86_64
64-bit integers | OFF
MPI | OFF
Fortran compiler | /usr/bin/gfortran
Fortran compiler version | GNU Fortran (Ubuntu/Linaro 4.6.3-1ubuntu5) 4.6.3
C compiler | /usr/bin/gcc
C compiler version | gcc (Ubuntu/Linaro 4.6.3-1ubuntu5) 4.6.3
C++ compiler | /usr/bin/g++
C++ compiler version | g++ (Ubuntu/Linaro 4.6.3-1ubuntu5) 4.6.3
BLAS | /usr/lib/libblas.so
LAPACK | /usr/lib/liblapack.so
Static linking | OFF
Last Git revision | 0be1cb1e4b1461124063216a63270ec99bfebef1
Configuration time | 2014-03-07 12:56:55.445014

Unfortunately, after having entered the ABACUS module and then computing for a short while
(for the same input as in my first post above), the job also crashed with the following error messages:

line 476: 27945 Segmentation fault (core dumped) ./dalton_copy.x
Error in /home/hh/Dalton13/DALTON-2013.2-Source/b_gftn/dalton.x, exit code 139

Once again many thanks for all your help, Hans

[erik]

Dear Hans

Unfortunately, I think you have hid an actual bug. My previous calculations did not invoke MCSCF geometry optimizations, but was only MCSCF response or wave function optimization calculations. Here everything has run fine, also with more than 255 basis functions...
However, for the specific case of geometry optimizations, I have tried your example (with gfortran compilers 4.4.6, running a serial calculation) - and I get the same error as you starting from when DALTON enters ABACUS.
This is to be reported as a bug...

Sorry that I am not of more help at the moment, but I am not in-dept familiar with the geometry optimization parts of the code.

With best regards
Erik

[kennethruud]

Hi!

I cannot quite recall the limitations, but it could very well be that Dalton cannot calculate geometrical derivatives for MCSCF wave functions when there are more than 255 basis functions. It is for sure such that it is not possible to do such calculations for molecular Hessians.

Peter Taylor may be able to give feedback as to whether a routine called PTRAN would work for more than 255 bfs. I suspect it does not.

You could do us a favour if you could try the following to calculations:

1) Run with a general print level of 4
**PROPERTIES
.PRINT
4
....

2) Run with the following input
**PROPERTIES
.....
*TWOEXP
.NOPV
.....


Best regards,

Kenneth

[taylor]

Yes, this is a known problem, and Kenneth points to the problem, although it's actually post PTRAN and in a subsequent sorting step that things are likely to go awry. Some background...

One contribution to the molecular gradient is from what is essentially a trace operation in which two-electron integrals differentiated with respect to the nuclear coordinates of the basis functions are multiplied by elements of the second-order reduced density matrix (2-matrix). For SCF wave functions the latter is trivially expressed in terms of products of first-order density matrix elements, but this is not in general true of correlated wave functions like MCSCF. What is done in Dalton is to take the trace operation
Sum_{pqrs} P_{pqrs} (pq|rs)'
in the MO basis, where P is the 2-matrix and the prime indicates differentiation of the integral, and rewrite this first in terms of symmetry orbitals a, b, etc. as
Sum_{abcd} P_{abcd} (ab|cd)',
where
P_{abcd} = Sum_{pqrs} P_{pqrs} C_{ap} C_{bq} C_{cr} C_{ds}
and C is the matrix of MO coefficients in the symmetry orbital basis. Incidentally, one reason for doing this rather than transforming the derivative integrals from the symmetry orbital basis into the MO basis is that there can be an order of magnitude more of these derivative integrals than 2-matrix elements in the symmetry orbital basis.

The routine PTRAN accomplishes this transformation. However, the transformed 2-matrix is not in the correct order to be "matched up" to the derivative integrals as they are calculated, and the 2-matrix needs to be sorted to do this. Further, the sorted 2-matrix is not used directly --- blocks of it are further "back-transformed" to the set of distinct atomic orbitals from which the derivative integrals are constructed, because this further reduces the operation count. All this sorting is quite messy, and was written a long time ago when 255 basis functions was a pretty big deal... And of course one can index up to 255 in 8 bits, which meant labels could be stored packed four to a 32-bit integer word. And unfortunately, perhaps, this has never been rewritten entirely so that either 64-bit integers, or pairs of 32-bit integers, are used for labels.

Kenneth was polite enough not to state that this was ever going to be rewritten it would likely be only me that could/would do it, but that is effectively the case. I would probably not do this even if I had the time (which I certainly don't --- my day job running a computing facility and life sciences research centre doesn't give me much opportunity for such things at the moment!) but because I think there are better ways to implement all this on modern machines, as I did several years ago for the conventional integral transformation problem. Whether any of that will find its way into Dalton2015 is an open question: it will \emph{not} be in Dalton2014.

Finally, and this may have been some sort of "dry run" before tackling a polyatomic system, but I notice you are running a diatomic. It is almost impossible, in terms of computer time, to make an analytical derivative approach to optimization followed by vibrational frequencies pay in this case. Unless one has absolutely no idea whatever of the equilibrium bond distance in a diatomic, a pretty good estimate of both r_e and omega_e can be obtained from a three-point fit around a guess r_e, especially if that fit is done in 1/r rather than r. Further, with six or seven points one could fit a quartic (r itself will do here as the independent variable) and get an even better estimate of r_e and omega_e as well as a very good estimate of omega_e x_e. Even with Dalton (and sparing the collective modesty of the Dalton authors) it is hard to see that doing two?, three? geometry steps and then analytical second derivatives can be competitive with doing six/seven single-point energies and then fitting.

Best regards
Pete

[hansh]

Dear Kenneth, dear Pete,

many thanks for your kind response and the clarifying and detailed remarks. While it is a pity that the current version of DALTON does not support geometry optimization for MCSCF with > 255 basis functions, on the other hand it certainly makes sense to ponder whether rewriting the code is worth the work and efforts needed. Perhaps, such work might be better invested into implementing a geometry optimization procedure for other multiconfiguration methods, like, e.g., for the LUCITA module???

And, yes, for PECs of diatomics one can do well without an automatic optimization. Actually, I first encountered the problem with > 255 basis functions when studying a polyatomic system whose proper description needs a multiconfiguration wavefunction. The diatomic example that I invoked in my post above was only meant to provide a simpler illustration of the issue.

Once again, warm thanks for your kind support and best wishes, Hans