whatIf? Wiki

This page is read only. You can view the source, but not change it. Ask your administrator if you think this is wrong.
====== Minimizing Root Mean Squared Difference with 5 point bracketing ======
A typical calibration problem is to get the model results to track on something while tweaking a variable.  The approach described here attempts to minimize the difference over time between the model variable and some given tracking variable.

Here is the code for the .t file for the feed back view.  The example implementation was to change the unitTA variable such that the difference between the model calculated space conditioning energy use and the energy use left over from the Stats Can total energy use minus the non-space conditioning energy use is minimized.

Instructions on how to customize the file for your situation is right in the top of the file.

<file>
! ** This is a generic bracketing feedback view that does a divide and conquer approach **
! ** REQUIREMENTS:**
! ** 		It requires a min and max starting value which bound the problem **
! ** ALGORITHM:**
! ** 		1. It divides the solution space into 5 evaluation positions, runs the model at each place and **
! ** 		chooses the position that minimizes the root mean squared difference between the model estimate **
! ** 		and the target **
! ** 		2. It then reduces the step size in half, chooses 5 more points centered on the best result from **
! **    the last evaluation and then runs the model (as needed) to get the result at these 5 new positions**
! ** 		3. Repeats the process of finding the min of the 5 points and sets up 5 more points until **
! ** 		the step size as a ratio of the best guess is less than the feedback view tolerance **
! ** HOW TO CUSTOMIZE FOR YOUR PROBLEM **
! **		1. Search through the file for BEGIN USER MUST ADJUST HERE and adjust the code for your variables

creview	unitTA[cr][]

view EUSpCTotal[cr,rfts,th][]
view SCEULeftForSC[cr,rfts,th][]

boolean $continue
if $numScns == 0
	say ("This view requires a starting scenario.")
	listen ("                              enter to continue >> ",$continue)
	quit (-1)
endif

! ** BEGIN USER MUST ADJUST HERE **
local modelEst[upperDim,otherDims] = EUSpCTotal[cr,rfts,th][1]
local target[upperDim,otherDims] = SCEULeftForSC[cr,rfts,th][1]

! sum up the "other dims"
local diff[upperDim] = sum (target[upperDim,otherDims] - modelEst[upperDim,otherDims]; dim1=resFuelTypeSim, dim2=time)

local outputVar[upperDim,auxDims] = copyshape (unitTA[cr][2])
! ** END USER MUST ADJUST HERE **

!integer	$fbViewIteration	readOnly
!integer	$fbViewMaxIters		readOnly
!real	$fbViewTolerance	readOnly

local tolerance[] = create (; data=$fbViewTolerance)

! Declare keep variables
keep xm[upperDim,xm]
keep currentGuess[upperDim,auxDims,xm]
keep currentMeshStep[upperDim,auxDims]
keep rmsq[upperDim,xm]
keep rmsqLastEval[upperDim,xm]
keep x_min[upperDim,auxDims]
keep x_max[upperDim,auxDims]
keep doEvaluation[]
keep lastGuessNumber[]
keep numEvalsCompleted[]

say ("")
say ("  processing for iteration ",$fbViewIteration)

if ($fbViewIteration == 1)
	! Setup shapes and initialize keep variables
	!set up the mesh points
	xm[upperDim,xm] = cumsum (create (; dimFrom1=diff[upperDim], dim=SEQ;xm:1:5:1;xmStep, data=1); multiplyEntity=off)

	!get initial guess from the backeting values which you must know
	currentGuess[upperDim,auxDims,xm] = copyshape (outputVar[upperDim,auxDims] * xm[upperDim,xm])
	currentMeshStep[upperDim,auxDims] = copyshape (sum (currentGuess[upperDim,auxDims,xm]))
	local rm[upperDim,xm] = changeshape (copyshape (diff[upperDim]); addDim2From2=xm[upperDim,xm])
	rmsq[upperDim,xm] = copyshape (rm[upperDim,xm] * rm[upperDim,xm])

! ** BEGIN USER MUST ADJUST HERE **
	x_min[upperDim,auxDims] = copyshape (unitTA[cr][1])
	x_min[upperDim,auxDims] = 0.000001
	x_max[upperDim,auxDims] = copyshape (unitTA[cr][1])
	x_max[upperDim,auxDims] = 0.000010
