gepmain.c

来自「GESPI 2.0动态系统模拟工具 」· C语言 代码 · 共 458 行

#include "copyleft.h"

/*
    GEPASI - a simulator of metabolic pathways and other dynamical systems
    Copyright (C) 1989, 1992, 1993  Pedro Mendes
*/

/*************************************/
/*                                   */
/*         main control block        */
/*                                   */
/*        Zortech C/C++ 3.0 r4       */
/*          MICROSOFT C 6.00         */
/*          Visual C/C++ 1.0         */
/*           QuickC/WIN 1.0          */
/*             ULTRIX cc             */
/*              GNU gcc              */
/*                                   */
/*   (include here compilers that    */
/*   compiled GEPASI successfully)   */
/*                                   */
/*************************************/


#ifndef MSDOS
#define COMMON
#endif

#ifdef _WINDOWS
#define COMMON
#endif

#if (_MSC_VER >= 610)                           /* QuickC/WIN or C/C++ 7.0                              */
 #define fcloseall _fcloseall
 #define unlink _unlink
 #define putenv _putenv
#endif
#if ((_MSC_VER == 600)||(_MSC_VER==700))        /* MICROSOFT C 6.00 & 7.00              */
 #pragma comment( exestr, "Compiled on " __DATE__ " at " __TIME__ )
#endif

#include <stdio.h>
#include <string.h>
#include <time.h>
#include <math.h>
#include <float.h>
#include <setjmp.h>
#include <signal.h>
#include <stdlib.h>
/*#include <conio.h>*/

#ifdef _WINDOWS
#include <io.h>
#include <fcntl.h>
#endif

#include "pmu.h"                                    /* several utilities     */
#include "globals.h"                                /* global symbols        */
#include "globvar.h"                                /* global variables      */

#ifndef _WINDOWS
#include "except.c"                                 /* exception error handl.*/
#endif

#include "r250.h"                                                                       /* pseudo-random numbers */
#include "matrix.h"                                 /* matrix algebra        */
#include "kinetics.h"                               /* rate eq. & deriv.     */
#include "rates.h"                                  /* net rates of met.     */
#include "gauss.h"                                  /* structural properties */
#include "newton.h"                                 /* steady state          */
#include "metconan.h"                               /* metabolic control     */
#include "lsoda.h"
#include "dynamic.h"                                    /* transient behaviour   */
#include "datab.h"
#include "getall.h"
#include "gepout.h"                                 /* data output           */
#include "heapchk.h"

#ifdef _ZTC
#define MEM_DEBUG 1
#include "mem.h"
#else
#define mem_malloc malloc
#define mem_free free
#define mem_realloc realloc
#endif

/* tidy up memory and disk                                                      */

void tidy( void )
{
 /* if -d switch delete the file with list of simulations       */
 if( delete ) unlink( listfile );

 /* free up all the dynamic memory blocks                                       */
 FreeMem();
}

/* terminate gepasi due to errors */

void term_gepasi( void )
{
 tidy();
 fprintf(stderr,errormsg[2]);
#ifdef _MSC_VER
 fcloseall();
#endif
#ifdef _ZTC
 mem_term();
#endif
 exit( 1 );                                         /* signal error          */
}

/* resets the state of the internal metabolites         */

void reset_state( double *keep )
{
 int i;

 for( i=0; i<indmet; i++ ) x[i] = x0[i];
 if( &endtime != keep ) endtime = options.endtime;
}

/* sets the values of the linked variables                      */

void set_linked( void )
{
 int i;

 for( i=0; i<nlinks; i++)
  switch( sparam[lindex[i]].operation )
  {
   case 1: *(sparam[lindex[i]].var) = sparam[lindex[i]].factor *
				      *(sparam[lindex[i]].linkedto);
	   break;
  }
}


/* main simulation procedure                                    */

