10 #ifndef IMPNPCTRANSPORT_SLAB_WITH_TOROIDAL_PORE_PAIR_SCORE_H
11 #define IMPNPCTRANSPORT_SLAB_WITH_TOROIDAL_PORE_PAIR_SCORE_H
13 #include "npctransport_config.h"
22 IMPNPCTRANSPORT_BEGIN_NAMESPACE
32 mutable double bottom_;
37 mutable bool is_pore_radius_optimized_;
69 unsigned int lower_bound,
70 unsigned int upper_bound)
const IMP_FINAL;
86 unsigned int lower_bound,
87 unsigned int upper_bound)
const
90 for (
unsigned int i = lower_bound; i < upper_bound; ++i) {
92 if (ret > max)
return std::numeric_limits<double>::max();
119 get_sphere_ellipsoid_penetration_depth
134 inline double get_sphere_penetration_depth
140 void update_cached_slab_params
145 SlabWithToroidalPorePairScore::update_cached_slab_params
151 bottom_= -0.5*thickness;
154 rh_= slab.get_horizontal_minor_radius();
155 rv_= slab.get_vertical_minor_radius();
156 is_pore_radius_optimized_= slab.get_pore_radius_is_optimized();
162 SlabWithToroidalPorePairScore::evaluate_index
169 update_cached_slab_params(slab);
175 double score=get_sphere_penetration_depth(d.
get_sphere(),
176 da ? &displacement :
nullptr);
182 if(is_pore_radius_optimized_){
183 double radial_displacement= displacement[0]*d_sphere[0]+displacement[1]*d_sphere[1];
193 SlabWithToroidalPorePairScore::evaluate_indexes
197 unsigned int lower_bound,
198 unsigned int upper_bound)
const
200 IMP_LOG_TERSE(
"SlabWithToroidalPore singleton - evaluate indexes"
202 if(upper_bound<lower_bound){
206 double radial_displacements(0.0);
208 m->access_spheres_data();
210 m->access_sphere_derivatives_data();
211 IMP::internal::BoolAttributeTableTraits::Container
const&
212 is_optimizable_table= m->access_optimizeds_data
213 (core::XYZ::get_coordinate_key(0));
216 update_cached_slab_params(slab);
219 for (
unsigned int i = lower_bound; i < upper_bound; ++i) {
221 int pi_index=pi.get_index();
225 "All particles are assumed to be evaluated against"
230 "optimable table inconsistent with d.get_coordinates_are_optimized for particle " << d);
232 "Different coords for particle " << d <<
" *** "
235 "Different radii for particle " << d <<
" *** "
236 << d.get_radius() <<
" vs. " << s.get_radius());
238 if(!is_optimizable_table[pi]) {
242 double cur_score = get_sphere_penetration_depth(spheres_table[pi_index],
243 da ? &displacement :
nullptr);
244 IMP_LOG_TERSE(
"SlabWithToroidalPore singleton score for sphere / displacement "
245 << spheres_table[pi_index] <<
" is " << cur_score <<
" / "
246 << displacement << std::endl);
248 if(cur_score>0.0 && da) {
251 for(
unsigned int j=0; j<3; j++) {
252 sphere_derivatives_table[pi_index][j] += (*da)(derivative_vector[j]);
255 radial_displacements+= displacement[0]*s[0] + displacement[1]*s[1];
259 if(da && is_pore_radius_optimized_){
268 SlabWithToroidalPorePairScore::
269 get_sphere_ellipsoid_penetration_depth
274 const double eps= 1e-9;
276 double dXY2= v[0]*v[0]+v[1]*v[1];
277 double dZ2= v[2]*v[2];
278 double dv2= dXY2 + dZ2 + eps;
280 double sinTheta2= dXY2/dv2;
281 double cosTheta2= dZ2/dv2;
282 double cur_r= std::sqrt(rv_*rv_*cosTheta2 + rh_*rh_*sinTheta2);
283 double dv=std::sqrt(dv2);
284 double surface_distance= dv - sphere.get_radius() - cur_r;
285 if(surface_distance>=0 ){
293 (*out_translation)= v.get_unit_vector();
295 return -surface_distance;
303 SlabWithToroidalPorePairScore::get_sphere_penetration_depth
307 double const x=sphere[0];
308 double const y=sphere[1];
309 double const z=sphere[2];
310 double const sr=sphere.get_radius();
311 double sphere_dz_top= (z - sr) - top_;
312 double sphere_dz_bottom = (z + sr) - bottom_;
313 bool is_above_top=sphere_dz_top>0;
314 bool is_below_bottom=sphere_dz_bottom<0;
315 if(is_above_top || is_below_bottom) {
322 double d_xy2 = x*x + y*y;
323 bool is_closer_to_top= (sphere_dz_bottom > -sphere_dz_top);
324 double abs_sphere_dz =
325 is_closer_to_top ? -sphere_dz_top : sphere_dz_bottom;
331 return abs_sphere_dz;
333 double d_xy = std::sqrt(d_xy2);
335 const double eps = 1e-9;
338 double scale= R_/d_xy;
340 ret= get_sphere_ellipsoid_penetration_depth
341 (sphere, origin, out_translation );
344 ret= get_sphere_ellipsoid_penetration_depth
345 (sphere, origin, out_translation );
352 IMPNPCTRANSPORT_END_NAMESPACE
void add_to_derivatives(const algebra::Vector3D &v, DerivativeAccumulator &d)
Add the vector v to the derivative vector of the x,y,z coordinates.
void add_to_pore_radius_derivative(double v, DerivativeAccumulator &d)
Abstract class for scoring object(s) of type ParticleIndexPair.
A decorator for a particle that's a slab with a toroidal pore.
SphereD< 3 > Sphere3D
Typedef for Python.
virtual double evaluate_if_good_indexes(Model *m, const ParticleIndexPairs &o, DerivativeAccumulator *da, double max, unsigned int lower_bound, unsigned int upper_bound) const
Macros for various classes.
#define IMP_FINAL
Have the compiler report an error if anything overrides this method.
#define IMP_OBJECT_METHODS(Name)
Define the basic things needed by any Object.
#define IMP_OBJECT_LOG
Set the log level to the object's log level.
Storage of a model, its restraints, constraints and particles.
bool get_coordinates_are_optimized() const
Get whether the coordinates are optimized.
A more IMP-like version of the std::vector.
#define IMP_LOG_TERSE(expr)
#define IMP_INTERNAL_CHECK(expr, message)
An assertion to check for internal errors in IMP. An IMP::ErrorException will be thrown.
Class for storing model, its restraints, constraints, and particles.
const algebra::Sphere3D & get_sphere() const
Return a sphere object.
Float get_pore_radius() const
get cylindrical pore radius
virtual ModelObjectsTemp do_get_inputs(Model *m, const ParticleIndexes &pis) const =0
Overload this method to specify the inputs.
Score for a slab with a toroidal pore.
const algebra::Vector3D & get_coordinates() const
Convert it to a vector.
#define IMP_CHECK_CODE(expr)
Only compile the code if checks are enabled.
A nullptr-initialized pointer to an IMP Object.
Float get_thickness() const
returns whether the particle last entered the transport moiety from its
Exception definitions and assertions.
virtual double evaluate_indexes(Model *m, const ParticleIndexPairs &o, DerivativeAccumulator *da, unsigned int lower_bound, unsigned int upper_bound) const
Compute the score and the derivative if needed over a set.
virtual double evaluate_if_good_index(Model *m, const ParticleIndexPair &vt, DerivativeAccumulator *da, double max) const
Compute the score and the derivative if needed, only if "good".
Decorator for a sphere-like particle.
virtual double evaluate_index(Model *m, const ParticleIndexPair &vt, DerivativeAccumulator *da) const =0
Compute the score and the derivative if needed.
#define IMP_OVERRIDE
Cause a compile error if this method does not override a parent method.
Class for adding derivatives from restraints to the model.
A decorator for a particle with x,y,z coordinates and a radius.