/*============================================================================*/
/* FILENAME: WRITE_MM2.C (WRITE_MM2.O)
/* PURPOSE: TO WRITE OUT THE CURRENT GEOMETRY AS AN MM2 INPUT FILE.
/* WRITTEN: M.V.GRIESHABER
/* MODIFICATIONS:
/*	22 JULY 1993 MVG (INCREASED RUN TIME TO MAX LIMIT).
/*	12 JULY 1993 MVG (CORRECTED SOME MM3 ONLY ATOM TYPES).
/*	8 FEBRUARY 1993 MVG (EXTENDED SYBYL ATOM TYPES, NEW MM3 TYPES).
/*============================================================================*/

#include "utility.h"
#include "newgeo.h"

#define MAX_RINGS 100			/* MAX RINGS IN A MOLECULE. */
#define MAX_CONNECTED_ATOM_LISTS 20	/* ONE LIST PER LINE. */
#define MAX_ATTACHED_ATOM_LISTS 20	/* REALLY ONLY ONE LIST, BUT SPREAD ACROSS THIS MANY LINES. */
#define MAX_ATOMS_PER_LIST 16		/* REALLY MAX ATOMS NUMBERS PER LINE. */
#define TERMINATION_ENERGY 0.00008	/* OPTIMIZATION TERMINATION (KCAL/MOLE). */
#define MM2_LONE_PAIR_TYPE 20		/* FOR AUTO ADD OF LONE PAIRS. */
#define MM2_LONE_PAIR_RADIUS 0.6	/* FOR AUTO ADD OF LONE PAIRS. */

typedef struct				/* DEFINE A LIST OF ATOMS. */
   {
   int natoms;				/* NUMBER OF ATOMS IN LIST. */
   int atoms[MAX_ATOMS];		/* ATOM INDEX. */
   } ATOM_LIST;

int nrings;				/* NUMBER OF RINGS IN LIST. */
ATOM_LIST rings[MAX_RINGS];		/* UNIQUE RINGS IN MOLECULE. */

int ring_candidate_atom[MAX_ATOMS];	/* T/F: CONSIDER ATOM IN RING SEARCH. */
int ring_check_atom[MAX_ATOMS];		/* T/F: CHECK ATOM FOR RING. */