int simulate( int p, int simple )
{
 int i, ccfs, diverge;
 time_t itime, ftime;

#ifndef _WINDOWS
  /* Save stack environment for return in case of error.                */
  /*  First time through, jmpret is 0, so true conditional              */
  /*  is executed. If an error occurs, jmpret will be set               */
  /*  to -1 and false conditional will be executed.                             */
  jmpret = setjmp( mark );
  if( jmpret == 0 )
  {
#endif

  /* reset concentration of external metabolites*/
  for( i=indmet; i<totmet; i++ ) x[i] = x0[i];

  /* set the values of the linked parameters    */
  if( nlinks ) set_linked();

  /* initialize the conserved moiety concentrations                                             */
  init_moiety();
  /* print the first part of the report file    */
  if( options.txt ) report_init( p );

  /* dynamic simulation */
  if( options.dyn )
  {
   time(&itime);
   dynamic();                           /* numeric integration */
   time(&ftime);
   ftime -= itime;

   if ( ! options.debug ) set_cursor(cx, cy);

   if( simple )
    printf( "\n  - integration completed in %ld s, ", ftime );
   else printf( "." );

   /* calculate the fluxes      */
   calcfluxes( x );

   /* print integration results in report file  */
   if( options.txt ) report_dyn();
   if( options.dat  && (!options.datss) ) dat_values();
  }
  else
   for( i=0; i<nmetab; i++) x[i] = 10 * options.hrcz;

  /* steady-state solution      */
  if ( options.ss )
  {
   if( simple )
    if( ! options.debug ) printf( "\n  - ss solution " );
   time( &itime );
   if( (diverge = newton()) )                               /* steady state                 */
   {
    endtime = 1e+10;                                        /* new endtime                  */
    int_to_ss();                                            /* new numeric integr.  		*/
    diverge = newton();										/* try again st. state			*/
   }
   time(&ftime);
   ftime -= itime;

   if ( options.debug ) printf( "\n\n  - ss solution in %ld s\n", ftime );
   else
    if( simple ) printf( "in %ld s", ftime );
    else printf( "." );

   if ( diverge == N_OK )
   {
    calcfluxes( xss );
    calc_tt();                                        /* calc. trans. times  */
    calc_Dxv();                                       /* calc. elasticities  */
    if ( ( ccfs = calc_Gamma() ) != MCA_NOCC )        /* calc. [] cnt.coeff. */
     calc_C();                                        /* calc. J cnt.coeff.  */
   }

   if( options.txt ) report_ss( diverge, ccfs );

   /* scale the control coefficients if necessary       */
   if( (options.dat) && (!options.datmca ) ) scale_coef( );
   /* copy xss to x                                                                     */
   memcpy( (void *) x, (void *) xss, sizeof( x ) );

   /* output data in columnar file      */
   if( options.dat && options.datss && (diverge==N_OK) ) dat_values();
  }

  if( options.txt && options.debug ) report_tech();
  return 0;
#ifndef _WINDOWS
  }
  else
  {
   fpcheck();
   if( options.txt ) report_error();
   if( options.dat ) dat_error();
  }
#endif
}

int hypergrid( int index, int p )
{
 int i;

 *(sparam[scindex[index]].var) = sparam[scindex[index]].low;
 if( options.dyn ) reset_state(sparam[scindex[index]].var);
 for( i=0; i<sparam[scindex[index]].dens+1; i++ )
 {
  if( index != 0 )
  {
   hypergrid( index-1, p );
  }
  else
  {
   printf("\n %lu ", ++iter );
   simulate( p, 0 );
  }
  if( sparam[scindex[index]].log )
   *(sparam[scindex[index]].var) = sparam[scindex[index]].low *
				   pow( 10, sparam[scindex[index]].ampl /
					    sparam[scindex[index]].dens * (i+1)
				      );
  else
   *(sparam[scindex[index]].var) = sparam[scindex[index]].low +
				   (i+1) *
				   sparam[scindex[index]].ampl /
				   sparam[scindex[index]].dens;
  if( options.dyn ) reset_state(sparam[scindex[index]].var);
  if( (totscan > 1) && (index == totscan - 1) ) dat_nl();
 }
 return 0;
}

int randomfill( int p )
{
 unsigned long i;
 int j;

 for( i=0; i< (unsigned long) options.scandens; i++ )
 {
  for( j=0; j<totscan; j++ )
  {
   if( sparam[scindex[j]].log )
    *(sparam[scindex[j]].var) = sparam[scindex[j]].low * pow( 10,
	   sparam[scindex[j]].ampl * dr250() );
   else
    *(sparam[scindex[j]].var) = sparam[scindex[j]].low +
	   sparam[scindex[j]].ampl * dr250();
   if( options.dyn ) reset_state(sparam[scindex[j]].var);
  }
  printf("\n %lu ", i );
  simulate( p, 0 );
 }
 return 0;
}

