#!/bin/sh  
# \
exec oagtclsh "$0" "$@"

#
# $Log: not supported by cvs2svn $
# Revision 1.15  2005/08/25 16:26:44  shang
# modified to work with any format rootname
#
# Revision 1.14  2005/06/21 15:36:54  shang
# removed the post-processing to avoid writing to the same proc file that other processor is writing to and replaced file copy by exec cp  and exec mv
#
# Revision 1.13  2005/03/29 22:01:50  shang
# added ReplaceHihgPoint proc to replace the high point of the simplex and removed
# the obselete codes
#
# Revision 1.12  2005/03/22 04:10:05  borland
# More tolerant of processing failures now.
#
# Revision 1.11  2005/03/18 06:06:25  borland
# Now uses rootname_ for filenames for 1D scan instead of rootname.  Makes it
# the same as the simplex filenames.
#
# Revision 1.10  2005/03/17 06:33:54  borland
# Did previous change correctly.
#
# Revision 1.9  2005/03/17 03:59:31  borland
# Use 1.0e200 when an invalid value is seen (second instance).
#
# Revision 1.8  2005/03/16 21:33:59  shang
# fixed a bug in createInitailSimplex proc and set the value to a big number when
# the returned value is not a numerical value.
#
# Revision 1.7  2005/01/07 16:00:10  shang
# added processors option to limit the maximum number of processors used at the same time.
#
# Revision 1.6  2004/11/22 16:45:44  shang
# modified CheckFileExistent proc to avoid the error from sdds2stream unitl correct
# results were obtained within the timeout limit.
#
# Revision 1.5  2004/11/18 17:28:35  shang
# added the columns of *.proc to the log file
#
# Revision 1.4  2004/11/17 18:41:46  shang
# added CheckFileExistent proc to check the existence of expected file, use
# sddscheck checking the existence and sdds2stream -rows to check if the expected
# file is empty.
#
# Revision 1.3  2004/11/10 16:34:22  shang
# It now skips the evaluation for the parameter values that has been evaluated before.
#
# Revision 1.2  2004/11/09 22:49:45  shang
# added logFile oneDScan options and added checking the existence of expected
# files, if no-existent after every 100 seconds, an message will be printed out.
#
# Revision 1.1  2004/10/26 19:31:49  shang
# first version
#
# parallize simplex optimization, Hairong Shang
# modified based on geneticoptimizer

set auto_path [linsert $auto_path 0  /usr/local/oag/apps/lib/$env(HOST_ARCH)]
set auto_path [linsert $auto_path 0 /usr/local/oag/lib_patch/$env(HOST_ARCH)]
APSStandardSetup

set usage {usage: parallelOptimizer -input <filename> [-errorFactor <number>] [-help 1] [-tolerance <value>] [-target <value>] [<-replaceHighPoint> 1] [<-reflThroughLowPoint> 1] [-restartFile <filename>] [-compareWithPrevBest <1>] [-logFile <filename>] [-processors <integer>] [-restart <number>] [-iterations <number>]}

