8 #ifndef IMPNPCTRANSPORT_SLAB_PAIR_SCORE_H
9 #define IMPNPCTRANSPORT_SLAB_PAIR_SCORE_H
11 #include "npctransport_config.h"
18 IMPNPCTRANSPORT_BEGIN_NAMESPACE
31 class IMPNPCTRANSPORTEXPORT
37 mutable double thickness_;
38 mutable double pore_radius_;
40 mutable double bottom_;
42 mutable bool is_pore_radius_optimized_;
56 double get_displacement_magnitude
81 unsigned int lower_bound,
82 unsigned int upper_bound)
const override final;
98 unsigned int lower_bound,
99 unsigned int upper_bound)
const override
102 for (
unsigned int i = lower_bound; i < upper_bound; ++i) {
104 if (ret > max)
return std::numeric_limits<double>::max();
123 inline double evaluate_sphere
134 std::pair<double, algebra::Vector3D> get_displacement_vector(
139 void update_cached_slab_params
144 SlabWithCylindricalPorePairScore::update_cached_slab_params
149 top_= 0.5*thickness_;
150 bottom_= -0.5*thickness_;
153 is_pore_radius_optimized_= slab.get_pore_radius_is_optimized();
158 SlabWithCylindricalPorePairScore::evaluate_index
165 "pip[0] is not a SlabWithCylindricalPore in evaluate_index()");
167 update_cached_slab_params(slab);
173 double score=evaluate_sphere(d_sphere,
174 da ? &displacement :
nullptr);
177 IMP_LOG(
PROGRESS,
"result in " << score <<
" and " << derivative_vector << std::endl);
179 if(is_pore_radius_optimized_){
181 double radial_displacement_magnitude=
182 std::sqrt(displacement[0]*displacement[0]+displacement[1]*displacement[1]);
191 SlabWithCylindricalPorePairScore::evaluate_indexes
195 unsigned int lower_bound,
196 unsigned int upper_bound)
const
198 if(upper_bound<lower_bound){
202 double radial_displacements_magnitude(0.0);
204 m->access_spheres_data();
206 m->access_sphere_derivatives_data();
207 IMP::internal::BoolAttributeTableTraits::Container
const& is_optimizable_table=
208 m->access_optimizeds_data(core::XYZ::get_coordinate_key(0));
211 update_cached_slab_params(slab);
213 for (
unsigned int i = lower_bound; i < upper_bound; ++i) {
215 int pi_index=pi.get_index();
219 "All particles are assumed to be evaluated against"
224 "optimable table inconsistent with d.get_coordinates_are_optimized for particle " << d);
226 "Different coords for particle " << d <<
" *** "
229 "Different radii for particle " << d <<
" *** "
230 << d.get_radius() <<
" vs. " << s.get_radius());
232 if(!is_optimizable_table[pi]) {
236 double cur_score = evaluate_sphere(spheres_table[pi_index],
237 da ? &displacement :
nullptr);
239 if(cur_score>0.0 && da) {
242 for(
unsigned int j=0; j<3; j++) {
243 sphere_derivatives_table[pi_index][j] += (*da)(derivative_vector[j]);
246 radial_displacements_magnitude+=
247 std::sqrt(displacement[0]*displacement[0] + displacement[1]*displacement[1]);
250 if(da && is_pore_radius_optimized_){
258 SlabWithCylindricalPorePairScore::evaluate_sphere
263 IMP_LOG(
VERBOSE,
"evaluate_sphere " << s << std::endl);
264 double const x= s[0];
265 double const y= s[1];
266 double const z= s[2];
267 double const sr= s.get_radius();
269 if ((z-sr > top_) || (z+sr < bottom_)) {
272 double const x2= x*x;
273 double const y2= y*y;
274 double const R= pore_radius_-sr;
275 double const R2= R*R;
280 std::pair<double, algebra::Vector3D> dp = get_displacement_vector(s.get_center());
282 "At point " << s.get_center() <<
" have distance " << dp.first
283 <<
" and direction " << dp.second << std::endl);
284 double const distance = dp.first;
288 double const score= k_ * (sr - distance);
289 if(out_displacement){
290 *out_displacement= dp.second;
292 "Not a unit vector");
302 inline std::pair<double, algebra::Vector3D>
303 SlabWithCylindricalPorePairScore::get_displacement_vector(
const algebra::Vector3D &v)
const {
304 double dXY2 = square(v[0]) + square(v[1]);
305 double dZ = v[2] - midZ_;
307 if (dXY2 > square(pore_radius_) ||
308 (v[2] <= top_ && v[2] >= bottom_)) {
309 double abs_dZ = std::abs(dZ);
310 double abs_dXY = std::sqrt(dXY2);
311 double dR = abs_dXY - pore_radius_;
312 if (dR + abs_dZ < .5 * thickness_) {
318 return std::make_pair(pore_radius_ - abs_dXY, rv.get_unit_vector());
347 return std::make_pair(diff.get_magnitude(), diff.get_unit_vector());
354 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_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.
Class for adding derivatives from restraints to the model.
A decorator for a particle with x,y,z coordinates and a radius.