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