Controlwriter.py

# This python file will take in a load of simulation data for a GCMC and mask it onto a dl_monte control template file
# Hopefully I'll be able to do it with classes, but we'll see
# To do this, I'm going to split the control file up into several blocks:
# Use block (including optional fed subblock)
# Physical data block
# Simulation data block
# GCMC block
import random
import attr
random.seed()


@attr.s
class dlmonte_control:
    title = None
    T = 77
    P = 0.001
    iterations = 10000
    eqlength = 0
    statsize = int(max(iterations / 1000, 1000))
    outfreq = int(max(iterations / 1000, 1000))

    statfreq = int(max(iterations / 1000, 1000))
    yamlfreq = int(max(iterations / 10000, 1000))
    archivefreq = int(max(iterations / 100, 10000))

    fed = True
    FEDprintfreq = iterations
    updatefreq = 1000
    low = 0
    high = 250
    win=True
    winlow = 0
    winhigh = 10
    res = False

    species = ''
    chempot = 1e-4

    simulation_data = [
        'title',
        'iterations',
        'equilibration_length',
        'statistics_block_size',
        'statistics_output_frequency',
        'output_frequency', # frequency of printing to OUTPUT.000
        'yaml_frequency', # frequency of writing to YAMLDATA.000
        'archive_frequency' # frequency of writing to configurations archive
    ]
    free_energy_data = [
        'fed',
        'fed_flavour',
        'bias_print_frequency',
        'bias_update_frequency',
        'collective_variable_minimum',
        'collective_variable_maximum',
        'collective_variable_window',
        'window_minimum',
        'window_maximum',
        'fed_restart'
    ]
    physics_data = [
        'T',
        'P'
    ]
    GCMC_data = [
        'species',
        'fugacity'
    ]
    all_data = simulation_data + free_energy_data + physics_data + GCMC_data
    @classmethod
    def updatefromdict(cls, dict):
        kwargs = {k: v for k, v in dict.items() if k in cls.all_data}
        cls(**kwargs)
    def printfedblock(self):
        span = self.high-self.low
        if not self.res:
            cont = 'new'
        else:
            cont = 'res'
        output = ''
        output += 'use fed generic\n'
        output += 'fed method tm {0:1.1e} {1} {2}\n'.format(self.FEDprintfreq, self.updatefreq, cont)
        output += 'fed order nmols {0} {1} {2} 1 '.format(span, self.low, self.high)
        if self.win ==True:
            output += 'win {0} {1}\n'.format(self.winlow, self.winhigh)
        else:
            output += '\n'
        output += 'fed done\n\n'
        return output


    def printuseblock(self):
        output = ''
        output += 'use ortho\n'
        output += 'use gaspressure\n\n'
        if self.fed is not None:
            output += self.printfedblock()
        output += 'finish\n\n'
        return output


    def printphysblock(self):
        ran0 = random.randint(1, 99)
        ran1 = random.randint(1, 99)
        ran2 = random.randint(1, 99)
        ran3 = random.randint(1, 99)
        output = ''
        output += 'seeds {0} {1} {2} {3} #4 randomseeds between 0 and 168\n'.format(ran0, ran1, ran2, ran3)
        output += 'temperature {0} #T in K\n'.format(self.T)
        output += '#pressure {0} #kilo atmospheres?\n'.format(self.P)
        output += 'ewald prec 1e-6 #no electrostatics\n\n'
        return output

    def printdatablock(self):
        output = ''
        output += 'steps {0:1.1e}  # number of MC iteractions \n'.format(self.iterations)
        output += 'equilibration {0}  # number of equil steps \n'.format(self.eqlength)
        output +='stack {0:1.1e}  # blocks for statistics \n'.format(self.statsize)
        output +='print {0:1.1e}  # frequency of writing for output files\n'.format(self.outfreq)
        output +='stats {0:1.1e}  # frequency of writing to a stats files\n'.format(self.statfreq)
        output += 'yamldata {0:1.1e}\n\n'.format(self.yamlfreq)
        output += 'sample coordinates {0:1.1e} # how often to archive coordinates\n'.format(self.archivefreq)
        output += '# archiveformat dcd # format of the trajectory file (ARCHIVE/HISTORY/TRAJECTORY)\n'
        output += 'revconformat dlmonte # format of REVCON file (to replace CONFIG if restarting)\n\n'
        return output

    def printGCMCblock(self):
        output = ''
        output += 'move molecule 1 25 \n {0}\n'.format(self.species)
        output += 'acceptmolmoveupdate 200 \n'
        output += 'move rotatemol 1 25 \n {0}\n'.format(self.species)
        output += 'acceptmolrotupdate 200 \n'
        output += 'move gcinsertmol 1 50 0.5 \n {0}  {1:1.4e}\n\n'.format(self.species, self.chempot)
        return output

    def printcontrol(self):
        output = ''
        if self.title:
            output += self.title
        else:
            output += 'An autogenerated DL_MONTE control file has no name\n'
        output += self.printuseblock()
        output += self.printphysblock()
        output += self.printdatablock()
        output += self.printGCMCblock()
        output += 'start simulation'
        return output