-
Notifications
You must be signed in to change notification settings - Fork 17
/
mkrhs.f90
123 lines (123 loc) · 5.1 KB
/
mkrhs.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
subroutine mkrhs(n,charge,x,y,z,ncav,ccav,phi,nylm,psi)
implicit none
!
! 888 888 .d8888b. .d88888b. .d8888b. 888b d888 .d88888b.
! 888 888 d88P Y88b d88P" "Y88b d88P Y88b 8888b d8888 d88P" "Y88b
! 888 888 888 888 888 888 Y88b. 88888b.d88888 888 888
! .d88888 .d88888 888 888 888 "Y888b. 888Y88888P888 888 888
! d88" 888 d88" 888 888 888 888 "Y88b. 888 Y888P 888 888 888
! 888 888 888 888 888 888 888 888 "888 888 Y8P 888 888 888
! Y88b 888 Y88b 888 Y88b d88P Y88b. .d88P Y88b d88P 888 " 888 Y88b. .d88P
! "Y88888 "Y88888 "Y8888P" "Y88888P" "Y8888P" 888 888 "Y88888P"
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! COPYRIGHT (C) 2015 by Filippo Lipparini, Benjamin Stamm, Paolo Gatto !
! Eric Cancès, Yvon Maday, Jean-Philip Piquemal, Louis Lagardère and !
! Benedetta Mennucci. !
! ALL RIGHT RESERVED. !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!
! A modular implementation of COSMO using a domain decomposition linear scaling
! strategy.
!
! This code is governed by the LGPL license and abiding by the rules of
! distribution of free software.
! This program is distributed in the hope that it will be useful, but
! WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
! or FITNESS FOR A PARTICULAR PURPOSE.
! See the GNU Lesser General Public License for more details.
!
! Users of this code are asked to include the following references in their
! publications:
!
! [1] E. Cancès, Y. Maday, B. Stamm
! "Domain decomposition for implicit solvation models"
! J. Chem. Phys. 139, 054111 (2013)
!
! [2] F. Lipparini, B. Stamm, E. Cancès, Y. Maday, B. Mennucci
! "Fast Domain Decomposition Algorithm for Continuum Solvation Models:
! Energy and First Derivatives"
! J. Chem. Theory Comput. 9, 3637–3648 (2013)
!
! Also, include one of the three following reference depending on whether you
! use this code in conjunction with a QM [3], Semiempirical [4] or Classical [5]
! description of the solute:
!
! [3] F. Lipparini, G. Scalmani, L. Lagardère, B. Stamm, E. Cancès, Y. Maday,
! J.-P. Piquemal, M. J. Frisch, B. Mennucci
! "Quantum, classical, and hybrid QM/MM calculations in solution: General
! implementation of the ddCOSMO linear scaling strategy"
! J. Chem. Phys. 141, 184108 (2014)
! (for quantum mechanical models)
!
! [4] F. Lipparini, L. Lagardère, G. Scalmani, B. Stamm, E. Cancès, Y. Maday,
! J.-P. Piquemal, M. J. Frisch, B. Mennucci
! "Quantum Calculations in Solution for Large to Very Large Molecules:
! A New Linear Scaling QM/Continuum Approach"
! J. Phys. Chem. Lett. 5, 953-958 (2014)
! (for semiempirical models)
!
! [5] F. Lipparini, L. Lagardère, C. Raynaud, B. Stamm, E. Cancès, B. Mennucci
! M. Schnieders, P. Ren, Y. Maday, J.-P. Piquemal
! "Polarizable Molecular Dynamics in a Polarizable Continuum Solvent"
! J. Chem. Theory Comput. 11, 623-634 (2015)
! (for classical models, including polarizable force fields
!
! The users of this code should also include the appropriate reference to the
! COSMO model. This distribution includes the routines to generate lebedev
! grids by D. Laikov and C. van Wuellen, as publicly available on CCL. If the routines
! are used, the following reference should also be included:
!
! [6] V.I. Lebedev, and D.N. Laikov
! "A quadrature formula for the sphere of the 131st
! algebraic order of accuracy"
! Doklady Mathematics, Vol. 59, No. 3, 1999, pp. 477-481.
!
! Written by Filippo Lipparini, October 2015.
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! !
! Silly routine to compute the potential and psi vector !
! !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
integer, intent(in) :: n, ncav, nylm
real*8, dimension(n), intent(in) :: x, y, z, charge
real*8, dimension(3,ncav), intent(in) :: ccav
real*8, dimension(ncav), intent(inout) :: phi
real*8, dimension(nylm,n), intent(inout) :: psi
!
integer :: isph, ic, j
real*8 :: v
real*8 :: dx, dy, dz, d2, d, pi, fac
real*8, parameter :: zero=0.0d0, one=1.0d0, four=4.0d0
!
pi = four*atan(one)
fac = sqrt(four*pi)
phi = zero
psi = zero
!
!$omp parallel do default(shared) private(ic,v,j,dx,dy,dz,d2,d)
do ic = 1, ncav
v = zero
do j = 1, n
dx = ccav(1,ic) - x(j)
dy = ccav(2,ic) - y(j)
dz = ccav(3,ic) - z(j)
d2 = dx*dx + dy*dy + dz*dz
d = sqrt(d2)
v = v + charge(j)/d
end do
phi(ic) = v
end do
!
! psi vector:
!
do isph = 1, n
psi(1,isph) = fac*charge(isph)
end do
!
return
end