diff --git a/mrdna/model/nbPot.py b/mrdna/model/nbPot.py index c26c3b07f116bc67121642c0a7dbf5c1d28df374..71a3c2a49ff6e93af7d2d10f13118714d7eead33 100644 --- a/mrdna/model/nbPot.py +++ b/mrdna/model/nbPot.py @@ -80,7 +80,9 @@ class AbstractNbDnaScheme(NonbondedScheme): if self.debye_length != self.debye_length0: ## units "e**2/(4 * pi * 80 epsilon0 AA)" kcal_mol A = 4.1507964 - du = (A/x) * (np.exp(-x/self.debye_length) - np.exp(-x/self.debye_length0)) + with np.errstate(divide='ignore',invalid='ignore'): + du = (A/x) * (np.exp(-x/self.debye_length) - np.exp(-x/self.debye_length0)) + du[x < 10] = du[x>=10][0] u = u + du diff --git a/mrdna/segmentmodel.py b/mrdna/segmentmodel.py index 8b85c3f0511f4e2c013ff4f4d3623d30f63f7718..22f205d26ee546b26dd2d596291222565c9d9137 100644 --- a/mrdna/segmentmodel.py +++ b/mrdna/segmentmodel.py @@ -453,6 +453,16 @@ class Segment(ConnectableElement, Group): tck, u = self.position_spline_params return np.mean(self.contour_to_position(u), axis=0) + def get_bounding_box( self, num_points=3 ): + positions = np.zeros( (num_points, 3) ) + i = 0 + for c in np.linspace(0,1,num_points): + positions[i] = (self.contour_to_position(c)) + i += 1 + min_ = np.array([np.min(positions[:,i]) for i in range(3)]) + max_ = np.array([np.max(positions[:,i]) for i in range(3)]) + return min_,max_ + def _get_location_positions(self): return [self.contour_to_nt_pos(l.address) for l in self.locations] @@ -3241,17 +3251,27 @@ proc calcforces {} { } """) - def dimensions_from_structure( self, padding_factor=1.5, isotropic=False ): - positions = [] + def get_bounding_box( self, num_points=3 ): + positions = np.zeros( (len(self.segments)*num_points, 3) ) + i = 0 for s in self.segments: - positions.append(s.contour_to_position(0)) - positions.append(s.contour_to_position(0.5)) - positions.append(s.contour_to_position(1)) - positions = np.array(positions) - dx,dy,dz = [(np.max(positions[:,i])-np.min(positions[:,i])+30)*padding_factor for i in range(3)] + for c in np.linspace(0,1,num_points): + positions[i] = (s.contour_to_position(c)) + i += 1 + min_ = np.array([np.min(positions[:,i]) for i in range(3)]) + max_ = np.array([np.max(positions[:,i]) for i in range(3)]) + return min_,max_ + + def get_bounding_box_center( self, num_points=3 ): + min_,max_ = self.get_bounding_box(num_points) + return 0.5*(max_+min_) + + def dimensions_from_structure( self, padding_factor=1.5, isotropic=False ): + min_,max_ = self.get_bounding_box() + dx,dy,dz = (max_-min_+30)*padding_factor if isotropic: dx = dy = dz = max((dx,dy,dz)) - return [dx,dy,dz] + return np.array([dx,dy,dz]) def add_grid_potential(self, grid_file, scale=1, per_nucleotide=True): grid_file = Path(grid_file)