int write_mm2_input(geometry_type,outfile)
/*============================================================================*/
/* PURPOSE: TO WRITE OUT THE CURRENT GEOMETRY AS AN MM2 OR MM3 INPUT FILE.
*/
   int geometry_type;		/* SPECIFIES MM2 OR MM3 TYPE FILE. */
   FILE *outfile;		/* MM2 FORMAT OUTPUT. */
   {
   int target_atom;		/* ATOM TO LOOK FOR. */
   int current_atom;		/* ATOM WE ARE PRESENTLY AT. */
   int previous_atom;		/* ATOM WE WERE PREVIOUSLY AT. */
   ATOM_LIST path;		/* PATH FOLLOWED DURING SEARCH. */
   int natom;			/* WHICH ATOM. */
   int aromatics_exist;		/* T/F: WHETHER AROMATIC ATOMS EXIST IN MOLECULE. */
   int nlone_pairs;		/* TOTAL ADD_LONE_PAIRS FLAGS. */
   int nconnected_atom_lists;	/* TOTAL CONNECTED ATOM LISTS. */
   int connected_atom_lists[MAX_CONNECTED_ATOM_LISTS][MAX_ATOMS_PER_LIST];
   int nattached_atoms;		/* TOTAL ATTACHED ATOMS IN MOLECULE. */
   int attached_atom_lists[MAX_ATTACHED_ATOM_LISTS][MAX_ATOMS_PER_LIST];
   char line[150];		/* MULTIPURPOSE STRING. */
   char default_title[100];	/* DEFAULT TITLE. */
   int nlist;			/* WHICH LIST. */
   int status;			/* RETURN STATUS FROM USER INTERFACE ROUTINES. */
   struct
      {
      int atom_type;		/* MM2 ATOM TYPE. */
      int aromatic;		/* T/F: AROMATIC ATOM TYPE. */
      int add_lone_pairs;	/* T/F: WHETHER TO AUTO ADD LONE PAIRS. */
      } mm2_data[MAX_ATOMS];	/* SAVED VALUES FROM MM2_ATOM_TYPE ROUTINE. */

   void find_rings();		/* FIND RINGS FROM CURRENT TO TARGET ATOM. */
   void mark_ring_atoms();	/* DETERMINE WHICH ATOMS TO CONSIDER FOR RINGS. */
   void filter_atom_type();	/* FILTER SPECIFIED ATOM TYPE. */

   /* FILTER OUT LONE PAIRS (AUTOMATICALLY ADDED AS NEEDED) AND DUMMY ATOMS. */
   /* WARN USER IF ANY DUMMY ATOMS WERE DELETED. */
   filter_atom_type("LP");
   for (natom=0; natom0) fprintf(outfile,"%5d",connected_atom_lists[nlist][natom]);
         else fprintf(outfile,"     ");
         }
      fprintf(outfile,"\n");
      }

   /* LINE 4 (ATTACHED ATOM LIST). */
   for (nlist=0; nlist<((nattached_atoms*2-1)/MAX_ATOMS_PER_LIST)+1; nlist++)
      {
      for (natom=0; natom0) fprintf(outfile,"%5d",attached_atom_lists[nlist][natom]);
         else fprintf(outfile,"     ");
         }
      fprintf(outfile,"\n");
      }

   /* LINE 5 (COORDINATES). */
   for (natom=0; natom=3 BONDS. */

   /* TO START, ASSUME ALL ATOMS ARE BOTH "RING" AND "RING CHECK" CANDIDATES. */
   for (natom=0; natom0)
      {
      /* HAVE A (POSSIBLY UNIQUE) RING, LOG IT. */
      log_ring(path);
      }
   else if (already_visited(current_atom,path))
      {
      /* ALREADY VISITED THIS ATOM, CAN'T CROSS PATH SO NO FURTHER ACTION. */
      }
   else
      {
      /* FOLLOW EACH BOND THAT IS AVAILABLE FROM THE CURRENT ATOM. */
      for (nbond=0; nbond3)
               mm2_type=60;	/* MM3 ONLY: P(V). */
            else
               mm2_type=25;	/* PHOSPHINE. */
            break;

         case 13:		/* SYBYL TYPE H (HYDROGEN). */
            if (connected_to("5",index,-1) || connected_to("6",index,-1) || connected_to("19",index,-1))
               mm2_type=23;	/* AMINE. */
            else if (connected_to("28",index,-1))
               mm2_type=28;	/* AMIDE. */
            else if (connected_to("31",index,-1))
               mm2_type=48;	/* AMMONIUM. */
            else if (connected_to("8 2 9",index,-1))
               mm2_type=24;	/* ACID. */
            else if (connected_to("8 2",index,-1))
               mm2_type=73;	/* ENOL. */
            else if (connected_to("8",index,-1))
               mm2_type=21;	/* ALCOHOL. */
            else if (geometry_type==TYPE_MM3_INPUT &&
               (connected_to("10",index,-1) || connected_to("18",index,-1) ||
               connected_to("29",index,-1) || connected_to("30",index,-1)))
               mm2_type=44;	/* MM3 ONLY: THIOL. */
            else
               mm2_type=5;	/* PLAIN HYDROGEN. */
            break;

         case 14:		/* SYBYL TYPE Br (BROMINE). */
            mm2_type=13;	/* BROMINE. */
            break;

         case 15:		/* SYBYL TYPE Cl (CHLORINE). */
            mm2_type=12;	/* CHLORINE. */
            break;

         case 16:		/* SYBYL TYPE F (FLUORINE). */
            mm2_type=11;	/* FLUORINE. */
            break;

         case 17:		/* SYBYL TYPE I (IODINE). */
            mm2_type=14;	/* IODINE. */
            break;

         case 18:		/* SYBYL TYPE S.2 (SP2 SULFUR). */
            if (in_non_c3_ring(5,index))
               {
               mm2_type=42;	/* THIOPHENE. */
               *set_aromatic_flag=TRUE;
               }
            else if (connected_to("2",index,-1) && intco[index].bond_cnt==1)
               mm2_type=17;	/* THIOCARBONYL SULFUR. */
            else
               mm2_type=15;	/* SP3 SULFUR. */
            break;

         case 19:		/* SYBYL TYPE N.pl3 (PLANAR TRIGONAL NITROGEN). */
            if (in_non_c3_ring(5,index) || in_non_c3_ring(6,index))
               {
               if (intco[index].bond_cnt==3)
                  {
                  mm2_type=40;	/* PYROLE NITROGEN. */
                  *set_aromatic_flag=TRUE;
                  }
               else
                  {
                  mm2_type=37;	/* PYRIDINE NITROGEN. */
                  *set_aromatic_flag=TRUE;
                  *add_lone_pairs=TRUE;
                  }
               }
            else if (connected_to("9",index,-1)==2)
               mm2_type=46;	/* NITRO. */
            else if (geometry_type==TYPE_MM3_INPUT &&
               (connected_to("6",index,-1) || connected_to("11",index,-1) || connected_to("19",index,-1)))
               {
               if (connected_to("8",index,-1))
                  mm2_type=109;	/* AZOXY NITROGEN. */
               else
                  mm2_type=107;	/* AZO NITROGEN. */
               }
            else
               {
               mm2_type=8;	/* SP3 NITROGEN. */
               *add_lone_pairs=TRUE;
               }
            break;

         case 20:		/* SYBYL TYPE LP (LONE PAIR). */
            /* ALL LONE PAIRS WERE INITIALLY DELETED AND REPLACED AS NEEDED. */
            mm2_type=99;	/* THIS CASE SHOULD NEVER BE EXECUTED. */
            break;

         case 21:		/* SYBYL TYPE Na (SODIUM). */
         case 22:		/* SYBYL TYPE K (POTASSIUM). */
         case 23:		/* SYBYL TYPE Ca (CALCIUM). */
         case 24:		/* SYBYL TYPE Li (LITHIUM). */
         case 25:		/* SYBYL TYPE Al (ALUMINUM). */
            /* NO CORRESPONDING MM2 TYPE. */
            printf("Warning: Atom #%d (SYBYL type %d) has no corresponding MM2 atom type.\n",index+1,intco[index].source_type);
            printf("         MM2 type 99 has been used; manual editing of file is required.\n");
            mm2_type=99;
            break;

         case 26:		/* SYBYL TYPE Du (DUMMY). */
            /* ALL DUMMIES WERE PREVIOUSLY DELETED. */
            mm2_type=99;	/* THIS CASE SHOULD NEVER BE EXECUTED. */
            break;

         case 27:		/* SYBYL TYPE Si (SILICON). */
            mm2_type=19;	/* SILANE. */
            break;

         case 28:		/* SYBYL TYPE N.am (AMIDE NITROGEN). */
            mm2_type=9;		/* AMIDE. */
            break;

         case 29:		/* SYBYL TYPE S.O (SULFOXIDE). */
            mm2_type=17;	/* SULFOXIDE. */
            break;

         case 30:		/* SYBYL TYPE S.O2 (SULFONE). */
            mm2_type=18;	/* SULFONE. */
            break;

         case 31:		/* SYBYL TYPE N.4 (AMMONIUM). */
            mm2_type=39;	/* AMMONIUM. */
            break;

         case 32:		/* SYBYL TYPE O.co2 (CARBOXYLATE). */
            mm2_type=47;	/* CARBOXYLATE. */
            break;

         case 33:		/* SYBYL TYPE C.cat (CARBOCATION). */
            mm2_type=30;	/* CARBOCATION. */
            break;

         case 906:		/* SYBYL TYPE Co.oh (COBALT). */
            mm2_type=66;	/* COBALT (III). */
            break;

         case 910:		/* SYBYL TYPE Fe (IRON). */
            mm2_type=62;	/* IRON (III). */
            break;

         case 911:		/* SYBYL TYPE Ni (NICKEL). */
            mm2_type=63;	/* NICKEL (II). */
            break;

         case 915:		/* SYBYL TYPE Ge (GERMANIUM). */
            mm2_type=31;	/* GERMANIUM. */
            break;

         case 917:		/* SYBYL TYPE Se (SELENIUM). */
            mm2_type=34;	/* SELENIUM. */
            break;

         case 918:		/* SYBYL TYPE Kr (KRYPTON). */
            mm2_type=54;	/* KRYPTON. */
            break;

         case 931:		/* SYBYL TYPE Sn (TIN). */
            mm2_type=32;	/* TIN. */
            break;

         case 933:		/* SYBYL TYPE Te (TELLURIUM). */
            mm2_type=35;	/* TELLURIUM. */
            break;

         case 934:		/* SYBYL TYPE Xe (XENON). */
            mm2_type=55;	/* XENON. */
            break;

         case 948:		/* SYBYL TYPE Pb (LEAD). */
            mm2_type=33;	/* LEAD. */
            break;

         case 984:		/* SYBYL TYPE He (HELIUM). */
            mm2_type=51;	/* HELIUM. */
            break;

         case 987:		/* SYBYL TYPE Ne (NEON). */
            mm2_type=52;	/* NEON. */
            break;

         case 988:		/* SYBYL TYPE Mg (MAGNESIUM). */
            mm2_type=59;	/* MAGNESIUM. */
            break;

         case 990:		/* SYBYL TYPE Ar (ARGON). */
            mm2_type=53;	/* ARGON. */
            break;

         default:		/* UNKNOWN SYBYL ATOM TYPE. */
            printf("Warning: Atom #%d (SYBYL type %d) is an unknown SYBYL atom type.\n",index+1,intco[index].source_type);
            printf("         MM2/MM3 type 99 has been used; manual editing of file is required.\n");
            mm2_type=99;
            break;
         }
         return(mm2_type);
      }
   else
      {
      /* DO A CRUDE CONVERSION FROM A NON SYBYL SOURCE TYPE TO MM2 ATOM TYPE. */
      printf("Error: Non SYBYL source type to MM2/MM3 conversion not yet implemented.\n");
      return(-1);
      }
   }

