13 namespace AmberParm7Reader {
17 size_t left_quotation_mark = str.find(
"'");
18 if (left_quotation_mark == std::string::npos) {
19 std::string num_fields, type, width,
dot, precision;
20 for (
size_t i = 0; i < str.size(); ++i) {
22 if (std::isdigit(str[i])) {
24 if (type.size() > 0) {
38 }
else if (std::isalpha(str[i])) {
41 }
else if (str[i] ==
'.') {
49 if (num_fields.size() > 0) {
50 specifier.
NumOfFields = std::stoul(num_fields.c_str());
56 specifier.
Type = type.front();
58 specifier.
Width = std::stoul(width.c_str());
59 if (precision.size() > 0) {
60 specifier.
Precision = std::stoul(precision.c_str());
66 size_t right_quotation_mark = str.find(
"'", left_quotation_mark + 1);
69 str.substr(left_quotation_mark + 1,
70 right_quotation_mark - (left_quotation_mark + 1));
77 const std::string& delim,
78 std::vector<std::string>& dest) {
79 size_t index = 0, new_index = 0;
81 while (index != data.length()) {
82 new_index = data.find(delim, index);
83 if (new_index != std::string::npos) tmpstr = data.substr(index, new_index - index);
84 else tmpstr = data.substr(index, data.length());
85 if (!tmpstr.empty()) {
86 dest.push_back(tmpstr);
88 if (new_index == std::string::npos)
break;
89 index = new_index + 1;
94 vector<FortranFormatSpecifier> &specifiers) {
96 vector<string> sub_format_strings;
98 specifiers.resize(sub_format_strings.size());
99 for (
size_t i = 0; i < sub_format_strings.size(); ++i) {
105 const string &source,
const vector<FortranFormatSpecifier> &specifiers,
106 vector<FortranData> &destination) {
107 size_t start_pos = 0;
108 for (
const auto &spec : specifiers) {
109 if (spec.IsPadding ==
true) {
113 for (
int i = 0; i < spec.NumOfFields; ++i) {
114 if (start_pos < source.size()) {
115 const string tmp = source.substr(start_pos, spec.Width);
116 switch (toupper(spec.Type)) {
118 destination.push_back(std::stoi(tmp));
121 destination.push_back(tmp);
129 destination.push_back(std::stod(tmp));
131 destination.push_back(std::stof(tmp));
135 start_pos += spec.Width;
144 using std::istringstream;
148 ifstream ifs_parm(filename);
149 if (!ifs_parm.is_open())
return false;
152 bool match_flag =
false;
153 bool match_format =
false;
155 vector<FortranFormatSpecifier> specifiers;
156 tuple<bool, vector<FortranData>> *p_map =
nullptr;
158 while (std::getline(ifs_parm, line)) {
164 if (line[0] ==
'%') {
165 if (line.find(
"%COMMENT") == 0) {
167 }
else if (line.find(
"%VERSION") == 0) {
169 char tmp_version_stamp[16];
172 if (sscanf(line.c_str(),
"%%VERSION VERSION_STAMP = %9s DATE = %8s %8s",
173 tmp_version_stamp, tmp_date, tmp_time) == 3) {
178 }
else if (line.find(
"%FLAG") == 0) {
180 if (match_flag && match_format) {
182 current_flag.clear();
185 match_format =
false;
187 istringstream iss(line.substr(strlen(
"%FLAG")));
189 p_map = &(toppar_map[current_flag]);
190 std::get<0>(*p_map) =
false;
193 }
else if (line.find(
"%FORMAT") == 0) {
195 if (match_flag ==
false) {
199 const size_t format_start = line.find(
"(", strlen(
"%FORMAT"));
200 const size_t format_end = line.find(
")", format_start + 1);
202 line.substr(format_start + 1, format_end - (format_start + 1)),
209 if (match_flag && match_format) {
225 auto current_section = toppar_map.find(
"TITLE");
226 if (current_section == toppar_map.end()) {
228 current_section = toppar_map.find(
"CTITLE");
229 if (current_section == toppar_map.end()) {
230 iout <<
iERROR <<
"TITLE/CTITLE section for AMBER/CHARMM force field is not "
236 iout <<
iWARN <<
"You are using CHARMM force field in AMBER format, which is not supported in NAMD!\n" <<
endi;
242 toppar_data.
ititl = get<1>(current_section->second).size() > 0 ? get<1>(current_section->second)[0].String :
string();
243 get<0>(current_section->second) =
true;
247 const string section_name{
"POINTERS"};
248 const auto ¤t_section = toppar_map.find(section_name);
249 if (current_section != toppar_map.end()) {
250 success =
parse_pointer(get<1>(current_section->second), toppar_data);
251 get<0>(current_section->second) =
true;
253 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
258 const string section_name{
"ATOM_NAME"};
259 const auto ¤t_section = toppar_map.find(section_name);
260 if (current_section != toppar_map.end()) {
263 get<0>(current_section->second) =
true;
265 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
270 const string section_name{
"CHARGE"};
271 const auto ¤t_section = toppar_map.find(section_name);
272 if (current_section != toppar_map.end()) {
274 toppar_data.
Charges, section_name);
275 get<0>(current_section->second) =
true;
277 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
282 const string section_name{
"ATOMIC_NUMBER"};
283 const auto ¤t_section = toppar_map.find(section_name);
284 if (current_section != toppar_map.end()) {
287 get<0>(current_section->second) =
true;
289 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
294 const string section_name{
"MASS"};
295 const auto ¤t_section = toppar_map.find(section_name);
296 if (current_section != toppar_map.end()) {
298 toppar_data.
Masses, section_name);
299 get<0>(current_section->second) =
true;
301 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
306 const string section_name{
"ATOM_TYPE_INDEX"};
307 const auto ¤t_section = toppar_map.find(section_name);
308 if (current_section != toppar_map.end()) {
310 toppar_data.
Iac, section_name);
311 get<0>(current_section->second) =
true;
313 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
318 const string section_name{
"NUMBER_EXCLUDED_ATOMS"};
319 const auto ¤t_section = toppar_map.find(section_name);
320 if (current_section != toppar_map.end()) {
322 toppar_data.
Iblo, section_name);
323 get<0>(current_section->second) =
true;
325 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
330 const int Ntype2d = toppar_data.
Ntypes * toppar_data.
Ntypes;
331 const string section_name{
"NONBONDED_PARM_INDEX"};
332 const auto ¤t_section = toppar_map.find(section_name);
333 if (current_section != toppar_map.end()) {
334 success =
parse_section(get<1>(current_section->second), Ntype2d, toppar_data.
Cno,
336 get<0>(current_section->second) =
true;
338 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
343 const string section_name{
"RESIDUE_LABEL"};
344 const auto ¤t_section = toppar_map.find(section_name);
345 if (current_section != toppar_map.end()) {
347 toppar_data.
ResNames, section_name);
348 get<0>(current_section->second) =
true;
350 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
355 const string section_name{
"RESIDUE_POINTER"};
356 const auto ¤t_section = toppar_map.find(section_name);
357 if (current_section != toppar_map.end()) {
359 toppar_data.
Ipres, section_name);
360 get<0>(current_section->second) =
true;
362 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
367 const string section_name{
"BOND_FORCE_CONSTANT"};
368 const auto ¤t_section = toppar_map.find(section_name);
369 if (current_section != toppar_map.end()) {
371 toppar_data.
Rk, section_name);
372 get<0>(current_section->second) =
true;
374 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
379 const string section_name{
"BOND_EQUIL_VALUE"};
380 const auto ¤t_section = toppar_map.find(section_name);
381 if (current_section != toppar_map.end()) {
383 toppar_data.
Req, section_name);
384 get<0>(current_section->second) =
true;
386 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
391 const string section_name{
"ANGLE_FORCE_CONSTANT"};
392 const auto ¤t_section = toppar_map.find(section_name);
393 if (current_section != toppar_map.end()) {
395 toppar_data.
Tk, section_name);
396 get<0>(current_section->second) =
true;
398 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
403 const string section_name{
"ANGLE_EQUIL_VALUE"};
404 const auto ¤t_section = toppar_map.find(section_name);
405 if (current_section != toppar_map.end()) {
407 toppar_data.
Teq, section_name);
408 get<0>(current_section->second) =
true;
410 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
415 const string section_name{
"DIHEDRAL_FORCE_CONSTANT"};
416 const auto ¤t_section = toppar_map.find(section_name);
417 if (current_section != toppar_map.end()) {
419 toppar_data.
Pk, section_name);
420 get<0>(current_section->second) =
true;
422 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
427 const string section_name{
"DIHEDRAL_PERIODICITY"};
428 const auto ¤t_section = toppar_map.find(section_name);
429 if (current_section != toppar_map.end()) {
431 toppar_data.
Pn, section_name);
432 get<0>(current_section->second) =
true;
434 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
439 const string section_name{
"DIHEDRAL_PHASE"};
440 const auto ¤t_section = toppar_map.find(section_name);
441 if (current_section != toppar_map.end()) {
443 toppar_data.
Phase, section_name);
444 get<0>(current_section->second) =
true;
446 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
451 const string section_name{
"SCEE_SCALE_FACTOR"};
452 const auto ¤t_section = toppar_map.find(section_name);
453 if (current_section != toppar_map.end()) {
456 get<0>(current_section->second) =
true;
458 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
463 const string section_name{
"SCNB_SCALE_FACTOR"};
464 const auto ¤t_section = toppar_map.find(section_name);
465 if (current_section != toppar_map.end()) {
468 get<0>(current_section->second) =
true;
470 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
475 const string section_name{
"SOLTY"};
476 const auto ¤t_section = toppar_map.find(section_name);
477 if (current_section != toppar_map.end()) {
481 get<0>(current_section->second) =
true;
483 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
488 const int Nttyp = toppar_data.
Ntypes * (toppar_data.
Ntypes + 1) / 2;
489 const string section_name{
"LENNARD_JONES_ACOEF"};
490 const auto ¤t_section = toppar_map.find(section_name);
491 if (current_section != toppar_map.end()) {
492 success =
parse_section(get<1>(current_section->second), Nttyp, toppar_data.
Cn1,
494 get<0>(current_section->second) =
true;
496 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
501 const int Nttyp = toppar_data.
Ntypes * (toppar_data.
Ntypes + 1) / 2;
502 const string section_name{
"LENNARD_JONES_BCOEF"};
503 const auto ¤t_section = toppar_map.find(section_name);
504 if (current_section != toppar_map.end()) {
505 success =
parse_section(get<1>(current_section->second), Nttyp, toppar_data.
Cn2,
507 get<0>(current_section->second) =
true;
509 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
514 const string section_name{
"BONDS_INC_HYDROGEN"};
515 const auto ¤t_section = toppar_map.find(section_name);
516 if (current_section != toppar_map.end()) {
521 get<0>(current_section->second) =
true;
523 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
528 const string section_name{
"BONDS_WITHOUT_HYDROGEN"};
529 const auto ¤t_section = toppar_map.find(section_name);
530 if (current_section != toppar_map.end()) {
531 success =
parse_section(get<1>(current_section->second), toppar_data.Nbona,
532 {ref(toppar_data.BondAt1), ref(toppar_data.BondAt2),
533 ref(toppar_data.BondNum)},
535 get<0>(current_section->second) =
true;
537 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
542 const string section_name{
"ANGLES_INC_HYDROGEN"};
543 const auto ¤t_section = toppar_map.find(section_name);
544 if (current_section != toppar_map.end()) {
545 success =
parse_section(get<1>(current_section->second), toppar_data.Ntheth,
546 {ref(toppar_data.AngleHAt1), ref(toppar_data.AngleHAt2),
547 ref(toppar_data.AngleHAt3), ref(toppar_data.AngleHNum)},
549 get<0>(current_section->second) =
true;
551 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
556 const string section_name{
"ANGLES_WITHOUT_HYDROGEN"};
557 const auto ¤t_section = toppar_map.find(section_name);
558 if (current_section != toppar_map.end()) {
559 success =
parse_section(get<1>(current_section->second), toppar_data.Ntheta,
560 {ref(toppar_data.AngleAt1), ref(toppar_data.AngleAt2),
561 ref(toppar_data.AngleAt3), ref(toppar_data.AngleNum)},
563 get<0>(current_section->second) =
true;
565 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
570 const string section_name{
"DIHEDRALS_INC_HYDROGEN"};
571 const auto ¤t_section = toppar_map.find(section_name);
572 if (current_section != toppar_map.end()) {
573 success =
parse_section(get<1>(current_section->second), toppar_data.Nphih,
574 {ref(toppar_data.DihHAt1), ref(toppar_data.DihHAt2),
575 ref(toppar_data.DihHAt3), ref(toppar_data.DihHAt4),
576 ref(toppar_data.DihHNum)},
578 get<0>(current_section->second) =
true;
580 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
585 const string section_name{
"DIHEDRALS_WITHOUT_HYDROGEN"};
586 const auto ¤t_section = toppar_map.find(section_name);
587 if (current_section != toppar_map.end()) {
588 success =
parse_section(get<1>(current_section->second), toppar_data.Nphia,
589 {ref(toppar_data.DihAt1), ref(toppar_data.DihAt2),
590 ref(toppar_data.DihAt3), ref(toppar_data.DihAt4),
591 ref(toppar_data.DihNum)},
593 get<0>(current_section->second) =
true;
595 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
600 const string section_name{
"EXCLUDED_ATOMS_LIST"};
601 const auto ¤t_section = toppar_map.find(section_name);
602 if (current_section != toppar_map.end()) {
603 success =
parse_section(get<1>(current_section->second), toppar_data.Nnb,
604 toppar_data.ExclAt, section_name);
605 get<0>(current_section->second) =
true;
607 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
612 const string section_name{
"HBOND_ACOEF"};
613 const auto ¤t_section = toppar_map.find(section_name);
614 if (current_section != toppar_map.end()) {
615 success =
parse_section(get<1>(current_section->second), toppar_data.Nphb,
616 toppar_data.HB12, section_name);
617 get<0>(current_section->second) =
true;
619 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
624 const string section_name{
"HBOND_BCOEF"};
625 const auto ¤t_section = toppar_map.find(section_name);
626 if (current_section != toppar_map.end()) {
627 success =
parse_section(get<1>(current_section->second), toppar_data.Nphb,
628 toppar_data.HB10, section_name);
629 get<0>(current_section->second) =
true;
631 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
636 const string section_name{
"HBCUT"};
637 const auto ¤t_section = toppar_map.find(section_name);
638 if (current_section != toppar_map.end()) {
640 get<0>(current_section->second) =
true;
642 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
647 const string section_name{
"AMBER_ATOM_TYPE"};
648 const auto ¤t_section = toppar_map.find(section_name);
649 if (current_section != toppar_map.end()) {
650 success =
parse_section(get<1>(current_section->second), toppar_data.Natom,
651 toppar_data.AtomSym, section_name);
652 get<0>(current_section->second) =
true;
654 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
659 const string section_name{
"TREE_CHAIN_CLASSIFICATION"};
660 const auto ¤t_section = toppar_map.find(section_name);
661 if (current_section != toppar_map.end()) {
662 success =
parse_section(get<1>(current_section->second), toppar_data.Natom,
663 toppar_data.AtomTree, section_name);
664 get<0>(current_section->second) =
true;
666 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
671 const string section_name{
"JOIN_ARRAY"};
672 const auto ¤t_section = toppar_map.find(section_name);
673 if (current_section != toppar_map.end()) {
677 get<0>(current_section->second) =
true;
679 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
684 const string section_name{
"IROTAT"};
685 const auto ¤t_section = toppar_map.find(section_name);
686 if (current_section != toppar_map.end()) {
690 get<0>(current_section->second) =
true;
692 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
696 if (toppar_data.Ifbox > 0) {
699 const string section_name{
"SOLVENT_POINTERS"};
700 const auto ¤t_section = toppar_map.find(section_name);
701 if (current_section != toppar_map.end()) {
702 if (get<1>(current_section->second).size() == 3) {
703 toppar_data.Iptres = get<1>(current_section->second)[0].Int;
704 toppar_data.Nspm = get<1>(current_section->second)[1].Int;
705 toppar_data.Nspsol = get<1>(current_section->second)[2].Int;
706 toppar_data.AtomsPerMolecule.reserve(toppar_data.Nspm);
707 toppar_data.BoxDimensions.reserve(4);
708 toppar_data.CapInfo.reserve(1);
709 toppar_data.CapInfo2.reserve(4);
710 get<0>(current_section->second) =
true;
712 iout <<
iERROR <<
"Expect " << 3 <<
" fields but only "
713 << get<1>(current_section->second).size() <<
" are in "
714 << section_name.c_str() <<
"\n";
718 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
723 const string section_name{
"ATOMS_PER_MOLECULE"};
724 const auto ¤t_section = toppar_map.find(section_name);
725 if (current_section != toppar_map.end()) {
726 success =
parse_section(get<1>(current_section->second), toppar_data.Nspm,
727 toppar_data.AtomsPerMolecule, section_name);
728 get<0>(current_section->second) =
true;
730 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
735 const string section_name{
"BOX_DIMENSIONS"};
736 const auto ¤t_section = toppar_map.find(section_name);
737 if (current_section != toppar_map.end()) {
739 toppar_data.BoxDimensions, section_name);
740 get<0>(current_section->second) =
true;
742 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
746 if (toppar_data.Ifcap > 0) {
749 const string section_name{
"CAP_INFO"};
750 const auto ¤t_section = toppar_map.find(section_name);
751 if (current_section != toppar_map.end()) {
752 success =
parse_section(get<1>(current_section->second), 1, toppar_data.CapInfo,
754 get<0>(current_section->second) =
true;
756 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
761 const string section_name{
"CAP_INFO2"};
762 const auto ¤t_section = toppar_map.find(section_name);
763 if (current_section != toppar_map.end()) {
764 success =
parse_section(get<1>(current_section->second), 4, toppar_data.CapInfo2,
766 get<0>(current_section->second) =
true;
768 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
774 const string section_name{
"RADIUS_SET"};
775 const auto ¤t_section = toppar_map.find(section_name);
776 if (current_section != toppar_map.end()) {
777 success =
parse_section(get<1>(current_section->second), 1, toppar_data.RadiusSet,
779 get<0>(current_section->second) =
true;
781 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
786 const string section_name{
"RADII"};
787 const auto ¤t_section = toppar_map.find(section_name);
788 if (current_section != toppar_map.end()) {
789 success =
parse_section(get<1>(current_section->second), toppar_data.Natom,
790 toppar_data.Radii, section_name);
791 get<0>(current_section->second) =
true;
793 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
798 const string section_name{
"SCREEN"};
799 const auto ¤t_section = toppar_map.find(section_name);
800 if (current_section != toppar_map.end()) {
801 success =
parse_section(get<1>(current_section->second), toppar_data.Natom,
802 toppar_data.Screen, section_name);
803 get<0>(current_section->second) =
true;
805 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
810 const string section_name{
"IPOL"};
811 const auto ¤t_section = toppar_map.find(section_name);
812 if (current_section != toppar_map.end()) {
813 toppar_data.Ipol = get<1>(current_section->second)[0].Int;
814 toppar_data.Polarizability.reserve(toppar_data.Natom);
815 get<0>(current_section->second) =
true;
817 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
820 if (toppar_data.Ipol > 0) {
823 const string section_name{
"POLARIZABILITY"};
824 const auto ¤t_section = toppar_map.find(section_name);
825 if (current_section != toppar_map.end()) {
826 success =
parse_section(get<1>(current_section->second), toppar_data.Natom,
827 toppar_data.Polarizability, section_name);
828 get<0>(current_section->second) =
true;
830 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
834 if (toppar_data.IsCharmmFF) {
840 const string section_name{
"CHARMM_UREY_BRADLEY_COUNT"};
841 const auto ¤t_section = toppar_map.find(section_name);
842 if (current_section != toppar_map.end()) {
843 if (get<1>(current_section->second).size() == 2) {
844 toppar_data.Nub = get<1>(current_section->second)[0].Int;
845 toppar_data.Nubtypes = get<1>(current_section->second)[1].Int;
846 toppar_data.UreyBradleyAt1.reserve(toppar_data.Nub);
847 toppar_data.UreyBradleyAt2.reserve(toppar_data.Nub);
848 toppar_data.UreyBradleyNum.reserve(toppar_data.Nub);
849 toppar_data.UreyBradleyForceConstants.reserve(toppar_data.Nubtypes);
850 toppar_data.UreyBradleyEquilValues.reserve(toppar_data.Nubtypes);
851 get<0>(current_section->second) =
true;
853 iout <<
iERROR <<
"Expect " << 2 <<
" fields but only "
854 << get<1>(current_section->second).size() <<
" are in "
855 << section_name.c_str() <<
"\n";
859 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
864 const string section_name{
"CHARMM_UREY_BRADLEY"};
865 const auto ¤t_section = toppar_map.find(section_name);
866 if (current_section != toppar_map.end()) {
867 success =
parse_section(get<1>(current_section->second), toppar_data.Nub,
868 {ref(toppar_data.UreyBradleyAt1),
869 ref(toppar_data.UreyBradleyAt2),
870 ref(toppar_data.UreyBradleyNum)},
872 get<0>(current_section->second) =
true;
874 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
879 const string section_name{
"CHARMM_UREY_BRADLEY_FORCE_CONSTANT"};
880 const auto ¤t_section = toppar_map.find(section_name);
881 if (current_section != toppar_map.end()) {
882 success =
parse_section(get<1>(current_section->second), toppar_data.Nubtypes,
883 toppar_data.UreyBradleyForceConstants, section_name);
884 get<0>(current_section->second) =
true;
886 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
891 const string section_name{
"CHARMM_UREY_BRADLEY_EQUIL_VALUE"};
892 const auto ¤t_section = toppar_map.find(section_name);
893 if (current_section != toppar_map.end()) {
894 success =
parse_section(get<1>(current_section->second), toppar_data.Nubtypes,
895 toppar_data.UreyBradleyEquilValues, section_name);
896 get<0>(current_section->second) =
true;
898 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
903 const string section_name{
"CHARMM_NUM_IMPROPERS"};
904 const auto ¤t_section = toppar_map.find(section_name);
905 if (current_section != toppar_map.end()) {
906 if (get<1>(current_section->second).size() == 1) {
907 toppar_data.Nimphi = get<1>(current_section->second)[0].Int;
908 toppar_data.ImproperAt1.reserve(toppar_data.Nimphi);
909 toppar_data.ImproperAt2.reserve(toppar_data.Nimphi);
910 toppar_data.ImproperAt3.reserve(toppar_data.Nimphi);
911 toppar_data.ImproperAt4.reserve(toppar_data.Nimphi);
912 toppar_data.ImproperNum.reserve(toppar_data.Nimphi);
913 get<0>(current_section->second) =
true;
915 iout <<
iERROR <<
"Expect " << 1 <<
" fields but only "
916 << get<1>(current_section->second).size() <<
" are in "
917 << section_name.c_str() <<
"\n";
921 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
926 const string section_name{
"CHARMM_IMPROPERS"};
927 const auto ¤t_section = toppar_map.find(section_name);
928 if (current_section != toppar_map.end()) {
930 get<1>(current_section->second), toppar_data.Nimphi,
931 {ref(toppar_data.ImproperAt1), ref(toppar_data.ImproperAt2),
932 ref(toppar_data.ImproperAt3), ref(toppar_data.ImproperAt4),
933 ref(toppar_data.ImproperNum)},
935 get<0>(current_section->second) =
true;
937 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
943 const string section_name{
"CHARMM_NUM_IMPR_TYPES"};
944 const auto ¤t_section = toppar_map.find(section_name);
945 if (current_section != toppar_map.end()) {
946 toppar_data.NimprTypes = get<1>(current_section->second)[0].Int;
947 toppar_data.ImproperForceConstants.reserve(toppar_data.NimprTypes);
948 toppar_data.ImproperPhases.reserve(toppar_data.NimprTypes);
949 get<0>(current_section->second) =
true;
951 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
956 const string section_name{
"CHARMM_IMPROPER_FORCE_CONSTANT"};
957 const auto ¤t_section = toppar_map.find(section_name);
958 if (current_section != toppar_map.end()) {
959 success =
parse_section(get<1>(current_section->second), toppar_data.NimprTypes,
960 toppar_data.ImproperForceConstants, section_name);
961 get<0>(current_section->second) =
true;
963 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
968 const string section_name{
"CHARMM_IMPROPER_PHASE"};
969 const auto ¤t_section = toppar_map.find(section_name);
970 if (current_section != toppar_map.end()) {
971 success =
parse_section(get<1>(current_section->second), toppar_data.NimprTypes,
972 toppar_data.ImproperPhases, section_name);
973 get<0>(current_section->second) =
true;
975 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
980 const int Nttyp = toppar_data.Ntypes * (toppar_data.Ntypes + 1) / 2;
981 const string section_name{
"LENNARD_JONES_14_ACOEF"};
982 const auto ¤t_section = toppar_map.find(section_name);
983 if (current_section != toppar_map.end()) {
985 toppar_data.LJ14ACoefficients.reserve(Nttyp);
986 success =
parse_section(get<1>(current_section->second), Nttyp,
987 toppar_data.LJ14ACoefficients, section_name);
988 get<0>(current_section->second) =
true;
990 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
995 const int Nttyp = toppar_data.Ntypes * (toppar_data.Ntypes + 1) / 2;
996 const string section_name{
"LENNARD_JONES_14_BCOEF"};
997 const auto ¤t_section = toppar_map.find(section_name);
998 if (current_section != toppar_map.end()) {
1000 toppar_data.LJ14BCoefficients.reserve(Nttyp);
1001 success =
parse_section(get<1>(current_section->second), Nttyp,
1002 toppar_data.LJ14BCoefficients, section_name);
1003 get<0>(current_section->second) =
true;
1005 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
1014 const string section_name = toppar_data.IsCharmmFF
1015 ? string(
"CHARMM_CMAP_COUNT")
1016 : string(
"CMAP_COUNT");
1017 const auto ¤t_section = toppar_map.find(section_name);
1018 if (current_section != toppar_map.end()) {
1019 if (get<1>(current_section->second).size() == 2) {
1020 toppar_data.CMAPTermCount = get<1>(current_section->second)[0].Int;
1021 toppar_data.CMAPTypeCount = get<1>(current_section->second)[1].Int;
1022 toppar_data.HasCMAP =
true;
1023 toppar_data.CMAPResolution.reserve(toppar_data.CMAPTypeCount);
1025 toppar_data.CMAPParameter.reserve(toppar_data.CMAPTypeCount);
1026 toppar_data.CMAPIndex.reserve(toppar_data.CMAPTermCount * 6);
1027 get<0>(current_section->second) =
true;
1029 iout <<
iERROR <<
"Expect " << 2 <<
" fields but only "
1030 << get<1>(current_section->second).size() <<
" are in "
1031 << section_name.c_str() <<
"\n";
1035 toppar_data.HasCMAP =
false;
1038 if (toppar_data.HasCMAP) {
1041 const string section_name = toppar_data.IsCharmmFF
1042 ? string(
"CHARMM_CMAP_RESOLUTION")
1043 : string(
"CMAP_RESOLUTION");
1044 const auto ¤t_section = toppar_map.find(section_name);
1045 if (current_section != toppar_map.end()) {
1047 toppar_data.CMAPTypeCount, toppar_data.CMAPResolution,
1049 for (
int i = 0; i < toppar_data.CMAPTypeCount; ++i) {
1050 toppar_data.CMAPParameter.push_back(vector<_REAL>());
1052 get<0>(current_section->second) =
true;
1054 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
1059 for (
int i = 0; i < toppar_data.CMAPTypeCount; ++i) {
1060 string name_index = std::to_string(i + 1);
1061 if (name_index.size() < 2)
1062 name_index.insert(0,
"0");
1063 const string section_name =
1064 (toppar_data.IsCharmmFF ? string(
"CHARMM_CMAP_PARAMETER_")
1065 : string(
"CMAP_PARAMETER_")) +
1067 const auto ¤t_section = toppar_map.find(section_name);
1068 if (current_section != toppar_map.end()) {
1069 const int resolution =
1070 toppar_data.CMAPResolution[i] * toppar_data.CMAPResolution[i];
1071 toppar_data.CMAPParameter[i].reserve(resolution);
1072 success =
parse_section(get<1>(current_section->second), resolution,
1073 toppar_data.CMAPParameter[i], section_name);
1074 get<0>(current_section->second) =
true;
1076 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
1082 const string section_name = toppar_data.IsCharmmFF
1083 ? string(
"CHARMM_CMAP_INDEX")
1084 : string(
"CMAP_INDEX");
1085 const auto ¤t_section = toppar_map.find(section_name);
1086 if (current_section != toppar_map.end()) {
1088 toppar_data.CMAPTermCount * 6, toppar_data.CMAPIndex,
1090 get<0>(current_section->second) =
true;
1092 iout <<
iWARN <<
"Missing " << section_name.c_str() <<
" section!\n" <<
endi;
1100 for (
int i = 0; i < toppar_data.Natom; ++i) {
1101 if (i + 1 == toppar_data.Ipres[res+1])
1103 toppar_data.AtomRes.push_back(res);
1107 for (
auto it = toppar_map.begin(); it != toppar_map.end(); ++it) {
1108 if (get<0>(it->second) ==
false) {
1109 iout <<
iWARN <<
"Unknown section " << (it->first).c_str()
1110 <<
" found in the AMBER format paramter file!\n" <<
endi;
1117 bool success =
true;
1118 if (source.size() > 30) {
1119 result.
Natom = source[0].Int;
1120 result.
Ntypes = source[1].Int;
1121 result.
Nbonh = source[2].Int;
1122 result.
Mbona = source[3].Int;
1123 result.
Ntheth = source[4].Int;
1124 result.
Mtheta = source[5].Int;
1125 result.
Nphih = source[6].Int;
1126 result.
Mphia = source[7].Int;
1127 result.
Nhparm = source[8].Int;
1128 result.
Nparm = source[9].Int;
1129 result.
Nnb = source[10].Int;
1130 result.
Nres = source[11].Int;
1131 result.
Nbona = source[12].Int;
1132 result.
Ntheta = source[13].Int;
1133 result.
Nphia = source[14].Int;
1134 result.
Numbnd = source[15].Int;
1135 result.
Numang = source[16].Int;
1136 result.
Nptra = source[17].Int;
1137 result.
Natyp = source[18].Int;
1138 result.
Nphb = source[19].Int;
1139 result.
Ifpert = source[20].Int;
1140 result.
Nbper = source[21].Int;
1141 result.
Ngper = source[22].Int;
1142 result.
Ndper = source[23].Int;
1143 result.
Mbper = source[24].Int;
1144 result.
Mgper = source[25].Int;
1145 result.
Mdper = source[26].Int;
1146 result.
Ifbox = source[27].Int;
1147 result.
Nmxrs = source[28].Int;
1148 result.
Ifcap = source[29].Int;
1151 result.
Ncopy = 31 < source.size() ? source[31].Int : 0;
1154 const int Nttyp = result.
Ntypes * (result.
Ntypes + 1) / 2;
1161 result.
Cno.reserve(Ntype2d);
1168 result.
Pk.reserve(result.
Nptra);
1169 result.
Pn.reserve(result.
Nptra);
1174 result.
Cn1.reserve(Nttyp);
1175 result.
Cn2.reserve(Nttyp);
1217 vector<string> &destination,
const string §ion_name) {
1218 bool success =
true;
1219 if (
int(source.size()) != count) {
1220 iout <<
iERROR <<
"Expect " << count <<
" fields but only "
1221 << source.size() <<
" are in "
1222 << section_name.c_str() <<
"\n";
1225 for (
int i = 0; i < count; ++i) {
1226 destination.push_back(source[i].String);
1228 iout <<
iINFO <<
"Complete reading " << section_name.c_str() <<
" with " << count <<
" fields.\n" <<
endi;
1234 vector<int> &destination,
const string §ion_name) {
1235 bool success =
true;
1236 if (
int(source.size()) != count) {
1237 iout <<
iERROR <<
"Expect " << count <<
" fields but only "
1238 << source.size() <<
" are in "
1239 << section_name.c_str() <<
"\n";
1242 for (
int i = 0; i < count; ++i) {
1243 destination.push_back(source[i].Int);
1245 iout <<
iINFO <<
"Complete reading " << section_name.c_str() <<
" with " << count <<
" fields.\n" <<
endi;
1251 vector<_REAL> &destination,
const string §ion_name) {
1252 bool success =
true;
1253 if (
int(source.size()) != count) {
1254 iout <<
iERROR <<
"Expect " << count <<
" fields but only "
1255 << source.size() <<
" are in "
1256 << section_name.c_str() <<
"\n";
1259 for (
int i = 0; i < count; ++i) {
1260 destination.push_back(source[i].
Real);
1262 iout <<
iINFO <<
"Complete reading " << section_name.c_str() <<
" with " << count <<
" fields.\n" <<
endi;
1268 const vector<FortranData> &source,
const int &count,
1269 std::initializer_list<std::reference_wrapper<vector<string>>> destination,
1270 const string §ion_name) {
1271 bool success =
true;
1272 const int total_size = destination.size() * count;
1273 if (
int(source.size()) != total_size) {
1274 iout <<
iERROR <<
"Expect " << total_size <<
" fields but only "
1275 << source.size() <<
" are in "
1276 << section_name.c_str() <<
"\n";
1279 const int stride = destination.size();
1280 for (
int i = 0; i < count; ++i) {
1282 for (
auto it_j = destination.begin(); it_j != destination.end(); ++it_j) {
1283 it_j->get().push_back(source[stride * i + j].String);
1287 iout <<
iINFO <<
"Complete reading " << section_name.c_str() <<
" with " << total_size <<
" fields.\n" <<
endi;
1293 const vector<FortranData> &source,
const int &count,
1294 std::initializer_list<std::reference_wrapper<vector<int>>> destination,
1295 const string §ion_name) {
1296 bool success =
true;
1297 const int total_size = destination.size() * count;
1298 if (
int(source.size()) != total_size) {
1299 iout <<
iERROR <<
"Expect " << total_size <<
" fields but only "
1300 << source.size() <<
" are in "
1301 << section_name.c_str() <<
"\n";
1304 const int stride = destination.size();
1305 for (
int i = 0; i < count; ++i) {
1307 for (
auto it_j = destination.begin(); it_j != destination.end(); ++it_j) {
1308 it_j->get().push_back(source[stride * i + j].Int);
1312 iout <<
iINFO <<
"Complete reading " << section_name.c_str() <<
" with " << total_size <<
" fields.\n" <<
endi;
1318 const vector<FortranData> &source,
const int &count,
1319 std::initializer_list<std::reference_wrapper<vector<_REAL>>> destination,
1320 const string §ion_name) {
1321 bool success =
true;
1322 const int total_size = destination.size() * count;
1323 if (
int(source.size()) != total_size) {
1324 iout <<
iERROR <<
"Expect " << total_size <<
" fields but only "
1325 << source.size() <<
" are in "
1326 << section_name.c_str() <<
"\n";
1329 const int stride = destination.size();
1330 for (
int i = 0; i < count; ++i) {
1332 for (
auto it_j = destination.begin(); it_j != destination.end(); ++it_j) {
1333 it_j->get().push_back(source[stride * i + j].
Real);
1337 iout <<
iINFO <<
"Complete reading " << section_name.c_str() <<
" with " << total_size <<
" fields.\n" <<
endi;
1345 iout <<
iINFO <<
"=== Start reading AMBER parm7 format file " << filename <<
" ===\n" <<
endi;
1350 iout <<
iERROR <<
"Failed to read AMBER parm7 file: " << filename <<
" at read_amber_parm_stage2\n" <<
endi;
1351 NAMD_die(
"Error in Ambertoppar readparm(const char *filename)\n");
1354 iout <<
iERROR <<
"Failed to read AMBER parm7 file: " << filename <<
" at read_amber_parm_stage1\n" <<
endi;
1355 NAMD_die(
"Error in Ambertoppar readparm(const char *filename)\n");
1357 iout <<
iINFO <<
"=== Finish reading AMBER parm7 format file. ===\n" <<
endi;
vector< string > RadiusSet
std::ostream & iINFO(std::ostream &s)
vector< string > ResNames
vector< string > AtomTree
vector< _REAL > ScnbScaleFactor
std::ostream & endi(std::ostream &s)
unordered_map< string, tuple< bool, vector< FortranData > > > AmberTopparMap
std::ostream & iWARN(std::ostream &s)
void split_string_by_specifiers(const string &source, const vector< FortranFormatSpecifier > &specifiers, vector< FortranData > &destination)
This function splits a string by a set of fortran format specifiers.
vector< int > AtomNumbers
bool read_amber_parm_stage2(AmberTopparMap &toppar_map, Ambertoppar &toppar_data)
Read an AmberTopparMap into an Ambertoppar struct for NAMD.
void parse_fortran_format(const std::string &str, FortranFormatSpecifier &specifier)
Parse a single fortran format specifier.
bool read_amber_parm_stage1(const char *filename, AmberTopparMap &toppar_map)
Read an AMBER topology file into an AmberTopparMap.
bool parse_pointer(const vector< FortranData > &source, Ambertoppar &result)
void split_string(const std::string &data, const std::string &delim, std::vector< std::string > &dest)
void NAMD_die(const char *err_msg)
vector< string > AtomNames
__device__ __forceinline__ float dot(const float3 v1, const float3 v2)
Ambertoppar readparm(const char *filename)
std::ostream & iERROR(std::ostream &s)
vector< _REAL > SceeScaleFactor
bool parse_section(const vector< FortranData > &source, const int &count, vector< string > &destination, const string §ion_name)
Copy data from source to destination.