input coordinates & charges

[Csilla Varnai]

Dear all,

I am trying to read a small peptide to run an MD calculation. I cannot
define the atomic charges in the FORCEFIELD section, because there is
no one-to-one mapping between the atomic types and the atomic charges
in case of peptides. Hence, I need to read in the charges. I have
found only 1 way to do it so far, namely reading it from a PDB file.
My problem is that in this case I loose the accuracy of the
coordinates, since this format is fairly restricted and I cannot read
the coordinates only at a +/-0.001 level.
Could anyone please help me, if there is a better way I've missed to
get both the charges and the coordinates at a given accuracy?
Thanks very much,

Csilla

[Teodoro Laino]

Dear Csilla,

you got exactly the point. 

In CP2K we take the CHARMM convention where a type has associated a charge.

This is different from AMBER where the same type can have different charges.

So the best thing to do is to redefine your types according the convention that a type

has a well defined charge (this most of the times is terribly easy.. just append a letter or a number to

the existing type).

Cheers

Teo

p.s.: even with PDB this won't necessarily work. I'm not TOTALLY sure but  in case you define two different

charges for the same type, only the first one may be assigned to that type. The result is that the

charges of your system will be wrong. I will add a warning for handling this exception, though reading 

the charges from the beta column of the PDB is not something really standard.

[Csilla Varnai]

Dear Teo,

Thanks very much!

Best regards,

Csilla

[Axel]

dear csilla,

On May 7, 11:55 am, Csilla Varnai <cv...@cam.ac.uk> wrote:
> Dear all,
>
> I am trying to read a small peptide to run an MD calculation. I cannot
> define the atomic charges in the FORCEFIELD section, because there is
> no one-to-one mapping between the atomic types and the atomic charges
> in case of peptides. Hence, I need to read in the charges. I have

i'm a bit confused as to why you want to manually enter
all the topology parameters? why don't you just use a
tool like psfgen from VMD/NAMD to build a .psf file that
contains all the information more or less automatically?

cheers,
axel.

> found only 1 way to do it so far, namely reading it from a PDB file.
> My problem is that in this case I loose the accuracy of the
> coordinates, since this format is fairly restricted and I cannot read
> the coordinates only at a +/-0.001 level.
> Could anyone please help me, if there is a better way I've missed to
> get both the charges and the coordinates at a given accuracy?
> Thanks very much,
>
> Csilla
>

> cp2k.inp
> 1KDownload
>
>

[Teodoro Laino]

well.. because she wants to exploit the very powerful and reliable topology constructor of CP2K.

I can understand her ;-) LOL

Teo

[Csilla Varnai]

Dear all,

I have redefined the atomic types, and regenerated the CHARMM22 potential file according to the new atomic types, at first sight in a very simple-minded way, substituting all the possible new atomic types (say, CC1, CC2) into the old ones (CC). But in this case, due to the large potential file, running the program is hopeless.

I must confess I have not tried the VMD psfgen, since I would like to invoke CP2K from another program, and that is the reason why I would like to avoid the use of other programs like the VMD psfgen. If I used the cp2k f77_interface, would it make any difference concerning to this problem? Should I also have a potential file containing atomic types with unambiguously defined charges? If that is the case then I run into the same memory problem. Does the psfgen seem to be the only solution?

Best regards,

Csilla

[Teodoro Laino]

Dear Csilla,

First of all let's say that there's no difference if you use cp2k standalone of you invoke it from another program.

Our convention is quite clear about it: a kind has a well defined charge. I.e. you can't provide an array of different charges

associated to the same kind.

The problem is in the definition of your problem.

So if you want some help, we need more details.

How did you generated the CHARMM parameter file? is it one of the standard publicly distributed?

or you generated it with external tools starting maybe from other brand topology files?

are you also providing a PSF or you rely on the internal topology builder?

I rarely used PSFGEN, so I cannot comment how helpful it would be in solving this problem.. But essentially,

with some lines of scripting, you can do the type conversion without getting mad (both for PDB, parameter file and so).

let us know withe above informations if you need further help!

Thanks,

Teo

[Axel]

On May 12, 6:19 am, Csilla Varnai <cv...@cam.ac.uk> wrote:
> Dear all,

dear csilla,

> I have redefined the atomic types, and regenerated the CHARMM22
> potential file according to the new atomic types, at first sight in a
> very simple-minded way, substituting all the possible new atomic types
> (say, CC1, CC2) into the old ones (CC). But in this case, due to the
> large potential file, running the program is hopeless.

what program? cp2k? hopeless in which sense?

it would be much easier to help you, if you'd explain
in more detail what you are actually trying to achieve
and what kind of systems you are planning to model.

> I must confess I have not tried the VMD psfgen, since I would like to
> invoke CP2K from another program, and that is the reason why I would
> like to avoid the use of other programs like the VMD psfgen. If I

whatever you do you _have_ to generate a "topology".
using a .psf and .par file is a very convenient way.
if you don't like using psfgen, you could consider
creating your topology directly. have a look at

http://www.ks.uiuc.edu/Training/Tutorials/namd/namd-tutorial-unix-html/node21.html

and the following pages for a minimal description of the format.
more details are in the parametrization tutorial on the same
server and - of course - on http://www.charmm.org/

unless you have to generate a large number of topologies
on the fly, i would guess that actually using psfgen with
some smart VMD/tcl scripting will be much faster than
writing your own tool.

> the cp2k f77_interface, would it make any difference concerning to
> this problem? Should I also have a potential file containing atomic
> types with unambiguously defined charges? If that is the case then I
> run into the same memory problem. Does the psfgen seem to be the only
> solution?

