From 0cd968bc7ec64efc0734a8eebe1d79aa733c47d3 Mon Sep 17 00:00:00 2001 From: Craig Warren Date: Mon, 25 Jan 2016 10:15:45 +0000 Subject: [PATCH] Adjustments to Taguchi optimisation structure. --- gprMax/optimisation_taguchi.py | 77 +++++++++++++++++++++++---------- user_libs/OA_18_7_3_2.npy | Bin 0 -> 206 bytes user_libs/OA_9_4_3_2.npy | Bin 0 -> 116 bytes 3 files changed, 53 insertions(+), 24 deletions(-) create mode 100644 user_libs/OA_18_7_3_2.npy create mode 100644 user_libs/OA_9_4_3_2.npy diff --git a/gprMax/optimisation_taguchi.py b/gprMax/optimisation_taguchi.py index cdb0a219..5f76280a 100644 --- a/gprMax/optimisation_taguchi.py +++ b/gprMax/optimisation_taguchi.py @@ -16,6 +16,7 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . +import os from collections import OrderedDict import numpy as np @@ -23,6 +24,8 @@ import numpy as np from gprMax.constants import floattype from gprMax.exceptions import CmdInputError +moduledirectory = os.path.dirname(os.path.abspath(__file__)) + def taguchi_code_blocks(inputfile, taguchinamespace): """Looks for and processes a Taguchi code block (containing Python code) in the input file. It will ignore any lines that are comments, i.e. begin with a double hash (##), and any blank lines. @@ -88,32 +91,58 @@ def construct_OA(optparams): # Number of parameters to optimise k = len(optparams) - p = int(np.ceil(np.log(k * (s - 1) + 1) / np.log(s))) - - # Number of experiments - N = s**p - - # Number of columns - cols = int((N - 1) / (s - 1)) - - # Algorithm to construct OA from: http://ieeexplore.ieee.org/xpl/articleDetails.jsp?reload=true&arnumber=6812898 - OA = np.zeros((N + 1, cols + 1), dtype=np.int8) - - # Construct basic columns - for ii in range(1, p + 1): - col = int((s**(ii - 1) - 1) / (s - 1) + 1) - for row in range(1, N + 1): - OA[row, col] = np.mod(np.floor((row - 1) / (s**(p - ii))), s) + # Load the appropriate OA + if k <= 4: + OA = np.load(os.path.join(moduledirectory, 'user_libs', 'OA_9_4_3_2.npy')) - # Construct non-basic columns - for ii in range(2, p + 1): - col = int((s**(ii - 1) - 1) / (s - 1) + 1) - for jj in range(1, col): - for kk in range(1, s): - OA[:, col + (jj - 1) * (s - 1) + kk] = np.mod(OA[:, jj] * kk + OA[:, col], s) + # Number of experiments + N = OA.shape[0] - # First row and first columns are unneccessary, only there to match algorithm, and cut down columns to number of parameters to optimise - OA = OA[1:, 1:k + 1] + # Number of columns of OA before cut down + cols = OA.shape[1] + + # Cut down OA columns to number of parameters to optimise + OA = OA[:, 0:k] + + elif k <= 7: + OA = np.load(os.path.join(moduledirectory, 'user_libs', 'OA_18_7_3_2.npy')) + + # Number of experiments + N = OA.shape[0] + + # Number of columns of OA before cut down + cols = OA.shape[1] + + # Cut down OA columns to number of parameters to optimise + OA = OA[:, 0:k] + + else: + p = int(np.ceil(np.log(k * (s - 1) + 1) / np.log(s))) + + # Number of experiments + N = s**p + + # Number of columns + cols = int((N - 1) / (s - 1)) + + # Algorithm to construct OA from: http://ieeexplore.ieee.org/xpl/articleDetails.jsp?reload=true&arnumber=6812898 + OA = np.zeros((N + 1, cols + 1), dtype=np.int8) + + # Construct basic columns + for ii in range(1, p + 1): + col = int((s**(ii - 1) - 1) / (s - 1) + 1) + for row in range(1, N + 1): + OA[row, col] = np.mod(np.floor((row - 1) / (s**(p - ii))), s) + + # Construct non-basic columns + for ii in range(2, p + 1): + col = int((s**(ii - 1) - 1) / (s - 1) + 1) + for jj in range(1, col): + for kk in range(1, s): + OA[:, col + (jj - 1) * (s - 1) + kk] = np.mod(OA[:, jj] * kk + OA[:, col], s) + + # First row and first columns are unneccessary, only there to match algorithm, and cut down columns to number of parameters to optimise + OA = OA[1:, 1:k + 1] return OA, N, cols, k, s, t diff --git a/user_libs/OA_18_7_3_2.npy b/user_libs/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