! ** END USER MUST ADJUST HERE **
	
	say ("    get initial 5 guesses and set independent variable to guess 1")
	local currentGuessx[upperDim,auxDims,xm] = copyshape (currentGuess[upperDim,auxDims,xm])
	currentGuessx[upperDim,auxDims,xm] = insert (x_min[upperDim,auxDims]; xm:1)
	currentGuessx[upperDim,auxDims,xm] = insert (x_max[upperDim,auxDims]; xm:5)
	currentGuess[upperDim,auxDims,xm] = linint (currentGuessx[upperDim,auxDims,xm])

	!get initial mesh step
	currentMeshStep[upperDim,auxDims] = x_max[upperDim,auxDims] - x_min[upperDim,auxDims] / 4

	outputVar[upperDim,auxDims] = extract (currentGuess[upperDim,auxDims,xm]; xm:1, shrink=on)
	lastGuessNumber[] = 1
	doEvaluation[] = 0
	numEvalsCompleted[] = 0
	skipto lastStepQuit1
endif

if ($fbViewIteration > 1)
	integer $lastGuessNumber
	getobjinfo (lastGuessNumber[], $lastGuessNumber; info=dataElement)
	say ("    get function for guess position ", $lastGuessNumber)
	rmsq[upperDim,xm] = insert (diff[upperDim] * diff[upperDim]; xm:$lastGuessNumber)
endif


if ($fbViewIteration > 1) and ($fbViewIteration <= 5)
	say ("    set independent variable to guess ",$fbViewIteration)
	outputVar[upperDim,auxDims] = extract (currentGuess[upperDim,auxDims,xm]; xm:$fbViewIteration, shrink=on)
	lastGuessNumber[] = $fbViewIteration

	if ($fbViewIteration==5)
		doEvaluation[] = 1
	endif	
	skipto lastStepQuit1
endif

integer $doEvaluation
getobjinfo (doEvaluation[], $doEvaluation; info=dataElement)
if 	$doEvaluation == 1
	say ("  processing for iteration ",$fbViewIteration)

	say ("    **do an evaluation**")
	!do evaluation here and get next guesses")
	!get best guess
	local minValue[upperDim], coordAtMin[upperDim,V], coordAtMin1[upperDim]
	minValue[upperDim], coordAtMin[upperDim,V] = coordmin (rmsq[upperDim,xm]; numDims=1)
	coordAtMin1[upperDim] = extract (coordAtMin[upperDim,V]; coord:+0, shrink=on)

	local isThisCoordCurrentBestGuess[upperDim,xm], currentBestGuess[upperDim,auxDims], rmsq_atCurrentGuess[upperDim]
	isThisCoordCurrentBestGuess[upperDim,xm] = booleq (xm[upperDim,xm], coordAtMin1[upperDim])
	currentBestGuess[upperDim,auxDims] = sum (isThisCoordCurrentBestGuess[upperDim,xm] * currentGuess[upperDim,auxDims,xm]; dim=xm)
	rmsq_atCurrentGuess[upperDim] = sum (isThisCoordCurrentBestGuess[upperDim,xm] * rmsq[upperDim,xm]; dim=xm)

	local numEvalsCompletedx[] = numEvalsCompleted[]
	numEvalsCompleted[] = numEvalsCompletedx[] + 1
	integer $numEvalsCompleted
	getobjinfo (numEvalsCompleted[], $numEvalsCompleted; info=dataElement)
	say ("      number of evaluations completed ",$fbViewIteration )
	if ($numEvalsCompleted > 1)