proc PrintHelp {} {
    global usage
    puts stdout "$usage\n"
    puts stdout "Description of commandline parameters:"
    puts stdout "input:  name of the input file, as described above."
    puts stdout "logFile: name of file for storing the result and parameter values of each simulation."
    puts stdout "tolerance:  convergence critia:"
    puts stdout "        if tolerance is <=0, fractional tolerance is used for evaluation "
    puts stdout "        if tolerance is >0, absolute tolerance (i.e., abs(HighResult-lowResult) ) is used for evaluation."
    puts stdout "target: if target is given and the best result obtained so far is less or equal to target, optimization is done."
    puts stdout "errorFactor: factor by which to multiply the errorLevels from the input file."
    puts stdout "replaceHighPoint: default is 0, if it is given as 1, then the high point of the"
    puts stdout "                  contraction simplex will be replaced by the overall minimum point"
    puts stdout "                  of reflection, shrinkage and expansion simplexes and only if"
    puts stdout "                  this minimum is higher than the low point but lower"
    puts stdout "                  than the high point in the contraction simplex. You may want to"
    puts stdout "                  do so if the high point is too big compared with the low point"
    puts stdout "reflThroughLowPoint: default reflecting point is through the centroid of original"
    puts stdout "                     simplex. If this value is 1, the reflecting point will be"
    puts stdout "                     the low point of the original simplex, which converges"
    puts stdout "                     if the low point is close to the desired position"
    puts stdout "restartFile: the best parameters will be saved into this file before program "
    puts stdout "             terminates, which can be used as input in next run"
    puts stdout "compareWithPrevBest: the default value is 0, which means that the difference between high point"
    puts stdout "                     and low point in the same simplex is used for convergence test; if it is 1,"
    puts stdout "                     the difference between the low point of new simplex and the previous best result"
    puts stdout "                     is used instead."
    puts stdout "sum:      default is 1, which uses smaller reflection factor (name as factor0."
    puts stdout "          if set sum as 0, the reflection factor is dimensions*factor0, where"
    puts stdout "          dimensions is the number of independent variables. The latter factor"
    puts stdout "          turns out to have better results and fast convergence for some mathmatical"
    puts stdout "          functions such as powell, trid-dixon functions."
    puts stdout "oneDScan: default is 0, no one dimensional scan will be perform, the initial simplex"
    puts stdout "          is created by adding the errorLevel to the initialValue" 
    puts stdout "          if it is 1, the initial simplex is created by one-dimensional scan, i.e.,"
    puts stdout "          a low point in each dimension is one of the vertexes in the initial simplex."
    puts stdout "processors: maximum number of processors to be used at the same time. If not given, the number of needed processors will be used."
    puts stdout "iterations: maximum allowed iterations allowed in one cycle."
    puts stdout "restart: if >0, optimizer restarts <restart> times if the iterations of one cycle exceeds the iterations."
    puts stdout ""
    puts stdout "Description of input file:"
    puts stdout "Parameters: "
    puts stdout "preProcessingScript: script to run prior to running the job. (May be blank.)"
    puts stdout "  The script must accept the following arguments and syntax: "
    puts stdout "    -rootname <string>       rootname for the job."
    puts stdout "    -tagList <string>        list of tags for varied quantities."
    puts stdout "    -valueList <string>      list of values associated with tags."
    puts stdout "    The tag names are the same as the parameter names given in the"
    puts stdout "    parameterName column (see below)."
    puts stdout "  E.g., the script will be called with arguments like these:"
    puts stdout "    -rootname run-000000 -tagList \"Var1 Var2\" -valueList \"1.7 -2.8\""
    puts stdout "postProcessingScript: script to run to postprocess a job.  Called with"
    puts stdout "  two arguments giving the rootname and runID of the job.  This script must produce"
    puts stdout "  a one-row SDDS file named <rootname>.proc that contains at least the "
    puts stdout "  following columns:"
    puts stdout "  runName --- the rootname of the run for this file"
    puts stdout "  runID --- the ID of this job."
    puts stdout "  penaltyValue --- the value of the penalty function"
    puts stdout "  This file must also have all the columns from <rootname>.inp sddsxref'd in."
    puts stdout "runScript: script to run the simulation job.  The script must accept the"
    puts stdout "  all the arguments as the preProcessingScript, plus runID."
    puts stdout "  it also accepts preProcessingScript and postProcessingScript arguments."
    puts stdout "  pre/post ProcessingScript arguments are optional."
    puts stdout "  The purpose of adding pre/post ProcessingScript to runScript"
    puts stdout "  is to parallel the pre/post processing instead of running in sequential"
    puts stdout "  If runScript is blank, then you"
    puts stdout "  must give inputFileTemplate, inputFileExtension, and programName."
    puts stdout "  The script must create a file <rootname>.run as soon as it starts, "
    puts stdout "   then run preProcessingScript (if given)."
    puts stdout "   If postProcessingscript is given, run it after all other tasks are done"
    puts stdout "      and create <rootname>.proc"
    puts stdout "   otherwise, create <rootname>.proc directly."
    puts stdout "   at last, create <rootname>.done right before finishing."
    puts stdout "inputFileTemplate: name of the template file for creating input files."
    puts stdout "  Sequences of the form <name> are substituted with values."
    puts stdout "inputFileExtension: extension to use for input files."
    puts stdout "programName: name of the program to run with the input files.  The"
    puts stdout "  input filename is given as the sole argument of the program."
    puts stdout ""
    puts stdout "Columns:"
    puts stdout "parameterName: names of the parameters to vary in the optimization."
    puts stdout "  These should appear in the input file template in the form <name>."
    puts stdout "initialValue: initial values of the parameters"
    puts stdout "errorLevel: rms width of the gaussian random errors to add to the"
    puts stdout "  parameter values as mutations."
    puts stdout "lowerLimit: smallest allowable values for the parameters"
    puts stdout "upperLimit: largest allowable values for the parameters."
}

set args $argv
set input ""
set logFile ""
set logFileID ""
set tolerance 0.00001
set sum 1
set iterations 1000
set restart 0
set target 0
set targetSupplied 0
set errorFactor 1
set reduce 0
set help 0
set replaceHighPoint 0
set reflThroughLowPoint 0
set compareWithPrevBest 0
set restartFile ""
set oneDScan 0
set processors 0
if {[APSStrictParseArguments {logFile input errorFactor reduce help reflThroughLowPoint \
                                tolerance target replaceHighPoint restart iterations \
                                restartFile compareWithPrevBest sum oneDScan processors}] || \
        ![string length  $input] || $errorFactor<=0 || $tolerance<=0} {
    if $help {
        PrintHelp
        exit 
    }    
    return -code error "$usage"
}

if {$restart<0} {
    puts stderr "Invalid restart value given: restart value is less than 0"
    exit
}
set totalIterations [expr ($restart+1)*$iterations]

if {$target} {
    set targetSupplied 1
}

if ![file exists $input] {
    return -code error "not found: $input"
}
if {[string length $logFile] && [file exists $logFile] } {
    puts stderr "log file $logFile already exists, exit."
    exit
}
set rootname [file rootname $input]

if [file exist $rootname.run] {
    return -code error "rootname in use: $rootname"
}

proc SubmitJob {args} {
    set code ""
    set input "" 
    APSStrictParseArguments {code input}

    eval file delete -force $input.log [file rootname $input].done 
    set tmpFile [file rootname $input].csh
    set fd [open $tmpFile w]
    puts $fd "#!/bin/csh "
    puts $fd "source ~/.cshrc"
    puts $fd "echo running $code $input on [exec uname -a]"
    puts $fd "$code $input >& $input.log"
    close $fd
    exec chmode +x $tmpFile
    catch {exec csub -noemail 1 -name [file root [file tail $input]] $tmpFile} result
    return "$result"
}

proc SubmitRunScript {args} {
    set script ""
    set tagList ""
    set valueList ""
    set rootname ""
    set runID ""
    set preProcessingScript ""
    set postProcessingScript ""
    APSStrictParseArguments {script rootname tagList valueList runID preProcessingScript postProcessingScript}
    global env
    eval file delete -force $rootname.log $rootname.done $rootname.run $rootname.proc
    set tmpFile [file rootname $rootname].csh
    APSAddToTempFileList $tmpFile
    set fd [open $tmpFile w]
    puts $fd "#!/bin/bash "
    puts $fd "source ~/.cshrc"
    puts $fd "echo running $script $rootname $tagList $valueList on [exec uname -a]"
    puts $fd "$script -rootname $rootname -tagList '$tagList' -valueList '$valueList' -runID $runID -preProcessingScript '$preProcessingScript' -postProcessingScript '$postProcessingScript' >& $rootname.log"
    close $fd
    catch {exec cat $tmpFile | qsub -V -cwd -j y -N [file root [file tail $rootname]] } result
    return "$result"
}

