*! 3.0.0  30Jul97 (Jeroen Weesie) STB-43 ip25
*     --  added support for parameters
*     --  replication counter is available to simulation program in scalar REPL  
*  2.0.1  by StataCorp

program define simul2
        version 5.0

      * parse command line
       
        parse "`*'", parse(" ,")
        local prog "`1'"
        mac shift
        if "`prog'" == "" {
                di in red "program name required"
                exit 198
        }
        local options "Args(str) Dots Par(str) Reps(int -1) report"
        parse "`*'"

      * scratch
      
        tempname hsim ic ncond nsim 
        tempfile result

      * check argument -reps-
                
        if `reps' < 1 {
                di in red "reps() required"
                exit 198
        }
        scalar `nsim' = `reps'

      * output -dots-     
      
        if "`dots'" == "" {
                local dots "*"
        }
        else    local dots ""
        
      * parse parameters

        scalar `ncond' = 1        
        local np 0
        if "`par'" != "" {
                * split -par- into -par1- -par2- etc
                parse "`par'", p(",")
                while "`1'" != "" {
                        if "`1'" != "," {
                                local np = `np' + 1
                                local par`np' "`1'"
                        }
                        mac shift
                }
                
                /* split par<i> into parameter name (pname<i>) and 
                     numlist-expanded parameter values (pvalues<i>) 
                   nvalue<i> is set to the number of levels of <i>.
                   ncond is set to the product of number of levels, i.e.,
                     the number of combinations of parameter values.
                */                   
                local ip 1
                while `ip' <= `np' {
                        parse "`par`ip''", p(" ")
                        local pname`ip' "`1'"      /* name of parameter */
                        mac shift
                        capt numlist "`*'", real min(1)
                        if _rc {
                                di in re "error in interpreting `par`ip''"
                                exit 198
                        }                                
                        local values`ip' "$S_1"    /* values for parameter */
                        local nvalue`ip' "$S_2"    /* number of values */
                        scalar `ncond' = `ncond' * `nvalue`ip''
                        local ip = `ip'+1
                } 
        } 
        
      * verify that enough memory is available
       
      * quietly desc
      * if _result(4)-20 < `nsim' * `ncond' {
      *         di in red "insufficient memory (observations)"
      *         exit 901
      * }

      * obtain and verify names of saved results
      
        global S_1
        capture `prog' ? `args'
        if _rc { 
                if _rc == 199 { 
                        di in red "program `prog' not found or"
                        di in red _quote "`prog' ?" _quote /* 
                          */ " attempted to execute an unrecognized command."
                        exit 199
                }
                di in red _quote "`prog' ?" _quote " returned:"
                error _rc
        }

        * verify that names of results-to-be-saved are defined
        
        if "$S_1" == "" {
                di in red _quote "`prog' ?" _quote " did not set \$S_1"
                exit 198
        }
                      
        * verify that the first names in S_1 coincide with the parameter names
        if `np' > 0 {
                local nS : word count $S_1
                if `nS' <= `np' {
                        di in re "S_1 should define parameters and names" /* 
                              */ " for saved results"
                        exit 198
                }                        
                local ip 1
                while `ip' <= `np' {
                        local p : word `ip' of $S_1
                        if "`p'" != "`pname`ip''" {
                                di in re "error in \$S_1: " 
                                di in re "name of result `ip' should equal "/*
                                    */ "the `ip'-th parameter (`pname`ip'')"
                                exit 198
                        }                                
                        local ip = `ip' + 1
                }
        }
       
      * display description of simulation
        
        if "`report'" != "" {
                di _n in gr "Simulation program: " in ye "`prog'"
                if `np' == 0 {
                        di in gr "no parameters specified"
                        di in gr "number of replication:" in ye `nsim'                       
                }
                else {
                        di _n in gr "Parameter  #values  values"
                        local ip 1
                        while `ip' <= `np' {
                                di in gr "`pname`ip''" in ye /*
                                  */ _col(14) %2.0f `nvalue`ip'' /*
                                  */ _col(21) "`values`ip''"
                                local ip = `ip' + 1
                        }  
                        di
                        di in gr "number of conditions:       " in ye `ncond'
                        di in gr "replications per condition: " in ye `nsim'
                        di
                }                        
        }
      
      * run simulation     
      
        * open file for simulation data with file handle hsim
        postfile `hsim' $S_1 using "`result'"
    
        if `np' == 0 {
                /* 
                   simple code to run simulation without parameters
                */                        
              
                scalar REPL = 1
                while REPL <= `nsim' {
                        `dots' di in gr "." _c
                        quietly `prog' `hsim' `args'
                        scalar REPL = REPL + 1
                }
                `dots' di _n
        }
        else {
              /* 
                   code to run simulation with "variable" number of parameters 
              */                        

              * initialize (I1,I2,I3,..) to (0,1,1,..)
              local I1 0
              local ip 2
              while `ip' <= `np' {
                        local I`ip' 1
                        local ip = `ip' + 1
              }                        
              
              * loop over all conditions
              scalar `ic' = 1
              while `ic' <= `ncond' {

                        * update condition (I1,I2,...)
                        local ip 1
                        while `I`ip'' == `nvalue`ip'' {
                                * note that this while loop cannot explode
                                local ip = `ip' + 1
                        }
                        local iip 1
                        while `iip' < `ip' {
                                local I`iip' 1
                                local iip = `iip' + 1
                        }
                        local I`ip' = `I`ip'' + 1                                                               
                        
                        * update parameters using scalars
                        local ip 1
                        while `ip' <= `np' {                      
                                local v : word `I`ip'' of `values`ip''
                                scalar `pname`ip'' = `v'
                                local ip = `ip' + 1
                        }                                
                                                                         
                        * display parameters associated with the condition
                        if "`dots'" == "" {
                                local ip 1
                                di in gr "[", _c
                                while `ip' < `np' {
                                        di in gr scalar(`pname`ip'') "," _c
                                        local ip = `ip'+1
                                }                                
                                di in gr scalar(`pname`np'') "] " _c 
                        }       
            
                        * loop for repetitions              
                        scalar REPL = 1
                        while REPL <= `nsim' {
                                `dots' di in gr "." _c
                                quietly `prog' `hsim' `args'
                                scalar REPL = REPL + 1
                        }
                        `dots' di 
                        
                        scalar `ic' = `ic' + 1
              }
        }
        
      * clean up screen and load file with simulation results
       
        postclose `hsim'
        qui use "`result'", clear
        global S_FN
        scalar drop REPL
end