packmol
packmol copied to clipboard
m3g
XYgauss surface
Added XYgauss function as constraint.
Works by taking the plane parsing of getinp, using over and below with 6 extra arguments:
over (or below) xygauss a b c d g h
based on h*exp(-(x-a)**2/2c**2 -(y-b)**2/2d**2)-(z-g)
To implement a new restraint type, there must an addition to the comprest.f90 file, where the function value relative to the constraints are computed. Also, after implementing the function and the gradient, run some examples with the flag
chkgrad
To check if the gradient is correct (this flag compares analytical and finite-difference gradients)
Hi Leandro, For the record, this function is working, and the gradient was checked as you requested. Best regards,
I would need a description, in the form of the user guide, to understand exactly what is this constraint, and some examples to test it.
A description, using similar variable names as the ones you have in webpage would be:
over (or below) xygauss
Gaussian planes are defined by the equation h exp(−(x-a₁)²/(2a₂²)−(y-b₁)²/(2b₂²))-(z-c)=0. Coordinates (a₁, b₁) define the x,y position, while c specifies the position in the z dimension. (a₂,b₂) set the width of the gaussian in x and y, respectively, while h specifies its height. As for conventional planes. it is possible to restrict atoms to be over or below the gaussian plane, although in this case it is restricted to be in the xy plane. The syntax is
over xygauss a₁ b₁ a₂ b₂ c h
Examples:
over xygauss 21.0 5.0 10.0 20.0 -23.0 15.0 below xygauss 21.0 5.0 10.0 20.0 -23.0 15.0
where the over keyword will make the atoms satisfy the condition
15.0 exp(− (x-21.0)²/(2(10.0)²) − (y-5.0)²/(2(20.0)²))-(z+23) >= 0
the below keyword will make the atoms satisfy
15.0 exp(− (x-21.0)²/(2(10.0)²) − (y-5.0)²/(2(20.0)²))-(z+23) <= 0
Let me know if that is ok or not.
Since there are too many changes since the first PR, I will probably try to merge by hand here.
Ok. I just did a pull of the upstream, and it seems to be fine here if that help.