set inputFileTemplate ""
set inputFileExtension ""
proc LoadAndCheckInputFile {} {
    global input inputParameterList inputColumnList inputFileTemplate
    if [catch {sdds load $input inputData} result] {
        return -code error "$result"
    }
    set inputParameterList $inputData(ParameterNames)
    foreach param $inputParameterList {
        global $param
        set $param [lindex $inputData(Parameter.$param) 0]
    }
    set inputColumnList $inputData(ColumnNames)
    foreach col $inputColumnList {
        global $col
        set $col [lindex $inputData(Column.$col) 0]
    }
    foreach newParam {runScript preProcessingScript} {
        if ![info exists $newParam] {
            global $newParam
            lappend inputParameterList $newParam
            set $newParam ""
        }
    }
    if [string length $runScript] {
        if [string length $programName] {
            return -code error "Give runScript or programName, but not both!"
        }
        set inputFileExtension run
    }
    if [string length $programName] {
        if ![string length $inputFileTemplate] {
            return -code error "inputFileTemplate is required for runProgram!"
        }
        if ![string length $inputFileExtension] {
            return -code error "inputFileExtension is required for runProgram!"
        }
    }
}


proc CreateInitialMatrix {args} {
    #for N dimension independent variables, create P0[N+1][N] array, while for i=1 to N
    # P[i]=P[0]+delta[i] delta[i] is the change of ith variable
    # P[0] is the initial value of N variables i.e. P[0][1],P[0][2], ..., P[0][N]
    # this array is named as pVector
    global inputParameterList inputColumnList simplexPoints
    eval global $inputParameterList 
    eval global $inputColumnList
    global rootname template errorFactor jobsStarted pVector pExist
    
    puts stderr "create initial simplex..."
    set simplexPoints 0
    if ![info exist pVector(parameterName)] {
        set pVector(parameterName) ""
        set pVector(Dimensions) [llength $parameterName]
        foreach pN $parameterName { lappend pVector(parameterName) <$pN> }
    }
    for {set i 0} {$i <= $pVector(Dimensions)} {incr i} {
        set pVector($i.Value) ""
    }
    for {set point 0} {$point <= $pVector(Dimensions)} {incr point} {
	set pVector($point.limit) 0
        set pExist($point) 0
        set index 1
        foreach pN $parameterName iV $initialValue eL $errorLevel lL $lowerLimit uL $upperLimit {
	    if {$iV<$lL || $iV>$uL || $uL<$lL} {
                puts stderr "Invalid initial value - $iV given for $pN, it exceeds the $lL - $uL."
                exit
            }
            if {$point==0 || $index != $point } {
                set value $iV
            } else {
                set value [expr $eL*$errorFactor+$iV]
                if {$value <= $lL} {
                    set el [expr ($uL-$lL)/12]
                    set value [expr $iL + $el]
                } elseif {$value >= $uL} {
                    set el [expr ($uL-$lL)/12]
                    set value [expr $uL-$el]
                }
            }
            lappend pVector($point.Value) $value
            lappend pVector($point.procOpt) -define=column,$pN,$value,type=double
            incr index
        }
    }
}

proc CreateInitialMatrix1 {args} {
    #for N dimension independent variables, create P0[N+1][N] array, while for i=1 to N
    # P[i]=P[0]+delta[i] delta[i] is the change of ith variable
    # P[0] is the initial value of N variables i.e. P[0][1],P[0][2], ..., P[0][N]
    # this array is named as pVector
    global inputParameterList inputColumnList 
    eval global $inputParameterList 
    eval global $inputColumnList
    global rootname template errorFactor jobsStarted pVector
    
    puts stderr "create initial simplex..."
    if ![info exist pVector(parameterName)] {
        set pVector(parameterName) ""
        set pVector(Dimensions) [llength $parameterName]
        foreach pN $parameterName { lappend pVector(parameterName) <$pN> }
    }
    for {set i 0} {$i <= $pVector(Dimensions)} {incr i} {
        set pVector($i.Value) ""
    }
    for {set point 0} {$point <= $pVector(Dimensions)} {incr point} {
        set index 1
        foreach pN $parameterName iV $initialValue eL $errorLevel lL $lowerLimit uL $upperLimit {
            if {$point==0 || $index != $point } {
                set value $iV
            } else {
                set value [expr $eL*$errorFactor+$iV]
                if [expr $uL==$lL] {set value $uL}
                if {$value < $lL} {set value $lL}
                if {$value > $uL} {set value $uL}
            }
            lappend pVector($point.Value) $value
            lappend pVector($point.procOpt) -define=column,$pN,$value,type=double
            incr index
        }
    }
}

proc SubmitOneJob {args} {
    set runName ""
    set procOpt ""
    set replOpt ""
    set valueList ""
    set id ""
    APSParseArguments {runName procOpt replOpt valueList id}
    global inputParameterList inputColumnList
    eval global $inputParameterList 
    eval global $inputColumnList
    
    set origOpt [concat <rootname> $parameterName]
    eval exec sddssequence -pipe -define=dummy,type=short -sequence=begin=0,end=0,number=1 \
      | sddsprocess -pipe=in $runName.inp \
      $procOpt \
      -print=column,runName,$runName \
      -define=column,runID,$id,type=long 
    if [string length $runScript] {
        catch {SubmitRunScript -script $runScript -rootname $runName -runID $id \
                 -tagList "$parameterName" -valueList "$valueList" \
                 -preProcessingScript "$preProcessingScript" \
                 -postProcessingScript "$postProcessingScript" } result
    } else {
        if [string length $preProcessingScript] {
            if [catch {exec $preProcessingScript \
                         -rootname $runName \
                         -tagList "$parameterName" \
                         -valueList "$valueList"} result] {
                return -code error "$result"
            }
            if [string length $result] {
                puts stderr "$result"
            }
        }
        exec replace -original=[join $origOpt ,] -replacement=[join $replOpt ,] \
          $inputFileTemplate $runName.$inputFileExtension
        catch {SubmitJob -code $programName -input $runName.$inputFileExtension \
             } result
    }
}

