/*ConvFail handle conergece failure */ #include #include "cddefines.h" #include "elecden.h" #include "zonecnt.h" #include "called.h" #include "phycon.h" #include "showme.h" #include "heat.h" #include "tcool.h" #include "dtemp.h" #include "limfal.h" #include "failz.h" #include "pressure.h" #include "mappar.h" #include "prtzone.h" #include "punt.h" #include "converge.h" void ConvFail(char chMode[5]) { double relerror ; # ifdef DEBUG_FUN fputs( "<+>ConvFail()\n", debug_fp ); # endif /* pressure failure */ if( strcmp( chMode , "pres" )==0 ) { /* record number of pressure failures */ conv.nPreFail += 1; if( called.lgTalk ) { fprintf( ioQQQ, " convfail: pressure not converged; zone%4ld Te:%11.2e Hden:%10.2e curr Pres:%10.2e Corr Pres%10.2e\n", ZoneCnt.nzone, phycon.te, phycon.hden, pressure.pnow, pressure.pcorec ); } } /* electron ensity failure */ else if( strcmp( chMode, "eden" ) == 0 ) { conv.nNeFail += 1; if( called.lgTalk ) { fprintf( ioQQQ, " convfail: electron density fails, correct=%11.3e " " assumed=%11.3e zone=%3ld.", ElecDen.EdenTrue, phycon.eden, ZoneCnt.nzone ); /* some extra information that may be printed */ /* heating cooling failure */ if( !conv.lgConvTemp ) { fprintf( ioQQQ, " Temperature failure also." ); } /* heating cooling failure */ if( !conv.lgConvIoniz ) { fprintf( ioQQQ, " Ionization failure also." ); } /* electron density is oscillating */ if( conv.lgEdenOscl ) { fprintf( ioQQQ, " Electron density oscillating." ); } /* heating cooling oscillating */ if( conv.lgCmHOsc ) { fprintf( ioQQQ, " Heating-cooling oscillating." ); } } fprintf( ioQQQ, " \n"); } else if( strcmp( chMode, "ioni" ) == 0 ) { /* ionization failure */ conv.nIonFail += 1; if( called.lgTalk ) { fprintf( ioQQQ, " convfail: ionization not converged %3ld Z=%4ld \n", conv.nIonFail, ZoneCnt.nzone ); } } /* rest of routine is temperature failure */ else if( strcmp( chMode, "temp" ) == 0 ) { ++conv.nTeFail; if( called.lgTalk ) { fprintf( ioQQQ, " convfail: temp fail %3ld Z=%4ld TE=%10.4e ERROR=%10.2e" "HTOT=%10.2e CTOT=%10.2e DERR=%10.2e DTe=%10.2e lgIonDone:%1c\n", conv.nTeFail, ZoneCnt.nzone, phycon.te, (heat.htot - tcool.ctot)/ heat.htot, heat.htot, tcool.ctot, tcool.dCoolmHeatdT, dtemp.dTemper, TorF(conv.lgConvIoniz) ); if( conv.lgTOscl && conv.lgCmHOsc ) { fprintf( ioQQQ, " Temperature and d(Cool-Heat)/dT were both oscillating.\n" ); } else if( conv.lgTOscl ) { fprintf( ioQQQ, " Temperature was oscillating.\n" ); } else if( conv.lgCmHOsc ) { fprintf( ioQQQ, " d(Cool-Heat)/dT was oscillating.\n" ); } /* not really a temperature failure, but something else */ if( !conv.lgConvIoniz ) { fprintf( ioQQQ, " Solution not converged due to %10.10s\n", conv.chConvIoniz ); } } /* remember which zone this is */ failz.ifailz[MIN2(conv.nTeFail,10)-1] = ZoneCnt.nzone; /* remember the relative error * convert to single precision for following max, abs (vax failed here) */ relerror = fabs((heat.htot-tcool.ctot)/ heat.htot); conv.failmx = (float)MAX2(conv.failmx,relerror); /*conv.failmx = (float)fmax(conv.failmx,fabs((float)((heat.htot-tcool.ctot)/ heat.htot)));*/ if( conv.nTeFail < limfal.LimFail ) { # ifdef DEBUG_FUN fputs( " <->ConvFail()\n", debug_fp ); # endif return; } /* punt */ if( called.lgTalk ) { fprintf( ioQQQ, " Stop after%3ld temperature failures.\n", limfal.LimFail ); fprintf( ioQQQ, " This limit can be reset with the FAILURES command.\n" ); /* adjust range of punting map */ MapPar.RangeMap[0] = (float)(phycon.te/100.); MapPar.RangeMap[1] = (float)MIN2(phycon.te*100.,9e9); if( !limfal.lgNoMap ) { /* need to make printout out now, before disturbing soln with map */ PrtZone(); punt("punt"); } } /* return out from here and let lgEndFun catch lgBusted set, * and generate normal output there */ called.lgBusted = TRUE; if( called.lgTalk ) { fprintf( ioQQQ, " ConvFail sets lgBusted since nTeFail>=LimFail%5ld%5ld\n", conv.nTeFail, limfal.LimFail ); } } else { fprintf( ioQQQ, " ConvFail called with insane mode %s \n", chMode ); ShowMe(); puts( "[Stop in ConvFail]" ); exit(1); } # ifdef DEBUG_FUN fputs( " <->ConvFail()\n", debug_fp ); # endif return; }