/home66/gary/public_html/cloudy/c08_branch/source/grid_do.cpp

Go to the documentation of this file.

/* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and
 * others.  For conditions of distribution and use see copyright notice in license.txt */
/*grid_do called by cdDrive, this returns 0 if things went ok, 1 for disaster */
#include "cddefines.h"
#include "conv.h"
#include "input.h"
#include "called.h"
#include "version.h"
#include "init.h"
#include "prt.h"
#include "trace.h"
#include "grains.h"
#include "parse.h"
#include "punch.h"
#include "optimize.h"
#include "grid.h"

/*grid_do called by cdDrive, calls gridXspec, returns false if ok, true for disaster */
bool grid_do(void)
{
        char chLine[INPUT_LINE_LENGTH], 
          chNote[8];
        long int i, 
          ii, 
          j;
        realnum ptem[LIMPAR]; 

        DEBUG_ENTRY( "grid_do()" );

        /* main driver for optimization runs
         * Drives cloudy to grid variables;*/

        /* code originally written by R.F. Carswell, IOA Cambridge */

        /* variables with optimizer */
        for( i=0; i < LIMPAR; i++ )
        {
                optimize.OptIncrm[i] = 0.;
                optimize.varang[i][0] = -FLT_MAX;
                optimize.varang[i][1] = FLT_MAX;
                /* this should be overwritten by format of vary line */
                strcpy( optimize.chVarFmt[i], "error - no optimizer line image was set" );
        }

        /* necessary to do this to keep all lines in */
        prt.lgFaintOn = false;
        conv.LimFail = 1000;

        /* this initializes variables at the start of each simulation
        * in a grid, before the parser is called - this must set any values
        * that may be changed by the command parser */
        InitDefaultsPreparse();

        /* call READR the first time to scan off all variable options */
        /* this is just an initial parsing to get the number of iterations and
         * the number of varied parameters.  The other Init* routines are not 
         * called after this because this is all done again later for each grid point */
        ParseCommands();

        /* done parsing, now turn off punch headers. */
        punch.lgPunHeader = false;

        /* >>chng 06 may 20, do not turn off printout if grid calculation is in place
         * rather than optimization */
        if( !grid.lgGrid )
        {
                called.lgTalk = false;
                /* this flag is needed to turn print on to have effect */
                called.lgTalkIsOK = false;
        }

        /* >>chng 00 aug 09, return memory allocated for grains, they are not used, PvH */
        ReturnGrainBins();

        optimize.nvary = optimize.nparm;

        /* check that more than 1 observed intensities or column densities were entered */
        if( ((optimize.nlobs + (long)optimize.lgOptLum + optimize.nTempObs + optimize.ncobs) < 1 ) && 
                !grid.lgGrid )
        {
                fprintf( ioQQQ, " The input stream has vary commands, but\n" );
                fprintf( ioQQQ, " no observed quantities were entered.  Whats up?\n" );
                fprintf( ioQQQ, " Use the NO VARY command to input vary options but not try to perform this.\n" );
                cdEXIT(EXIT_FAILURE);
        }

        /* check that the total number of parameters to vary is greater than 1 */
        if( optimize.nvary < 1 )
        {
                fprintf( ioQQQ, " No parameters to vary were entered. Whats up?\n" );
                cdEXIT(EXIT_FAILURE);
        }

        if( ( strcmp(optimize.chOptRtn,"XSPE") == 0 ) && ( optimize.nRangeSet != optimize.nvary ) )
        {
                fprintf( ioQQQ, " Every parameter with a VARY option must have a GRID specified,\n" );
                fprintf( ioQQQ, " and the GRID must be specified after the VARY option.\n" );
                fprintf( ioQQQ, " These requirements were not satisfied for %ld parameter(s).\n",  abs(optimize.nvary - optimize.nRangeSet) );
                cdEXIT(EXIT_FAILURE);
        }

        /* lgTrOptm set with trace grid command */
        if( trace.lgTrOptm )
        {
                for( i=0; i < optimize.nvary; i++ )
                {
                        /*print the command format as debugging aid */
                        fprintf( ioQQQ, "%s\n", optimize.chVarFmt[i]);

                        /* now generate the actual command with parameter,
                         * there will be from 1 to 3 numbers on the line */
                        if( optimize.nvarxt[i] == 1 )
                        {
                                /* case with 1 parameter */
                                sprintf( chLine , optimize.chVarFmt[i], optimize.vparm[0][i] );
                        }

                        else if( optimize.nvarxt[i] == 2 )
                        {
                                /* case with 2 parameter */
                                sprintf( chLine , optimize.chVarFmt[i], optimize.vparm[0][i], optimize.vparm[1][i]);
                        }

                        else if( optimize.nvarxt[i] == 3 )
                        {
                                /* case with 3 parameter */
                                sprintf( chLine , optimize.chVarFmt[i], 
                                        optimize.vparm[0][i], optimize.vparm[1][i] , optimize.vparm[2][i] );
                        }

                        else if( optimize.nvarxt[i] == 4 )
                        {
                                /* case with 4 parameter */
                                sprintf( chLine , optimize.chVarFmt[i], 
                                        optimize.vparm[0][i], optimize.vparm[1][i] , optimize.vparm[2][i], optimize.vparm[3][i] );
                        }

                        else if( optimize.nvarxt[i] == 5 )
                        {
                                /* case with 5 parameter */
                                sprintf( chLine , optimize.chVarFmt[i], 
                                        optimize.vparm[0][i], optimize.vparm[1][i] , optimize.vparm[2][i],
                                        optimize.vparm[3][i], optimize.vparm[4][i]);
                        }

                        else
                        {
                                fprintf(ioQQQ,"The number of variable options on this line makes no sense to me1\n");
                                cdEXIT(EXIT_FAILURE);
                        }

                        /* print the resulting command line*/
                        fprintf( ioQQQ, "%s\n", chLine );
                }
        }

        /* option to change default increments; if zero then leave as is */
        for( i=0; i < LIMPAR; i++ )
        {
                if( optimize.OptIncrm[i] != 0. )
                {
                        optimize.vincr[i] = optimize.OptIncrm[i];
                }
        }

        /* say who we are */
        if( strcmp(optimize.chOptRtn,"XSPE") == 0 )
        {
                fprintf( ioQQQ, "                                                          Grid  Driver\n" );
        }
        else
        {
                fprintf( ioQQQ, "                                                      Optimization  Driver\n" );
        }
        int indent = (int)((122 - strlen(t_version::Inst().chVersion))/2);
        fprintf( ioQQQ, "%*cCloudy %s\n\n",indent,' ',t_version::Inst().chVersion);
        fprintf( ioQQQ, "%23c**************************************%7.7s**************************************\n", ' ', t_version::Inst().chDate );
        fprintf( ioQQQ, "%23c*%81c*\n", ' ', ' ' );

        /* now echo initial input quantities with flag for vary */
        /* first loop steps over all command lines entered */
        for( i=0; i <= input.nSave; i++ )
        {
                /* put space to start line, overwrite if vary found */
                strcpy( chNote, "       " );
                /* loop over all vary commands, see if this is one */
                for( j=0; j < optimize.nvary; j++ )
                {
                        /* input.nSave is on C array counting, rest are on fortran */
                        if( i == optimize.nvfpnt[j] )
                        {
                                /* this is a vary command, put keyword at start */
                                strcpy( chNote, "VARY>>>" );
                        }
                }

                fprintf( ioQQQ, "%22.7s * %-80s*\n", chNote, input.chCardSav[i] );
        }
        fprintf( ioQQQ, "%23c*%81c*\n", ' ', ' ' );
        fprintf( ioQQQ, "%23c***********************************************************************************\n\n\n", ' ' );

        /* option to trace logical flow within this sub */
        if( optimize.lgOptimFlow )
        {
                for( j=0; j < optimize.nvary; j++ )
                {
                        i = optimize.nvfpnt[j];
                        fprintf( ioQQQ, " trace:%80.80s\n", input.chCardSav[i]);
                        fprintf( ioQQQ, "%80.80s\n", optimize.chVarFmt[j]);
                        fprintf( ioQQQ, " number of variables on line:%4ld\n", 
                          optimize.nvarxt[j] );
                        fprintf( ioQQQ, " Values:" );
                        for( ii=1; ii <= optimize.nvarxt[j]; ii++ )
                        {
                                fprintf( ioQQQ, "%10.2e", optimize.vparm[ii-1][j] );
                        }
                        fprintf( ioQQQ, "\n" );
                }
        }


        if( strcmp(optimize.chOptRtn,"PHYM") == 0 )
        {
                fprintf( ioQQQ, " Up to%5ld iterations will be performed,\n", 
                  optimize.nIterOptim );
                fprintf( ioQQQ, " and the final version of the input file will be written to the file %s\n", 
                  chOptimFileName );

                fprintf( ioQQQ, " The optimize_phymir method will be used" );
                if( optimize.lgParallel ) {
                        fprintf( ioQQQ, " in parallel mode.\n The maximum no. of CPU's to be used is %1ld.\n",optimize.useCPU );
                }
                else {
                        fprintf( ioQQQ, " in sequential mode.\n" );
                }
        }

        else if( strcmp(optimize.chOptRtn,"SUBP") == 0 )
        {
                fprintf( ioQQQ, " Up to%5ld iterations will be performed,\n", 
                  optimize.nIterOptim );
                fprintf( ioQQQ, " and the final version of the input file will be written to the file %s\n", 
                  chOptimFileName );

                fprintf( ioQQQ, " Subplex method will be used.\n" );
        }

        else if( strcmp(optimize.chOptRtn,"XSPE") == 0 )
        {
                fprintf( ioQQQ, " Producing grid output.\n" );
        }

        else
        {
                fprintf( ioQQQ, " I do not understand what method to use.\n" );
                fprintf( ioQQQ, " Sorry.\n" );
                cdEXIT(EXIT_FAILURE);
        }

        fprintf( ioQQQ, "\n%4ld parameter(s) will be varied.  The first lines, and the increments are:\n", 
          optimize.nvary );

        for( i=0; i < optimize.nvary; i++ )
        {
                optimize.varmax[i] = -FLT_MAX;
                optimize.varmin[i] = FLT_MAX;
                /*  write formatted to output using the format held in chVarFmt(np) */

                /* now generate the actual command with parameter,
                 * there will be from 1 to 3 numbers on the line */
                if( optimize.nvarxt[i] == 1 )
                {
                        /* case with 1 parameter */
                        sprintf( chLine , optimize.chVarFmt[i], optimize.vparm[0][i] );
                }

                else if( optimize.nvarxt[i] == 2 )
                {
                        /* case with 2 parameter */
                        sprintf( chLine , optimize.chVarFmt[i], optimize.vparm[0][i], optimize.vparm[1][i]);
                }

                else if( optimize.nvarxt[i] == 3 )
                {
                        /* case with 3 parameter */
                        sprintf( chLine , optimize.chVarFmt[i], 
                                optimize.vparm[0][i], optimize.vparm[1][i] , optimize.vparm[2][i] );
                }

                else if( optimize.nvarxt[i] == 4 )
                {
                        /* case with 4 parameter */
                        sprintf( chLine , optimize.chVarFmt[i], 
                                optimize.vparm[0][i], optimize.vparm[1][i] , optimize.vparm[2][i] , optimize.vparm[3][i] );
                }

                else if( optimize.nvarxt[i] == 5 )
                {
                        /* case with 5 parameter */
                        sprintf( chLine , optimize.chVarFmt[i], 
                                optimize.vparm[0][i], optimize.vparm[1][i] , optimize.vparm[2][i],
                                 optimize.vparm[3][i] , optimize.vparm[4][i]);
                }

                else
                {
                        fprintf(ioQQQ,"The number of variable options on this line makes no sense to me2\n");
                        cdEXIT(EXIT_FAILURE);
                }

                fprintf( ioQQQ, "\n %s\n", chLine );
                fprintf( ioQQQ, " Initial increment is%6.3f, the limits are%10.2e to %10.2e\n", 
                  optimize.vincr[i], optimize.varang[i][0], optimize.varang[i][1] );
        }

        /* this will be number of times grid calls cloudy */
        optimize.nOptimiz = 0;
        /* >>chng 06 jan 26, moved this to an option on the command line */
        /* grid.numParamValues = 10; */

        if( strcmp(optimize.chOptRtn,"XSPE") == 0 )
        {
                for( j=0; j < optimize.nvary; j++ )
                {
                        /* ptem[j] = optimize.vparm[0][j]; */
                        ptem[j] = optimize.varang[j][0]; 
                        /* >> chng 06 sep 4, delta is now read at command line. */
                        /* delta[j] = optimize.vincr[j]; */
                        /* delta[j] = ( optimize.varang[j][1] - optimize.varang[j][0] )/
                                ((realnum)grid.numParamValues[j] - 1.f); */
                }
                /* >>chng 06 aug 23, fill rest of array with zeros. */
                for( j=optimize.nvary; j < LIMPAR; j++ )
                {
                        ptem[j] = 0.f; 
                        grid.paramIncrements[j] = 0.f;
                }

                gridXspec(ptem,optimize.nvary);
                for( j=0; j < optimize.nvary; j++ )
                {
                        optimize.vparm[0][j] = ptem[j];
                }

                fprintf( ioQQQ, " **************************************************\n" );
                fprintf( ioQQQ, " **************************************************\n" );
                fprintf( ioQQQ, " **************************************************\n" );
                fprintf( ioQQQ, "\n Cloudy was called %4ld times.\n\n", optimize.nOptimiz );

                /* though a page eject */
                fprintf( ioQQQ, "\f" );
        }
        else
        {
                lgAbort = lgOptimize_do();
        }

        if( lgAbort )
        {
                /* busted set means there were serious problems somewhere */
                return 1;
        }
        else
        {
                /* return 0 if everything is ok */
                return 0;
        }
}

---