00001
00002
00003
00004
00005 #include "cddefines.h"
00006 #include "taulines.h"
00007 #include "iso.h"
00008 #include "phycon.h"
00009 #include "mole.h"
00010 #include "elementnames.h"
00011 #include "dynamics.h"
00012 #include "stopcalc.h"
00013 #include "dense.h"
00014 #include "iterations.h"
00015 #include "colden.h"
00016 #include "punch.h"
00017 #include "rt.h"
00018 #include "conv.h"
00019
00020
00021
00022 void ConvIterCheck( void )
00023 {
00024 bool lgConverged;
00025 long int nelem,
00026 i,
00027 ipISO,
00028 ipHi, ipLo;
00029
00030 DEBUG_ENTRY( "ConvIterCheck()" );
00031
00032
00033
00034
00035
00036
00037
00038
00039 lgConverged = true;
00040 strcpy( conv.chNotConverged, "Converged!" );
00041 if( iteration > 1 && conv.lgAutoIt )
00042 {
00043 for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
00044 {
00045 for( nelem=ipISO; nelem < LIMELM; nelem++ )
00046 {
00047 if( dense.lgElmtOn[nelem] )
00048 {
00049
00050
00051
00052 if(ipISO==ipH_LIKE )
00053 {
00054 ipHi = ipH3p;
00055 ipLo = ipH2s;
00056 }
00057 else if( ipISO==ipHE_LIKE )
00058 {
00059 ipHi = ipHe2p3P2;
00060 ipLo = ipHe2s3S;
00061 }
00062 else
00063
00064 TotalInsanity();
00065
00066
00067
00068
00069 if( Transitions[ipISO][nelem][ipHi][ipLo].Emis->TauIn > 0.5 )
00070 {
00071
00072 if( fabs(Transitions[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Emis->TauTot -
00073 Transitions[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Emis->TauIn*rt.DoubleTau) >
00074 conv.autocv*fabs(Transitions[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Emis->TauIn*rt.DoubleTau) )
00075 {
00076
00077 lgConverged = false;
00078
00079
00080
00081 sprintf( conv.chNotConverged, "%s-like Lya",elementnames.chElementSym[ipISO] );
00082
00083 if( punch.lgPunConv )
00084 {
00085 fprintf( punch.ipPunConv, " %s-like Lya thick, "
00086 "nelem= %s iteration %li old %.3e new %.3e \n",
00087 elementnames.chElementSym[ipISO] ,
00088 elementnames.chElementSym[nelem],
00089 iteration,
00090 Transitions[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Emis->TauTot ,
00091 Transitions[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Emis->TauIn);
00092 }
00093 }
00094
00095 if( fabs(Transitions[ipISO][nelem][ipHi][ipLo].Emis->TauTot -
00096 Transitions[ipISO][nelem][ipHi][ipLo].Emis->TauIn*rt.DoubleTau) >
00097 conv.autocv*fabs(Transitions[ipISO][nelem][ipHi][ipLo].Emis->TauIn*rt.DoubleTau) )
00098 {
00099
00100 lgConverged = false;
00101
00102
00103
00104 sprintf( conv.chNotConverged, "%s-like subord",elementnames.chElementSym[ipISO] );
00105
00106 if( punch.lgPunConv )
00107 {
00108 fprintf( punch.ipPunConv, " %s-like subord, nelem= %s iteration %li old %.3e new %.3e \n" ,
00109 elementnames.chElementSym[ipISO],
00110 elementnames.chElementSym[nelem],
00111 iteration,
00112 Transitions[ipISO][nelem][ipHi][ipLo].Emis->TauTot ,
00113 Transitions[ipISO][nelem][ipHi][ipLo].Emis->TauIn );
00114 }
00115 }
00116 }
00117 }
00118 }
00119 }
00120
00121
00122
00123 for( i=0; i<NCOLD; ++i )
00124 {
00125
00126
00127 if( colden.colden[i]/colden.colden[ipCOL_HTOT] > 1e-5 &&
00128 fabs(colden.colden_old[i]-colden.colden[i]) > conv.autocv*colden.colden[i] )
00129 {
00130
00131 lgConverged = false;
00132
00133
00134
00135 strcpy( conv.chNotConverged, "H mole col" );
00136
00137 if( punch.lgPunConv )
00138 {
00139 fprintf( punch.ipPunConv, " H mole col species %li iteration %li old %.2e new %.2e H col den %.2e\n",
00140 i,iteration,
00141 colden.colden_old[i],
00142 colden.colden[i],
00143 colden.colden[ipCOL_HTOT] );
00144 }
00145 }
00146 }
00147
00148
00149
00150 for( i=0; i<mole.num_comole_calc; ++i )
00151 {
00152 if(COmole[i]->n_nuclei == 1)
00153 continue;
00154
00155
00156 if( COmole[i]->hevcol/colden.colden[ipCOL_HTOT] > 1e-5 &&
00157 fabs(COmole[i]->hevcol_old-COmole[i]->hevcol) > conv.autocv*COmole[i]->hevcol )
00158 {
00159
00160 lgConverged = false;
00161
00162
00163
00164 strcpy( conv.chNotConverged, "CO mol col" );
00165
00166
00167
00168 if( punch.lgPunConv )
00169 {
00170 fprintf( punch.ipPunConv, "CO mol col, old:%.3e new:%.3e\n" ,
00171 COmole[i]->hevcol_old ,
00172 COmole[i]->hevcol );
00173 }
00174 }
00175 }
00176
00177
00178 if( dynamics.lgAdvection )
00179 {
00180
00181 if( dynamics.convergence_error > conv.autocv*dynamics.error_scale2*dynamics.convergence_tolerance ||
00182 dynamics.discretization_error > conv.autocv*dynamics.error_scale2 )
00183 {
00184 lgConverged = false;
00185
00186
00187 strcpy( conv.chNotConverged, "Dynamics " );
00188 if( punch.lgPunConv )
00189 {
00190 fprintf( punch.ipPunConv, " Dynamics\n" );
00191 }
00192 }
00193 }
00194
00195 if( punch.lgPunConv && lgConverged )
00196 {
00197 fprintf( punch.ipPunConv, " converged\n" );
00198 }
00199
00200
00201 if( lgConverged )
00202 iterations.itermx = MIN2(iterations.itermx,iteration);
00203
00204
00205
00206 if( phycon.te < StopCalc.tend && nzone == 1 )
00207 {
00208 lgConverged = true;
00209 strcpy( conv.chNotConverged, " " );
00210 iterations.itermx = MIN2(iterations.itermx,iteration);
00211 }
00212 }
00213 return;
00214 }