set jobsrunning 0
set jobList ""
set doneList ""
proc DistributeJobs {args} {
    global inputParameterList inputColumnList processors jobList doneList
    eval global $inputParameterList 
    eval global $inputColumnList
    global template errorFactor jobsStarted
    global pVector rootname reflVector expanVector contractVector
    global pExist reflExist expanExist contractExist

    set option ""
    #option could be {p refl expan and contract}
    APSParseArguments {rootname option}
    set jobs [expr $pVector(Dimensions) +1]
    for {set job 0} {$job< $jobs} {incr job} {
        if [set ${option}Exist($job)] {
            puts stderr "$option: id $job already exist, skip"
            continue
        }
        set runName [format ${rootname}_${option}-%06ld $job]
        set fileList [glob -nocomplain ${runName}.*]
        eval file delete -force file $fileList
	after 1000
        set origOpt [concat <rootname> $pVector(parameterName)]
        set valueList [set ${option}Vector($job.Value)]
        set replOpt [concat $runName $valueList]
        set procOpt ""
        foreach pN $parameterName val $valueList {
            lappend procOpt -redefine=column,$pN,$val,type=double
        }
	if [set ${option}Vector($job.limit)] {
            set options -defaultType=double
            foreach pN $parameterName val $valueList {
                lappend options "-col=$pN"
                lappend options "-data=$val"
            }
            lappend options "-col=penaltyValue"
            lappend options "-data=1.0e200"
	    exec touch $runName.done
            if [catch {eval exec sddsmakedataset $runName.proc $options} result] {
                puts stderr "Unable to create $runName.proc file for beyond limit point"
                exit
            }
        } else {
	    lappend jobList "SubmitOneJob -runName $runName -procOpt \"$procOpt\" \
                  -replOpt \"$replOpt\" -valueList \"$valueList\" -id $job"
	    lappend doneList $runName.done
	}
    }
}

proc SubmitParallelJobs {args} {
    global jobList doneList processors
    if ![llength $jobList] {
        puts stderr "No jobs for submitting!"
        exit
    }
    set jobs [llength $jobList]
    set submitted 0
    set jobsrunning 0
    foreach done $doneList {
        set $done 0
        set $done.checked 0
    }
    if {$processors>0} {
        while {$submitted<$jobs} {
            foreach job $jobList done $doneList {
                if {$jobsrunning<$processors && ![set $done]} {
                    eval $job
                    incr submitted
                    incr jobsrunning
                    set $done 1
                    puts stderr "job [file root $done] submitted."
                }
                if {[file exist $done] && ![set $done.checked]} {
                    set jobsrunning [expr $jobsrunning -1]
                    set $done.checked 1
                }
            }
            after 1000
        }
    } else {
        foreach job $jobList {
            eval $job
        }
    }
    set jobList ""
    set doneList ""
}

set reflectionFactors {-0.5 -1 -2}
set reflectionJobs [llength $reflectionFactors]
set steps 0
proc StartOptimize {} {
    global inputParameterList inputColumnList reflectionJobs reflectionFactors steps 
    eval global $inputParameterList 
    eval global $inputColumnList
    global rootname template errorFactor jobsStarted nTotalJobs pLowResult 
    global bestResult pVector target tolerance targetSupplied pHighResult bestValues pLowID
    global contractLowID contractHighID expanLowID expanHighID reflLowID reflHighID
    global contractVector expanVector reflVector pHighID replaceHighPoint oneDScan
    global iterations restart noprogress

    puts stderr "[clock format [clock seconds] -format %H:%M:%S]: First simplex run ...: "
    if $oneDScan {
        puts stderr "run one-d scan to create initial simplex..."
        OneDScan
    } else {
        CreateInitialMatrix
        set jobs [expr $pVector(Dimensions)+1]
        puts stderr "evaluating initial simplex ..."
        if [catch {DistributeJobs -option p} result] {
            puts stderr $result
            exit
        }
        SubmitParallelJobs
        AnalyzeResults1 -jobs $jobs -option p -step [expr $steps + 1] -operation initial_simplex
    }
    set jobs [expr $pVector(Dimensions)+1]
    incr steps
    set result [CheckTolerance -tolerance $tolerance -target $target -targetSupplied $targetSupplied]
  
    if {$result} {
        PrintBestResult
        exit
    }
    set factorList {-0.5 -1 -2}
    set noprogress 0
    set iter 0
    while {1} {
	incr iter
        set replace 0
        CalculatePSum
        #reflection
        puts stderr "reflecting simplex ..."
        foreach lamda $factorList option {contract refl expan} {
            RunParallelOptions -lamda $lamda -option $option
        }
        puts stderr "shrink toward the best point..."
        ContractAroundBestPoint -jobs $jobs
        SubmitParallelJobs
        catch {eval file delete -force file ${rootname}_p.results ${rootname}_contract.results ${rootname}_refl.results ${rootname}_expan.results}
    
        foreach option {contract refl expan p} {
            set result_$option [AnalyzeResults1 -option $option -jobs $jobs -step [expr $steps + 1]]
        }
        set min $result_contract
        set min_opt contract
        foreach opt {refl expan} {
            set result [set result_$opt]
            if {$result<$min} {
                set min $result
                set min_opt $opt
            }
        }
        puts stderr "min=$min, min_opt=$min_opt, pHigh=$pHighResult"
        if {$min<$bestResult && $min<$result_p} {
            puts stderr "replace previous simplex with $min_opt simplex.."
            ReplaceSimplex -repl $min_opt
            set replace 1
        } elseif {$min!=$result_p && $min<$pHighResult} {
            if {$replaceHighPoint} {
                puts stderr "shrinkage holds, replace the high point"
                ReplaceHighPoint -repl $min_opt
            }
        }
        incr steps
        puts stderr "checking convergency..."
        set final [CheckTolerance -tolerance $tolerance -target $target -targetSupplied $targetSupplied -replace $replace]
        if {$final} {
            PrintBestResult -result 1
            exit
        } else {
            if {$iterations && !$restart && $iter>$iterations} {
                PrintBestResult -result 2
                exit
            } 
            if {$restart && $iterations && $iter>$iterations} {
		set restart [expr $restart-1]
                global initialValue bestValues
                set initialValue $bestValues
                puts stderr "[exec date +%H:%M:%S] No progress after $noprogress iterations: restart optimization."
                CreateInitialMatrix
                set jobs [expr $pVector(Dimensions)+1]
                puts stderr "evaluating initial simplex ..."
                if [catch {DistributeJobs -option p} result] {
                    puts stderr $result
                    exit
                }
                set noprogress 0
		set iter 0
                SubmitParallelJobs
                AnalyzeResults1 -jobs $jobs -option p -step [expr $steps + 1] -operation restart_simplex
                incr steps
                set result [CheckTolerance -tolerance $tolerance -target $target -targetSupplied $targetSupplied]
                if $result {
                    PrintBestResult -result 1
                    exit
                }
                
            }
        }
    }
}