what memory problem?

cheers,
axel.

[Csilla Varnai]

Dear All,

I am trying to run a program which invokes CP2K at each time step.
Basically, I will need to run CP2K to calculate the forces during my
simulations, so I am trying to find the best way to do that.
I have the atomic coordinates and connectivities stored during the
run, from which I can create CHARMM format to pass it to CP2K.

I have generated a modified forcefield file of the original
distribution of CHARMM from the www.charmm.org home page to have a 1-
to-1 mapping between the atomic types and atomic charges. Having
revised the way of generating the modified forcefield I have a
forcefield of reasonable size, and the forcefield can be read by CP2K
fine. I seem to have no problem with the forcefield any longer, thanks
for the help!

I thought of relying on the internal topology builder of CP2K.
But I am getting the following error message with the default
BONDLENGTH_MAX = 3 Angstroms:

GENERATE| ERROR in connectivity generation!
GENERATE| The THRESHOLD to select possible bonds is larger than
the max. bondlength
GENERATE| used to build the neighbors lists. Increase the
BONDLENGTH_MAX parameter
GENERATE| Present THRESHOLD ( 7.233872 ). Present
BONDLENGTH_MAX ( 5.669178 )

Why does it say 5.669178 instead of 6? Is its unit bohr?
If I increase the BONDLENGTH_MAX, it does not help. If I set 4.5
Angstroms then I am getting the following:

GENERATE| Preliminary Number of Bonds
generated: 15757

GENERATE| WARNING in connectivity generation!
GENERATE| Two molecules/residues named (HSE) have different
number of atoms.
GENERATE| Molecule starting at position (238) has Nr. <4> of atoms.
GENERATE| while the other same molecules have Nr. <11> of atoms!
GENERATE| Increasing bondparm_factor by 1.05.. An error was found
in the generated
GENERATE| connectivity. Retry...
GENERATE| Present value of BONDPARM_FACTOR ( 1.155000 ).
...

GENERATE| ERROR in connectivity generation!
GENERATE| The THRESHOLD to select possible bonds is bigger than
the MAX bondlength
GENERATE| used to build the neighbors lists. Increase the
BONDLENGTH_MAX patameter
GENERATE| Present THRESHOLD ( 8.792816 ). Present
BONDLENGTH_MAX ( 8.503768 )

Can the hydrogen bonds cause problems in this case?
Thanks,

Csilla

[Teodoro Laino]

Dear Csilla,

3/0.529177
5.66918063332306581729
It's in bohr..

> If I increase the BONDLENGTH_MAX, it does not help. If I set 4.5
> Angstroms then I am getting the following:
>
>
>
> GENERATE| Preliminary Number of Bonds
> generated: 15757
>
> GENERATE| WARNING in connectivity generation!
> GENERATE| Two molecules/residues named (HSE) have different
> number of atoms.
> GENERATE| Molecule starting at position (238) has Nr. <4> of
> atoms.
> GENERATE| while the other same molecules have Nr. <11> of atoms!
> GENERATE| Increasing bondparm_factor by 1.05.. An error was found
> in the generated
> GENERATE| connectivity. Retry...
> GENERATE| Present value of BONDPARM_FACTOR ( 1.155000 ).
> ...
>
> GENERATE| ERROR in connectivity generation!
> GENERATE| The THRESHOLD to select possible bonds is bigger than
> the MAX bondlength
> GENERATE| used to build the neighbors lists. Increase the
> BONDLENGTH_MAX patameter
> GENERATE| Present THRESHOLD ( 8.792816 ). Present
> BONDLENGTH_MAX ( 8.503768 )
>
>
> Can the hydrogen bonds cause problems in this case?

Most of the times yes but these errors may appear also because the
values of the parameters
in the input files are totally wrong.. The construction of the
topology can be very hard and difficult and requires some
idea of the bond presents and how to catch them properly...
If you want to use the internal topology generator you've to play
with the parameter of the &GENERATE section.
(In this case the calue of 4.5 it's extremely large!!)

Otherwise I would suggest you to provide a PSF.

Cheers
Teo

[Fawzi Mohamed]

On May 13, 2008, at 2:31 PM, Csilla Varnai wrote:

>
> Dear All,
>
> I am trying to run a program which invokes CP2K at each time step.
> Basically, I will need to run CP2K to calculate the forces during my
> simulations, so I am trying to find the best way to do that.
> I have the atomic coordinates and connectivities stored during the
> run, from which I can create CHARMM format to pass it to CP2K.

To repetively invoke cp2k, especially if just the coordinates change,
you can look at building cpshell (make cpshell VERSION=sopt) that
produces a cp2k_shell.sopt executable that can be called interactively.
For it there is a python wrapping script in tool/cp2kController.
You can also use the f77_interface, it isn't so difficult.
Any of these solutions let you calculate the topology just once (so
maybe even the internal topology builder can cope with it), and then
change coordinates and get updated energy/forces.

ciao
Fawzi

> [...]

[Teodoro Laino]

Anyway if you succeeded to build the topology (with the internal topology creator) at least

once, you can always dump a PSF (see section DUMP_PSF) and use the PSF for the whole set

of other calculations..

Teo

[Csilla Varnai]

Thanks very much! Great. That sounds really useful.

If I would like to run CP2K via the f77_interface, how should I compile CP2K as a library?

Regards,

Csilla

[Fawzi Mohamed]

On May 13, 2008, at 3:26 PM, Csilla Varnai wrote:

Thanks very much! Great. That sounds really useful.

If I would like to run CP2K via the f77_interface, how should I compile CP2K as a library?

Regards,

make libs VERSION=sopt

Fawzi