[gmx-users] Re: Indexing problem when using genconf

[Ryan S Davis (rsdavis1)]

Justin A. Lemkul wrote:
>
>
> Ryan S Davis (rsdavis1) wrote:
>> I wanted to copy a bilayer into a grid of 2x2x1 replicas. I used
>> genconf and everything seemed to work fine exept that annoying feature
>> that the command does not reorder the molecule types, so I end up with
>> a .top file looking like this...
>>
>>   1 #include "martini_v2.1.itp"
>>   2 #include "martini_v2.0_lipids.itp"
>>   3 #include "martini_v2.0_cholesterol.itp"
>>   4   5 [ system ]
>>   6 CHOL
>>   7   8 [ molecules ]
>>   9 DPPC 832
>>  10 CHOL 208
>>  11 W 8320
>>  12 DPPC 832
>>  13 CHOL 208
>>  14 W 8320
>>  15 DPPC 832
>>  16 CHOL 208
>>  17 W 8320
>>  18 DPPC 832
>>  19 CHOL 208
>>  20 W 8320
>>
>>
>> Anyway, I run the simulation...no errors. I make an ndx file using
>> make_ndx...indices look fine despite the repetitive order. HOWEVER,
>> when I try to run commands such as
>> trjconv with the index file as input, it reads all the way up to the
>> first block of Waters and quits with the error
>>
>> """
>> Program trjconv, VERSION 4.0.7
>> Source code file: gmx_trjconv.c, line: 1037
>>
>> Fatal error:
>> Index[29952] 46593 is larger than the number of atoms in the
>> trajectory file (46592)
>> """
>>
>> which I didnt expect, but makes perfect sense knowing that I specified
>> in the .mdp file to not output water to the xtc file...
>>
>> """
>>  xtc-grps                 = dppc chol
>> """
>>
>> Normally this isnt an issue because waters are typically last in the
>> topology. But, I still need access to this data. How can I force the
>> post-processing commands to read past the absent water blocks?
>>
>> The only options I see at the moments is to
>> 1) scrap genconf, make new topology somehow, and rerun
>> 2) reset to output water, and rerun
>> 3) limit my analysis to the very sparse output from the .trr file
>>
>
> I see two viable options, one of which is to use your option 1 via a few
> standard Unix tools to create a proper coordinate file, i.e.:
>
> grep DPPC conf.gro > dppc
> grep CHOL conf.gro > chol
> grep W conf.gro > w
>
> cat dppc chol w > system.gro
>
> Add a title, number of atoms, and box vectors at the end, and you're
> done.  Then sum up the entries in [molecules] and the system should run
> just fine.
>
> The second option is to just create a dummy system that has only DPPC
> and CHOL molecules.  Take your input coordinate file, strip out the
> waters, and renumber with genconf.  Remove all W blocks from the .top
> and create a new .tpr file. This should now match the number of atoms in
> the trajectory and the order in which they appear.
>

...or just do this all in one step with tpbconv.  I always manage to quickly
think of the hardest way to do things :)

-Justin


Thanks, got it working with tpbconv.
It would be nice to have an option on genconf to reorder the molecule types in the output...
Not that I am willing to make that modification myself, ha :)
Ryan