!		Can look at last evaluations results compared to this one then store this one for future comparisons	
!		table (rmsqLastEval[upperDim,xm], rmsq[upperDim,xm])
	endif
	rmsqLastEval[upperDim,xm] = rmsq[upperDim,xm]

	!check for stoping
	local numberNotConverged[] = sum (boolgt (abs (2 * currentMeshStep[upperDim,auxDims] / currentBestGuess[upperDim,auxDims]), tolerance[]))
	integer $numberNotConverged
	getobjinfo (numberNotConverged[], $numberNotConverged; info=dataElement, default=1)

	if ($numberNotConverged == 0) or ($fbViewIteration == $fbViewMaxIters)
		if ($numberNotConverged == 0)
			say ("      stopping by convergence at iteration ",$fbViewIteration )
		endif
		if ($fbViewIteration == $fbViewMaxIters)
			say ("      stopping by reaching maximum iterations ",$fbViewIteration )
			say ("        number not converged is ",$numberNotConverged)
		endif
		outputVar[upperDim,auxDims]	= currentBestGuess[upperDim,auxDims]
		skipto lastStepQuit0
	else
		say ("    set new guesses for points 1 3 and 5")
		local currentGuessx[upperDim,auxDims,xm] = copyshape (currentGuess[upperDim,auxDims,xm])
		currentGuessx[upperDim,auxDims,xm] = insert (currentBestGuess[upperDim,auxDims]; xm:3)
		currentGuessx[upperDim,auxDims,xm] = insert (max (currentBestGuess[upperDim,auxDims] - currentMeshStep[upperDim,auxDims], x_min[upperDim,auxDims]); xm:1)
		currentGuessx[upperDim,auxDims,xm] = insert (min (currentBestGuess[upperDim,auxDims] + currentMeshStep[upperDim,auxDims], x_max[upperDim,auxDims]); xm:5)
		currentGuess[upperDim,auxDims,xm] = linint (currentGuessx[upperDim,auxDims,xm])
		local currentMeshStepx[upperDim,auxDims] = currentMeshStep[upperDim,auxDims] / 2
		currentMeshStep[upperDim,auxDims] = currentMeshStepx[upperDim,auxDims]

		!get new rmsq by putting the best guess in the middle, the guess one below as the new min, and 
		!the one above as the new max
		local rmsqx[upperDim,xm] = copyshape (rmsq[upperDim,xm])
		rmsqx[upperDim,xm] = insert (rmsq_atCurrentGuess[upperDim]; xm:3)
		
		local isThisCoord[upperDim,xm], guess[upperDim,auxDims], rmsq_atGuess[upperDim]
		isThisCoord[upperDim,xm] = booleq (max (coordAtMin1[upperDim] - 1, 1), xm[upperDim,xm])
		guess[upperDim,auxDims] = sum (isThisCoord[upperDim,xm] * currentGuess[upperDim,auxDims,xm]; dim=xm)
		rmsq_atGuess[upperDim] = sum (isThisCoord[upperDim,xm] * rmsq[upperDim,xm]; dim=xm)
		rmsqx[upperDim,xm] = insert (rmsq_atGuess[upperDim]; xm:1)
	
		isThisCoord[upperDim,xm] = booleq (min (coordAtMin1[upperDim] + 1, 5), xm[upperDim,xm])
		guess[upperDim,auxDims] = sum (isThisCoord[upperDim,xm] * currentGuess[upperDim,auxDims,xm]; dim=xm)
		rmsq_atGuess[upperDim] = sum (isThisCoord[upperDim,xm] * rmsq[upperDim,xm]; dim=xm)
		rmsqx[upperDim,xm] = insert (rmsq_atGuess[upperDim]; xm:5)
		
		rmsq[upperDim,xm] = rmsqx[upperDim,xm]

		!this leaves only xm=2,4 to be calculated before getting next minimum
		say ("    set independent variable to guess 2")
		outputVar[upperDim,auxDims] = extract (currentGuess[upperDim,auxDims,xm]; xm:2, shrink=on)		
		lastGuessNumber[] = 2
		doEvaluation[] = 0
		skipto lastStepQuit1
	endif
else
	say ("    set independent variable to guess 4")
	outputVar[upperDim,auxDims] = extract (currentGuess[upperDim,auxDims,xm]; xm:4, shrink=on)		
	lastGuessNumber[] = 4
	doEvaluation[] = 1
	skipto lastStepQuit1
endif

! ** BEGIN USER MUST ADJUST HERE **
lastStepQuit0:
	unitTA[cr][2] = outputVar[upperDim,auxDims]
	quit (0)
	
lastStepQuit1:
	unitTA[cr][2] = outputVar[upperDim,auxDims]
	quit (1)
! ** END USER MUST ADJUST HERE **
</file>
whatIf? Wiki

User Tools

Site Tools

Page Tools
