00001
00002
00003
00004 #include "cddefines.h"
00005 #define NPLXMX (LIMPAR*(LIMPAR+6)+1)
00006 #include "input.h"
00007 #include "called.h"
00008 #include "prt.h"
00009 #include "save.h"
00010 #include "optimize.h"
00011
00012
00013 bool lgOptimize_do(void)
00014 {
00015 long int i,
00016 iflag,
00017 ii,
00018 iworke[NPLXMX],
00019 j,
00020 mode,
00021 need,
00022 nfe,
00023 np;
00024 chi2_type fret;
00025 realnum fx,
00026 ptem[LIMPAR],
00027 delta[LIMPAR],
00028 toler,
00029 worke[NPLXMX];
00030
00031 DEBUG_ENTRY( "lgOptimize_do()" );
00032
00033
00034
00035
00036
00037
00038 toler = (realnum)log10(1. + optimize.OptGlobalErr);
00039
00040 if( strcmp(optimize.chOptRtn,"PHYM") == 0 )
00041 {
00042
00043 for( j=0; j < optimize.nvary; j++ )
00044 {
00045 ptem[j] = optimize.vparm[0][j];
00046 delta[j] = optimize.vincr[j];
00047 }
00048
00049 for( j=optimize.nvary; j < LIMPAR; j++ )
00050 {
00051 ptem[j] = -FLT_MAX;
00052 delta[j] = -FLT_MAX;
00053 }
00054 optimize_phymir(ptem,delta,optimize.nvary,&fret,toler);
00055 for( j=0; j < optimize.nvary; j++ )
00056 optimize.vparm[0][j] = ptem[j];
00057 }
00058
00059 else if( strcmp(optimize.chOptRtn,"SUBP") == 0 )
00060 {
00061 fprintf( ioQQQ, " Begin optimization with SUBPLEX\n" );
00062 need = 2*optimize.nvary + optimize.nvary*(optimize.nvary + 4) + 1;
00063 if( need > NPLXMX )
00064 {
00065 fprintf( ioQQQ, " Increase size of NPLXMX in parameter statements to handle this many variables.\n" );
00066 fprintf( ioQQQ, " I need at least %5ld\n", need );
00067 cdEXIT(EXIT_FAILURE);
00068 }
00069 for( j=0; j < optimize.nvary; j++ )
00070 ptem[j] = optimize.vparm[0][j];
00071
00072
00073
00074
00075
00076 mode = 0;
00077
00078
00079 optimize_subplex(
00080
00081 optimize.nvary,
00082
00083 toler,
00084
00085 optimize.nIterOptim,
00086
00087 mode,
00088
00089 optimize.vincr,
00090
00091 ptem,
00092
00093 &fx,
00094
00095 &nfe,
00096
00097
00098 worke,
00099
00100
00101 iworke,
00102
00103 &iflag);
00104
00105 if( iflag == -1 )
00106 {
00107 fprintf( ioQQQ, " SUBPLEX exceeding maximum iterations.\n"
00108 " This can be reset with the OPTIMZE ITERATIONS command.\n" );
00109 }
00110
00111 for( j=0; j < optimize.nvary; j++ )
00112 optimize.vparm[0][j] = ptem[j];
00113
00114 if( optimize.lgOptimFlow )
00115 {
00116 fprintf( ioQQQ, " trace return optimize_subplex:\n" );
00117 for( j=0; j < optimize.nvary; j++ )
00118 {
00119 fprintf( ioQQQ, " Values:" );
00120 for( ii=1; ii <= optimize.nvarxt[j]; ii++ )
00121 fprintf( ioQQQ, " %.2e", optimize.vparm[ii-1][j] );
00122 fprintf( ioQQQ, "\n" );
00123 }
00124 }
00125 }
00126 else
00127 TotalInsanity();
00128
00129
00130 called.lgTalk = cpu.i().lgMPI_talk();
00131 called.lgTalkIsOK = cpu.i().lgMPI_talk();
00132 prt.lgFaintOn = true;
00133
00134 if( called.lgTalk )
00135 {
00136 fprintf( ioQQQ, " **************************************************\n" );
00137 fprintf( ioQQQ, " **************************************************\n" );
00138 fprintf( ioQQQ, " **************************************************\n" );
00139 fprintf( ioQQQ, "\n Cloudy was called %4ld times.\n\n", optimize.nOptimiz );
00140
00141 for( i=0; i < optimize.nvary; i++ )
00142 {
00143 np = optimize.nvfpnt[i];
00144
00145
00146 if( optimize.nvarxt[i] == 1 )
00147 {
00148 sprintf( input.chCardSav[np], optimize.chVarFmt[i], optimize.vparm[0][i] );
00149 }
00150
00151 else if( optimize.nvarxt[i] == 2 )
00152 {
00153 sprintf( input.chCardSav[np], optimize.chVarFmt[i], optimize.vparm[0][i],
00154 optimize.vparm[1][i]);
00155 }
00156
00157 else if( optimize.nvarxt[i] == 3 )
00158 {
00159 sprintf( input.chCardSav[np], optimize.chVarFmt[i],
00160 optimize.vparm[0][i], optimize.vparm[1][i], optimize.vparm[2][i] );
00161 }
00162
00163 else if( optimize.nvarxt[i] == 4 )
00164 {
00165 sprintf( input.chCardSav[np], optimize.chVarFmt[i],
00166 optimize.vparm[0][i], optimize.vparm[1][i], optimize.vparm[2][i],
00167 optimize.vparm[3][i] );
00168 }
00169
00170 else if( optimize.nvarxt[i] == 5 )
00171 {
00172 sprintf( input.chCardSav[np], optimize.chVarFmt[i],
00173 optimize.vparm[0][i], optimize.vparm[1][i], optimize.vparm[2][i],
00174 optimize.vparm[3][i], optimize.vparm[4][i]);
00175 }
00176
00177 else
00178 {
00179 fprintf(ioQQQ,"The number of variable options on this line makes no sense to me.\n");
00180 cdEXIT(EXIT_FAILURE);
00181 }
00182
00183 fprintf( ioQQQ, " Optimal command: %s\n", input.chCardSav[np]);
00184 fprintf( ioQQQ, " Smallest value:%10.2e Largest value:%10.2e Allowed range %10.2e to %10.2e\n",
00185 optimize.varmin[i], optimize.varmax[i], optimize.varang[i][0],
00186 optimize.varang[i][1] );
00187 }
00188 }
00189
00190 if( cpu.i().lgMaster() )
00191 {
00192
00193 FILE *ioOptim = open_data( chOptimFileName, "w", AS_LOCAL_ONLY );
00194 for( i=0; i <= input.nSave; i++ )
00195 fprintf( ioOptim, "%s\n", input.chCardSav[i]);
00196 fclose( ioOptim );
00197
00198
00199
00200 fprintf( ioQQQ, "\f" );
00201
00202 for( i=0; i < optimize.nvary; i++ )
00203 ptem[i] = optimize.vparm[0][i];
00204
00205 (void)optimize_func(ptem);
00206 }
00207
00208 return lgAbort;
00209 }