cloudy
trunk
|
00001 /* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and 00002 * others. For conditions of distribution and use see copyright notice in license.txt */ 00003 /*ParseHDEN parse the hden command */ 00004 #include "cddefines.h" 00005 #include "input.h" 00006 #include "dense.h" 00007 #include "optimize.h" 00008 #include "parse.h" 00009 00010 void ParseHDEN(char *chCard ) 00011 { 00012 bool lgEOL; 00013 long int i; 00014 00015 DEBUG_ENTRY( "ParseHDEN()" ); 00016 00017 if( dense.gas_phase[ipHYDROGEN] > 0. ) 00018 { 00019 fprintf( ioQQQ, " PROBLEM DISASTER More than one density command was entered.\n" ); 00020 cdEXIT(EXIT_FAILURE); 00021 } 00022 00023 /* log of hydrogen density */ 00024 i = 5; 00025 dense.gas_phase[ipHYDROGEN] = (realnum)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL); 00026 if( lgEOL ) 00027 { 00028 fprintf( ioQQQ, " DISASTER The density MUST be entered with this command. STOP\n" ); 00029 cdEXIT(EXIT_FAILURE); 00030 } 00031 00032 /* check for further options */ 00033 if( ! nMatch("LINE",chCard) ) 00034 { 00035 /* check size of density - will we crash? */ 00036 if( dense.gas_phase[ipHYDROGEN] > log10(FLT_MAX) || 00037 dense.gas_phase[ipHYDROGEN] < log10(FLT_MIN) ) 00038 { 00039 fprintf(ioQQQ, 00040 " DISASTER - the log of the entered hydrogen density is %.3f - too extreme for this processor.\n", 00041 dense.gas_phase[ipHYDROGEN]); 00042 if( dense.gas_phase[ipHYDROGEN] > 0. ) 00043 fprintf(ioQQQ, 00044 " DISASTER - the log of the largest hydrogen density this processor can do is %.3f.\n", 00045 log10(FLT_MAX) ); 00046 else 00047 fprintf(ioQQQ, 00048 " DISASTER - the log of the smallest hydrogen density this processor can do is %.3f.\n", 00049 log10(FLT_MIN) ); 00050 fprintf(ioQQQ," Sorry.\n" ); 00051 00052 cdEXIT(EXIT_FAILURE); 00053 } 00054 00055 dense.gas_phase[ipHYDROGEN] = (realnum)pow((realnum)10.f,dense.gas_phase[ipHYDROGEN]); 00056 } 00057 00058 if( dense.gas_phase[ipHYDROGEN] <= 0. ) 00059 { 00060 fprintf( ioQQQ, " PROBLEM DISASTER Hydrogen density must be > 0.\n" ); 00061 cdEXIT(EXIT_FAILURE); 00062 } 00063 00064 /* this is the linear initial density */ 00065 dense.den0 = dense.gas_phase[ipHYDROGEN]; 00066 00067 /* check if exponent entered */ 00068 dense.DensityPower = (realnum)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL); 00069 /* this branch when exponent was entered - do something with it */ 00070 if( !lgEOL ) 00071 { 00072 /* not constant density 00073 * some sort of power law density distribution */ 00074 if( nMatch("COLU",chCard) ) 00075 { 00076 /* density will depend on column density to a power 00077 * number entered is col den, convert to scale radius 00078 * at this point HDEN is LOG10 of hydrogen density */ 00079 dense.rscale = (realnum)pow(10.,FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL)); 00080 if( lgEOL ) 00081 { 00082 fprintf( ioQQQ, " The column density MUST be set if the col den option is to be used.\n" ); 00083 cdEXIT(EXIT_FAILURE); 00084 } 00085 strcpy( dense.chDenseLaw, "POWC" ); 00086 } 00087 else if( nMatch("DEPT",chCard) ) 00088 { 00089 /* depth option, sets scale radius, log cm */ 00090 dense.rscale = (realnum)pow(10.,FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL)); 00091 if( lgEOL ) 00092 { 00093 fprintf( ioQQQ, " The scale depth MUST be set if the depth option is to be used.\n" ); 00094 cdEXIT(EXIT_FAILURE); 00095 } 00096 strcpy( dense.chDenseLaw, "POWD" ); 00097 } 00098 else 00099 { 00100 /* radius option, will be relative to inner radius */ 00101 strcpy( dense.chDenseLaw, "POWR" ); 00102 } 00103 } 00104 00105 /* vary option */ 00106 if( optimize.lgVarOn ) 00107 { 00108 /* pointer to where to write */ 00109 optimize.nvfpnt[optimize.nparm] = input.nRead; 00110 optimize.vparm[0][optimize.nparm] = (realnum)log10(dense.gas_phase[ipHYDROGEN]); 00111 optimize.vincr[optimize.nparm] = 1.; 00112 00113 /* these are varios options for density laws, 00114 * first is constant density or pressre*/ 00115 if( strcmp(dense.chDenseLaw ,"CDEN") == 0 || 00116 strcmp(dense.chDenseLaw ,"CPRE") == 0 ) 00117 { 00118 strcpy( optimize.chVarFmt[optimize.nparm], "HDEN=%f" ); 00119 optimize.nvarxt[optimize.nparm] = 1; 00120 } 00121 00122 /* power law density distrution */ 00123 else if( strcmp(dense.chDenseLaw,"POWR") == 0 ) 00124 { 00125 strcpy( optimize.chVarFmt[optimize.nparm], "HDEN=%f, power=%f" ); 00126 optimize.vparm[1][optimize.nparm] = dense.DensityPower; 00127 optimize.nvarxt[optimize.nparm] = 2; 00128 } 00129 00130 /* power law with density scale depending on column density */ 00131 else if( strcmp(dense.chDenseLaw,"POWC") == 0 ) 00132 { 00133 strcpy( optimize.chVarFmt[optimize.nparm], "HDEN=%f, power=%f, column=%f" ); 00134 optimize.nvarxt[optimize.nparm] = 3; 00135 optimize.vparm[1][optimize.nparm] = dense.DensityPower; 00136 optimize.vparm[2][optimize.nparm] = (realnum)log10(dense.rscale); 00137 } 00138 00139 /* power law with density scale depending on depth */ 00140 else if( strcmp(dense.chDenseLaw,"POWD") == 0 ) 00141 { 00142 strcpy( optimize.chVarFmt[optimize.nparm], "HDEN=%f, power=%f, depth=%f" ); 00143 optimize.nvarxt[optimize.nparm] = 3; 00144 optimize.vparm[1][optimize.nparm] = dense.DensityPower; 00145 optimize.vparm[2][optimize.nparm] = (realnum)log10(dense.rscale); 00146 } 00147 00148 /* could not idensity an option */ 00149 else 00150 { 00151 fprintf( ioQQQ, " Internal error in HDEN\n" ); 00152 cdEXIT(EXIT_FAILURE); 00153 } 00154 ++optimize.nparm; 00155 } 00156 return; 00157 }