1! This file is part of xtb.
2!
3! Copyright (C) 2019-2020 Sebastian Ehlert
4!
5! xtb is free software: you can redistribute it and/or modify it under
6! the terms of the GNU Lesser General Public License as published by
7! the Free Software Foundation, either version 3 of the License, or
8! (at your option) any later version.
9!
10! xtb is distributed in the hope that it will be useful,
11! but WITHOUT ANY WARRANTY; without even the implied warranty of
12! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
13! GNU Lesser General Public License for more details.
14!
15! You should have received a copy of the GNU Lesser General Public License
16! along with xtb.  If not, see <https://www.gnu.org/licenses/>.
17
18!> API for dealing with the calculation results
19module xtb_api_results
20   use, intrinsic :: iso_c_binding
21   use xtb_mctc_accuracy, only : wp
22   use xtb_mctc_convert, only : evtoau
23   use xtb_api_environment
24   use xtb_api_utils
25   use xtb_type_restart
26   implicit none
27   private
28
29   public :: VResults
30   public :: newResults_api, delResults_api, getEnergy_api, getGradient_api, &
31      & getVirial_api, getDipole_api, getCharges_api, getBondOrders_api, &
32      & getNao_api, getOrbitalEigenvalues_api, getOrbitalOccupations_api, &
33      & getOrbitalCoefficients_api, copyResults_api
34
35
36   !> Void pointer to wavefunction and result objects
37   type :: VResults
38      type(TRestart), allocatable :: chk
39      real(wp), allocatable :: energy
40      real(wp), allocatable :: gradient(:, :)
41      real(wp), allocatable :: sigma(:, :)
42      real(wp), allocatable :: dipole(:)
43      real(wp), allocatable :: egap
44   end type VResults
45
46
47contains
48
49
50!> Create new singlepoint results object
51function newResults_api() result(vres) &
52      & bind(C, name="xtb_newResults")
53   type(VResults), pointer :: res
54   type(c_ptr) :: vres
55
56   call checkGlobalEnv
57
58   allocate(res)
59   vres = c_loc(res)
60
61end function newResults_api
62
63
64!> Delete singlepoint results object
65subroutine delResults_api(vres) &
66      & bind(C, name="xtb_delResults")
67   type(c_ptr), intent(inout) :: vres
68   type(VResults), pointer :: res
69
70   call checkGlobalEnv
71
72   if (c_associated(vres)) then
73      call c_f_pointer(vres, res)
74      deallocate(res)
75      vres = c_null_ptr
76   end if
77
78end subroutine delResults_api
79
80
81!> Create copy from a singlepoint results object
82function copyResults_api(vold) result(vres) &
83      & bind(C, name="xtb_copyResults")
84   type(VResults), pointer :: res
85   type(c_ptr) :: vres
86   type(VResults), pointer :: old
87   type(c_ptr), value :: vold
88
89   vres = c_null_ptr
90
91   if (c_associated(vold)) then
92      call c_f_pointer(vold, old)
93      call checkGlobalEnv
94
95      allocate(res, source=old)
96      vres = c_loc(res)
97   end if
98
99end function copyResults_api
100
101
102!> Query singlepoint results object for energy
103subroutine getEnergy_api(venv, vres, dptr) &
104      & bind(C, name="xtb_getEnergy")
105   character(len=*), parameter :: source = "xtb_api_getEnergy"
106   type(c_ptr), value :: venv
107   type(VEnvironment), pointer :: env
108   type(c_ptr), value :: vres
109   type(VResults), pointer :: res
110   real(c_double), intent(inout) :: dptr
111
112   if (c_associated(venv)) then
113      call c_f_pointer(venv, env)
114      call checkGlobalEnv
115
116      if (.not.c_associated(vres)) then
117         call env%ptr%error("Results object is not allocated", source)
118         return
119      end if
120      call c_f_pointer(vres, res)
121
122      if (.not.allocated(res%energy)) then
123         call env%ptr%error("Energy is not available in results", source)
124         return
125      end if
126
127      dptr = res%energy
128
129   end if
130
131end subroutine getEnergy_api
132
133
134!> Query singlepoint results object for gradient
135subroutine getGradient_api(venv, vres, dptr) &
136      & bind(C, name="xtb_getGradient")
137   character(len=*), parameter :: source = "xtb_api_getGradient"
138   type(c_ptr), value :: venv
139   type(VEnvironment), pointer :: env
140   type(c_ptr), value :: vres
141   type(VResults), pointer :: res
142   real(c_double), intent(inout) :: dptr(3, *)
143
144   if (c_associated(venv)) then
145      call c_f_pointer(venv, env)
146      call checkGlobalEnv
147
148      if (.not.c_associated(vres)) then
149         call env%ptr%error("Results object is not allocated", source)
150         return
151      end if
152      call c_f_pointer(vres, res)
153
154      if (.not.allocated(res%gradient)) then
155         call env%ptr%error("Gradient is not available in results", source)
156         return
157      end if
158
159      dptr(1:3, 1:size(res%gradient, 2)) = res%gradient
160
161   end if
162
163end subroutine getGradient_api
164
165
166!> Query singlepoint results object for virial
167subroutine getVirial_api(venv, vres, dptr) &
168      & bind(C, name="xtb_getVirial")
169   character(len=*), parameter :: source = "xtb_api_getVirial"
170   type(c_ptr), value :: venv
171   type(VEnvironment), pointer :: env
172   type(c_ptr), value :: vres
173   type(VResults), pointer :: res
174   real(c_double), intent(inout) :: dptr(3, 3)
175
176   if (c_associated(venv)) then
177      call c_f_pointer(venv, env)
178      call checkGlobalEnv
179
180      if (.not.c_associated(vres)) then
181         call env%ptr%error("Results object is not allocated", source)
182         return
183      end if
184      call c_f_pointer(vres, res)
185
186      if (.not.allocated(res%sigma)) then
187         call env%ptr%error("Virial is not available in results", source)
188         return
189      end if
190
191      dptr(1:3, 1:3) = res%sigma
192
193   end if
194
195end subroutine getVirial_api
196
197
198!> Query singlepoint results object for dipole moment
199subroutine getDipole_api(venv, vres, dptr) &
200      & bind(C, name="xtb_getDipole")
201   character(len=*), parameter :: source = "xtb_api_getDipole"
202   type(c_ptr), value :: venv
203   type(VEnvironment), pointer :: env
204   type(c_ptr), value :: vres
205   type(VResults), pointer :: res
206   real(c_double), intent(inout) :: dptr(3)
207
208   if (c_associated(venv)) then
209      call c_f_pointer(venv, env)
210      call checkGlobalEnv
211
212      if (.not.c_associated(vres)) then
213         call env%ptr%error("Results object is not allocated", source)
214         return
215      end if
216      call c_f_pointer(vres, res)
217
218      if (.not.allocated(res%dipole)) then
219         call env%ptr%error("Dipole moment is not available in results", source)
220         return
221      end if
222
223      dptr(1:3) = res%dipole
224
225   end if
226
227end subroutine getDipole_api
228
229
230!> Query singlepoint results object for partial charges
231subroutine getCharges_api(venv, vres, dptr) &
232      & bind(C, name="xtb_getCharges")
233   character(len=*), parameter :: source = "xtb_api_getCharges"
234   type(c_ptr), value :: venv
235   type(VEnvironment), pointer :: env
236   type(c_ptr), value :: vres
237   type(VResults), pointer :: res
238   real(c_double), intent(inout) :: dptr(*)
239
240   if (c_associated(venv)) then
241      call c_f_pointer(venv, env)
242      call checkGlobalEnv
243
244      if (.not.c_associated(vres)) then
245         call env%ptr%error("Results object is not allocated", source)
246         return
247      end if
248      call c_f_pointer(vres, res)
249
250      if (.not.allocated(res%chk)) then
251         call env%ptr%error("Partial charges are not available in results", source)
252         return
253      end if
254
255      dptr(1:size(res%chk%wfn%q)) = res%chk%wfn%q
256
257   end if
258
259end subroutine getCharges_api
260
261
262!> Query singlepoint results object for bond orders
263subroutine getBondOrders_api(venv, vres, dptr) &
264      & bind(C, name="xtb_getBondOrders")
265   character(len=*), parameter :: source = "xtb_api_getBondOrders"
266   type(c_ptr), value :: venv
267   type(VEnvironment), pointer :: env
268   type(c_ptr), value :: vres
269   type(VResults), pointer :: res
270   real(c_double), intent(inout) :: dptr(*)
271
272   if (c_associated(venv)) then
273      call c_f_pointer(venv, env)
274      call checkGlobalEnv
275
276      if (.not.c_associated(vres)) then
277         call env%ptr%error("Results object is not allocated", source)
278         return
279      end if
280      call c_f_pointer(vres, res)
281
282      if (.not.allocated(res%chk)) then
283         call env%ptr%error("Bond orders are not available in results", source)
284         return
285      end if
286
287      dptr(1:size(res%chk%wfn%wbo)) = reshape(res%chk%wfn%wbo, [size(res%chk%wfn%wbo)])
288
289   end if
290
291end subroutine getBondOrders_api
292
293
294!> Query singlepoint results object for bond orders
295subroutine getNao_api(venv, vres, iptr) &
296      & bind(C, name="xtb_getNao")
297   character(len=*), parameter :: source = "xtb_api_getNao"
298   type(c_ptr), value :: venv
299   type(VEnvironment), pointer :: env
300   type(c_ptr), value :: vres
301   type(VResults), pointer :: res
302   integer(c_int), intent(inout) :: iptr
303
304   if (c_associated(venv)) then
305      call c_f_pointer(venv, env)
306      call checkGlobalEnv
307
308      if (.not.c_associated(vres)) then
309         call env%ptr%error("Results object is not allocated", source)
310         return
311      end if
312      call c_f_pointer(vres, res)
313
314      if (allocated(res%chk)) then
315         iptr = res%chk%wfn%nao
316      else
317         iptr = 0
318      end if
319
320   end if
321
322end subroutine getNao_api
323
324
325!> Query singlepoint results object for orbital energies
326subroutine getOrbitalEigenvalues_api(venv, vres, dptr) &
327      & bind(C, name="xtb_getOrbitalEigenvalues")
328   character(len=*), parameter :: source = "xtb_api_getOrbitalEigenvalues"
329   type(c_ptr), value :: venv
330   type(VEnvironment), pointer :: env
331   type(c_ptr), value :: vres
332   type(VResults), pointer :: res
333   real(c_double), intent(inout) :: dptr(*)
334
335   if (c_associated(venv)) then
336      call c_f_pointer(venv, env)
337      call checkGlobalEnv
338
339      if (.not.c_associated(vres)) then
340         call env%ptr%error("Results object is not allocated", source)
341         return
342      end if
343      call c_f_pointer(vres, res)
344
345      if (.not.allocated(res%chk)) then
346         call env%ptr%error("Orbital eigenvalues are not available in results", &
347            & source)
348         return
349      end if
350
351      dptr(1:size(res%chk%wfn%emo)) = res%chk%wfn%emo * evtoau
352
353   end if
354
355end subroutine getOrbitalEigenvalues_api
356
357
358!> Query singlepoint results object for occupation numbers
359subroutine getOrbitalOccupations_api(venv, vres, dptr) &
360      & bind(C, name="xtb_getOrbitalOccupations")
361   character(len=*), parameter :: source = "xtb_api_getOrbitalOccupations"
362   type(c_ptr), value :: venv
363   type(VEnvironment), pointer :: env
364   type(c_ptr), value :: vres
365   type(VResults), pointer :: res
366   real(c_double), intent(inout) :: dptr(*)
367
368   if (c_associated(venv)) then
369      call c_f_pointer(venv, env)
370      call checkGlobalEnv
371
372      if (.not.c_associated(vres)) then
373         call env%ptr%error("Results object is not allocated", source)
374         return
375      end if
376      call c_f_pointer(vres, res)
377
378      if (.not.allocated(res%chk)) then
379         call env%ptr%error("Occupation numbers are not available in results", &
380            & source)
381         return
382      end if
383
384      dptr(1:size(res%chk%wfn%focc)) = res%chk%wfn%focc
385
386   end if
387
388end subroutine getOrbitalOccupations_api
389
390
391!> Query singlepoint results object for orbital coefficients
392subroutine getOrbitalCoefficients_api(venv, vres, dptr) &
393      & bind(C, name="xtb_getOrbitalCoefficients")
394   character(len=*), parameter :: source = "xtb_api_getOrbitalCoefficients"
395   type(c_ptr), value :: venv
396   type(VEnvironment), pointer :: env
397   type(c_ptr), value :: vres
398   type(VResults), pointer :: res
399   real(c_double), intent(inout) :: dptr(*)
400
401   if (c_associated(venv)) then
402      call c_f_pointer(venv, env)
403      call checkGlobalEnv
404
405      if (.not.c_associated(vres)) then
406         call env%ptr%error("Results object is not allocated", source)
407         return
408      end if
409      call c_f_pointer(vres, res)
410
411      if (.not.allocated(res%chk)) then
412         call env%ptr%error("Orbital coefficients are not available in results", &
413            & source)
414         return
415      end if
416
417      dptr(1:size(res%chk%wfn%C)) = reshape(res%chk%wfn%C, [size(res%chk%wfn%C)])
418
419   end if
420
421end subroutine getOrbitalCoefficients_api
422
423
424end module xtb_api_results
425