From fa508cbcc6bab30c50e88098e38364d1f43cc499 Mon Sep 17 00:00:00 2001 From: Craig Warren Date: Tue, 15 Dec 2015 12:08:18 +0000 Subject: [PATCH] Added Taguchi optimisation module. --- user_libs/optimisations/OA_18_7_3_2.npy | Bin 0 -> 206 bytes user_libs/optimisations/OA_9_4_3_2.npy | Bin 0 -> 116 bytes user_libs/optimisations/taguchi.py | 212 ++++++++++++++++++++++++ 3 files changed, 212 insertions(+) create mode 100644 user_libs/optimisations/OA_18_7_3_2.npy create mode 100644 user_libs/optimisations/OA_9_4_3_2.npy create mode 100644 user_libs/optimisations/taguchi.py diff --git a/user_libs/optimisations/OA_18_7_3_2.npy b/user_libs/optimisations/OA_18_7_3_2.npy new file mode 100644 index 0000000000000000000000000000000000000000..8c40ded14881d03d7e62baa96ce7153b94c06a5b GIT binary patch literal 206 zcmbR27wQ`j$;jZwP_3SlTAW;@Zl$1JQ);NLqoAIaUsO_*m=~X4l#&V(4=E~51qv5u zBo?Fsxf+HRItu2RItsN4aKOdD$b=OEg@Fl2*gZ`3;@BP6^Q@< literal 0 HcmV?d00001 diff --git a/user_libs/optimisations/taguchi.py b/user_libs/optimisations/taguchi.py new file mode 100644 index 00000000..8c64b288 --- /dev/null +++ b/user_libs/optimisations/taguchi.py @@ -0,0 +1,212 @@ +# Copyright (C) 2015, Craig Warren +# +# This module is licensed under the Creative Commons Attribution-ShareAlike 4.0 International License. +# To view a copy of this license, visit http://creativecommons.org/licenses/by-sa/4.0/. +# +# Please use the attribution at http://dx.doi.org/10.1190/1.3548506 + +import os +import numpy as np +from collections import OrderedDict + +import h5py + +from gprMax.constants import floattype +from gprMax.exceptions import CmdInputError + +moduledirectory = os.path.dirname(os.path.abspath(__file__)) + + +def select_OA(optparams): + """Load an orthogonal array (OA) from a numpy file. Configure and return OA and properties of OA. + + Args: + optparams (dict): Dictionary containing name of parameters to optimise and their initial ranges + + Returns: + OA (array): Orthogonal array + N (int): Number of experiments in OA + k (int): Number of parameters to optimise in OA + s (int): Number of levels in OA + t (int): Strength of OA + """ + + # Load the appropriate OA + if len(optparams) <= 4: + OA = np.load(os.path.join(moduledirectory, 'OA_9_4_3_2.npy')) + elif len(optparams) <= 7: + OA = np.load(os.path.join(moduledirectory, 'OA_18_7_3_2.npy')) + else: + raise CmdInputError('Too many parameters to optimise for the available orthogonal arrays (OA). Please find and load a bigger, suitable OA.') + + # Cut down OA columns to number of parameters to optimise + OA = OA[:, 0:len(optparams)] + + # Number of experiments + N = OA.shape[0] + + # Number of parameters to optimise + k = OA.shape[1] + + # Number of levels + s = 3 + + # Strength + t = 2 + + return OA, N, k, s + + +def calculate_ranges_experiments(optparams, optparamsinit, levels, levelsopt, levelsdiff, OA, N, k, s, i): + """Calculate values for parameters to optimise for a set of experiments. + + Args: + optparams (dict): Ordered dictionary containing name of parameters to optimise and their values + optparamsinit (list): Initial ranges for parameters to optimise + levels (array): Lower, central, and upper values for each parameter + levelsopt (array): Optimal level for each parameter from previous iteration + levelsdiff (array): Difference used to set values in levels array + OA (array): Orthogonal array + N (int): Number of experiments in OA + k (int): Number of parameters to optimise in OA + s (int): Number of levels in OA + i (int): Iteration number + + Returns: + optparams (dict): Ordered dictionary containing name of parameters to optimise and their values + levels (array): Lower, central, and upper values for each parameter + levelsdiff (array): Difference used to set values in levels array + """ + + # Reducing function used for calculating levels + RR = np.exp(-(i/18)**2) + + # Calculate levels for each parameter + for p in range(0, k): + # Central levels - for first iteration set to midpoint of initial range and don't use RR + if i == 0: + levels[1, p] = ((optparamsinit[p][1][1] - optparamsinit[p][1][0]) / 2) + optparamsinit[p][1][0] + levelsdiff[p] = (optparamsinit[p][1][1] - optparamsinit[p][1][0]) / (s + 1) + # Central levels - set to optimum from previous iteration + else: + levels[1, p] = levels[levelsopt[p], p] + levelsdiff[p] = RR * levelsdiff[p] + + # Lower levels set using central level and level differences values; and check they are not outwith initial ranges + if levels[1, p] - levelsdiff[p] < optparamsinit[p][1][0]: + levels[0, p] = optparamslist[p][1][0] + else: + levels[0, p] = levels[1, p] - levelsdiff[p] + + # Upper levels set using central level and level differences values; and check they are not outwith initial ranges + if levels[1, p] + levelsdiff[p] > optparamsinit[p][1][1]: + levels[2, p] = optparamslist[p][1][1] + else: + levels[2, p] = levels[1, p] + levelsdiff[p] + + # Update dictionary of parameters to optimise with lists of new values; clear dictionary first + optparams = OrderedDict((key, list()) for key in optparams) + p = 0 + for key, value in optparams.items(): + for exp in range(0, N): + if OA[exp, p] == 0: + optparams[key].append(levels[0, p]) + elif OA[exp, p] == 1: + optparams[key].append(levels[1, p]) + elif OA[exp, p] == 2: + optparams[key].append(levels[2, p]) + p += 1 + + return optparams, levels, levelsdiff + + + +def calculate_optimal_levels(optparams, levels, levelsopt, fitness, OA, N, k): + """Calculate optimal levels from results of fitness metric by building a response table. + + Args: + optparams (dict): Ordered dictionary containing name of parameters to optimise and their values + levels (array): Lower, central, and upper values for each parameter + levelsopt (array): Optimal level for each parameter from previous iteration + fitness (array): Values from results of fitness metric + OA (array): Orthogonal array + N (int): Number of experiments in OA + k (int): Number of parameters to optimise in OA + + Returns: + optparams (dict): Ordered dictionary containing name of parameters to optimise and their values + levelsopt (array): Optimal level for each parameter from previous iteration + """ + + # Build a table of responses based on the results of the fitness metric + for p in range(0, k): + responses = np.zeros(3, dtype=floattype) + + cnt1 = 0 + cnt2 = 0 + cnt3 = 0 + + for exp in range(1, N): + if OA[exp, p] == 0: + responses[0] += fitness[exp] + cnt1 += 1 + elif OA[exp, p] == 1: + responses[1] += fitness[exp] + cnt2 += 1 + elif OA[exp, p] == 2: + responses[2] += fitness[exp] + cnt3 += 1 + + responses[0] /= cnt1 + responses[1] /= cnt2 + responses[2] /= cnt3 + + # Calculate optimal level from table of responses + tmp = np.where(responses == np.amax(responses))[0] + + # If there is more than one level found use the first + if len(tmp) > 1: + tmp = tmp[0] + + levelsopt[p] = tmp + + # Update dictionary of parameters to optimise with lists of new values; clear dictionary first + optparams = OrderedDict((key, list()) for key in optparams) + p = 0 + for key, value in optparams.items(): + optparams[key].append(levels[levelsopt[p], p]) + p += 1 + + return optparams, levelsopt + + +def fitness_max(filename, outputnames): + """Return the maximum value from specific outputs in a file. + + Args: + filename (dict): Ordered dictionary containing name of parameters to optimise and their values + outputnames (list): Names (IDs) of outputs (rxs) + + Returns: + maxvalue (array): Maximum value(s) from specific outputs + """ + + maxvalue = np.zeros(len(outputnames), dtype=floattype) + f = h5py.File(filename, 'r') + nrx = f.attrs['nrx'] + + i = 0 + for rx in range(1, nrx + 1): + tmp = f['/rxs/rx' + str(rx) + '/'] + if tmp.attrs['Name'] in outputnames: + fieldname = list(tmp.keys())[0] + maxvalue[i] = np.amax(tmp[fieldname]) + i += 1 + + return maxvalue + + + + + + \ No newline at end of file