forester/archive/RIO/others/hmmer/src/hmmcalibrate-pvm.c

   1 /************************************************************
   2  * HMMER - Biological sequence analysis with profile HMMs
   3  * Copyright (C) 1992-1999 Washington University School of Medicine
   4  * All Rights Reserved
   5  *
   6  *     This source code is distributed under the terms of the
   7  *     GNU General Public License. See the files COPYING and LICENSE
   8  *     for details.
   9  ************************************************************/
  10
  11 #ifdef HMMER_PVM
  12
  13 /* hmmcalibrate-pvm.c
  14  * SRE, Tue Aug 18 15:19:28 1998
  15  * Redesigned for better parallelization: SRE, Wed Dec  1 09:48:58 1999
  16  *
  17  * Design:
  18  *   Initialization:
  19  *       receive parameters of random sequence synthesis, and an HMM.
  20  *       send an OK signal to the master.
  21  *
  22  *   Main loop:
  23  *       receive work packet: # of seqs to make
  24  *       Synthesize and score # seqs
  25  *       send results: # raw scores.
  26  *
  27  *   Termination:
  28  *       master sends a shutdown signal instead of a work packet.
  29  *
  30  * PVM slave for hmmcalibrate.
  31  * RCS $Id: hmmcalibrate-pvm.c,v 1.1.1.1 2005/03/22 08:34:12 cmzmasek Exp $
  32  */
  33
  34 #include <stdio.h>
  35 #include <stdlib.h>
  36 #include <float.h>
  37 #include <pvm3.h>
  38
  39 #include "version.h"
  40 #include "structs.h"            /* data structures, macros, #define's   */
  41 #include "config.h"             /* compile-time configuration constants */
  42 #include "funcs.h"              /* function declarations                */
  43 #include "globals.h"            /* alphabet global variables            */
  44 #include "squid.h"              /* general sequence analysis library    */
  45 #include "stopwatch.h"          /* CPU timing routines                  */
  46
  47 static void leave_pvm(void);
  48
  49 int
  50 main(void)
  51 {
  52   int      master_tid;          /* PVM TID of our master */
  53   int      slaveidx;            /* my slave index (0..nslaves-1) */
  54   struct plan7_s *hmm;          /* HMM to calibrate, sent from master */
  55   char    *seq;                 /* synthetic random sequence */
  56   char    *dsq;                 /* digitized seq */
  57   int      len;                 /* length of seq */
  58   float   *sc;                  /* scores of seqs */
  59   int      seed;                /* random number seed */
  60   int      nsample;             /* number of seqs to sample */
  61   int      fixedlen;            /* if nonzero, fixed length of seq */
  62   float    lenmean;             /* Gaussian mean length of seq */
  63   float    lensd;               /* Gaussian length std. dev. for seq */
  64   float    randomseq[MAXABET];  /* iid frequencies of residues */
  65   float    p1;
  66   int      alphatype;           /* alphabet type, hmmAMINO or hmmNUCLEIC    */
  67   int      idx;
  68   int      code;
  69   Stopwatch_t stopwatch;       /* CPU timings */
  70
  71   /* Register leave_pvm() cleanup function so any exit() call
  72    * first calls pvm_exit().
  73    */
  74   if (atexit(leave_pvm) != 0) {
  75     pvm_exit(); Die("slave couldn't register leave_pvm()");
  76   }
  77
  78   /*****************************************************************
  79    * initialization.
  80    * Master broadcasts the problem to us:
  81    *    an HMM;
  82    *    parameters of the HMM calibration.
  83    * We send back:
  84    *    an OK flag, and our RELEASE, for some sanity checking.
  85    ******************************************************************/
  86
  87   StopwatchStart(&stopwatch);
  88
  89   master_tid = pvm_parent();    /* who's our master? */
  90
  91   pvm_recv(master_tid, HMMPVM_INIT);
  92   pvm_upkfloat(&lenmean,  1, 1); /* mean length of random seqs    */
  93   pvm_upkfloat(&lensd,    1, 1); /* std. dev. of random seq len   */
  94   pvm_upkint(&fixedlen,   1, 1); /* if non-zero, override lenmean */
  95   pvm_upkint(&alphatype,  1, 1); /* alphabet type, hmmAMINO or hmmNUCLEIC */
  96   pvm_upkint(&seed,       1, 1); /* random number seed */
  97   SetAlphabet(alphatype);        /* must set alphabet before reading HMM! */
  98   hmm = PVMUnpackHMM();
  99   if (hmm == NULL) Die("oh no, the HMM never arrived");
 100
 101   P7DefaultNullModel(randomseq, &p1);
 102   P7Logoddsify(hmm, TRUE);
 103
 104   /* tell the master we're OK and ready to go (or not)
 105    */
 106   code = HMMPVM_OK;
 107   pvm_initsend(PvmDataDefault);
 108   pvm_pkint(&code, 1, 1);
 109   PVMPackString(RELEASE);
 110   pvm_send(master_tid, HMMPVM_RESULTS);
 111
 112   /*****************************************************************
 113    * Main loop.
 114    * Receive: a number of sequences we're supposed to do.
 115    *          If we receive a 0, we have no work, so wait for shutdown;
 116    *          if we receive a -1, shut down.
 117    *****************************************************************/
 118   slaveidx = -1;
 119   for (;;)
 120     {
 121       pvm_recv(master_tid, HMMPVM_WORK);
 122       pvm_upkint(&nsample,  1, 1);
 123       pvm_upkint(&idx,      1, 1);
 124
 125       if (nsample == 0)  continue;  /* go into stasis */
 126       if (nsample == -1) break;     /* shut down      */
 127
 128       if (slaveidx == -1) {     /* first time: set id, seed sre_random  */
 129         slaveidx = idx;
 130         sre_srandom(seed+idx);  /* unique seed in current PVM   */
 131       }
 132
 133       sc = MallocOrDie(sizeof(float) * nsample);
 134       for (idx = 0; idx < nsample; idx++)
 135         {
 136                                 /* choose length of random sequence */
 137           if (fixedlen) len = fixedlen;
 138           else do len = (int) Gaussrandom(lenmean, lensd); while (len < 1);
 139                                 /* generate it */
 140           seq = RandomSequence(Alphabet, randomseq, Alphabet_size, len);
 141           dsq = DigitizeSequence(seq, len);
 142           SQD_DPRINTF2(("slave %d seq: %d : %20.20s...\n", slaveidx, len, seq));
 143
 144           if (P7ViterbiSize(len, hmm->M) <= RAMLIMIT)
 145             sc[idx] = P7Viterbi(dsq, len, hmm, NULL);
 146           else
 147             sc[idx] = P7SmallViterbi(dsq, len, hmm, NULL);
 148
 149           free(seq);
 150           free(dsq);
 151         }
 152
 153       /* Return output to master, some of which is sanity checking.
 154        *   1. our slave index.
 155        *   2. how many seqs we simulated.
 156        *   3. the array of scores we got, so the master can stuff
 157        *      them into a histogram.
 158        */
 159       pvm_initsend(PvmDataDefault);
 160       pvm_pkint(&slaveidx, 1, 1);
 161       pvm_pkint(&nsample,  1, 1);
 162       pvm_pkfloat(sc, nsample,1);
 163       pvm_send(master_tid, HMMPVM_RESULTS);
 164
 165       /* cleanup
 166        */
 167       free(sc);
 168     }
 169
 170   /***********************************************
 171    * Cleanup, return.
 172    ***********************************************/
 173
 174   FreePlan7(hmm);
 175   StopwatchStop(&stopwatch);
 176
 177   /* tell the master we heard his shutdown signal, and
 178    * give him our CPU times; then exit.
 179    */
 180   pvm_initsend(PvmDataDefault);
 181   pvm_pkint(&slaveidx, 1, 1);
 182   StopwatchPVMPack(&stopwatch);
 183   pvm_send(master_tid, HMMPVM_RESULTS);
 184
 185   return 0;     /* pvm_exit() is called by atexit() registration. */
 186 }
 187
 188 /* Function: leave_pvm()
 189  *
 190  * Purpose:  Cleanup function, to deal with crashes. We register
 191  *           this function using atexit() so it gets called before
 192  *           the slave dies.
 193  */
 194 void leave_pvm(void)
 195 {
 196   SQD_DPRINTF1(("slave leaving PVM.\n"));
 197   pvm_exit();
 198 }
 199
 200 #else /* if HMMER_PVM not defined: include a dummy */
 201
 202 #include <stdio.h>
 203 int main(void)
 204 {
 205   printf("hmmcalibrate-pvm disabled. PVM support was not compiled into HMMER.\n");
 206   exit(0);
 207 }
 208
 209 #endif