set bestResult -10
set noprogress 0
proc CheckTolerance {args} {
    set tolerance 0
    set targetSupplied 0
    set target 0
    set replace 0
    APSParseArguments {tolerance targetSupplied target replace}
     
    global pLowResult pHighResult bestResult steps pLowID pVector bestValues prevBest compareWithPrevBest pHighID
    global noprogress
    if {$bestResult<0} {
        set bestResult $pVector($pLowID.fValue)
        set prevBest $pVector($pHighID.fValue)
        set bestValues $pVector($pLowID.Value)
    } elseif {$pLowResult < $bestResult} {
        set prevBest $bestResult
        set bestResult $pVector($pLowID.fValue)
        set bestValues $pVector($pLowID.Value)
        set noprogress 0
    } else {
        incr noprogress
    }
    puts stderr "noprogress: $noprogress"
    puts stderr "steps: $steps, p high result: $pHighResult; p low result: $pLowResult, prev best=$prevBest"
    if {$targetSupplied && $pLowResult<=$target} {
            return 1
    }
    if {$compareWithPrevBest} {
        set compareWith $prevBest
    } else {
        set compareWith $pHighResult
    }
    if {$compareWith==$pLowResult} {
	return 0
    }
    if {$tolerance<=0} {
        set denominator [expr (abs($compareWith) + abs($bestResult))/2]
        if {$denominator} {
            set merit [expr abs($bestResult-$compareWith)/$denominator]
        } else {
            puts stderr "Error: divide-by-zero in fractional tolerance evaluation (simplexMin)\n"
            exit
        }
    } else {
        set merit [expr abs($compareWith-$bestResult)]
    }
    if {$merit<=$tolerance} {
        return 1
    } else {
        return 0
    }
}


proc AnalyzeResults1 {args} {
    set jobs 0
    set option p
    set step ""
    set operation ""
    APSParseArguments {jobs option step operation}
    
    global pLowID pHighID reflLowID reflHighID expanLowID expanHighID contractLowID contractHighID
    global pLowResult pHighResult reflLowResult reflHighResult 
    global bestResult pVector reflVector rootname expanVector contractVector bestValues
    global postProcessingScript programName parameterName logFile

    set jobs [expr $pVector(Dimensions)+1]
    for {set i 0} {$i<$jobs} {incr i} {
        set root [format ${rootname}_${option}-%06ld $i]
        lappend rootList $root
	lappend idList $i
        set $root 0
    }
    set count 0
    while {1} {
        set notdone 0
        foreach root $rootList ID $idList {
	    if [file exists ${root}.done] {
		if [file exist $root.proc] {
		    if ![CheckFileExistent -filename $root.proc] {
			set notdone 1
		    }
		} else {
		    set notdone 1
		}
            } else {
                set notdone 1
                if {$count>100} {
                    puts stderr "[clock format [clock seconds] -format %H:%M:%S]: job $root is not done yet."
                    set count 0
                }
            }
        }
        if {!$notdone} {
            break
        }
        after 1000
        incr count
    }
    after 1000
    foreach root $rootList ID $idList {
        if [catch {CheckFileExistent -filename $root.proc} result] {
            puts stderr $result
            exit
        }
	if [catch {exec sdds2stream $root.proc -col=penaltyValue} results] {
            puts stderr "unable to get the penalty value of $root.proc, $results, [exec sddscheck $root.proc]"
            exit
        }
        set ${option}Vector($ID.fValue) [lindex $results 0]
        if [catch {expr [set ${option}Vector($ID.fValue)] / 2} result] {
	    set ${option}Vector($ID.fValue) 1.0e200
        }
        puts stderr "$option ID=$ID done, penaltyValue=[set ${option}Vector($ID.fValue)], parameter values are:"
        set index 0
        foreach pN $parameterName {
            puts stderr "        $pN:  [lindex [set ${option}Vector($ID.Value)] $index]"
            incr index
        }
    }
    
    set lowID 0
    set highID 0
    set lowResult [set ${option}Vector(0.fValue)]
    set highResult $lowResult
    for {set i 1} {$i < $jobs} {incr i} {
        set val [set ${option}Vector($i.fValue)]
        if {$val<$lowResult} {
            set lowID $i
            set lowResult $val
        }
        if {$val>$highResult} {
            set highID $i
            set highResult $val
        }
    }
    set ${option}LowID $lowID
    set ${option}HighID $highID
    set ${option}LowResult  $lowResult
    set ${option}HighResult $highResult
    puts stderr "option=$option: HighID=$highID, LowID=$lowID"
    puts stderr "HighResult=[set ${option}HighResult], LowResult=[set ${option}LowResult]"
    SaveToLogFile -step $step -option $option -jobs $jobs -operation $operation
    return [set ${option}LowResult]
}

