1 /*-
2  * Copyright (c) 2012-2017 Ilya Kaliman
3  *
4  * Redistribution and use in source and binary forms, with or without
5  * modification, are permitted provided that the following conditions
6  * are met:
7  *
8  * 1. Redistributions of source code must retain the above copyright
9  *    notice, this list of conditions and the following disclaimer.
10  * 2. Redistributions in binary form must reproduce the above copyright
11  *    notice, this list of conditions and the following disclaimer in the
12  *    documentation and/or other materials provided with the distribution.
13  *
14  * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS "AS IS" AND
15  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
16  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
17  * ARE DISCLAIMED. IN NO EVENT SHALL AUTHOR OR CONTRIBUTORS BE LIABLE
18  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
19  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
20  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
21  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
22  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
23  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
24  * SUCH DAMAGE.
25  */
26 
27 #ifndef LIBEFP_EFP_H
28 #define LIBEFP_EFP_H
29 
30 #include <stddef.h>
31 
32 /** \file efp.h
33  * Public libefp interface.
34  *
35  * A note on units: masses are in AMU, everything else is in atomic units.
36  */
37 
38 #ifdef __cplusplus
39 extern "C" {
40 #endif
41 
42 /** Version string. */
43 #define LIBEFP_VERSION_STRING "1.5.0"
44 
45 /** Result of an operation. */
46 enum efp_result {
47 	/** Operation was successful. */
48 	EFP_RESULT_SUCCESS = 0,
49 	/** Fatal error has occurred. */
50 	EFP_RESULT_FATAL,
51 	/** Insufficient memory. */
52 	EFP_RESULT_NO_MEMORY,
53 	/** File not found. */
54 	EFP_RESULT_FILE_NOT_FOUND,
55 	/** Syntax error. */
56 	EFP_RESULT_SYNTAX_ERROR,
57 	/** Unknown EFP fragment. */
58 	EFP_RESULT_UNKNOWN_FRAGMENT,
59 	/** Polarization SCF procedure did not converge. */
60 	EFP_RESULT_POL_NOT_CONVERGED
61 };
62 
63 /** Flags to specify EFP energy terms. */
64 enum efp_term {
65 	/** EFP/EFP electrostatics. */
66 	EFP_TERM_ELEC = 1 << 0,
67 	/** EFP/EFP polarization. */
68 	EFP_TERM_POL = 1 << 1,
69 	/** EFP/EFP dispersion. */
70 	EFP_TERM_DISP = 1 << 2,
71 	/** EFP/EFP exchange repulsion. */
72 	EFP_TERM_XR = 1 << 3,
73 	/** EFP/EFP charge transfer, reserved for future use. */
74 	EFP_TERM_CHTR = 1 << 4,
75 	/** Ab initio/EFP electrostatics. */
76 	EFP_TERM_AI_ELEC = 1 << 5,
77 	/** Ab initio/EFP polarization. */
78 	EFP_TERM_AI_POL = 1 << 6,
79 	/** Ab initio/EFP dispersion, reserved for future use. */
80 	EFP_TERM_AI_DISP = 1 << 7,
81 	/** Ab initio/EFP exchange repulsion, reserved for future use. */
82 	EFP_TERM_AI_XR = 1 << 8,
83 	/** Ab initio/EFP charge transfer, reserved for future use. */
84 	EFP_TERM_AI_CHTR = 1 << 9
85 };
86 
87 /** Fragment-fragment dispersion damping type. */
88 enum efp_disp_damp {
89 	EFP_DISP_DAMP_OVERLAP = 0, /**< Overlap-based damping (default). */
90 	EFP_DISP_DAMP_TT,          /**< Tang-Toennies damping. */
91 	EFP_DISP_DAMP_OFF          /**< No dispersion damping. */
92 };
93 
94 /** Fragment-fragment electrostatic damping type. */
95 enum efp_elec_damp {
96 	EFP_ELEC_DAMP_SCREEN = 0,  /**< SCREEN-controlled damping (default). */
97 	EFP_ELEC_DAMP_OVERLAP,     /**< Overlap-based damping. */
98 	EFP_ELEC_DAMP_OFF          /**< No electrostatic damping. */
99 };
100 
101 /** Fragment-fragment polarization damping type. */
102 enum efp_pol_damp {
103 	EFP_POL_DAMP_TT = 0,       /**< Tang-Toennies like damping (default). */
104 	EFP_POL_DAMP_OFF           /**< No polarization damping. */
105 };
106 
107 /** Describes the way fragment coordinates are specified. */
108 enum efp_coord_type {
109 	/** Coordinates of center of mass of a fragment and Euler angles. */
110 	EFP_COORD_TYPE_XYZABC = 0,
111 	/** Coordinates of three points belonging to a fragment. */
112 	EFP_COORD_TYPE_POINTS,
113 	/** Coordinates of fragment center of mass and its rotation matrix. */
114 	EFP_COORD_TYPE_ROTMAT
115 };
116 
117 /** Driver used for solving polarization equations. */
118 enum efp_pol_driver {
119 	/** Iterative solution of polarization equations. */
120 	EFP_POL_DRIVER_ITERATIVE = 0,
121 	/** Direct solution of polarization equations. */
122 	EFP_POL_DRIVER_DIRECT
123 };
124 
125 /** \struct efp
126  * Main EFP opaque structure.
127  */
128 struct efp;
129 
130 /** Options controlling EFP computation. */
131 struct efp_opts {
132 	/** Specifies which energy terms should be computed.
133 	 * This field is a collection of bit flags where each bit specifies
134 	 * whether the term is enabled (bit is set to 1) or disabled (bit is
135 	 * set to 0). To enable the term, use bitwise OR with the corresponding
136 	 * efp_term constant (e.g., terms |= EFP_TERM_ELEC). To disable the
137 	 * term, use bitwise AND NOT (e.g., terms &= ~EFP_TERM_POL). */
138 	unsigned terms;
139 	/** Dispersion damping type (see #efp_disp_damp). */
140 	enum efp_disp_damp disp_damp;
141 	/** Electrostatic damping type (see #efp_elec_damp). */
142 	enum efp_elec_damp elec_damp;
143 	/** Polarization damping type (see #efp_pol_damp). */
144 	enum efp_pol_damp pol_damp;
145 	/** Driver used to find polarization induced dipoles. */
146 	enum efp_pol_driver pol_driver;
147 	/** Enable periodic boundary conditions if nonzero. */
148 	int enable_pbc;
149 	/** Enable fragment-fragment interaction cutoff if nonzero. */
150 	int enable_cutoff;
151 	/** Cutoff distance for fragment-fragment interactions. */
152 	double swf_cutoff;
153 };
154 
155 /** EFP energy terms. */
156 struct efp_energy {
157 	/**
158 	 * EFP/EFP electrostatic energy. */
159 	double electrostatic;
160 	/**
161 	 * Charge penetration energy from overlap-based electrostatic
162 	 * damping. Zero if overlap-based damping is turned off. */
163 	double charge_penetration;
164 	/**
165 	 * Interaction energy of EFP electrostatics with point charges. */
166 	double electrostatic_point_charges;
167 	/**
168 	 * All polarization energy goes here. Polarization is computed
169 	 * self-consistently so it can't be separated into EFP/EFP and AI/EFP
170 	 * parts. */
171 	double polarization;
172 	/**
173 	 * EFP/EFP dispersion energy. */
174 	double dispersion;
175 	/**
176 	 * AI/EFP dispersion energy. */
177 	double ai_dispersion;
178 	/**
179 	 * EFP/EFP exchange-repulsion energy. */
180 	double exchange_repulsion;
181 	/**
182 	 * Sum of all the above energy terms. */
183 	double total;
184 };
185 
186 /** EFP atom info. */
187 struct efp_atom {
188 	char label[32];   /**< Atom label. */
189 	double x;         /**< X coordinate of atom position. */
190 	double y;         /**< Y coordinate of atom position. */
191 	double z;         /**< Z coordinate of atom position. */
192 	double mass;      /**< Atom mass. */
193 	double znuc;      /**< Nuclear charge. */
194 };
195 
196 /**
197  * Callback function which is called by libefp to obtain electric field in the
198  * specified points.
199  *
200  * This function is used to obtain the electric field from electrons
201  * in the \a ab \a initio part. This callback is called by libefp during
202  * polarization calculation. Libefp supplies the \p xyz array with
203  * coordinates of the points where the field should be computed.
204  *
205  * \param[in] n_pt Number of points in \p xyz array.
206  *
207  * \param[in] xyz Coordinates of points where electric field should be
208  * computed. The size of this array must be [3 * \p n_pt] elements.
209  *
210  * \param[out] field Computed \a x \a y \a z components of electric field. The
211  * size of this array must be at least [3 * \p n_pt] elements.
212  *
213  * \param[in] user_data User data which was specified during initialization.
214  *
215  * \return The implemented function should return ::EFP_RESULT_FATAL on error
216  * and ::EFP_RESULT_SUCCESS if the calculation has succeeded.
217  */
218 typedef enum efp_result (*efp_electron_density_field_fn)(size_t n_pt,
219     const double *xyz, double *field, void *user_data);
220 
221 /**
222  * Get a human readable banner string with information about the library.
223  *
224  * \return Banner string, zero-terminated.
225  */
226 const char *efp_banner(void);
227 
228 /**
229  * Print libefp banner to stdout.
230  */
231 void efp_print_banner(void);
232 
233 /**
234  * Create a new efp object.
235  *
236  * \return A new efp object or NULL on error.
237  */
238 struct efp *efp_create(void);
239 
240 /**
241  * Get default values of simulation options.
242  *
243  * \param[out] opts Structure to store the defaults. See ::efp_opts.
244  */
245 void efp_opts_default(struct efp_opts *opts);
246 
247 /**
248  * Set the error log callback function.
249  *
250  * The callback function can be used to print verbose diagnostic messages from
251  * libefp. By default libefp prints to stderr using the function shown below.
252  * Logging can be disabled by setting \a cb to NULL.
253  *
254  * \code
255  * void log_cb(const char *msg)
256  * {
257  *     fprintf(stderr, "LIBEFP: %s\n", msg);
258  * }
259  * \endcode
260  *
261  * \param[in] cb Error log callback function or NULL if none.
262  */
263 void efp_set_error_log(void (*cb)(const char *));
264 
265 /**
266  * Set computation options.
267  *
268  * \param[in] efp The efp structure.
269  *
270  * \param[in] opts New options for EFP computation. See ::efp_opts.
271  *
272  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
273  */
274 enum efp_result efp_set_opts(struct efp *efp, const struct efp_opts *opts);
275 
276 /**
277  * Get currently set computation options.
278  *
279  * \param[in] efp The efp structure.
280  *
281  * \param[out] opts Current options for EFP computation. See ::efp_opts.
282  *
283  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
284  */
285 enum efp_result efp_get_opts(struct efp *efp, struct efp_opts *opts);
286 
287 /**
288  * Add EFP potential from a file.
289  *
290  * \param[in] efp The efp structure.
291  *
292  * \param[in] path Path to the EFP potential file, zero terminated string.
293  *
294  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
295  */
296 enum efp_result efp_add_potential(struct efp *efp, const char *path);
297 
298 /**
299  * Add a new fragment to the EFP subsystem.
300  *
301  * \param[in] efp The efp structure.
302  *
303  * \param[in] name Fragment name, zero terminated string.
304  *
305  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
306  */
307 enum efp_result efp_add_fragment(struct efp *efp, const char *name);
308 
309 /**
310  * Prepare the calculation.
311  *
312  * New fragments must NOT be added after a call to this function.
313  *
314  * \param[in] efp The efp structure.
315  *
316  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
317  */
318 enum efp_result efp_prepare(struct efp *efp);
319 
320 /**
321  * Skip interactions between the fragments.
322  *
323  * \param[in] efp The efp structure.
324  *
325  * \param[in] i Index of the first fragment.
326  *
327  * \param[in] j Index of the second fragment.
328  *
329  * \param[in] value Specifies whether to skip i-j interactions (true/false).
330  *
331  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
332  */
333 enum efp_result efp_skip_fragments(struct efp *efp, size_t i, size_t j,
334     int value);
335 
336 /**
337  * Set the callback function which computes electric field from electrons
338  * in \a ab \a initio subsystem.
339  *
340  * \param[in] efp The efp structure.
341  *
342  * \param[in] fn The callback function. See ::efp_electron_density_field_fn.
343  *
344  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
345  */
346 enum efp_result efp_set_electron_density_field_fn(struct efp *efp,
347     efp_electron_density_field_fn fn);
348 
349 /**
350  * Set user data to be passed to ::efp_electron_density_field_fn.
351  *
352  * \param[in] efp The efp structure.
353  *
354  * \param[in] user_data User data which will be passed as a last parameter to
355  * ::efp_electron_density_field_fn.
356  *
357  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
358  */
359 enum efp_result efp_set_electron_density_field_user_data(struct efp *efp,
360     void *user_data);
361 
362 /**
363  * Setup arbitrary point charges interacting with EFP subsystem.
364  *
365  * This can be used to compute contributions from \a ab \a initio nuclei.
366  *
367  * \param[in] efp The efp structure.
368  *
369  * \param[in] n_ptc Number of point charges.
370  *
371  * \param[in] ptc Array of \p n_ptc elements with charge values.
372  *
373  * \param[in] xyz Array of [3 * \p n_ptc] elements with \a x \a y \a z
374  * coordinates of charge positions.
375  *
376  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
377  */
378 enum efp_result efp_set_point_charges(struct efp *efp, size_t n_ptc,
379     const double *ptc, const double *xyz);
380 
381 /**
382  * Get the number of currently set point charges.
383  *
384  * \param[in] efp The efp structure.
385  *
386  * \param[out] n_ptc Number of point charges.
387  *
388  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
389  */
390 enum efp_result efp_get_point_charge_count(struct efp *efp, size_t *n_ptc);
391 
392 /**
393  * Get values of currently set point charges.
394  *
395  * \param[in] efp The efp structure.
396  *
397  * \param[out] ptc Array of \p n_ptc elements where charges will be stored.
398  *
399  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
400  */
401 enum efp_result efp_get_point_charge_values(struct efp *efp, double *ptc);
402 
403 /**
404  * Set values of point charges.
405  *
406  * \param[in] efp The efp structure.
407  *
408  * \param[in] ptc Array of \p n_ptc elements with charge values.
409  *
410  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
411  */
412 enum efp_result efp_set_point_charge_values(struct efp *efp, const double *ptc);
413 
414 /**
415  * Get coordinates of currently set point charges.
416  *
417  * \param[in] efp The efp structure.
418  *
419  * \param[out] xyz Array where \a x \a y \a z coordinates of point charges will
420  * be stored. The size of the array must be at least [3 * \a n] elements, where
421  * \a n is the total number of point charges.
422  *
423  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
424  */
425 enum efp_result efp_get_point_charge_coordinates(struct efp *efp, double *xyz);
426 
427 /**
428  * Set coordinates of point charges.
429  *
430  * \param[in] efp The efp structure.
431  *
432  * \param[in] xyz Array with \a x \a y \a z coordinates of point charges. The
433  * size of the array must be at least [3 * \a n] elements, where \a n is the
434  * total number of point charges.
435  *
436  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
437  */
438 enum efp_result efp_set_point_charge_coordinates(struct efp *efp,
439     const double *xyz);
440 
441 /**
442  * Get gradient on point charges from EFP subsystem.
443  *
444  * \param[in] efp The efp structure.
445  *
446  * \param[out] grad For each point charge \a x \a y \a z components of energy
447  * gradient are stored. The size of this array must be at least [3 * \a n]
448  * elements, where \a n is the total number of point charges.
449  *
450  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
451  */
452 enum efp_result efp_get_point_charge_gradient(struct efp *efp, double *grad);
453 
454 /**
455  * Update positions and orientations of effective fragments.
456  *
457  * \param[in] efp The efp structure.
458  *
459  * \param[in] coord_type Specifies the type of coordinates in the \a coord
460  * array (see #efp_coord_type).
461  *
462  * \param[in] coord Array of fragment coordinates.
463  *
464  * If \p coord_type is \a EFP_COORD_TYPE_XYZABC then for each fragment the \p
465  * coord array should contain \a x \a y \a z components of the center of mass
466  * position and three Euler rotation angles representing orientation of a
467  * fragment. The size of the \p coord array must be at least [6 * \a n]
468  * elements, where \a n is the number of fragments.
469  *
470  * If \p coord_type is \a EFP_COORD_TYPE_POINTS then for each fragment the \p
471  * coord array should contain the coordinates of 3 points in space. For each
472  * fragment point 1 and first atom of fragment are made to coincide. The vector
473  * connecting points 1 and 2 is aligned with the corresponding vector
474  * connecting fragment atoms. The plane defined by points 1, 2, and 3 is made
475  * to coincide with the corresponding fragment plane. The size of the \p coord
476  * array must be at least [9 * \a n] elements, where \a n is the number of
477  * fragments.
478  *
479  * If \p coord_type is \a EFP_COORD_TYPE_ROTMAT then for each fragment the \p
480  * coord array should contain \a x \a y \a z components of the center of mass
481  * position and nine elements of the rotation matrix representing orientation
482  * of a fragment. The size of the \p coord array must be at least [12 * \a n]
483  * elements, where \a n is the number of fragments.
484  *
485  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
486  */
487 enum efp_result efp_set_coordinates(struct efp *efp,
488     enum efp_coord_type coord_type, const double *coord);
489 
490 /**
491  * Update position and orientation of the specified effective fragment.
492  *
493  * \param[in] efp The efp structure.
494  *
495  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
496  * the total number of fragments minus one.
497  *
498  * \param[in] coord_type Specifies the type of coordinates in the \a coord
499  * array (see #efp_coord_type).
500  *
501  * \param[in] coord Array of coordinates specifying fragment position and
502  * orientation.
503  *
504  * If \p coord_type is \a EFP_COORD_TYPE_XYZABC then the \p coord array should
505  * contain \a x \a y \a z components of the center of mass position and three
506  * Euler rotation angles representing orientation of a fragment. The \p coord
507  * array must contain a total of 6 elements.
508  *
509  * If \p coord_type is \a EFP_COORD_TYPE_POINTS then the \p coord array should
510  * contain the coordinates of 3 points in space. Point 1 and first atom of
511  * fragment are made to coincide. The vector connecting points 1 and 2 is
512  * aligned with the corresponding vector connecting fragment atoms. The plane
513  * defined by points 1, 2, and 3 is made to coincide with the corresponding
514  * fragment plane. The \p coord array must contain a total of 9 elements.
515  *
516  * If \p coord_type is \a EFP_COORD_TYPE_ROTMAT then the \p coord array should
517  * contain \a x \a y \a z components of the center of mass position and nine
518  * elements of the rotation matrix representing orientation of a fragment. The
519  * \p coord array must contain a total of 12 elements.
520  *
521  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
522  */
523 enum efp_result efp_set_frag_coordinates(struct efp *efp, size_t frag_idx,
524     enum efp_coord_type coord_type, const double *coord);
525 
526 /**
527  * Get center of mass positions and Euler angles of the effective fragments.
528  *
529  * \param[in] efp The efp structure.
530  *
531  * \param[out] xyzabc Upon return the coordinates of the center of mass and
532  * Euler rotation angles for each fragment will be written to this array. The
533  * size of the \p xyzabc array must be at least [6 * \a n] elements, where \a n
534  * is the total number of fragments.
535  *
536  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
537  */
538 enum efp_result efp_get_coordinates(struct efp *efp, double *xyzabc);
539 
540 /**
541  * Get center of mass position and Euler angles of a fragment.
542  *
543  * \param[in] efp The efp structure.
544  *
545  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
546  * the total number of fragments minus one.
547  *
548  * \param[out] xyzabc Upon return the coordinates of the center of mass and
549  * Euler rotation angles for the fragment will be written to this array. The
550  * size of the \p xyzabc array must be at least [6] elements.
551  *
552  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
553  */
554 enum efp_result efp_get_frag_xyzabc(struct efp *efp, size_t frag_idx,
555     double *xyzabc);
556 
557 /**
558  * Setup periodic box size.
559  *
560  * \param[in] efp The efp structure.
561  * \param[in] x Box size in x dimension.
562  * \param[in] y Box size in y dimension.
563  * \param[in] z Box size in z dimension.
564  *
565  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
566  */
567 enum efp_result efp_set_periodic_box(struct efp *efp, double x, double y,
568     double z);
569 
570 /**
571  * Get periodic box size.
572  *
573  * \param[in] efp The efp structure.
574  *
575  * \param[out] xyz Array of 3 elements where 3 dimensions of box size will be
576  * stored.
577  *
578  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
579  */
580 enum efp_result efp_get_periodic_box(struct efp *efp, double *xyz);
581 
582 /**
583  * Get the stress tensor.
584  *
585  * \param[in] efp The efp structure.
586  *
587  * \param[out] stress Array of 9 elements where the stress tensor will be
588  * stored.
589  *
590  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
591  */
592 enum efp_result efp_get_stress_tensor(struct efp *efp, double *stress);
593 
594 /**
595  * Get the ab initio screening parameters.
596  *
597  * \param[in] efp The efp structure.
598  *
599  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
600  * the total number of fragments minus one.
601  *
602  * \param[out] screen Array of N elements where screening parameters will be
603  * stored. N is the total number of multipole points.
604  *
605  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
606  */
607 enum efp_result efp_get_ai_screen(struct efp *efp, size_t frag_idx,
608     double *screen);
609 
610 /**
611  * Set ab initio orbital energies.
612  *
613  * \param[in] efp The efp structure.
614  *
615  * \param[in] n_core Number of core orbitals.
616  *
617  * \param[in] n_act Number of active orbitals.
618  *
619  * \param[in] n_vir Number of virtual orbitals.
620  *
621  * \param[in] oe Array of orbital energies. The size of this array must be
622  * (n_core + n_act + n_vir) elements.
623  *
624  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
625  */
626 enum efp_result efp_set_orbital_energies(struct efp *efp, size_t n_core,
627     size_t n_act, size_t n_vir, const double *oe);
628 
629 /**
630  * Set ab initio dipole integrals.
631  *
632  * \param[in] efp The efp structure.
633  *
634  * \param[in] n_core Number of core orbitals.
635  *
636  * \param[in] n_act Number of active orbitals.
637  *
638  * \param[in] n_vir Number of virtual orbitals.
639  *
640  * \param[in] dipint Dipole integral matrices for x,y,z axes. The total size of
641  * this array must be 3 * (n_core + n_act + n_vir) ^ 2 elements.
642  *
643  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
644  */
645 enum efp_result efp_set_dipole_integrals(struct efp *efp, size_t n_core,
646     size_t n_act, size_t n_vir, const double *dipint);
647 
648 /**
649  * Update wave function dependent energy terms.
650  *
651  * This function must be called during \a ab \a initio SCF.
652  *
653  * \param[in] efp The efp structure.
654  *
655  * \param[out] energy Wave function dependent EFP energy.
656  *
657  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
658  */
659 enum efp_result efp_get_wavefunction_dependent_energy(struct efp *efp,
660     double *energy);
661 
662 /**
663  * Perform the EFP computation.
664  *
665  * \param[in] efp The efp structure.
666  *
667  * \param[in] do_gradient If nonzero value is specified in addition to energy
668  * compute the gradient.
669  *
670  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
671  */
672 enum efp_result efp_compute(struct efp *efp, int do_gradient);
673 
674 /**
675  * Get total charge of a fragment.
676  *
677  * \param[in] efp The efp structure.
678  *
679  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
680  * the total number of fragments minus one.
681  *
682  * \param[out] charge Total charge of a fragment.
683  *
684  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
685  */
686 enum efp_result efp_get_frag_charge(struct efp *efp, size_t frag_idx,
687     double *charge);
688 
689 /**
690  * Get spin multiplicity of a fragment.
691  *
692  * \param[in] efp The efp structure.
693  *
694  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
695  * the total number of fragments minus one.
696  *
697  * \param[out] mult Spin multiplicity of a fragment.
698  *
699  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
700  */
701 enum efp_result efp_get_frag_multiplicity(struct efp *efp, size_t frag_idx,
702     int *mult);
703 
704 /**
705  * Get number of electrostatic multipole points for a particular fragment.
706  *
707  * \param[in] efp The efp structure.
708  *
709  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
710  * the total number of fragments minus one.
711  *
712  * \param[out] n_mult Number of electrostatic multipole points.
713  *
714  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
715  */
716 enum efp_result efp_get_frag_multipole_count(struct efp *efp, size_t frag_idx,
717     size_t *n_mult);
718 
719 /**
720  * Get total number of multipoles from EFP electrostatics.
721  *
722  * \param[in] efp The efp structure.
723  *
724  * \param[out] n_mult Number of electrostatics multipoles.
725  *
726  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
727  */
728 enum efp_result efp_get_multipole_count(struct efp *efp, size_t *n_mult);
729 
730 /**
731  * Get coordinates of electrostatics multipoles.
732  *
733  * \param[in] efp The efp structure.
734  *
735  * \param[out] xyz Array where coordinates of EFP electrostatics multipoles
736  * will be stored. Size of the \p xyz array must be at least [3 * \p n_mult]
737  * elements, where \p n_mult is the value returned by the
738  * ::efp_get_multipole_count function.
739  *
740  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
741  */
742 enum efp_result efp_get_multipole_coordinates(struct efp *efp, double *xyz);
743 
744 /**
745  * Get electrostatics multipoles from EFP fragments.
746  *
747  * \param[in] efp The efp structure.
748  *
749  * \param[out] mult Array where charges, dipoles, quadrupoles, and octupoles
750  * for each point will be stored.
751  *
752  * The size of the \p mult array must be at least [(1 + 3 + 6 + 10) * \p
753  * n_mult] elements (charges + dipoles + quadrupoles + octupoles), where \p
754  * n_mult is the value returned by the ::efp_get_multipole_count function.
755  *
756  * Quadrupoles are stored in the following order:
757  *    \a xx, \a yy, \a zz, \a xy, \a xz, \a yz
758  *
759  * Octupoles are stored in the following order:
760  *    \a xxx, \a yyy, \a zzz, \a xxy, \a xxz,
761  *    \a xyy, \a yyz, \a xzz, \a yzz, \a xyz
762  *
763  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
764  */
765 enum efp_result efp_get_multipole_values(struct efp *efp, double *mult);
766 
767 /**
768  *  Get the number of polarization induced dipoles.
769  *
770  * \param[in] efp The efp structure.
771  *
772  * \param[out] n_dip Number of polarization induced dipoles.
773  *
774  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
775  */
776 enum efp_result efp_get_induced_dipole_count(struct efp *efp, size_t *n_dip);
777 
778 /**
779  * Get coordinates of induced dipoles.
780  *
781  * \param[in] efp The efp structure.
782  *
783  * \param[out] xyz Array where the coordinates of polarizable points will be
784  * stored. The size of the \p xyz array must be at least [3 * \p n_dip]
785  * elements.
786  *
787  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
788  */
789 enum efp_result efp_get_induced_dipole_coordinates(struct efp *efp,
790     double *xyz);
791 
792 /**
793  * Get values of polarization induced dipoles.
794  *
795  * \param[in] efp The efp structure.
796  *
797  * \param[out] dip Array where induced dipoles will be stored. The size of the
798  * array must be at least [3 * \p n_dip] elements.
799  *
800  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
801  */
802 enum efp_result efp_get_induced_dipole_values(struct efp *efp, double *dip);
803 
804 /**
805  * Get values of polarization conjugated induced dipoles.
806  *
807  * \param[in] efp The efp structure.
808  *
809  * \param[out] dip Array where induced dipoles will be stored. The size of the
810  * array must be at least [3 * \p n_dip] elements.
811  *
812  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
813  */
814 enum efp_result efp_get_induced_dipole_conj_values(struct efp *efp,
815     double *dip);
816 
817 /**
818  * Get the number of LMOs in a fragment.
819  *
820  * \param[in] efp The efp structure.
821  *
822  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
823  * the total number of fragments minus one.
824  *
825  * \param[out] n_lmo Number of LMOs.
826  *
827  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
828  */
829 enum efp_result efp_get_lmo_count(struct efp *efp, size_t frag_idx,
830     size_t *n_lmo);
831 
832 /**
833  * Get coordinates of LMO centroids.
834  *
835  * \param[in] efp The efp structure.
836  *
837  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
838  * the total number of fragments minus one.
839  *
840  * \param[out] xyz Array where the coordinates of LMO centroids will be
841  * stored. The size of the \p xyz array must be at least [3 * \p n_lmo]
842  * elements.
843  *
844  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
845  */
846 enum efp_result efp_get_lmo_coordinates(struct efp *efp, size_t frag_idx,
847     double *xyz);
848 
849 /**
850  * Get parameters of fitted exchange-repulsion.
851  *
852  * \param[in] efp The efp structure.
853  *
854  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
855  * the total number of fragments minus one.
856  *
857  * \param[out] xrfit Array where the parameters will be stored. The size of the
858  * \p xrfit array must be at least [4 * \p n_lmo] elements.
859  *
860  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
861  */
862 enum efp_result efp_get_xrfit(struct efp *efp, size_t frag_idx, double *xrfit);
863 
864 /**
865  * Get computed energy components.
866  *
867  * \param[in] efp The efp structure.
868  * \param[out] energy Computed EFP energy components (see efp_energy).
869  *
870  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
871  */
872 enum efp_result efp_get_energy(struct efp *efp, struct efp_energy *energy);
873 
874 /**
875  * Get computed EFP energy gradient.
876  *
877  * \param[in] efp The efp structure.
878  *
879  * \param[out] grad For each fragment \a x \a y \a z components of negative
880  * force and torque will be written to this array. The size of this array must
881  * be at least [6 * \a n] elements, where \a n is the total number of
882  * fragments.
883  *
884  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
885  */
886 enum efp_result efp_get_gradient(struct efp *efp, double *grad);
887 
888 /**
889  * Get computed EFP energy gradient on individual atoms.
890  *
891  * \param[in] efp The efp structure.
892  *
893  * \param[out] grad For each atom, \a x \a y \a z components of negative force
894  * will be added to this array. The size of this array must be at least
895  * [3 * \a n] elements, where \a n is the total number of atoms from all
896  * fragments. An atom is a point with non-zero mass inside a fragment.
897  * Any initial gradient from this array will be gathered on fragments at the
898  * beginning and then redistributed back to the atoms. This can be used to
899  * account for other interactions, e.g., bonded forces from MM forcefield.
900  *
901  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
902  */
903 enum efp_result efp_get_atomic_gradient(struct efp *efp, double *grad);
904 
905 /**
906  * Get the number of fragments in this computation.
907  *
908  * \param[in] efp The efp structure.
909  *
910  * \param[out] n_frag Total number of fragments in this simulation.
911  *
912  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
913  */
914 enum efp_result efp_get_frag_count(struct efp *efp, size_t *n_frag);
915 
916 /**
917  * Get the name of the specified effective fragment.
918  *
919  * \param[in] efp The efp structure.
920  *
921  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
922  * the total number of fragments minus one.
923  *
924  * \param[in] size Size of a \p frag_name buffer.
925  *
926  * \param[out] frag_name A buffer where name of the fragment will be stored.
927  *
928  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
929  */
930 enum efp_result efp_get_frag_name(struct efp *efp, size_t frag_idx, size_t size,
931     char *frag_name);
932 
933 /**
934  * Get total mass of a fragment.
935  *
936  * \param[in] efp The efp structure.
937  *
938  * \param[in] frag_idx Index of a fragment. This must be a value between zero
939  * and the total number of fragments minus one.
940  *
941  * \param[out] mass Output mass value.
942  *
943  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
944  */
945 enum efp_result efp_get_frag_mass(struct efp *efp, size_t frag_idx,
946     double *mass);
947 
948 /**
949  * Get fragment principal moments of inertia.
950  *
951  * \param[in] efp The efp structure.
952  *
953  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
954  * the total number of fragments minus one.
955  *
956  * \param[out] inertia Array of 3 elements where principal moments of
957  * inertia of a fragment will be stored.
958  *
959  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
960  */
961 enum efp_result efp_get_frag_inertia(struct efp *efp, size_t frag_idx,
962     double *inertia);
963 
964 /**
965  * Get the number of atoms in the specified fragment.
966  *
967  * \param[in] efp The efp structure.
968  *
969  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
970  * the total number of fragments minus one.
971  *
972  * \param[out] n_atoms Total number of atoms in the fragment.
973  *
974  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
975  */
976 enum efp_result efp_get_frag_atom_count(struct efp *efp, size_t frag_idx,
977     size_t *n_atoms);
978 
979 /**
980  * Get atoms comprising the specified fragment.
981  *
982  * \param[in] efp The efp structure.
983  *
984  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
985  * the total number of fragments minus one.
986  *
987  * \param[in] size Size of the \p atoms array. Must be greater than or equal to
988  * the size returned by the ::efp_get_frag_atom_count function.
989  *
990  * \param[out] atoms Array where atom information will be stored.
991  *
992  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
993  */
994 enum efp_result efp_get_frag_atoms(struct efp *efp, size_t frag_idx,
995     size_t size, struct efp_atom *atoms);
996 
997 /**
998  * Get electric field for a point on a fragment.
999  *
1000  * \param[in] efp The efp structure.
1001  *
1002  * \param[in] frag_idx Index of a fragment. Must be a value between zero and
1003  * the total number of fragments minus one.
1004  *
1005  * \param[in] xyz Coordinates of a point for which electric field should be
1006  * computed.
1007  *
1008  * \param[out] field Electric field \a x \a y \a z components in atomic units.
1009  *
1010  * \return ::EFP_RESULT_SUCCESS on success or error code otherwise.
1011  */
1012 enum efp_result efp_get_electric_field(struct efp *efp, size_t frag_idx,
1013     const double *xyz, double *field);
1014 
1015 /**
1016  * Convert rigid body torque to derivatives of energy by Euler angles.
1017  *
1018  * \param[in] euler Array of 3 elements containing the values of Euler angles
1019  * for the current orientation of the rigid body.
1020  *
1021  * \param[in] torque Array of 3 elements containing torque components.
1022  *
1023  * \param[out] deriv Array of 3 elements where derivatives of energy by
1024  * Euler angles will be stored. This can point to the same memory location as
1025  * the \p torque argument.
1026  */
1027 void efp_torque_to_derivative(const double *euler, const double *torque,
1028     double *deriv);
1029 
1030 /**
1031  * Release all resources used by this \a efp.
1032  *
1033  * \param[in] efp The efp structure to be released.
1034  */
1035 void efp_shutdown(struct efp *efp);
1036 
1037 /**
1038  * Convert #efp_result to a human readable message.
1039  *
1040  * \param res Result value to be converted to string.
1041  *
1042  * \return Human readable string with the description of the result,
1043  * zero-terminated.
1044  */
1045 const char *efp_result_to_string(enum efp_result res);
1046 
1047 #ifdef __cplusplus
1048 } /* extern "C" */
1049 #endif
1050 
1051 #endif /* LIBEFP_EFP_H */
1052