1

I want to read the data from two different files. In my code I have two polymer objects and I want to read their initial coordinates from the csv file. At this moment I can read only one file. Could someone suggest how I can do this for the second?

Here is the complete code and the method used for reading the first file. 

#Define an object called 'monomer' which consists of a position (list with 3 entries)
class monomer:
    def __init__(self, pos):
    # assign the given position to the 'pos' attribute
        self.pos = pos
        #self.m = np.loadtxt('traj_T_4.csv')

    @classmethod
    def from_file(cls, filename):
        instances = []
        with open(filename, 'r') as file:
            for line in file:
                pos = str(map(float, line.split(',')))
                instances.append(cls(pos))
        return instances

    def LJ(self, other):
    # Define the LJ interaction between one and another monomer

        # Calculate the distance between position of one and the other monomer
        rij =self.distance(other)

        sig_by_r6 = np.power(sigma**2 / rij, 6)
        sig_by_r12 = np.power(sigma**2 / rij, 12)

        return 4*epsilon*(sig_by_r12 - sig_by_r6)

    def FENE(self, other):
        rij =self.distance(other)


        fene = (-0.5 * K * np.power(R,2) * np.log(1 - ((np.sqrt(rij) - r0) / R)**2))


        if not np.isfinite(fene): # Make sure, that NaN (or - infinity) values (because of the log) are never passed to the energy calculation, as the programme behaves uncontrollably then
            return float("inf")   # Rather return infinity, since this will always fail the energy comparison
        else:
            return fene

    def distance(self, other):

        return np.linalg.norm(list(np.array(self.pos) - np.array(other.pos)))

    def shift(self, vector):
    # Move a monomer by 'vector'
         self.pos += vector

# Define an object called 'polymer', which is a chain (list) of monomer objects
class polymer:
    def __init__(self, base, N):

        self.chain = []
        # append each monomer to the polymer's chain
        # base (list of 2) represents the initial x and y value of all monomers
        for i in range(N):
            self.chain.append( monomer([base[0], base[1], i*0.7]) )

    def FENE_Energy(self):
        # Define polymer FENE energy as sum of FENE interactions between bonded monomers
        fene  = 0.0
        for i in range(N-1):
            fene += self.chain[i].FENE(self.chain[i+1])
        return fene

    def LJ_Energy(self):
        LJ = 0.0

        for i in range(N):
            for j in range(N):
                if np.absolute(i-j) > 1:
                    if self.chain[i].distance(self.chain[j]) < rcutoff: # LJ only for non-bonded monomers, i.e. distance is greater than 1
                        LJ += self.chain[i].LJ(self.chain[j])

        return LJ/2.

    def Energy(self):
        return self.LJ_Energy() + self.FENE_Energy()

    def move(self): # Method returns a 'polymer' object, which is either the same as put in or a shifted one
        backup = copy.deepcopy(self) # "backup = self" does not work as expected, as 'backup' becomes a mere synonym of 'self', but not an independent (immutable) object, as intended

        monomer_index = np.random.randint(low=1, high=N) # choose a random monomer but 0, since it is grafted

        direction = np.random.rand(3) - [0.5, 0.5, 0.5] # random vector consists of movement between [-0.5, -0.5, -0.5] and [0.5, 0.5, 0.5]

        self.chain[monomer_index].shift(direction)
        if self.chain[monomer_index].pos[2] < 0.0:
            return backup

        rnd = np.random.rand(1)
        deltaE = self.Energy() - backup.Energy()

        if np.exp(-(1/k_B*T)*deltaE) < rnd:
            return backup
        else:
            return self # return the shifted polymer

if __name__ == "__main__":

    list_of_monomers = monomer.from_file('monomer.csv')
    list_of_monomers = monomer.from_file('monomer1.csv')

    K = 40.0
    R = 0.3
    r0 = 0.7
    sigma = r0/2**6
    rcutoff = 2.5*sigma
    epsilon = 1.0
    max_delta = 0.21 # -sigma/10 to sigma/10
    N = 20
    k_B = 1
    T = 4.0
    shifts = 100

# Declare two variables to be 'polymers'
    pol1 = polymer([-r0/2, 0.0], N)
    pol2 = polymer([ r0/2, 0.0], N)
Improve this question
3
  • 1
    Just open the file and set appropriate values in __init__. Commented Feb 12, 2020 at 12:37
  • or add a classmethod monomer.from_file(cls, filename) and reuse the current init there. Commented Feb 12, 2020 at 12:56
  • @Lars Can you please show that with an example? Commented Feb 12, 2020 at 14:57

1 Answer 1

Reset to default
1

I dont know the format of the file, but lets say it is comma separated with a 3D coordinate on each line:

class monomer:
    @classmethod
    def from_file(cls, filename):
        instances = []
        with open(filename, ‘r’) as file:
            for line in file:
                pos = list(map(float, line.split(‘,’)))
                instances.append(cls(pos))
        return instances

Which can be called as:

list_of_monomers = monomer.from_file(‘monomers.csv’)
Improve this answer
Sign up to request clarification or add additional context in comments.

3 Comments

I have a follow up on this. I tried to do this but I have a problem where I have two different polymer chains. The method you suggested above works well for the first chain, but how can I also include the second chain from the file. I edit the code above. Thanks!
@Mahesh, then you need some way to encode the separate the positions for the different polymers in the file (or use different files), say separate them with ‘|’ and split on ‘|’ first, before creating the list of monomers from the positions as in the code above. (So it depends on how you separate the positions for the different polymers in the file).
Then again it might be easier to use separate files for each polymer.