proc SaveToLogFile {args} {
    set option p
    set step ""
    set jobs 0
    set operation ""
    APSParseArguments {jobs option step jobs operation}
    global logFile parameterName logFileID pVector contractVector reflVector expanVector
    global rootname logDataArray procColumns logColumns

    if ![string length $logFile] {
        return
    }
    set runName ${rootname}_${option}-
    set procFiles [lsort [glob ${runName}??????.proc]]
    set tmpRoot /tmp/[APSTmpString]
    if [catch {eval exec sddscombine $procFiles $tmpRoot.1 -merge -overwrite} result] {
	puts stderr "sddscombine $procFiles $tmpRoot.1 -merge -overwrite"
        return -code error $result
    }
    
    if ![string length $operation] {
        switch $option {
            p {
                set operation shrinkage
                #note that it is actually shrinkage operation -- shrink toward low point
            }
            refl {
                set operation reflection
            }
            contract {
                set operation contraction
                #note that it is actually contraction operation  -- contract the simplex
            }
            expan {
                set operation expansion
            }
        }
    }
    set options ""
    set index 0
    for {set i 0} {$i<$jobs} {incr i} {
        lappend idList $i
        lappend stepList $step
        lappend operList $operation
    }
    foreach name $parameterName {
        for {set i 0} {$i<$jobs} {incr i} {
            lappend $name [lindex [set ${option}Vector($i.Value)] $index]
        }
        incr index
        lappend options "-col=$name,type=double"
        lappend options "-data=[join [set $name] ,]"
    }
    if [catch {eval exec sddsmakedataset $tmpRoot.2 \
                 -col=step,type=long -data=[join $stepList ,] \
                 -col=operation,type=string -data=[join $operList ,]\
                 -col=runID,type=long -data=[join $idList ,] $options} result] {
        return -code error $result
    }
    if [catch {exec sddsxref $tmpRoot.2 $tmpRoot.1 $tmpRoot.4 -fillIn -nowarnings \
                 -leave=[join $parameterName ,],step,operation,runID } result] {
        return -code error $result
    }
    if {$step==1} {
        exec mv $tmpRoot.4 $logFile
    } else {
        if [catch {exec sddscombine $logFile $tmpRoot.4 -merge $tmpRoot.5 } result] {
	    puts stderr "sddscombine $logFile $tmpRoot.4 -merge $tmpRoot.5"
            return -code error $result
        }
        exec mv $tmpRoot.5 $logFile
    }
    set files [glob ${tmpRoot}.*]
    eval file delete -force $files
}

proc CalculatePSum {args} {
    global pVector 
    global inputParameterList inputColumnList 
    eval global $inputParameterList 
    eval global $inputColumnList 

    set index 0
    foreach pN $parameterName {
        set pVector($index.sum) 0
        for {set point 0} {$point <=$pVector(Dimensions)} {incr point} {
            set pVector($index.sum) [expr $pVector($index.sum) + [lindex $pVector($point.Value) $index]]
        }
        set pVector($index.center) [expr $pVector($index.sum)/($pVector(Dimensions)+1)]
        incr index
    }
}

proc ContractAroundBestPoint {args} {
    global inputParameterList inputColumnList 
    eval global $inputParameterList 
    eval global $inputColumnList
    global pVector pHighID pLowID bestResult parameterList lowerLimit upperLimit rootname pExist

    set jobs 0
    APSParseArguments {jobs} 
    
    set bestValues $pVector($pLowID.Value)
    for {set job 0} {$job < $jobs} {incr job} {
        if {$job !=$pLowID} {
            set index 0
            set valueList ""
            foreach val $pVector($job.Value) lL $lowerLimit uL $upperLimit {
                set val [expr 0.5 * ($val + [lindex $bestValues $index])]
                if {$val < $lL} {set $val $lL}
                if {$val > $uL} {set val $uL}
                lappend valueList $val
                incr index
            }
            set pVector($job.Value) $valueList
            set pExist($job) 0
        }
    }
    set pExist($pLowID) 1
    DistributeJobs -option p
}

proc PrintBestResult {args} {
    global parameterName pVector pLowID pHighID rootname steps bestResult bestValues pLowResult
    global restartFile input logFileID startTime restart iterations totalIterations
    set result 1
    APSParseArguments {result}
    
    if [string length $logFileID] {close $logFileID}
    puts stderr "Optimization started at $startTime"
    switch $result {
        1 {
            puts stderr "[exec date +%H:%M:%S] Optimization done after $steps steps."
            puts stderr "The best result is: $bestResult"
        }
        2 {
            puts stderr "[exec date +%H:%M:%S] the number of total iterations exceeds $totalIterations, exit."
            puts stderr "The best result obtained so far is: $bestResult"
        }
    }
    puts stderr "The parameter values are:"
    foreach pN $parameterName val $bestValues {
        puts stderr "$pN: $val"
    }
    set files [glob -nocomplain ${rootname}_*]
    catch {eval file delete -force file $files }
    set files [glob -nocomplain ${rootname}-*]
    catch {eval file delete -force file $files}
    set tmpfile /tmp/[APSTmpString]
    if [string length $restartFile] {
        if [catch {exec sddsmakedataset $tmpfile -col=parameterName,type=string \
                     -data=[join $parameterName ,] \
                     -col=initialValue1,type=double -data=[join $bestValues ,] \
                     -par=penaltyValue,type=double -data=$bestResult } result] {
            return -code error $result
        }
        if [catch {exec sddsxref $input $tmpfile -pipe=out -match=parameterName \
                     -take=initialValue1 \
                     | sddsconvert -pipe=in -delete=col,initialValue \
                     -rename=col,initialValue1=initialValue $restartFile } result] {
            return -code error $result
        }
        file delete -force $tmpfile
    }
}

