8 #ifndef IMPNPCTRANSPORT_SLAB_PAIR_SCORE_H
9 #define IMPNPCTRANSPORT_SLAB_PAIR_SCORE_H
11 #include "npctransport_config.h"
18 IMPNPCTRANSPORT_BEGIN_NAMESPACE
27 class IMPNPCTRANSPORTEXPORT
33 mutable double thickness_;
34 mutable double pore_radius_;
36 mutable double bottom_;
38 mutable bool is_pore_radius_optimized_;
50 double get_displacement_magnitude
75 unsigned int lower_bound,
92 unsigned int lower_bound,
96 for (
unsigned int i = lower_bound; i < upper_bound; ++i) {
98 if (ret > max)
return std::numeric_limits<double>::max();
117 inline double evaluate_sphere
128 std::pair<double, algebra::Vector3D> get_displacement_vector(
133 void update_cached_slab_params
138 SlabWithCylindricalPorePairScore::update_cached_slab_params
143 top_= 0.5*thickness_;
144 bottom_= -0.5*thickness_;
147 is_pore_radius_optimized_= slab.get_pore_radius_is_optimized();
152 SlabWithCylindricalPorePairScore::evaluate_index
159 "pip[0] is not a SlabWithCylindricalPore in evaluate_index()");
161 update_cached_slab_params(slab);
167 double score=evaluate_sphere(d_sphere,
168 da ? &displacement :
nullptr);
171 IMP_LOG(
PROGRESS,
"result in " << score <<
" and " << derivative_vector << std::endl);
173 if(is_pore_radius_optimized_){
175 double radial_displacement_magnitude=
176 std::sqrt(displacement[0]*displacement[0]+displacement[1]*displacement[1]);
185 SlabWithCylindricalPorePairScore::evaluate_indexes
189 unsigned int lower_bound,
190 unsigned int upper_bound)
const
192 if(upper_bound<lower_bound){
196 double radial_displacements_magnitude(0.0);
198 m->access_spheres_data();
200 m->access_sphere_derivatives_data();
201 IMP::internal::BoolAttributeTableTraits::Container
const& is_optimizable_table=
202 m->access_optimizeds_data(core::XYZ::get_coordinate_key(0));
205 update_cached_slab_params(slab);
207 for (
unsigned int i = lower_bound; i < upper_bound; ++i) {
209 int pi_index=pi.get_index();
213 "All particles are assumed to be evaluated against"
218 "optimable table inconsistent with d.get_coordinates_are_optimized for particle " << d);
220 "Different coords for particle " << d <<
" *** "
223 "Different radii for particle " << d <<
" *** "
224 << d.get_radius() <<
" vs. " << s.get_radius());
226 if(!is_optimizable_table[pi]) {
230 double cur_score = evaluate_sphere(spheres_table[pi_index],
231 da ? &displacement :
nullptr);
233 if(cur_score>0.0 && da) {
236 for(
unsigned int j=0; j<3; j++) {
237 sphere_derivatives_table[pi_index][j] += (*da)(derivative_vector[j]);
240 radial_displacements_magnitude+=
241 std::sqrt(displacement[0]*displacement[0] + displacement[1]*displacement[1]);
244 if(da && is_pore_radius_optimized_){
252 SlabWithCylindricalPorePairScore::evaluate_sphere
257 IMP_LOG(
VERBOSE,
"evaluate_sphere " << s << std::endl);
258 double const x= s[0];
259 double const y= s[1];
260 double const z= s[2];
261 double const sr= s.get_radius();
263 if ((z-sr > top_) || (z+sr < bottom_)) {
266 double const x2= x*x;
267 double const y2= y*y;
268 double const R= pore_radius_-sr;
269 double const R2= R*R;
274 std::pair<double, algebra::Vector3D> dp = get_displacement_vector(s.get_center());
276 "At point " << s.get_center() <<
" have distance " << dp.first
277 <<
" and direction " << dp.second << std::endl);
278 double const distance = dp.first;
282 double const score= k_ * (sr - distance);
283 if(out_displacement){
284 *out_displacement= dp.second;
286 "Not a unit vector");
296 inline std::pair<double, algebra::Vector3D>
297 SlabWithCylindricalPorePairScore::get_displacement_vector(
const algebra::Vector3D &v)
const {
298 double dXY2 = square(v[0]) + square(v[1]);
299 double dZ = v[2] - midZ_;
301 if (dXY2 > square(pore_radius_) ||
302 (v[2] <= top_ && v[2] >= bottom_)) {
303 double abs_dZ = std::abs(dZ);
304 double abs_dXY = std::sqrt(dXY2);
305 double dR = abs_dXY - pore_radius_;
306 if (dR + abs_dZ < .5 * thickness_) {
312 return std::make_pair(pore_radius_ - abs_dXY, rv.get_unit_vector());
341 return std::make_pair(diff.get_magnitude(), diff.get_unit_vector());
348 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.
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.
bool get_coordinates_are_optimized() const
Get whether the coordinates are optimized.
#define IMP_LOG_PROGRESS(expr)
A more IMP-like version of the std::vector.
#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.
A decorator for a particle that's a slab with a cylindrical pore.
Float get_pore_radius() const
get cylindrical pore radius
Produce copious output to allow someone to trace through the computation.
virtual ModelObjectsTemp do_get_inputs(Model *m, const ParticleIndexes &pis) const =0
Overload this method to specify the inputs.
const algebra::Vector3D & get_coordinates() const
Convert it to a vector.
#define IMP_CHECK_CODE(expr)
Only compile the code if checks are enabled.
Float get_thickness() const
returns whether the particle last entered the transport moiety from its
Helper macros for throwing and handling exceptions.
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.
#define IMP_USAGE_CHECK(expr, message)
A runtime test for incorrect usage of a class or method.
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.