int in_three_member_ring(index)
/*============================================================================*/
/* PURPOSE: DETERMINE WHETHER OR NOT THE SPECIFIED ATOM IS CONTAINED IN A
/*	THREE MEMBER RING.
*/
   int index;	/* ATOM INDEX. */
   {
   int nring;	/* WHICH RING. */
   int natom;	/* WHICH ATOM. */

   /* EXAMINE EACH THREE MEMBER RING FOR INCLUSION OF INDEX ATOM. */
   for (nring=0; nring0) nconnections+=connected_to(chain,new_atom,current_atom);
         else nconnections++;
         }
      }

   /* RETURN NUMBER OF DIFFERENT SUBPATHS FROM THE CURRENT ATOM THAT MEET */
   /* THE REST OF THE CHAIN ORDER. */
   return(nconnections);
   }

int make_connected_atom_lists(connected_atom_lists)
/*============================================================================*/
/* PURPOSE: CREATE THE "CONNECTED" ATOM LISTS FOR MM2 FROM THE CURRENT MOLECULE.
/*	A "CONNECTED" ATOM IS ANY ATOM WITH AT LEAST TWO BONDS TO IT.  ALL BONDS
/*	IN THE MOLECULE MUST SHOW UP IN THE CONNECTED LISTS (EXCEPT THOSE TO
/*	"ATTACHED" ATOMS).
*/
   int connected_atom_lists[MAX_CONNECTED_ATOM_LISTS][MAX_ATOMS_PER_LIST];
   {
   int natom;			/* WHICH ATOM. */
   int nbond;			/* WHICH BOND. */
   int nlist;			/* WHICH LIST. */
   int first_unused_atom;	/* FIRST ATOM WITH UNUSED BONDS. */
   int bonds_used[MAX_ATOMS][MAX_BONDS_PER_ATOM];	/* T/F; BOND USED IN LIST. */

   void mark_all_attached_atom_bonds();	/* MARK ATTACHED ATOM BONDS AS USED. */
   void make_connected_atom_list();	/* FORM SINGLE CHAIN FROM ROOT ATOM. */

   for (nlist=0; nlist1)
            {
            nlist++;
            connected_atom_lists[nlist][0]=natom+1;	/* INDEX FROM ONE FOR MM2. */
            break;	/* THERE CAN ONLY BE ONE SUCH ATOM. */
            }
         }
      }

   return(nlist+1);	/* RETURN TOTAL NUMBER OF CONNECTION LISTS USED. */
   }