proc RunParallelOptions {args} {
    set lamda 0
    set option p
    APSParseArguments {lamda option}

    global pVector reflVector expanVector contractVector pLowID lowerLimit upperLimit reflThroughLowPoint
    global pExist reflExist expanExist contractExist sum errorFactor pLowID rootname
    
    # CalculatePSum
    set fac1 [expr (1.0 - $lamda) / $pVector(Dimensions)]
    set fac2 [expr $fac1 - $lamda]
    if {$reflThroughLowPoint} {
        set valueList1 $pVector($pLowID.Value)
    } else {
        for {set vetex 0} {$vetex < $pVector(Dimensions)} {incr vetex} {
            lappend valueList1 $pVector($vetex.center)
        }
    }
    set vetexes [expr $pVector(Dimensions) + 1]
    if $sum {
        set factor $vetexes
    } else {
        set factor 1
    }
    for {set vetex 0} {$vetex<$vetexes} {incr vetex} {
	set ${option}Vector($vetex.limit) 0
        set valueList ""
        set index 0
        foreach val $pVector($vetex.Value) lL $lowerLimit uL $upperLimit {
	    set val1 [lindex $valueList1 $index]
            set value [expr $val1*$fac1*$factor - $val * $fac2]
            if {$value<=$lL} {
                set value $lL
		set ${option}Vector($vetex.limit) 1
            }
            if {$value>=$uL} {
                set value $uL
		set ${option}Vector($vetex.limit) 1
            }
            incr index
            lappend valueList $value
        }
        set ${option}Vector($vetex.Value) $valueList
        set ${option}Exist($vetex) 0
    }
    if {$reflThroughLowPoint && $sum} {
        if {$option!="p"} {
            set root [format ${rootname}_p-%06ld $pLowID]
            set files [glob ${root}.*]
            set root1 [format ${rootname}_${option}-%06ld $pLowID]
            foreach file $files {
                set ext [file extension $file]
                file copy -force $file ${root1}$ext
            }
            set ${option}Exist($pLowID) 1
        }
    }
    DistributeJobs -option $option
}

proc ReplaceHighPoint {args} {
    global pVector reflVector expanVector contractVector bestResult bestValues  pLowResult pHighResult
    global pLowID pHighID reflLowID reflHighID expanLowID expanHighID contractLowID contractHighID rootname
    set repl ""
    APSParseArguments {repl LowResult}
    
    set pVector($pHighID.fValue) [set ${repl}Vector([set ${repl}LowID].fValue)]
    set pVector($pHighID.Value) [set ${repl}Vector([set ${repl}LowID].Value)]
    set max $pVector($pHighID.fValue)
    for {set i 0} {$i<=$pVector(Dimensions)} {incr i} {
        if {$pVector($i.fValue)>$max} {
            set pHighID $i
            set max $pVector($i.fValue)
        }
    }
    set pHighResult $max
}

proc ReplaceSimplex {args} {
    global pVector reflVector expanVector contractVector bestResult bestValues  pLowResult pHighResult
    global pLowID pHighID reflLowID reflHighID expanLowID expanHighID contractLowID contractHighID rootname
    set repl ""
    APSParseArguments {repl LowResult}
    
    set dim $pVector(Dimensions)
    for {set point 0} {$point <=$dim} {incr point} {
        set pVector($point.Value) [set ${repl}Vector($point.Value)]
        set pVector($point.fValue) [set ${repl}Vector($point.fValue)]
    }
    
    set pLowID [set ${repl}LowID]
    set pHighID [set ${repl}HighID]
    set pLowResult $pVector($pLowID.fValue)
    set pHighResult $pVector($pHighID.fValue)
    set runname1 [format ${rootname}_$repl-%06ld $pLowID]
    set runname2 [format ${rootname}_p-%06ld $pLowID]
    set fileList [glob ${runname1}.*]
    foreach file $fileList {
        set ext [file extension $file]
        set replFile ${runname2}$ext
        file copy -force $file $replFile
    }
    puts stderr "lowResult after replace=$pLowResult, lowID=$pLowID"   
}