void process_sim( int p )
{
 int nRc, err;

 get_cursor(&cx, &cy);

 printf("\n\nprocessing %s: ",filename[p]);

 /* read data for this simulation       */
 if ((nRc = ReadSim(filename[p])))
 {
  out_error( nRc );
  term_gepasi();
 }
 if( options.dyn )
 {
  reset_state( NULL );
  strcpy( dyname, filename[p] );                         /* setup .dyn filename    */
  fixext( dyname, '.', ".dyn" );
  unlink( dyname );
 }

 /* set up indexes to effectors and pointers to rate equations  */
 eff_rateq();

 if( options.dat  && (!options.append) ) dat_fileovr();
 if( (options.dat) && (options.dattit == 0) ) dat_titles();

 err = 0;

 switch( options.scan )
 {
  case 0:  if( options.dyn ) reset_state( NULL );
	   simulate( p, 1 );
	   break;
  case 1:  iter=0;
	   hypergrid( totscan-1, p );
	   break;
  case 2:  randomfill( p );
	   break;
  default: if( options.dyn ) reset_state( NULL );
	   simulate( p, 1 );
 }

 if( (options.dat) & (options.dattit == 1) ) dat_titles();
}

/* the main body */

int main(int argc, char *argv[])
{
 int p;

#ifdef _WINDOWS
 struct _wsizeinfo ws;
#if _MSC_VER == 800
 ws._version = _winver;
#else
 ws._version = _WINVER;
#endif
 ws._type = _WINSIZEMAX;
 _wsetsize( fileno( stdout ), &ws );
 _wsetscreenbuf( fileno( stdout ), _WINBUFINF );
 _wabout( "GEPASI 2.0 for MS-Windows\n       (C) P. Mendes 1992/93" );
 _wsetexit( _WINEXITPERSIST );
#endif

#ifdef _ZTC
mem_init();
#endif

 if( getenv( "TZ" ) == NULL )
  if( putenv( "TZ=GMT0BST" ) == -1 )
   printf( "\n  * unable to set GMT, using Pacific time zone (PST)\n" );

 /* randomize the rand seed with the time function      */
 r250_init( (int) time( NULL ) );

 delete=0;

 /* check command line arguments and act accordingly    */
 if ( argc > 1 )
 {
  if ( ! strcmp( "-h", argv[1] ) )
  {
   fprintf( stderr, "\n%s (%s)\n\n", gepasi, version_no );
   fprintf( stderr, "This program is free software; you can redistribute it and/or modify\n");
   fprintf( stderr, "it under the terms of the GNU General Public License as published by\n");
   fprintf( stderr, "the Free Software Foundation; either version 1, or (at your option)\n");
   fprintf( stderr, "any later version.\n\n");
   fprintf( stderr, "This program is distributed in the hope that it will be useful,\n");
   fprintf( stderr, "but WITHOUT ANY WARRANTY; without even the implied warranty of\n");
   fprintf( stderr, "MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the\n");
   fprintf( stderr, "GNU General Public License for more details.\n\n");
   fprintf( stderr, "You should have received a copy of the GNU General Public License\n");
   fprintf( stderr, "along with this program, in the file COPYING; if not, write to the\n");
   fprintf( stderr, "Free Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.\n\n");
   return 0;
  }
  else
  {
   strcpy( listfile, argv[1] );
   getlist( listfile );              /* file with list of simulations   */
   if( argc > 2 )
    if( !strcmp( argv[2], "-d" ) )   /* -d switch ?                                             */
     delete = 1;
  }
 }
 else getlist( (char *) NULL );      /* get list of simuls. from console*/

 /* allocate memory for several data structures         */
 if( MemAlloc() == -1 )
 {
  fprintf( stderr, "\nnot enough memory for GEPASI!" );
  term_gepasi();
 }

 /* setup ktype database                                        */
 InitDataBase();

#ifndef _WINDOWS
 /* set up floating point error handler         */
 /* this line generates a warning in Microsoft C compilers, just ignore it! */
 if( signal( SIGFPE, fphandler ) == SIG_ERR )
 {
  printf("  * ERROR - Couldn't set floating point error handler\n" );
  term_gepasi();
 }
#endif

cls();

 /* print screen banner                                     */
 printf("%s\n%s\n",gepasi,version_no);

 /* print no.of simulations*/
 printf("\n\n%d simulation%s",nrunfiles, nrunfiles==1 ? "" : "s ");

 /* execute all simulations     */
 for(p=0;p<nrunfiles;p++)
 {
  process_sim( p );
 }
 /* clean memory and delete file if necessary                           */
 tidy();

#ifdef _ZTC
  mem_term();
#endif

 /* make sure a carriage return is sent to stdout                       */
 printf("\n");

 /* sound signal and screen message of program termination      */
 beep();
 fprintf(stderr,"\nAll done.\n\n");
}