void mark_all_attached_atom_bonds(bonds_used)
/*============================================================================*/
/* PURPOSE: MARK ALL BONDS TO OR FROM AN ATTACHED ATOM (ANY ATOM WITH ONLY ONE
/*	BOND IS AN ATTACHED ATOM) AS USED SO THEY WILL NOT BE CONSIDERED IN
/*	THE SEARCH FOR CONNECTED ATOMS.
*/
   int bonds_used[MAX_ATOMS][MAX_BONDS_PER_ATOM];	/* T/F; BOND USED IN LIST. */
   {
   int natom;			/* WHICH ATOM. */

   void mark_bond_used();	/* MARK BOND AND RECIPROCAL BOND AS USED. */

   for (natom=0; natom= 1 UNUSED BOND. */
   int bonds_used[MAX_ATOMS][MAX_BONDS_PER_ATOM];	/* T/F; BOND USED IN LIST. */
   {
   int natom;			/* WHICH ATOM. */
   int nbond;			/* WHICH BOND. */

   /* SCAN EACH BOND OF EACH ATOM FOR THE FIRST UNUSED BOND. */
   for (natom=0; natom=MAX_ATOMS_PER_LIST) return;	/* LIST IS FULL, STOP HERE. */

   for (nbond=0; nbond MAX_ATOMS_PER_LIST)
            {
            nlist++;	/* FILLED UP ONE LIST, ROLL OVER TO THE NEXT ONE. */
            natoms_in_list=1;
            if (nlist+1 > MAX_ATTACHED_ATOM_LISTS)
               {
               printf("Error: Maximum number of attached atom lists exceeded, file is incorrect.\n");
               return(nlist);	/* TRUNCATE AT THE MAXIMUM. */
               }
            }
         /* ADD ROOT AND ATTACHED ATOM (ADJUST TO INDEX FROM ONE FOR MM2). */
         attached_atom_lists[nlist][(natoms_in_list-1)*2]=intco[natom].bonds[0]+1;
         attached_atom_lists[nlist][(natoms_in_list-1)*2+1]=natom+1;
         }
      }

   return(nattached_atoms);	/* RETURN TOTAL NUMBER OF ATTACHED ATOMS. */
   }