#scan one dimension (defined by index) by +- delta and +- 2* delta
proc OneDScan {args} {
    global inputParameterList inputColumnList processors
    eval global $inputParameterList 
    eval global $inputColumnList
    global rootname errorFactor pVector

    set rootname1 ${rootname}_
    set submitted 0
    set pVector(Dimensions) [llength $initialValue]
    set pVector(parameterName) $parameterName
    set pVector(0.Value) $initialValue
    set procOpt0 ""
    foreach pN $parameterName val $initialValue {
        lappend procOpt0 -define=column,$pN,$val,type=double
    }
    set replOpt [concat ${rootname1}center-0 $initialValue]
    #submit the job for initial point
    set runName ${rootname1}center-0
    set fileList [glob -nocomplain ${runName}.*]
    eval file delete -force file $fileList
    lappend jobList "SubmitOneJob -runName $runName -procOpt \"$procOpt0\" \
      -replOpt \"$replOpt\" -valueList \"$initialValue\" -id 0"
    incr submitted
    set ${runName}.Value $initialValue
    set $runName 0
    set doneFiles $runName.done
    set procFiles $runName.proc
    set idList 0
    for {set index 0} {$index<$pVector(Dimensions)} {incr index} {
        # scan one dimension with four points 
        set runname ${rootname1}${index}
        set i 0
        set lL [lindex $lowerLimit $index]
        set uL [lindex $upperLimit $index]
        set eL [lindex $errorLevel $index]
        set iV [lindex $initialValue $index]
        set pN [lindex $parameterName $index]
        foreach factor {-2 -1 1 2} id {0 1 2 3} {
            set val [expr $iV + $factor * $eL * $errorFactor]
            if {$val<$lL} {set val $lL}
            if {$val>$uL} {set val $uL}
            set runName ${runname}-$id
            set valueList [lreplace $initialValue $index $index $val]
            set ${runName}.Value $valueList
            set $runName 0
            set fileList [glob -nocomplain ${runName}.*]
            eval file delete -force file $fileList
            set replOpt [concat $runName $valueList]
            set procOpt [lreplace $procOpt0 $index $index -define=column,$pN,$val,type=double]
            lappend jobList "SubmitOneJob -runName $runName -procOpt $procOpt -replOpt \"$replOpt\" \
              -id $id -valueList \"$valueList\""
            lappend doneFiles $runName.done
            lappend procFiles $runName.proc
	    lappend idList $id
        }
    }
    foreach done $doneFiles {
        set $done 0
        set $done.checked 0
    }
    if {$processors>0} {
        set submitted 0
        set jobsrunning 0
        set jobs [llength $doneFiles]
        while {$submitted<$jobs} {
            foreach job $jobList done $doneFiles {
                if {$jobsrunning<$processors && ![set $done]} {
                    incr submitted
                    incr jobsrunning
                    set $done 1
                    puts stderr "job [file root $done] submitted."
                    eval $job
                }
                if {[file exist $done] && ![set $done.checked]} {
                    set jobsrunning [expr $jobsrunning -1]
                    set $done.checked 1
                }
            }
            after 1000
        } 
    } else {
        foreach job $jobList {
            eval $job
        }
    }
    puts stderr "Jobs for one-d scan submitted."
    while {1} {
        set notdone 0
        foreach file $doneFiles ID $idList {
            set runName [file rootname $file]
            if [file exists $file] {
                if ![file exists ${runName}.proc] {
                    if ![string length $postProcessingScript] {
                        puts stderr "postProcessingScript is not provided for creating *.proc file."
                        exit
                    }
                    if [catch {exec $postProcessingScript -tagList "$parameterName" \
                                 -valueList "[set ${runName}.Value]" \
                                 -runID $ID -rootname $runName } result] {
                        puts stderr "$result"
                        exit
                    }
                }
                if {![set $runName]} {
                    if [catch {CheckFileExistent -filename $runName.proc} result] {
                        puts stderr $result
                        exit
                    }
                    if [catch {exec sdds2stream $runName.proc -col=penaltyValue} results] {
                        puts stderr "unable to read $root.proc: $results"
                        exit
                    }
                    if ![string length $results] {
                        puts stderr "unable to get the penalty value of $runName.proc, [exec sddscheck $runName.proc]"
                        exit
                    }
                    set ${runName}.fValue  [lindex $results 0]
                    if [catch {expr [set ${runName}.fValue] / 2} result] {
                        set ${runName}.fValue 1.0e200
                    }
                    set $runName 1
                    puts stderr "$runName done, f=[set ${runName}.fValue], parameter values are:"
                    puts stderr "               [set ${runName}.Value]"
                }
            } else {
                set notdone 1
            }
        }
        if {!$notdone} {
            break
        }
        after 1000
    }
    global pLowID pHighID pLowResult pHighResult
    set pVector(0.limit) 0
    set pVector(0.fValue) [set ${rootname1}center-0.fValue]
    set pVector(0.Value) $initialValue
    set lowID 0
    set highID 0
    set lowResult $pVector(0.fValue)
    set highResult $pVector(0.fValue)
    set minRoot ${rootname1}center-0
    set fileList [glob ${minRoot}.*]
    foreach file $fileList {
        set ext [file extension $file]
        file copy -force $file ${rootname}_p-000000$ext
    }
    for {set index 0} {$index < $pVector(Dimensions)} {incr index} {
        set runname ${rootname1}${index}
        set min [set ${runname}-0.fValue]
        set minIndex 0
        for {set i 1} {$i<=3} {incr i} {
            set fValue [set ${runname}-$i.fValue]
            if {$fValue < $min} {
                set minIndex $i
                set min $fValue
            }
        }    
        set point [expr $index +1]
	set pVector($point.limit) 0
        set pVector($point.Value) [set ${runname}-${minIndex}.Value]
        set pVector($point.fValue) $min
        set fileroot [format ${rootname}_p-%06ld $point]
        set minRoot ${rootname1}${index}-$minIndex
        set fileList [glob ${minRoot}.*]
        foreach file $fileList {
            set ext [file extension $file]
            file copy -force $file ${fileroot}$ext
        }
        if {$min<$lowResult} {
            set lowID $point
            set lowResult $min
        }
        if {$min>$highResult} {
            set highID $point
            set highResult $min
        }
    }
    set pHighID $highID
    set pLowID $lowID
    set pHighResult $highResult
    set pLowResult $lowResult
  
    puts stderr "one-d scan done, lowID=$pLowID, highID=$pHighID,lowResult=$lowResult, highResult=$highResult"
    set fileList [glob ${rootname1}*-?.*]
    eval file delete -force $fileList
    SaveToLogFile -step 1 -option p -jobs [expr $pVector(Dimensions) + 1]\
      -operation oneDScan 
}


proc CheckFileExistent {args} {
    set filename ""
    set timeout 500
    APSParseArguments {filename timeout}
    
    set timeout [expr [clock seconds] + $timeout]
    set rows 0
    while {[clock seconds]<$timeout} {
        if {[exec sddscheck $filename]=="ok"} {
            if [catch {exec sdds2stream $filename -rows=bar} rows] {
                set rows 0
            }
            if {$rows} {
                break
            }
        }
        after 1000
    }
    return $rows
}

LoadAndCheckInputFile 
set files [glob -nocomplain ${rootname}*.run]
if [llength $files] {
    puts stderr "There are jobs for $rootname running, can not start new run!"
    exit
}
set jobsStarted 0
set startTime "[exec date +%H:%M:%S]"
StartOptimize


