Differences

This shows you the differences between two versions of the page.

--- howtos:feedback:tracking_method_-_minimizing_root_mean_squared_difference_with_5_point_bracketing [2010/12/14 19:28]
shona.weldon
+++ howtos:feedback:tracking_method_-_minimizing_root_mean_squared_difference_with_5_point_bracketing [2013/02/10 17:10] (current)
marcus.williams
@@ Line 1: / Line 1: @@
-====== Minimizing Root Mean Squared Difference with 5 point bracketing ======
+====== Minimizing Root Mean Squared Difference with 5 point bracketing (bi-section variant) ======
 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.
@@ Line 5: / Line 5: @@
 Instructions on how to customize the file for your situation is right in the top of the file.
+===== Example code - simple =====
 <file>
@@ Line 65: / Line 67: @@
 if ($fbViewIteration == 1)
+! ** NOTE: If there is more than one upper dimension the following block of code must be adjusted accordingly ***
 	! Setup shapes and initialize keep variables
 	!set up the mesh points
@@ Line 76: / Line 81: @@
 ! ** BEGIN USER MUST ADJUST HERE **
+	! The [1] object indices assume calibration; for simulation views increment accordingly
 	x_min[upperDim,auxDims] = copyshape (unitTA[cr][1])
 	x_min[upperDim,auxDims] = 0.000001
@@ Line 103: / Line 109: @@
 	say ("    get function for guess position ", $lastGuessNumber)
 	rmsq[upperDim,xm] = insert (diff[upperDim] * diff[upperDim]; xm:$lastGuessNumber)
+!	display (rmsq[upperDim,xm]) !helpful if debugging
 endif
@@ Line 133: / Line 140: @@
 	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)
+!display (coordAtMin1[upperDim], minValue[upperDim]) !helpful when debugging
 	local numEvalsCompletedx[] = numEvalsCompleted[]
@@ Line 138: / Line 146: @@
 	integer $numEvalsCompleted
 	getobjinfo (numEvalsCompleted[], $numEvalsCompleted; info=dataElement)
-	say ("      number of evaluations completed ",$fbViewIteration )
+	say ("      number of evaluations completed ",$numEvalsCompleted)
 	if ($numEvalsCompleted > 1)
 !		Can look at last evaluations results compared to this one then store this one for future comparisons
@@ Line 146: / Line 154: @@
 	!check for stoping
-	local numberNotConverged[] = sum (boolgt (abs (2 * currentMeshStep[upperDim,auxDims] / currentBestGuess[upperDim,auxDims]), tolerance[]))
+	local numberNotConverged[] = sum (boolgt (abs (2 * currentMeshStep[upperDim,auxDims] / currentBestGuess[upperDim,auxDims]), \
+             tolerance[]); dim=all)
 	integer $numberNotConverged
 	getobjinfo (numberNotConverged[], $numberNotConverged; info=dataElement, default=1)
@@ Line 230: / Line 239: @@
 	quit (1)
 ! ** END USER MUST ADJUST HERE **
+</file>
+===== Example code for if you want several modes your view =====
+<file>
+creview	themCondC[uss,usdt,rsta][]
+creview	useExoThemCond[][]
+view themSurfAreaCond[uss,usdt,rsta][]
+view 	detEUSpC[uss,usdt,usrft->usreud,th][]
+view	resEnrEndUseSEDSandRECSEst[uss,usrft->usreu,t19912008][]
+boolean $continue
+if $numScns == 0
+	say ("This view requires a starting scenario.")
+	listen ("                              enter to continue >> ",$continue)
+	quit (-1)
+endif
+keep mode[]
+keep firstIterForCurMode[]
+keep fbViewMaxItersCurMode[]
+if ($fbViewIteration == 1)
+	mode[] = 1
+	firstIterForCurMode[] = 1
+	fbViewMaxItersCurMode[] = 40
+endif
+integer $mode
+getobjinfo (mode[], $mode; info=dataElement)
+integer $firstIterForCurMode
+getobjinfo (firstIterForCurMode[], $firstIterForCurMode; info=dataElement)
+integer $firstIterForCurModeP4
+$firstIterForCurModeP4 = $firstIterForCurMode + 4
+integer $fbViewMaxItersCurMode
+getobjinfo (fbViewMaxItersCurMode[], $fbViewMaxItersCurMode; info=dataElement)
+if (mode == 1)
+	! ** BEGIN USER MUST ADJUST HERE **
+	local modelEst[upperDim,otherDims] = sum (extract (mapcat (detEUSpC[uss,usdt,usrft->usreud,th][1]; \
+		dim=usrftCOusreud); \
+		USResEndUseDet:spaceHeating, time:$firstCompYr..$lastCompYr, shrink=on); \
+		dim=USDwellType)
+	local target[upperDim,otherDims] = extract (mapcat (resEnrEndUseSEDSandRECSEst[uss,usrft->usreu,t19912008][1]; \
+		dim=usrftCOusreu); \
+		USResEndUse:spaceHeating, shrink=on)
+	! sum up the "other dims"
+	local diff[upperDim] = sum (target[upperDim,otherDims] - modelEst[upperDim,otherDims]; dim1=USResFuelType, dim2=time)
+	local outputVar[upperDim,auxDims] = copyshape (themCondC[uss,usdt,rsta][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 == $firstIterForCurMode)
+		! 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 **
+		local HOURPerDay[] = create (; entity=HOUR/day, data=24)
+		x_min[upperDim,auxDims] = themSurfAreaCond[upperDim,auxDims][1] / HOURPerDay[] / 2
+		x_max[upperDim,auxDims] = themSurfAreaCond[upperDim,auxDims][1] / HOURPerDay[] * 5
+	! ** 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 > $firstIterForCurMode)
+		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 > $firstIterForCurMode) and ($fbViewIteration <= $firstIterForCurModeP4)
+		integer $iter = $fbViewIteration - $firstIterForCurMode + 1
+		say ("    set independent variable to guess ",$iter)
+		outputVar[upperDim,auxDims] = extract (currentGuess[upperDim,auxDims,xm]; xm:$iter, shrink=on)
+		lastGuessNumber[] = $iter
+	if ($fbViewIteration == $firstIterForCurModeP4)
+			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 == $fbViewMaxItersCurMode)
+			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 ("    Num not converged: ", $numberNotConverged)
+			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:
+		themCondC[uss,usdt,rsta][2] = outputVar[upperDim,auxDims]
+		useExoThemCond[][2] = 1
+!NOTE IF THERE IS ANOTHER MODE YOU WOULD WANT TO SWITHCH THE MODE HERE AND RESET THE FOLLOWING:
+!	mode[] = 2
+!       getobjinfo (mode[], $mode; info=dataElement)
+!	firstIterForCurMode[] = $fbViewIteration + 1
+!       getobjinfo (firstIterForCurMode[], $firstIterForCurMode; info=dataElement)
+!       $firstIterForCurModeP4 = $firstIterForCurMode + 4
+!	fbViewMaxItersCurMode[] = $fbViewIteration + 10
+!       getobjinfo (fbViewMaxItersCurMode[], $fbViewMaxItersCurMode; info=dataElement)
+	lastStepQuit1:
+		themCondC[uss,usdt,rsta][2] = outputVar[upperDim,auxDims]
+		useExoThemCond[][2] = 1
+		quit (1)
+	! ** END USER MUST ADJUST HERE **
+endif !$mode == 1
 </file>

whatIf? Wiki

User Tools

Site Tools

Differences

Page Tools