tools/gen_irr_diag.c

   1 #include <stdio.h>
   2 #include <stdlib.h>
   3 #include <sys/time.h>
   4 #include <sys/resource.h>
   5 #include <unistd.h>
   6 #include <string.h>
   7 #include <strings.h>
   8
   9 int  nChmOpts = 0;
  10 char rxt_tbl[5][1000][128];
  11 char chm_scheme[5][128];
  12 int  rxt_cnt[5];
  13
  14 void strip_blanks( char *instring, char *outstring )
  15 {
  16   int i, j, slen;
  17   char *c;
  18
  19   slen = strlen( instring );
  20   strcpy( outstring,instring );
  21   c = instring;
  22   for( i=0,j=0; i < slen; i++ )
  23   {
  24     if( strncmp( c, " ", 1 ) )
  25       strncpy( (outstring+(j++)), c, 1 );
  26     c++;
  27   }
  28   c = index( outstring,'\n' );
  29   *c = '\0';
  30 }
  31
  32 int AppendReg( char *chem_opt, int ndx )
  33 {
  34    int Nrxt;
  35    int i,j;
  36    int slen;
  37    char *stradj = "abcdef";
  38    char *strt, *end;
  39    char *token;
  40    char *wstrg1;
  41    char path[256];
  42    char fname[256];
  43    char inln[1024],winln[1024],s[1024];
  44    char rxtstr[128];
  45    char rxtstr_tbl[1000][128];
  46    char buffer[128];
  47    char rxtsym[128];
  48    FILE *fp_eqn, *fp_reg;
  49
  50    strcpy( fname,chem_opt );
  51    slen = sprintf( path, "chem/KPP/mechanisms/%s/%s.eqn",fname,fname );
  52    fprintf(stderr,"Using file:\n");
  53    fprintf(stderr,"%s\n", path);
  54
  55 //
  56 //  open input *eqn file
  57 //
  58    if( (fp_eqn = fopen( path, "r" )) == NULL )
  59    {
  60       fprintf(stderr,"Can not open %s for reading\n", path);
  61       return(-1);
  62    }
  63    fprintf(stderr,"\n");
  64 //
  65 //  open output registry file
  66 //
  67    if( access( "Registry/registry.irr_diag",F_OK ) )
  68      fp_reg = fopen( "Registry/registry.irr_diag", "w" );
  69    else
  70      fp_reg = fopen( "Registry/registry.irr_diag", "a" );
  71
  72    if( fp_reg == NULL ) {
  73      fprintf(stderr,"Can not open registry.irr_diag for writing\n");
  74      return(-2);
  75    }
  76    strcpy( buffer,"\"Integrated Reaction Rate\"  \"\"");
  77    fprintf(fp_reg,"state  real  -  ikjf  irr_diag_%s  -  -  -  - %s\n",fname,buffer);
  78
  79    Nrxt = 0;
  80 //
  81 //  loop over input lines
  82 //
  83    while( fgets( inln, 1024, fp_eqn ) != NULL )
  84    {
  85       if( strncmp( inln, "//", 2 ) && strncmp( inln, "#", 1 ) )
  86       {
  87 //
  88 //  strip blanks from reaction string
  89 //
  90         strip_blanks( inln, inln );
  91 //
  92 //  concatentate lines?
  93 //
  94         if( rindex( inln,';' ) == NULL ) {
  95           do {
  96             fgets( winln, 1024, fp_eqn );
  97             strip_blanks( winln, winln );
  98             strcat( inln,winln );
  99           } while( rindex( inln,';' ) == NULL );
 100         }
 101         strt = strchr( inln,'{' );
 102         while( strt != NULL )
 103         {
 104           end = strchr( strt,'}' );
 105           slen = 0;
 106           if( end != NULL )
 107             slen = end - strt;
 108           if( slen > 0 )
 109           {
 110             for( i = 0; i < slen; i++ )
 111               if( !strncmp( strt+i,"=",1 ) ) *(strt+i) = '$';
 112           }
 113           strt += slen;
 114           strt = strchr( strt,'{' );
 115         }
 116
 117         strt = strchr( inln, '}' );
 118         if( strt == NULL ) continue;
 119         strt++;
 120         end  = rindex( strt, '=' );
 121         if( end != NULL )
 122         {
 123           *end= '\0';
 124           slen = strlen( strt );
 125           wstrg1 = strt;
 126 //
 127 //  string to upper case
 128 //
 129           slen = strlen( wstrg1 );
 130           for( i=0; i <= slen; i++ )
 131             wstrg1[i] = toupper( wstrg1[i] );
 132 //
 133 //  remove all text between {} pair including delimiters
 134 //
 135           strt = strchr( wstrg1,'{' );
 136           if( strt != NULL )
 137           {
 138             end = index( wstrg1,'}' );
 139             if( end != NULL )
 140             {
 141               char c[128];
 142               *strt = '\0';
 143               strcpy( c,wstrg1 );
 144               end++;
 145               strcat( c,end );
 146               strcpy( wstrg1,c );
 147             }
 148           }
 149
 150           strcpy( rxtstr_tbl[Nrxt],wstrg1 );
 151 //
 152 //  check for unique reaction string
 153 //
 154           if( Nrxt > 0 )
 155           {
 156             int Nmatch = 0;
 157             for( i = 0; i < Nrxt; i++)
 158               if( !strcmp( wstrg1,rxtstr_tbl[i] ) )
 159                 Nmatch++;
 160
 161             if( Nmatch > 0 )
 162             {
 163               Nmatch--;
 164               strcat( wstrg1,"_" );
 165               strncat( wstrg1,stradj+Nmatch,1 );
 166             }
 167           }
 168
 169           strcpy( rxtstr,wstrg1 );
 170           strcpy( rxtsym,wstrg1 );
 171           strcpy( rxt_tbl[ndx][Nrxt],wstrg1 );
 172 //
 173 //  change + to _
 174 //
 175           for( i=0; i < slen; i++ )
 176           {
 177             if( ! strncmp( rxtsym+i, "+", 1 ) )
 178               strncpy( rxtsym+i, "_", 1 );
 179           }
 180           strcat( rxtsym,"_IRR" );
 181 //
 182 //  form output line
 183 //
 184 //        fprintf(fp_reg,"state  real  %s  ikjf  irr_diag_%s  1  -  rh9  \"%s\"  \"%s Integrated Reaction Rate\"  \"molecules/cm^3/s\"\n",rxtsym,fname,rxtstr,rxtstr);
 185           fprintf(fp_reg,"state  real  %s  ikjf  irr_diag_%s  1  -  rh9  \"%s\"  \"%s Integrated Reaction Rate\"  \"ppmv\"\n",rxtsym,fname,rxtsym,rxtstr);
 186           Nrxt++;
 187         }
 188       }
 189    }
 190
 191    nChmOpts++;
 192    rxt_cnt[ndx] = Nrxt;
 193    strcpy( chm_scheme[ndx],chem_opt );
 194
 195    fclose(fp_eqn);
 196    fclose(fp_reg);
 197
 198    return(0);
 199 }
 200
 201 int irr_diag_scalar_indices( char *dirname )
 202 {
 203    int Nrxt;
 204    int i, j;
 205    int first, flush, s1;
 206    char fname[256];
 207    char line[132];
 208    char piece[132];
 209    char *blank = "                                                           ";
 210    FILE *fp_inc;
 211
 212 //
 213 //  open output inc file
 214 //
 215    sprintf( fname, "inc/scalar_indices_irr_diag_decls.inc");
 216    fp_inc = fopen( fname, "w" );
 217
 218    if( fp_inc == NULL ) {
 219      fprintf(stderr,"Can not open %s for writing\n",fname);
 220      return(-2);
 221    }
 222    fprintf( fp_inc," \n");
 223    fprintf( fp_inc,"  INTEGER, PARAMETER :: nchm_opts = %d\n",nChmOpts);
 224    fprintf( fp_inc,"  INTEGER            :: chm_opt_ndx\n");
 225    fprintf( fp_inc," \n");
 226    fprintf( fp_inc,"  INTEGER            :: chm_opts_cnt(nchm_opts)\n");
 227    fprintf( fp_inc,"  INTEGER            :: chm_opts_ndx(nchm_opts)\n");
 228    fprintf( fp_inc,"  CHARACTER(len=32)  :: chm_opts_name(nchm_opts)\n");
 229    for( i = 0,j = 0; i < nChmOpts; i++ ) {
 230      if( rxt_cnt[i] > j )
 231        j = rxt_cnt[i];
 232    }
 233    fprintf( fp_inc,"  CHARACTER(len=64), TARGET  :: rxtsym(%d,nchm_opts)\n",j);
 234    fprintf( fp_inc," \n");
 235    fclose(fp_inc);
 236 //
 237 //  open output inc file
 238 //
 239    sprintf( fname, "inc/scalar_indices_irr_diag.inc");
 240    fp_inc = fopen( fname, "w" );
 241
 242    if( fp_inc == NULL ) {
 243      fprintf(stderr,"Can not open %s for writing\n",fname);
 244      return(-2);
 245    }
 246
 247    fprintf( fp_inc," \n");
 248    sprintf( line,"  chm_opts_cnt(:nchm_opts) = (/ ");
 249    for( i = 0; i < nChmOpts; i++ ) {
 250      if( i == 0 )
 251        sprintf( piece," %d",rxt_cnt[i]);
 252      else
 253        sprintf( piece," ,%d",rxt_cnt[i]);
 254      strcat( line,piece );
 255    }
 256    strcat( line," /)\n" );
 257    fprintf( fp_inc,line );
 258    fprintf( fp_inc," \n");
 259
 260    for( i = 0; i < nChmOpts; i++ ) {
 261      sprintf( line,"  chm_opts_name(%d) = '%s'\n",i+1,chm_scheme[i]);
 262      fprintf( fp_inc,line );
 263    }
 264    fprintf( fp_inc," \n");
 265
 266    sprintf( line,"  chm_opts_ndx(:nchm_opts) = (/ ");
 267    for( i = 0; i < nChmOpts; i++ ) {
 268      if( i == 0 )
 269        sprintf( piece,"%s_kpp",chm_scheme[i]);
 270      else
 271        sprintf( piece," ,%s_kpp",chm_scheme[i]);
 272      strcat( line,piece );
 273    }
 274    strcat( line," /)\n" );
 275    fprintf( fp_inc,line );
 276    fprintf( fp_inc," \n");
 277
 278    for( i = 0; i < nChmOpts,rxt_cnt[i] > 0; i++ ) {
 279      for( j = 0; j < rxt_cnt[i]; j++ ) {
 280        sprintf( line,"     rxtsym(%d,%d) = '%s'\n",j+1,i+1,rxt_tbl[i][j]);
 281        fprintf( fp_inc,"%s",line);
 282      }
 283      fprintf( fp_inc," \n");
 284    }
 285 /*
 286    fprintf( fp_inc,"  IF( model_config_rec%%irr_opt(idomain) == 1 ) THEN\n");
 287    fprintf( fp_inc,"    CALL nl_get_chem_opt( idomain,chem_opt )\n");
 288    fprintf( fp_inc,"    DO chm_opt_ndx = 1,nchm_opts\n");
 289    fprintf( fp_inc,"      IF( chem_opt == chm_opts_ndx(chm_opt_ndx) ) THEN\n");
 290    fprintf( fp_inc,"        EXIT\n");
 291    fprintf( fp_inc,"      ENDIF\n");
 292    fprintf( fp_inc,"    ENDDO\n");
 293    fprintf( fp_inc,"    IF( chm_opt_ndx > nchm_opts ) THEN\n");
 294    fprintf( fp_inc,"      write(err_mes,*) 'IRR not supported for chem option ',chem_opt\n");
 295    fprintf( fp_inc,"      CALL wrf_error_fatal( trim(err_mes) )\n");
 296    fprintf( fp_inc,"    ENDIF\n");
 297    fprintf( fp_inc,"  ENDIF\n");
 298 */
 299    fclose(fp_inc);
 300
 301    return(0);
 302 }