Initialize the non-bonded object allocating the parameters vectors
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(ommp_nonbonded_type), | intent(inout) | :: | vdw | |||
type(ommp_topology_type), | intent(in), | target | :: | top | ||
character(len=*) | :: | vdw_type | ||||
character(len=*) | :: | radius_rule | ||||
character(len=*) | :: | radius_size | ||||
character(len=*) | :: | radius_type | ||||
character(len=*) | :: | epsrule | ||||
real(kind=rp) | :: | cutoff |
subroutine vdw_init(vdw, top, vdw_type, radius_rule, radius_size, &
radius_type, epsrule, cutoff)
!! Initialize the non-bonded object allocating the parameters vectors
use mod_memory, only: mallocate
use mod_io, only: fatal_error
use mod_neighbor_list, only: nl_init
use mod_constants, only: OMMP_DEFAULT_NL_SUB
implicit none
type(ommp_nonbonded_type), intent(inout) :: vdw
type(ommp_topology_type), intent(in), target :: top
character(len=*) :: vdw_type, radius_rule, radius_size, radius_type, &
epsrule
real(rp) :: cutoff
select case(trim(vdw_type))
case("lennard-jones")
vdw%vdwtype = OMMP_VDWTYPE_LJ
continue
case("buckingham")
call fatal_error("VdW type buckingham is not implemented")
case("buffered-14-7")
vdw%vdwtype = OMMP_VDWTYPE_BUF714
continue
case("mm3-hbond")
call fatal_error("VdW type mm3-hbond is not implemented")
case("gaussian")
call fatal_error("VdW type gaussian is not implemented")
case default
call fatal_error("VdW type specified is not understood")
end select
select case(trim(radius_rule))
case("arithmetic")
vdw%radrule = OMMP_RADRULE_ARITHMETIC
continue
case("geometric")
call fatal_error("radiusrule geometric is not implemented")
case("cubic-mean")
vdw%radrule = OMMP_RADRULE_CUBIC
continue
case default
call fatal_error("radiusrule specified is not understood")
end select
select case(trim(radius_size))
case("radius")
vdw%radf = 2.0
continue
case("diameter")
vdw%radf = 1.0
continue
case default
call fatal_error("radiussize specified is not understood")
end select
select case(trim(radius_type))
case("sigma")
call fatal_error("radiustype sigma is not implemented")
case("r-min")
vdw%radtype = OMMP_RADTYPE_RMIN
continue
case default
call fatal_error("radiustype specified is not understood")
end select
select case(trim(epsrule))
case("geometric")
vdw%epsrule = OMMP_EPSRULE_GEOMETRIC
continue
case("arithmetic")
call fatal_error("epsilonrule arithmetic is not implemented")
case("harmonic")
call fatal_error("epsilonrule harmonic is not implemented")
case("w-h")
call fatal_error("epsilonrule w-h is not implemented")
case("hhg")
vdw%epsrule = OMMP_EPSRULE_HHG
continue
case default
call fatal_error("epsilonrule specified is not understood")
end select
vdw%top => top
call mallocate('vdw_init [vdw_r]', top%mm_atoms, vdw%vdw_r)
call mallocate('vdw_init [vdw_e]', top%mm_atoms, vdw%vdw_e)
call mallocate('vdw_init [vdw_f]', top%mm_atoms, vdw%vdw_f)
call mallocate('vdw_init [vdw_pair_mask_a]', &
top%mm_atoms, pair_allocation_chunk,vdw%vdw_pair_mask_a)
call mallocate('vdw_init [vdw_pair_mask_b]', &
top%mm_atoms, pair_allocation_chunk,vdw%vdw_pair_mask_b)
call mallocate('vdw_init [vdw_pair_r]', pair_allocation_chunk, &
vdw%vdw_pair_r)
call mallocate('vdw_init [vdw_pair_e]', pair_allocation_chunk, &
vdw%vdw_pair_e)
vdw%vdw_f = 1.0_rp
if(cutoff > 0.0) then
vdw%use_nl = .true.
if(vdw%use_nl) call nl_init(vdw%nl, top%cmm, cutoff, OMMP_DEFAULT_NL_SUB)
else
vdw%use_nl = .false.
end if
end subroutine vdw_init