00001
00002
00003
00004 #include "cddefines.h"
00005 #include "physconst.h"
00006 #include "optimize.h"
00007 #include "doppvel.h"
00008 #include "stopcalc.h"
00009 #include "abund.h"
00010 #include "geometry.h"
00011 #include "dense.h"
00012 #include "plot.h"
00013 #include "grid.h"
00014 #include "rfield.h"
00015 #include "dynamics.h"
00016 #include "magnetic.h"
00017 #include "trace.h"
00018 #include "h2.h"
00019 #include "mole.h"
00020 #include "hmi.h"
00021 #include "rt.h"
00022 #include "thermal.h"
00023 #include "opacity.h"
00024 #include "atomfeii.h"
00025 #include "called.h"
00026 #include "extinc.h"
00027 #include "wind.h"
00028 #include "hextra.h"
00029 #include "iterations.h"
00030 #include "radius.h"
00031 #include "input.h"
00032 #include "assertresults.h"
00033 #include "phycon.h"
00034 #include "fudgec.h"
00035 #include "version.h"
00036 #include "neutrn.h"
00037 #include "conv.h"
00038 #include "parse.h"
00039
00040
00041 void ParseCommands(void)
00042 {
00043 char chCard[INPUT_LINE_LENGTH],
00044 chKey2[3],
00045 chKey3[4],
00046 chKey4[5],
00047 chKey5[6];
00048
00049 bool lgDSet,
00050 lgEOF,
00051 lgEOL,
00052 lgNu2;
00053 bool lgStop ,
00054 lgStop_not_enough_info;
00055
00056 long int i,
00057 j,
00058 nqh;
00059 long nInitFile=0;
00060
00061 float a,
00062 ar1,
00063 teset,
00064 z;
00065 double thickness;
00066
00067 DEBUG_ENTRY( "ParseCommands()" );
00068
00069
00070 abund.lgAbnSolar = true;
00071
00072
00073 plotCom.lgPlotON = false;
00074 plotCom.nplot = 0;
00075
00076
00077
00078
00079
00080
00081 lgDSet = false;
00082
00083 radius.Radius = 0.;
00084 radius.rinner = 0.;
00085 nqh = 0;
00086 rfield.nspec = 0;
00087
00088
00089 input.lgSetNoBuffering = false;
00090
00091
00092 InitAssertResults();
00093
00094 for( i=0; i < LIMSPC; i++ )
00095 {
00096 strcpy( rfield.chRSpec[i], "UNKN" );
00097 }
00098 optimize.nparm = 0;
00099
00100
00101
00102 if( optimize.lgTrOpt )
00103 {
00104
00105
00106 if( optimize.nTrOpt == optimize.nOptimiz )
00107 {
00108 trace.lgTrace = true;
00109 called.lgTalk = true;
00110 trace.lgTrOvrd = true;
00111 fprintf( ioQQQ, " READR turns on trace from optimize option.\n" );
00112 }
00113 }
00114
00115
00116 if( version.nBetaVer > 0 && called.lgTalk )
00117 {
00118 fprintf( ioQQQ,
00119 "\n This is %s (beta %ld) of Cloudy, and is intended for testing only.\n\n",
00120 version.chVersion, version.nBetaVer );
00121 }
00122
00123 if( called.lgTalk )
00124 {
00125
00126 for( i=0; i<57; i++ )
00127 {
00128 fprintf(ioQQQ,"%c",' ');
00129 }
00130 fprintf( ioQQQ, "Cloudy %s\n", version.chVersion);
00131
00132 for( i=0; i<57; i++ )
00133 {
00134 fprintf(ioQQQ,"%c",' ');
00135 }
00136 fprintf( ioQQQ, "www.nublado.org\n\n");
00137
00138 for( i=0; i<23; i++ )
00139 {
00140 fprintf(ioQQQ,"%c",' ');
00141 }
00142
00143
00144 fprintf( ioQQQ, "**************************************");
00145 fprintf( ioQQQ, "%8.8s", version.chDate);
00146 fprintf( ioQQQ, "*************************************\n");
00147
00148 for( i=0; i<23; i++ )
00149 {
00150 fprintf(ioQQQ,"%c",' ');
00151 }
00152 fprintf( ioQQQ, "*");
00153
00154 for( i=0; i<81; i++ )
00155 {
00156 fprintf(ioQQQ,"%c",' ');
00157 }
00158 fprintf( ioQQQ, "*\n");
00159 }
00160
00161
00162
00163
00164
00165
00166 input.iReadWay = 1;
00167
00168
00169
00170 input_init();
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181
00182
00183
00184 input_readarray(chCard,&lgEOF);
00185
00186
00187 while( !lgEOF && chCard[0] != ' ' )
00188 {
00189
00190
00191
00192 if( called.lgTalk && !nMatch("HIDE",input.chCARDCAPS) )
00193 {
00194 fprintf( ioQQQ, " * ");
00195 i=0;
00196 while( chCard[i]!='\0' )
00197 {
00198 fprintf(ioQQQ,"%c",chCard[i]);
00199 ++i;
00200 }
00201
00202 while( i<80 )
00203 {
00204 fprintf(ioQQQ,"%c",' ');
00205 ++i;
00206 }
00207 fprintf(ioQQQ,"*\n");
00208 }
00209
00210
00211 caps( chCard );
00212
00213
00214
00215 strncpy( chKey2 , chCard , 2 );
00216 chKey2[2] = '\0';
00217
00218
00219 strncpy( chKey3 , chCard , 3 );
00220 chKey3[3] = '\0';
00221
00222
00223 strncpy( chKey4 , chCard , 4 );
00224 chKey4[4] = '\0';
00225
00226
00227 strncpy( chKey5 , chCard , 5 );
00228 chKey5[5] = '\0';
00229
00230
00231 if( nMatch("VARY",chCard) )
00232 {
00233 optimize.lgVarOn = true;
00234 if( optimize.nparm + 1 > LIMPAR )
00235 {
00236 fprintf( ioQQQ, " Too many VARY lines entered; the limit is%4ld\n",
00237 LIMPAR );
00238 puts( "[Stop in ParseCommands]" );
00239 cdEXIT(EXIT_FAILURE);
00240 }
00241 }
00242
00243 else
00244 {
00245 optimize.lgVarOn = false;
00246 }
00247
00248
00249
00250
00251 if( chCard[0]=='C' && (chCard[1]==' ' || chCard[1]== 0) )
00252 {
00253
00254
00255 i = 1;
00256 }
00257
00258
00259 else if( strcmp(chKey4,"ABSO") == 0 )
00260 {
00261
00262 ParseAbsMag(chCard,&nqh);
00263 }
00264
00265 else if( strcmp(chKey3,"AGE") == 0 )
00266 {
00267
00268 ParseAge(chCard);
00269 }
00270
00271 else if( strcmp(chKey4,"AGN ") == 0 )
00272 {
00273
00274 ParseAgn(chCard);
00275 }
00276
00277 else if( strcmp(chKey4,"ABUN") == 0 )
00278 {
00279
00280
00281 ParseAbundances(chCard,lgDSet);
00282
00283 abund.lgAbnSolar = false;
00284 }
00285
00286 else if( strcmp(chKey4,"APER") == 0 )
00287 {
00288
00289
00290
00291
00292
00293
00294
00295
00296
00297
00298
00299 if( nMatch("SLIT",chCard) )
00300 {
00301
00302
00303 geometry.iEmissPower = 1;
00304 }
00305 else if( nMatch("BEAM",chCard) )
00306 {
00307
00308
00309
00310 geometry.iEmissPower = 0;
00311 }
00312 else
00313 {
00314 fprintf( ioQQQ, " One of the keywords SLIT or BEAM must appear.\n" );
00315 fprintf( ioQQQ, " Sorry.\n" );
00316 puts( "[Stop in ParseCommands]" );
00317 cdEXIT(EXIT_FAILURE);
00318 }
00319 }
00320
00321 else if( strcmp(chKey4,"ASSE") == 0 )
00322 {
00323
00324 ParseAssertResults();
00325 }
00326
00327 else if( strcmp(chKey4,"ATOM") == 0 )
00328 {
00329
00330 if( nMatch("FEII",chCard) || nMatch("FE II",chCard) )
00331 {
00332
00333 ParseAtomFeII(chCard);
00334 }
00335
00336 else if( nMatch("H-LI",chCard) )
00337 {
00338
00339 ParseAtomHLike(chCard);
00340 }
00341
00342 else if( nMatch("HE-L",chCard) )
00343 {
00344
00345 ParseAtomHeLike(chCard);
00346 }
00347
00348 else if( nMatch("ROTO",chCard) )
00349 {
00350
00351 fprintf(ioQQQ,"PROBLEM - the atom rotor command is now the ATOM CO command. "
00352 "Please use that instead. I will accept this for now but may not in future versions.\n");
00353 ParseAtomCO(chCard);
00354 }
00355
00356 else if( nMatch(" CO ",chCard) )
00357 {
00358
00359
00360 ParseAtomCO(chCard);
00361 }
00362
00363 else if( nMatch(" H2 ",chCard) )
00364 {
00365
00366
00367 ParseAtomH2(chCard);
00368 }
00369
00370 else
00371 {
00372 fprintf( ioQQQ, " I could not recognize a keyword on this atom command.\n");
00373 fprintf( ioQQQ, " The available keys are FeII, H-Like, He-like, rotor and H2.\n");
00374 fprintf( ioQQQ, " Sorry.\n" );
00375 puts( "[Stop in ParseCommands]" );
00376 cdEXIT(EXIT_FAILURE);
00377 }
00378 }
00379
00380 else if( strcmp(chKey4,"BACK") == 0 )
00381 {
00382
00383 ParseBackgrd(&nqh,chCard,&ar1);
00384 }
00385
00386 else if( strcmp(chKey4,"BLAC") == 0 )
00387 {
00388
00389 ParseBlackbody(chCard,&nqh,&ar1);
00390
00391
00392 if( optimize.lgVarOn )
00393 {
00394
00395 optimize.nvarxt[optimize.nparm] = 1;
00396 strcpy( optimize.chVarFmt[optimize.nparm], "BLACKbody=%f" );
00397
00398 optimize.nvfpnt[optimize.nparm] = input.nRead;
00399
00400
00401 optimize.vparm[0][optimize.nparm] = (float)log10(rfield.slope[rfield.nspec-1]);
00402
00403 optimize.vincr[optimize.nparm] = 0.5;
00404 optimize.nparm += 1;
00405 }
00406 }
00407
00408 else if( strcmp(chKey4,"BREM") == 0 )
00409 {
00410
00411 strcpy( rfield.chSpType[rfield.nspec], "BREMS" );
00412 i = 5;
00413 rfield.slope[rfield.nspec] =
00414 (float)FFmtRead(chCard,&i, INPUT_LINE_LENGTH,&lgEOL);
00415 if( lgEOL )
00416 {
00417 NoNumb(chCard);
00418 }
00419
00420
00421 if( rfield.slope[rfield.nspec] <= 10. )
00422 {
00423 rfield.slope[rfield.nspec] =
00424 pow(10.,rfield.slope[rfield.nspec]);
00425 }
00426 rfield.cutoff[rfield.nspec][0] = 0.;
00427
00428
00429 if( optimize.lgVarOn )
00430 {
00431
00432 optimize.nvarxt[optimize.nparm] = 1;
00433 strcpy( optimize.chVarFmt[optimize.nparm], "BREMS, T=%f" );
00434
00435 optimize.nvfpnt[optimize.nparm] = input.nRead;
00436
00437 optimize.vparm[0][optimize.nparm] = (float)rfield.slope[rfield.nspec];
00438 optimize.vincr[optimize.nparm] = 0.5;
00439 ++optimize.nparm;
00440 }
00441 ++rfield.nspec;
00442 if( rfield.nspec >= LIMSPC )
00443 {
00444
00445 fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
00446 puts( "[Stop in ParseCommands]" );
00447 cdEXIT(EXIT_FAILURE);
00448 }
00449 }
00450
00451 else if( strcmp(chKey4,"CASE") == 0 )
00452
00453 {
00454
00455 opac.lgCaseB = true;
00456
00457
00458 i = 5;
00459 opac.tlamin = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00460 if( lgEOL )
00461 {
00462
00463
00464
00465
00466
00467 opac.tlamin = 1e5;
00468 }
00469 else
00470 {
00471 opac.tlamin = (float)pow(10.f,opac.tlamin);
00472 }
00473
00474
00475 if( nMatch("HUMM",chCard) )
00476 {
00477 opac.lgCaseB_HummerStorey = true;
00478 }
00479
00480
00481 if( nMatch("O PH",chCard) )
00482 {
00483 opac.lgCaseB_no_photo = true;
00484 }
00485
00486
00487 if( nMatch("O PDE",chCard) )
00488 {
00489
00490 opac.lgCaseB_no_pdest = true;
00491 }
00492 }
00493
00494 else if( strcmp(chKey4,"CEXT") == 0 )
00495 {
00496
00497 thermal.lgCExtraOn = true;
00498 i = 5;
00499 thermal.CoolExtra = (float)pow(10.,FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL));
00500 if( lgEOL )
00501 {
00502 NoNumb(chCard);
00503 }
00504 thermal.cextpw = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00505 }
00506
00507
00508 else if( (strcmp(chKey4,"CMB ") == 0) || (strcmp(chKey4,"FIRE") == 0) )
00509 {
00510
00511 i = 5;
00512
00513 z = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00514
00515 ParseCMB(z,&nqh,&ar1);
00516 }
00517
00518 else if( strcmp(chKey4,"COMP") == 0 )
00519 {
00520
00521 ParseCompile(chCard);
00522 }
00523
00524 else if( strcmp(chKey4,"CONS") == 0 )
00525 {
00526
00527
00528 ParseConstant(chCard);
00529 }
00530
00531 else if( strcmp(chKey4,"CORO") == 0 )
00532 {
00533
00534
00535 ParseCoronal( chCard,&nqh,&ar1);
00536 }
00537
00538 else if( strcmp(chKey4,"COSM") == 0 )
00539 {
00540
00541
00542 ParseCosmicRays( chCard );
00543 }
00544
00545 else if( strcmp(chKey4,"COVE") == 0 )
00546 {
00547
00548 i = 5;
00549
00550
00551 geometry.covgeo = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00552
00553 if( lgEOL )
00554 {
00555 NoNumb(chCard);
00556 }
00557
00558
00559 if( geometry.covgeo <= 0. )
00560 {
00561 geometry.covgeo = (float)pow(10.f,geometry.covgeo);
00562 }
00563
00564 if( geometry.covgeo > 1. )
00565 {
00566 fprintf( ioQQQ, " A covering factor greater than 1 makes no physical sense. Sorry.\n" );
00567 puts( "[Stop in ParseCommands]" );
00568 cdEXIT(EXIT_FAILURE);
00569 }
00570
00571
00572 geometry.covrt = geometry.covgeo;
00573 }
00574
00575 else if( strcmp(chKey4,"CRAS") == 0 )
00576 {
00577
00578 ParseCrashDo(chCard);
00579 }
00580
00581 else if( strcmp(chKey4,"CYLI") == 0 )
00582 {
00583
00584 i = 5;
00585 radius.lgCylnOn = true;
00586 radius.CylindHigh = pow(10.,FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL));
00587 if( lgEOL )
00588 {
00589 NoNumb(chCard);
00590 }
00591 }
00592
00593 else if( strcmp(chKey4,"DIEL") == 0 )
00594 {
00595
00596 fprintf( ioQQQ, " The DIELectronic command has been replaced with the SET DIELectronic recombination command.\n" );
00597 fprintf( ioQQQ, " Please have a look at Hazy.\n Sorry.\n\n" );
00598 puts( "[Stop in ParseCommands]" );
00599 cdEXIT(EXIT_FAILURE);
00600 }
00601
00602 else if( strcmp(chKey4,"DIFF") == 0 )
00603 {
00604
00605
00606 if( nMatch(" OTS",chCard) )
00607 {
00608 if( nMatch("SIMP",chCard) )
00609 {
00610
00611 strcpy( rfield.chDffTrns, "OSS" );
00612 }
00613 else
00614 {
00615
00616 strcpy( rfield.chDffTrns, "OTS" );
00617 }
00618 rfield.lgOutOnly = false;
00619 }
00620 else if( nMatch(" OUT",chCard) )
00621 {
00622 rfield.lgOutOnly = true;
00623 i = 4;
00624 j = (long int)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00625 if( lgEOL )
00626 {
00627
00628 strcpy( rfield.chDffTrns, "OU2" );
00629 }
00630 else
00631 {
00632 if( j > 0 && j < 10 )
00633 {
00634 sprintf( rfield.chDffTrns, "OU%1ld", j );
00635 }
00636 else
00637 {
00638 fprintf( ioQQQ, " must be between 1 and 9 \n" );
00639 puts( "[Stop in ParseCommands]" );
00640 cdEXIT(EXIT_FAILURE);
00641 }
00642 }
00643 }
00644
00645 else
00646 {
00647 fprintf( ioQQQ, " There should have been OUTward or OTS on this line. Sorry.\n" );
00648 puts( "[Stop in ParseCommands]" );
00649 cdEXIT(EXIT_FAILURE);
00650 }
00651 }
00652
00653 else if( strcmp(chKey4,"DIST") == 0 )
00654 {
00655
00656 i = 5;
00657 radius.distance = FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00658 if( lgEOL )
00659 {
00660 NoNumb(chCard);
00661 }
00662
00663
00664 if( !nMatch("LINE",chCard ) )
00665 {
00666 radius.distance = pow(10., radius.distance );
00667 }
00668
00669
00670 if( nMatch("PARS",chCard ) )
00671 {
00672 radius.distance *= PARSEC;
00673 }
00674 }
00675
00676 else if( strcmp(chKey4,"DLAW") == 0 )
00677 {
00678
00679 ParseDLaw(chCard);
00680 }
00681
00682 else if( strcmp(chKey4,"DOUB") == 0 )
00683 {
00684
00685 rt.DoubleTau = 2.;
00686 }
00687
00688 else if( strcmp(chKey4,"DRIV") == 0 )
00689 {
00690
00691 ParseDrive(chCard);
00692 }
00693
00694 else if( strcmp(chKey4,"EDEN") == 0 )
00695 {
00696
00697
00698 i = 5;
00699 dense.EdenExtra = (float)pow(10.,FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL));
00700 if( lgEOL )
00701 {
00702 NoNumb(chCard);
00703 }
00704
00705 phycon.lgPhysOK = false;
00706 }
00707
00708 else if( strcmp(chKey4,"ELEM") == 0 )
00709 {
00710
00711
00712
00713
00714 ParseElement(chCard);
00715 }
00716
00717 else if( strcmp(chKey4,"ENER") == 0 )
00718 {
00719
00720 i = 5;
00721 if( nqh >= LIMSPC )
00722 {
00723
00724 fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
00725 puts( "[Stop in ParseCommands]" );
00726 cdEXIT(EXIT_FAILURE);
00727 }
00728
00729 ASSERT( nqh < LIMSPC );
00730 strcpy( rfield.chRSpec[nqh], "SQCM" );
00731 teset = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00732 if( lgEOL )
00733 {
00734 NoNumb(chCard);
00735 }
00736
00737
00738 if( !radius.lgRadiusKnown )
00739 {
00740
00741 ar1 = (float)radius.rdfalt;
00742 radius.Radius = pow(10.,radius.rdfalt);
00743 }
00744
00745
00746 if( nMatch("LINE",chCard) || teset > 10. )
00747 {
00748
00749 teset = (float)log10(teset);
00750 }
00751
00752 if( teset > 5. )
00753 {
00754 fprintf( ioQQQ, " This intensity may be too large. The code may crash due to overflow. Was log intended?\n" );
00755 }
00756
00757
00758 strcpy( rfield.chSpNorm[nqh], "LUMI" );
00759
00760
00761 rfield.range[nqh][0] = rfield.emm;
00762 rfield.range[nqh][1] = rfield.egamry;
00763 rfield.totpow[nqh] = (4.*teset - 4.2464476 + 0.60206);
00764
00765
00766 if( nMatch( "TIME" , chCard ) )
00767 rfield.lgTimeVary[nqh] = true;
00768
00769
00770 if( optimize.lgVarOn )
00771 {
00772 strcpy( optimize.chVarFmt[optimize.nparm], "ENERGY DENSITY %f log " );
00773
00774 optimize.nvfpnt[optimize.nparm] = input.nRead;
00775 optimize.vparm[0][optimize.nparm] = (float)log10(rfield.totpow[nqh]);
00776 optimize.vincr[optimize.nparm] = 0.1f;
00777 optimize.nvarxt[optimize.nparm] = 1;
00778 optimize.nparm += 1;
00779 }
00780
00781
00782 ++nqh;
00783 }
00784
00785 else if( strcmp(chKey4,"EXTI") == 0 )
00786 {
00787
00788
00789
00790
00791
00792
00793
00794
00795 i = 5;
00796 extinc.excolm = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00797 if( lgEOL )
00798 {
00799 NoNumb(chCard);
00800 }
00801
00802
00803
00804
00805 if( !nMatch("LINE" , chCard ) )
00806 {
00807
00808 if( extinc.excolm>35. )
00809 {
00810 fprintf(ioQQQ,
00811 " The first parameter on this command line is the log of either the column density or optical depth.\n");
00812 fprintf(ioQQQ,
00813 " The value seems pretty big to me - please check it.\n");
00814
00815 fflush(ioQQQ);
00816 }
00817 extinc.excolm = (float)pow(10.f,extinc.excolm);
00818 }
00819
00820
00821 extinc.exleak = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00822
00823 if( lgEOL )
00824 {
00825 extinc.exleak = 0.;
00826 }
00827
00828 if( extinc.exleak < 0. )
00829 {
00830 extinc.exleak = (float)pow(10.f,extinc.exleak);
00831 }
00832 if( extinc.exleak > 1. )
00833 {
00834
00835 if( called.lgTalk )
00836 {
00837 fprintf( ioQQQ, " A leakage of%9.0f%% was entered - this must be less than 100%%\n",
00838 extinc.exleak*100. );
00839 }
00840 puts( "[Stop in ParseCommands]" );
00841 cdEXIT(EXIT_FAILURE);
00842 }
00843
00844
00845 extinc.exlow = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00846 if( lgEOL )
00847 {
00848 extinc.exlow = 0.99946f;
00849 }
00850
00851 else
00852 {
00853 if( extinc.exlow <= 0. )
00854 extinc.exlow = (float)pow(10.f,extinc.exlow);
00855 if( extinc.exlow < 0.99946 )
00856 {
00857 fprintf( ioQQQ, " Energy less than 1 Ryd!!\n" );
00858 }
00859 }
00860
00861
00862 if( nMatch("OPTI" , chCard ) )
00863 {
00864
00865 extinc.excolm /= (float)(extinc.cnst_col2optdepth*
00866 pow(extinc.exlow,extinc.cnst_power) );
00867 }
00868 }
00869
00870 else if( strcmp(chKey4,"FAIL") == 0 )
00871 {
00872
00873 i = 5;
00874
00875
00876 j = conv.LimFail;
00877 conv.LimFail = (long int)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00878 if( lgEOL )
00879 {
00880 NoNumb(chCard);
00881 }
00882
00883
00884
00885
00886 if( nMatch(" MAP",chCard) && !nMatch(" NO ",chCard) )
00887 {
00888 conv.lgMap = true;
00889 }
00890
00891
00892 if( conv.LimFail > j )
00893 {
00894 fprintf( ioQQQ, " This command should not be necessary.\n" );
00895 fprintf( ioQQQ, " Please show this input stream to Gary Ferland if this command is really needed for this simulation.\n" );
00896 }
00897 }
00898
00899 else if( strcmp(chKey4,"FILL") == 0 )
00900 {
00901
00902 i = 5;
00903 a = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00904 if( lgEOL )
00905 {
00906 NoNumb(chCard);
00907 }
00908
00909 if( a <= 0. )
00910 {
00911
00912 geometry.FillFac = (float)pow(10.f,a);
00913 }
00914 else
00915 {
00916
00917 geometry.FillFac = a;
00918 if( geometry.FillFac > 1.0 )
00919 {
00920 if( called.lgTalk )
00921 {
00922 fprintf( ioQQQ, " Filling factor > 1, reset to 1\n" );
00923 }
00924 geometry.FillFac = 1.;
00925 }
00926 }
00927 geometry.fiscal = geometry.FillFac;
00928
00929
00930 geometry.filpow = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00931
00932
00933 if( optimize.lgVarOn )
00934 {
00935 strcpy( optimize.chVarFmt[optimize.nparm], "FILLING FACTOR= %f power=%f" );
00936
00937
00938 optimize.nvfpnt[optimize.nparm] = input.nRead;
00939 optimize.vparm[0][optimize.nparm] = (float)log10(geometry.FillFac);
00940 optimize.vincr[optimize.nparm] = 0.5;
00941
00942
00943 optimize.vparm[1][optimize.nparm] = geometry.filpow;
00944 optimize.nvarxt[optimize.nparm] = 2;
00945
00946
00947 optimize.varang[optimize.nparm][0] = -1e38f;
00948 optimize.varang[optimize.nparm][1] = 0.f;
00949 optimize.nparm += 1;
00950 }
00951 }
00952
00953 else if( strcmp(chKey4,"FLUC") == 0 )
00954 {
00955
00956 ParseFluc(chCard);
00957 }
00958
00959 else if( strcmp(chKey4,"F(NU") == 0 )
00960 {
00961
00962
00963 lgNu2 = false;
00964
00965 ParseF_nu(chCard,&nqh,&ar1,"SQCM",lgNu2);
00966 }
00967
00968 else if( strcmp(chKey4,"FORC") == 0 )
00969 {
00970
00971
00972
00973 i = 5;
00974 thermal.ConstTemp = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00975 if( lgEOL )
00976 {
00977 NoNumb(chCard);
00978 }
00979
00980 if( nMatch(" LOG",chCard) || (thermal.ConstTemp <= 10. &&
00981 !nMatch("LINE",chCard)) )
00982 {
00983 thermal.ConstTemp = (float)pow(10.f,thermal.ConstTemp);
00984 }
00985
00986
00987 if( thermal.ConstTemp < 3. )
00988 {
00989 fprintf( ioQQQ, " TE reset to 3K: entered number too small.\n" );
00990 thermal.ConstTemp = 3.;
00991 }
00992 }
00993
00994 else if( strcmp(chKey4,"FUDG") == 0 )
00995 {
00996
00997 i = 5;
00998 for( j=0; j < NFUDGC; j++ )
00999 {
01000 fudgec.fudgea[j] = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01001 if( !lgEOL )
01002 fudgec.nfudge = j+1;
01003 }
01004 }
01005
01006 else if( strcmp(chKey4,"GLOB") == 0 )
01007 {
01008
01009
01010 ParseGlobule(chCard);
01011 }
01012
01013 else if( (strcmp(chKey4,"GRAI") == 0 )||( strcmp(chKey4,"PGRA") == 0 ) )
01014 {
01015
01016 ParseGrain(chCard,&lgDSet);
01017 }
01018
01019 else if( strcmp(chKey4,"GRID") == 0 )
01020 {
01021
01022
01023 ParseOptimize(chCard);
01024 }
01025
01026 else if( strcmp(chKey4,"HDEN") == 0 )
01027 {
01028
01029 ParseHDEN(chCard);
01030 }
01031
01032 else if( strcmp(chKey4,"HELI") == 0 )
01033 {
01034
01035 ParseAtomHeLike(chCard);
01036 }
01037
01038 else if( strcmp(chKey4,"HEXT") == 0 )
01039 {
01040
01041
01042
01043
01044 i = 5;
01045 hextra.TurbHeat = (float)pow(10.,FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL));
01046 a = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01047 if( lgEOL )
01048 {
01049 hextra.turrad = 0.;
01050 }
01051 else
01052 {
01053 hextra.turrad = (float)pow(10.f,a);
01054 }
01055 a = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01056 if( lgEOL )
01057 {
01058 hextra.turback = 0.;
01059 }
01060 else
01061 {
01062 hextra.turback = (float)pow(10.f,a);
01063 }
01064 }
01065
01066 else if( strcmp(chKey4,"HIGH") == 0 )
01067 {
01068
01069 thermal.lgTeHigh = true;
01070 }
01071
01072 else if( strncmp( chCard ,"HYDROGEN",8) == 0 )
01073 {
01074 fprintf(ioQQQ," This command has been replaced with the ATOM H-LIKE command.\n");
01075 fprintf(ioQQQ," I will parse the command for now, but may not in the future.\n");
01076
01077 ParseAtomHLike(chCard);
01078 }
01079
01080 else if( strcmp(chKey4,"ILLU") == 0 )
01081 {
01082
01083 i = 5;
01084 geometry.AngleIllumRadian = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01085 if( lgEOL )
01086 {
01087 NoNumb(chCard);
01088 }
01089
01090 if( nMatch("RADI",chCard) )
01091 {
01092
01093 geometry.AngleIllumRadian /= (float)(PI/180.);
01094 }
01095
01096 if( geometry.AngleIllumRadian < 0. || geometry.AngleIllumRadian >= 90. )
01097 {
01098 fprintf( ioQQQ, " Angle of illumination must be between 0 and 90.\n" );
01099 puts( "[Stop in ParseCommands]" );
01100 cdEXIT(EXIT_FAILURE);
01101 }
01102
01103
01104 geometry.AngleIllumRadian = (float)(geometry.AngleIllumRadian*PI/180.);
01105
01106
01107 geometry.DirectionalCosin = (float)(1./cos(geometry.AngleIllumRadian));
01108
01109
01110 if( optimize.lgVarOn )
01111 {
01112
01113 optimize.nvarxt[optimize.nparm] = 1;
01114 strcpy( optimize.chVarFmt[optimize.nparm], "ILLUminate %f radians " );
01115
01116 optimize.nvfpnt[optimize.nparm] = input.nRead;
01117 optimize.vparm[0][optimize.nparm] = geometry.AngleIllumRadian;
01118
01119 optimize.vincr[optimize.nparm] = 0.1f;
01120 optimize.nparm += 1;
01121 }
01122 }
01123
01124 else if( strcmp(chKey4,"INIT") == 0 )
01125 {
01126
01127
01128
01129
01130 ParseInit(chCard);
01131
01132
01133
01134 ++nInitFile;
01135 if( nInitFile > 1 )
01136 {
01137 fprintf( ioQQQ,
01138 " This is the second init file, I can only handle one.\nSorry.\n" );
01139 puts( "[Stop in ParseCommands]" );
01140 cdEXIT(EXIT_FAILURE);
01141 }
01142
01143
01144
01145 input.iReadWay = -1;
01146
01147
01148
01149 input_init();
01150 }
01151
01152 else if( strcmp(chKey5,"INTEN") == 0 )
01153 {
01154
01155 i = 5;
01156 if( nqh >= LIMSPC )
01157 {
01158
01159 fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
01160 puts( "[Stop in ParseCommands]" );
01161 cdEXIT(EXIT_FAILURE);
01162 }
01163
01164
01165 ASSERT( nqh < LIMSPC );
01166 strcpy( rfield.chRSpec[nqh], "SQCM" );
01167 rfield.totpow[nqh] = FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01168 if( lgEOL )
01169 {
01170 NoNumb(chCard);
01171 }
01172
01173
01174
01175 if( !radius.lgRadiusKnown )
01176 {
01177
01178 ar1 = (float)radius.rdfalt;
01179 radius.Radius = pow(10.,radius.rdfalt);
01180 }
01181 if( nMatch("LINE",chCard) )
01182 {
01183
01184 ASSERT( nqh < LIMSPC );
01185
01186 rfield.totpow[nqh] = log10(rfield.totpow[nqh]);
01187 }
01188 strcpy( rfield.chSpNorm[nqh], "LUMI" );
01189
01190 ParseRangeOption(nqh,chCard);
01191
01192
01193 if( nMatch( "TIME" , chCard ) )
01194 rfield.lgTimeVary[nqh] = true;
01195
01196
01197 if( optimize.lgVarOn )
01198 {
01199
01200
01201
01202
01203
01204 strcpy( optimize.chVarFmt[optimize.nparm], "INTENSITY %f range %f %f " );
01205
01206 optimize.nvfpnt[optimize.nparm] = input.nRead;
01207 optimize.vparm[0][optimize.nparm] = (float)rfield.totpow[nqh];
01208 optimize.vincr[optimize.nparm] = 0.5;
01209
01210
01211
01212
01213 optimize.vparm[1][optimize.nparm] = (float)rfield.range[nqh][0];
01214 optimize.vparm[2][optimize.nparm] = (float)rfield.range[nqh][1];
01215 optimize.nvarxt[optimize.nparm] = 3;
01216 ++optimize.nparm;
01217 }
01218 ++nqh;
01219 }
01220
01221 else if( strcmp(chKey5,"INTER") == 0 )
01222 {
01223
01224
01225
01226
01227
01228 ParseInterp(chCard,&lgEOF);
01229 }
01230
01231 else if( strcmp(chKey4,"IONI") == 0 )
01232 {
01233
01234
01235
01236 ParseIonPar(&nqh,chCard,&ar1);
01237 }
01238
01239 else if( strcmp(chKey4,"ITER") == 0 )
01240 {
01241
01242 i = 5;
01243 iterations.itermx = (long int)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL) - 1;
01244 iterations.itermx = MAX2(iterations.itermx,1);
01245
01246
01247
01248 if( iterations.itermx > iterations.iter_malloc - 1 )
01249 {
01250 long int iter_malloc_save = iterations.iter_malloc;
01251
01252
01253 iterations.iter_malloc = iterations.itermx+3;
01254 iterations.IterPrnt = (long int*)REALLOC(iterations.IterPrnt ,
01255 (size_t)iterations.iter_malloc*sizeof(long int) );
01256 geometry.nend = (long int*)REALLOC(geometry.nend ,
01257 (size_t)iterations.iter_malloc*sizeof(long int) );
01258 radius.router = (double*)REALLOC(radius.router ,
01259 (size_t)iterations.iter_malloc*sizeof(double) );
01260
01261 for(j=iter_malloc_save; j<iterations.iter_malloc; ++j )
01262 {
01263 iterations.IterPrnt[j] = iterations.IterPrnt[iter_malloc_save-1];
01264 geometry.nend[j] = geometry.nend[iter_malloc_save-1];
01265 radius.router[j] = radius.router[iter_malloc_save-1];
01266 }
01267 }
01268
01269 if( nMatch("CONV",chCard) )
01270 {
01271
01272
01273 conv.lgAutoIt = true;
01274
01275 if( lgEOL )
01276 {
01277 iterations.itermx = 10 - 1;
01278 }
01279 a = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01280
01281 if( !lgEOL )
01282 {
01283 conv.autocv = a;
01284 }
01285 }
01286 }
01287
01288 else if( strcmp(chKey4,"L(NU") == 0 )
01289 {
01290
01291
01292 lgNu2 = false;
01293
01294 ParseF_nu(chCard,&nqh,&ar1,"4 PI",lgNu2);
01295 }
01296
01297 else if( strcmp(chKey4,"LASE") == 0 )
01298 {
01299
01300
01301 strcpy( rfield.chSpType[rfield.nspec], "LASER" );
01302
01303
01304 i = 5;
01305 rfield.slope[rfield.nspec] = FFmtRead(chCard,&i, INPUT_LINE_LENGTH,&lgEOL);
01306
01307
01308 if( rfield.slope[rfield.nspec] <= 0. )
01309 {
01310 rfield.slope[rfield.nspec] =
01311 pow(10.,rfield.slope[rfield.nspec]);
01312 }
01313 if( lgEOL )
01314 {
01315 NoNumb(chCard);
01316 }
01317
01318
01319 rfield.cutoff[rfield.nspec][0] = FFmtRead(chCard,&i, INPUT_LINE_LENGTH,&lgEOL);
01320 if( lgEOL )
01321 {
01322
01323 rfield.cutoff[rfield.nspec][0] = 0.05;
01324 }
01325
01326
01327 ++rfield.nspec;
01328 if( rfield.nspec >= LIMSPC )
01329 {
01330
01331 fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
01332 puts( "[Stop in ParseCommands]" );
01333 cdEXIT(EXIT_FAILURE);
01334 }
01335 }
01336
01337 else if( strcmp(chKey4,"LUMI") == 0 )
01338 {
01339
01340 i = 5;
01341 if( nqh >= LIMSPC )
01342 {
01343
01344 fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
01345 puts( "[Stop in ParseCommands]" );
01346 cdEXIT(EXIT_FAILURE);
01347 }
01348
01349 strcpy( rfield.chRSpec[nqh], "4 PI" );
01350 rfield.totpow[nqh] = FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01351 if( lgEOL )
01352 {
01353 NoNumb(chCard);
01354 }
01355 if( nMatch("LINE",chCard) )
01356 {
01357
01358 rfield.totpow[nqh] = log10(rfield.totpow[nqh]);
01359 }
01360
01361 strcpy( rfield.chSpNorm[nqh], "LUMI" );
01362
01363
01364 if( nMatch("SOLA",chCard) )
01365 {
01366
01367 rfield.range[nqh][0] = rfield.emm;
01368 rfield.range[nqh][1] = rfield.egamry;
01369
01370 rfield.totpow[nqh] += 33.5827f;
01371 }
01372 else
01373 {
01374
01375
01376 ParseRangeOption(nqh,chCard);
01377 }
01378
01379
01380 if( nMatch( "TIME" , chCard ) )
01381 rfield.lgTimeVary[nqh] = true;
01382
01383
01384 if( optimize.lgVarOn )
01385 {
01386
01387
01388
01389 strcpy( optimize.chVarFmt[optimize.nparm], "LUMINOSITY %f range %f %f " );
01390
01391 optimize.nvfpnt[optimize.nparm] = input.nRead;
01392 optimize.vparm[0][optimize.nparm] = (float)rfield.totpow[nqh];
01393 optimize.vincr[optimize.nparm] = 0.5;
01394
01395 optimize.vparm[1][optimize.nparm] = (float)rfield.range[nqh][0];
01396 optimize.vparm[2][optimize.nparm] = (float)rfield.range[nqh][1];
01397 optimize.nvarxt[optimize.nparm] = 3;
01398 optimize.nparm += 1;
01399 }
01400 ++nqh;
01401 }
01402
01403 else if( strcmp(chKey4,"MAGN") == 0 )
01404 {
01405
01406 ParseMagnet( chCard );
01407 }
01408
01409 else if( strcmp(chKey4,"MAP ") == 0 )
01410 {
01411
01412 ParseMap(chCard);
01413 }
01414
01415 else if( strcmp(chKey4,"META") == 0 )
01416 {
01417
01418
01419 ParseMetal(chCard);
01420
01421 abund.lgAbnSolar = false;
01422 }
01423
01424 else if( strcmp(chKey4,"NEUT") == 0 )
01425 {
01426
01427
01428 i = 5;
01429 neutrn.frcneu = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01430 if( lgEOL )
01431 {
01432 NoNumb(chCard);
01433 }
01434 if( neutrn.frcneu > 0. )
01435 {
01436 neutrn.frcneu = (float)log10(neutrn.frcneu);
01437 }
01438
01439 neutrn.lgNeutrnHeatOn = true;
01440 neutrn.effneu = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01441 if( lgEOL )
01442 {
01443 neutrn.effneu = 1.0;
01444 }
01445 else
01446 {
01447 if( neutrn.effneu <= 0. )
01448 neutrn.effneu = (float)pow(10.f,neutrn.effneu);
01449 }
01450 }
01451
01452 else if( strcmp(chKey3,"NO ") == 0 )
01453 {
01454
01455 ParseDont(chCard);
01456 }
01457
01458 else if( strcmp(chKey4,"NORM") == 0 )
01459 {
01460
01461 ParseNorm(chCard);
01462 }
01463
01464 else if( strcmp(chKey4,"NUF(") == 0 )
01465 {
01466
01467
01468 lgNu2 = true;
01469
01470 ParseF_nu(chCard,&nqh,&ar1,"SQCM",lgNu2);
01471 }
01472
01473 else if( strcmp(chKey4,"NUL(") == 0 )
01474 {
01475
01476
01477 lgNu2 = true;
01478
01479 ParseF_nu(chCard,&nqh,&ar1,"4 PI",lgNu2);
01480 }
01481
01482 else if( strcmp(chKey4,"OPTI") == 0 )
01483 {
01484
01485
01486 ParseOptimize(chCard);
01487 }
01488
01489 else if( strcmp(chKey4,"PHI(") == 0 )
01490 {
01491
01492 i = 5;
01493 if( nqh >= LIMSPC )
01494 {
01495
01496 fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
01497 puts( "[Stop in ParseCommands]" );
01498 cdEXIT(EXIT_FAILURE);
01499 }
01500
01501 ASSERT( nqh < LIMSPC );
01502 strcpy( rfield.chRSpec[nqh], "SQCM" );
01503 strcpy( rfield.chSpNorm[nqh], "PHI " );
01504 rfield.totpow[nqh] = FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01505 if( lgEOL )
01506 {
01507 NoNumb(chCard);
01508 }
01509
01510
01511 if( !radius.lgRadiusKnown )
01512 {
01513
01514 ar1 = (float)radius.rdfalt;
01515 radius.Radius = pow(10.,radius.rdfalt);
01516 }
01517
01518 if( rfield.totpow[nqh] > 29. )
01519 {
01520 fprintf( ioQQQ, " Is the flux for this continuum correct?\n" );
01521 fprintf( ioQQQ, " It appears too bright to me.\n" );
01522 }
01523
01524 ParseRangeOption(nqh,chCard);
01525
01526
01527 if( nMatch( "TIME" , chCard ) )
01528 rfield.lgTimeVary[nqh] = true;
01529
01530
01531 if( optimize.lgVarOn )
01532 {
01533
01534
01535
01536 strcpy( optimize.chVarFmt[optimize.nparm], "phi(h) %f range %f %f" );
01537
01538 optimize.nvfpnt[optimize.nparm] = input.nRead;
01539 optimize.vparm[0][optimize.nparm] = (float)rfield.totpow[nqh];
01540 optimize.vincr[optimize.nparm] = 0.5;
01541
01542 optimize.vparm[1][optimize.nparm] = (float)rfield.range[nqh][0];
01543 optimize.vparm[2][optimize.nparm] = (float)rfield.range[nqh][1];
01544 optimize.nvarxt[optimize.nparm] = 3;
01545 optimize.nparm += 1;
01546 }
01547 ++nqh;
01548 }
01549
01550 else if( strcmp(chKey4,"PLOT") == 0 )
01551 {
01552
01553
01554 ParsePlot(chCard);
01555 }
01556
01557 else if( strcmp(chKey4,"POWE") == 0 )
01558 {
01559
01560 ParsePowerlawContinuum(chCard);
01561 }
01562
01563 else if( strcmp(chKey4,"PRIN") == 0 )
01564 {
01565
01566 ParsePrint(chCard);
01567 }
01568
01569 else if( strcmp(chKey4,"PUNC") == 0 )
01570 {
01571
01572 ParsePunch(chCard);
01573 }
01574
01575 else if( strcmp(chKey4,"Q(H)") == 0 )
01576 {
01577
01578 i = 5;
01579 if( nqh >= LIMSPC )
01580 {
01581
01582 fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
01583 puts( "[Stop in ParseCommands]" );
01584 cdEXIT(EXIT_FAILURE);
01585 }
01586
01587
01588 ASSERT( nqh < LIMSPC );
01589 strcpy( rfield.chRSpec[nqh], "4 PI" );
01590 strcpy( rfield.chSpNorm[nqh], "Q(H)" );
01591 rfield.totpow[nqh] = FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01592 if( rfield.totpow[nqh] > 100. && called.lgTalk )
01593 {
01594 fprintf( ioQQQ, " Is this reasonable?\n" );
01595 }
01596 if( lgEOL )
01597 {
01598 NoNumb(chCard);
01599 }
01600
01601 ParseRangeOption(nqh,chCard);
01602
01603
01604 if( nMatch( "TIME" , chCard ) )
01605 rfield.lgTimeVary[nqh] = true;
01606
01607
01608 if( optimize.lgVarOn )
01609 {
01610
01611
01612
01613 strcpy( optimize.chVarFmt[optimize.nparm], "Q(H) %f range %f %f" );
01614
01615 optimize.nvfpnt[optimize.nparm] = input.nRead;
01616 optimize.vparm[0][optimize.nparm] = (float)rfield.totpow[nqh];
01617 optimize.vincr[optimize.nparm] = 0.5;
01618
01619 optimize.vparm[1][optimize.nparm] = (float)rfield.range[nqh][0];
01620 optimize.vparm[2][optimize.nparm] = (float)rfield.range[nqh][1];
01621 optimize.nvarxt[optimize.nparm] = 3;
01622 ++optimize.nparm;
01623 }
01624
01625 ++nqh;
01626 }
01627
01628 else if( strcmp(chKey4,"RATI") == 0 )
01629 {
01630
01631
01632
01633
01634
01635
01636
01637 ParseRatio(chCard,&nqh);
01638 }
01639
01640 else if( strcmp(chKey4,"RADI") == 0 )
01641 {
01642
01643
01644
01645
01646 ParseRadius(chCard,&ar1);
01647 }
01648
01649 else if( strcmp(chKey4,"ROBE") == 0 )
01650 {
01651
01652 radius.dRadSign = -1.;
01653 }
01654
01655 else if( strcmp(chKey4,"SET ") == 0 )
01656 {
01657
01658 ParseSet(chCard);
01659 }
01660
01661 else if( strcmp(chKey4,"SPEC") == 0 )
01662 {
01663
01664 cdEXIT(EXIT_FAILURE);
01665 }
01666
01667 else if( strcmp(chKey4,"SPHE") == 0 )
01668 {
01669
01670
01671 ParseSphere(chCard);
01672 }
01673
01674 else if( strcmp(chKey4,"STAT") == 0 )
01675 {
01676
01677 ParseState(chCard);
01678 }
01679
01680 else if( strcmp(chKey4,"STOP") == 0 )
01681 {
01682
01683
01684 ParseStop(chCard);
01685 }
01686
01687 else if( strcmp(chKey4,"TABL") == 0 )
01688 {
01689
01690
01691
01692
01693 ParseTable(&nqh,chCard,&ar1);
01694 }
01695
01696 else if( strcmp(chKey4,"TAUM") == 0 )
01697 {
01698
01699
01700 i = 5;
01701 opac.taumin = (float)pow(10.,FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL));
01702 if( lgEOL )
01703 {
01704 NoNumb(chCard);
01705 }
01706 }
01707
01708 else if( strcmp(chKey4 , "TEST" ) == 0 )
01709 {
01710 char chStuff[INPUT_LINE_LENGTH];
01711
01712
01713
01714
01715 int nPrintTest = nMatch("PRIN" , chCard );
01716
01717
01718 sprintf( chStuff , "HDEN 4 " );
01719 if( nPrintTest )
01720 fprintf(ioQQQ , "%s\n" , chStuff );
01721 ParseHDEN(chStuff);
01722
01723
01724 sprintf( chStuff , "CONSTANT TEMPER 4 " );
01725 if( nPrintTest )
01726 fprintf(ioQQQ , "%s\n" , chStuff );
01727 ParseConstant(chStuff);
01728
01729
01730 sprintf( chStuff , "TABLE AGN " );
01731 if( nPrintTest )
01732 fprintf(ioQQQ , "%s\n" , chStuff );
01733 ParseTable(&nqh,chStuff,&ar1);
01734
01735
01736 sprintf( chStuff , "IONIZATION PARAMETER -2 " );
01737 if( nPrintTest )
01738 fprintf(ioQQQ , "%s\n" , chStuff );
01739 ParseIonPar(&nqh,chStuff,&ar1);
01740
01741
01742 sprintf( chStuff , "ABUNDANCES OLD SOLAR 84 " );
01743 if( nPrintTest )
01744 fprintf(ioQQQ , "%s\n" , chStuff );
01745 ParseAbundances(chStuff,lgDSet);
01746
01747
01748
01749
01750 sprintf( chStuff , "STOP LYMAN OPTICAL -4 " );
01751 if( nPrintTest )
01752 fprintf(ioQQQ , "%s\n" , chStuff );
01753 ParseStop(chStuff);
01754
01755
01756 if( nMatch("LARG",chCard) )
01757 {
01758 sprintf( chStuff , "ATOM H-LIKE ELEMENT HYDROGEN LEVELS LIMIT " );
01759 if( nPrintTest )
01760 fprintf(ioQQQ , "%s\n" , chStuff );
01761 ParseAtomHLike(chStuff);
01762 }
01763
01764
01765 if( nMatch(" H2 ",chCard) )
01766 {
01767
01768
01769 sprintf( chStuff , "ATOM H2 LIMIT -20 " );
01770 if( nPrintTest )
01771 fprintf(ioQQQ , "%s\n" , chStuff );
01772 ParseAtomH2(chStuff);
01773 }
01774
01775
01776 sprintf( chStuff , "STOP ZONE 2 " );
01777 if( nPrintTest )
01778 fprintf(ioQQQ , "%s\n" , chStuff );
01779 ParseStop(chStuff);
01780
01781
01782 sprintf( chStuff , "SET DR 0 " );
01783 if( nPrintTest )
01784 fprintf(ioQQQ , "%s\n" , chStuff );
01785 ParseSet(chStuff);
01786
01787
01788 if( nMatch("FEII",chCard) || nMatch("FE II",chCard) )
01789 {
01790 sprintf( chStuff , "ATOM FEII VERNER " );
01791 if( nPrintTest )
01792 fprintf(ioQQQ , "%s\n" , chStuff );
01793 ParseAtomFeII(chStuff);
01794 }
01795
01796
01797 sprintf( input.chCARDCAPS , "ASSERT HYDROGEN 1 IONIZATION -3.040 " );
01798 if( nPrintTest )
01799 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01800 ParseAssertResults();
01801
01802 sprintf( input.chCARDCAPS , "ASSERT HELIUM 2 IONIZATION -1.067 " );
01803 if( nPrintTest )
01804 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01805 ParseAssertResults();
01806
01807 sprintf( input.chCARDCAPS , "ASSERT CARBON 2 IONIZATION -2.301 " );
01808 if( nPrintTest )
01809 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01810 ParseAssertResults();
01811
01812
01813 sprintf( input.chCARDCAPS , "ASSERT CARBON 3 IONIZATION -0.560 " );
01814 if( nPrintTest )
01815 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01816 ParseAssertResults();
01817
01818
01819 sprintf( input.chCARDCAPS , "ASSERT CARBON 4 IONIZATION -0.373 " );
01820 if( nPrintTest )
01821 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01822 ParseAssertResults();
01823
01824
01825 sprintf( input.chCARDCAPS , "ASSERT CARBON 5 IONIZATION -0.530 " );
01826 if( nPrintTest )
01827 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01828 ParseAssertResults();
01829
01830
01831 sprintf( input.chCARDCAPS , "ASSERT OXYGEN 3 IONIZATION -0.861 " );
01832 if( nPrintTest )
01833 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01834 ParseAssertResults();
01835
01836
01837 sprintf( input.chCARDCAPS , "ASSERT OXYGEN 4 IONIZATION -0.157 " );
01838 if( nPrintTest )
01839 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01840 ParseAssertResults();
01841
01842
01843 sprintf( input.chCARDCAPS , "ASSERT OXYGEN 5 IONIZATION -0.808 " );
01844 if( nPrintTest )
01845 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01846 ParseAssertResults();
01847
01848
01849 sprintf( input.chCARDCAPS , "ASSERT LINE \"CA B\" 4861 0.946 " );
01850 if( nPrintTest )
01851 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01852
01853 strcpy( input.chOrgCard , input.chCARDCAPS);
01854 ParseAssertResults();
01855
01856
01857
01858
01859 sprintf( input.chCARDCAPS , "ASSERT LINE \"O 3\" 5007 2.72 " );
01860 if( nPrintTest )
01861 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01862
01863 strcpy( input.chOrgCard , input.chCARDCAPS);
01864 ParseAssertResults();
01865
01866 sprintf( input.chCARDCAPS , "ASSERT HTOT -15.011" );
01867 if( nPrintTest )
01868 fprintf(ioQQQ , "%s\n" , input.chCARDCAPS );
01869 ParseAssertResults();
01870
01871 }
01872
01873 else if( strcmp(chKey4,"TIME") == 0 )
01874 {
01875
01876 ParseDynaTime(chCard);
01877 }
01878
01879 else if( strcmp(chKey4,"TITL") == 0 )
01880 {
01881
01882 strcpy( input.chTitle , &input.chOrgCard[5] );
01883 }
01884
01885 else if( strcmp(chKey4,"TLAW") == 0 )
01886 {
01887
01888 ParseTLaw(chCard);
01889 }
01890
01891 else if( strcmp(chKey4,"TOLE") == 0 )
01892 {
01893 fprintf(ioQQQ,
01894 "This command has been replaced with the SET TEMPERATURE TOLERANCE command.\n");
01895 fprintf(ioQQQ,
01896 "I will parse the command for now, but may not in a future version.\n");
01897
01898 i = 5;
01899 conv.HeatCoolRelErrorAllowed = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01900 if( lgEOL )
01901 {
01902 NoNumb(chCard);
01903 }
01904 if( conv.HeatCoolRelErrorAllowed <= 0. )
01905 {
01906 conv.HeatCoolRelErrorAllowed = (float)pow(10.f,conv.HeatCoolRelErrorAllowed);
01907 }
01908 }
01909
01910 else if( strcmp(chKey4,"TRAC") == 0 )
01911 {
01912
01913 ParseTrace(chCard);
01914 }
01915
01916 else if( strcmp(chKey4,"TURB") == 0 )
01917 {
01918
01919 i = 5;
01920 DoppVel.TurbVel = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01921
01922
01923
01924
01925 if( nMatch(" NO ",chCard) && nMatch("PRES",chCard) )
01926 {
01927 DoppVel.lgTurb_pressure = false;
01928 }
01929 else
01930 {
01931
01932 DoppVel.lgTurb_pressure = true;
01933 }
01934
01935
01936 DoppVel.lgTurbEquiMag = false;
01937
01938 if( nMatch("EQUI",chCard) && nMatch("PART",chCard) )
01939 {
01940
01941 DoppVel.lgTurbEquiMag = true;
01942 }
01943
01944 if( nMatch(" LOG",chCard) )
01945 {
01946 DoppVel.TurbVel = (float)pow(10.f,DoppVel.TurbVel);
01947 }
01948
01949 DoppVel.TurbVel *= 1e5;
01950
01951
01952
01953
01954 DoppVel.Heiles_Troland_F = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
01955 if( lgEOL )
01956 {
01957
01958 DoppVel.Heiles_Troland_F = 3.f;
01959 }
01960
01961
01962
01963 if( nMatch("DISS",chCard) )
01964 {
01965 DoppVel.DispScale = (float)pow(10., FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL) );
01966 if( lgEOL )
01967 {
01968 NoNumb(chCard);
01969 }
01970 }
01971
01972
01973 if( optimize.lgVarOn )
01974 {
01975
01976 optimize.nvarxt[optimize.nparm] = 1;
01977 strcpy( optimize.chVarFmt[optimize.nparm], "TURBULENCE= %f log" );
01978
01979 optimize.nvfpnt[optimize.nparm] = input.nRead;
01980
01981 optimize.vparm[0][optimize.nparm] = (float)log10(DoppVel.TurbVel/
01982 1e5);
01983 optimize.vincr[optimize.nparm] = 0.5;
01984 optimize.nparm += 1;
01985 }
01986 }
01987
01988 else if( strcmp(chKey4,"WIND") == 0 )
01989 {
01990
01991
01992 ParseDynaWind(chCard);
01993 }
01994
01995 else if( strcmp(chKey2,"XI") == 0 )
01996 {
01997
01998
01999
02000 ParseIonPar(&nqh,chCard,&ar1);
02001 }
02002
02003
02004
02005 else if( !lgInputComment(chCard) )
02006 {
02007
02008 fprintf( ioQQQ, " Unrecognized command. Key=\"%4.4s\". This is routine ParseCommands.\n",
02009 chKey4 );
02010 fprintf( ioQQQ, " The line image was==%s==\n",
02011 chCard );
02012 fprintf( ioQQQ, " Sorry.\n" );
02013 puts( "[Stop in ParseCommands]" );
02014 cdEXIT(EXIT_FAILURE);
02015 }
02016
02017
02018
02019 input_readarray(chCard,&lgEOF);
02020 }
02021
02022
02023
02024
02025
02026
02027 for( i=0; i < INPUT_LINE_LENGTH; i++ )
02028 {
02029 chCard[i] = ' ';
02030 }
02031 chCard[INPUT_LINE_LENGTH-1] = '\0';
02032
02033 if( called.lgTalk )
02034 {
02035 fprintf( ioQQQ, " * %80.80s*\n",
02036 chCard );
02037 fprintf( ioQQQ, " ***********************************************************************************\n\n\n\n" );
02038 }
02039
02040
02041
02042
02043 if( optimize.lgTrOpt )
02044 {
02045
02046
02047 if( optimize.nTrOpt != optimize.nOptimiz )
02048 {
02049 trace.lgTrace = false;
02050
02051 trace.lgTrOvrd = false;
02052 }
02053 else
02054 {
02055 trace.lgTrace = true;
02056 called.lgTalk = true;
02057 trace.lgTrOvrd = true;
02058 fprintf( ioQQQ, " READR turns on trace from optimize option.\n" );
02059 }
02060 }
02061
02062
02063
02064
02066 if( dense.gas_phase[ipHYDROGEN] == -99. )
02067 {
02068 fprintf( ioQQQ, " Hydrogen density MUST be specified.\n" );
02069 lgStop_not_enough_info = true;
02070 lgStop = true;
02071 }
02072
02073
02074 dense.gas_phase[ipHYDROGEN] = (float)pow(10.f,dense.gas_phase[ipHYDROGEN]);
02075
02076 radius.rinner = radius.Radius;
02077
02078
02079 wind.emdot = dense.gas_phase[ipHYDROGEN]*wind.windv0;
02080
02081
02082 if( iterations.lgConverge_set)
02083 {
02084 iterations.itermx = MIN2( iterations.itermx , iterations.lim_iter );
02085 for( j=0; j < iterations.iter_malloc; j++ )
02086 {
02087 geometry.nend[j] = MIN2( geometry.nend[j] , iterations.lim_zone );
02088 }
02089 }
02090
02091
02092
02093 if( strcmp(dense.chDenseLaw,"DLW1") == 0 )
02094 {
02095 dense.gas_phase[ipHYDROGEN] = (float)dense_fabden(radius.Radius,radius.depth);
02096 }
02097 else if( strcmp(dense.chDenseLaw,"DLW2") == 0 )
02098 {
02099
02100 dense.gas_phase[ipHYDROGEN] =(float)dense_tabden(radius.Radius,radius.depth);
02101 }
02102
02103
02104
02105
02106
02107 lgStop_not_enough_info = false;
02108 lgStop = false;
02109
02110
02111 if( input.nSave < 0 )
02112 {
02113 fprintf( ioQQQ, " No commands were entered - whats up?\n" );
02114 puts( "[Stop in ParseCommands]" );
02115 cdEXIT(EXIT_FAILURE);
02116 }
02117
02118
02119 if( wind.windv0 > 0. && conv.lgAutoIt )
02120 {
02121 if( called.lgTalk )
02122 {
02123 fprintf( ioQQQ, " >>> PROBLEM Due to the nature of the Sobolev approximation, it makes no sense to converge a windy model.\n" );
02124 fprintf( ioQQQ, " >>> Iterate to convergence is turned off\n\n\n" );
02125 }
02126 conv.lgAutoIt = false;
02127 iterations.itermx = 0;
02128 }
02129
02130
02131
02132 if( opac.lgCaseB && conv.lgAutoIt && wind.windv0 >= 0. )
02133 {
02134 if( called.lgTalk )
02135 {
02136 fprintf( ioQQQ, " >>> Case B is an artificial test, it makes no sense to converge this model.\n" );
02137 fprintf( ioQQQ, " >>> Iterate to convergence is turned off.\n\n\n" );
02138 }
02139 conv.lgAutoIt = false;
02140 iterations.itermx = 0;
02141 }
02142
02143
02144 if( dense.DensityPower!=0. && strcmp( dense.chDenseLaw, "CPRE" )==0 )
02145 {
02146 if( called.lgTalk )
02147 {
02148 fprintf( ioQQQ, " Specifying both a density power law and constant pressure is impossible.\n" );
02149 }
02150 lgStop = true;
02151 }
02152
02153 if( !rfield.lgIonizReevaluate && strcmp( dense.chDenseLaw, "CDEN" )!=0 )
02154 {
02155 if( called.lgTalk )
02156 {
02157 fprintf( ioQQQ, " NO REEVALUATE IONIZATION can only be used with constant density.\n" );
02158 fprintf( ioQQQ, " Resetting to reevaluate ionization.\n\n" );
02159 }
02160 rfield.lgIonizReevaluate = true;
02161 }
02162
02163
02164 thickness = min( MIN3( StopCalc.tauend, StopCalc.colpls, StopCalc.colnut ),
02165 MIN3( StopCalc.col_h2, StopCalc.col_h2_nut, StopCalc.HColStop ) );
02166 if( thickness < COLUMN_INIT )
02167 {
02168
02169 thickness /= (dense.gas_phase[ipHYDROGEN]*geometry.FillFac);
02170
02171 if( thickness > 1e25 && radius.router[0] > 1e25 )
02172 {
02173 fprintf( ioQQQ,
02174 ">>> The specified column density and hydrogen density correspond to a thickness of %.2e cm.\n",
02175 thickness);
02176 fprintf( ioQQQ,
02177 ">>> This seems large to me.\n");
02178 fprintf(ioQQQ,">>> a very large radius may cause overflow.\n\n");
02179 }
02180 }
02181
02182
02183 if( !rfield.lgDoLineTrans && rfield.lgOpacityFine )
02184 {
02185 if( called.lgTalk )
02186 {
02187 fprintf( ioQQQ, " NO LINE TRANSER set but fine opacities still computed.\n" );
02188 fprintf( ioQQQ, " Turning off fine opacities.\n\n" );
02189 }
02190 rfield.lgOpacityFine = false;
02191 }
02192
02193
02194 if( h2.lgH2ON && (!rfield.lgDoLineTrans || !rfield.lgOpacityFine) )
02195 {
02196 if( called.lgTalk )
02197 {
02198 fprintf( ioQQQ, " Large H2 molecule turned on but line transfer and fine opacities are not.\n" );
02199 fprintf( ioQQQ, " Turning on line transfer and fine opacities.\n\n" );
02200 }
02201 rfield.lgOpacityFine = true;
02202 rfield.lgDoLineTrans = true;
02203 }
02204
02205
02206
02207 if( called.lgTalk && (StopCalc.tend < TEND ||
02208
02209 (thermal.ConstTemp > 0. && thermal.ConstTemp < TEND ) ) )
02210 {
02211
02212
02213 if( (hextra.cryden == 0.) && !co.lgNoCOMole && !hmi.lgNoH2Mole )
02214 {
02215 fprintf( ioQQQ, "\n >>>\n"
02216 " >>> The simulation is going into neutral gas but cosmic rays are not included.\n" );
02217 fprintf( ioQQQ, " >>> Ion-molecule chemistry will not occur without a source of ionization.\n" );
02218 fprintf( ioQQQ, " >>> The chemistry network may collapse deep in molecular regions.\n" );
02219 fprintf( ioQQQ, " >>> Add galactic background cosmic rays with the COSMIC RAYS BACKGROUND command.\n" );
02220 fprintf( ioQQQ, " >>> YOU ARE MAKING A BIG MISTAKE!!\n >>>\n\n\n\n" );
02221 }
02222 }
02223
02224
02228 if( called.lgTalk && dense.gas_phase[ipHYDROGEN] == 0 )
02229 {
02230 fprintf( ioQQQ, " PROBLEM DISASTER The hydrogen density has not been specified - use the HDEN command to do this.\n");
02231 lgStop_not_enough_info = true;
02232 }
02233 else if( called.lgTalk && dense.gas_phase[ipHYDROGEN] < 1e-4 )
02234 {
02235 fprintf( ioQQQ, " >>> Is the entered value of the hydrogen density (%.2e) reasonable?\n",
02236 dense.gas_phase[ipHYDROGEN]);
02237 fprintf( ioQQQ, " >>> It seems pretty low to me.\n\n\n" );
02238 }
02239 else if( called.lgTalk && dense.gas_phase[ipHYDROGEN] > 1e15 )
02240 {
02241 fprintf( ioQQQ, " >>> Is this value of the hydrogen density reasonable?\n" );
02242 fprintf( ioQQQ, " >>> It seems pretty high to me.\n\n\n" );
02243 }
02244
02245
02246 if( called.lgTalk && !lgStop && !lgStop_not_enough_info )
02247 {
02248 if( dense.gas_phase[ipHYDROGEN] < 1e-6 || dense.gas_phase[ipHYDROGEN] > 1e19 )
02249 {
02250 fprintf( ioQQQ, " >>> Simulation may crash because of extreme "
02251 "density. The value was %.2e\n\n" ,
02252 dense.gas_phase[ipHYDROGEN] );
02253 }
02254 }
02255
02256 if( rfield.nspec == 0 )
02257 {
02258 fprintf( ioQQQ, " PROBLEM DISASTER Type of continuum MUST be specified.\n" );
02259 lgStop = true;
02260 lgStop_not_enough_info = true;
02261 }
02262
02263 if( nqh == 0 )
02264 {
02265 fprintf( ioQQQ, " PROBLEM DISASTER Luminosity of continuum MUST be specified.\n" );
02266 lgStop = true;
02267 lgStop_not_enough_info = true;
02268 }
02269
02270
02271
02272
02273
02274 if( radius.Radius == 0. )
02275 {
02276 fprintf( ioQQQ, " PROBLEM DISASTER Starting radius MUST be specified.\n" );
02277 lgStop = true;
02278 lgStop_not_enough_info = true;
02279 }
02280
02281 if( rfield.nspec != nqh )
02282 {
02283 fprintf( ioQQQ, " PROBLEM DISASTER There were not the same number of continuum shapes and luminosities entered.\n" );
02284 lgStop = true;
02285 }
02286
02287 if( (optimize.nparm==0 && optimize.nOptimiz==0) && grid.lgGrid )
02288 {
02289
02290 fprintf( ioQQQ, " PROBLEM DISASTER The GRID command was entered but no parameters include a VARY option.\n" );
02291 fprintf( ioQQQ, " There is nothing for grid to vary.\n" );
02292 lgStop = true;
02293 }
02294
02295 if( lgStop_not_enough_info )
02296 {
02297 fprintf( ioQQQ, " PROBLEM DISASTER I do not have enough information to do the simulation, I cannot go on.\n" );
02298 fprintf( ioQQQ, "\n\n Sorry.\n\n\n" );
02299 puts( "[Stop in ParseCommands]" );
02300 cdEXIT(EXIT_FAILURE);
02301 }
02302
02303 if( lgStop )
02304 {
02305 puts( "[Stop in ParseCommands]" );
02306 cdEXIT(EXIT_FAILURE);
02307 }
02308
02309
02310
02311 DEBUG_EXIT( "ParseCommands()" );
02312
02313 return;
02314 }
02315