Dear Jia,
Sorry for the delay... Yes, using the last version of the R.E.D. Tools
is better. Once the R.E.D. Tools III.4 will be available you will have
to apply the same corrections as these below.
-I- Modifications of Ante_R.E.D. from the III.3 distribution
Yes you need to add "AU" (capital letters) in An_R.E.D. as you started below
*1* line 23 add in %Elements (authorize Gold to be used & provide its
atomic number)
,AU =>"79" just after BR =>"35"
*2* line 26 add in %Valence (define the maximum valency for Gold:
useful for the topology of the FF library: the MOL2 file format
generated at the end by R.E.D.)
,AU =>"4" (or Au =>"6", I do not know) just after BR =>"1"
*3* line 26 add in %Radius (it is useful for the definition of the
atom connectivities - it is different from the radius used in MEP
computation!!!)
,AU =>"1.4" (1.5 should be OK; adapt if not) just after BR =>"1.210"
=> this will generate a P2N file & inputs for geometry optimization
using Firefly, GAMESS-US and Gaussian.
-1- check chemical equivalencing in the P2N file (in the first column
of atom names)
see
http://q4md-forcefieldtools.org/Tutorial/Tutorial-1.php#ATOM-NAME
(Ante_R.E.D. 2.0 solves all here).
-2- check the atom connectivity for Gold in the P2N file: you should
have 4 connections with Gold in your bio-inorganic complex (affected
by the radius defined in *3* above: if not correct increase 1.5 to 1.6
etc...)
see
http://www.wwpdb.org/documentation/format30/index.html
-3- check the total charge & spin multiplicity of your complex in the
P2N file & QM inputs
-4- adapt the QM input you use (GAMESS, Firefly or Gaussian) so that
you select the theory level of your choice; as a starting point I
would use B3LYP/6-31G* & you have to define an Effective Core
Potential (ECP) for Gold (I think).
=> run the geometry optimization output without R.E.D. first to see
how it goes with ECP...
I would start using a small model such as [AuH4]+/-x before to switch
on your bio-inorganic complex (quite big) to set up all the charge
derivation process. This will allow you to quickly get results.
=======
-II- Modifications of R.E.D.-III.3
*1* line 3872 add in %Elements (authorize Gold to be used & provide
its atomic number)
,'AU' =>"79" just after 'BR' =>"35"
*2* in the function "Readpdb"
line 656 just before %L4atoms, add the following new code Line:
%atomAU = ('AU'=>"79"); # Gold
*3* in the function "Readpdb"
line 658 add in "Twoletterelt"
,'AU' =>"AU" just after 'BR'=>"BR"
*4* just before the ending bracket "}" line 870 of the function "Readpdb", add
$TestAU[$NM]=0;
for ($i=0; $i<$nbatoms[$NM]; $i++){
if(exists $atomAU{$tab[1][$i][$NM]}){ $TestAU[$NM]=1; last; }
}
=> Here, you create a variable $TestLAU[$NM] for each molecule (you
will be able to derive charge for different Au complexes if you need);
this variable = 1 if the Gold element is found in the P2N file.
*5* line 1206-1230 you will find the lines to modify the geometry
optimization input ("file handle JOB1_FILE"). With the variable you
defined in *4* you can now modify the input according to the presence
of the Gold element
if ($TestAU[$NM]==1)) { print JOB1_FILE " blabla for GAMESS/Firefly \n"; }
same idea with Gaussian input line 1306
if ($TestAU[$NM]==1)) { print JOB1_FILE " blabla for Gaussian... \n"; }
*6* line 1757-1793 you will find the lines to modify the MEP
computation input ("file handle JOB2_FILE"). With the variable you
defined in *4* you can again modify the input according to the
presence of the Gold element
if ($TestAU[$NM]==1)) { print JOB2_FILE " blabla for GAMESS/Firefly \n"; }
Same idea with Gaussian input line 1874
if ($TestAU[$NM]==1)) { print JOB2_FILE " blabla for Gaussian \n"; }
*7* You need to find out about the radius of Gold to be used in MEP
computation. For GAMESS/Firefly it is hard coded in the source code.
For Gaussian, you need to provide it at the end of the input (after
having jumped a line): See 1905-1911 (this part has been re-written is
III.4)
I would start with the 1.8 value & see how this cuisine works first
;-) Then you might decide to fine tune this value...
This should be OK. It is not that difficult; you only need to proceed
step by step starting from a basic Gold model and understand what you
do. Then, when you will switch to your bio-inorganic complex you will
only have to handle additional intra-molecular charge constraints for
the two dipeptides involved in the charge derivation.
Fore sure, here, a tutorial should be written... I hope this helps.
regards, Francois
> Dear Francois,
>
> I really appreciate your kind help!
>
> If possible, would you please to tell me how to modify the source
> code according to step 1 to step 3?
>
> On Fri, Jun 25, 2010 at 1:25 PM, FyD <fyd.q4md-forcefieldtools.org> wrote:
>
>> Dear Jia,
>>
>> I am trying to simulate a peptide with a metal-containing complex. As
>>> shown in the attached picture, the gold atom (yellow ball) of
>>> phenanthroline-gold is bound to N and S atoms of HIS and CYS,
>>> respectively.
>>> I used to calculate the resp charge of non-standard residue by R. E. D.
>>> III.
>>> It works fine. But I have no experience to work on small molecules.
>>>
>>> I would like to know if it possible to obtain the RESP charge of such a
>>> heavy-metal contained complex which is bound to two residues?
>>>
>>
>> You need to solve several problems in R.E.D.-III.x to get charges for your
>> Gold complex:
>>
>> 1) Add the "Au" atom in the source codes to authorize R.E.D. & Ante_R.E.D.
>> to use "Au"
>> 2) Add which "radius" value you want to use for "Au"
>>
>
> I have checked the perl scripts and found lines as follows. I guess the
> type and radius of atom should be added in these lines, is it right? However
> I am not sure which one of the three kinds of radii(# Radii used by
> Ante_R.E.D.; # Empirical radii; # Calculed radii ) is default, would you
> please to tell me?
>
> # Information:
> # One can modify the "$sort" & "$notlig" variables in the main section of
> the program (line 591)
> #
> ######## Chemical elements recognized by Ante_R.E.D.: up to Bromine (Z =
> 35)...
> # 'Terminal' chemical elements that may generate vdW bumps in a bad initial
> structure
> my %Terminal = (H =>"1",F =>"9",CL =>"17",BR =>"35");
>
> # Elements: Warning: Calcium is XX !!!
> my %Elements = (H =>"1",LI =>"3",BE =>"4",B =>"5",C =>"6",N =>"7",O =>"8",F
> =>"9",NA =>"11",MG =>"12",AL =>"13",SI =>"14",P =>"15",S =>"16",CL =>"17",K
> =>"19",XX =>"20",SC=>"21",TI => "22",V=>"23",CR =>"24",MN=>"25",FE =>"26",CO
> =>"27",NI =>"28",CU =>"29",ZN =>"30",GA =>"31","GE" =>"32",AS =>"33",SE
> =>"34",BR =>"35");
>
> # Maximal number of bonds for each chemical element
> my %Valence = (H =>"1",LI =>"1",BE =>"2",B =>"3",C =>"4",N =>"4",O =>"2",F
> =>"1",NA =>"1",MG =>"2",AL =>"6",SI =>"6",P =>"5",S =>"6",CL =>"1",K
> =>"1",XX =>"2",SC=>"6",TI => "6",V=>"6",CR =>"6",MN=>"8",FE =>"6",CO
> =>"6",NI =>"6",CU =>"6",ZN =>"6",GA =>"3","GE" =>"4",AS =>"3",SE =>"2",BR
> =>"1");
>
> # Radii used by Ante_R.E.D.
> my %Radius = (H =>"0.230",LI =>"0.680",BE =>"0.350",B =>"0.830",C
> =>"0.680",N =>"0.680",O =>"0.680",F =>"0.640",NA =>"0.970",MG =>"1.100",AL
> =>"1.350",SI =>"1.200",P =>"0.750",S =>"1.020",CL =>"0.990",K =>"1.330",CA
> =>"0.990",SC=>"1.440",TI => "1.470",V=>"1.330",CR =>"1.350",MN=>"1.350",FE
> =>"1.340",CO =>"1.330",NI =>"1.500",CU =>"1.520",ZN =>"1.450",GA
> =>"1.220","GE" =>"1.170",AS =>"1.210",SE =>"1.220",BR =>"1.210");
>
> # Empirical radii: You might try the "Empirical radii" instead of the "Radii
> used by Ante_R.E.D."...
> #my %Radius = (H =>"0.250",LI =>"1.450",BE =>"1.050",B =>"0.850",C
> =>"0.700",N =>"0.650",O =>"0.600",F =>"0.500",NA =>"1.800",MG =>"1.500",AL
> =>"1.250",SI =>"1.100",P =>"1.000",S =>"1.000",CL =>"1.000",K =>"2.200",CA
> =>"1.800",SC=>"1.600",TI => "1.400",V=>"1.350",CR =>"1.400",MN=>"1.400",FE
> =>"1.400",CO =>"1.350",NI =>"1.350",CU =>"1.350",ZN =>"1.350",GA
> =>"1.300","GE" =>"1.250",AS =>"1.150",SE =>"1.150",BR =>"1.150");
>
> # Calculed radii: You might try the "Calculated radii" instead of the "Radii
> used by Ante_R.E.D."...
> # my %Radius = (H =>"0.530",LI =>"1.670",BE =>"1.120",B =>"0.870",C
> =>"0.670",N =>"0.560",O =>"0.480",F =>"0.420",NA =>"1.900",MG =>"1.450",AL
> =>"1.180",SI =>"1.110",P =>"0.980",S =>"0.880",CL =>"0.790",K =>"2.430",CA
> =>"1.940",SC=>"1.840",TI => "1.760",V=>"1.710",CR =>"1.660",MN=>"1.610",FE
> =>"1.560",CO =>"1.520",NI =>"1.490",CU =>"1.450",ZN =>"1.420",GA
> =>"1.360","GE" =>"1.250",AS =>"1.140",SE =>"1.030",BR =>"0.94");
>
>
>
>> 3) Add which "theory level" you want to use during the geometry
>> optimization & MEP computation step
>>
>
> I have little knowledge about quantum chemistry software such as
> GUASSIAN and GAMESS, possibly as a result I failed to find the "theory
> level" in the pl scripts. Would you please to help me?
>
>
>> 4) Select a target complex to be involved in charge derivation using the
>> RESP program. Here, you need to imagine how you could include your complex
>> in your whole molecule. You select a model & you use constraint(s) during
>> the RESP fit to make your model compatible with the whole system: Once you
>> understood the use of intra-molecular charge constraint in dipeptide for
>> instance you can easily apply the same strategy to your Au complex (using
>> two amino acid central fragment).
>> See http://q4md-forcefieldtools.org/Tutorial/Tutorial-1.php#10
>>
>> I can help you setting these steps if you need...
>>
>> The next version of the R.E.D. Tools source code is going to be distributed
>> under the GNU General Public License. We will see if we can release III.4
>> next week.
>>
>> I am now using R.E.D-III.2&GAMESS-US, is it better to wait for the
> upcoming III.4 version?
>
>
>> regards, Francois
>>
>>
> Thank you so much!
>
> Best regards,
> Sincerely,
> Jia
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Tue Jun 29 2010 - 23:30:03 PDT
